import matplotlib.pyplot as plt
import numpy as np
class diffusion:
    def __init__(self,constant=0.3,l=500,update_times=50):
        self.l=l
        self.c=constant
        self.ss=[[0]*self.l]
        self.n=update_times
    def diff(self):
        self.ss[-1][int(self.l/2)]=1
        temp=[0]*self.l
        counter=1
        while(1):
            for i in range(len(self.ss[-1])-2):
                temp[i+1]=self.ss[-1][i+1]+self.c*(self.ss[-1][i+2]+self.ss[-1][i]-2*self.ss[-1][i+1])
            self.ss.append(temp)
            counter+=1
            if counter>self.n:
                break
    def show(self):
        x=np.arange(0,self.l,1)
        plt.plot(x,self.ss[-1],'.',label='time=%.2f'%(len(self.ss)-1))

        plt.legend(frameon=True)
        plt.grid(True)
        plt.xlabel('x')
        plt.ylabel('density')
        plt.title('diffusion in one dimension')  

a=diffusion(update_times=100)
a.diff()
a.show()
a=diffusion(update_times=500)
a.diff()
a.show()
a=diffusion(update_times=1000)
a.diff()
a.show()
a=diffusion(update_times=10000)
a.diff()
a.show()
a=diffusion(update_times=50000)
a.diff()
a.show()