program FortranFigure4_3 implicit none save integer, parameter :: long = selected_real_kind(15,50) integer, parameter :: numbermesh=40 integer, parameter :: nahbon=100000 real(kind=long), dimension(nahbon) :: RanNum real(kind=long), dimension(numbermesh) :: A,B real(kind=long), dimension(numbermesh) :: k1S,k2S,k3S,k4S,DhS real(kind=long), parameter :: DiffConst=0.0001_long ! diffusion constant real(kind=long), parameter :: k1=0.001_long ! degradation rate real(kind=long), parameter :: k2=0.01_long ! creation rate real(kind=long), parameter :: k3=1.2_long ! degradation rate real(kind=long), parameter :: k4=1.0_long ! creation rate real(kind=long) :: time,Dh,h,ss,a0,finaltime,ate,tau integer :: ii,ir,ns integer size, seed(2) data seed /10,10/ size = 2 h=0.025_long ! h=1/real(numbermesh) Dh=DiffConst/(h*h) time=0_long ir=0 a0=0.0_long finaltime=30.0_long*60.0_long call random_seed(SIZE=size) call random_seed(PUT=seed(1:size)) call random_number(RanNum) do ii=1,numbermesh A(ii)=0.0_long end do do ii=2,numbermesh-1 DhS(ii)=2.0_long*Dh end do ii=1 DhS(ii)=Dh ii=numbermesh DhS(ii)=Dh do ii=1,numbermesh k1S(ii)=k1*A(ii)*(A(ii)-1) k2S(ii)=k2*A(ii)*B(ii) k3S(ii)=0.0_long k4S(ii)=0.0_long end do do ii=1,9*numbermesh/10 k3S(ii)=k3 end do do ii=2*numbermesh/5+1,numbermesh k4S(ii)=k4 end do do ii=1,numbermesh a0=a0+DhS(ii)*A(ii)+DhS(ii)*B(ii)+k1S(ii)+k2S(ii)+k3S(ii)+k4S(ii) end do do while (timenahbon) then call random_number(RanNum) ir=1 end if tau=(1/a0)*log(1/RanNum(ir)) ir=ir+1 if (ir>nahbon) then call random_number(RanNum) ir=1 end if ss=0_long ii=0 do while ((ssRanNum(ir)*a0) then a0=a0-k1S(ii)-k2S(ii) A(ii)=A(ii)-2_long k1S(ii)=k1*A(ii)*(A(ii)-1) k2S(ii)=k2*A(ii)*B(ii) a0=a0+k1S(ii)+k2S(ii)-2.0_long*DhS(ii) else ii=0 do while ((ssRanNum(ir)*a0) then a0=a0-k1S(ii)-k2S(ii) A(ii)=A(ii)-1_long B(ii)=B(ii)-1_long k1S(ii)=k1*A(ii)*(A(ii)-1) k2S(ii)=k2*A(ii)*B(ii) a0=a0+k1S(ii)+k2S(ii)-2.0_long*DhS(ii) else ii=0 do while ((ssRanNum(ir)*a0) then a0=a0-k1S(ii)-k2S(ii) A(ii)=A(ii)+1_long k1S(ii)=k1*A(ii)*(A(ii)-1) k2S(ii)=k2*A(ii)*B(ii) a0=a0+k1S(ii)+k2S(ii)+DhS(ii) else ii=0 do while ((ssRanNum(ir)*a0) then a0=a0-k2S(ii) B(ii)=B(ii)+1_long k2S(ii)=k2*A(ii)*B(ii) a0=a0+k2S(ii)+DhS(ii) else ii=0 do while ((ssRanNum(ir)*a0) then a0=a0-k1S(ii)-k2S(ii)-k1S(ii+1)-k2S(ii+1) A(ii)=A(ii)-1_long A(ii+1)=A(ii+1)+1_long k1S(ii)=k1*A(ii)*(A(ii)-1) k2S(ii)=k2*A(ii)*B(ii) k1S(ii+1)=k1*A(ii+1)*(A(ii+1)-1) k2S(ii+1)=k2*A(ii+1)*B(ii+1) a0=a0+k1S(ii)+k2S(ii)+k1S(ii+1)+k2S(ii+1) a0=a0+DhS(ii+1)-DhS(ii) else ii=1 do while ((ssRanNum(ir)*a0) then a0=a0-k1S(ii)-k2S(ii)-k1S(ii-1)-k2S(ii-1) A(ii)=A(ii)-1_long A(ii-1)=A(ii-1)+1_long k1S(ii)=k1*A(ii)*(A(ii)-1) k2S(ii)=k2*A(ii)*B(ii) k1S(ii-1)=k1*A(ii-1)*(A(ii-1)-1) k2S(ii-1)=k2*A(ii-1)*B(ii-1) a0=a0+k1S(ii)+k2S(ii)+k1S(ii-1)+k2S(ii-1) a0=a0+DhS(ii-1)-DhS(ii) else ii=0 do while ((ssRanNum(ir)*a0) then a0=a0-k2S(ii)-k2S(ii+1) B(ii)=B(ii)-1_long B(ii+1)=B(ii+1)+1_long k2S(ii)=k2*A(ii)*B(ii) k2S(ii+1)=k2*A(ii+1)*B(ii+1) a0=a0+k2S(ii)+k2S(ii+1)+DhS(ii+1)-DhS(ii) else ii=1 do while ((ssRanNum(ir)*a0) then a0=a0-k2S(ii)-k2S(ii-1) B(ii)=B(ii)-1_long B(ii-1)=B(ii-1)+1_long k2S(ii)=k2*A(ii)*B(ii) k2S(ii-1)=k2*A(ii-1)*B(ii-1) a0=a0+k2S(ii)+k2S(ii-1)+DhS(ii-1)-DhS(ii) else write(*,*) 'Error.' end if end if end if end if end if end if end if end if time=time+tau if (int(0.01_long*time)-int(0.01_long*(time-tau))>0.5_long) then write(*,*) 'Time computed:',nint(time),'seconds (final time 1800 sec)' end if end do open(unit=1, file='dataFigure4_3.dat', status='replace') 100 format(i6,X,i6) do ii=1,numbermesh write(1,100) nint(A(ii)),nint(B(ii)) end do close(unit=1) end program FortranFigure4_3