[gghey]

Hello ,i want to ask this :

If i had to create a program where a particle will execute a random walk for N=1000 and must calculate the number of grid positions where the particle visited at least one time.And if i wanted to make 10000 runs and find 10 points (one for every 100 steps ,from 0 to 1000) , which will be the means of 10000 runs

I am doing this:

import scipy as sc
import matplotlib.pyplot as plt
import random

plegma=1000
grid=sc.ones(plegma)   # grid full of available positions(ones)    

for p in range(10000):
    #-------------------Initialize problem-------------------------------------------------
    his_pos=[]                  # list which holds the position of the particle in the grid
    in_pos = int(sc.random.randint(0,len(grid),1))            #initial position of particle
    means=[]                                                    #list which holds the means 
    #--------------------------------------------------------------------------------------
    for i in range(0,1000,100):
        step=2*sc.random.random_integers(0,1)-1        #the step of the particle can be -1 or 1
        # Check position for edges and fix if required
        # Move by step
        in_pos += step
        #Correct according to periodic boundaries
        in_pos = in_pos % len(grid)  

        #Keep track of random walk
        his_pos.append(in_pos)
        history=sc.array(his_pos)
        mean_his=sc.mean(history) 
        means.append(mean_his)


plt.plot(means,'bo')
plt.show()

I have these problems :

1) The 10000 runs is the accuracy we want,right?But ,for example, from my code i have this:

('his_pos=', [243]) 

('history=', array([243])) 

('mean_his=', 243.0) 

('means=', [243.0]) 

('his_pos=', [243, 244]) 

('history=',array([243, 244]))

('mean_his=', 243.5) 

('means=', [243.0, 243.5]) 

('his_pos=', [243, 244, 243]) 

('history=', array([243, 244, 243])) 

('mean_his=', 243.33333333333334) 

('means=', [243.0, 243.5,243.33333333333334])

and this goes on all over again from the beginning 10000 times. But ,like this ,the information gets lost!(am i right?).But if i put out of the loops the “mean=[]”,then this will keep 10000 values,but i want only 10!I am a bit confused!
Anyway,my plot is not the desired curve ,so somewhere i am making a mistake.

I'll be grateful if someone explains me some things.
Thanks!

(It's a double post but i can't figure the solution..)

[Motoma]

Please verify that I understand what you are trying to do:

• Run 10000 tests.
  
• Each test takes 1000 random steps.
  
• You want to average the first 100 steps across all 10000 runs, the second 100 steps across all 10000 runs, etc...

Is this correct, or am I misunderstanding?

[gghey]

Hello, yes you are right.Here is the exercise:

Quote

Create a program where a particle will execute a random walk for N=1000 steps for these two cases: i) in a 1D system ii) in a 2D system. The program must calculate the where S is the number of grid positions where the particle visited at least one time.You will make 10000 runs and find 10 points (one for every 100 steps ,from 0 to 1000) , which will be the means of 10000 runs.Do the plot of in relation to time t.

[Motoma]

Your loops shouldn't be nested, then:

runs = [[], [], [], [], [], [], [], [], [], []]

for i in range(10000):
    position = starting_position
    for j in range(10):
        position_slice = []
        for k in range(100):
            # take random step
            position_slice.append(position)
        runs[j].append(average(position_slice))

for i in range(10):
    slice_average = average(runs[i])
    print("Run %i: %r" % (i, slice_average))

As you can see above, I'm taking the shortcut of averaging the each 100 steps in order to save memory.

[gghey]

I did this ,but it gives me only 1 value!

plegma=1000
grid=sc.ones(plegma)   # grid full of available positions(ones)    

for p in range(10000):
    #-------------------Initialize problem-------------------------------------------------
    in_pos = int(sc.random.randint(0,len(grid),1))            #initial position of particle
    means=[]                                                    #list which holds the means 
    #--------------------------------------------------------------------------------------
    for k in range(10):
        his_pos=[]                  # list which holds the position of the particle in the grid
        for i in range(100):
            step=2*sc.random.random_integers(0,1)-1        #the step of the particle can be -1 or 1
            # Check position for edges and fix if required
            # Move by step
            in_pos += step
            #Correct according to periodic boundaries
            in_pos = in_pos % len(grid)  
                
            #Keep track of random walk
            his_pos.append(in_pos)
            history=sc.array(his_pos)
        mean_his=sc.mean(history) 
        means.append(mean_his)

for l in range(10):
    means_average = sc.mean(means[l])

This post has been edited by gghey: 16 November 2011 - 11:37 AM

[Motoma]

You are overwriting means_average ten times. In my code, it printed out each result. If you want, you can store it to a list:

final_list = []
for i in range(10):
    final_list.append(sc.mean(means[i])

[gghey]

I am still getting 1 value.
I can't understand!

[Motoma]

Perhaps I'm not understanding you: My code returns exactly ten values. Please post your code again.

[gghey]

........................
for i in range(100):
            step=2*sc.random.random_integers(0,1)-1        #the step of the particle can be -1 or 1
            # Check position for edges and fix if required
            # Move by step
            in_pos += step
            #Correct according to periodic boundaries
            in_pos = in_pos % len(grid)  
                
            #Keep track of random walk
            his_pos.append(in_pos)
            history=sc.array(his_pos)
            mean_his=sc.mean(history) 
        means.append(mean_his)
  
plt.plot(means,'bo')
plt.show()

It gives me 10 points but it supposes to form a curve but the don't

[Motoma]

How much of a curve do you think a random direction will create? For ever step you take in one direction, you'll take one in the opposite.

[gghey]

I must compare the points i am taking from the program with expression sn:

steps=sc.arange(1000)
sn=(8*steps/sc.pi)**(1/2.0)*(1+1/4*steps-3/64*(steps**2))
plt.plot(steps,sn,'r')


So,the points must be as close to this curve.

[gghey]

I am still having difficulties to understand the "make 10000 runs".

Like i said in my first post "..and this goes on all over again from the beginning 10000 times. But ,like this ,the information gets lost!"

I have a table "means" which has 10 values for every 100 steps.Until now,ok.
When i am putting it to the loop "for p in range(10000)" the table changes and this has as a result for me,that to the end i am taking another table with different values and the old ones are lost.So,what's the meaning ?
What is the purpose of making 10000 runs?It has some kind to do with the accuracy ?

[Motoma]

I'm quite sure that understanding and answering these two vital questions are the real purpose of the exercise.

The information doesn't get lost, it gets averaged. Aberrations get smoothed out, and you begin to see a statistical pattern happening.

If the process was run once, your randomness may leave you on the far edge of your map. The likelihood of this happening 10000 times is highly improbable. The more times you run the test, the closer you will get to having the probabilistic average, that is, the most likely scenario.

[gghey]

I understand what you are saying but from the other point of view ,we only get the last values from the "means" array.
Even if i run the code for 10 times or for 10000 times ,i still get the last value (one only value) from the "means".
I can't understand how this smoothness applies.


(Also,somewhere i am wrong because i can't get points close to the curve i want to compare)

And thanks for the help.

[Motoma]

In the pseudo-code example I posted above, I am creating 10 lists, each containing the mean position of each of the 10000 trials. Each of these 10000 element lists is then averaged to get ten, single average points for the entire session. If your code is written similarly (which I can't tell because you haven't posted the code in full), it is not overwriting anything, rather it is averaging.