! module ox_module ! ! Calculate matrix coefficients in FS for partitioning of OX into O and O3. ! This must be called before the composition module (comp_major.F). ! implicit none contains !----------------------------------------------------------------------- subroutine comp_ox(tn,o2,ox,op,h2o,h2o2,oh,ho2,h2,h,o1,o1d,o3,no, | no2,n4s,n2d,o2p,n2p,nplus,hplus,nop,ne,co,ch4,cl,clo,xiop2p, | xiop2d,barm,xnmbar,lev0,lev1,lon0,lon1,lat) ! ! This routine calculates partitioning coefficients for OX in ! the array fs(lon0:lon1,nlevp1,2,0:2,lat0:lat1). The fs array ! is allocated in the chemrates module (chemrates.F). ! use cons_module,only: | rmassinv_o1 ,rmassinv_oh ,rmassinv_ho2,rmassinv_h , | rmassinv_no ,rmassinv_no2,rmassinv_n4s,rmassinv_n2d, | rmassinv_h2o2,rmassinv_ch4,rmassinv_h2 ,rmassinv_co , | rmassinv_o3 ,rmassinv_o2 ,rmassinv_n2 , | rmass_o1 ,rmass_o2 ,rmass_o3 use qrj_module,only: | pdo3d, ! photodissociation of o3(d) (F(NJO3D)) | pdo3p, ! photodissociation of o3(p) (F(NJO3P)) | pdo2, ! photodissociation of o2 (F(NRJ)) | qo2p, ! o2+ ionization (F(NQO2P)) | qo2j, ! (F(NQO2J)) is zero | qop, ! o+ ionization (F(NQOP)) | pdh2ol, ! lyman-alpha photodissociation of h2o (F(NJH2OL)) | pdnosrb,! photodissociation of no by srb (XJNO) | pdno2 ! XJNO2 use comp_meta_module,only: ! (nlevp1,lon0:lon1,lat0:lat1) | oh_h, ! oh/h (ratio1) | ho2_h, ! ho2/h (ratio2) | h_hox, ! h/hox (ratio3) | rmtru_hox, ! "true" mass of hox | rmtru_ox, ! "true" mass of ox (defined by this routine) | ch3,ch2o,ch3o2,ch3o,cho use chemrates_module,only: beta1,beta2,beta3,beta6,beta9,beta11, | beta12,beta13,rkm7a,rkm7b,rkm20,rkm21,rkm22,rkm23,rkm24,rkm25, | rkm26,rkm27,rkm28,rkm29,rkm31,rkm34,rkm35,rkm36,rkm37,rkm38, | rkm30,rkm40, rk1,rk3,rk4,rk5,rk6,rk7,rk8,rk9,rk10,rk11,rk12, | rk16,rk23,rk24,rk27,gam3,gam4,gam6,gam7,gam9,gam11,gam12,gam13, | gam14,gam15,del1,del2,ra1,ra2, | fs ! (fs(lon0:lon1,nlevp1,2,0:2,lat0:lat1) (defined by this routine) use init_module,only: istep ! for debug only implicit none ! ! Input args: integer,intent(in) :: lev0,lev1,lon0,lon1,lat real,dimension(lev0:lev1,lon0-2:lon1+2),intent(in) :: | tn, ! neutral temperature | o2, ! molecular oxygen | ox, ! ox from previous step | op, ! o+ ion | o1, ! atomic oxygen | o1d, ! o(1d) | o3, ! ozone | no, ! no | n4s, ! n(4s) | n2d, ! n(2d) | o2p, ! o2+ from elden | n2p, ! n2+ from elden | nplus, ! n+ from oplus | hplus, ! h+ from oplus (was xnhp) | nop, ! no+ from elden | ne, ! ne from elden | xiop2p, ! O+(2p) from oplus module | xiop2d, ! O+(2d) from oplus module | barm, ! mbar (interfaces) | xnmbar ! p0e(-z), etc. ! ! Local: integer :: i,k,i0,i1,nk,nkm1,ibndbot,ibndtop real,dimension(lev0:lev1,lon0:lon1) :: | o3_o, ! o3/o ratio s7 | o_ox, ! o/ox ratio s8 | barmi, ! barm at midpoints | n2, ! n2 (mmr) n2=1-o2-o ! Photodissociation: | pdo3i, ! pdo3d+pdo3p at midpoints s13 (NJO3D+NJO3P) | pdo2i, ! pdo2 at midpoints s12 (NRJ) ! Ionization: | qo2pi, ! qo2p at midpoints s11 (NQO2P) | qo2ji, ! qo2j at midpoints s10 (NQO2J) | qopi, ! qop at midpoints s9 (NQOP) ! Production and loss rates: | ox_prod1, ! ox production s1 | ox_prod2, ! ox production s2 | ox_loss1, ! ox loss s3 | ox_loss2, ! ox loss s4 | o2_prod1, ! o2 production s6 | o2_prod2, ! o2 production s9 | o2_prod3, ! o2 production s10 | o2_loss ! o2 loss s13 ! ! Diagnostics for secondary history (see addfsech calls): real,dimension(lev0:lev1,lon0:lon1) :: | fs11, ! fs(i,k,1,1,lat) | fs12, ! fs(i,k,1,2,lat) | fs21, ! fs(i,k,2,1,lat) | fs22, ! fs(i,k,2,2,lat) | fs10, ! fs(i,k,1,0,lat) | fs20 ! fs(i,k,2,0,lat) ! i0=lon0 ; i1=lon1 ; nk=lev1-lev0+1 ; nkm1=nk-1 ! ! Convert several fields from full-levels (interfaces) to half-levels ! (midpoints). Most of these are from the qrj module. ! do i=lon0,lon1 do k=lev0,lev1-1 pdo3i(k,i) = 0.5*(pdo3d (k,i,lat)+pdo3d(k+1,i,lat))+ ! s13 | 0.5*(pdo3p (k,i,lat)+pdo3p(k+1,i,lat)) pdo2i(k,i) = 0.5*(pdo2 (k,i,lat)+pdo2 (k+1,i,lat)) ! s12 qo2pi(k,i) = 0.5*(qo2p (k,i,lat)+qo2p (k+1,i,lat)) ! s11 qo2ji(k,i) = 0.5*(qo2j (k,i,lat)+qo2j (k+1,i,lat)) ! s10 qopi (k,i) = 0.5*(qop (k,i,lat)+qop (k+1,i,lat)) ! s9 barmi(k,i) = 0.5*(barm(k,i)+barm(k+1,i)) enddo ! k=lev0,lev1-1 enddo ! i=lon0,lon1 ! call addfsech('PDO3I' ,' ',' ',pdo3i ,i0,i1,nk,nkm1,lat) ! call addfsech('PDO2I' ,' ',' ',pdo2i ,i0,i1,nk,nkm1,lat) ! call addfsech('QO2PI' ,' ',' ',qo2pi ,i0,i1,nk,nkm1,lat) ! call addfsech('QO2JI' ,' ',' ',qo2ji ,i0,i1,nk,nkm1,lat) ! call addfsech('QOPI' ,' ',' ',qopi ,i0,i1,nk,nkm1,lat) do i=lon0,lon1 do k=lev0,lev1-1 ! ! Ratio O3/O: ! o3_o(k,i) = xnmbar(k,i)*rkm21(k,i,lat)*o2(k,i)*rmassinv_o2* ! s7 | xnmbar(k,i)/barmi(k,i) / (xnmbar(k,i)* | (rkm24(k,i,lat)*o1 (k,i)*rmassinv_o1+ | (rkm7a+rkm7b) *o1d(k,i)*rmassinv_o1+ | rkm29(k,i,lat) *oh (k,i)*rmassinv_oh+ | rkm34(k,i,lat) *ho2(k,i)*rmassinv_ho2+ | rkm37(k,i,lat) *h (k,i)*rmassinv_h+ | beta9(k,i,lat) *no (k,i)*rmassinv_no+ | beta12(k,i,lat)*no2(k,i)*rmassinv_no2)+ | pdo3i(k,i)+del1(k,i,lat)*cl(k,i)*xnmbar(k,i)/barmi(k,i)) ! ! Ratio O/OX: o_ox(k,i) = 1./(1.+o3_o(k,i)) ! s8 ! ! OX production terms: ! n2(k,i) = 1.-o2(k,i)-ox(k,i) ox_prod1(k,i) = xnmbar(k,i)**2* ! s1 | (beta3(k,i,lat)*n4s(k,i)*rmassinv_n4s*no (k,i)*rmassinv_no + | beta6 *n2d(k,i)*rmassinv_n2d*no (k,i)*rmassinv_no) | +0.5*(pdnosrb(k,i,lat)+pdnosrb(k+1,i,lat))* | no(k,i)*rmassinv_no*xnmbar(k,i)+ | xnmbar(k,i)*(rk4*o2p(k,i)*n4s(k,i)*rmassinv_n4s | +rk10*op(k,i)*n2d(k,i)*rmassinv_n2d+ | (rk16*xiop2p(k,i)+rk24*xiop2d(k,i))*n2(k,i)*rmassinv_n2)+ | (ra1(k,i,lat)*nop(k,i)+2.*ra2(k,i,lat)*o2p(k,i))* | sqrt(ne(k,i)*ne(k+1,i)) ! ox_prod2(k,i) = xnmbar(k,i)* ! s2 | (beta1(k,i,lat)*n4s(k,i)*rmassinv_n4s+ | beta2*n2d(k,i)*rmassinv_n2d)+ | rk1(k,i,lat)*op(k,i)+rk7*nplus(k,i)+ | rk23*xiop2p(k,i)+rk27*xiop2d(k,i)+ | 2.*pdo2i(k,i) ! ! OX loss terms: ! ox_loss1(k,i) = ! s3 | (2.*rkm20(k,i,lat)*xnmbar(k,i)/barmi(k,i)+ | 2.*rkm24(k,i,lat)*o3_o(k,i))*o_ox(k,i)**2 ! ox_loss2(k,i) = (xnmbar(k,i)* ! s4 | (beta9 (k,i,lat)*no (k,i)*rmassinv_no *o3_o(k,i)+ | rk3(k,i,lat)*n2p(k,i) +rk8 *nplus(k,i)+ | xnmbar(k,i)/barmi(k,i) ! ! O2 production: ! o2_prod1(k,i) = ! s6 | (rkm20(k,i,lat)*xnmbar(k,i)/barmi(k,i)+ | 2.*rkm24(k,i,lat)*o3_o(k,i))*o_ox(k,i)**2 ! o2_prod2(k,i) = ! s9 | ((beta9(k,i,lat)*no (k,i)*rmassinv_no *o3_o(k,i)+ | xnmbar(k,i)+pdo3i(k,i)*o3_o(k,i))*o_ox(k,i)+ ! o2_prod3(k,i) = xnmbar(k,i)**2* ! s10 | (rkm31(k,i,lat)*oh(k,i)*rmassinv_oh*ho2(k,i)*rmassinv_ho2+ | rkm38 *h (k,i)*rmassinv_h *ho2(k,i)*rmassinv_ho2+ | rkm35(k,i,lat)*(ho2(k,i)*rmassinv_ho2)**2)+ | rk5*no(k,i)*rmassinv_no*o2p(k,i)*xnmbar(k,i)+ | gam6(k,i,lat)*ch3o2(k,i,lat)*ho2(k,i)*rmassinv_ho2* | xnmbar(k,i)+gam7(k,i,lat)*ch3o2(k,i,lat)**2+ | gam15*co(k,i)*rmassinv_co*o1d(k,i)*rmassinv_o1* | xnmbar(k,i)**2 ! ! O2 loss: ! o2_loss(k,i) = xnmbar(k,i)* ! s13 | (rkm36(k,i,lat)*xnmbar(k,i)/barmi(k,i)*h(k,i)*rmassinv_h+ | beta1(k,i,lat)*n4s(k,i)*rmassinv_n4s+ | beta2 *n2d(k,i)*rmassinv_n2d)+ | rk1(k,i,lat)*op(k,i)+(rk6+rk7)*nplus(k,i)+ | rk9*n2p(k,i)+pdo2i(k,i)+rk23*xiop2p(k,i)+rk27*xiop2d(k,i)+ | gam3(k,i,lat)*ch3(k,i,lat)*xnmbar(k,i)/barmi(k,i)+ | gam9(k,i,lat)*ch3o(k,i,lat)+gam12(k,i,lat)*cho(k,i,lat) ! ! "True" OX mass: ! rmtru_ox(k,i,lat) = o_ox(k,i)*(rmass_o1+o3_o(k,i)*rmass_o3) ! s7 ! ! Matrix coefficients (fs(lon0:lon1,nlevp1,2,0:2,lat0:lat1)) ! fs(i,k,1,1,lat) = -o2_loss(k,i)-xnmbar(k,i)**2*o_ox(k,i)* | ox(k,i)/rmtru_ox(k,i,lat)*(1./3.*rkm21(k,i,lat)*o_ox(k,i)* | ox(k,i)/rmtru_ox(k,i,lat)+ 2./3.*rkm22(k,i,lat)* | o2(k,i)*rmassinv_o2 + 1./2.*rkm23(k,i,lat)* | n2(k,i)*rmassinv_n2) fs(i,k,1,2,lat) = (o2_prod2(k,i)+xnmbar(k,i)*o2_prod1(k,i)* | ox(k,i)/rmtru_ox(k,i,lat)-xnmbar(k,i)**2*o_ox(k,i)* | o2(k,i)*rmassinv_o2*(2./3.*rkm21(k,i,lat)*o_ox(k,i)* | ox(k,i)/rmtru_ox(k,i,lat)+1./3.*rkm22(k,i,lat)*o2(k,i)* | rmassinv_o2+1./2.*rkm23(k,i,lat)*n2(k,i)*rmassinv_n2))* | rmass_o2/rmtru_ox(k,i,lat) fs(i,k,2,1,lat) = ox_prod2(k,i)*rmtru_ox(k,i,lat)*rmassinv_o2 fs(i,k,2,2,lat) = -ox_loss2(k,i)-ox_loss1(k,i)*ox(k,i)/ | rmtru_ox(k,i,lat)*xnmbar(k,i) fs(i,k,1,0,lat) = (o2_prod3(k,i)-qo2pi(k,i))*rmass_o2/ | xnmbar(k,i) fs(i,k,2,0,lat) = (ox_prod1(k,i)-qopi(k,i))*rmtru_ox(k,i,lat)/ | xnmbar(k,i) ! ! Save fs to secondary history: fs11(k,i) = fs(i,k,1,1,lat) fs12(k,i) = fs(i,k,1,2,lat) fs21(k,i) = fs(i,k,2,1,lat) fs22(k,i) = fs(i,k,2,2,lat) fs10(k,i) = fs(i,k,1,0,lat) fs20(k,i) = fs(i,k,2,0,lat) enddo ! k=lev0,lev1-1 enddo ! i=lon0,lon1 ! 5/21/04 btf: NaNQ at step 5 (80,4,20), lats 8-12 (-52.5->-32.5) ! Both xiop2p and xiop2d have NaNQ at lat=10, k=31-44, ! i=30-32, among 1.e-20 values. ! These come out of oplus.F, so am retreating to there.. ! if (istep==5) then ! do k=lev0,lev1-1 ! write(6,"(/,'comp_ox: istep=',i2,' lat=',i3,' k=',i3, ! | ' xiop2p(k,lon0:lon1)=',/,(6e12.4))") istep,lat,k, ! | xiop2p(k,lon0:lon1) ! write(6,"(/,'xiop2d(k,lon0:lon1)=',/,(6e12.4))") ! | xiop2d(k,lon0:lon1) ! enddo ! endif ! 5/21/04 btf: NaNQ at step 5 (80,4,20), lats 8-12 (-52.5->-32.5) ! if (istep==5) then ! do k=lev0,lev1-1 ! write(6,"(/,'comp_ox: istep=',i2,' lat=',i3,' k=',i3, ! | ' ox_prod1(k,:)=',/,(6e12.4))") istep,lat,k,ox_prod1(k,:) ! NaN ! write(6,"(/,'ox_prod2(k,:)=',/,(6e12.4))") ox_prod2(k,:) ! NaN ! write(6,"(/,'ox_loss1(k,:)=',/,(6e12.4))") ox_loss1(k,:) ! ok ! write(6,"(/,'ox_loss2(k,:)=',/,(6e12.4))") ox_loss2(k,:) ! ok ! write(6,"(/,'o2_prod1(k,:)=',/,(6e12.4))") o2_prod1(k,:) ! ok ! write(6,"(/,'o2_prod2(k,:)=',/,(6e12.4))") o2_prod2(k,:) ! ok ! write(6,"(/,'o2_prod3(k,:)=',/,(6e12.4))") o2_prod3(k,:) ! ok ! write(6,"(/,'o2_loss (k,:)=',/,(6e12.4))") o2_loss (k,:) ! NaNQ ! enddo ! endif ! 5/21/04 btf: NaNQ at step 5 (80,4,20), lats 8-12 (-52.5->-32.5) ! fs(:,k,iox,0,lat) was found to be NaNs in comp_major. ! Here, fs(:,k,1,1,lat) are NaNs ! fs(:,k,2,1,lat) are NaNs ! fs(:,k,2,0,lat) are NaNs ! if (istep==5) then ! do k=lev0,lev1-1 ! write(6,"(/,'comp_ox: istep=',i2,' lat=',i3,' k=',i3, ! | ' fs(lon0:lon1,k,1,1,lat)=',/,(6e12.4))") istep,lat,k, ! NaN ! | fs(lon0:lon1,k,1,1,lat) ! write(6,"(/,'fs(lon0:lon1,k,1,2,lat)=',/,(6e12.4))") ! ok ! | fs(lon0:lon1,k,1,2,lat) ! write(6,"(/,'fs(lon0:lon1,k,2,1,lat)=',/,(6e12.4))") ! NaN ! | fs(lon0:lon1,k,2,1,lat) ! write(6,"(/,'fs(lon0:lon1,k,2,2,lat)=',/,(6e12.4))") ! ok ! | fs(lon0:lon1,k,2,2,lat) ! write(6,"(/,'fs(lon0:lon1,k,1,0,lat)=',/,(6e12.4))") ! ok ! | fs(lon0:lon1,k,1,0,lat) ! write(6,"(/,'fs(lon0:lon1,k,2,0,lat)=',/,(6e12.4))") ! NaN ! | fs(lon0:lon1,k,2,0,lat) ! enddo ! endif ! ! 12/03: Diffs vs tgcm24 in ox_prod1, ox_prod2 and o2_loss near zp -5 ! around -60 and +60 deg latitude. ! ! call addfsech('O3_O' ,' ',' ',o3_o ,i0,i1,nk,nkm1,lat) ! call addfsech('OX_PROD1',' ',' ',ox_prod1,i0,i1,nk,nkm1,lat) ! call addfsech('OX_PROD2',' ',' ',ox_prod2,i0,i1,nk,nkm1,lat) ! call addfsech('OX_LOSS1',' ',' ',ox_loss1,i0,i1,nk,nkm1,lat) ! call addfsech('OX_LOSS2',' ',' ',ox_loss2,i0,i1,nk,nkm1,lat) ! call addfsech('O2_PROD1',' ',' ',o2_prod1,i0,i1,nk,nkm1,lat) ! call addfsech('O2_PROD2',' ',' ',o2_prod2,i0,i1,nk,nkm1,lat) ! call addfsech('O2_PROD3',' ',' ',o2_prod3,i0,i1,nk,nkm1,lat) ! call addfsech('O2_LOSS' ,' ',' ',o2_loss ,i0,i1,nk,nkm1,lat) ! call addfsech('RMTRU_OX',' ',' ',rmtru_ox(:,lon0:lon1,lat), ! | i0,i1,nk,nkm1,lat) ! ! 12/03: Small diffs vs tgcm24 in fs11, fs21, fs10, and fs20, probably ! related to the diffs noted above in ox production and o2 loss. ! ! call addfsech('FS11' ,' ',' ',fs11,i0,i1,nk,nkm1,lat) ! call addfsech('FS12' ,' ',' ',fs12,i0,i1,nk,nkm1,lat) ! call addfsech('FS21' ,' ',' ',fs21,i0,i1,nk,nkm1,lat) ! call addfsech('FS22' ,' ',' ',fs22,i0,i1,nk,nkm1,lat) ! call addfsech('FS10' ,' ',' ',fs10,i0,i1,nk,nkm1,lat) ! call addfsech('FS20' ,' ',' ',fs20,i0,i1,nk,nkm1,lat) ! end subroutine comp_ox !----------------------------------------------------------------------- end module ox_module