! module qrj_module ! ! Calculate ionization and heating rates. ! 10/10/06 btf: new qrj from Liying ! use params_module,only: nlevp1,nlonp4,zibot,dz use addfld_module,only: addfld implicit none ! ! Coefficients for wavelength scans in qrj. ! These are referenced in routines qrj, init_sflux, init_sigmas, and ! set_qrj_coeff (all in qrj_module.F). ! integer,parameter :: lmax=37 ! number of binsfrom 0.5A-1750A integer,parameter :: l1=15 ! number of bins in (1050-1750A) real :: | sigeuv(8,lmax), ! absorption coefficients | euveff(nlevp1), ! | rlmeuv(lmax), ! | feuv(lmax), ! | fsrc(l1), ! | sigsrc(l1), ! | rlmsrc(l1), ! | sigin4s(lmax), ! photoionization cross section if N | quench(4), ! | wave1(lmax), ! short bound of wave bins | wave2(lmax), ! long bound of wave bins | sfmin(lmax), ! reference solar minimum flux of EUVAC model | afac(lmax), ! The A factor of EUVAC model | sflux(lmax) ! Solar flux for each time step ! ! Branching ratios for photon are branching ratios from photon absorption rate. ! Branching ratios for photolectron are branching ratios from photoionization rate. real:: | BPhotonI(3,lmax), ! Photoionization branching ratio for major species | BElectronI(3,lmax), ! Photoelectron ionization branching ratio ! for major species | brop2pPh(lmax), ! photoionization branching ratio for O+(2p) | brop2dPh(lmax), ! photoionization branching ratio for O+(2d) | brop4sPh(lmax), ! photoionization branching ratio for O+(4s) | bro2DPh(lmax), ! photodissociation branching ratio for O2 | brn2DPh(lmax), ! photodissociation branching ratio for N2 | bro2DIPh(lmax), ! Photon dissociative ionization branching ratio for O2 | brn2DIPh(lmax), ! Photon dissociative ionization branching ratio for N2 | brop2pEl(lmax), ! electron impact ionization branching ratio for O+(2p) | brop2dEl(lmax), ! electron impact ionization branching ratio for O+(2d) | brop4sEl(lmax), ! electron impact ionization branching ratio for O+(4s) | bro2DIEl(lmax), ! Photoelectron dissociative ionization branching ratio for O2 | brn2DIEl(lmax), ! Photoelectron dissociative ionization branching ratio for N2 | brn2DEl(lmax), ! Photoelectron dissociation branching ratio for N2 | bro2DEl(lmax) ! Photoelectron dissociation branching ratio for O2 ! ! Heating and ionization terms set by qrj, and used by other routines. ! These are allocated for task subdomains by alloc_q (called from allocdata) ! e.g.: allocate(qtotal(nlevp1,lon0:lon1,lat0:lat1),stat=istat) ! real,dimension(:,:,:),allocatable,save :: | qtotal, ! total heating ! F(NQ) | qop2p, ! o+(2p) ionization ! F(NQOP2P) | qop2d, ! o+(2d) ionization ! F(NQOP2D) | qo2p, ! o2+ ionization ! F(NQO2P) | qo2j, ! ! F(NQO2J) zeroed out | qop, ! o+ ionization ! F(NQOP) | qn2p, ! n2+ ionization ! F(NQN2P) | qnp, ! n+ ionization ! F(NQNP) | qnop, ! no+ ionization ! F(NQNOP) | qnoplya,! no+ nighttime ionization by lyman-alpha (XJNOPN) | qnolya ! no ionization by lyman-alpha (XJNOP) ! ! Photodissociation terms set by qrj, and used by other routines. real,dimension(:,:,:),allocatable :: | pdo2, ! total o2 dissociation (photon and electron), ! in frequency (s^-1) (F(NRJ)) | pdo2d, ! o2(d) photodissociation (F(NJO2D)) | pdn2, ! total N2 dissociation (photon and electron), ! in N production rate (times 2) (cm^-3 s^-1) (F(NQTEF)) | pdco2t, ! total photodissociation of co2 (F(NJCO2T)) | pdco2d, ! total photodissociation of co2 (XJCO2D) | pdh2ol, ! lyman-alfa dissociation of h2o (F(NJH2OL)) | pdh2ot, ! total photodissociation of h2o (F(NJH2OT)) | pdo3d, ! photodissociation of o3(d) (F(NJO3D)) | pdo3p, ! photodissociation of o3(p) (F(NJO3P)) | pdch4a, ! partial photodissociation of ch4 (XJCH4A) | pdch4b, ! partial photodissociation of ch4 (XJCH4B) | pdch4t, ! total photodissociation of ch4 (F(NJCH4T)) | pdnoeuv,! photodissociation of no by euv (DNOEUV) | pdnosrb,! photodissociation of no by srb (XJNO) | pdh2o2, ! XJH2O2 | pdch2oa,! XJCH2OA | pdch2ob,! XJCH2OB | pdn2o, ! XJN2O | pdho2, ! XJHO2 | pdno2, ! XJNO2 | pdch3oo,! XJCH3OO | pd762 ! XJ762 ! ! Attenuation of lyman-alfa at night (s.a., qinite.F) ! real,dimension(:,:,:),allocatable :: attlya ! SJNMLYA ! real,external :: expo contains !----------------------------------------------------------------------- subroutine qrj(sco2,sco1,scn2,sco3,scno,vo2,tn,no,o2,o1,he,o3,n4s, | o21d,xnmbari,cp,lev0,lev1,lon0,lon1,lat) ! ! Calculate heating, ionization, and dissociation rates. ! (Note that if check_exp is set true (cons.F), then qrj_exp is called ! instead of this routine) ! Also note input arg o1 is actually ox. ! use input_module,only: f107,f107a use init_module,only: sfeps,istep use cons_module,only: t0,expz,avo,rmassinv_n4s,rmassinv_no, | rmassinv_o2,rmassinv_o1,rmassinv_o3,rmassinv_n2,expzmid,gask, | amu,expzmid_inv,p0 use lbc,only: fb,b use chemrates_module,only: disn2p,beta9 use fields_module,only: tlbc use o2srbc ! module in o2srbc.f for O2 heating/dissoc in Shumann-Runge use chapman_module,only: idn ! ! Args: integer,intent(in) :: lev0,lev1,lon0,lon1,lat real,dimension(lev0:lev1,lon0-2:lon1+2),intent(in) :: | sco2,sco1,scn2,sco3,scno, ! chapman slant column integrals | vo2, ! chapman vertical integral | tn,no,o2,o1,he,o3,n4s, ! tn and species mass mixing ratios | o21d, | xnmbari, ! p0*e(-z)*barm/kT at interfaces | cp ! specific heat ! ! VT vampir tracing: ! #ifdef VT #include #endif ! ! Local: integer :: k,i,l,nlevs,ier integer :: levmin ! levmin=-7, It is the level above which ionization ! are calculated. real,parameter :: | do2=8.203E-12 , | do22=1.1407E-11 , | aband=0.143 , ! shumann-runge | bband=9.64E8 , ! shumann-runge | cband=9.03E-19 , ! shumann-runge | e3=0.33 , | hc = 1.9845E-16 , ! C(60) | zpmin = -7.0 ! was zkmin real :: rlmeuvinv(lmax),factor,fmin,fmax real,dimension(lev0:lev1,lon0:lon1) :: | o2i, ! o2 at interfaces (s1) | o1i, ! o at interfaces (s2) | o3i, ! o3 at interfaces (s1) | hei, ! he at interfaces | n2i, ! n2 at interfaces (s3) | n4si, ! n4s at interfaces (s4) | o21di, ! o21d at interfaces (s3) | tni ! tn at interfaces (s6) real,dimension(lev0:lev1,lon0:lon1) :: | quenchfac, ! (s8) | sigchap, ! (s9) | p3f, ! (s7) | pdnolya, ! photodissociation of NO by lyman-alpha (SJMLYA) | pdo2lya, ! photodissociation of O2 by lyman-alpha (SJO2LYA) | sco3t ! temporary chapman integral ! ! temporary loop variables real:: | absorp_o, ! photoabsorption frequency of O | absorp_o2, ! photoabsoption frequency of o2 | absorp_n2, ! photoabsoption frequency of n2 | ioniz_o, ! photoionization frequency of o | ioniz_o2, ! photoionization frequency of o2 | ioniz_n2, ! photoionization frequency of n2 | htfac ! hc/wavelength, for calculation of heating real,dimension(lev0:lev1,lon0:lon1) :: | di_o2, ! total dissociative ionization frequency of o2 | di_n2, ! total dissociative ionization frequency of n2 | mn_o2, ! transfer mass density to number density (O2) | mn_o1, ! transfer mass density to number density (O) | mn_n2, ! transfer mass density to number density (N2) | mn_n ! transfer mass density to number density (N) ! real,dimension(lev0:lev1,lon0:lon1) :: | sum1, ! sum(o2,o,n2)(sigma*chapman) (s5) | sum2 ! sum(o2,o,n2)(sigma*psi/rmass) (s6) ! ! Declarations for O2 photolysis and heating in Shumann-Runge Continuum ! (SRC) and Shumann-Runge Bands (SRB) (see o2srbc.f). ! 12/11/00: declarations for o2srbc as per kibo12: ! 8/5/04: Sub jo2h2osrb returns xdo2srb and dh2osrb, but only dh2osrb ! is being used (do2srb from sub mkdo2srb is used instead of xdo2srb) real :: | xno2(lev0:lev1), ! o2 column density (cm-2) | rho(lev0:lev1), ! total density | do2src(lev0:lev1,lon0:lon1), ! o2 dissoc in SRC (mkdo2src output) | ho2src(lev0:lev1,lon0:lon1), ! o2 heating in SRC (mkho2src output) | do2srb(lev0:lev1,lon0:lon1), ! o2 dissoc in SRB (mkdo2srb output) | ho2srb(lev0:lev1,lon0:lon1), ! o2 heating in SRB (mkdo2srb output) | xdo2srb(lev0,lev1,lon0:lon1), ! alternate o2 srb dissoc (sub jo2h2osrb) | dh2osrb(lev0:lev1,lon0:lon1) ! h2o dissoc in SRB (jo2h2osrb output) real :: sfac,temp1,temp2 ! ! Heating diagnostics: real,dimension(lev0:lev1,lon0:lon1) :: | qtotal1,qtotal2,qtotal3,qtotal4,qtotal5,qtotal6,qtotal7 integer :: i0,i1,nk,nkm1 ! ! If mks > 0, ho2src or ho2srb are returned in deg K/sec (mks) ! If mks <= 0, ho2src or ho2srb are returned in ergs/gm-1/s-1 (cgs) integer,parameter :: mks=0 ! units flag for ho2src ! logical,parameter :: | debug=.false. ! insert print statements ! if (debug) write(6,"('Enter qrj: lat=',i3,' lon0,1=',2i3)") | lat,lon0,lon1 ! #ifdef VT ! code = 118 ; state = 'qrj' ; activity='ModelCode' call vtbegin(118,ier) #endif ! ! Exec: nlevs = lev1-lev0+1 i0=lon0 ; i1=lon1 ; nk=lev1-lev0+1 ; nkm1=nk-1 ! write(6,"('qrj: lat=',i2,' lon0,1=',2i3,' idn(lon0:lon1,lat)=',/, ! | (15i3))") lat,lon0,lon1,idn(lon0:lon1,lat) ! call addfld('SCO2',' ',' ',sco2(:,lon0:lon1), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('SCO1',' ',' ',sco1(:,lon0:lon1), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('SCN2',' ',' ',scn2(:,lon0:lon1), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('XNMBARI',' ',' ',xnmbari(:,lon0:lon1), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! ! calculate inverse of wave length ! do i=1,lmax rlmeuvinv(i) = 1./rlmeuv(i) enddo ! ! O2,O,N4S at interface levels: do i=lon0,lon1 do k=lev0,lev1-1 o2i (k+1,i) = 0.5*(o2(k,i)+o2(k+1,i)) o1i (k+1,i) = 0.5*(o1(k,i)+o1(k+1,i)) hei (k+1,i) = 0.5*(he(k,i)+he(k+1,i)) n4si(k+1,i) = 0. enddo enddo ! ! Bottom boundary: do i=lon0,lon1 o2i(1,i) = .5*((b(i,1,1)+1.)*o2(1,i)+b(i,1,2)*o1(1,i)+ | fb(i,1,lat)) o1i(1,i) = .5*(b(i,2,1)*o2(1,i)+(b(i,2,2)+1.)*o1(1,i)+ | fb(i,2,lat)) hei(1,i) = .5*(b(i,3,1)*o2(1,i)+(b(i,3,2)+1.)*o1(1,i)+ | fb(i,3,lat)) n4si(1,i) = 0. enddo ! ! N2 at interfaces: do k=lev0,lev1 n2i(k,:) = 1.-o2i(k,:)-o1i(k,:)-hei(k,:) enddo ! ! Zero diffs, 7/30/03: ! call addfld('O2I' ,' ',' ',o2i, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('O1I' ,' ',' ',o1i, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('N2I' ,' ',' ',n2i, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! ! calculate variables for transferring mass density to number density do i=lon0,lon1 do k=lev0,lev1 mn_o2(k,i)=o2i(k,i)*rmassinv_o2 mn_o1(k,i)=o1i(k,i)*rmassinv_o1 mn_n2(k,i)=n2i(k,i)*rmassinv_n2 mn_n(k,i)=n4si(k,i)*rmassinv_n4s enddo enddo ! ! Initialize terms on current processor subdomain: ! (global module data for use by other routines) ! (in tgcm24, some of these inits occur in sub clearnq) ! qtotal(lev0:lev1,lon0:lon1,lat) = 0. ! NQ qop2p (lev0:lev1,lon0:lon1,lat) = 0. ! NQOP2P qop2d (lev0:lev1,lon0:lon1,lat) = 0. ! NQOP2D qo2p (lev0:lev1,lon0:lon1,lat) = 0. ! NQO2P qo2j (lev0:lev1,lon0:lon1,lat) = 0. ! NQO2J qop (lev0:lev1,lon0:lon1,lat) = 0. ! NQOP qn2p (lev0:lev1,lon0:lon1,lat) = 0. ! NQN2P qnp (lev0:lev1,lon0:lon1,lat) = 0. ! NQNP qnop (lev0:lev1,lon0:lon1,lat) = 0. ! NQNOP ! pdo2 (lev0:lev1,lon0:lon1,lat) = 0. ! NRJ pdn2 (lev0:lev1,lon0:lon1,lat) = 0. ! NQTEF disn2p (lev0:lev1,lon0:lon1,lat) = 0. ! pdco2t (lev0:lev1,lon0:lon1,lat) = 0. ! NJCO2T pdh2ol (lev0:lev1,lon0:lon1,lat) = 0. ! NJH2OL pdh2ot (lev0:lev1,lon0:lon1,lat) = 0. ! NJH2OT pdo3d (lev0:lev1,lon0:lon1,lat) = 0. ! NJO3D pdo3p (lev0:lev1,lon0:lon1,lat) = 0. ! NJO3P pdo2d (lev0:lev1,lon0:lon1,lat) = 0. ! NJO2D (added 5/25/04) pdch4t (lev0:lev1,lon0:lon1,lat) = 0. ! NJCH4T pdnoeuv(lev0:lev1,lon0:lon1,lat) = 0. ! DNOEUV ! ! Initialize local arrays di_o2=0. di_n2=0. ! ! Summation over wavelength. Define parameters at task subdomains, and ! vertical column only from zpmin to top. ! levmin = int((zpmin-zibot)/dz)+1 ! zpmin == -17.0 for timegcm ! write(6,"('qrj: lat=',i3,' zpmin=',f6.2,' levmin=',i3)") ! | lat,zpmin,levmin ! ! Summation over wavelength: do l=l1+1,lmax ! from 0.5A to 1050A sum1(:,:) = 0. ! sum(o2,o,n2)(sigma*chapman) sum2(:,:) = 0. ! sum(o2,o,n2)(sigma*psi/rmass) htfac=hc*rlmeuvinv(l) do i=lon0,lon1 do k=levmin,lev1 sum1(k,i) = sum1(k,i)+sigeuv(1,l)*sco2(k,i)+ | sigeuv(2,l)*sco1(k,i)+ | sigeuv(3,l)*scn2(k,i) ! sum1(k,i) = feuv(l)*exp(-sum1(k,i)) sum2(k,i) = sum2(k,i)+sigeuv(1,l)*mn_o2(k,i)+ | sigeuv(2,l)*mn_o1(k,i)+ | sigeuv(3,l)*mn_n2(k,i) enddo enddo ! ! Longitude and column domain of current process: do i=lon0,lon1 do k=levmin,lev1 ! absorption/ionization frequency for the three major species (O2, O, and N2) absorp_o2= sum1(k,i)*sigeuv(1,l) absorp_o= sum1(k,i)*sigeuv(2,l) absorp_n2= sum1(k,i)*sigeuv(3,l) ioniz_o2=absorp_o2*BPhotonI(1,l) ioniz_o=absorp_o*BPhotonI(2,l) ioniz_n2=absorp_n2*BPhotonI(3,l) ! ! ionization/dissociative ionization frequency (s^-1) di_o2(k,i)=di_o2(k,i)+absorp_o2*bro2DIPh(l)+ | ioniz_o2*bro2DIEl(l) di_n2(k,i)=di_n2(k,i)+absorp_n2*brn2DIPh(l)+ | ioniz_n2*brn2DIEl(l) qnp(k,i,lat) = qnp(k,i,lat)+sigin4s(l)*sum1(k,i) qn2p(k,i,lat) = qn2p(k,i,lat)+ioniz_n2+ | ioniz_n2*BElectronI(3,l) disn2p(k,i,lat) = disn2p(k,i,lat)+ioniz_n2*BElectronI(3,l) qo2p(k,i,lat) = qo2p(k,i,lat)+ioniz_o2+ | ioniz_o2*BElectronI(1,l) qop2p(k,i,lat) = qop2p(k,i,lat)+absorp_o*brop2pPh(l)+ | ioniz_o*brop2pEl(l) qop2d(k,i,lat) = qop2d(k,i,lat)+absorp_o*brop2dPh(l)+ | ioniz_o*brop2dEl(l) qop(k,i,lat) = qop(k,i,lat)+absorp_o*brop4sPh(l)+ | ioniz_o*brop4sEl(l) ! ! total dissociation and EUV heating pdo2(k,i,lat) = pdo2(k,i,lat)+(absorp_o2*bro2DPh(l) | +ioniz_o2*bro2DEl(l)) pdco2t (k,i,lat) = pdco2t (k,i,lat)+sum1(k,i)*sigeuv(5,l) ! f(njco2t) pdnoeuv(k,i,lat) = pdnoeuv(k,i,lat)+sum1(k,i)*sigeuv(7,l) ! dnoeuv pdo3d (k,i,lat) = pdo3d (k,i,lat)+sum1(k,i)*sigeuv(4,l) ! f(njo3d) pdch4t (k,i,lat) = pdch4t (k,i,lat)+sum1(k,i)*sigeuv(8,l) ! f(njch4t) pdn2(k,i,lat) = pdn2(k,i,lat)+absorp_n2*brn2DPh(l) | +ioniz_n2*brn2DEl(l) qtotal(k,i,lat) = qtotal(k,i,lat)+htfac* | sum1(k,i)*sum2(k,i) enddo enddo enddo ! l=l1+1,lmax ! ! call addfld('QOP2P' ,' ',' ',qop2p(lev0:lev1,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('QOP2D' ,' ',' ',qop2d(lev0:lev1,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('SUM1' ,' ',' ',sum1, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('Q0' ,' ',' ',qtotal(:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('R_QRJ',' ',' ',r,'lev',lev0,lev1,'lon',lon0,lon1,lat) ! ! Multiply Q by efficiency factor: do i=lon0,lon1 do k=lev0,lev1 qtotal(k,i,lat) = qtotal(k,i,lat)*euveff(k)*avo enddo enddo ! call addfld('Q1' ,' ',' ',qtotal(:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! ! transfer frequency (s^-1) to production rates (cm^-3 s^-1) do i=lon0,lon1 do k=lev0,lev1 mn_o2(k,i)=mn_o2(k,i)*xnmbari(k,i) mn_n2(k,i)=mn_n2(k,i)*xnmbari(k,i) mn_o1(k,i)=mn_o1(k,i)*xnmbari(k,i) mn_n(k,i)=mn_n(k,i)*xnmbari(k,i) di_o2(k,i)=di_o2(k,i)*mn_o2(k,i) di_n2(k,i)=di_n2(k,i)*mn_n2(k,i) qo2p(k,i,lat) = qo2p(k,i,lat)*mn_o2(k,i)-di_o2(k,i) qn2p(k,i,lat) = qn2p(k,i,lat)*mn_n2(k,i)-di_n2(k,i) disn2p(k,i,lat)=disn2p(k,i,lat)*mn_n2(k,i) qnp(k,i,lat) = qnp(k,i,lat)*mn_n(k,i)+di_n2(k,i) qop(k,i,lat) = qop(k,i,lat)*mn_o1(k,i)+0.54*di_o2(k,i) qop2p(k,i,lat) = qop2p(k,i,lat)*mn_o1(k,i)+0.22*di_o2(k,i) qop2d(k,i,lat) = qop2d(k,i,lat)*mn_o1(k,i)+0.24*di_o2(k,i) pdn2(k,i,lat) = 2.*pdn2(k,i,lat)*mn_n2(k,i)+di_n2(k,i) enddo enddo ! ! lyman_alpha ! do i=lon0,lon1 do k=lev0,lev1 ! pdnolya: photodissociation of NO by lyman-alpha (local) (SJMLYA) pdnolya(k,i) = (0.68431 *exp(-8.22114e-21*sco2(k,i))+ | 0.229841 *exp(-1.77556e-20*sco2(k,i))+ | 0.0865412*exp(-8.22112e-21*sco2(k,i)))* | sflux(12) ! ! pdo2lya: photodissociation of O2 by lyman-alpha (local) (SJO2LYA) pdo2lya(k,i) = (6.0073e-21*exp(-8.2166e-21 *sco2(k,i))+ | 4.28569e-21*exp(-1.63296e-20*sco2(k,i))+ | 1.28059e-20*exp(-4.85121e-17*sco2(k,i))) | *sflux(12) ! ! qnolya: no ionization by lyman-alpha (XJNOP) qnolya(k,i,lat) = 2.02e-18*pdnolya(k,i) ! ! Calc of qnoplya moved to qinite. ! ! attlya: attenuation of lyman-alpha at night (SJNMLYA) ! attlya(k,i,lat) = (0.68431 *exp(-8.22114e-21*vo2(k,i))+ ! | 0.229841 *exp(-1.77556e-20*vo2(k,i))+ ! | 0.0865412*exp(-8.22112e-21*vo2(k,i))) ! | *sfeps ! ! qnoplya: no+ nighttime ionization by lyman-alpha (XJNOPN) ! qnoplya(k,i,lat) = 2.02e-18*sflux(12)/100.*attlya(k,i,lat) ! enddo ! k=lev0,lev1 enddo ! i=lon0,lon1 ! ! Total NO dissociation pdnosrb (XJNO) is returned by sub jno. ! allocate(pdnosrb(nlevp1,lon0:lon1,lat0:lat1),stat=istat) ! call jno( | pdnosrb(:,lon0:lon1,lat), o2(:,lon0:lon1), o1(:,lon0:lon1), | he(:,lon0:lon1), xnmbari(:,lon0:lon1), | sco2(:,lon0:lon1),sco3(:,lon0:lon1), | scno (:,lon0:lon1),lev0,lev1,lon0,lon1,lat) ! call addfld('pdnosrb',' ',' ',pdnosrb(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! ! Add no ionization to qnop: do i=lon0,lon1 qnop(lev0,i,lat) = qnop(lev0,i,lat)+(qnolya(lev0,i,lat)+ | pdnoeuv(lev0,i,lat))*no(lev0,i)*xnmbari(lev0,i)*rmassinv_no enddo do i=lon0,lon1 do k=lev0+1,lev1 qnop(k,i,lat) = qnop(k,i,lat)+(qnolya(k,i,lat)+ | pdnoeuv(k,i,lat))*.5*(no(k,i)+no(k-1,i))*xnmbari(k,i)* | rmassinv_no enddo enddo ! ! Minimum values for ionization rates: do i=lon0,lon1 do k=lev0,lev1 if (qo2p (k,i,lat) < 1.e-20) qo2p (k,i,lat) = 1.e-20 ! (f(nqo2p)) if (qop (k,i,lat) < 1.e-20) qop (k,i,lat) = 1.e-20 ! (f(nqop)) if (qop2p(k,i,lat) < 1.e-20) qop2p(k,i,lat) = 1.e-20 ! (f(nqop2p)) if (qop2d(k,i,lat) < 1.e-20) qop2d(k,i,lat) = 1.e-20 ! (f(nqop2d)) if (qn2p (k,i,lat) < 1.e-20) qn2p (k,i,lat) = 1.e-20 ! (f(nqn2p)) if (qnp (k,i,lat) < 1.e-20) qnp (k,i,lat) = 1.e-20 ! (f(nqnp)) if (pdn2 (k,i,lat) < 1.e-20) pdn2 (k,i,lat) = 1.e-20 ! (f(nqtef)) if (qnop (k,i,lat) < 1.e-20) qnop (k,i,lat) = 1.e-20 ! (f(nqnop)) enddo enddo ! ! tn at interfaces (S6): do i=lon0,lon1 ! tni(1,i) = tn(lev1,i) ! tn bottom boundary is stored in top slot tni(lev0,i) = tlbc(i,lat) do k=lev0+1,lev1-1 tni(k,i) = .5*(tn(k-1,i)+tn(k,i)) enddo tni(lev1,i) = tn(lev1-1,i) ! nlevp1 <- nlev enddo ! ! Quench: do i=lon0,lon1 do k=lev0,lev1 factor = avo*p0/gask*expz(1)*expzmid**(2*k-3) quenchfac(k,i) = factor/((o2i(k,i)*rmassinv_o2+ ! s8 | o1i(k,i)*rmassinv_o1+ | n2i(k,i)*rmassinv_n2)*tni(k,i)) quenchfac(k,i) = quenchfac(k,i)* | (quench(1)*n2i(k,i)*rmassinv_n2+ | quench(2)*o2i(k,i)*rmassinv_o2) quenchfac(k,i) = quench(3)*quenchfac(k,i)/ | (quench(4)+quenchfac(k,i)) enddo enddo ! ! Zero diffs, 8/03: ! call addfld('qo2p',' ',' ',qo2p(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qn2p',' ',' ',qn2p(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qnp',' ',' ',qnp(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdn2',' ',' ',pdn2(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qop2pd',' ',' ',qop2pd, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qop2p',' ',' ',qop2p(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qop2d',' ',' ',qop2d(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qop',' ',' ',qop(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdnosrb',' ',' ',pdnosrb(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdnolya',' ',' ',pdnolya, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdo2lya',' ',' ',pdo2lya, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qnolya',' ',' ',qnolya(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('attlya',' ',' ',attlya(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qnoplya',' ',' ',qnoplya(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('quench',' ',' ',quenchfac, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qnop',' ',' ',qnop(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! ! Summation over wave length: do l=1,l1 ! from 1050A to 1750A do i=lon0,lon1 do k=lev0,lev1 sigchap(k,i) = fsrc(l)*exp(-sigsrc(l)*sco2(k,i)) pdo3d (k,i,lat) = pdo3d (k,i,lat)+sigeuv(4,l)*sigchap(k,i) pdco2t(k,i,lat) = pdco2t(k,i,lat)+sigeuv(5,l)*sigchap(k,i) pdh2ot(k,i,lat) = pdh2ot(k,i,lat)+sigeuv(6,l)*sigchap(k,i) qnop (k,i,lat) = qnop (k,i,lat)+sigeuv(7,l)*sigchap(k,i) enddo ! k=lev0,lev1 enddo ! i=lon0,lon1 enddo ! l=1,l1 ! ! call addfld('pdo3d',' ',' ',pdo3d(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdco2t',' ',' ',pdco2t(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdh2ot',' ',' ',pdh2ot(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qnop',' ',' ',qnop(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! ! Most code from here to end of qrj is time-gcm only (not in tiegcm): ! ! Calculate O2 photolysis and heating in Shumann-Runge Continuum (src) ! Shumann-Runge Bands (srb): ! do i=lon0,lon1 ho2src(:,i) = 0. do2src(:,i) = 0. do k=lev0,lev1 rho(k) = (o2i(k,i)*xnmbari(k,i)+ | o1i(k,i)*xnmbari(k,i)+ | n2i(k,i)*xnmbari(k,i))*amu ! gm/cm3 xno2(k) = o2i(k,i)*rmassinv_o2*xnmbari(k,i) enddo ! lev0,lev1 ! ! Pass columns at each grid point to the o2src routines. ! mkdo2src returns do2src (dissociation), mkho2src returns ho2src (heating) ! subroutine mkdo2src(sco2,f107d,do2src,nlev) ! subroutine mkho2src(sco2,xno2,rho,cp,f107d,ho2src,nlev,mks,lat) ! call mkdo2src(sco2(:,i),f107,do2src(:,i),nlevs) call mkho2src(sco2(:,i),xno2,rho,cp(:,i),f107,ho2src(:,i), | nlevs,mks) ! ! mkdo2srb returns do2srb and ho2srb columns: call mkdo2srb(sco2(:,i),xno2,rho,cp(:,i),f107,sfeps, | do2srb(:,i),ho2srb(:,i),nlevs,mks) enddo ! i=lon0,lon1 ! ! Add src and srb to total heating and o2 heating and dissociation: do i=lon0,lon1 do k=lev0,lev1 qtotal(k,i,lat) = qtotal(k,i,lat)+ho2src(k,i)+ho2srb(k,i) ! f(nq) pdo2(k,i,lat) = pdo2(k,i,lat)+do2src(k,i)+do2srb(k,i) ! f(nrj) pdo2d(k,i,lat) = pdo2d(k,i,lat)+do2src(k,i) ! f(njo2d) enddo ! k=lev0,lev1 enddo ! i=lon0,lon1 ! qtotal3(:,:) = qtotal(:,:,lat) ! call addfld('QTOTAL3',' ',' ',qtotal3, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('do2src',' ',' ',do2src, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('ho2src',' ',' ',ho2src, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('do2srb',' ',' ',do2srb, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('ho2srb',' ',' ',ho2srb, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qtotal',' ',' ',qtotal(:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdo2' ,' ',' ',pdo2 (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdo2d' ,' ',' ',pdo2d (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! ! Contributions from ozone dissociation from the Herbzberg, Hartley, ! Huggins and chapius bands. ! Contributions from solar lyman-alpha, SRB and Herzberg to O2 ! dissociation and heating. ! sfac = 1.+0.11*(f107-65.)/165. do i=lon0,lon1 do k=lev0,lev1-1 o3i (k+1,i) = 0.5*(o3 (k,i)+o3 (k+1,i)) ! s1 o2i (k+1,i) = 0.5*(o2 (k,i)+o2 (k+1,i)) ! s2 o21di(k+1,i) = 0.5*(o21d(k,i)+o21d(k+1,i)) ! s3 enddo ! k=lev0,lev1-1 o3i (lev0,i) = o3 (lev0,i)**1.5/sqrt(o3 (lev0+1,i)) o2i (lev0,i) = o2 (lev0,i)**1.5/sqrt(o2 (lev0+1,i)) o21di(lev0,i) = o21d(lev0,i)**1.5/sqrt(o21d(lev0+1,i)) enddo ! i=lon0,lon1 ! ! Sub jo2h2osrb calculates o2 and h2o dissociation at SRB: ! dh2osrb is dissociation of h2o at SRB. ! xdo2srb is an alternate calculation of o2 dissociation at SRB, ! and is not used in this version. (see sub mkdo2srb for o2 ! dissociation, called above) ! call jo2h2osrb(xdo2srb,dh2osrb,sco2,sco3,lev0,lev1,lon0,lon1,lat) ! call addfld('JO2SRB' ,' ',' ',xdo2srb, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('JH2OSRB',' ',' ',dh2osrb, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! ! Scan subdomain: do i=lon0,lon1 do k=lev0,lev1 sco3t(k,i) = sco3(k,i) ! s9 if (sco3t(k,i) < 1.e+5) sco3t(k,i) = 1.e+5 ! ! Dissociation and heating of O3 from lyman-alpha: pdo3d (k,i,lat) = pdo3d (k,i,lat)+2.27e-17*pdnolya(k,i) qtotal(k,i,lat) = qtotal(k,i,lat)+2.27e-17*pdnolya(k,i)* | o3i(k,i)*rmassinv_o3*avo*9.944e-12 ! qtotal4(k,i) = qtotal(k,i,lat) ! diagnostic temp1 = 1.0e-3*exp(-1.5577e-13*sco3t(k,i)**0.6932) if (sco3t(k,i) < 1.6e+20) | temp1 = 1.04e-2*exp(-1.0217e-6*sco3t(k,i)**0.3587) ! ! Hartley bands of O3: pdo3d(k,i,lat) = pdo3d(k,i,lat)+ | (temp1*0.68+exp(-1.4912e-14*sco2 (k,i)**0.5298)* | (4.053e-4 *exp(-8.1381e-16*sco3t(k,i)**0.8856)+ | 4.7e-6 *exp(-1.5871e-14*sco3t(k,i)**0.7665))* | 1.085 *exp(-1.4655e-25*sco2 (k,i)**1.0743)*0.68)* | sfeps ! ! Chappius and Huggins bands: pdo3p(k,i,lat) = pdo3p(k,i,lat)+ | (4.5e-4*exp(-3.4786e-11*sco2(k,i) | **0.3366)*exp(-1.0061e-20*sco3t(k,i)**0.9719) | +(7.5e-4 *exp(-2.7663e-19*sco3t(k,i)**1.0801) | +2.5e-4/(1.+1.5772e-18 *sco3t(k,i)**0.9516)) | *exp(-1.0719e-10*sco2(k,i)**0.3172))*sfeps ! ! Add ozone heating to total heat: temp1 = exp(-5.50e-24*sco2(k,i)-6.22E-18*sco3(k,i)) temp2 = exp(-1.34e-26*sco2(k,i)-1.66E-18*sco3(k,i)) qtotal(k,i,lat) = qtotal(k,i,lat)+ | ((5.512e-16*exp(-3.16E-21*sco3(k,i)) | +(41.21e-16*exp(-8.94E-19*sco3(k,i)) | + 10.90e-16*exp(-1.09E-19*sco3(k,i)) | + 52.80e-16*exp(-3.07E-18*sco3(k,i))) | +(69.60e-16*exp(-9.65E-18*sco3(k,i)) ! | + 9.39e-16*exp(-4.79E-18*sco3(k,i))) | + 9.39e-16*exp(-4.79E-18*sco3(k,i)))*1.6 | +(14.94e-16*temp1+2.76E-16*temp2)) | *o3i(k,i)*rmassinv_o3*avo)*sfeps ! qtotal5(k,i) = qtotal(k,i,lat) ! diagnostic ! ! Dissociation and heating of O2 from lyman-alpha: pdo2 (k,i,lat) = pdo2 (k,i,lat)+pdo2lya(k,i) pdo2d (k,i,lat) = pdo2d (k,i,lat)+pdo2lya(k,i)*0.53 qtotal(k,i,lat) = qtotal(k,i,lat)+pdo2lya(k,i)* | o2i(k,i)*rmassinv_o2*avo*8.13816E-12 ! qtotal6(k,i) = qtotal(k,i,lat) ! diagnostic ! ! Dissociation and heating of O2 from Herzberg: temp1 = | 7.4e-10*exp(-2.5352e-15*sco2 (k,i)**0.6288)* | exp(-4.6661E-18*sco3t(k,i)**0.9538)* | 0.9*exp(-1.4110E-25*sco2 (k,i)**1.0667) pdo2(k,i,lat) = pdo2(k,i,lat)+temp1*sfeps qtotal(k,i,lat) = qtotal(k,i,lat)+temp1* | o2i(k,i)*rmassinv_o2*avo*0.46458e-12*sfeps ! qtotal7(k,i) = qtotal(k,i,lat) ! diagnostic ! ! O(1D) photodixxociation production from CO2. ! From lyman-alpha: pdco2t(k,i,lat) = pdco2t(k,i,lat)+8.14e-20*pdnolya(k,i) pdco2d(k,i,lat) = pdco2t(k,i,lat) ! ! Total co2 photodissociation: if (sco2(k,i) >= 6.3e+22) then pdco2t(k,i,lat) = pdco2t(k,i,lat)+(2.0e-11* | exp(-1.9087e-20*sco2 (k,i)**0.8675)* | exp(-2.9570e-14*sco3t(k,i)**0.7582)* | exp(-5.9648e-15*sco2 (k,i)**0.6172))*sfeps else pdco2t(k,i,lat) = pdco2t(k,i,lat)+(8.5E-9* | exp(-3.4368E-3 *sco2 (k,i)**0.1456)* | exp(-2.9570E-14*sco3t(k,i)**0.7582)* | exp(-5.9648E-15*sco2 (k,i)**0.6172))*sfeps endif ! ! Total h2o photodissociation: pdh2ol(k,i,lat) = 1.53e-17*pdnolya(k,i) ! ! 8/5/04 btf: Use dh2osrb from sub jo2h2osrb. ! temp1=exp(-1.e-7*sco2(k,i)**0.35)*sfeps ! pdh2ot(k,i,lat) = pdh2ot(k,i,lat)+pdh2ol(k,i,lat)+ ! | sfac*1.2e-6*temp1 pdh2ot(k,i,lat) = (pdh2ot(k,i,lat)+pdh2ol(k,i,lat)+ | dh2osrb(k,i))*sfeps ! ! Total o2(1d) ionization: qo2p(k,i,lat) = qo2p(k,i,lat)+ | (0.549e-9*exp(-2.406e-20*sco2(k,i))+2.6e-9*exp(-8.508E-20* | sco2(k,i)))*o21di(k,i)*rmassinv_o2*xnmbari(k,i)*sfeps ! ! Total CH4 photodissociation: pdch4a(k,i,lat) = 1.3e-6*exp(-8.4899e-21*sco2(k,i)**1.0034)* | exp(-3.1398e-17*sco3t(k,i)**1.0031)* | exp(-0.2776/sco3t(k,i)**0.8240) pdch4b(k,i,lat) = 3.888e-6*exp(-4.7152e-21*sco2(k,i)**1.0153)* | exp(-4.0876e-16*sco3t(k,i)**0.9347)* | exp(-0.4976/sco3t(k,i)**0.0916) pdch4t(k,i,lat) = pdch4t(k,i,lat)+(pdch4a(k,i,lat)+ | pdch4b(k,i,lat))*sfeps ! ! New rates for H2O2, CH2O, N2O AND HO2: ! pdh2o2(k,i,lat) = ! XJH2O2 | ((1.0E-4/(1.+1.6951E-17*sco3t(k,i)**0.8573))+ | 1.0E-4*exp(-2.0818E-23*sco2(k,i)**0.9415)* | exp(-8.5266E-14*sco3t(k,i)**0.7466)*0.8721* | exp(-1.4871E-20*sco2 (k,i)**0.8573))*sfeps pdch2oa(k,i,lat) = ! XJCH2OA | (1.2E-4*exp(-2.3481E-40*sco2(k,i)**1.4962)* | exp(-2.1444E-10*sco3t(k,i)**0.5043)* | 0.95*exp(-7.1534E-54*sco2(k,i)**2.0170))*sfeps pdch2ob(k,i,lat) = ! XJCH2OB | (1.1E-4*exp(-6.0858E-70*sco2(k,i)**2.6383)* | exp(-1.8189E-10*sco3t(k,i)**0.4812)* | exp(-2.4759E-7 *sco2 (k,i)**0.1970))*sfeps pdn2o(k,i,lat) = ! XJN2O | ((3.5E-6*exp(-1.1232E-5 *sco2(k,i)**0.2638)+ | 3.5E-7*exp(-3.4971E-18*sco2(k,i)**0.7601))* | exp(-7.5897E-14*sco3t(k,i)**0.7283)* | exp(-1.8121E-33*sco3t(k,i)**1.7512))*sfeps pdho2(k,i,lat) = ! XJHO2 | (8.2E-4*exp(-6.4971E-16*sco2(k,i)**0.6354)* | exp(-3.2108E-12*sco3t(k,i)**0.6469)* | exp(-7.2877E-21*sco3t(k,i)**1.1077))*sfeps ! pdno2(k,i,lat) = 1.3e-2*float(idn(i,lat))*sfeps ! XJNO2 pdch3oo(k,i,lat) = 2.71e-5*float(idn(i,lat))*sfeps ! XJCH3OO pd762(k,i,lat) = 0. ! XJ762 ! enddo ! k=lev0,lev1-1 enddo ! i=lon0,lon1 ! call addfld('QTOTAL4',' ',' ',qtotal4, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('QTOTAL5',' ',' ',qtotal5, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('QTOTAL6',' ',' ',qtotal6, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('QTOTAL7',' ',' ',qtotal7, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! ! call addfld('QTOTAL2',' ',' ',qtotal (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qop2p' ,' ',' ',qop2p (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qop2d' ,' ',' ',qop2d (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qo2p' ,' ',' ',qo2p (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qop' ,' ',' ',qop (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qn2p' ,' ',' ',qn2p (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qnp' ,' ',' ',qnp (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qnop' ,' ',' ',qnop (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qnoplya',' ',' ',qnoplya(:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qnolya',' ',' ',qnolya (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! ! call addfld('pdo2',' ',' ',pdo2 (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdo2d',' ',' ',pdo2d (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdn2',' ',' ',pdn2 (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdco2t',' ',' ',pdco2t (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdh2ol',' ',' ',pdh2ol (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdh2ot',' ',' ',pdh2ot (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdo3d',' ',' ',pdo3d (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdo3p',' ',' ',pdo3p (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdch4t',' ',' ',pdch4t (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdnoeuv',' ',' ',pdnoeuv(:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdnosrb',' ',' ',pdnosrb(:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdh2o2',' ',' ',pdh2o2 (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdch2oa',' ',' ',pdch2oa(:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdch2ob',' ',' ',pdch2ob(:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdn2o',' ',' ',pdn2o (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdho2',' ',' ',pdho2 (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdno2',' ',' ',pdno2 (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pdch3oo',' ',' ',pdch3oo(:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('pd762',' ',' ',pd762 (:,:,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! #ifdef VT ! code = 118 ; state = 'qrj' ; activity='ModelCode' call vtend(118,ier) #endif end subroutine qrj !----------------------------------------------------------------------- subroutine init_sflux(iprint) use input_module,only: f107,f107a use init_module,only: sfeps ! ! Flux initialization once per time step, called from advance. ! ! Args: integer,intent(in) :: iprint ! Local: integer :: n ! sflux: scaled solar flux returned by subroutine ssflux() (photons cm-2 s-1) ! call ssflux(f107,f107a) do n = l1+1,lmax feuv(n) = sflux(n)*sfeps enddo do n = 1,l1 fsrc(n) = sflux(n)*sfeps enddo end subroutine init_sflux !----------------------------------------------------------------------- subroutine init_qrj integer :: m,n ! ! Called once per run, from tgcm.F ! ! Initialize bins (37 bins) ! wave1 = (/1700.00, 1650.00, 1600.00, 1550.00, 1500.00, | 1450.00, 1400.00, 1350.00, 1300.00, 1250.00, | 1200.00, 1215.67, 1150.00, 1100.00, 1050.00, | 1027.00, 987.00, 975.00, 913.00, 913.00, | 913.00, 798.00, 798.00, 798.00, 650.00, | 650.00, 540.00, 320.00, 290.00, 224.00, | 155.00, 70.00, 32.00, 18.00, 8.00, | 4.00, 0.50/) wave2 = (/1750.00, 1700.00, 1650.00, 1600.00, 1550.00, | 1500.00, 1450.00, 1400.00, 1350.00, 1300.00, | 1250.00, 1215.67, 1200.00, 1150.00, 1100.00, | 1050.00, 1027.00, 987.00, 975.00, 975.00, | 975.00, 913.00, 913.00, 913.00, 798.00, | 798.00, 650.00, 540.00, 320.00, 290.00, | 224.00, 155.00, 70.00, 32.00, 18.00, | 8.00, 4.00/) ! ! Solar spectrum based on EUVAC and glow for wave length less than 1050 A ! and Woods for wavelength greater than 1050 A ! ! solar minimum flux (when P_index=80, unit:photon cm^-2 S^-1) ! sfmin=(/3.397e+11, 1.998e+11, 1.055e+11, 7.260e+10, | 5.080e+10, 2.802e+10, 1.824e+10, 1.387e+10, | 2.659e+10, 7.790e+09, 1.509e+10, 3.940e+11, | 8.399e+09, 3.200e+09, 3.298e+09, 4.235e+09, | 4.419e+09, 4.482e+09, 7.156e+08, 1.028e+09, | 3.818e+08, 8.448e+08, 3.655e+09, 2.364e+09, | 1.142e+09, 1.459e+09, 4.830e+09, 2.861e+09, | 8.380e+09, 4.342e+09, 5.612e+09, 1.270e+09, | 5.326e+08, 2.850e+07, 2.000e+06, 1.000e+04, | 5.010e+01/) ! ! scaling factor A as defined in EUVAC model ! afac=(/5.937e-04, 6.089e-04, 1.043e-03, 1.125e-03, | 1.531e-03, 1.202e-03, 1.873e-03, 2.632e-03, | 2.877e-03, 2.610e-03, 3.739e-03, 4.230e-03, | 2.541e-03, 2.099e-03, 3.007e-03, 4.825e-03, | 5.021e-03, 3.950e-03, 4.422e-03, 4.955e-03, | 4.915e-03, 5.437e-03, 5.261e-03, 5.310e-03, | 3.680e-03, 5.719e-03, 5.857e-03, 1.458e-02, | 7.059e-03, 2.575e-02, 1.433e-02, 9.182e-03, | 1.343e-02, 6.247e-02, 2.000e-01, 3.710e-01, | 6.240e-01/) ! ! 0.5*(wave1+wave2) in centimeter ! rlmeuv=(/1.725e-05, 1.675e-05, 1.625e-05, 1.575e-05, | 1.525e-05, 1.475e-05, 1.425e-05, 1.375e-05, | 1.325e-05, 1.275e-05, 1.225e-05, 1.216e-05, | 1.175e-05, 1.125e-05, 1.075e-05, 1.038e-05, | 1.007e-05, 9.810e-06, 9.440e-06, 9.440e-06, | 9.440e-06, 8.555e-06, 8.555e-06, 8.555e-06, | 7.240e-06, 7.240e-06, 5.950e-06, 4.300e-06, | 3.050e-06, 2.570e-06, 1.895e-06, 1.125e-06, | 5.100e-07, 2.500e-07, 1.300e-07, 6.000e-08, | 2.250e-08/) ! ! O2 absorption cross section: ! sigeuv(1,:) = (/ | 5.00e-01, 1.50e+00, 3.40e+00, 6.00e+00, 1.00e+01, | 1.30e+01, 1.50e+01, 1.20e+01, 2.20e+00, 4.00e-01, | 1.30e+01, 1.00e-02, 1.40e+00, 4.00e-01, 1.00e+00, | 1.15e+00, 1.63e+00, 1.87e+01, 3.25e+01, 1.44e+01, | 1.34e+01, 1.33e+01, 1.09e+01, 1.05e+01, 2.49e+01, | 2.36e+01, 2.70e+01, 2.03e+01, 1.68e+01, 1.32e+01, | 7.63e+00, 2.63e+00, 6.46e-01, 2.10e-01, 2.25e-01, | 3.40e-02, 4.54e-03/) ! ! O absorption cross section: ! sigeuv(2,:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 3.79e+00, 4.10e+00, 3.00e+00, 4.79e+00, | 8.52e+00, 1.31e+01, 1.07e+01, 7.72e+00, 6.02e+00, | 3.78e+00, 1.32e+00, 3.25e-01, 1.05e-01, 1.13e-01, | 1.70e-02, 2.27e-03/) ! ! N2 absorption cross section: ! sigeuv(3,:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 2.55e+00, 1.15e+02, 1.44e+01, | 2.18e+00, 7.17e+01, 1.31e+01, 2.14e+00, 5.45e+01, | 2.30e+01, 2.31e+01, 1.97e+01, 1.17e+01, 9.94e+00, | 5.09e+00, 1.53e+00, 3.46e-01, 1.14e+00, 1.41e-01, | 2.01e-02, 2.53e-03/) ! ! O3 absorption cross section ! sigeuv(4,:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 1.25e+01, 9.20e+00, | 9.20e+00, 9.20e+00, 9.16e+00, 9.50e+00, 9.50e+00, | 9.50e+00, 1.47e+01, 1.47e+01, 1.47e+01, 2.74e+01, | 2.02e+01, 3.33e+01, 7.74e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00/) ! ! CO2 absorption cross section ! sigeuv(5,:) = (/ | 5.00e-02, 1.00e-01, 1.50e-01, 3.00e-01, 4.00e-01, | 5.50e-01, 5.00e-01, 5.00e-01, 8.00e-01, 5.00e-01, | 2.00e-01, 0.00e+00, 0.00e+00, 1.85e+01, 1.48e+01, | 1.42e+01, 1.66e+01, 3.98e+01, 7.41e+01, 7.41e+01, | 7.41e+01, 1.74e+01, 1.74e+01, 1.74e+01, 3.23e+01, | 3.31e+01, 3.33e+01, 2.50e+01, 2.34e+01, 2.33e+01, | 1.11e+01, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00/) ! ! H2O absorption cross section ! sigeuv(6,:) = (/ | 5.00e+00, 5.00e+00, 5.00e+00, 3.00e+00, 1.50e+00, | 8.00e-01, 8.00e-01, 1.10e+00, 5.00e+00, 8.00e+00, | 8.00e+00, 0.00e+00, 4.44e+00, 4.44e+00, 4.44e+00, | 1.41e+01, 2.46e+01, 1.10e+01, 1.85e+01, 1.85e+01, | 1.85e+01, 2.36e+01, 2.36e+01, 2.36e+01, 3.98e+01, | 2.33e+01, 2.28e+01, 2.45e+01, 2.22e+01, 1.91e+01, | 1.11e+01, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00/) ! ! NO absorption cross section ! sigeuv(7,:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 1.85e+00, 2.04e+00, | 2.04e+00, 0.00e+00, 2.41e+00, 3.70e+00, 6.48e+00, | 3.70e+00, 7.82e+00, 1.98e+01, 2.41e+01, 2.41e+01, | 2.41e+01, 1.55e+01, 1.55e+01, 1.55e+01, 1.36e+01, | 1.40e+01, 2.00e+01, 2.02e+01, 2.30e+01, 2.40e+01, | 2.22e+01, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00/) ! ! CH4 absorption cross section ! sigeuv(8,:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 2.97e+01, 1.68e+01, 3.74e+01, 4.23e+01, 4.23e+01, | 4.23e+01, 4.90e+01, 4.90e+01, 4.90e+01, 4.21e+01, | 3.89e+01, 3.00e+01, 1.81e+01, 1.03e+01, 1.00e+01, | 1.00e+01, 1.00e+01, 1.00e+01, 1.00e+01, 1.00e+01, | 1.00e+01, 1.00e+01/) ! ! The three major species' ionization branching ratio (off absorption): ! ! O2 ! BPhotonI(1,:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 6.13e-01, 8.30e-01, 6.20e-01, 7.86e-01, | 7.56e-01, 5.34e-01, 5.74e-01, 5.49e-01, 4.76e-01, | 6.73e-01, 9.83e-01, 1.00e+00, 1.00e+00, 1.00e+00, | 1.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, | 1.00e+00, 1.00e+00/) ! ! O ! BPhotonI(2,:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, | 1.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, | 1.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, | 1.00e+00, 1.00e+00/) ! ! N2 ! BPhotonI(3,:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 4.29e-01, | 6.80e-01, 1.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, | 1.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, | 1.00e+00, 1.00e+00/) ! ! photon ionization branching ratio for O+(2p),o+(2d),O+(4s) ! (off O photon ionization) ! ! O+(2p) ! brop2pPh(:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 8.56e-03, 2.52e-01, 2.60e-01, 2.46e-01, 2.41e-01, | 2.33e-01, 2.27e-01, 2.26e-01, 2.24e-01, 2.24e-01, | 2.24e-01, 2.24e-01/) ! ! O+(2d) ! brop2dPh(:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 6.98e-02, | 3.37e-01, 4.51e-01, 4.24e-01, 4.03e-01, 4.02e-01, | 3.92e-01, 3.77e-01, 3.74e-01, 3.78e-01, 3.78e-01, | 3.78e-01, 3.78e-01/) ! ! O+(4s) ! brop4sPh(:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, 9.30e-01, | 6.55e-01, 2.98e-01, 3.17e-01, 3.46e-01, 3.50e-01, | 3.67e-01, 3.89e-01, 3.93e-01, 3.90e-01, 3.90e-01, | 3.90e-01, 3.90e-01/) ! ! O2 photon dissociation braching ratio ! bro2DPh(:) = (/ | 1.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, | 1.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, | 1.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, | 1.00e+00, 3.87e-01, 1.70e-01, 3.80e-01, 2.14e-01, | 2.44e-01, 4.66e-01, 4.26e-01, 4.51e-01, 5.24e-01, | 3.27e-01, 1.74e-02, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00/) ! ! n2 photon dissociation braching ratio for n(2d) ! brn2DPh(:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, | 1.00e+00, 1.00e+00, 1.00e+00, 1.00e+00, 5.71e-01, | 3.20e-01, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00/) ! ! O2 photon dissociative ionization braching ratio ! bro2DIPh(:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 5.35e-04, 1.08e-01, 2.40e-01, 3.51e-01, 3.76e-01, | 4.47e-01, 6.53e-01, 8.92e-01, 1.00e+00, 1.00e+00, | 1.00e+00, 1.00e+00/) ! ! n2 photon dissociative ionization braching ratio ! brn2DIPh(:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 4.72e-03, 9.27e-02, 2.46e-01, | 2.53e-01, 2.49e-01, 2.82e-01, 9.60e-01, 9.60e-01, | 9.60e-01, 9.60e-01/) ! ! n(4s) photoionization cross section ! sigin4s(:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 3.24e+00, 2.48e+00, 2.11e+00, 1.12e+01, | 1.18e+01, 1.16e+01, 9.35e+00, 6.43e+00, 4.80e+00, | 2.55e+00, 6.81e-01, 1.66e-01, 5.68e-01, 7.05e-02, | 1.00e-02, 1.27e-03/) ! ! The three major species' electron impact ionization branching ratio ! off its photon ionization rate ! ! O2 ! BElectronI(1,:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 2.38e-02, 1.05e-01, 2.42e-01, | 5.79e-01, 1.61e+00, 4.27e+00, 6.00e+01, 2.03e+01, | 5.02e+01, 2.11e+02/) ! ! O ! BElectronI(2,:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 1.27e-01, 4.18e-01, 6.94e-01, | 1.09e+00, 2.19e+00, 4.99e+00, 7.14e+01, 2.36e+01, | 5.06e+01, 2.17e+02/) ! ! N2 ! BElectronI(3,:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 3.07e-02, 1.78e-01, 3.61e-01, | 9.33e-01, 2.86e+00, 7.79e+00, 1.08e+01, 3.22e+01, | 8.09e+01, 3.43e+02/) ! ! electron impact ionization branching ratio for O+(2p),o+(2d),O+(4s) ! (off O photon ionization) ! ! O+(2p) ! brop2pEl(:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 9.15e-03, 6.12e-02, 1.16e-01, | 2.03e-01, 4.36e-01, 1.01e+00, 1.46e+01, 4.77e+00, | 1.10e+01, 4.74e+01/) ! ! O+(2d) ! brop2dEl(:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 3.39e-02, 1.48e-01, 2.53e-01, | 4.18e-01, 8.53e-01, 1.96e+00, 2.82e+01, 9.36e+00, | 2.07e+01, 8.85e+01/) ! ! O+(4s) ! brop4sEl(:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 8.36e-02, 2.09e-01, 3.25e-01, | 4.70e-01, 9.02e-01, 2.02e+00, 2.86e+01, 9.42e+00, | 1.89e+01, 8.12e+01/) ! ! photoelectron dissociative ionization branching ratio of O2 ! bro2DIEl(:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 9.60e-04, 1.43e-02, 5.23e-02, | 1.63e-01, 5.21e-01, 1.44e+00, 2.03e+01, 6.98e+00, | 1.79e+01, 7.61e+01/) ! ! photoelectron dissociative ionization branching ratio of N2 ! brn2DIEl(:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 1.84e-04, 8.49e-03, 3.66e-02, | 1.46e-01, 5.71e-01, 1.65e+00, 2.29e+00, 6.95e+00, | 1.83e+01, 7.87e+01/) ! ! photoelectron dissociation branching ratio of N2 ! brn2DEl(:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 1.57e-01, 5.15e-01, 7.64e-01, | 1.37e+00, 2.91e+00, 6.53e+00, 9.05e+00, 2.53e+01, | 5.21e+01, 2.45e+02/) ! ! photoelectron dissociation branching ratio of o2 ! bro2DEl(:) = (/ | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, | 0.00e+00, 1.10e-02, 6.53e-01, 7.62e-01, 9.96e-01, | 1.27e+00, 2.04e+00, 4.11e+00, 5.70e+01, 1.78e+01, | 2.03e+01, 8.79e+01/) ! ! transfer units of cross section to cm^2 do m=1,8 do n=1,lmax sigeuv(m,n)=sigeuv(m,n)*1.E-18 enddo enddo do n = 1,lmax sigin4s(n) = sigin4s(n)*1.E-18 enddo ! do n = 1,l1 rlmsrc(n) = rlmeuv(n) sigsrc(n) = sigeuv(1,n) enddo ! euveff(:) = 0.05 quench = (/7.E-11,5.E-11,3.1401E-12,9.1E-3/) end subroutine init_qrj !--------------------------------------------------------------- subroutine ssflux (f107, f107a) ! ! Args: real,intent(in) :: f107,f107a ! ! Local: real :: pind integer :: l ! ! solar model is the same as EUVAC, i.e., solar flux at ! (f107,f107a) is scale as: sflux=sfmin(1+afac(P-80.)) ! where P=0.5(f107+f107a) pind=0.5*(f107+f107a) do l=1,lmax sflux(l)=sfmin(l)*(1+afac(l)*(pind-80.)) ! set solar flux to be 80% of the value when pind=80 ! if it becomes negative if (sflux(l) .le. 0.1*sfmin(l)) sflux(l) = 0.1*sfmin(l) enddo ! end subroutine ssflux !----------------------------------------------------------------------- subroutine alloc_q(lon0,lon1,lat0,lat1) ! ! Args: integer,intent(in) :: lon0,lon1,lat0,lat1 ! ! Local: integer :: istat ! ! Allocate heating and ionization rates: ! allocate(qtotal(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' qtotal: stat=',i3)") istat allocate(qop2p(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' qop2p: stat=',i3)") istat allocate(qop2d(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' qop2d: stat=',i3)") istat allocate(qo2p(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' qo2p: stat=',i3)") istat allocate(qo2j(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' qo2j: stat=',i3)") istat allocate(qop(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' qop: stat=',i3)") istat allocate(qn2p(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' qn2p: stat=',i3)") istat allocate(qnp(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' qnp: stat=',i3)") istat allocate(qnop(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' qnop: stat=',i3)") istat allocate(qnoplya(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' qnoplya: stat=',i3)") istat allocate(qnolya(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' qnolya: stat=',i3)") istat ! ! Allocate photodissociation rates: ! allocate(pdo2(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdo2: stat=',i3)") istat allocate(pdo2d(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdo2d: stat=',i3)") istat allocate(pdn2(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdn2: stat=',i3)") istat allocate(pdco2t(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdco2t: stat=',i3)") istat allocate(pdco2d(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdco2d: stat=',i3)") istat allocate(pdh2ol(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdh2ol: stat=',i3)") istat allocate(pdh2ot(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdh2ot: stat=',i3)") istat allocate(pdo3d(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdo3d: stat=',i3)") istat allocate(pdo3p(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdo3p: stat=',i3)") istat allocate(pdch4a(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdch4a: stat=',i3)") istat allocate(pdch4b(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdch4b: stat=',i3)") istat allocate(pdch4t(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdch4t: stat=',i3)") istat allocate(pdnoeuv(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdnoeuv: stat=',i3)") istat allocate(pdnosrb(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdnosrb: stat=',i3)") istat ! allocate(attlya(nlevp1,lon0:lon1,lat0:lat1),stat=istat) ! if (istat /= 0) write(6,"('>>> alloc_q: error allocating', ! | ' attlya: stat=',i3)") istat allocate(pdh2o2(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdh2o2: stat=',i3)") istat allocate(pdch2oa(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdch2oa: stat=',i3)") istat allocate(pdch2ob(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdch2ob: stat=',i3)") istat allocate(pdn2o(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdn2o: stat=',i3)") istat allocate(pdho2(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdho2: stat=',i3)") istat allocate(pdno2(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdno2: stat=',i3)") istat allocate(pdch3oo(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pdch3oo: stat=',i3)") istat allocate(pd762(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_q: error allocating', | ' pd762: stat=',i3)") istat end subroutine alloc_q !----------------------------------------------------------------------- subroutine jno(xjno,o2mmr,o1mmr,hemmr,xnmbar,sco2,sco3,scno, | lev0,lev1,lon0,lon1,lat) ! ! Return NO dissociation in xjno. This is called from QRJ, and replaces ! former single-statement calculation of XJNO in QRJ. ! 8/4/04 btf: adapted for timegcm1. ! use cons_module,only: rmassinv_n2 implicit none ! ! Args: integer,intent(in) :: lev0,lev1,lon0,lon1,lat real,intent(out) :: xjno(lev0:lev1,lon0:lon1) real,intent(in),dimension(lev0:lev1,lon0:lon1) :: | o2mmr,o1mmr, ! o2,o densities (mmr) | hemmr, ! he (mmr) | xnmbar, ! mbar (for n2 number density) | sco2,sco3,scno ! slant column integrations from chapman ! ! Local: integer :: k,i,il real :: prob,xnn2,tauo3,sum,xjno1,xjno2,xjno3 ! ! O2(5,0) BAND, NO D(0-0) BAND 190.2-192.5NM ! real,dimension(6) :: so2,w1no,s1no,w2no,s2no data so2 /1.12e-23,2.45e-23,7.19e-23,3.04e-22,1.75e-21,1.11e-20/ data w1no /0. ,5.12e-2 ,1.36e-1 ,1.65e-1 ,1.41e-1 ,4.5e-2 / data s1no /0. ,1.32e-18,6.35e-19,7.09e-19,2.18e-19,4.67e-19/ data w2no /0. ,5.68e-3 ,1.52e-2 ,1.83e-2 ,1.57e-2 ,5.0e-3 / data s2no /0. ,4.41e-17,4.45e-17,4.5e-17 ,2.94e-17,4.35e-17/ ! ! O2(9,0) BAND, NO D(1-0) BAND 183.1-184.6NM ! real,dimension(6) :: s1o2,w11no,s11no,w12no,s12no data s1o2 /1.35e-22,2.99e-22,7.33e-22,3.07e-21,1.69e-20,1.66e-19/ data w11no /0. ,0. ,1.93e-03,9.73e-2 ,9.75e-2 ,3.48e-2 / data s11no /0. ,0. ,3.05e-21,5.76e-19,2.29e-18,2.21e-18/ data w12no /0. ,0. ,2.14e-4 ,1.08e-2 ,1.08e-2 ,3.86e-3 / data s12no /0.,0.,3.2e-21,5.71e-17,9.09e-17,6.0e-17/ ! ! O2(10,0) BAND, NO D(1-0) BAND 181.6-183.1NM ! real,dimension(6) :: s2o2,w21no,s21no,w22no,s22no data s2o2 /2.97e-22,5.83e-22,2.05e-21,8.19e-21,4.8e-20 ,2.66e-19/ data w21no /4.5e-2 ,1.8e-1 ,2.25e-1 ,2.25e-1 ,1.8e-1 ,4.5e-2 / data s21no /1.8e-18 ,1.5e-18 ,5.01e-19,7.20e-20,6.72e-20,1.49e-21/ data w22no /5.e-3 ,2.e-2 ,2.5e-2 ,2.5e-2 ,2.0e-2 ,5.0e-3 / data s22no /1.4e-16 ,1.52e-16,7.e-17 ,2.83e-17,2.73e-17,6.57e-18/ ! ! O3 ABSORPTION CROSS SECTION FOR EACH BAND ! real :: so3a1,so3a2,so3a3 data so3a1,so3a2,so3a3 /42.8e-20, 62.2e-20, 68.8e-20/ real :: dlam,dlam1,dlam2 data dlam,dlam1,dlam2 /2.3, 1.5, 1.5/ ! ! Loop over zp and longitude: ! do i=lon0,lon1 do k=lev0,lev1 xnn2 = (1.-o2mmr(k,i)-o1mmr(k,i)-hemmr(k,i))*xnmbar(k,i)* | rmassinv_n2 prob = 1.65e+9 / (5.1e+7 + 1.65e+9 + 1.5e-9 * xnn2) tauo3 = exp(-so3a1*sco3(k,i)) sum = 0. do il = 1,6 sum = sum + | (exp(-so2(il) *sco2(k,i))*(w1no(il)*s1no(il)* | exp(-s1no(il)*scno(k,i))+w2no(il)*s2no(il)* | exp(-s2no(il)*scno(k,i)))) enddo xjno1 = dlam*3.98e+11*tauo3*sum*prob ! tauo3 = exp(-so3a2*sco3(k,i)) sum = 0. do il = 1,6 sum = sum + | (exp(-s1o2(il) *sco2(k,i))*(w11no(il)*s11no(il)* | exp(-s11no(il)*scno(k,i))+w12no(il)*s12no(il)* | exp(-s12no(il)*scno(k,i)))) enddo xjno2 = dlam1*2.21e+11*tauo3*sum*prob ! tauo3 = exp(-so3a3*sco3(k,i)) sum = 0. do il = 1,6 sum = sum + | (exp(-s2o2(il) *sco2(k,i))*(w21no(il)*s21no(il)* | exp(-s21no(il)*scno(k,i))+w22no(il)*s22no(il)* | exp(-s22no(il)*scno(k,i)))) enddo xjno3 = dlam2*2.30e+11*tauo3*sum*prob ! ! Total NO dissociation: ! xjno(k,i) = xjno1+xjno2+xjno3 enddo enddo end subroutine jno !----------------------------------------------------------------------- subroutine jno_exp(xjno,o2mmr,o1mmr,hemmr,xnmbar,sco2,sco3,scno, | lev0,lev1,lon0,lon1,lat) ! ! Return NO dissociation in xjno. This is called from QRJ, and replaces ! former single-statement calculation of XJNO in QRJ. ! 8/4/04 btf: adapted for timegcm1. ! use cons_module,only: rmassinv_n2 implicit none ! ! Args: integer,intent(in) :: lev0,lev1,lon0,lon1,lat real,intent(out) :: xjno(lev0:lev1,lon0:lon1) real,intent(in),dimension(lev0:lev1,lon0:lon1) :: | o2mmr,o1mmr, ! o2,o densities (mmr) | hemmr, ! he (mmr) | xnmbar, ! mbar (for n2 number density) | sco2,sco3,scno ! slant column integrations from chapman ! ! Local: integer :: k,i,il real :: prob,xnn2,tauo3,sum,xjno1,xjno2,xjno3 ! ! O2(5,0) BAND, NO D(0-0) BAND 190.2-192.5NM ! real,dimension(6) :: so2,w1no,s1no,w2no,s2no data so2 /1.12e-23,2.45e-23,7.19e-23,3.04e-22,1.75e-21,1.11e-20/ data w1no /0. ,5.12e-2 ,1.36e-1 ,1.65e-1 ,1.41e-1 ,4.5e-2 / data s1no /0. ,1.32e-18,6.35e-19,7.09e-19,2.18e-19,4.67e-19/ data w2no /0. ,5.68e-3 ,1.52e-2 ,1.83e-2 ,1.57e-2 ,5.0e-3 / data s2no /0. ,4.41e-17,4.45e-17,4.5e-17 ,2.94e-17,4.35e-17/ ! ! O2(9,0) BAND, NO D(1-0) BAND 183.1-184.6NM ! real,dimension(6) :: s1o2,w11no,s11no,w12no,s12no data s1o2 /1.35e-22,2.99e-22,7.33e-22,3.07e-21,1.69e-20,1.66e-19/ data w11no /0. ,0. ,1.93e-03,9.73e-2 ,9.75e-2 ,3.48e-2 / data s11no /0. ,0. ,3.05e-21,5.76e-19,2.29e-18,2.21e-18/ data w12no /0. ,0. ,2.14e-4 ,1.08e-2 ,1.08e-2 ,3.86e-3 / data s12no /0.,0.,3.2e-21,5.71e-17,9.09e-17,6.0e-17/ ! ! O2(10,0) BAND, NO D(1-0) BAND 181.6-183.1NM ! real,dimension(6) :: s2o2,w21no,s21no,w22no,s22no data s2o2 /2.97e-22,5.83e-22,2.05e-21,8.19e-21,4.8e-20 ,2.66e-19/ data w21no /4.5e-2 ,1.8e-1 ,2.25e-1 ,2.25e-1 ,1.8e-1 ,4.5e-2 / data s21no /1.8e-18 ,1.5e-18 ,5.01e-19,7.20e-20,6.72e-20,1.49e-21/ data w22no /5.e-3 ,2.e-2 ,2.5e-2 ,2.5e-2 ,2.0e-2 ,5.0e-3 / data s22no /1.4e-16 ,1.52e-16,7.e-17 ,2.83e-17,2.73e-17,6.57e-18/ ! ! O3 ABSORPTION CROSS SECTION FOR EACH BAND ! real :: so3a1,so3a2,so3a3 data so3a1,so3a2,so3a3 /42.8e-20, 62.2e-20, 68.8e-20/ real :: dlam,dlam1,dlam2 data dlam,dlam1,dlam2 /2.3, 1.5, 1.5/ ! ! Loop over zp and longitude: ! do i=lon0,lon1 do k=lev0,lev1 xnn2 = (1.-o2mmr(k,i)-o1mmr(k,i)-hemmr(k,i))*xnmbar(k,i)* | rmassinv_n2 prob = 1.65e+9 / (5.1e+7 + 1.65e+9 + 1.5e-9 * xnn2) tauo3 = expo(-so3a1*sco3(k,i)) sum = 0. do il = 1,6 sum = sum + | (expo(-so2(il) *sco2(k,i))*(w1no(il)*s1no(il)* | expo(-s1no(il)*scno(k,i))+w2no(il)*s2no(il)* | expo(-s2no(il)*scno(k,i)))) enddo xjno1 = dlam*3.98e+11*tauo3*sum*prob ! tauo3 = expo(-so3a2*sco3(k,i)) sum = 0. do il = 1,6 sum = sum + | (expo(-s1o2(il) *sco2(k,i))*(w11no(il)*s11no(il)* | expo(-s11no(il)*scno(k,i))+w12no(il)*s12no(il)* | expo(-s12no(il)*scno(k,i)))) enddo xjno2 = dlam1*2.21e+11*tauo3*sum*prob ! tauo3 = expo(-so3a3*sco3(k,i)) sum = 0. do il = 1,6 sum = sum + | (expo(-s2o2(il) *sco2(k,i))*(w21no(il)*s21no(il)* | expo(-s21no(il)*scno(k,i))+w22no(il)*s22no(il)* | expo(-s22no(il)*scno(k,i)))) enddo xjno3 = dlam2*2.30e+11*tauo3*sum*prob ! ! Total NO dissociation: ! xjno(k,i) = xjno1+xjno2+xjno3 enddo enddo end subroutine jno_exp !----------------------------------------------------------------------- subroutine jo2h2osrb(xjo2srb,xjh2osrb,sco2,sco3, | lev0,lev1,lon0,lon1,lat) ! ! Dissociation of O2 and H2O from Shumann-Runge bands. ! implicit none ! ! Args: integer,intent(in) :: lev0,lev1,lon0,lon1,lat real,intent(out),dimension(lev0:lev1,lon0:lon1) :: | xjo2srb,xjh2osrb ! ! 11/29/04 btf: fix dimension bug in sco2, sco3: real,intent(in),dimension(lev0:lev1,lon0-2:lon1+2) :: | sco2,sco3 ! ! Local: integer,parameter :: | nwt=6, ! number of weights | nwl=16 ! number of wavelength bins real :: sumo2,sumh2o,trano3(nwl),exptmp integer :: i,k,il,jl real :: w(nwt) = (/0.05, 0.20, 0.25, 0.25, 0.20, 0.05/) real :: sg(nwt,nwl) = reshape(source=(/ | 3.770E-21, 1.280E-20, 3.600E-20, 7.690E-20, 2.640E-19, 1.050E-18, | 1.033E-21, 1.750E-21, 4.588E-21, 1.714E-20, 1.012E-19, 1.097E-18, | 7.680E-22, 1.244E-21, 2.670E-21, 1.080E-20, 8.363E-20, 7.697E-19, | 1.127E-21, 2.015E-21, 4.660E-21, 1.650E-20, 9.295E-20, 5.017E-19, | 5.560E-22, 1.579E-21, 3.715E-21, 1.377E-20, 7.221E-20, 3.439E-19, | 2.968E-22, 5.831E-22, 2.053E-21, 8.192E-21, 4.802E-20, 2.655E-19, | 1.350E-22, 2.991E-22, 7.334E-22, 3.074E-21, 1.689E-20, 1.658E-19, | 1.079E-22, 2.092E-22, 5.875E-22, 2.590E-21, 1.578E-20, 1.028E-19, | 4.494E-23, 6.684E-23, 2.384E-22, 1.126E-21, 6.988E-21, 5.549E-20, | 1.910E-23, 3.438E-23, 1.168E-22, 4.736E-22, 2.653E-21, 2.529E-20, | 1.117E-23, 2.447E-23, 7.188E-23, 3.042E-22, 1.748E-21, 1.112E-20, | 9.600E-24, 1.264E-23, 2.553E-23, 1.049E-22, 5.190E-22, 2.824E-21, | 6.815E-24, 7.499E-24, 1.030E-23, 2.478E-23, 1.517E-22, 1.254E-21, | 6.977E-24, 7.493E-24, 9.308E-24, 1.653E-23, 7.869E-23, 4.630E-22, | 6.739E-24, 6.774E-24, 6.931E-24, 7.357E-24, 1.039E-23, 5.183E-23, | 6.842E-24, 6.850E-24, 6.878E-24, 7.113E-24, 8.407E-24, 2.869E-23 | /),shape=(/nwt,nwl/)) real :: si(nwl) = (/ 0.825E+10, | 1.334E+10, 1.535E+10, 1.531E+10, 1.835E+10, 2.301E+10, | 2.208E+10, 2.200E+10, 2.879E+10, 3.370E+10, 3.980E+10, | 4.176E+10, 5.823E+10, 6.270E+10, 6.554E+10, 7.928E+10/) real :: dl(nwl) = (/ 20., | 8., 10., 11., 12., 15., 15., 17., 19., 20., 23., 22., 25., 13., | 15.,25./) real :: sigo3(nwl) = (/ 81.1E-20, | 79.9E-20,78.6E-20,78.6E-20,76.3E-20,76.3E-20,68.8E-20,62.2E-20, | 57.6E-20,52.6E-20,47.6E-20,42.8E-20,38.3E-20,34.7E-20,32.3E-20, | 31.4E-20/) real :: sgh2o(nwt,nwl) = reshape(source=(/ | 1.870E-18, 2.050E-18, 2.180E-18, 2.270E-18, 2.190E-18, 2.190E-18, | 1.600E-18, 1.610E-18, 1.670E-18, 1.690E-18, 1.720E-18, 1.720E-18, | 1.250E-18, 1.270E-18, 1.330E-18, 1.370E-18, 1.340E-18, 1.397E-18, | 8.450E-19, 8.520E-19, 9.150E-19, 9.780E-19, 9.690E-19, 1.050E-18, | 4.300E-19, 4.510E-19, 5.010E-19, 5.730E-19, 6.110E-19, 6.640E-19, | 1.890E-19, 2.050E-19, 2.550E-19, 2.900E-19, 3.020E-19, 3.400E-19, | 8.160E-20, 8.980E-20, 9.940E-20, 1.270E-19, 1.370E-19, 1.500E-19, | 4.720E-20, 4.980E-20, 5.360E-20, 6.070E-20, 6.430E-20, 6.760E-20, | 2.960E-20, 2.950E-20, 3.340E-20, 3.710E-20, 3.970E-20, 4.110E-20, | 1.190E-20, 1.270E-20, 1.430E-20, 1.890E-20, 2.180E-20, 2.300E-20, | 3.490E-21, 3.520E-21, 4.240E-21, 5.840E-21, 6.610E-21, 8.170E-21, | 1.600E-21, 1.810E-21, 2.040E-21, 2.380E-21, 2.780E-21, 2.940E-21, | 3.840E-22, 4.440E-22, 5.790E-22, 7.330E-22, 9.430E-22, 1.100E-21, | 1.840E-22, 1.920E-22, 1.990E-22, 2.280E-22, 2.430E-22, 2.620E-22, | 7.300E-23, 6.700E-23, 7.360E-23, 9.580E-23, 7.270E-23, 2.320E-23, | 0., 0., 0., 0., 0., 0. | /),shape=(/nwt,nwl/)) ! ! o2 and h2o photodissociation ! do i=lon0,lon1 do k=lev0,lev1 xjo2srb (k,i) = 0. xjh2osrb(k,i) = 0. do il=1,16 trano3(il) = exp(-sigo3(il)*sco3(k,i)) sumo2 = 0. sumh2o = 0. do jl=1,6 exptmp = exp(-sg(jl,il)*sco2(k,i)) sumo2 = sumo2 +sg (jl,il)*w(jl)*exptmp sumh2o = sumh2o+sgh2o(jl,il)*w(jl)*exptmp enddo ! jl=1,6 xjo2srb(k,i) = xjo2srb(k,i) +si(il)*dl(il)*sumo2 *trano3(il) xjh2osrb(k,i)= xjh2osrb(k,i)+si(il)*dl(il)*sumh2o*trano3(il) enddo ! il=1,16 enddo ! k=lev0,lev1 enddo ! i=lon0,lon1 end subroutine jo2h2osrb !----------------------------------------------------------------------- end module qrj_module