c*********************************************************************** c * c This a version that is more general and should be applicable to * c other latitudes. There are stIJl 0/0 problems, in particular what * c to do if the entire profIJe does not go to zero at once. * c * c*********************************************************************** subroutine hox(IJ,IK) #include "com2d.h" c hox is defined as oh+ho2+hco+h c this loop iterates to get self-consistent hox if(tfd(IJ) .gt. .10) then do 2334 ih1=1,5 c make sure non-zero (after polar night) if(c(13,IJ,IK) .gt. 0.0) then RAT(IJ,IK)=C(14,IJ,IK)/C(13,IJ,IK) RHOH(IJ,IK)=C(12,IJ,IK)/C(13,IJ,IK) else c(13,IJ,IK)=1.e4 c(14,IJ,IK)=c(13,IJ,IK) rat(IJ,IK)=1. rhoh(IJ,IK)=.01 c(12,IJ,IK)=.01*c(13,IJ,IK) endif c set up quadratic terms for HOx; solve for OH by eliminating HO2 and c H using ratios HO2/OH (RAT) and H/OH (RHOH) c loss proportional to oh*oh G7=-2.*((K(40,IJ,IK))*RAT(IJ,IK) c +K(59,IJ,IK)+K(12,IJ,IK)*M(IJ,IK)+(K(19,IJ,IK)+k(64,IJ,IK)* c m(IJ,IK))*RAT(IJ,IK)*RAT(IJ,IK)+(k(71,IJ,IK)+k(72,IJ,IK))* c rhoh(IJ,IK)*rat(IJ,IK)) c G7=-2.*(K(59,IJ,IK)+K(12,IJ,IK)*M(IJ,IK) c 1 +(K(40,IJ,IK)+K(65,IJ,IK)*M(IJ,IK) c 2 +(K(19,IJ,IK)+K(64,IJ,IK)*M(IJ,IK)) c 3 *RAT(IJ,IK) c 4 +(K(71,IJ,IK)+K(72,IJ,IK))*RHOH(IJ,IK)) c 5 *RAT(IJ,IK)) c loss proportional to oh G8=-(K(14,IJ,IK)*C(18,IJ,IK) 1 +K(58,IJ,IK)*C(26,IJ,IK) c +K(18,IJ,IK)*C(6,IJ,IK)*M(IJ,IK) c +K(37,IJ,IK)*C(10,IJ,IK) C + K(27,IJ,IK)*C(29,IJ,IK) 1 +K(55,IJ,IK)*C(32,IJ,IK) 2 +K(13,IJ,IK)*C(30,IJ,IK) C +K(61,IJ,IK)*C(28,IJ,IK) 2 +K(62,IJ,IK)*C(25,IJ,IK) 3 +(K(22,IJ,IK)*C(21,IJ,IK) C +K(53,IJ,IK)*C(27,IJ,IK) 4 +K(6,IJ,IK)*C(28,IJ,IK) C +K(108,ij,ik)*C(44,ij,ik) c Add BrO + HO2 --> HOBr + O2 K(108) C +K(34,IJ,IK)*C(6,IJ,IK)*M(IJ,IK)) 5 *RAT(IJ,IK)) c production by dissociation G9a=2.*J(4,IJ,IK)*C(15,IJ,IK) 1 +2.*J(8,IJ,IK)*C(16,IJ,IK) 2 +J(13,IJ,IK)*C(26,IJ,IK) 3 +J(6,IJ,IK)*C(10,IJ,IK) 4 +J(18,IJ,IK)*C(29,IJ,IK) 5 +J(10,IJ,IK)*C(23,IJ,IK) 6 +J(23,IJ,IK)*C(25,IJ,IK) 7 +J(24,IJ,IK)*C(32,IJ,IK) c production by reaction with o(1d) and o, and decomposition of HO2NO2, and photolysis of HOBr --> OH + Br g9b=(2.*(K(39,IJ,IK)*C(15,IJ,IK) 1 +K(48,IJ,IK)*C(17,IJ,IK)) 2 +(K(49,IJ,IK)+K(110,IJ,IK))*C(18,IJ,IK))*C(2,IJ,IK) 1 +J(62,ij,ik)*C(68,ij,ik) g9c=(K(21,IJ,IK)*C(23,IJ,IK)+K(33,IJ,IK)* 1 C(16,IJ,IK) 1 +K(54,IJ,IK)*C(10,IJ,IK))*C(1,IJ,IK) g9c=(K(21,IJ,IK)*C(23,IJ,IK)+K(33,IJ,IK)* 1 C(16,IJ,IK))*C(1,IJ,IK) g9d=(K(23,IJ,IK)*C(17,IJ,IK) 1 +K(7,IJ,IK)*C(16,IJ,IK))*C(27,IJ,IK) 2 +K(17,IJ,IK)*C(3,IJ,IK)*C(22,IJ,IK) 3 +K(70,IJ,IK)*M(IJ,IK)*C(32,IJ,IK) g9=g9a+g9b+g9c+g9d c solve a quadratic equation for OH c add half to half the old value for stabIJity C(13,IJ,IK) = (c(13,IJ,IK)+QUAD(G7,G8,G9))/2. c if(ij .eq. 9 .and. ik .le. 2)type *,ik,c(13,ij,ik),g7,rat(ij,ik),rhoh(ij,ik) c if(ij .eq. 1)type *,ik,rat(ij,ik), rhoh(ij,ik) c if(ij .eq. 1)type *,K(40,IJ,IK),k(65,IJ,IK),m(IJ,IK) c c ,K(59,IJ,IK),K(12,IJ,IK),K(19,IJ,IK),k(64,IJ,IK), c c k(71,IJ,IK),k(72,IJ,IK) c if(ij .eq. 1)type *,ik,g9a,g9b,g9c,g9d c if(ij .eq. 1)type *,J(4,IJ,IK),C(15,IJ,IK),J(8,IJ,IK),C(16,IJ,IK) c 2 ,J(13,IJ,IK),C(26,IJ,IK) c 3 ,J(6,IJ,IK),C(10,IJ,IK) c 4 ,J(18,IJ,IK),C(29,IJ,IK) c 5 ,J(10,IJ,IK),C(23,IJ,IK) c 6 ,J(23,IJ,IK),C(25,IJ,IK) c 7 ,J(24,IJ,IK),C(32,IJ,IK) c solve a single equation for H PROD=(K(41,IJ,IK)*C(1,IJ,IK)+K(36,IJ,IK)* 1 C(19,IJ,IK) 1 +K(30,IJ,IK)*C(17,IJ,IK))*C(13,IJ,IK) 2 +(K(48,IJ,IK)*C(2,IJ,IK) 3 +K(23,IJ,IK)*C(27,IJ,IK))*C(17,IJ,IK) 4 +J(4,IJ,IK)*C(15,IJ,IK) 5 +J(10,IJ,IK)*C(23,IJ,IK) 6 +J(18,IJ,IK)*C(29,IJ,IK) LOSS=K(3,IJ,IK)*C(3,IJ,IK)*M(IJ,IK) 1 +K(11,IJ,IK)*C(4,IJ,IK) 2 +(k(71,IJ,IK)+k(72,IJ,IK)+k(73,IJ,IK)) 6 *RAT(IJ,IK)*C(13,IJ,IK) C(12,IJ,IK)=PROD/LOSS c if(ij .eq. 9 .and. ik .le. 2)type *,c(12,ij,ik),prod,K(41,IJ,IK) c 1 ,C(1,IJ,IK),K(36,IJ,IK), c 1 C(19,IJ,IK), c 1 K(30,IJ,IK),C(17,IJ,IK),C(13,IJ,IK) c 2 ,K(48,IJ,IK),C(2,IJ,IK) c 3 ,K(23,IJ,IK),C(27,IJ,IK),C(17,IJ,IK) c 4 ,J(4,IJ,IK),C(15,IJ,IK) c 5 ,J(10,IJ,IK),C(23,IJ,IK) c 6 +J(18,IJ,IK)*C(29,IJ,IK) c set up terms for HO2 c loss proportional to HO2**2 G4=-2.*(K(19,IJ,IK)+k(64,IJ,IK)*m(IJ,IK)) c loss proportional to HO2 G5=-((K(40,IJ,IK))* 1 C(13,IJ,IK) 1 +K(38,IJ,IK)*C(5,IJ,IK) 2 +K(5,IJ,IK)*C(4,IJ,IK) 3 +K(42,IJ,IK)*C(1,IJ,IK) 4 +K(22,IJ,IK)*C(21,IJ,IK) 5 +K(53,IJ,IK)*C(27,IJ,IK) 5 +K(6,IJ,IK)*C(28,IJ,IK) 6 +K(34,IJ,IK)*C(6,IJ,IK)*M(IJ,IK) 7 +K(67,IJ,IK)*C(27,IJ,IK) 8 +(k(71,IJ,IK)+k(72,IJ,IK) 9 +k(73,IJ,IK))*c(12,IJ,IK) C +K(108,ij,ik)*C(44,ij,ik)) c Add BrO + HO2 --> HOBr + O2 K(108) c production of HO2 G6=K(3,IJ,IK)*C(12,IJ,IK)*C(3,IJ,IK)*M(IJ,IK) 1 +K(4,IJ,IK)*C(13,IJ,IK)*C(4,IJ,IK) 2 +K(7,IJ,IK)* C(27,IJ,IK)*C(16,IJ,IK) 3 +K(17,IJ,IK)*C(22,IJ,IK)* C(3,IJ,IK) 4 +K(60,IJ,IK)*C(13,IJ,IK)*C(28,IJ,IK) 5 +K(29,IJ,IK)*C(16,IJ,IK)*C(13,IJ,IK) 6 +K(33,IJ,IK)*C(16,IJ,IK)*C(1,IJ,IK) 7 +K(52,IJ,IK)*C(24,IJ,IK)*C(3,IJ,IK) 8 +K(70,IJ,IK)*C(32,IJ,IK)*M(IJ,IK) c solve a quadratic equation for HO2 c if(ij .eq. 9 .and. ik .le. 2)type *,c(14,ij,ik),g6 c 1 ,K(3,IJ,IK),C(12,IJ,IK),C(3,IJ,IK),M(IJ,IK) c 1 ,K(4,IJ,IK),C(13,IJ,IK),C(4,IJ,IK) c 2 ,K(7,IJ,IK), C(27,IJ,IK),C(16,IJ,IK) c 3 ,K(17,IJ,IK),C(22,IJ,IK),C(3,IJ,IK) c 4 ,K(60,IJ,IK),C(13,IJ,IK),C(28,IJ,IK) c 5 ,K(29,IJ,IK),C(16,IJ,IK),C(13,IJ,IK) c 6 ,K(33,IJ,IK),C(16,IJ,IK),C(1,IJ,IK) c 7 ,K(52,IJ,IK),C(24,IJ,IK),C(3,IJ,IK) c 8 ,K(70,IJ,IK),C(32,IJ,IK),M(IJ,IK) C(14,IJ,IK)=QUAD(G4,G5,G6) 2334 CONTINUE cn(13,IJ,IK)=c(13,IJ,IK) cn(12,IJ,IK)=c(12,IJ,IK) cn(14,IJ,IK)=c(14,IJ,IK) else c polar night! -- put all HOx into H2 or H2O before setting to zero, for H conservation C make it dependent on the relative amounts of H2 and H2O hoxall = (c(12,ij,ik) + c(13,ij,ik) + c(14,ij,ik))/2. wvrat = c(15,ij,ik)/(c(15,ij,ik) + c(17,ij,ik)) cn(15,ij,ik) = cn(15,ij,ik) + wvrat*hoxall h2rat = c(17,ij,ik)/(c(15,ij,ik) + c(17,ij,ik)) cn(17,ij,ik) = cn(17,ij,ik) + h2rat*hoxall cn(13,IJ,IK)=0.0 cn(12,IJ,IK)=0.0 cn(14,IJ,IK)=0.0 endif C H2O2 IF(TFD(IJ) .GT. 0.10) THEN CN(16,IJ,IK)=(K(19,IJ,IK)*C(14,IJ,IK)**2+K(12,IJ,IK)* C C(13,IJ,IK)**2*M(IJ,IK))/(J(8,IJ,IK)+K(7,IJ,IK) C *C(27,IJ,IK)+K(29,IJ,IK)*C(13,IJ,IK)+K(33,IJ,IK) C *C(1,IJ,IK)+2.4E-7*C(38,IJ,IK)) ELSE CN(16,IJ,IK)=C(16,IJ,IK)*EXP(-2.4E-7*C(38,IJ,IK)*86400.) ENDIF SAVE return END