mardi 7 juin 2016

Simulation and Histogram of Random Hopping in Python

I explain you my problem:

Imagine you have a bar, with say s positions. Each position can be counted as position 0, position 1, .. , position s-1. Now what I want to do is simulate the following : At a certain point in time, a number of particles, say n particles, start in a state of the bar (assume a position in the middle). At this point with random probabilities pr and pl (pr + pl =1) this particles go right or left respectively. So basically the probability reflects the proportion of particles swapping right or left.

I want to repeat this a number of times and see what is the final position of the particles and plot an histogram of it. This is my function hop, which I made in order to simulate the hopping of the particles.

def hop(n):
'''
n is the number of particles starting in the middle position.
'''

s = [0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,n,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0]
global ls
ls = len(s)
i = 0

while i < 100:
    for j in range(0,int(len(s))):
        if s[j] != 0 :
            pr = random.random()
            pl = 1 - pr
            if j - 1 < 0:
                s[j+1] = s[j+1]+int(s[j]*pr)
                s[j] = s[j] - int(s[j]*pr)
            elif len(s) <= j+1:
                s[j-1] = s[j-1] + int(s[j]*pl)
                s[j] = s[j] - int(s[j]*pl)
            else:
                s[j-1] = s[j-1] + int(s[j]*pl)
                s[j+1] = s[j+1] + int(s[j]*pr)
                s[j] = s[j] - int(s[j]*pr) - int(s[j]*pl)
            j+=1
        elif s[j] == 0:
            s[j] = 0
            j+=1
    i+=1
return s

And here is the rest, that I used to plot the histogram:

x = hop(100)
y = sum(x) #This is for debugging purposes, I want to check that I'm left 
           with the right number of particles
list1 = []
for k in range(0,ls):
    list1.append(k)
plt.hist(x,list1)
plt.show()

Where I have imported mathplotlib and specifically I've imported

import matplotlib.pyplot as plt
import random

My problem is that from the histograms that I obtain this is doing something very wrong. Indeed the histograms are all skewed to the left, which wouldn't be possible If the probabilities are taken at random. Furthermore the histogram doesn't show me the right amounts of particles.

Does anyone understand what is going wrong?

Thanks




Aucun commentaire:

Enregistrer un commentaire