| 1 | !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ |
|---|
| 2 | ! User-defined Rate Law functions |
|---|
| 3 | ! Note: the default argument type for rate laws, as read from the equations file, is single precision |
|---|
| 4 | ! but all the internal calculations are performed in double precision |
|---|
| 5 | !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ |
|---|
| 6 | ! |
|---|
| 7 | ! RFo: IMPORTANT: Comments must ALWAYS stand BEFORE the code to be processed |
|---|
| 8 | ! This is just a savety copy. The UserRateLaws.f90 which is actually used |
|---|
| 9 | ! is copied from GASPASE_PREPROC/mechanisms/UserRateLaws.f90 to kpp/util/UserRateLaws.f90 |
|---|
| 10 | ! each time when run_kpp4palm.ksh is run |
|---|
| 11 | |
|---|
| 12 | !~~~> Arrhenius |
|---|
| 13 | KPP_REAL FUNCTION ARR( A0,B0,C0 ) |
|---|
| 14 | REAL A0,B0,C0 |
|---|
| 15 | ARR = DBLE(A0) * EXP(-DBLE(B0)/TEMP) * (TEMP/300.0_dp)**DBLE(C0) |
|---|
| 16 | END FUNCTION ARR |
|---|
| 17 | |
|---|
| 18 | |
|---|
| 19 | !~~~> Simplified Arrhenius, with two arguments, fixed TEMP=20 deg C |
|---|
| 20 | !~~~> Note: The argument B0 has a changed sign when compared to ARR |
|---|
| 21 | ! |
|---|
| 22 | !~~~> Simplified Arrhenius, with two arguments |
|---|
| 23 | ! RFo: same as above, but with TEMP in parameter list and _dp |
|---|
| 24 | ! KPP_REAL (kind=dp) FUNCTION ARR2( A0, B0 ) |
|---|
| 25 | ! REAL (kind=dp) A0,B0 |
|---|
| 26 | ! ARR2 = A0 * EXP( B0 / 20.0_dp ) |
|---|
| 27 | ! END FUNCTION ARR2 |
|---|
| 28 | REAL(kind=dp) FUNCTION ARR2( a0, b0, temp ) |
|---|
| 29 | REAL(kind=dp) :: temp |
|---|
| 30 | REAL(kind=dp) :: a0,b0 |
|---|
| 31 | ARR2 = a0 * EXP( -b0 / temp ) |
|---|
| 32 | END FUNCTION ARR2 |
|---|
| 33 | |
|---|
| 34 | |
|---|
| 35 | KPP_REAL FUNCTION EP2(A0,C0,A2,C2,A3,C3) |
|---|
| 36 | REAL A0,C0,A2,C2,A3,C3 |
|---|
| 37 | REAL(kind=dp) K0,K2,K3 |
|---|
| 38 | K0 = DBLE(A0) * EXP(-DBLE(C0)/TEMP) |
|---|
| 39 | K2 = DBLE(A2) * EXP(-DBLE(C2)/TEMP) |
|---|
| 40 | K3 = DBLE(A3) * EXP(-DBLE(C3)/TEMP) |
|---|
| 41 | K3 = K3*CFACTOR*1.0E6_dp |
|---|
| 42 | EP2 = K0 + K3/(1.0_dp+K3/K2 ) |
|---|
| 43 | END FUNCTION EP2 |
|---|
| 44 | |
|---|
| 45 | KPP_REAL FUNCTION EP3(A1,C1,A2,C2) |
|---|
| 46 | REAL A1, C1, A2, C2 |
|---|
| 47 | REAL(kind=dp) K1, K2 |
|---|
| 48 | K1 = DBLE(A1) * EXP(-DBLE(C1)/TEMP) |
|---|
| 49 | K2 = DBLE(A2) * EXP(-DBLE(C2)/TEMP) |
|---|
| 50 | EP3 = K1 + K2*(1.0E6_dp*CFACTOR) |
|---|
| 51 | END FUNCTION EP3 |
|---|
| 52 | |
|---|
| 53 | KPP_REAL FUNCTION FALL ( A0,B0,C0,A1,B1,C1,CF) |
|---|
| 54 | REAL A0,B0,C0,A1,B1,C1,CF |
|---|
| 55 | REAL(kind=dp) K0, K1 |
|---|
| 56 | K0 = DBLE(A0) * EXP(-DBLE(B0)/TEMP)* (TEMP/300.0_dp)**DBLE(C0) |
|---|
| 57 | K1 = DBLE(A1) * EXP(-DBLE(B1)/TEMP)* (TEMP/300.0_dp)**DBLE(C1) |
|---|
| 58 | K0 = K0*CFACTOR*1.0E6_dp |
|---|
| 59 | K1 = K0/K1 |
|---|
| 60 | FALL = (K0/(1.0_dp+K1))* & |
|---|
| 61 | DBLE(CF)**(1.0_dp/(1.0_dp+(LOG10(K1))**2)) |
|---|
| 62 | END FUNCTION FALL |
|---|
| 63 | |
|---|
| 64 | !--------------------------------------------------------------------------- |
|---|
| 65 | |
|---|
| 66 | ELEMENTAL REAL(kind=dp) FUNCTION k_3rd(temp,cair,k0_300K,n,kinf_300K,m,fc) |
|---|
| 67 | |
|---|
| 68 | INTRINSIC LOG10 |
|---|
| 69 | |
|---|
| 70 | REAL(kind=dp), INTENT(IN) :: temp ! temperature [K] |
|---|
| 71 | REAL(kind=dp), INTENT(IN) :: cair ! air concentration [molecules/cm3] |
|---|
| 72 | REAL, INTENT(IN) :: k0_300K ! low pressure limit at 300 K |
|---|
| 73 | REAL, INTENT(IN) :: n ! exponent for low pressure limit |
|---|
| 74 | REAL, INTENT(IN) :: kinf_300K ! high pressure limit at 300 K |
|---|
| 75 | REAL, INTENT(IN) :: m ! exponent for high pressure limit |
|---|
| 76 | REAL, INTENT(IN) :: fc ! broadening factor (usually fc=0.6) |
|---|
| 77 | REAL(kind=dp) :: zt_help, k0_T, kinf_T, k_ratio |
|---|
| 78 | |
|---|
| 79 | zt_help = 300._dp/temp |
|---|
| 80 | k0_T = k0_300K * zt_help**(n) * cair ! k_0 at current T |
|---|
| 81 | kinf_T = kinf_300K * zt_help**(m) ! k_inf at current T |
|---|
| 82 | k_ratio = k0_T/kinf_T |
|---|
| 83 | k_3rd = k0_T/(1._dp+k_ratio)*fc**(1._dp/(1._dp+LOG10(k_ratio)**2)) |
|---|
| 84 | |
|---|
| 85 | END FUNCTION k_3rd |
|---|
| 86 | |
|---|
| 87 | !--------------------------------------------------------------------------- |
|---|
| 88 | |
|---|
| 89 | ELEMENTAL REAL(kind=dp) FUNCTION k_arr (k_298,tdep,temp) |
|---|
| 90 | ! Arrhenius function |
|---|
| 91 | |
|---|
| 92 | REAL, INTENT(IN) :: k_298 ! k at T = 298.15K |
|---|
| 93 | REAL, INTENT(IN) :: tdep ! temperature dependence |
|---|
| 94 | REAL(kind=dp), INTENT(IN) :: temp ! temperature |
|---|
| 95 | |
|---|
| 96 | INTRINSIC EXP |
|---|
| 97 | |
|---|
| 98 | k_arr = k_298 * EXP(tdep*(1._dp/temp-3.3540E-3_dp)) ! 1/298.15=3.3540e-3 |
|---|
| 99 | |
|---|
| 100 | END FUNCTION k_arr |
|---|
| 101 | |
|---|
| 102 | !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ |
|---|
| 103 | ! End of User-defined Rate Law functions |
|---|
| 104 | !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ |
|---|
