mercredi 19 mai 2021

How to simulate First passage time in python for a Random Walk?

I have the following code for the calculation of the number of particles that hit a trap in x=-10 in 2D random walk, but the plot is very different from the expected curve:

import matplotlib.pyplot as plt 
import matplotlib
import numpy as np
import pylab
import random

n = 1000

n_simulations=1000

x = numpy.zeros((n_simulations,n))
y = numpy.zeros((n_simulations,n))


steps = np.arange(0, n, 1)
for i in range (n_simulations):
    for j in range (1,n):
        val=random.randint(1,5)
        if val == 1:
            x[i,j] = x[i,j - 1] + 1
            y[i,j] = y[i,j - 1] 
        elif val == 2:
            x[i,j] = x[i,j - 1] - 1
            y[i,j] = y[i,j - 1] 
        elif val == 3:
            x[i,j] = x[i,j - 1]
            y[i,j] = y[i,j - 1] + 1
        elif val == 4:
            x[i,j] = x[i,j - 1]
            y[i,j] = y[i,j - 1] - 1
        else:
            x[i,j] = x[i,j - 1]
            y[i,j] = y[i,j - 1]
        if x[i,j] == -10:
            break
fp = numpy.zeros((n_simulations,n)) # number of paricles that hit the trap for each simulation. 
for i in range(n_simulations):
    for j in range (1,n):
        if x[i,j] == -30:
            fp[i,j] = fp[i,j-1] +1
        else:
            fp[i,j] = fp[i,j-1]
        

plt.xlim(0,1000)

plt.plot(steps,fp[:,999])

plt.show()

The 2d random trajectory is correct but the probability given by the fp(number os particles that hit the trap for the first time) is wrong. The plot Nfpt/N with Nfpt, the number of particles that hit the trap should have a peak around t~70-100 but in this case the peak is in zero, wich doesn't make sense.




Aucun commentaire:

Enregistrer un commentaire