chem/KPP/kpp/kpp-2.1/test/pan.k

   1 #MODEL pan
   2 #INTEGRATOR ros3
   3 #DRIVER ./pan_drv
   4
   5 #INLINE F_DECL_INT
   6       REAL PRESS
   7 #ENDINLINE
   8
   9 #INLINE F95_UTIL_INT
  10 !************** SPECIAL RATE FUNCTIONS **********************
  11       DOUBLE PRECISION FUNCTION RJPL( K0300, Q, KU300, R, M, T )
  12       IMPLICIT NONE
  13       DOUBLE PRECISION k0300,q,ku300,r,m,t
  14       DOUBLE PRECISION tt,k0,ku,k0m,kk,lgkk,e,f
  15 ! JPL standard three body reaction rate format extended
  16       TT= T / 3.D2
  17       K0= K0300 * exp(-1.D0*Q*log(TT))
  18       KU= KU300 * exp(-1.D0*R*log(TT))
  19       K0M= K0 * M
  20       KK= K0M / KU
  21       LGKK=0.43429448190324926D0 * LOG(KK) ! = log10(KK)
  22       E=1.D0 / ( 1.D0 + LGKK*LGKK )
  23       F=exp(-0.5108256237659887D0*E)       ! -0.51=log(0.6)
  24       RJPL = F * K0M / ( 1.D0 + KK )
  25       END FUNCTION RJPL
  26 !---------------------------------------------------------------------
  27 !TROE(FIX(I_M),TEMP,.933,2.85E-30,-2.67,3.13E-11,363.)!Dransfield et al.'99(GRL) this call in XYZ.eqn
  28 !---------------------------------------------------------------------
  29       DOUBLE PRECISION FUNCTION TROE(M,T,beta,k0,k0e,kinf,Tc)
  30       IMPLICIT NONE
  31       DOUBLE PRECISION M,T,beta,k0,k0e,kinf,Tc
  32       DOUBLE PRECISION k0t,bcrit,Trat,dN,N,Bx,F
  33 !  real Troe rate constants: for OH + NO2 -> HNO3, Dransfield et al. 1999 (GRL)
  34       k0t = k0 * exp(k0e*log(T/3.D2))
  35       bcrit = beta*M*k0t/kinf
  36       Trat = T/Tc
  37       dN=sign(0.1D0-0.2605766891419492D0*Trat,1.D0-bcrit) ! 0.26=0.6*.434;log-->log10
  38       N = 0.75D0 + 0.5515539920171264D0*Trat           ! 0.55=1.27*.434;log-->log10
  39       Bx = (0.43429448190324926D0*log(bcrit)-0.12D0) / (N+dN)
  40       F = exp(-1.D0*Trat/(1. + Bx*Bx))
  41       TROE = k0t * (beta*M/(1.+bcrit)) * F
  42       END FUNCTION TROE
  43 !---------------------------------------------------------------------
  44       DOUBLE PRECISION FUNCTION RHNO3(M,T)
  45       IMPLICIT NONE
  46       DOUBLE PRECISION M,T
  47       DOUBLE PRECISION K0,K2,K3
  48 !     SPECIAL RATE CONSTANTS:    OH + HNO3 {+M} -->  NO3
  49 !     taken from S. Brown et al. 1999 GRL, JPL 2000
  50       K0=2.4D-14*EXP(460.D0/T)
  51       K2=2.7D-17*EXP(2199.D0/T)
  52       K3=M*6.5D-34*EXP(1335.D0/T)
  53       RHNO3 = K0 + K2 / ( 1 + K2/K3 )
  54       END FUNCTION RHNO3
  55 !---------------------------------------------------------------------
  56       DOUBLE PRECISION FUNCTION RHO2HO2(M,H2O,T)
  57       IMPLICIT NONE
  58       DOUBLE PRECISION M,H2O,T
  59       DOUBLE PRECISION RX1,RX2,RX3
  60 ! rate constant of the HO2 + HO2 --> H2O2 + O2 reaction
  61       RX1= 2.3D-13 *EXP(600.D0/T)
  62       RX2= 1.7D-33 *EXP(1000.D0/T) * M
  63       RX3= 1.4D-21 *EXP(2200.D0/T) * H2O
  64       RHO2HO2 = (RX1 + RX2)*(1.D0 + RX3)
  65       END FUNCTION RHO2HO2
  66 !---------------------------------------------------------------------
  67       DOUBLE PRECISION FUNCTION PHUX(X,Y,Z,CHI)
  68 !     BERECHNUNG VON PHOTOLYSERATEN MIT EINEM ALGORITHMUS AUS
  69 !     ROETHS FLUX-PROGRAMM
  70 !     CHI IN RADIANT(BOGENMASS)
  71 !     X,Y,Z WERDEN VON ROETH UEBERNOMMEN
  72 !     X  IST EINE MAXIMALPHOTOLYSERATE FUER CHI=0
  73 !     KUHN 07.09.93
  74 !rvk: no minimal photolysis rate (use zero instead, since KPP has no problems with that)
  75       IMPLICIT NONE
  76       DOUBLE PRECISION X,Y,Z,CHI,CHIZ,YCHIZ,MINYZ,EYCHIZ,EMINYZ
  77       PARAMETER (MINYZ = -30.D0, EMINYZ = 9.357623D-14 ) !EMINYZ=EXP(MINYZ)
  78       CHIZ   = CHI * Z
  79 ! BERECHNUNG DES AUSDRUCKES NUR FUER CHIZ KLEINER PI/2  (COS > 0)
  80       IF (CHIZ.LT.1.57079632679489D0) THEN
  81          YCHIZ = Y * (1. - (1./ COS(CHIZ) ) )
  82 ! SKALIERUNGSFAKTOR GROESSER EXP(-MINYZ)
  83          IF (YCHIZ.GT.MINYZ) THEN
  84             EYCHIZ =  EXP (YCHIZ)
  85          ELSE
  86 !            EYCHIZ =  EMINYZ
  87             EYCHIZ = 0.D0
  88          ENDIF
  89       ELSE
  90 !         EYCHIZ = EMINYZ
  91          EYCHIZ = 0.D0
  92       ENDIF
  93       PHUX = X * EYCHIZ
  94       if (PHUX.lt.1.D-10) PHUX = 0.D0
  95       END FUNCTION PHUX
  96 #ENDINLINE
  97