MODULE module_ctrans_aqchem CONTAINS !*********************************************************************** ! Portions of Models-3/CMAQ software were developed or based on * ! information from various groups: Federal Government employees, * ! contractors working on a United States Government contract, and * ! non-Federal sources (including research institutions). These * ! research institutions have given the Government permission to * ! use, prepare derivative works, and distribute copies of their * ! work in Models-3/CMAQ to the public and to permit others to do * ! so. EPA therefore grants similar permissions for use of the * ! Models-3/CMAQ software, but users are requested to provide copies * ! of derivative works to the Government without restrictions as to * ! use by others. Users are responsible for acquiring their own * ! copies of commercial software associated with Models-3/CMAQ and * ! for complying with vendor requirements. Software copyrights by * ! the MCNC Environmental Modeling Center are used with their * ! permissions subject to the above restrictions. * !*********************************************************************** ! RCS file, release, date & time of last delta, author, state, [and locker] ! $Header: /project/work/rep/CCTM/src/cloud/cloud_acm/aqchem.F,v 1.32 2008/09/10 19:40:39 sjr Exp $ ! what(1) key, module and SID; SCCS file; date and time of last delta: SUBROUTINE AQCHEM ( TEMP, PRES_PA, TAUCLD, PRCRATE, & WCAVG, WTAVG, AIRM, ALFA0, ALFA2, ALFA3, GAS, & AEROSOL, LIQUID, GASWDEP, AERWDEP, HPWDEP ) !----------------------------------------------------------------------- ! ! DESCRIPTION: ! Compute concentration changes in cloud due to aqueous chemistry, ! scavenging and wet deposition amounts. ! ! Revision History: ! No Date Who What ! -- -------- --- ----------------------------------------- ! 0 / /86 CW BEGIN PROGRAM - Walceks's Original Code ! 1 / /86 RB INCORPORATE INTO RADM ! 2 03/23/87 DH REFORMAT ! 3 04/11/88 SJR STREAMLINED CODE - ADDED COMMENTS ! 4 08/27/88 SJR COMMENTS, MODIFIED FOR RPM ! 4a 03/15/96 FSB Scanned hard copy to develop Models3 ! Version. ! 5 04/24/96 FSB Made into Models3 Format ! 6 02/18/97 SJR Revisions to link with Models3 ! 7 08/12/97 SJR Revised for new concentration units (moles/mole) ! and new treatment of nitrate and nitric acid ! 8 01/15/98 sjr revised to add new aitken mode scavenging ! and aerosol number scavenging ! 9 12/15/98 David Wong at LM: ! -- change division of XL, TEMP to multiplication of XL, TEMP ! reciprocal, respectively ! -- change / TOTOX / TSIV to / ( TOTOX * TSIV ) ! 10 03/18/99 David Wong at LM: ! -- removed "* 1.0" redundant calculation at TEMP1 calculation ! 11 04/27/00 sjr Added aerosol surface area as modeled species ! 12 12/02 sjr changed calls to HLCONST and updated the dissociation ! constants ! 13 06/26/03 sjr revised calculations of DTW based on CMAS website ! discussions ! 14 08/05/03 sjr revision made to the coarse aerosol number washout ! 15 04/20/05 us revisions to add sea salt species in the fine and ! coarse aerosol modes, and HCl dissolution/dissociation ! 16 10/13/05 sjr fixed bug in the integration time step calculation ! (reported by Bonyoung Koo) ! 17 03/01/06 sjr added elemental carbon aerosol; organic aerosols ! replaced with primary, secondary biogenic, and ! secondary anthropogenic; fixed 3rd moment calc to ! include EC and primary organics (not secondary); ! re-arranged logic for setting Cl & Na ending conc; ! added pointers/indirect addressing for arrays WETDEP ! and LIQUID ! 16 03/30/07 sjr Limit integration timestep by cloud washout time ! 17 04/10/07 sjr increased loop limits as follows: I20C <10000, ! I7777C <10000, I30C <10000, ICNTAQ <60000 ! ! ! Reference: ! Walcek & Taylor, 1986, A theoretical Method for computing ! vertical distributions of acidity and sulfate within cumulus ! clouds, J. Atmos Sci., Vol. 43, no. 4 pp 339 - 355 ! ! Called by: AQMAP ! ! Calls the following subroutines: none ! ! Calls the following functions: HLCONST ! ! ARGUMENTS TYPE I/O DESCRIPTION ! --------- ---- ------------ -------------------------------- ! GAS(ngas) real input&output Concentration for species i=1,12 ! GASWDEP(ngas) real output wet deposition for species ! (01)= SO2 conc (mol/mol) ! (02)= HNO3 conc (mol/mol) ! (03)= N2O5 conc (mol/mol) ! (04)= CO2 conc (mol/mol) ! (05)= NH3 conc (mol/mol) ! (06)= H2O2 conc (mol/mol) ! (07)= O3 conc (mol/mol) ! (08)= FOA conc (mol/mol) ! (09)= MHP conc (mol/mol) ! (10)= PAA conc (mol/mol) ! (11)= H2SO4 conc (mol/mol) ! (12)= HCL conc (mol/mol) ! ! AEROSOL(naer) real input&output Concentration for species i=1,36 ! AERWDEP(naer) real output wet deposition for species ! (01)= SO4AKN conc (mol/mol) ! (02)= SO4ACC conc (mol/mol) ! (03)= SO4COR conc (mol/mol) ! (04)= NH4AKN conc (mol/mol) ! (05)= NH4ACC conc (mol/mol) ! (06)= NO3AKN conc (mol/mol) ! (07)= NO3ACC conc (mol/mol) ! (08)= NO3COR conc (mol/mol) ! (09)= ORGAAKN conc (mol/mol) ! (10)= ORGAACC conc (mol/mol) ! (11)= ORGPAKN conc (mol/mol) ! (12)= ORGPACC conc (mol/mol) ! (13)= ORGBAKN conc (mol/mol) ! (14)= ORGBACC conc (mol/mol) ! (15)= ECAKN conc (mol/mol) ! (16)= ECACC conc (mol/mol) ! (17)= PRIAKN conc (mol/mol) ! (18)= PRIACC conc (mol/mol) ! (19)= PRICOR conc (mol/mol) ! (20)= NAAKN conc (mol/mol) ! (21)= NAACC conc (mol/mol) ! (22)= NACOR conc (mol/mol) ! (23)= CLAKN conc (mol/mol) ! (24)= CLACC conc (mol/mol) ! (25)= CLCOR conc (mol/mol) ! (26)= NUMAKN conc ( #/mol ) ! (27)= NUMACC conc ( #/mol ) ! (28)= NUMCOR conc ( #/mol ) ! (29)= SRFAKN conc (m2/mol ) ! (30)= SRFACC conc (m2/mol ) ! (31)= NACL conc (mol/mol) ! (32)= CACO3 conc (mol/mol) ! (33)= MGCO3 conc (mol/mol) ! (34)= A3FE conc (mol/mol) ! (35)= B2MN conc (mol/mol) ! (36)= K conc (mol/mol) !----------------------------------------------------------------------- IMPLICIT NONE ! INCLUDE SUBST_IOPARMS ! I/O parameters definitions ! INCLUDE SUBST_RXCMMN ! Mechanism reaction common block !....................................................................... ! INCLUDE FILE CONST.EXT ! Contains: Fundamental constants for air quality modeling ! Dependent Upon: none ! Revision History: ! Adapted 6/92 by CJC from ROM's PI.EXT. ! 3/1/93 John McHenry - include constants needed by LCM aqueous chemistry ! 9/93 by John McHenry - include additional constants needed for FMEM clouds ! and aqueous chemistry ! 3/4/96 Dr. Francis S. Binkowski - reflect current Models3 view that MKS ! units should be used wherever possible and that sources be documented. ! Some variables have been added, names changed, and values revised. ! 3/7/96 - add universal gas constant and compute gas constant in chemical ! form. TWOPI is now calculated rather than input. ! 3/13/96 - group declarations and parameter statements ! 9/13/96 - include more physical constants ! 12/24/96 - eliminate silly EPSILON, AMISS ! 1/06/97 - eliminate most derived constants - YOJ ! 1/17/97 (comments only) to provide numerical values as reference - DWB ! 4/30/08 - Changed REARTH to match default value in MM5 and WRF - TLO ! FSB References: ! CRC76, "CRC Handbook of Chemistry and Physics (76th Ed)", ! CRC Press, 1995 ! Hobbs, P.V. "Basic Physical Chemistry for the Atmospheric Sciences", ! Cambridge Univ. Press, 206 pp, 1995. ! Snyder, J.P., "Map Projections-A Working Manual, U.S. Geological Survey ! Paper 1395 U.S.GPO, Washington, DC, 1987. ! Stull, R. B., "An Introduction to Bounday Layer Meteorology", Kluwer, ! Dordrecht, 1988 !....................................................................... ! Geometric Constants: REAL PI ! pi (single precision 3.141593) PARAMETER ( PI = 3.14159265358979324 ) REAL PI180 ! pi/180 [ rad/deg ] PARAMETER ( PI180 = PI / 180.0 ) ! Geodetic Constants: REAL REARTH ! radius of the earth [ m ] ! FSB: radius of sphere having same surface area as ! Clarke ellipsoid of 1866 ( Source: Snyder, 1987) ! PARAMETER ( REARTH = 6370997.0 ) PARAMETER ( REARTH = 6370000.0 ) ! default Re in MM5 and WRF REAL SIDAY ! length of a sidereal day [ sec ] ! FSB: Source: CRC76 pp. 14-6 PARAMETER ( SIDAY = 86164.09 ) REAL GRAV ! mean gravitational acceleration [ m/sec**2 ] ! FSB: Value is mean of polar and equatorial values. ! Source: CRC Handbook (76th Ed) pp. 14-6 PARAMETER ( GRAV = 9.80622 ) REAL DG2M ! latitude degrees to meters PARAMETER ( DG2M = REARTH * PI180 ) ! Solar Constant: REAL SOLCNST ! Solar constant [ W/m**2 ], p14-2 CRC76 PARAMETER ( SOLCNST = 1373.0 ) ! Fundamental Constants: ( Source: CRC76, pp. 1-1 to 1-6) REAL AVO ! Avogadro's Constant [ number/mol ] PARAMETER ( AVO = 6.0221367E23 ) REAL RGASUNIV ! universal gas constant [ J/mol-K ] PARAMETER ( RGASUNIV = 8.314510 ) REAL STDATMPA ! standard atmosphere [ Pa ] PARAMETER ( STDATMPA = 101325.0 ) REAL STDTEMP ! Standard Temperature [ K ] PARAMETER ( STDTEMP = 273.15 ) REAL STFBLZ ! Stefan-Boltzmann [ W/(m**2 K**4) ] PARAMETER ( STFBLZ = 5.67051E-8 ) ! FSB Non-MKS REAL MOLVOL ! Molar volume at STP [ L/mol ] Non MKS units PARAMETER ( MOLVOL = 22.41410 ) ! Atmospheric Constants: REAL MWAIR ! mean molecular weight for dry air [ g/mol ] ! FSB: 78.06% N2, 21% O2, and 0.943% A on a mole ! fraction basis ( Source : Hobbs, 1995) pp. 69-70 PARAMETER ( MWAIR = 28.9628 ) REAL RDGAS ! dry-air gas constant [ J / kg-K ] PARAMETER ( RDGAS = 1.0E3 * RGASUNIV / MWAIR ) ! 287.07548994 REAL MWWAT ! mean molecular weight for water vapor [ g/mol ] PARAMETER ( MWWAT = 18.0153 ) REAL RWVAP ! gas constant for water vapor [ J/kg-K ] PARAMETER ( RWVAP = 1.0E3 * RGASUNIV / MWWAT ) ! 461.52492604 ! FSB NOTE: CPD, CVD, CPWVAP and CVWVAP are calculated assuming dry air and ! water vapor are classical ideal gases, i.e. vibration does not contribute ! to internal energy. REAL CPD ! specific heat of dry air at constant pressure [ J/kg-K ] PARAMETER ( CPD = 7.0 * RDGAS / 2.0 ) ! 1004.7642148 REAL CVD ! specific heat of dry air at constant volume [ J/kg-K ] PARAMETER ( CVD = 5.0 * RDGAS / 2.0 ) ! 717.68872485 REAL CPWVAP ! specific heat for water vapor at constant pressure [ J/kg-K ] PARAMETER ( CPWVAP = 4.0 * RWVAP ) ! 1846.0997042 REAL CVWVAP ! specific heat for water vapor at constant volume [ J/kg-K ] PARAMETER ( CVWVAP = 3.0 * RWVAP ) ! 1384.5747781 REAL VP0 ! vapor press of water at 0 C [ Pa ] Source: CRC76 pp. 6-15 PARAMETER ( VP0 = 611.29 ) ! FSB The following values are taken from p. 641 of Stull (1988): REAL LV0 ! latent heat of vaporization of water at 0 C [ J/kg ] PARAMETER ( LV0 = 2.501E6 ) REAL DLVDT ! Rate of change of latent heat of vaporization with ! respect to temperature [ J/kg-K ] PARAMETER ( DLVDT = 2370.0 ) REAL LF0 ! latent heat of fusion of water at 0 C [ J/kg ] PARAMETER ( LF0 = 3.34E5 ) !...Aqueous species pointers INCLUDE File !...........PARAMETERS and their descriptions: INTEGER, PARAMETER :: NGAS = 12 ! number of gas-phase species for AQCHEM INTEGER, PARAMETER :: NAER = 36 ! number of aerosol species for AQCHEM INTEGER, PARAMETER :: NLIQS = 41 ! number of liquid-phase species in AQCHEM !...pointers for the AQCHEM array GAS INTEGER, PARAMETER :: LSO2 = 1 ! Sulfur Dioxide INTEGER, PARAMETER :: LHNO3 = 2 ! Nitric Acid INTEGER, PARAMETER :: LN2O5 = 3 ! Dinitrogen Pentoxide INTEGER, PARAMETER :: LCO2 = 4 ! Carbon Dioxide INTEGER, PARAMETER :: LNH3 = 5 ! Ammonia INTEGER, PARAMETER :: LH2O2 = 6 ! Hydrogen Perioxide INTEGER, PARAMETER :: LO3 = 7 ! Ozone INTEGER, PARAMETER :: LFOA = 8 ! Formic Acid INTEGER, PARAMETER :: LMHP = 9 ! Methyl Hydrogen Peroxide INTEGER, PARAMETER :: LPAA = 10 ! Peroxyacidic Acid INTEGER, PARAMETER :: LH2SO4 = 11 ! Sulfuric Acid INTEGER, PARAMETER :: LHCL = 12 ! Hydrogen Chloride !...pointers for the AQCHEM array AEROSOL INTEGER, PARAMETER :: LSO4AKN = 1 ! Aitken-mode Sulfate INTEGER, PARAMETER :: LSO4ACC = 2 ! Accumulation-mode Sulfate INTEGER, PARAMETER :: LSO4COR = 3 ! Coarse-mode Sulfate INTEGER, PARAMETER :: LNH4AKN = 4 ! Aitken-mode Ammonium INTEGER, PARAMETER :: LNH4ACC = 5 ! Accumulation-mode Ammonium INTEGER, PARAMETER :: LNO3AKN = 6 ! Aitken-mode Nitrate INTEGER, PARAMETER :: LNO3ACC = 7 ! Accumulation-mode Nitrate INTEGER, PARAMETER :: LNO3COR = 8 ! Coarse-mode Nitrate INTEGER, PARAMETER :: LORGAAKN = 9 ! Aitken-mode anthropogenic SOA INTEGER, PARAMETER :: LORGAACC = 10 ! Accumulation-mode anthropogenic SOA INTEGER, PARAMETER :: LORGPAKN = 11 ! Aitken-mode primary organic aerosol INTEGER, PARAMETER :: LORGPACC = 12 ! Accumulation-mode primary organic aerosol INTEGER, PARAMETER :: LORGBAKN = 13 ! Aitken-mode biogenic SOA INTEGER, PARAMETER :: LORGBACC = 14 ! Accumulation-mode biogenic SOA INTEGER, PARAMETER :: LECAKN = 15 ! Aitken-mode elemental carbon INTEGER, PARAMETER :: LECACC = 16 ! Accumulation-mode elemental carbon INTEGER, PARAMETER :: LPRIAKN = 17 ! Aitken-mode primary aerosol INTEGER, PARAMETER :: LPRIACC = 18 ! Accumulation-mode primary aerosol INTEGER, PARAMETER :: LPRICOR = 19 ! Coarse-mode primary aerosol INTEGER, PARAMETER :: LNAAKN = 20 ! Aitken-mode Sodium INTEGER, PARAMETER :: LNAACC = 21 ! Accumulation-mode Sodium INTEGER, PARAMETER :: LNACOR = 22 ! Coarse-mode Sodium INTEGER, PARAMETER :: LCLAKN = 23 ! Aitken-mode Chloride ion INTEGER, PARAMETER :: LCLACC = 24 ! Accumulation-mode Chloride ion INTEGER, PARAMETER :: LCLCOR = 25 ! Coarse-mode Chloride ion INTEGER, PARAMETER :: LNUMAKN = 26 ! Aitken-mode number INTEGER, PARAMETER :: LNUMACC = 27 ! Accumulation-mode number INTEGER, PARAMETER :: LNUMCOR = 28 ! Coarse-mode number INTEGER, PARAMETER :: LSRFAKN = 29 ! Aitken-mode surface area INTEGER, PARAMETER :: LSRFACC = 30 ! Accumulation-mode surface area INTEGER, PARAMETER :: LNACL = 31 ! Sodium Chloride aerosol for AE3 only {depreciated in AE4} INTEGER, PARAMETER :: LCACO3 = 32 ! Calcium Carbonate aerosol (place holder) INTEGER, PARAMETER :: LMGCO3 = 33 ! Magnesium Carbonate aerosol (place holder) INTEGER, PARAMETER :: LA3FE = 34 ! Iron aerosol (place holder) INTEGER, PARAMETER :: LB2MN = 35 ! Manganese aerosol (place holder) INTEGER, PARAMETER :: LK = 36 ! Potassium aerosol (Cl- tracked separately) (place holder) !...pointers for the AQCHEM arrays LIQUID and WETDEP INTEGER, PARAMETER :: LACL = 1 ! Hydrogen ion INTEGER, PARAMETER :: LNH4L = 2 ! Ammonium INTEGER, PARAMETER :: LCAL = 3 ! Calcium INTEGER, PARAMETER :: LNAACCL = 4 ! Sodium INTEGER, PARAMETER :: LOHL = 5 ! Hydroxyl radical ion INTEGER, PARAMETER :: LSO4ACCL = 6 ! Sulfate (attributed to accumulation mode) INTEGER, PARAMETER :: LHSO4ACCL = 7 ! bisulfate (attributed to accumulation mode) INTEGER, PARAMETER :: LSO3L = 8 ! sulfite INTEGER, PARAMETER :: LHSO3L = 9 ! bisulfite INTEGER, PARAMETER :: LSO2L = 10 ! sulfur dioxide INTEGER, PARAMETER :: LCO3L = 11 ! carbonate INTEGER, PARAMETER :: LHCO3L = 12 ! bicarbonate INTEGER, PARAMETER :: LCO2L = 13 ! carbon dioxide INTEGER, PARAMETER :: LNO3ACCL = 14 ! nitrate(attributed to accumulation mode) INTEGER, PARAMETER :: LNH3L = 15 ! ammonia INTEGER, PARAMETER :: LCLACCL = 16 ! chloride ion (attributed to accumulation mode) INTEGER, PARAMETER :: LH2O2L = 17 ! hydrogen peroxide INTEGER, PARAMETER :: LO3L = 18 ! ozone INTEGER, PARAMETER :: LFEL = 19 ! iron INTEGER, PARAMETER :: LMNL = 20 ! Manganese INTEGER, PARAMETER :: LAL = 21 ! generalized anion associated with iron INTEGER, PARAMETER :: LFOAL = 22 ! Formic acid INTEGER, PARAMETER :: LHCO2L = 23 ! HCOO- ion INTEGER, PARAMETER :: LMHPL = 24 ! Methyl hydrogen peroxide INTEGER, PARAMETER :: LPAAL = 25 ! Peroxyacidic acid INTEGER, PARAMETER :: LHCLL = 26 ! Hydrogen chloride INTEGER, PARAMETER :: LPRIML = 27 ! primary aerosol INTEGER, PARAMETER :: LMGL = 28 ! Magnesium INTEGER, PARAMETER :: LKL = 29 ! potassium INTEGER, PARAMETER :: LBL = 30 ! generalized anion associated with manganese INTEGER, PARAMETER :: LHNO3L = 31 ! nitric acid INTEGER, PARAMETER :: LPRIMCORL = 32 ! coarse-mode primary aerosol INTEGER, PARAMETER :: LNUMCORL = 33 ! coarse-mode number INTEGER, PARAMETER :: LTS6CORL = 34 ! sulfate (attributed to coarse mode) INTEGER, PARAMETER :: LNACORL = 35 ! sodium (attributed to coarse mode) INTEGER, PARAMETER :: LCLCORL = 36 ! chloride ion (attributed to coarse mode) INTEGER, PARAMETER :: LNO3CORL = 37 ! nitrate (attributed to coarse mode) INTEGER, PARAMETER :: LORGAL = 38 ! anthropogenic SOA INTEGER, PARAMETER :: LORGPL = 39 ! primary organic aerosols INTEGER, PARAMETER :: LORGBL = 40 ! biogenic SOA INTEGER, PARAMETER :: LECL = 41 ! elemental carbon !...surrogate names, their background values, and units !... for AQCHEM's GAS species CHARACTER*16, SAVE :: SGRGAS ( NGAS ) ! surrogate name for gases CHARACTER*16, SAVE :: BUNTSGAS( NGAS ) ! units of bkgnd values REAL, SAVE :: BGNDGAS( NGAS ) ! background values for each gas DATA SGRGAS( LSO2 ), BGNDGAS( LSO2 ) /'SO2 ', 0.0 / DATA SGRGAS( LHNO3 ), BGNDGAS( LHNO3 ) /'HNO3 ', 0.0 / DATA SGRGAS( LN2O5 ), BGNDGAS( LN2O5 ) /'N2O5 ', 0.0 / DATA SGRGAS( LCO2 ), BGNDGAS( LCO2 ) /'CO2 ', 340.0 / DATA SGRGAS( LNH3 ), BGNDGAS( LNH3 ) /'NH3 ', 0.0 / DATA SGRGAS( LH2O2 ), BGNDGAS( LH2O2 ) /'H2O2 ', 0.0 / DATA SGRGAS( LO3 ), BGNDGAS( LO3 ) /'O3 ', 0.0 / DATA SGRGAS( LFOA ), BGNDGAS( LFOA ) /'FOA ', 0.0 / DATA SGRGAS( LMHP ), BGNDGAS( LMHP ) /'MHP ', 0.0 / DATA SGRGAS( LPAA ), BGNDGAS( LPAA ) /'PAA ', 0.0 / DATA SGRGAS( LH2SO4 ), BGNDGAS( LH2SO4 ) /'H2SO4 ', 0.0 / DATA SGRGAS( LHCL ), BGNDGAS( LHCL ) /'HCL ', 0.0 / DATA BUNTSGAS( LSO2 ) / 'ppm' / DATA BUNTSGAS( LHNO3 ) / 'ppm' / DATA BUNTSGAS( LN2O5 ) / 'ppm' / DATA BUNTSGAS( LCO2 ) / 'ppm' / DATA BUNTSGAS( LNH3 ) / 'ppm' / DATA BUNTSGAS( LH2O2 ) / 'ppm' / DATA BUNTSGAS( LO3 ) / 'ppm' / DATA BUNTSGAS( LFOA ) / 'ppm' / DATA BUNTSGAS( LMHP ) / 'ppm' / DATA BUNTSGAS( LPAA ) / 'ppm' / DATA BUNTSGAS( LH2SO4 ) / 'ppm' / DATA BUNTSGAS( LHCL ) / 'ppm' / !...surrogate names, their background values, units, and molecular weights !... for AQCHEM's AEROSOL species CHARACTER*16, SAVE :: SGRAER ( NAER ) ! surrogate name for aerosols CHARACTER*16, SAVE :: BUNTSAER( NAER ) ! units of bkgnd values REAL, SAVE :: SGRAERMW( NAER ) ! molecular weight for aerosol species REAL, SAVE :: BGNDAER ( NAER ) ! bkground vals each aerosols DATA SGRAER( LSO4AKN ), SGRAERMW( LSO4AKN ) / 'SO4_AITKEN ' , 96.0 / DATA SGRAER( LSO4ACC ), SGRAERMW( LSO4ACC ) / 'SO4_ACCUM ' , 96.0 / DATA SGRAER( LSO4COR ), SGRAERMW( LSO4COR ) / 'SO4_COARSE ' , 96.0 / DATA SGRAER( LNH4AKN ), SGRAERMW( LNH4AKN ) / 'NH4_AITKEN ' , 18.0 / DATA SGRAER( LNH4ACC ), SGRAERMW( LNH4ACC ) / 'NH4_ACCUM ' , 18.0 / DATA SGRAER( LNO3AKN ), SGRAERMW( LNO3AKN ) / 'NO3_AITKEN ' , 62.0 / DATA SGRAER( LNO3ACC ), SGRAERMW( LNO3ACC ) / 'NO3_ACCUM ' , 62.0 / DATA SGRAER( LNO3COR ), SGRAERMW( LNO3COR ) / 'NO3_COARSE ' , 62.0 / DATA SGRAER( LORGAAKN ), SGRAERMW( LORGAAKN ) / 'ORGA_AITKEN ' , 150.0 / DATA SGRAER( LORGAACC ), SGRAERMW( LORGAACC ) / 'ORGA_ACCUM ' , 150.0 / DATA SGRAER( LORGPAKN ), SGRAERMW( LORGPAKN ) / 'ORGP_AITKEN ' , 220.0 / DATA SGRAER( LORGPACC ), SGRAERMW( LORGPACC ) / 'ORGP_ACCUM ' , 220.0 / DATA SGRAER( LORGBAKN ), SGRAERMW( LORGBAKN ) / 'ORGB_AITKEN ' , 177.0 / DATA SGRAER( LORGBACC ), SGRAERMW( LORGBACC ) / 'ORGB_ACCUM ' , 177.0 / DATA SGRAER( LECAKN ), SGRAERMW( LECAKN ) / 'EC_AITKEN ' , 12.0 / DATA SGRAER( LECACC ), SGRAERMW( LECACC ) / 'EC_ACCUM ' , 12.0 / DATA SGRAER( LPRIAKN ), SGRAERMW( LPRIAKN ) / 'PRI_AITKEN ' , 200.0 / DATA SGRAER( LPRIACC ), SGRAERMW( LPRIACC ) / 'PRI_ACCUM ' , 200.0 / DATA SGRAER( LPRICOR ), SGRAERMW( LPRICOR ) / 'PRI_COARSE ' , 100.0 / DATA SGRAER( LNAAKN ), SGRAERMW( LNAAKN ) / 'NA_AITKEN ' , 23.0 / DATA SGRAER( LNAACC ), SGRAERMW( LNAACC ) / 'NA_ACCUM ' , 23.0 / DATA SGRAER( LNACOR ), SGRAERMW( LNACOR ) / 'NA_COARSE ' , 23.0 / DATA SGRAER( LCLAKN ), SGRAERMW( LCLAKN ) / 'CL_AITKEN ' , 35.5 / DATA SGRAER( LCLACC ), SGRAERMW( LCLACC ) / 'CL_ACCUM ' , 35.5 / DATA SGRAER( LCLCOR ), SGRAERMW( LCLCOR ) / 'CL_COARSE ' , 35.5 / DATA SGRAER( LNUMAKN ), SGRAERMW( LNUMAKN ) / 'NUM_AITKEN ' , 1.0 / DATA SGRAER( LNUMACC ), SGRAERMW( LNUMACC ) / 'NUM_ACCUM ' , 1.0 / DATA SGRAER( LNUMCOR ), SGRAERMW( LNUMCOR ) / 'NUM_COARSE ' , 1.0 / DATA SGRAER( LSRFAKN ), SGRAERMW( LSRFAKN ) / 'SRF_AITKEN ' , 1.0 / DATA SGRAER( LSRFACC ), SGRAERMW( LSRFACC ) / 'SRF_ACCUM ' , 1.0 / DATA SGRAER( LNACL ), SGRAERMW( LNACL ) / 'NACL ' , 58.4 / ! AE3 NaCl aerosol {depreciated in AE4} DATA SGRAER( LCACO3 ), SGRAERMW( LCACO3 ) / 'CACO3 ' , 100.1 / DATA SGRAER( LMGCO3 ), SGRAERMW( LMGCO3 ) / 'MGCO3 ' , 84.3 / DATA SGRAER( LA3FE ), SGRAERMW( LA3FE ) / 'A3FE ' , 55.8 / DATA SGRAER( LB2MN ), SGRAERMW( LB2MN ) / 'B2MN ' , 54.9 / DATA SGRAER( LK ), SGRAERMW( LK ) / 'K ' , 39.1 / DATA BGNDAER( LSO4AKN ), BUNTSAER( LSO4AKN ) / 0.0, 'ug/m3' / DATA BGNDAER( LSO4ACC ), BUNTSAER( LSO4ACC ) / 0.0, 'ug/m3' / DATA BGNDAER( LSO4COR ), BUNTSAER( LSO4COR ) / 0.0, 'ug/m3' / DATA BGNDAER( LNH4AKN ), BUNTSAER( LNH4AKN ) / 0.0, 'ug/m3' / DATA BGNDAER( LNH4ACC ), BUNTSAER( LNH4ACC ) / 0.0, 'ug/m3' / DATA BGNDAER( LNO3AKN ), BUNTSAER( LNO3AKN ) / 0.0, 'ug/m3' / DATA BGNDAER( LNO3ACC ), BUNTSAER( LNO3ACC ) / 0.0, 'ug/m3' / DATA BGNDAER( LNO3COR ), BUNTSAER( LNO3COR ) / 0.0, 'ug/m3' / DATA BGNDAER( LORGAAKN ), BUNTSAER( LORGAAKN ) / 0.0, 'ug/m3' / DATA BGNDAER( LORGAACC ), BUNTSAER( LORGAACC ) / 0.0, 'ug/m3' / DATA BGNDAER( LORGPAKN ), BUNTSAER( LORGPAKN ) / 0.0, 'ug/m3' / DATA BGNDAER( LORGPACC ), BUNTSAER( LORGPACC ) / 0.0, 'ug/m3' / DATA BGNDAER( LORGBAKN ), BUNTSAER( LORGBAKN ) / 0.0, 'ug/m3' / DATA BGNDAER( LORGBACC ), BUNTSAER( LORGBACC ) / 0.0, 'ug/m3' / DATA BGNDAER( LECAKN ), BUNTSAER( LECAKN ) / 0.0, 'ug/m3' / DATA BGNDAER( LECACC ), BUNTSAER( LECACC ) / 0.0, 'ug/m3' / DATA BGNDAER( LPRIAKN ), BUNTSAER( LPRIAKN ) / 0.0, 'ug/m3' / DATA BGNDAER( LPRIACC ), BUNTSAER( LPRIACC ) / 0.0, 'ug/m3' / DATA BGNDAER( LPRICOR ), BUNTSAER( LPRICOR ) / 0.0, 'ug/m3' / DATA BGNDAER( LNAAKN ), BUNTSAER( LNAAKN ) / 0.0, 'ug/m3' / DATA BGNDAER( LNAACC ), BUNTSAER( LNAACC ) / 0.0, 'ug/m3' / DATA BGNDAER( LNACOR ), BUNTSAER( LNACOR ) / 0.0, 'ug/m3' / DATA BGNDAER( LCLAKN ), BUNTSAER( LCLAKN ) / 0.0, 'ug/m3' / DATA BGNDAER( LCLACC ), BUNTSAER( LCLACC ) / 0.0, 'ug/m3' / DATA BGNDAER( LCLCOR ), BUNTSAER( LCLCOR ) / 0.0, 'ug/m3' / DATA BGNDAER( LNUMAKN ), BUNTSAER( LNUMAKN ) / 0.0, ' #/m3' / DATA BGNDAER( LNUMACC ), BUNTSAER( LNUMACC ) / 0.0, ' #/m3' / DATA BGNDAER( LNUMCOR ), BUNTSAER( LNUMCOR ) / 0.0, ' #/m3' / DATA BGNDAER( LSRFAKN ), BUNTSAER( LSRFAKN ) / 0.0, 'm2/m3' / DATA BGNDAER( LSRFACC ), BUNTSAER( LSRFACC ) / 0.0, 'm2/m3' / DATA BGNDAER( LNACL ), BUNTSAER( LNACL ) / 0.0, 'ug/m3' / ! AE3 NaCl aerosol {depreciated in AE4} DATA BGNDAER( LCACO3 ), BUNTSAER( LCACO3 ) / 0.0, 'ug/m3' / DATA BGNDAER( LMGCO3 ), BUNTSAER( LMGCO3 ) / 0.0, 'ug/m3' / DATA BGNDAER( LA3FE ), BUNTSAER( LA3FE ) / 0.010, 'ug/m3' / DATA BGNDAER( LB2MN ), BUNTSAER( LB2MN ) / 0.005, 'ug/m3' / DATA BGNDAER( LK ), BUNTSAER( LK ) / 0.0, 'ug/m3' / CHARACTER*120 XMSG ! Exit status message DATA XMSG / ' ' / !...........PARAMETERS and their descriptions: INTEGER NUMOX ! number of oxidizing reactions PARAMETER ( NUMOX = 5 ) REAL H2ODENS ! density of water at 20 C and 1 ATM PARAMETER ( H2ODENS = 1000.0 ) ! (kg/m3) REAL ONETHIRD ! 1/3 PARAMETER ( ONETHIRD = 1.0 / 3.0 ) REAL TWOTHIRDS ! 2/3 PARAMETER ( TWOTHIRDS = 2.0 / 3.0 ) REAL CONCMIN ! minimum concentration PARAMETER ( CONCMIN = 1.0E-30 ) REAL SEC2HR ! convert seconds to hours PARAMETER ( SEC2HR = 1.0 / 3600.0 ) !...........ARGUMENTS and their descriptions ! INTEGER JDATE ! current model date, coded YYYYDDD ! INTEGER JTIME ! current model time, coded HHMMSS REAL AIRM ! total air mass in cloudy layers (mol/m2) REAL ALFA0 ! scav coef for aitken aerosol number REAL ALFA2 ! scav coef for aitken aerosol sfc area REAL ALFA3 ! scav coef for aitken aerosol mass REAL HPWDEP ! hydrogen wet deposition (mm mol/liter) REAL PRCRATE ! precip rate (mm/hr) REAL PRES_PA ! pressure (Pa) REAL TAUCLD ! timestep for cloud (s) REAL TEMP ! temperature (K) REAL WCAVG ! liquid water content (kg/m3) REAL WTAVG ! total water content (kg/m3) REAL, INTENT(INOUT) :: GAS ( NGAS ) ! gas phase concentrations (mol/molV) REAL, INTENT(INOUT) :: AEROSOL( NAER ) ! aerosol concentrations (mol/molV) REAL, INTENT(OUT) :: LIQUID( NLIQS ) ! liquid concentrations (moles/liter) REAL, INTENT(OUT) :: GASWDEP( NGAS ) ! gas phase wet deposition array (mm mol/liter) REAL, INTENT(OUT) :: AERWDEP( NAER ) ! aerosol wet deposition array (mm mol/liter) !...........LOCAL VARIABLES (scalars) and their descriptions: LOGICAL, SAVE :: FIRSTIME = .TRUE. ! flag for first pass thru CHARACTER*6 PNAME ! driver program name DATA PNAME / 'AQCHEM' / SAVE PNAME CHARACTER( 16 ), SAVE :: AE_VRSN ! Aerosol version name INTEGER I20C ! loop counter for do loop 20 INTEGER I30C ! loop counter for do loop 30 INTEGER ITERAT ! # iterations of aqueous chemistry solver INTEGER I7777C ! aqueous chem iteration counter INTEGER ICNTAQ ! aqueous chem iteration counter INTEGER LIQ ! loop counter for liquid species INTEGER IOX ! index over oxidation reactions REAL DEPSUM REAL BETASO4 REAL A ! iron's anion concentration REAL AC ! H+ concentration in cloudwater (mol/liter) REAL ACT1 ! activity corretion factor!single ions REAL ACT2 ! activity factor correction!double ions REAL ACTB ! REAL AE ! guess for H+ conc in cloudwater (mol/liter) REAL B ! manganese's anion concentration REAL PRES_ATM ! pressure (Atm) REAL BB ! lower limit guess of cloudwater pH REAL CA ! Calcium conc in cloudwater (mol/liter) REAL CAA ! inital Calcium in cloudwater (mol/liter) REAL CL ! total Cl- conc in cloudwater (mol/liter) REAL CLACC ! fine Cl- in cloudwater (mol/liter) REAL CLACCA ! initial fine Cl in cloudwater (mol/liter) REAL CLAKNA ! initial interstitial aero Cl (mol/liter) REAL CLCOR ! coarse Cl- conc in cloudwater (mol/liter) REAL CLCORA ! init coarse Cl- in cloudwater (mol/liter) REAL CO2H ! Henry's Law constant for CO2 REAL CO21 ! First dissociation constant for CO2 REAL CO22 ! Second dissociation constant for CO2 REAL CO212 ! CO21*CO22 REAL CO212H ! CO2H*CO21*CO22 REAL CO21H ! CO2H*CO21 REAL CO2L ! CO2 conc in cloudwater (mol/liter) REAL CO3 ! CO3= conc in cloudwater (mol/liter) REAL CO3A ! initial CO3 in cloudwater (mol/liter) REAL CTHK1 ! cloud thickness (m) REAL DTRMV ! REAL DTS6 ! REAL EBETASO4T ! EXP( -BETASO4 * TAUCLD ) REAL EALFA0T ! EXP( -ALFA0 * TAUCLD ) REAL EALFA2T ! EXP( -ALFA2 * TAUCLD ) REAL EALFA3T ! EXP( -ALFA3 * TAUCLD ) REAL EC ! elemental carbon acc+akn aerosol in cloudwater (mol/liter) REAL ECACCA ! init EC ACC aerosol in cloudwater (mol/liter) REAL ECAKNA ! init EC AKN aerosol in cloudwater (mol/liter) REAL FA ! functional value ?? REAL FB ! functional value ?? REAL FE ! Fe+++ conc in cloudwater (mol/liter) REAL FEA ! initial Fe in cloudwater (mol/liter) REAL FNH3 ! frac weight of NH3 to total ammonia REAL FNH4ACC ! frac weight of NH4 acc to total ammonia REAL FHNO3 ! frac weight of HNO3 to total NO3 REAL FNO3ACC ! frac weight of NO3 acc to total NO3 REAL FRACLIQ ! fraction of water in liquid form REAL FOA1 ! First dissociation constant for FOA REAL FOAH ! Henry's Law constant for FOA REAL FOA1H ! FOAH*FOA1 REAL FOAL ! FOA conc in cloudwater (mol/liter) REAL FTST ! REAL GM ! REAL GM1 ! REAL GM1LOG ! REAL GM2 ! activity correction factor REAL GM2LOG ! REAL HA ! REAL HB ! REAL H2OW ! REAL H2O2H ! Henry's Law Constant for H2O2 REAL H2O2L ! H2O2 conc in cloudwater (mol/liter) REAL HCLH ! Henry's Law Constant for HCL REAL HCL1 ! First dissociation constant for HCL REAL HCL1H ! HCL1*HCLH REAL HCLL ! HCl conc in cloudwater (mol/liter) REAL HCO2 ! HCO2 conc in cloudwater (mol/liter) REAL HCO3 ! HCO3 conc in cloudwater (mol/liter) REAL HNO3H ! Henry's Law Constant for HNO3 REAL HNO31 ! First dissociation constant for HNO3 REAL HNO31H ! REAL HNO3L ! HNO3 conc in cloudwater (mol/liter) REAL HSO3 ! HSO3 conc in cloudwater (mol/liter) REAL HSO4 ! HSO4 concn in cloudwater (mol/liter) REAL HSO4ACC ! accumulation mode HSO4 concn in cloudwater (mol/liter) REAL HSO4COR ! coarse HSO4 concn in cloudwater (mol/liter) REAL HTST ! REAL K ! K conc in cloudwater (mol/liter) REAL KA ! initial K in cloudwater (mol/liter) REAL LGTEMP ! log of TEMP REAL M3NEW ! accumulation mode mass at time t REAL M3OLD ! accumulation mode mass at time 0 REAL MG ! REAL MGA ! inital Mg in cloudwater (mol/liter) REAL MHPH ! Henry's Law Constant for MHP REAL MHPL ! MHP conc in cloudwater (mol/liter) REAL MN ! Mn++ conc in cloudwater (mol/liter) REAL MNA ! initial Mn in cloudwater (mol/liter) REAL NA ! Na conc in cloudwater (mol/liter) REAL NAACC ! Na in cloudwater (mol/liter) REAL NAACCA ! initial Na in cloudwater (mol/liter) REAL NAAKNA ! init Aitken mode aer conc (mol/liter) REAL NACOR ! coarse Na in cloudwater (mol/liter) REAL NACORA ! init Coarse Na in cloudwater (mol/liter) REAL NH31 ! First dissociation constant for NH3 REAL NH3H ! Henry's Law Constant for NH3 REAL NH3DH20 ! REAL NH31HDH ! REAL NH3L ! NH3 conc in cloudwater (mol/liter) REAL NH4 ! NH4+ conc in cloudwater (mol/liter) REAL NH4AKNA ! init NH4 akn conc in cloudwater (mol/liter) REAL NH4ACCA ! init NH4 acc conc in cloudwater (mol/liter) REAL NITAER ! total aerosol nitrate REAL NO3 ! NO3 conc in cloudwater (mol/liter) REAL NO3ACC ! NO3 acc conc in cloudwater (mol/liter) REAL NO3ACCA ! init NO3 acc conc in cloudwater (mol/liter) REAL NO3AKNA ! init NO3 akn conc in cloudwater (mol/liter) REAL NO3CORA ! init NO3 coa conc in cloudwater (mol/liter) REAL NO3COR ! NO3 coarse conc in cloudwater (mol/liter) REAL NUMCOR ! coarse aerosol number in cloudwater (mol/liter) REAL NUMCORA ! initial coarse aerosol number in cloudwater (mol/liter) REAL O3H ! Henry's Law Constant for O3 REAL O3L ! O3 conc in cloudwater (mol/liter) REAL OH ! OH conc in cloudwater (mol/liter) REAL ORGA ! anthro SOA in cloudwater (mol/liter) REAL ORGAACCA ! init anthro ACC SOA in cloudwater (mol/liter) REAL ORGAAKNA ! init anthro AKN SOA in cloudwater (mol/liter) REAL ORGP ! primary ORGANIC aerosol in cloudwater (mol/liter) REAL ORGPACCA ! init primary ORG ACC aerosol in cloudwater (mol/liter) REAL ORGPAKNA ! init primary ORG AKN aerosol in cloudwater (mol/liter) REAL ORGB ! biogenic SOA in cloudwater (mol/liter) REAL ORGBACCA ! init biogenic ACC SOA in cloudwater (mol/liter) REAL ORGBAKNA ! init biogenic AKN SOA in cloudwater (mol/liter) REAL PAAH ! Henry's Law Constant for PAA REAL PAAL ! PAA conc in cloudwater (mol/liter) REAL PCO20 ! total CO2 partial pressure (atm) REAL PCO2F ! gas only CO2 partial pressure (atm) REAL PFOA0 ! total ORGANIC acid partial pressure (atm) REAL PFOAF ! gas only ORGANIC ACID partial press (atm) REAL PH2O20 ! total H2O2 partial pressure (atm) REAL PH2O2F ! gas only H2O2 partial pressure (atm) REAL PHCL0 ! total HCL partial pressure (atm) REAL PHCLF ! gas only HCL partial pressure (atm) REAL PHNO30 ! total HNO3 partial pressure (atm) REAL PHNO3F ! gas only HNO3 partial pressure (atm) REAL PMHP0 ! total MHP partial pressure (atm) REAL PMHPF ! gas only MHP partial pressure (atm) REAL PNH30 ! total NH3 partial pressure (atm) REAL PNH3F ! gas only NH3 partial pressure (atm) REAL PO30 ! total O3 partial pressure (atm) REAL PO3F ! gas only O3 partial pressure (atm) REAL PPAA0 ! total PAA partial pressure (atm) REAL PPAAF ! gas only PAA partial pressure (atm) REAL PRIM ! PRIMARY acc+akn aerosol in cloudwater (mol/liter) REAL PRIMCOR ! PRIMARY coarse aerosol in cloudwater (mol/liter) REAL PRIACCA ! init PRI ACC aerosol in cloudwater (mol/liter) REAL PRIAKNA ! init PRI AKN aerosol in cloudwater (mol/liter) REAL PRICORA ! init PRI COR aerosol in cloudwater (mol/liter) REAL PSO20 ! total SO2 partial pressure (atm) REAL PSO2F ! gas only SO2 partial pressure (atm) REAL RATE ! REAL RECIPA1 ! REAL RECIPA2 ! REAL RECIPAP1 ! one over pressure (/atm) REAL RH2O2 ! REAL RMHP ! REAL RPAA ! REAL RT ! gas const * temperature (liter atm/mol) REAL SCVEFF ! Scavenging efficiency (%) DATA SCVEFF / 100.0 / ! currently set to 100% SAVE SCVEFF REAL SIV ! dissolved so2 in cloudwater (mol/liter) REAL SK6 ! REAL SK6TS6 ! REAL SO21 ! First dissociation constant for SO2 REAL SO22 ! Second dissociation constant for SO2 REAL SO2H ! Henry's Law Constant for SO2 REAL SO212 ! SO21*SO22 REAL SO212H ! SO21*SO22*SO2H REAL SO21H ! SO21*SO2H REAL SO2L ! SO2 conc in cloudwater (mol/liter) REAL SO3 ! SO3= conc in cloudwater (mol/liter) REAL SO4 ! SO4= conc in cloudwater (mol/liter) REAL SO4ACC ! accumulation mode SO4= conc in cloudwater (mol/liter) REAL SO4COR ! coarse SO4= conc in cloudwater (mol/liter) REAL STION ! ionic strength REAL TAC ! REAL TEMP1 ! REAL TIMEW ! cloud chemistry clock (sec) REAL TOTOX ! REAL TOTAMM ! total ammonium REAL TOTNIT ! total nitrate (excluding coarse mode) REAL TS6 ! SO4 conc in cloudwater (mol/liter) REAL TS6AKNA ! init SO4 akn conc in cloudwater (mol/liter) REAL TS6ACC ! SO4 acc conc in cloudwater (mol/liter) REAL TS6ACCA ! init SO4 acc conc in cloudwater (mol/liter) REAL TS6COR ! coarse SO4 conc in cloudwater (mol/liter) REAL TS6CORA ! init SO4 coa conc in cloudwater (mol/liter) REAL TSIV ! REAL TST ! REAL TWASH ! washout time for clouds (sec) REAL WETFAC ! converts mol/l to mm-mol/l based on precip REAL XC1 ! (/mm) REAL XC2 ! (liter-atm/mol/mm) REAL XL ! conversion factor (liter-atm/mol) REAL ONE_OVER_XL ! 1.0 / XL REAL PRES_ATM_OVER_XL ! PRES_ATM / XL REAL XLCO2 ! REAL XLH2O2 ! REAL XLHCL ! const in calc of HCL final partial pres REAL XLHNO3 ! REAL XLMHP ! REAL XLNH3 ! REAL XLO3 ! REAL XLPAA ! REAL XLSO2 ! !...........LOCAL VARIABLES (arrays) and their descriptions: REAL WETDEP( NLIQS ) ! wet deposition array (mm mol/liter) REAL DSIVDT( 0:NUMOX ) ! rate of so2 oxid incloud (mol/liter/sec) REAL DS4 ( 0:NUMOX ) ! S(IV) oxidized over timestep DTW(0) REAL DTW ( 0:NUMOX ) ! cloud chemistry timestep (sec) REAL ONE_OVER_TEMP ! 1.0 / TEMP !...........EXTERNAL FUNCTIONS and their descriptions: ! REAL HLCONST ! EXTERNAL HLCONST !********************************************************************* ! begin body of subroutine AQCHEM ONE_OVER_TEMP = 1.0 / TEMP !...check for bad temperature, cloud air mass, or pressure IF ( TEMP .LE. 0.0 .OR. AIRM .LE. 0.0 .OR. PRES_PA .LE. 0.0 ) THEN XMSG = 'MET DATA ERROR, EXITING ROUTINE.' ! CALL M3EXIT ( PNAME, JDATE, JTIME, XMSG, XSTAT2 ) write(0,*) '' write(0,*) PNAME,' : ',XMSG write(0,*) '' write(0,*) 'TEMP :' write(0,*) TEMP write(0,*) 'PRES_PA :' write(0,*) PRES_PA write(0,*) 'TAUCLD :' write(0,*) TAUCLD write(0,*) 'PRCRATE :' write(0,*) PRCRATE write(0,*) 'WCAVG :' write(0,*) WCAVG write(0,*) 'WTAVG :' write(0,*) WTAVG write(0,*) 'AIRM :' write(0,*) AIRM write(0,*) 'ALFA0 :' write(0,*) ALFA0 write(0,*) 'ALFA2 :' write(0,*) ALFA2 write(0,*) 'ALFA3 :' write(0,*) ALFA3 write(0,*) 'GAS :' write(0,*) GAS write(0,*) 'AEROSOL :' write(0,*) AEROSOL write(0,*) 'GASWDEP :' write(0,*) GASWDEP write(0,*) 'AERWDEP :' write(0,*) AERWDEP write(0,*) 'HPWDEP :' write(0,*) HPWDEP write(0,*) '' return END IF !...initialize counters and compute several conversion factors ICNTAQ = 0 ITERAT = 0 RT = ( MOLVOL / STDTEMP ) * TEMP ! R * T (liter atm / mol) PRES_ATM = PRES_PA / STDATMPA ! pressure (atm) CTHK1 = AIRM * RT / ( PRES_ATM * 1000.0 ) ! cloud thickness (m) XL = WCAVG * RT / H2ODENS ! conversion factor (l-atm/mol) ONE_OVER_XL = 1.0 / XL PRES_ATM_OVER_XL = PRES_ATM / XL TST = 0.999 GM = SCVEFF / 100.0 ACT1 = 1.0 ACT2 = 1.0 GM2 = 1.0 TIMEW = 0.0 RECIPAP1 = 1.0 / PRES_ATM XC1 = 1.0 / ( WCAVG * CTHK1 ) XC2 = RT / ( 1000.0 * CTHK1 ) FRACLIQ = WCAVG / WTAVG TWASH = WTAVG * 1000.0 * CTHK1 * 3600.0 & / ( H2ODENS * MAX( 1.0E-20, PRCRATE ) ) !...set equilibrium constants as a function of temperature !... Henry's law constants SO2H = HLCONST( 'SO2 ', TEMP, .FALSE., 0.0 ) CO2H = HLCONST( 'CO2 ', TEMP, .FALSE., 0.0 ) NH3H = HLCONST( 'NH3 ', TEMP, .FALSE., 0.0 ) H2O2H = HLCONST( 'H2O2 ', TEMP, .FALSE., 0.0 ) O3H = HLCONST( 'O3 ', TEMP, .FALSE., 0.0 ) HCLH = HLCONST( 'HCL ', TEMP, .FALSE., 0.0 ) HNO3H = HLCONST( 'HNO3 ', TEMP, .FALSE., 0.0 ) MHPH = HLCONST( 'METHYLHYDROPEROX', TEMP, .FALSE., 0.0 ) PAAH = HLCONST( 'PEROXYACETIC_ACI', TEMP, .FALSE., 0.0 ) FOAH = HLCONST( 'FORMIC_ACID ', TEMP, .FALSE., 0.0 ) TEMP1 = ONE_OVER_TEMP - 1.0 / 298.0 !...dissociation constants FOA1 = 1.80E-04 * EXP( -2.00E+01 * TEMP1 ) ! Martell and Smith (1977) SK6 = 1.02E-02 * EXP( 2.72E+03 * TEMP1 ) ! Smith and Martell (1976) SO21 = 1.30E-02 * EXP( 1.96E+03 * TEMP1 ) ! Smith and Martell (1976) SO22 = 6.60E-08 * EXP( 1.50E+03 * TEMP1 ) ! Smith and Martell (1976) CO21 = 4.30E-07 * EXP( -1.00E+03 * TEMP1 ) ! Smith and Martell (1976) CO22 = 4.68E-11 * EXP( -1.76E+03 * TEMP1 ) ! Smith and Martell (1976) H2OW = 1.00E-14 * EXP( -6.71E+03 * TEMP1 ) ! Smith and Martell (1976) NH31 = 1.70E-05 * EXP( -4.50E+02 * TEMP1 ) ! Smith and Martell (1976) HCL1 = 1.74E+06 * EXP( 6.90E+03 * TEMP1 ) ! Marsh and McElroy (1985) HNO31 = 1.54E+01 * EXP( 8.70E+03 * TEMP1 ) ! Schwartz (1984) !...Kinetic oxidation rates !... From Chamedies (1982) ! RH2O2 = 8.0E+04 * EXP( -3650.0 * TEMP1 ) !KW based on CMAQv5.0 From Jacobson (1997) RH2O2 = 7.45E+07 * EXP( -15.96E0 * ( ( 298.0E0 / TEMP ) - 1.0E0 ) ) !...From Kok ! RMHP = 1.75E+07 * EXP( -3801.0 * TEMP1 ) ! RPAA = 3.64E+07 * EXP( -3994.0 * TEMP1 ) !KW based on CMAQv5.0 From Jacobson (1997) RMHP = 1.90E+07 * EXP( -12.75E0 * ( ( 298.0E0 / TEMP ) - 1.0E0 ) ) RPAA = 3.67E+07 * EXP( -13.42E0 * ( ( 298.0E0 / TEMP ) - 1.0E0 ) ) !...make initializations DO LIQ = 1, NLIQS WETDEP( LIQ ) = 0.0 END DO DO IOX = 0, NUMOX DSIVDT( IOX ) = 0.0 DTW ( IOX ) = 0.0 DS4 ( IOX ) = 0.0 END DO !...compute the initial accumulation aerosol 3rd moment !... secondary organic aerosol and water are not included M3OLD = ( AEROSOL( LSO4ACC ) * SGRAERMW( LSO4ACC ) / 1.8e6 & + AEROSOL( LNH4ACC ) * SGRAERMW( LNH4ACC ) / 1.8e6 & + AEROSOL( LNO3ACC ) * SGRAERMW( LNO3ACC ) / 1.8e6 & + AEROSOL( LORGPACC ) * SGRAERMW( LORGPACC ) / 2.0e6 & + AEROSOL( LECACC ) * SGRAERMW( LECACC ) / 2.2e6 & + AEROSOL( LPRIACC ) * SGRAERMW( LPRIACC ) / 2.2e6 & + AEROSOL( LNAACC ) * SGRAERMW( LNAACC ) / 2.2e6 & + AEROSOL( LCLACC ) * SGRAERMW( LCLACC ) / 2.2e6 ) !cc & * 6.0 / PI ! cancels out in division at end of subroutine !...compute fractional weights for several species TOTNIT = GAS( LHNO3 ) + AEROSOL( LNO3ACC ) IF ( TOTNIT .GT. 0.0 ) THEN FHNO3 = GAS( LHNO3 ) / TOTNIT FNO3ACC = AEROSOL( LNO3ACC ) / TOTNIT ELSE FHNO3 = 1.0 FNO3ACC = 0.0 END IF TOTAMM = GAS( LNH3 ) + AEROSOL( LNH4ACC ) IF ( TOTAMM .GT. 0.0 ) THEN FNH3 = GAS( LNH3 ) / TOTAMM FNH4ACC = AEROSOL( LNH4ACC ) / TOTAMM ELSE FNH3 = 1.0 FNH4ACC = 0.0 END IF !...initial concentration from accumulation-mode aerosol loading (mol/liter) !... an assumption is made that all of the accumulation-mode !... aerosol mass in incorporated into the cloud droplets TS6ACCA = ( AEROSOL( LSO4ACC ) & + GAS ( LH2SO4 ) ) * PRES_ATM_OVER_XL NO3ACCA = AEROSOL( LNO3ACC ) * PRES_ATM_OVER_XL NH4ACCA = AEROSOL( LNH4ACC ) * PRES_ATM_OVER_XL ORGAACCA = AEROSOL( LORGAACC ) * PRES_ATM_OVER_XL ORGPACCA = AEROSOL( LORGPACC ) * PRES_ATM_OVER_XL ORGBACCA = AEROSOL( LORGBACC ) * PRES_ATM_OVER_XL ECACCA = AEROSOL( LECACC ) * PRES_ATM_OVER_XL PRIACCA = AEROSOL( LPRIACC ) * PRES_ATM_OVER_XL NAACCA = AEROSOL( LNAACC ) * PRES_ATM_OVER_XL CLACCA = AEROSOL( LCLACC ) * PRES_ATM_OVER_XL !...initial concentration from coarse-mode aerosol loading (mol/liter) !... an assumption is made that all of the coarse-mode !... aerosol mass in incorporated into the cloud droplets TS6CORA = AEROSOL( LSO4COR ) * PRES_ATM_OVER_XL NO3CORA = AEROSOL( LNO3COR ) * PRES_ATM_OVER_XL IF ( AE_VRSN .EQ. 'AE3' ) THEN CLCORA = AEROSOL( LNACL ) * PRES_ATM_OVER_XL NACORA = AEROSOL( LNACL ) * PRES_ATM_OVER_XL ELSE CLCORA = AEROSOL( LCLCOR ) * PRES_ATM_OVER_XL NACORA = AEROSOL( LNACOR ) * PRES_ATM_OVER_XL END IF KA = AEROSOL( LK ) * PRES_ATM_OVER_XL CAA = AEROSOL( LCACO3 ) * PRES_ATM_OVER_XL MGA = AEROSOL( LMGCO3 ) * PRES_ATM_OVER_XL FEA = AEROSOL( LA3FE ) * PRES_ATM_OVER_XL MNA = AEROSOL( LB2MN ) * PRES_ATM_OVER_XL CO3A = ( AEROSOL( LCACO3 ) & + AEROSOL( LMGCO3 ) ) * PRES_ATM_OVER_XL PRICORA = AEROSOL( LPRICOR ) * PRES_ATM_OVER_XL NUMCORA = AEROSOL( LNUMCOR ) * PRES_ATM_OVER_XL !...set constant factors that will be used in later multiplications (moles/atm) XLH2O2 = H2O2H * XL XLO3 = O3H * XL XLMHP = MHPH * XL XLPAA = PAAH * XL XLSO2 = SO2H * XL XLNH3 = NH3H * XL XLHCL = HCLH * XL XLHNO3 = HNO3H * XL XLCO2 = CO2H * XL SO212 = SO21 * SO22 SO21H = SO21 * SO2H SO212H = SO212 * SO2H CO212 = CO21 * CO22 CO21H = CO21 * CO2H CO212H = CO22 * CO21H NH3DH20 = NH31 / H2OW NH31HDH = NH3H * NH3DH20 FOA1H = FOA1 * FOAH HCL1H = HCL1 * HCLH HNO31H = HNO31 * HNO3H !...If kinetic calculations are made, return to this point I20C = 0 20 CONTINUE I20C = I20C + 1 IF ( I20C .GE. 10000 ) THEN XMSG = 'EXCESSIVE LOOPING AT I20C, EXITING ROUTINE.' ! CALL M3EXIT ( PNAME, JDATE, JTIME, XMSG, XSTAT2 ) write(0,*) '' write(0,*) PNAME,' : ',XMSG write(0,*) '' write(0,*) 'TEMP :' write(0,*) TEMP write(0,*) 'PRES_PA :' write(0,*) PRES_PA write(0,*) 'TAUCLD :' write(0,*) TAUCLD write(0,*) 'PRCRATE :' write(0,*) PRCRATE write(0,*) 'WCAVG :' write(0,*) WCAVG write(0,*) 'WTAVG :' write(0,*) WTAVG write(0,*) 'AIRM :' write(0,*) AIRM write(0,*) 'ALFA0 :' write(0,*) ALFA0 write(0,*) 'ALFA2 :' write(0,*) ALFA2 write(0,*) 'ALFA3 :' write(0,*) ALFA3 write(0,*) 'GAS :' write(0,*) GAS write(0,*) 'AEROSOL :' write(0,*) AEROSOL write(0,*) 'GASWDEP :' write(0,*) GASWDEP write(0,*) 'AERWDEP :' write(0,*) AERWDEP write(0,*) 'HPWDEP :' write(0,*) HPWDEP write(0,*) '' return END IF !...set aitken-mode aerosol loading (mol/liter) NO3AKNA = AEROSOL( LNO3AKN ) * PRES_ATM_OVER_XL & * ( 1.0 - EXP( -ALFA3 * TIMEW ) ) NH4AKNA = AEROSOL( LNH4AKN ) * PRES_ATM_OVER_XL & * ( 1.0 - EXP( -ALFA3 * TIMEW ) ) TS6AKNA = AEROSOL( LSO4AKN ) * PRES_ATM_OVER_XL & * ( 1.0 - EXP( -ALFA3 * TIMEW ) ) ORGAAKNA = AEROSOL( LORGAAKN ) * PRES_ATM_OVER_XL & * ( 1.0 - EXP( -ALFA3 * TIMEW ) ) ORGPAKNA = AEROSOL( LORGPAKN ) * PRES_ATM_OVER_XL & * ( 1.0 - EXP( -ALFA3 * TIMEW ) ) ORGBAKNA = AEROSOL( LORGBAKN ) * PRES_ATM_OVER_XL & * ( 1.0 - EXP( -ALFA3 * TIMEW ) ) ECAKNA = AEROSOL( LECAKN ) * PRES_ATM_OVER_XL & * ( 1.0 - EXP( -ALFA3 * TIMEW ) ) PRIAKNA = AEROSOL( LPRIAKN ) * PRES_ATM_OVER_XL & * ( 1.0 - EXP( -ALFA3 * TIMEW ) ) NAAKNA = AEROSOL( LNAAKN ) * PRES_ATM_OVER_XL & * ( 1.0 - EXP( -ALFA3 * TIMEW ) ) CLAKNA = AEROSOL( LCLAKN ) * PRES_ATM_OVER_XL & * ( 1.0 - EXP( -ALFA3 * TIMEW ) ) !...Initial gas phase partial pressures (atm) !... = initial partial pressure - amount deposited partial pressure PSO20 = GAS( LSO2 ) * PRES_ATM & + DS4( 0 ) * XL & - ( WETDEP( LSO3L ) + WETDEP( LHSO3L ) + WETDEP( LSO2L ) ) * XC2 PNH30 = GAS( LNH3 ) * PRES_ATM & + ( NH4ACCA + NH4AKNA ) * XL & - ( WETDEP( LNH4L ) + WETDEP( LNH3L ) ) * XC2 PHNO30 = ( GAS( LHNO3 ) + 2.0 * GAS( LN2O5 ) ) * PRES_ATM & + ( NO3ACCA + NO3CORA + NO3AKNA ) * XL & - ( WETDEP( LNO3ACCL ) + WETDEP( LHNO3L ) + WETDEP( LNO3CORL ) ) * XC2 PHCL0 = GAS( LHCL ) * PRES_ATM & + ( CLACCA + CLCORA + CLAKNA ) * XL & ! new for sea salt - ( WETDEP( LCLACCL ) + WETDEP( LHCLL ) + WETDEP( LCLCORL ) ) * XC2 PH2O20 = GAS( LH2O2 ) * PRES_ATM - WETDEP( LH2O2L ) * XC2 PO30 = GAS( LO3 ) * PRES_ATM - WETDEP( LO3L ) * XC2 PFOA0 = GAS( LFOA ) * PRES_ATM & - ( WETDEP( LFOAL ) + WETDEP( LHCO2L ) ) * XC2 PMHP0 = GAS( LMHP ) * PRES_ATM - WETDEP( LMHPL ) * XC2 PPAA0 = GAS( LPAA ) * PRES_ATM - WETDEP( LPAAL ) * XC2 PCO20 = GAS( LCO2 ) * PRES_ATM & + CO3A * XL & - ( WETDEP( LCO3L ) + WETDEP( LHCO3L ) + WETDEP( LCO2L ) ) * XC2 !...don't allow gas concentrations to go below zero PSO20 = MAX( PSO20, 0.0 ) PNH30 = MAX( PNH30, 0.0 ) PH2O20 = MAX( PH2O20, 0.0 ) PO30 = MAX( PO30, 0.0 ) PFOA0 = MAX( PFOA0, 0.0 ) PMHP0 = MAX( PMHP0, 0.0 ) PPAA0 = MAX( PPAA0, 0.0 ) PCO20 = MAX( PCO20, 0.0 ) PHCL0 = MAX( PHCL0, 0.0 ) PHNO30 = MAX( PHNO30, 0.0 ) !...Molar concentrations of soluble aerosols !... = Initial amount - amount deposited (mol/liter) TS6COR = MAX( TS6CORA - WETDEP( LTS6CORL ) * XC1, 0.0 ) NO3COR = MAX( NO3CORA - WETDEP( LNO3CORL ) * XC1, 0.0 ) NACOR = MAX( NACORA - WETDEP( LNACORL ) * XC1, 0.0 ) CLCOR = MAX( CLCORA - WETDEP( LCLCORL ) * XC1, 0.0 ) TS6 = TS6ACCA + TS6AKNA + TS6COR & - ( WETDEP( LSO4ACCL ) + WETDEP( LHSO4ACCL ) ) * XC1 & - DS4( 0 ) NA = NAACCA + NAAKNA + NACOR & - WETDEP( LNAACCL ) * XC1 CA = CAA - WETDEP( LCAL ) * XC1 MG = MGA - WETDEP( LMGL ) * XC1 K = KA - WETDEP( LKL ) * XC1 FE = FEA - WETDEP( LFEL ) * XC1 MN = MNA - WETDEP( LMNL ) * XC1 ORGA = ORGAACCA + ORGAAKNA - WETDEP( LORGAL ) * XC1 ORGP = ORGPACCA + ORGPAKNA - WETDEP( LORGPL ) * XC1 ORGB = ORGBACCA + ORGBAKNA - WETDEP( LORGBL ) * XC1 EC = ECACCA + ECAKNA - WETDEP( LECL ) * XC1 PRIM = PRIACCA + PRIAKNA - WETDEP( LPRIML ) * XC1 PRIMCOR = PRICORA - WETDEP( LPRIMCORL ) * XC1 NUMCOR = NUMCORA - WETDEP( LNUMCORL ) * XC1 A = 3.0 * FE B = 2.0 * MN !...don't allow aerosol concentrations to go below zero TS6 = MAX( TS6, 0.0 ) NA = MAX( NA, 0.0 ) CA = MAX( CA, 0.0 ) MG = MAX( MG, 0.0 ) K = MAX( K, 0.0 ) FE = MAX( FE, 0.0 ) MN = MAX( MN, 0.0 ) ORGA = MAX( ORGA, 0.0 ) ORGP = MAX( ORGP, 0.0 ) ORGB = MAX( ORGB, 0.0 ) EC = MAX( EC, 0.0 ) PRIM = MAX( PRIM, 0.0 ) PRIMCOR = MAX( PRIMCOR, 0.0 ) NUMCOR = MAX( NUMCOR, 0.0 ) A = MAX( A, 0.0 ) B = MAX( B, 0.0 ) SK6TS6 = SK6 * TS6 !...find solution of the equation using a method of reiterative !... bisections Make initial guesses for pH: between .01 to 10. HA = 0.01 HB = 10.0 I7777C = 0 7777 CONTINUE I7777C = I7777C + 1 IF ( I7777C .GE. 10000 ) THEN XMSG = 'EXCESSIVE LOOPING AT I7777C, EXITING ROUTINE.' ! CALL M3EXIT ( PNAME, JDATE, JTIME, XMSG, XSTAT2 ) write(0,*) '' write(0,*) PNAME,' : ',XMSG write(0,*) '' write(0,*) 'TEMP :' write(0,*) TEMP write(0,*) 'PRES_PA :' write(0,*) PRES_PA write(0,*) 'TAUCLD :' write(0,*) TAUCLD write(0,*) 'PRCRATE :' write(0,*) PRCRATE write(0,*) 'WCAVG :' write(0,*) WCAVG write(0,*) 'WTAVG :' write(0,*) WTAVG write(0,*) 'AIRM :' write(0,*) AIRM write(0,*) 'ALFA0 :' write(0,*) ALFA0 write(0,*) 'ALFA2 :' write(0,*) ALFA2 write(0,*) 'ALFA3 :' write(0,*) ALFA3 write(0,*) 'GAS :' write(0,*) GAS write(0,*) 'AEROSOL :' write(0,*) AEROSOL write(0,*) 'GASWDEP :' write(0,*) GASWDEP write(0,*) 'AERWDEP :' write(0,*) AERWDEP write(0,*) 'HPWDEP :' write(0,*) HPWDEP write(0,*) '' return END IF HA = MAX( HA - 0.8, 0.1 ) HB = MIN( HB + 0.8, 9.9 ) AE = 10.0**( -HA ) RECIPA1 = 1.0 / ( AE * ACT1 ) RECIPA2 = 1.0 / ( AE * AE * ACT2 ) !...calculate final gas phase partial pressure of SO2, NH3, HNO3 !... HCOOH, and CO2 (atm) PSO2F = PSO20 / ( 1.0 + XLSO2 * ( 1.0 + SO21 * RECIPA1 & + SO212 * RECIPA2 ) ) PNH3F = PNH30 / ( 1.0 + XLNH3 * ( 1.0 + NH3DH20 * AE ) ) PHCLF = PHCL0 / ( 1.0 + XLHCL * ( 1.0 + HCL1 * RECIPA1 ) ) PFOAF = PFOA0 / ( 1.0 + XL * ( FOAH + FOA1H * RECIPA1 ) ) PHNO3F = PHNO30 / ( 1.0 + XLHNO3 * ( 1.0 + HNO31 * RECIPA1 ) ) PCO2F = PCO20 / ( 1.0 + XLCO2 * ( 1.0 + CO21 * RECIPA1 & + CO212 * RECIPA2 ) ) !...calculate liquid phase concentrations (moles/liter) SO4 = SK6TS6 / ( AE * GM2 + SK6 ) HSO4 = TS6 - SO4 SO3 = SO212H * PSO2F * RECIPA2 HSO3 = SO21H * PSO2F * RECIPA1 CO3 = CO212H * PCO2F * RECIPA2 HCO3 = CO21H * PCO2F * RECIPA1 OH = H2OW * RECIPA1 NH4 = NH31HDH * PNH3F * AE HCO2 = FOA1H * PFOAF * RECIPA1 NO3 = HNO31H * PHNO3F * RECIPA1 CL = HCL1H * PHCLF * RECIPA1 ! new for sea salt !...compute functional value FA = AE + NH4 + NA + 2.0 * ( CA + MG - CO3 - SO3 - SO4 ) & - OH - HCO3 - HSO3 - NO3 - HSO4 - HCO2 - CL !...Start iteration and bisection ****************<<<<<<< I30C = 0 30 CONTINUE I30C = I30C + 1 IF ( I30C .GE. 10000 ) THEN XMSG = 'EXCESSIVE LOOPING AT I30C, EXITING ROUTINE.' ! CALL M3EXIT ( PNAME, JDATE, JTIME, XMSG, XSTAT2 ) write(0,*) '' write(0,*) PNAME,' : ',XMSG write(0,*) '' write(0,*) 'TEMP :' write(0,*) TEMP write(0,*) 'PRES_PA :' write(0,*) PRES_PA write(0,*) 'TAUCLD :' write(0,*) TAUCLD write(0,*) 'PRCRATE :' write(0,*) PRCRATE write(0,*) 'WCAVG :' write(0,*) WCAVG write(0,*) 'WTAVG :' write(0,*) WTAVG write(0,*) 'AIRM :' write(0,*) AIRM write(0,*) 'ALFA0 :' write(0,*) ALFA0 write(0,*) 'ALFA2 :' write(0,*) ALFA2 write(0,*) 'ALFA3 :' write(0,*) ALFA3 write(0,*) 'GAS :' write(0,*) GAS write(0,*) 'AEROSOL :' write(0,*) AEROSOL write(0,*) 'GASWDEP :' write(0,*) GASWDEP write(0,*) 'AERWDEP :' write(0,*) AERWDEP write(0,*) 'HPWDEP :' write(0,*) HPWDEP write(0,*) '' return END IF BB = ( HA + HB ) / 2.0 AE = 10.0**( -BB ) ICNTAQ = ICNTAQ + 1 IF ( ICNTAQ .GE. 60000 ) THEN XMSG = 'Maximum AQCHEM total iterations exceeded, EXITING ROUTINE.' ! CALL M3EXIT ( PNAME, JDATE, JTIME, XMSG, XSTAT2 ) write(0,*) '' write(0,*) PNAME,' : ',XMSG write(0,*) '' write(0,*) 'TEMP :' write(0,*) TEMP write(0,*) 'PRES_PA :' write(0,*) PRES_PA write(0,*) 'TAUCLD :' write(0,*) TAUCLD write(0,*) 'PRCRATE :' write(0,*) PRCRATE write(0,*) 'WCAVG :' write(0,*) WCAVG write(0,*) 'WTAVG :' write(0,*) WTAVG write(0,*) 'AIRM :' write(0,*) AIRM write(0,*) 'ALFA0 :' write(0,*) ALFA0 write(0,*) 'ALFA2 :' write(0,*) ALFA2 write(0,*) 'ALFA3 :' write(0,*) ALFA3 write(0,*) 'GAS :' write(0,*) GAS write(0,*) 'AEROSOL :' write(0,*) AEROSOL write(0,*) 'GASWDEP :' write(0,*) GASWDEP write(0,*) 'AERWDEP :' write(0,*) AERWDEP write(0,*) 'HPWDEP :' write(0,*) HPWDEP write(0,*) '' return END IF RECIPA1 = 1.0 / ( AE * ACT1 ) RECIPA2 = 1.0 / ( AE * AE * ACT2 ) !...calculate final gas phase partial pressure of SO2, NH3, HCL, HNO3 !... HCOOH, and CO2 (atm) PSO2F = PSO20 / ( 1.0 + XLSO2 & * ( 1.0 + SO21 * RECIPA1 + SO212 * RECIPA2 ) ) PNH3F = PNH30 / ( 1.0 + XLNH3 * ( 1.0 + NH3DH20 * AE ) ) PHCLF = PHCL0 / ( 1.0 + XLHCL * ( 1.0 + HCL1 * RECIPA1 ) ) PHNO3F = PHNO30 / ( 1.0 + XLHNO3 * ( 1.0 + HNO31 * RECIPA1 ) ) PFOAF = PFOA0 / ( 1.0 + XL * ( FOAH + FOA1H * RECIPA1 ) ) PCO2F = PCO20 / ( 1.0 + XLCO2 * ( 1.0 + CO21 * RECIPA1 & + CO212 * RECIPA2 ) ) !...calculate liquid phase concentrations (moles/liter) SO4 = SK6TS6 / ( AE * GM2 + SK6 ) HSO4 = TS6 - SO4 SO3 = SO212H * PSO2F * RECIPA2 HSO3 = SO21H * PSO2F * RECIPA1 CO3 = CO212H * PCO2F * RECIPA2 HCO3 = CO21H * PCO2F * RECIPA1 OH = H2OW * RECIPA1 NH4 = NH31HDH * PNH3F * AE HCO2 = FOA1H * PFOAF * RECIPA1 NO3 = HNO31H * PHNO3F * RECIPA1 CL = HCL1H * PHCLF * RECIPA1 ! new for sea salt !...compute functional value FB = AE + NH4 + NA + 2.0 * ( CA + MG - CO3 - SO3 - SO4 ) & - OH - HCO3 - HSO3 - NO3 - HSO4 - HCO2 - CL !...Calculate and check the sign of the product of the two functional values FTST = FA * FB IF ( FTST .LE. 0.0 ) THEN HB = BB ELSE HA = BB FA = FB END IF !...Check convergence of solutions HTST = HA / HB IF ( HTST .LE. TST ) GO TO 30 !...end of zero-finding routine ****************<<<<<<<<<<<< !...compute Ionic strength and activity coefficient by the Davies equation STION = 0.5 * (AE + NH4 + OH + HCO3 + HSO3 & + 4.0 * (SO4 + CO3 + SO3 + CA + MG + MN) & + NO3 + HSO4 + 9.0 * FE + NA + K + CL + A + B + HCO2) GM1LOG = -0.509 * ( SQRT( STION ) & / ( 1.0 + SQRT( STION ) ) - 0.2 * STION ) GM2LOG = GM1LOG * 4.0 GM1 = 10.0**GM1LOG GM2 = MAX( 10.0**GM2LOG, 1.0E-30 ) ACTB = ACT1 ACT1 = MAX( GM1 * GM1, 1.0E-30 ) ACT2 = MAX( GM1 * GM1 * GM2, 1.0E-30 ) !...check for convergence and possibly go to 7777, to recompute !... Gas and liquid phase concentrations TAC = ABS( ACTB - ACT1 ) / ACTB IF ( TAC .GE. 1.0E-2 ) GO TO 7777 !...return an error if the pH is not in range !cc IF ( ( HA .LT. 0.02 ) .OR. ( HA .GT. 9.49 ) ) THEN IF ( ( HA .LT. 0.1 ) .OR. ( HA .GT. 9.9 ) ) THEN print *, ha XMSG = 'PH VALUE OUT OF RANGE, EXITING ROUTINE.' ! CALL M3EXIT ( PNAME, JDATE, JTIME, XMSG, XSTAT2 ) write(0,*) '' write(0,*) PNAME,' : ',XMSG write(0,*) '' write(0,*) 'TEMP :' write(0,*) TEMP write(0,*) 'PRES_PA :' write(0,*) PRES_PA write(0,*) 'TAUCLD :' write(0,*) TAUCLD write(0,*) 'PRCRATE :' write(0,*) PRCRATE write(0,*) 'WCAVG :' write(0,*) WCAVG write(0,*) 'WTAVG :' write(0,*) WTAVG write(0,*) 'AIRM :' write(0,*) AIRM write(0,*) 'ALFA0 :' write(0,*) ALFA0 write(0,*) 'ALFA2 :' write(0,*) ALFA2 write(0,*) 'ALFA3 :' write(0,*) ALFA3 write(0,*) 'GAS :' write(0,*) GAS write(0,*) 'AEROSOL :' write(0,*) AEROSOL write(0,*) 'GASWDEP :' write(0,*) GASWDEP write(0,*) 'AERWDEP :' write(0,*) AERWDEP write(0,*) 'HPWDEP :' write(0,*) HPWDEP write(0,*) '' return END IF !...Make those concentration calculations which can be made outside !... of the function. SO2L = SO2H * PSO2F AC = 10.0**( -BB ) SIV = SO3 + HSO3 + SO2L !...Calculate final gas phase concentrations of oxidants (atm) PH2O2F = ( PH2O20 + XL * DS4( 1 ) ) / ( 1.0 + XLH2O2 ) PO3F = ( PO30 + XL * DS4( 2 ) ) / ( 1.0 + XLO3 ) PMHPF = ( PMHP0 + XL * DS4( 4 ) ) / ( 1.0 + XLMHP ) PPAAF = ( PPAA0 + XL * DS4( 5 ) ) / ( 1.0 + XLPAA ) PH2O2F = MAX( PH2O2F, 0.0 ) PO3F = MAX( PO3F, 0.0 ) PMHPF = MAX( PMHPF, 0.0 ) PPAAF = MAX( PPAAF, 0.0 ) !...Calculate liquid phase concentrations of oxidants (moles/liter) H2O2L = PH2O2F * H2O2H O3L = PO3F * O3H MHPL = PMHPF * MHPH PAAL = PPAAF * PAAH FOAL = PFOAF * FOAH NH3L = PNH3F * NH3H CO2L = PCO2F * CO2H HCLL = PHCLF * HCLH HNO3L = PHNO3F * HNO3H !...compute modal concentrations SO4COR = SK6 * TS6COR / ( AE * GM2 + SK6 ) HSO4COR = MAX( TS6COR - SO4COR, 0.0 ) TS6ACC = MAX( TS6 - TS6COR, 0.0 ) SO4ACC = MAX( SO4 - SO4COR, 0.0 ) HSO4ACC = MAX( HSO4 - HSO4COR, 0.0 ) NO3ACC = MAX( NO3 - NO3COR, 0.0 ) NAACC = MAX( NA - NACOR, 0.0 ) CLACC = MAX( CL - CLCOR, 0.0 ) !...load the liquid concentration array with current values LIQUID( LACL ) = AC LIQUID( LNH4L ) = NH4 LIQUID( LCAL ) = CA LIQUID( LNAACCL ) = NAACC LIQUID( LOHL ) = OH LIQUID( LSO4ACCL ) = SO4ACC LIQUID( LHSO4ACCL ) = HSO4ACC LIQUID( LSO3L ) = SO3 LIQUID( LHSO3L ) = HSO3 LIQUID( LSO2L ) = SO2L LIQUID( LCO3L ) = CO3 LIQUID( LHCO3L ) = HCO3 LIQUID( LCO2L ) = CO2L LIQUID( LNO3ACCL ) = NO3ACC LIQUID( LNH3L ) = NH3L LIQUID( LCLACCL ) = CLACC LIQUID( LH2O2L ) = H2O2L LIQUID( LO3L ) = O3L LIQUID( LFEL ) = FE LIQUID( LMNL ) = MN LIQUID( LAL ) = A LIQUID( LFOAL ) = FOAL LIQUID( LHCO2L ) = HCO2 LIQUID( LMHPL ) = MHPL LIQUID( LPAAL ) = PAAL LIQUID( LHCLL ) = HCLL LIQUID( LORGAL ) = ORGA LIQUID( LPRIML ) = PRIM LIQUID( LMGL ) = MG LIQUID( LKL ) = K LIQUID( LBL ) = B LIQUID( LHNO3L ) = HNO3L LIQUID( LPRIMCORL ) = PRIMCOR LIQUID( LNUMCORL ) = NUMCOR LIQUID( LTS6CORL ) = TS6COR LIQUID( LNACORL ) = NACOR LIQUID( LCLCORL ) = CLCOR LIQUID( LNO3CORL ) = NO3COR LIQUID( LORGPL ) = ORGP LIQUID( LORGBL ) = ORGB LIQUID( LECL ) = EC !...if the maximum cloud lifetime has not been reached, then compute !... the next timestep. IF ( TIMEW .LT. TAUCLD ) THEN !...make kinetics calculations !... note: DS4(i) and DSIV(I) are negative numbers! DTRMV = 300.0 IF ( ( CTHK1 .GT. 1.0E-10 ) .AND. ( PRCRATE .GT. 1.0E-10 ) ) & DTRMV = 3.6 * WTAVG * 1000.0 * CTHK1 / PRCRATE ! <<HSO3+H : Smith and Martell (1976) DATA B( LHSO3 ), D( LHSO3 ) / 6.60E-08, 1.50E+03 / ! HSO3<=>SO3+H : Smith and Martell (1976) DATA B( LHNO2 ), D( LHNO2 ) / 5.10E-04, -1.26E+03 / ! HNO2(aq)<=>NO2+H : Schwartz and White (1981) DATA B( LHNO3 ), D( LHNO3 ) / 1.54E+01, 8.70E+03 / ! HNO3(aq)<=>NO3+H : Schwartz (1984) DATA B( LCO2 ), D( LCO2 ) / 4.30E-07, -1.00E+03 / ! CO2*H2O<=>HCO3+H : Smith and Martell (1976) DATA B( LHCO3 ), D( LHCO3 ) / 4.68E-11, -1.76E+03 / ! HCO3<=>CO3+H : Smith and Martell (1976) DATA B( LH2O2 ), D( LH2O2 ) / 2.20E-12, -3.73E+03 / ! H2O2(aq)<=>HO2+H : Smith and Martell (1976) DATA B( LHCHO ), D( LHCHO ) / 2.53E+03, 4.02E+03 / ! HCHO(aq)<=>H2C(OH)2 : Le Hanaf (1968) DATA B( LHCOOH ), D( LHCOOH ) / 1.80E-04, -2.00E+01 / ! HCOOH(aq)<=>HCOO+H : Martell and Smith (1977) DATA B( LHO2 ), D( LHO2 ) / 3.50E-05, 0.00E+00 / ! HO2(aq)<=>H+O2 : Perrin (1982) DATA B( LNH4OH ), D( LNH4OH ) / 1.70E-05, -4.50E+02 / ! NH4*OH<=>NH4+OH : Smith and Martell (1976) DATA B( LH2O ), D( LH2O ) / 1.00E-14, -6.71E+03 / ! H2O<=>H+OH : Smith and Martell (1976) DATA B( LATRA ), D( LATRA ) / 2.09E-02, 0.00E+00 / ! C8H14ClN5<=>C8H13ClN5+H : Weber (1970) DATA B( LCL2 ), D( LCL2 ) / 5.01E-04, 0.00E+00 / ! CL2*H2O <=> HOCL + H + CL : LIN AND PEHKONEN, JGR, 103, D21, 28093-28102, NOVEMBER 20, 1998. ALSO SEE NOTE BELOW DATA B( LHOCL ), D( LHOCL ) / 3.16E-08, 0.00E+00 / ! HOCL <=>H + OCL : LIN AND PEHKONEN, JGR, 103, D21, 28093-28102, NOVEMBER 20, 1998 DATA B( LHCL ), D( LHCL ) / 1.74E+06, 6.90E+03 / ! HCL <=> H + CL : Marsh and McElroy (1985) DATA B( LHYDRAZINE), D( LHYDRAZINE) / 1.11E-08, 0.00E+00 / ! HYDRAZINE <=> HYDRAZINE+ + OH- : Moliner and Street (1989) !------------------------------------------------------------------------------- ! Note for dissociation constant for equation 14: CL2*H2O <=> HOCL + H + CL ! Need aqueous [CL-] concentration to calculate effective henry's law coefficient ! Used a value of 2.0 mM following Lin and Pehkonen, JGR, 103, D21, 28093-28102, November 20, 1998 !------------------------------------------------------------------------------- !...........EXTERNAL FUNCTIONS and their descriptions: ! INTEGER, EXTERNAL :: INDEX1 ! array position for string matching ! INTEGER, EXTERNAL :: TRIMLEN ! string length, excl. trailing blanks !----------------------------------------------------------------------- ! begin body of subroutine HLCONST SPC = INDEX1( CNAME, MXSPCS, SUBNAME ) !...error if species not found in table IF ( SPC .LE. 0 ) THEN XMSG = CNAME( 1:TRIMLEN( CNAME ) ) // ' not found in Henry''s Law Constant table, aborting.' ! CALL M3EXIT ( PNAME, 0, 0, XMSG, XSTAT2 ) write(0,*) '' write(0,*) PNAME,' : ',XMSG stop END IF !...compute the Henry's Law Constant TFAC = ( 298.0 - TEMP) / ( 298.0 * TEMP ) KH = A( SPC ) * EXP( E( SPC ) * TFAC ) HLCONST = KH !...compute the effective Henry's law constants IF ( EFFECTIVE ) THEN IF ( HPLUS .LE. 0.0 ) THEN XMSG = 'Negative or Zero [H+] concentration specified, aborting.' ! CALL M3EXIT ( PNAME, JDATE, JTIME, XMSG, XSTAT2 ) write(0,*) '' write(0,*) PNAME,' : ',XMSG stop END IF HPLUSI = 1.0 / HPLUS HPLUS2I = HPLUSI * HPLUSI !...assign a value for clminus. use 2.0 mM based on Lin and Pehkonene, 1998, JGR CLMINUS = 2.0E-03 ! chlorine ion conc [CL-] CLMINUSI = 1.0 / CLMINUS ! 1 / CLMINUS CHECK_NAME: SELECT CASE ( CNAME( 1:TRIMLEN( CNAME ) ) ) CASE ('SO2') ! SO2H2O <=> HSO3- + H+ ! & HSO3- <=> SO3= + H+ AKEQ1 = B( LSO2 ) * EXP( D( LSO2 ) * TFAC ) AKEQ2 = B( LHSO3 ) * EXP( D( LHSO3 ) * TFAC ) HLCONST = KH * ( 1.0 + AKEQ1 * HPLUSI + AKEQ1 * AKEQ2 * HPLUS2I ) CASE ('HNO2') ! HNO2(aq) <=> NO2- + H+ AKEQ1 = B( LHNO2 ) * EXP( D( LHNO2 ) * TFAC ) HLCONST = KH * ( 1.0 + AKEQ1 * HPLUSI ) CASE ('HNO3') ! HNO3(aq) <=> NO3- + H+ AKEQ1 = B( LHNO3 ) * EXP( D( LHNO3 ) * TFAC ) HLCONST = KH * ( 1.0 + AKEQ1 * HPLUSI ) CASE ('CO2') ! CO2H2O <=> HCO3- + H+ ! & HCO3- <=> CO3= + H+ AKEQ1 = B( LCO2 ) * EXP( D( LCO2 ) * TFAC ) AKEQ2 = B( LHCO3 ) * EXP( D( LHCO3 ) * TFAC ) HLCONST = KH & * ( 1.0 + AKEQ1 * HPLUSI + AKEQ1 * AKEQ2 * HPLUS2I ) CASE ('H2O2') ! H2O2(aq) <=> HO2- + H+ AKEQ1 = B( LH2O2 ) * EXP( D( LH2O2 ) * TFAC ) HLCONST = KH * ( 1.0 + AKEQ1 * HPLUSI ) CASE ('FORMALDEHYDE') ! HCHO(aq) <=> H2C(OH)2(aq) AKEQ1 = B( LHCHO ) * EXP( D( LHCHO ) * TFAC ) HLCONST = KH * ( 1.0 + AKEQ1 ) CASE ('FORMIC_ACID') ! HCOOH(aq) <=> HCOO- + H+ AKEQ1 = B( LHCOOH ) * EXP( D( LHCOOH ) * TFAC ) HLCONST = KH * ( 1.0 + AKEQ1 * HPLUSI ) CASE ('HO2') ! HO2(aq) <=> H+ + O2- AKEQ1 = B( LHO2 ) * EXP( D( LHO2 ) * TFAC ) HLCONST = KH * ( 1.0 + AKEQ1 * HPLUSI ) CASE ('NH3') ! NH4OH <=> NH4+ + OH- AKEQ1 = B( LNH4OH ) * EXP( D( LNH4OH ) * TFAC ) AKEQ2 = B( LH2O ) * EXP( D( LH2O ) * TFAC ) OHION = AKEQ2 * HPLUSI HLCONST = KH * ( 1.0 + AKEQ1 / OHION ) CASE ('HYDRAZINE') ! HYDRAZINE <=> HYDRAZINE+ + OH- AKEQ1 = B( LHYDRAZINE ) * EXP( D( LHYDRAZINE ) * TFAC ) AKEQ2 = B( LH2O ) * EXP( D( LH2O ) * TFAC ) OHION = AKEQ2 * HPLUSI HLCONST = KH * ( 1.0 + AKEQ1 / OHION ) CASE ('ATRA', 'DATRA') ! ATRA(aq) <=> ATRA- + H ! or DATRA(aq) <=> DATRA- + H AKEQ1 = B( LATRA ) * EXP( D( LATRA ) * TFAC ) HLCONST = KH * ( 1.0 + AKEQ1 * HPLUSI ) CASE ( 'CL2' ) ! CL2*H2O <=> HOCL + H + CL ! HOCL <=>H + OCL AKEQ1 = B( LCL2 ) * EXP( D( LCL2 ) * TFAC ) AKEQ2 = B( LHOCL ) * EXP( D( LHOCL ) * TFAC ) HLCONST = KH * ( 1.0 + AKEQ1 * HPLUSI * CLMINUSI & + AKEQ1 * AKEQ2 * HPLUS2I * CLMINUSI ) CASE ( 'HCL' ) ! HCL <=> H+ + CL- AKEQ1 = B( LHCL ) * EXP( D( LHCL ) * TFAC ) HLCONST = KH * ( 1.0 + AKEQ1 * HPLUSI ) CASE ( 'HOCL' ) ! HOCL <=> H+ + OCL- AKEQ1 = B( LHOCL ) * EXP( D( LHOCL ) * TFAC ) HLCONST = KH * ( 1.0 + AKEQ1 * HPLUSI ) END SELECT CHECK_NAME END IF RETURN END FUNCTION HLCONST !......................................................................... ! Version "@(#)$Header: /env/proj/archive/cvs/ioapi/./ioapi/src/index1.f,v 1.2 2000/11/28 21:22:49 smith_w Exp $" ! EDSS/Models-3 I/O API. Copyright (C) 1992-1999 MCNC ! Distributed under the GNU LESSER GENERAL PUBLIC LICENSE version 2.1 ! See file "LGPL.txt" for conditions of use. !......................................................................... INTEGER FUNCTION INDEX1 (NAME, N, NLIST) !*********************************************************************** ! subroutine body starts at line 46 ! ! FUNCTION: ! ! Searches for NAME in list NLIST and returns the subscript ! (1...N) at which it is found, or returns 0 when NAME not ! found in NLIST ! ! PRECONDITIONS REQUIRED: none ! ! SUBROUTINES AND FUNCTIONS CALLED: none ! ! REVISION HISTORY: ! ! 5/88 Modified for ROMNET ! 9/94 Modified for Models-3 by CJC ! !*********************************************************************** IMPLICIT NONE !....... Arguments and their descriptions: CHARACTER*(*) NAME ! Character string being searched for INTEGER N ! Length of array to be searched CHARACTER*(*) NLIST(*) ! array to be searched !....... Local variable: INTEGER I ! loop counter !..................................................................... !....... begin body of INDEX1() DO 100 I = 1, N IF ( NAME .EQ. NLIST( I ) ) THEN ! Found NAME in NLIST INDEX1 = I RETURN ENDIF 100 CONTINUE INDEX1 = 0 ! not found RETURN END FUNCTION INDEX1 END MODULE module_ctrans_aqchem