source: palm/trunk/UTIL/chemistry/gasphase_preproc/kpp/drv/general_stochastic.f90 @ 3136

Last change on this file since 3136 was 2696, checked in by kanani, 7 years ago

Merge of branch palm4u into trunk

File size: 1.6 KB
Line 
1PROGRAM KPP_ROOT_Driver
2
3  USE KPP_ROOT_Model
4  USE KPP_ROOT_Initialize, ONLY: Initialize
5
6  IMPLICIT NONE
7     
8!~~~> No. of Stochastic events per simulation snaphot     
9      INTEGER :: Nevents
10!~~~> No. of Stochastic tau-steps simulation snaphot     
11      INTEGER :: Nsteps
12!~~~> No. of Molecules     
13      INTEGER :: NmlcV(NVAR), NmlcF(NFIX)
14!~~~> Random numbers     
15      REAL :: r1, r2
16!~~~> Local variables     
17      INTEGER :: i
18      KPP_REAL :: T, Tau, SCT(NREACT)
19
20!~~~~> Initialize and prescribe volume and no. of events
21      CALL Initialize()
22      Volume  = 1000.0d0
23      Nevents = 5000
24      Nsteps = 10
25      Tau = (TEND-TSTART)/100.0
26!~~~~~~~~~~~~~~~~     
27     
28!~~~> Translate initial values from conc. to NmlcVules
29      NmlcV(1:NVAR) = INT(Volume*VAR(1:NVAR))
30      NmlcF(1:NFIX) = INT(Volume*FIX(1:NFIX))
31!~~~> Compute the stochastic reaction rates
32      CALL Update_RCONST()
33      CALL StochasticRates( RCONST, Volume, SCT )   
34
35!~~~> Save initial data
36      OPEN(10, file='KPP_ROOT_stochastic.dat')
37      WRITE(10,992) T, (NmlcV(i),i=1,NVAR)
38           
39!~~~> TIME loop starts
40      T = TSTART
41kron: DO WHILE (T < TEND)
42
43        WRITE(6,991) T,(SPC_NAMES(i),NmlcV(i), i=1,NVAR)
44
45!~~~> Choose here one of the following time-stepping routines
46        CALL Gillespie(Nevents, T, SCT, NmlcV, NmlcF)
47!        CALL TauLeap(Nsteps, Tau, T, SCT, NmlcV, NmlcF)
48       
49        WRITE(10,992) T, (NmlcV(i),i=1,NVAR)
50       
51      END DO kron
52!~~~> TIME loop ends
53
54      WRITE(6,991) T,(SPC_NAMES(i), NmlcV(i), i=1,NVAR)
55
56      CLOSE(10)
57     
58991   FORMAT('T=',F12.3,20(A,'=',I5,'; '))
59992   FORMAT(E24.14,200(1X,I12))
60
61
62END PROGRAM KPP_ROOT_Driver
63
Note: See TracBrowser for help on using the repository browser.