! module chemrates_module use params_module,only: nlev,nlevp1,nlonp4 ! ! Chemical reaction rates. ! ! Temperature independent rates are real parameters (constants), ! set at compile time. ! Temperature dependent rates are dimensioned (nlev,nlon) and are ! updated by sub rates_tdep at every latitude of every time step. ! See refs: ! Roble,R.G., 1995. Energetics of the Mesosphere and Thermosphere, ! Geophysical Monograph 87, American Geophysical Union ! R.G. Roble, E.C. Ridley ! An auroral model for the NCAR thermospheric general circulation model, ! Annales Geophysicae,5A, (6), 369-382, 1987. ! implicit none ! ! Temperature independent reaction rate constants: real,parameter :: ! ! Ion chemistry: | rk4 = 1.0E-10, ! O2+ + N4S -> NO+ + O + 4.21 eV | rk5 = 4.4E-10, ! O2+ + NO -> NO+ + O2 + 2.813 eV | rk6 = 4.0E-10, ! N+ + O2 -> O2+ + N4S + 2.486 eV | rk7 = 2.0E-10, ! N+ + O2 -> NO+ + O + 6.669 eV | rk8 = 1.0E-12, ! N+ + O -> O+ + N + 0.98 eV | rk9 = 6.0E-11, ! N2+ + O2 -> O2+ + N2 + 3.52 eV | rk10 = 5.0E-11, ! O+ + N2D -> N+ + O + 1.45 eV | rk12 = 6.0E-10, | rk13 = 9.4E-10, | rk14 = 2.0E-9, | rk15 = 2.4E-9, | rk16 = 3.4E-10, ! O+(2P) + N2 -> N2+ + O + 3.02 eV | rk17 = 1.0E-10, ! O+(2P) + N2 -> N+ + NO + 0.70 eV | rk18 = 4.0E-10, ! O+(2P) + O -> O+ + O + 5.2 eV | rk21 = 0.047, ! O+(2P) -> O+ + hv (2470A) | rk22 = 0.172, ! O+(2P) -> O+ + hv (7320A) | rk23 = 4.8E-10, ! O+(2D) -> N2+O + 1.33 eV | rk24 = 8.0E-10, ! O+(2D) + O -> O+(4S) +e+ 3.31 eV | rk25 = 5.0E-12, ! O+(2D) + e -> O+(4S) + e + 3.31 eV | rk27 = 7.0E-10, | rk28 = 7.7E-5, ! ! Odd nitrogen reaction rates. | beta2 = 5.0E-12, ! N2D + O2 -> NO + O1D + 1.84 eV ! | beta4 = 5.0E-13, ! N2D + O -> N4S + O + 2.38 eV | beta4 = 7.0E-13, ! N2D + O -> N4S + O + 2.38 eV | beta6 = 7.0E-11, ! N2D + NO -> N2 + O + 5.63 eV | beta7 = 1.06E-5, ! N2D -> N4S + hv | beta8 = 5.0E-11, ! NO + hv -> N4S + O | beta11= 9.3E-12, | beta13= 3.0E-12, | beta14= 6.7E-11, ! ! Neutral chemical reaction rates for mesospheric extension. | rkm3 = 2.2e-10, | rkm4 = 1.0E-10, | rkm5a = 1.4e-10, | rkm5b = 1.4E-11, | rkm7a = 1.2E-10, | rkm7b = 1.2E-10, | rkm8 = 8.0E-12, ! | rkm8 = 2.0E-11, | rkm10 = 1.0E-20, | rkm11 = 1.3E-16, | rkm13 = 2.1E-15, | rkm14 = 4.2E-13, | rkm15 = 2.2E-11, | rkm16 = 8.0E-14, | rkm17 = 3.9E-17 ! ! OXYGEN-HYDROGEN ! real,parameter :: RKM38 = 8.1E-12 ! cmph2, cmphox, cmpox, qjchem real,parameter :: RKM45 = 1.4E-10 ! cmpch4, cmpco, cmph2o, cmphox ! ! Carbon reactions ! real,parameter :: GAM2A = 0. ! cmphox, cmpmeta, qjchem real,parameter :: GAM2B = 0. ! cmph2, cmpmeta, qjchem real,parameter :: GAM4 = 1.0E-10 ! cmphox, cmpmeta, cmpox, qjchem real,parameter :: GAM10 = 1.0E-11 ! cmph2o, cmphox, cmpmeta, qjchem real,parameter :: GAM15 = 1.E-13 ! cmpco, cmpco2, cmpmeta, cmpox,qjchem ! ! CO2 number density mixing ratio: ! real,parameter :: co2mix = 180.e-6 real,parameter :: co2mix = 360.e-6 ! real,parameter :: co2mix = 720.e-6 ! ! a1d,asg,adl moved from /NMBTRF/ (crates_tdep.h) real,parameter :: a1d = 6.81E-3 real,parameter :: asg = 0.0758 real,parameter :: adl = 2.58E-4 ! ! Temperature dependent chemical reaction rate coefficients. ! These are set by rates_tdep for every latitude at every timestep. ! real,dimension(:,:,:),allocatable,save :: ! (nlevp1,lon0:lon1,lat0:lat1) ! ! Ion chemistry: | rk1, ! O+ + O2 -> O + O2+ + 1.555 eV | rk2, ! O+ + N2 -> NO+ + N4s + 1.0888 eV | rk3, ! N2+ + O -> NO+ + N2D + 0.70 eV | ra1, ! NO+ + e -> (20% N4S + O + 2.75 eV), (80% N2D + O + 0.38 eV) | ra2, ! O2+ + e -> (15% O + O + 6.95 eV), (85% O + O + 4.98 eV) | ra3, ! N2+ + e -> (10% N4S + N4S + 5.82 eV), (90% N4S + N2D + 3.44 eV) ! ! Neutral chemistry: | beta1, ! N4S + O2 -> NO + O + 1.4 eV | beta3, ! N4S + NO -> N2 + O + 3.25 eV | beta5, ! N2D + e -> N4S + e + 2.38 eV | beta9, ! NO + hv/Ly-a -> NO+ + e | beta10, | beta12, ! | rk11, | rk19, ! O+(2P) + e -> O+(4S) + e + 5.0 eV | rk20, ! O+(2P) + e -> O+(2O) + e + 1.69 eV | rk26, ! | rkm12, ! O + O + N2 -> O2 + N2 | tvib, | disn2p,! photodissociation of N2 by EUV (move to qrj?) ! ! Time-gcm only: ! ! Metastables: | rkm1 , ! | rkm2a, ! | rkm2b, ! | rkm6 , ! | rkm9 , ! ! ! Oxygen-hydrogen interactions | rkm20, | rkm21, | rkm22, | rkm23, | rkm24, | rkm25, | rkm26, | rkm27, | rkm28, | rkm29, | rkm30, | rkm31, | rkm32, | rkm33, | rkm34, | rkm35, | rkm36, | rkm37, | rkm39, | rkm40, | rkm41, | rkm42, | rkm43, | rkm44, ! ! D-region ion chemistry | rin9 , | rin10, | rin11, | rin12, | rin13, | rin17, | rin18, | rin19, | rin20, | rin21, | rin40, ! ! Carbon reactions | gam1 , | gam3 , | gam5 , | gam6 , | gam7 , | gam8 , | gam9 , | gam11, | gam12, | gam13, | gam14, ! ! Chlorine reactions | del1, | del2, | del3 ! ! Matrices for partitioning of OX (set by comp_ox.F, used by comp_major.F) real,allocatable :: fs(:,:,:,:,:) ! (i,k,3,0:3,j) ! contains !----------------------------------------------------------------------- subroutine alloc_tdep ! ! Allocate temperature-dependent reaction rates for task subdomain: ! Called once per run from init_fields. ! use mpi_module,only: lon0,lon1,lat0,lat1 integer :: istat ! allocate(rk1(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' rk1: stat=',i3)") istat allocate(rk2(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' rk2: stat=',i3)") istat allocate(rk3(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' rk3: stat=',i3)") istat allocate(ra1(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' ra1: stat=',i3)") istat allocate(ra2(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' ra2: stat=',i3)") istat allocate(ra3(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' ra3: stat=',i3)") istat allocate(beta1(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' beta1: stat=',i3)") istat allocate(beta3(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' beta3: stat=',i3)") istat allocate(beta5(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' beta5: stat=',i3)") istat allocate(beta9(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' beta9: stat=',i3)") istat allocate(beta10(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' beta10: stat=',i3)") istat allocate(beta12(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' beta12: stat=',i3)") istat allocate(rk11(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' rk11: stat=',i3)") istat allocate(rk19(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' rk19: stat=',i3)") istat allocate(rk20(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' rk20: stat=',i3)") istat allocate(rk26(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' rk26: stat=',i3)") istat allocate(rkm12(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' rkm12: stat=',i3)") istat allocate(tvib(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' tvib: stat=',i3)") istat allocate(disn2p(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' disn2p: stat=',i3)") istat ! allocate(fs(lon0:lon1,nlevp1,3,0:3,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_tdep: error allocating', | ' fs: stat=',i3)") istat ! ! Time-gcm only: ! ! Metastables: allocate(rkm1 (nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm2a(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm2b(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm6 (nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm9 (nlevp1,lon0:lon1,lat0:lat1),stat=istat) ! ! Oxygen-hydrogen interactions allocate(rkm20(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm21(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm22(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm23(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm24(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm25(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm26(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm27(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm28(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm29(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm30(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm31(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm32(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm33(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm34(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm35(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm36(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm37(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm39(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm40(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm41(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm42(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm43(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rkm44(nlevp1,lon0:lon1,lat0:lat1),stat=istat) ! ! D-region ion chemistry allocate(rin9 (nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rin10(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rin11(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rin12(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rin13(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rin17(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rin18(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rin19(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rin20(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rin21(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(rin40(nlevp1,lon0:lon1,lat0:lat1),stat=istat) ! ! Carbon reactions allocate(gam1 (nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(gam3 (nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(gam5 (nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(gam6 (nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(gam7 (nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(gam8 (nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(gam9 (nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(gam11(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(gam12(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(gam13(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(gam14(nlevp1,lon0:lon1,lat0:lat1),stat=istat) ! ! Chlorine reactions allocate(del1(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(del2(nlevp1,lon0:lon1,lat0:lat1),stat=istat) allocate(del3(nlevp1,lon0:lon1,lat0:lat1),stat=istat) ! end subroutine alloc_tdep !----------------------------------------------------------------------- subroutine chemrates_tdep(tn,te,ti,fno2,fnvo2,lev0,lev1,lon0,lon1, | lat) ! ! Calculate temperature-dependent reaction rates (called at each latitude) ! use input_module,only: f107 ! 10.7 cm flux (from input and/or gpi) use init_module,only: sfeps ! flux variation from orbital excentricity use init_module,only: istep ! for debug only use cons_module,only: check_exp ! ! Args: integer,intent(in) :: | lev0,lev1, ! first and last level indices, this task | lon0,lon1, ! first and last longitude indices, this task | lat ! latitude index real,dimension(lev0:lev1,lon0-2:lon1+2),intent(in) :: | tn, ! neutral temperature (deg K) | te, ! electron temperature (deg K) | ti, ! ion temperature (deg K) | fno2, ! o2 line integral (see chapman.F) | fnvo2 ! o2 column number density (see chapman.F) ! ! Local: integer :: k,i real :: ti1(lev0:lev1,lon0:lon1), | ti2(lev0:lev1,lon0:lon1), | ti3(lev0:lev1,lon0:lon1), | etvib(lev0:lev1,lon0:lon1) ! ! expo() (util.F) is used only if check_exp is true. This will avoid ! NaNS fpe, but will degrade performance. Check_exp is in cons.F. real,external :: expo ! ! write(6,"('Enter chemrates_tdep: lat=',i3,' lon0,1=',2i3, ! | ' lev0,1=',2i3)") lat,lon0,lon1,lev0,lev1 do i=lon0,lon1 do k=lev0,lev1-1 ! ! ti1 = T1/300. (S15) ti2 = T2/300. (S14) ti3 = TR/300. (S13) ! ti1(k,i) = (0.667* ti(k,i)+0.333 *tn(k,i))/300. ti2(k,i) = (0.6363*ti(k,i)+0.3637*tn(k,i))/300. ti3(k,i) = .5*(ti(k,i)+tn(k,i))/300. ! ! Rate coefficient for O+ + O2 -> O + O2+ + 1.555eV (HIERL) rk1(k,i,lat)=1.6E-11*ti1(k,i)**(-0.52)+5.5E-11* | exp(-22.85/ti1(k,i)) ! ! k2: O+ + N2 -> NO+ + N4s + 1.0888 eV ! (assumes t2/300 in ti2) rk2(k,i,lat)=(8.6E-14*ti2(k,i)-5.92E-13)*ti2(k,i)+1.533E-12 ! ! Rate coefficient for O+ + N2 -> N + NO+ ! See also sub altv (altv.f). Tvib is output of sub altv, if it ! is called (see dynamics.f). ! tvib(k,i,lat) = tn(k,i) ! tvib in crates_tdep.h etvib(k,i) = exp(-3353./tvib(k,i,lat)) ! etvib is local ! ! rk2: O+ + N2 -> NO+ + N4s + 1.0888 eV rk2(k,i,lat) = (((((270.*etvib(k,i)+220.)*etvib(k,i)+85.)* | etvib(k,i)+38.)*etvib(k,i)+1.)*rk2(k,i,lat)*etvib(k,i)+ | rk2(k,i,lat))*(1.-etvib(k,i)) ! ! rk3: N2+ + O -> NO+ + N2D + 0.70 eV if (300.*ti3(k,i) >= 1500.) then rk3(k,i,lat)=5.2E-11*ti3(k,i)**0.2 else rk3(k,i,lat)=1.4E-10*ti3(k,i)**(-0.44) endif ! ! rk11: rk11(k,i,lat) = rk12*8./9.*sqrt((ti(k,i)+tn(k,i)/16.)/ | (tn(k,i)+ti(k,i)/16.)) ! ! rk19: O+(2P) + e -> O+(4S) + e + 5.0 eV rk19(k,i,lat) = 4.0E-8*sqrt(300./te(k,i)) ! ! rk20: O+(2P) + e -> O+(2O) + e + 1.69 eV rk20(k,i,lat) = 1.5E-7*sqrt(300./te(k,i)) ! ! rd26: rk26(k,i,lat) = 6.6e-8*sqrt(300./te(k,i)) ! ! rkm12: O + O + N2 -> O2 + N2 (loss of O to O2) ! (same for O + O + O2 -> O2 + O2, Dickinson, et.al., 1984) rkm12(k,i,lat) = 9.59E-34*exp(480./tn(k,i)) ! ! ra1: NO+ + e -> (20% N4S + O + 2.75 eV), (80% N2D + O + 0.38 eV) ra1(k,i,lat)=4.2E-7*(300./te(k,i))**0.85 ! ! ra2: O2+ + e -> (15% O + O + 6.95 eV), (85% O + O + 4.98 eV) if (te(k,i) >= 1200.) then ra2(k,i,lat)=1.6E-7*(300./te(k,i))**0.55 else ! ra2: notes replace 2.7e-7 with 1.95e-7 ra2(k,i,lat)=2.7E-7*(300./te(k,i))**0.7 endif ! ! ra3: N2+ + e -> (10% N4S + N4S + 5.82 eV), (90% N4S + N2D + 3.44 eV) ra3(k,i,lat)=1.8E-7*(300./te(k,i))**0.39 ! ! beta1: N4S + O2 -> NO + O + 1.4 eV beta1(k,i,lat) = 1.5E-11*exp(-3600./tn(k,i)) ! ! beta3: N4S + NO -> N2 + O + 3.25 eV ! beta3(k,i,lat) = 2.5E-10*sqrt(tn(k,i)/300.)*exp(-600./tn(k,i)) beta3(k,i,lat) = 3.4E-11*sqrt(tn(k,i)/300.) ! ! beta5: N2D + e -> N4S + e + 2.38 eV beta5(k,i,lat) = 3.6E-10*sqrt(te(k,i)/300.) ! ! beta9: NO + hv/Ly-a -> NO+ + e (use JPL-94 for timegcm) beta9(k,i,lat) = 1.8e-12*exp(-1370./tn(k,i)) enddo ! i=lon0,lon1 enddo ! k=lev0,lev1-1 ! 5/21/04 btf: NaNQ at step 5 (80,4,20): ! Seeing NaNs in te at lat=10, k=31-44, i=30-32 ! Am retreating to settei at step 4.. ! if (istep==5) then ! do k=lev0,lev1 ! write(6,"(/,'chemrates: istep=',i2,' lat=',i2,' k=',i2, ! | ' rk19(k,:,lat)=',/,(6e12.4))") istep,lat,k,rk19(k,:,lat) ! write(6,"('rk20(k,:,lat)=',/,(6e12.4))") rk20(k,:,lat) ! write(6,"('te(k,lon0:lon1)=',/,(6e12.4))") te(k,lon0:lon1) ! enddo ! endif ! ! Additional reaction rates for time-gcm (not in tiegcm): ! (Neutral chemical reaction rates for mesospheric extension) ! do i=lon0,lon1 do k=lev0,lev1-1 ! ! Odd nitrogen beta10, 12: beta10(k,i,lat) = 3.5E-12*exp(250./tn(k,i)) beta12(k,i,lat) = 0. ! ! Metastables ! rkm1 (k,i,lat) = 1.8E-11*exp(110./tn(k,i)) rkm2a(k,i,lat) = 3.2E-11*exp( 70./tn(k,i)) rkm2b(k,i,lat) = rkm2a(k,i,lat) rkm6 (k,i,lat) = 7.4E-11*exp( 120./tn(k,i)) rkm9 (k,i,lat) = 3.6E-18*exp(-220./tn(k,i)) ! ! Oxygen-hydrogen interactions ! rkm20(k,i,lat) = 9.59E-34*exp(480./tn(k,i)) rkm21(k,i,lat) = 6.0E-34*(300./tn(k,i))**2.3 rkm22(k,i,lat) = rkm21(k,i,lat) rkm23(k,i,lat) = rkm21(k,i,lat) rkm24(k,i,lat) = 8.0E-12*exp(-2060./tn(k,i)) rkm25(k,i,lat) = 2.2E-11*exp( 120./tn(k,i)) rkm26(k,i,lat) = 1.5E-11*exp( 200./tn(k,i)) rkm27(k,i,lat) = 1.4E-12*exp(-2000./tn(k,i)) rkm28(k,i,lat) = 1.6E-11*exp(-4570 /tn(k,i)) rkm29(k,i,lat) = 1.6E-12*exp( -940./tn(k,i)) rkm30(k,i,lat) = 4.2E-12*exp( -240./tn(k,i)) rkm31(k,i,lat) = 4.8E-11*exp( 250./tn(k,i)) rkm32(k,i,lat) = 2.9E-12*exp( -160./tn(k,i)) rkm33(k,i,lat) = 5.5E-12*exp(-2000./tn(k,i)) rkm34(k,i,lat) = 1.1E-14*exp( -500./tn(k,i)) rkm35(k,i,lat) = 2.3E-13*exp( 600./tn(k,i)) rkm36(k,i,lat) = 5.7E-32*(300./tn(k,i))**1.6 rkm37(k,i,lat) = 1.4E-10*exp( -470./tn(k,i)) rkm39(k,i,lat) = 4.2E-10*exp( -950./tn(k,i)) rkm40(k,i,lat) = 8.3E-11*exp( -500./tn(k,i)) rkm41(k,i,lat) = 5.7E-32*(300./tn(k,i))**1.6 rkm42(k,i,lat) = 1.35E-13*exp( -100./tn(k,i)) rkm43(k,i,lat) = 2.4 E-12*exp(-1710./tn(k,i)) rkm44(k,i,lat) = 3.5 E-11*exp(-4550./tn(k,i)) ! ! D-region ion chemistry ! rin9 (k,i,lat) = 2.E-31*(300./tn(k,i))**4.4 rin10(k,i,lat) = 1.5E+6/(tn(k,i)**5.4*exp(2150./tn(k,i))) rin11(k,i,lat) = 7.E-30*(300./tn(k,i))**3. rin12(k,i,lat) = 3.11E+4/(tn(k,i)**4.*exp(4025./tn(k,i))) rin13(k,i,lat) = 1.E-27*(308./tn(k,i))**4.7 rin17(k,i,lat) = 2.E-31*(300./tn(k,i))**4.4 rin18(k,i,lat) = 1.5E+6/(tn(k,i)**5.4*exp(1800./tn(k,i))) rin19(k,i,lat) = 7.E-30*(300./tn(k,i))**3. rin20(k,i,lat) = 3.11E+4/(tn(k,i)**4.*exp(3335./tn(k,i))) rin21(k,i,lat) = 1.0E-27*(308./tn(k,i))**4.7 rin40(k,i,lat) = 2.6E-30*(300./tn(k,i))**3.2 ! ! Carbon reactions ! gam1(k,i,lat) = 2.3E-12*exp(-1700./tn(k,i)) gam3(k,i,lat) = 4.5E-31*(300./tn(k,i))**3 gam5(k,i,lat) = 4.2E-12*exp( 180./tn(k,i)) gam6(k,i,lat) = 3.3E-13*exp( 800./tn(k,i)) gam7(k,i,lat) = 2.2E-13*exp( 200./tn(k,i)) gam8(k,i,lat) = 3.8E-12*exp( 200./tn(k,i)) gam9(k,i,lat) = 3.9E-14*exp( -900./tn(k,i)) gam11(k,i,lat) = 3.4E-11*exp(-1600./tn(k,i)) gam12(k,i,lat) = 3.5E-12*exp( 140./tn(k,i)) gam13(k,i,lat) = 3.5E-13*sqrt(tn(k,i))*exp(-3650./tn(k,i)) gam14(k,i,lat) = 1.7E-33*exp(-1510./tn(k,i)) ! ! Chlorine reactions ! del1(k,i,lat) = 2.9E-11*exp( -260./tn(k,i)) del2(k,i,lat) = 3.0E-11*exp( +70./tn(k,i)) del3(k,i,lat) = 1.1E-11*exp(-1400./tn(k,i)) enddo ! k=lev0,lev1-1 enddo ! i=lon0,lon1 end subroutine chemrates_tdep end module chemrates_module