Plotting the mean square displacement of a 2D random walk as a function of δt

491 Views Asked by At

I've already created a code for random walk of 10000 steps and then repeated it 12 times and stored each run in a separate text file (which was required in the question). I then calculated the mean square displacement of it(not sure if it's done correct). I now need to 'plot my Mean Square Displacement as a function of δt, including errorbars σ = std(MSD)/√N, where std(MSD) is the standard deviation among the different runs and N is the number of runs.' and then compute the diffusion constant D from the curve and check that D = 2 (∆/dt) where dt = 1.

Here is my code so far:

import numpy as np 
import matplotlib.pyplot as plt
import random as rd
import math

a = (np.zeros((10000, 2), dtype=np.float))  
     
def randwalk(x,y):
     theta= 2*math.pi*rd.random() 
     x+=math.cos(theta); # This uses the equation given, since we are told the spatial unit = 1
     y+=math.sin(theta); 
     return (x,y)
     
x, y = 0.,0.  
for i in range(10000): # Using for loop and range function to initialize the array
     x, y = randwalk(x,y)
     a[i,:] = x,y
     

fn_base = "random_walk_%i.txt" # Saves each run in a numbered text file, fn_base is a varaible to hold format

N = 12
for j in range(N): 
    rd.seed(j) # seed(j) explicitly sets the seed to random numbers
    x , y = 0., 0.
    for i in range(10000):
        x, y = randwalk(x,y)
        a[i,:] = x, y
    fn = fn_base % j
    np.savetxt(fn, a)

destinations = np.zeros((12, 2), dtype=np.float)
for j in range(12):
    x, y = 0., 0.
    for i in range(10000):
        x, y = randwalk(x, y)
    destinations[j] = x, y

square_distances = destinations[:,0] ** 2 + destinations[:,1] ** 2
m_s_d = np.mean(square_distances)

I think that to do it I just have to plot the msd against the number of steps? But I'm not sure how to do this. I saw a similar question on stackoverflow but the code for it is different than mine and I don't understand how to use that for my code.

I tried to do next

plt.figure()
t = 10000
plt.plot(m_s_d, t) 
plt,show()

But this gives an error as the dimensions are not equal.

Edit ** I think my issue is that I am trying to plot it against number of steps when I should be plotting it against the change in time. However I can’t work out how to calculate the change in time dt?

Apologies in advance is question isn't formulated well, I am fairly new to computing. Thank you.

0

There are 0 best solutions below