#include "dims.h" SUBROUTINE CMPMETA use cons_module,only: len1,len2,kmax,kmaxp1,brn2d, | rmass_o2, rmass_o, rmassinv_o2, rmassinv_o, rmassinv_n2 implicit none C **** C **** THIS SUBROUTINE CALCULATES PSI(N(2D)),PSI(O(1D)), C **** PSI(O2(1DELTA)),PSI(O2(1SIGMA)) C **** ALL ASSUMING PHOTOCHEMICAL EQUILIBRIUM C **** #include "params.h" #include "fcom.h" #include "vscr.h" #include "index.h" #include "buff.h" #include "crates_const.h" #include "crates_tdep.h" #include "compcom.h" #include "mwt.h" #include "phys.h" ! ! Local: integer :: ntk,npohk,npho2k,nphk,nphoxk,npcok,npn4sk,npnok,npn2dk, | nps1k,nps2k,npo1k,nek,nopk,nn2pk,nnopk,nqtefk,no2pk,npo3k, | nph2ok,nph2k,npch4k,nmsk,npo21dk,njo3dk,njo2dk,npco2k,i,k real :: fmin,fmax real :: RO1DN2D(ZIMXP,ZKMXP),E6300(ZIMXP,ZKMXP),E5200(ZIMXP,ZKMXP) real,parameter :: a1=.079, a2=.083, rk2a=4.e-17, rk2b=2.2e-15, + rk2c=8.e-14, dp=6.6, dpp=19. real :: xo2(zimxp,zkmxp),xo1(zimxp,zkmxp),xn2(zimxp,zkmxp),rkk1, | eo200a(zimxp,zkmxp),xno21s_cm3(zimxp,zkmxp) ! NPNOK = NJ+NPNO NPN2DK = NJNP+NPN2D NPS1K = NJ+NPS NPO1K = NJ+NPO1 NEK = NJ+NE NOPK = NJ+NOP NN2PK = NN2P NNOPK = NNOP NQTEFK = NQTEF-1 DO 1 K=1,KMAX NQTEFK = NQTEFK+1 DO 1 I=1,LEN1 C **** C **** S1 = EXPS(K+1/2), S2 = PN2D(K+1/2) C **** S2(I,K) = .5*(F(I,NQTEFK)+F(I,NQTEFK+1))*BRN2D 1 CONTINUE ! DO 3 I=1,LEN2 C **** S1 = N*MBAR (K+1/2) S1(I,1) = XNMBAR(I,1) C **** C **** S2 = TOTAL N2D PRODUCTION C **** S2(I,1)=S2(I,1)+RK3(I,1)*F(I,NN2PK)*S1(I,1)*F(I,NPO1K)* | rmassinv_o+(RA1(I,1)*F(I,NNOPK)*0.85+RA3(I,1)*F(I,NN2PK)*0.9)* | SQRT(F(I,NEK)*F(I,NEK+1)) C **** C **** S3 = (TOTAL N2D LOSS)/N(N2D) C **** S3(I,1) = S1(I,1)*(BETA2*F(I,NPS1K)*rmassinv_o2+BETA4* 1 F(I,NPO1K)*rmassinv_o+BETA6*F(I,NPNOK)/RMNO) 2 +BETA7+BETA5(I,1)*SQRT(F(I,NEK)*F(I,NEK+1)) 3 +RK10*F(I,NOPK) C **** C **** F(NPN2DK) = PSI(N2D) C **** F(I,NPN2DK)=RMN2D*S2(I,1)/(S3(I,1)*S1(I,1)) E5200(I,1) = BETA7*F(I,NPN2DK) E5200(I,1) = BETA7*S2(I,1)/S3(I,1) ! ! 12/11/00: this merge commented out as per kibo12: ! F(I,NPN2DK)=merge(F(I,NPN2DK),1.E-60,F(I,NPN2DK)-1.E-60>=0.) ! if (f(i,npn2dk) < 1.e-60) f(i,npn2dk) = 1.e-60 3 CONTINUE C **** C **** METASTABLE FOR O(1D),O2(1DELTA),AND O2(1SIGMA) C **** NO2PK = NJ+NO2P NPO1K = NJ+NPO1 NPO3K = NJ+NPO3 NPH2OK = NJ+NPH2O NPH2K = NJ+NPH2 NPCH4K = NJ+NPCH4 NMSK = NJ+NMS NJO3DK = NJO3D-1 NJO2DK = NJO2D-1 NPCO2K = NJ+NPCO2 NPN2DK = NPN2D+NJNP NPS1K = NJ+NPS ! O2 NPS2K = NPS1K+KMAXP1 ! OX NEK = NJ+NE DO 11 K=1,KMAX NJO3DK = NJO3DK+1 NJO2DK = NJO2DK+1 C **** C **** S1 = EXPS(K+1/2), S2 = J(O3 (1D)) ,S3 = J(O2 (1D)) C **** DO 11 I=1,LEN1 S2(I,K)=.5*(F(I,NJO3DK)+F(I,NJO3DK+1)) S3(I,K)=.5*(F(I,NJO2DK)+F(I,NJO2DK+1)) S10(I,K) = .5*(XJ762(I,K)+XJ762(I,K+1)) ! xj762 set to 0 in qrj S11(I,K) = .5*(XJN2O(I,K)+XJN2O(I,K+1)) ! xjn2o set in qrj S12(I,K) = .5*(XJCO2D(I,K)+XJCO2D(I,K+1)) ! xjco2d set in qrj 11 CONTINUE ntk = nt+nj DO 33 I=1,LEN2 C **** C **** TOTAL O(1D) PRODUCTION C **** S4(I,1) = S1(I,1)*(S3(I,1)*F(I,NPS1K)*rmassinv_o2+S2(I,1)* 1 F(I,NPO3K)/RMO3+BETA2*F(I,NPN2DK)/RMN2D 2 *F(I,NPS1K)*rmassinv_o2*0.1*S1(I,1))+0.85*RA2(I,1)* 3 F(I,NO2PK)*SQRT(F(I,NEK)*F(I,NEK+1))+ 4 RMN2O(I,1)*S1(I,1)/(.5*(F(I,NMSK)+F(I,NMSK+1))) 5 *S11(I,1)+S12(I,1)*F(I,NPCO2K)*S1(I,1)/RMCO2 C **** C **** TOTAL O(1D) LOSS C **** S5(I,1) = S1(I,1)*(RKM2A(I,1)*F(I,NPS1K)*rmassinv_o2+RKM1(I,1)* 1 (1.-F(I,NPS1K)-F(I,NPS2K))*rmassinv_n2+RKM3 2 *F(I,NPH2OK)/RMH2O+RKM4*F(I,NPH2K)/RMH2 3 +(RKM5A+RKM5B)*F(I,NPCH4K) 4 /RMCH4+RKM6(I,1)*F(I,NPCO2K)/RMCO2+RKM7A 5 *F(I,NPO3K)/RMO3+RKM8*F(I,NPO1K)*rmassinv_o 6 +GAM15*F(I,NPCO2K)/RMCO2)+A1D C **** C **** CALCULATE O(1D) IN PHOTCHEMICAL EQUILIBRIUM C **** XNO1D(I,1) = PSI(O(1D)) C **** (xno1d is in crates_tdep.h) XNO1D(I,1) = rmass_o*S4(I,1)/(S5(I,1)*S1(I,1)) ! 12/11/00: as per kibo12: ! XNO1D(I,1) = merge(XNO1D(I,1),1.E-60,XNO1D(I,1)-1.E-60>=0.) if (xno1d(i,1) < 1.e-60) xno1d(i,1) = 1.e-60 C **** C **** C **** TOTAL O2(1SIGMA) PRODUCTION C **** S6(I,1) = (RKM2A(I,1)*XNO1D(I,1)*rmassinv_o*F(I,NPS1K)*rmassinv_o2 1 *S1(I,1)*0.75+S10(I,1)*F(I,NPS1K)*rmassinv_o2) 2 *S1(I,1) C **** C **** TOTAL O2(1SIGMA) LOSS C **** S13(I,1) = S1(I,1)*(RKM13*(1.-F(I,NPS1K)-F(I,NPS2K))* 1 rmassinv_n2+RKM14*F(I,NPCO2K)/RMCO2+RKM15 2 *F(I,NPO3K)/RMO3+RKM16*F(I,NPO1K)*rmassinv_o 3 +RKM17*F(I,NPS1K)*rmassinv_o2) S7(I,1) = S13(I,1)+ASG C **** C **** CALCULATE O2(1S) IN PHOTCHEMICAL EQUILIBRIUM C **** XNO21S(I,1) = PSI(O2(1SIGMA)) C **** ! 10/23/01: Save o21s in number density rather than mmr (see addfsech below) XNO21S(I,1) = rmass_o2*S6(I,1)/(S7(I,1)*S1(I,1)) ! MMR xno21s_cm3(i,1) = s6(i,1)/s7(i,1) ! CM3 ! 12/11/00: as per kibo12: if (xno21s(i,1) < 1.e-60) xno21s(i,1) = 1.e-60 if (xno21s_cm3(i,1) < 1.e-60) xno21s_cm3(i,1) = 1.e-60 ! ! subroutine mkeo200(tn,xo2,xo,xn2,id,fout) ! rk1 = 4.7e-33*(300./tn(k))**2 ! fout(k) = (a1*rk1*xo(k)**2 * (xn2(k)+xo2(k)) * xo2(k)) / ! + ((a2 + rk2a*xo2(k) + rk2b*xn2(k) + rk2c*xo(k)) * ! + (dp*xo2(k)+dpp*xo(k))) ! ! 10/12/01 bf: use O1 (NPO1K) rather than OX (NPS2K) in eo200 calculation: ! xo2(I,1) = xnmbar(I,1)*F(I,NPS1K)*rmassinv_o2 ! xo1(I,1) = xnmbar(I,1)*F(I,NPS2K)*rmassinv_o xo1(I,1) = xnmbar(I,1)*F(I,NPO1K)*rmassinv_o xn2(I,1) = xnmbar(I,1)*(1.-f(i,nps1k)-f(i,nps2k))*rmassinv_n2 rkk1 = 4.7e-33*(300./f(i,ntk))**2 eo200a(i,1) = (a1*rkk1*xo1(i,1)**2 * (xn2(i,1)+xo2(i,1)) * | xo2(i,1)) / ((a2 + rk2a*xo2(i,1) + rk2b*xn2(i,1) + rk2c* | xo1(i,1)) * (dp*xo2(i,1)+dpp*xo1(i,1))) C **** C **** C **** C **** TOTAL O2(1DELTA) PRODUCTION C **** S8(I,1) = S1(I,1)*0.9*S2(I,1)*F(I,NPO3K)/RMO3+S13(I,1)* 1 XNO21S(I,1)*rmassinv_o2*S1(I,1) C **** C **** TOTAL O2(1DELTA) LOSS C **** S9(I,1) = S1(I,1)*(RKM9(I,1)*F(I,NPS1K)*rmassinv_o2+RKM10* 1 (1.-F(I,NPS1K)-F(I,NPS2K))*rmassinv_n2+RKM11* 2 F(I,NPO1K)*rmassinv_o)+ADL C **** C **** CALCULATE O2(1D) IN PHOTCHEMICAL EQUILIBRIUM C **** XNO21D(I,1) = PSI(O2(1DLETA)) C **** XNO21D(I,1) = rmass_o2*S8(I,1)/(S9(I,1)*S1(I,1)) ! 12/11/00: as per kibo12: if (xno21d(i,1) < 1.e-60) xno21d(i,1) = 1.e-60 C RO1DN2D(I,1) = A1D*XNO1D(I,1)/(BETA7*F(I,NPN2DK)) E6300(I,1) = A1D*S4(I,1)/S5(I,1) C 33 CONTINUE C NPN2DK = NJNP+NPN2D call addfsech('EO200A' ,' ',' ',eo200a,zimxp,zkmxp,zkmx,j) ! call addfsech('O21SIG' ,' ',' ',XNO21S ,zimxp,zkmxp,zkmx,j) call addfsech('O21SIG' ,' ',' ',xno21s_cm3,zimxp,zkmxp,zkmx,j) call addfsech('E6300' ,' ',' ',E6300 ,zimxp,zkmxp,zkmx,j) call addfsech('E5200' ,' ',' ',E5200 ,zimxp,zkmxp,zkmx,j) call addfsech('RO1DN2D' ,' ',' ',RO1DN2D,zimxp,zkmxp,zkmx,j) C NPO21DK = NJNP+NPO21D-1 DO 333 K=1,KMAX NPO21DK=NPO21DK+1 DO 333 I=1,LEN1 F(I,NPO21DK) = XNO21D(I,K) 333 CONTINUE C **** C **** SET LEVEL KMAX+1 USING LOGARITHMIC INTERPOLATION ! 6/26/98 btf: do same for npn2d. C **** NPO21DK = NJNP+NPO21D+KMAX npn2dk = njnp+npn2d+kmax DO I = 1,LEN1 F(I,NPO21DK) = F(I,NPO21DK-1)**2/F(I,NPO21DK-2) f(i,npn2dk) = f(i,npn2dk-1)**2/f(i,npn2dk-2) ENDDO NTK = NJ+NT-1 NMSK = NJ+NMS-1 NPS1K = NJ+NPS-1 NPS2K = NJ+NPS2-1 NOPK = NJ+NOP-1 NPH2OK = NJ+NPH2O-1 NPOHK = NJ+NPOH-1 NPHO2K = NJ+NPHO2-1 NPH2K = NJ+NPH2-1 NPHK = NJ+NPH-1 NPHOXK = NJ+NPHOX-1 NPO1K = NJ+NPO1-1 NPO3K = NJ+NPO3-1 NPNOK = NJ+NPNO-1 NPCOK = NJ+NPCO-1 NPN4SK = NJ+NPN4S-1 NPCH4K = NJ+NPCH4-1 DO 7 K=1,KMAX NTK = NTK+1 NMSK = NMSK+1 NPS1K = NPS1K+1 NPS2K = NPS2K+1 NOPK = NOPK+1 NPH2OK = NPH2OK+1 NPOHK = NPOHK+1 NPHO2K = NPHO2K+1 NPH2K = NPH2K+1 NPHK = NPHK+1 NPHOXK = NPHOXK+1 NPO1K = NPO1K+1 NPO3K = NPO3K+1 NPNOK = NPNOK+1 NPCOK = NPCOK+1 NPN4SK = NPN4SK+1 NPCH4K = NPCH4K+1 DO 7 I=1,LEN1 C **** C **** S3 = N*MBAR (K+1/2) C **** S3(I,K) = XNMBAR(I,K) C **** C **** PARTITIONING OF HOX: S7=OH/H RATIO, C **** S8=HO2/H RATIO AND S9=H/HOX NUMBER DENSITY RATIOS C **** S11(I,K) = .5*(F(I,NMSK)+F(I,NMSK+1)) C **** S13(I,K) = (RKM26(I,K)*F(I,NPO1K)/RMO1+RKM34(I,K)*F(I,NPO3K) 1 /RMO3+BETA10(I,K)*F(I,NPNOK)/RMNO)*S3(I,K) 2 +(XJHO2(I,K)+XJHO2(I,K+1))*0.5 S7(I,K) = (RKM36(I,K)*F(I,NPS1K)*rmassinv_o2*S3(I,K)*S3(I,K)/ 1 S11(I,K)+RKM37(I,K)*F(I,NPO3K)/RMO3*S3(I,K))/ 2 ((RKM25(I,K)*F(I,NPO1K)/RMO1+RKM33(I,K)* 3 F(I,NPH2K)/RMH2+BETA8*F(I,NPN4SK)/RMN4S 4 +RKM42(I,K)*F(I,NPCOK)/RMCO)*S3(I,K)) S8(I,K) = (RKM36(I,K)*F(I,NPS1K)*rmassinv_o2*S3(I,K)*S3(I,K)/ 1 S11(I,K))/S13(I,K)+(RKM29(I,K)*F(I,NPO3K)/RMO3) 3 *S3(I,K)*S7(I,K)/S13(I,K) S9(I,K) = 1./(1.+S7(I,K)+S8(I,K)) S10(I,K) = S9(I,K)*(RMH+S7(I,K)*RMOH+S8(I,K)*RMHO2) ! RATIO1(I,K) = merge(S7(I,K),1.E-6,S7(I,K)-1.E-6>=0.) ! RATIO2(I,K) = merge(S8(I,K),1.E-6,S8(I,K)-1.E-6>=0.) ! RATIO3(I,K) = merge(S9(I,K),1.E-6,S9(I,K)-1.E-6>=0.) ! RMTRU(I,K) = merge(S10(I,K),1.,S10(I,K)-1.>=0.) ratio1(i,k) = s7(i,k) if (s7(i,k) < 1.e-6) ratio1(i,k) = 1.e-6 ratio2(i,k) = s8(i,k) if (s8(i,k) < 1.e-6) ratio2(i,k) = 1.e-6 ratio3(i,k) = s9(i,k) if (s9(i,k) < 1.e-6) ratio3(i,k) = 1.e-6 rmtru(i,k) = s10(i,k) if (s10(i,k) < 1.) rmtru(i,k) = 1. C **** C **** CALCULATE N(H2O2) IN PHOTOCHEMICAL EQUILIBRIUM C **** C S5(I,K) = 0.5*(XJH2O2(I,K)+XJH2O2(I,K+1)) C XNH2O2(I,K) = RKM35(I,K)*(RATIO2(I,K)*RATIO3(I,K)* C 1 F(I,NPHOXK)/RMTRU(I,K)*S3(I,K))**2/ C 2 (S5(I,K)+S3(I,K)*(RKM27(I,K)*F(I,NPO1K)/RMO1+ C 3 RKM32(I,K)*F(I,NPHOXK)/RMTRU(I,K)*RATIO1(I,K) C 4 *RATIO3(I,K)))*RMH2O2/S3(I,K) C XNH2O2(I,K) = merge(XNH2O2(I,K),1.E-60,XNH2O2(I,K)-1.E-60 C 1 >=0.) ! if (xnh2o2(i,k) < 1.e-60) xnh2o2(i,k) = 1.e-60 XNH2O2(I,K) = 1.E-60 7 CONTINUE C NMSK = NJ+NMS-1 C NPS1K = NJ+NPS-1 C NPHOXK = NJ+NPHOX-1 C NPO1K = NJ+NPO1-1 C NPNOK = NJ+NPNO-1 C NPCH4K = NJ+NPCH4-1 C NJCH4TK = NJCH4T-1 DO 8 K=1,KMAX C NMSK = NMSK+1 C NPS1K = NPS1K+1 C NPHOXK = NPHOXK+1 C NPO1K = NPO1K+1 C NPNOK = NPNOK+1 C NPCH4K = NPCH4K+1 C NJCH4TK = NJCH4TK+1 DO 8 I=1,LEN1 C **** C **** METHANE OXIDATION CHEMISTRY C **** C **** S15 = OH , S14 = HO2 C S15(I,K) = F(I,NPHOXK)/RMTRU(I,K)*RATIO1(I,K)*RATIO3(I,K) C 1 *S3(I,K) C S14(I,K) = F(I,NPHOXK)/RMTRU(I,K)*RATIO2(I,K)*RATIO3(I,K) C 1 *S3(I,K) C S4(I,K) = 0.5*(F(I,NJCH4TK)+F(I,NJCH4TK+1)) C **** C **** C CH3(I,K) = (GAM1(I,K)*S15(I,K)+S4(I,K)+GAM13(I,K)* C 1 F(I,NPO1K)/RMO1*S3(I,K)+GAM2A*XNO1D(I,K) C 2 /RMO1*S3(I,K))*F(I,NPCH4K)/RMCH4*S3(I,K)/(( C 3 GAM3(I,K)*F(I,NPS1K)*rmassinv_o2*S3(I,K)/ C 4 (.5*(F(I,NMSK)+F(I,NMSK+1))) C 5 +GAM4*F(I,NPO1K)/RMO1)*S3(I,K)) C CH3O2(I,K) = GAM3(I,K)*F(I,NPS1K)*rmassinv_o2*S3(I,K)/ C 1 (.5*(F(I,NMSK)+F(I,NMSK+1))) C 2 *S3(I,K)*CH3(I,K)/(GAM5(I,K)*F(I,NPNOK)/RMNO C 3 *S3(I,K)+GAM6(I,K)*S14(I,K)-GAM8(I,K)*S15(I,K) C 4 *GAM6(I,K)*S14(I,K)/(XJCH3OO(I,K)+GAM8(I,K)* C 5 S15(I,K))) C CH3OOH(I,K) = GAM6(I,K)*S14(I,K)*CH3O2(I,K)/(XJCH3OO(I,K) C 1 +GAM8(I,K)*S15(I,K)) C CH3O(I,K) = (GAM5(I,K)*F(I,NPNOK)/RMNO*S3(I,K)*CH3O2(I,K)+ C 1 XJCH3OO(I,K)*CH3OOH(I,K))/(GAM9(I,K)* C 2 F(I,NPS1K)*rmassinv_o2*S3(I,K)) C CH2O(I,K) = (GAM2B*XNO1D(I,K)/RMO1*F(I,NPCH4K)/RMCH4 C 1 *S3(I,K)*S3(I,K)+GAM4*F(I,NPO1K)/RMO1*S3(I,K) C 2 *CH3(I,K)+GAM9(I,K)*F(I,NPS1K)*rmassinv_o2*S3(I,K)* C 3 CH3O(I,K))/(GAM10*S15(I,K)+XJCH2OA(I,K) C 4 +XJCH2OB(I,K)+GAM11(I,K)*F(I,NPO1K)/RMO1*S3(I,K)) C CHO(I,K) = (XJCH2OA(I,K)+GAM11(I,K)*F(I,NPO1K)/RMO1*S3(I,K) C 1 +GAM10*S15(I,K))*CH2O(I,K)/(GAM12(I,K)* C 2 F(I,NPS1K)*rmassinv_o2*S3(I,K)) CH3(I,K) = 1.E-60 CH3O2(I,K) = 1.E-60 CH3OOH(I,K) = 1.E-60 CH3O(I,K) = 1.E-60 CH2O(I,K) = 1.E-60 CHO(I,K) = 1.E-60 8 CONTINUE ! call fminmax(f(1,njnp+npn2d),zimxp*zkmxp,fmin,fmax) ! write(6,"('cmpmeta j=',i2,' f(1,njnp+npn2d) min,max=',2e12.4)") ! | j,fmin,fmax ! call fminmax(f(1,njnp+npo21d),zimxp*zkmxp,fmin,fmax) ! write(6,"('cmpmeta j=',i2,' f(1,njnp+npo21d) min,max=',2e12.4)") ! | j,fmin,fmax RETURN END