SUBROUTINE HGAQSCHEM ( CHG0, CHG2, CO3, CSO2, COH, CHO2, CHCL, & CCL2, CPM, CPH, LWC, TEMP, PRES, DELT, & SO2H, SO2TMP, O3H, O3TMP ) C c----CAMx v7Beta6 190902 C C HGAQSCHEM performs the Hg aqueous-phase chemistry calculations. C C Copyright 2003 C AER, Inc. C C Revision History: C*********************************************************************** C Version 1.0 written March 2003 by Prakash Karamchandani, AER * C Code updated August 2003 to use Henry's Law and SO2 dissociation * C constants consistent with those used in the CAMx aqueous-phase * C chemistry module * C*********************************************************************** C C Called by: CHEMDRIV C C*********************************************************************** IMPLICIT NONE C Includes: INCLUDE 'hgaqschm.prm' ! Hg aqueous-phase chemistry parameters C Arguments: REAL CHG0 ! Hg(0) concentration, ppm (Input/Output) REAL CHG2 ! Hg(2) concentration, ppm (Input/Output) REAL CO3 ! O3 concentration, ppm (Input only) REAL CSO2 ! SO2 concentration, ppm (Input only) REAL COH ! OH concentration, ppm (Input only) REAL CHO2 ! HO2 concentration, ppm (Input only) REAL CHCL ! HCL concentration, ppm (Input only) REAL CCL2 ! CL2 concentration, ppm (Input only) REAL CPM ! Total PM concentration, ug/m3 (Input only) REAL CPH ! Cloudwater pH (Input only) REAL LWC ! Liquid water content, g/m3 (Input only) REAL TEMP ! Temperature, K (Input only) REAL PRES ! Pressure, mb (Input only) REAL DELT ! Chemistry time step, seconds (Input only) REAL SO2H ! Henry's Law constant for SO2, M/atm (Input only) REAL O3H ! Henry's Law constant for O3, M/atm (Input only) REAL SO2TMP ! SO2 Henry's Law temperature dependence (Input only) REAL O3TMP ! O3 Henry's Law temperature dependence (Input only) ! Local variables ! Standard atmosphere in mb REAL STDATMMB PARAMETER ( STDATMMB = 1013.25 ) ! Pressure in atmospheres REAL PRES_ATM ! Conv. factor from ppm to atmospheres REAL CFACT ! Conversion factor (M to atm), liter-atm/mole REAL XL ! Henry's Law constant for OH and HO2, M/atm REAL OHH, HO2H ! Temperature in K and Centigrade and inverse temperature (K-1) REAL TEMPK, TEMPC, TEMPI ! Double precision local variables (for accuracy): DOUBLE PRECISION Y00 ! Initial total Hg(0) partial pressure, atm DOUBLE PRECISION Y20 ! Initial total Hg(2) partial pressure, atm DOUBLE PRECISION Z0T ! Final total Hg(0) partial pressure, atm DOUBLE PRECISION Z2T ! Final total Hg(2) partial pressure, atm ! --- Local concentration array DOUBLE PRECISION Y( NUMSP ) ! Y(1) [Hg(0)]-Gaseous (variable), atm ! Y(2) [Hg(2)]-Gaseous (variable), atm ! Y(3) [O3]-Gaseous (constant), atm ! Y(4) [Cl2]-Gaseous (constant), atm ! Y(5) [HCl]-Gaseous (constant), atm ! Y(6) [SO2]-Gaseous (constant), atm ! Y(7) [OH]-Aqueous (constant), M ! Y(8) [HO2]-Aqueous (constant), M DOUBLE PRECISION PM ! PM10 conc. (g/L) ! --- Intermediate variables for calculating Hg(2) complexes with SO2 ! --- in solution and their reaction rates DOUBLE PRECISION CON1, CON2, CON3, CON4, CON51, CON52 ! --- Intermediate variables for calculating Hg(2) concentrations in ! --- solution as Hg(2)2+, HgCl2 and Hg(OH)2 and total Hg(2) DOUBLE PRECISION CON5, CON6, CON7, TOTHG ! --- Aqueous concentrations (M) of dissociation products of Cl2, ! --- HOCL and OCL- DOUBLE PRECISION CONHOCL, CONOCL ! --- Aqueous concentrations (M) of HCl and Cl2 DOUBLE PRECISION CONHCL, CONCL ! --- Intermediate variables for calculating Cl species distribution ! --- in solution DOUBLE PRECISION THCL, TCL2, CA, CB, CC, CD, AA, AB, AC, AD ! --- Aqueous equilibrium constants for SO2 and HSO3 ! --- SO2 <=> H+ + HSO3- (M) ! --- HSO3- <=> SO3-- + H+ (M) DOUBLE PRECISION EKSO2, EKHSO3 ! --- Henry's Law constants (M/atm) DOUBLE PRECISION HNRHG0, HNRHGCL2, HNRHGOH2 DOUBLE PRECISION HNRSO2, HNRO3, HNRHCL, HNRCL2 ! --- H+ concentration, OH- concentration DOUBLE PRECISION PHCON, POHCON DOUBLE PRECISION RLCONV ! Conversion factor (M to atm) ! --- Chemistry parameters DOUBLE PRECISION DT ! Chemistry time step, sec DOUBLE PRECISION ALPHA1, ALPHA2, ALPH12 ! Intermediate variables DOUBLE PRECISION DENOM1, DENOM2 ! --- Temperature-dependent rate for HGSO3 reaction DOUBLE PRECISION RKHGSO3I ! --- Assign values to local variables ! --- Local copy of temperature (don't allow temperature to go below ! --- freezing) TEMPK = MAX( TEMP , REAL(TZERO) ) ! --- Inverse temperature TEMPI = 1.0 / TEMPK ! --- Calculate partial pressures in atm. PRES_ATM = PRES / STDATMMB CFACT = 1.E-6 * PRES_ATM ! --- XL is the conversion factor from M (moles/liter of water) to atm XL = ( LWC / H2ODENS ) * ( MOLVOL * TEMPK / TZERO ) Y00 = CHG0 * CFACT Y20 = CHG2 * CFACT Y( KO3G ) = CO3 * CFACT Y( KSO2G ) = CSO2 * CFACT Y( KHCLG ) = CHCL * CFACT Y( KCL2G ) = CCL2 * CFACT ! --- Calculate liquid-phase concs of OH and HO2 radicals ! --- Henry's Law constants for OH and HO2 ! --- Values below are from Jacobson's book (1999) OHH = 25. * EXP( 17.72 * ( 298.15 * TEMPI - 1.0 ) ) HO2H = 2000. * EXP( 22.28 * ( 298.15 * TEMPI - 1.0 ) ) ! --- Since liquid-phase chemistry of OH and HO2 radicals is not explicitly ! --- simulated in CAMx, assume that in-cloud processes remove 50% of OH and ! --- 90% of HO2 Y( KOHAQ ) = COH * CFACT * OHH * 0.5 / ( 1.0 + OHH * XL ) Y( KHO2AQ ) = CHO2 * CFACT * HO2H * 0.1 / ( 1.0 + HO2H * XL ) ! --- H+ and OH- aqueous concentrations in moles/liter PHCON = 10. ** ( -CPH ) POHCON = 10. ** ( CPH - 14 ) ! --- concentration of PM in g/L of water PM = CPM * H2ODENS * 1.E-9 / LWC ! No need to check if LWC = 0. since ! this routine is not called if ! LWC < LWMIN ! --- Assign some double precision variables ! --- Conversion factor from M to atm. RLCONV = XL ! --- Chemistry time step DT = DELT ! --- Henry's Law constants for SO2 and O3 HNRO3 = O3H*exp(O3TMP*(1./298. - 1./tempk)) HNRSO2 = SO2H*exp(SO2TMP*(1./298. - 1./tempk)) EKSO2 = 10.**(853./tempk)/54950. EKHSO3 = 10.**(621.9/tempk)/1.897e+9 ! --- Henry's Law constants for Hg(0), HgCl2, Hg(OH)2, HCl and Cl2 ! --- (Not available from CAMx) ! --- Special expression for Hg(0) Henry's Law constant HNRHG0 = 54.78775913129319 * EXP( -55.7339 + 9540.36 * TEMPI + & 16.0477 * ALOG( TEMPK / 100. ) ) HNRHGCL2 = 1.4E6 HNRHGOH2 = 1.2E4 HNRHCL = 1.1 ! --- Special expression for Cl2 Henry's Law constant TEMPC = TEMPK - TZERO HNRCL2 = 0.149 + TEMPC * ( -3.59E-03 + & TEMPC * ( 2.5E-05 + TEMPC * 5.67E-08 ) ) ! --- Gas-liquid partitioning of reactant species Y( KO3G ) = Y( KO3G ) / ( 1.0 + RLCONV * HNRO3 ) ! --- Calculate distribution of chlorine species IF ( RLCONV .GT. 0. ) THEN ! Check for RLCONV not really required ! since this routine is only called for ! liquid water >= lwmin THCL = Y( KHCLG ) / RLCONV TCL2 = 2. * Y( KCL2G ) / RLCONV C C Calculate subcomponents CA = 2. * ( 1.0 + PHCON / EKHOCL ) CB = 2. * ( 1.0 + 1.0 / ( HNRCL2 * RLCONV ) ) * PHCON**2. / & EKCL2 / EKHOCL CC = 1.0 + PHCON / EKHCL * ( 1.0 + 1.0 / ( HNRHCL * RLCONV ) ) CD = TCL2 * (1.0 + PHCON / EKHOCL ) C C Solve quadratic equation AA = CB * CC AB = CA * CC - CB * THCL AC = -( CD + CA * THCL ) AD = AB * AB - 4. * AA * AC CONCL = ( -AB + SQRT( AD ) ) / (2. * AA ) CONHCL = CONCL * PHCON / EKHCL Y( KHCLG ) = CONHCL / HNRHCL Y( KCL2G ) = 0.5 * TCL2 / (1. / ( RLCONV ) + HNRCL2 * & ( 1.0 + EKCL2 / ( EKHCL * HNRHCL * Y( KHCLG ) ) * & ( 1.0 + EKHOCL / PHCON ) ) ) END IF ! ( RLCONV > 0. ) ! --- Calculate partitioning of S(IV), assuming that mercury sulfites ! --- are negligible compared to other sulfite ions. Y( KSO2G ) = Y( KSO2G ) / ( 1.0 + RLCONV * HNRSO2 * ( 1.0 + & EKSO2 / PHCON * ( 1.0 + EKHSO3 / PHCON ) ) ) ! --- Calculate constants CON2 = EKHGSO3 * EKSO2 * EKHSO3 * HNRSO2 / ( PHCON * PHCON ) CON3 = EKHGSO3 * EKHGSO32 * ( EKSO2 * EKSO2 ) * & ( EKHSO3 * EKHSO3 ) * ( HNRSO2 * HNRSO2 ) / ( PHCON**4 ) CON5 = ( 1.0 + Y( KSO2G ) * ( CON2 + CON3 * Y( KSO2G ) ) ) CON6 = ( ( EKHCL * HNRHCL * Y( KHCLG ) / PHCON )**2 ) / EKHGCL2 CON7 = ( POHCON * POHCON ) / EKHGOH2 CONHOCL = HNRCL2 * Y( KCL2G ) * EKCL2 / & ( EKHCL * HNRHCL * Y( KHCLG ) ) CONOCL = CONHOCL * EKHOCL / PHCON CON51 = CON2 * Y( KSO2G ) CON52 = CON3 * Y( KSO2G ) * Y( KSO2G ) DENOM1 = ( 1.0 + RLCONV * HNRHG0 ) TOTHG = CON5 + ( CON51 + CON52 ) * EKP * PM + ( CON6 + CON7 ) * 2 ( 1.0 + EKP * PM ) DENOM2 = CON6 / HNRHGCL2 + CON7 / HNRHGOH2 + RLCONV * TOTHG ! --- Hg aqueous chemistry calculations for time-step DT ! --- Calculate HGSO3 reaction rate (Van Loon et al., 2000) RKHGSO3I = RKHGSO3 * EXP( -105000. / 8.3 * ( TEMPI - TREFI ) ) CON1 = RKHGSO32 * CON3 CON4 = RKHGSO3I * CON2 ALPHA1 = ( RLCONV * HNRHG0 * ( RKO3 * HNRO3 * Y( KO3G ) + & RKOH * Y( KOHAQ ) + & RKHOCL * CONHOCL + & RKOCL * CONOCL ) ) / DENOM1 ALPHA2 = RLCONV / DENOM2 * ( CON1 * ( Y( KSO2G ) * Y( KSO2G ) ) + & CON4 * Y( KSO2G ) + RKHO2 * Y( KHO2AQ ) * & ( CON5 + CON6 + CON7 ) ) ALPH12 = ALPHA1 + ALPHA2 ! --- Calculate analytical solution Z0T = ALPHA2 * ( Y00 + Y20 ) / ALPH12 * & (1.0 - EXP( -ALPH12 * DT ) ) + Y00 * EXP( -ALPH12 * DT ) Z2T = Y00 + Y20 - Z0T CHG0 = Z0T / CFACT CHG2 = Z2T / CFACT RETURN END