Line data    Source code

       1             : ! Include shortname defintions, so that the F77 code does not have to be modified to
       2             : ! reference the CARMA structure.
       3             : #include "carma_globaer.h"
       4             : 
       5             : !!  Calculates the vapor pressure for sulfuric acid.
       6             : !!
       7             : !!  <pvap_liq> and <pvap_ice> are vapor pressures in units of [dyne/cm^2]
       8             : !!
       9             : !!  Created   Dec-1995  (Ackerman)
      10             : !!  Modified  Sep-1997  (McKie)
      11             : !!  Modified Jul-2001 (Mills)
      12             : !!
      13             : !!  NOTE: To calculate vapor pressure of H2SO4 water vapor pressure (pvapl(iz, igash2o))
      14             : !!  should be calculated before this calculation.
      15             : !!
      16             : !! @author Mike Mills, Tianyi Fan
      17             : !! @version Feb-2011
      18  2914821276 : subroutine vaporp_H2SO4_Ayers1980(carma, cstate, iz, rc, pvap_liq, pvap_ice)
      19             : !     types
      20             :   use carma_precision_mod
      21             :   use carma_enums_mod
      22             :   use carma_constants_mod
      23             :   use carma_types_mod
      24             :   use carmastate_mod
      25             :   use carma_mod
      26             :   use sulfate_utils
      27             : 
      28             :   implicit none
      29             : 
      30             :   type(carma_type), intent(in)         :: carma     !! the carma object
      31             :   type(carmastate_type), intent(inout) :: cstate    !! the carma state object
      32             :   integer, intent(in)                  :: iz        !! z index
      33             :   real(kind=f), intent(out)            :: pvap_liq  !! vapor pressure wrt liquid [dyne/cm2]
      34             :   real(kind=f), intent(out)            :: pvap_ice  !! vapor pressure wrt ice [dyne[cm2]
      35             :   integer, intent(inout)               :: rc        !! return code, negative indicates failure
      36             : 
      37             :   ! Local declarations
      38             :   real(kind=f)            :: gc_cgs                          ! water vapor mass concentration [g/cm3]
      39             :   real(kind=f)            :: fk1, fk4, fk4_1, fk4_2
      40             :   real(kind=f)            :: factor_kulm                     ! Kulmala correction terms
      41             :   real(kind=f)            :: en, temp
      42             :   real(kind=f)            :: sulfeq
      43             :   real(kind=f), parameter :: t0_kulm     = 340._f            !  T0 set in the low end of the Ayers measurement range (338-445K)
      44             :   real(kind=f), parameter :: t_crit_kulm = 905._f            !  Critical temperature = 1.5 * Boiling point
      45             :   real(kind=f), parameter :: fk0         = -10156._f / t0_kulm + 16.259_f   !  Log(Kulmala correction factor)
      46             :   real(kind=f), parameter :: fk2         = 1._f / t0_kulm
      47             :   real(kind=f), parameter :: fk3         = 0.38_f / (t_crit_kulm - t0_kulm)
      48             : 
      49             : 
      50             :   ! Saturation vapor pressure of sulfuric acid
      51             :   !
      52             :   ! Don't allow saturation vapor pressure to underflow at very low temperatures
      53  2914821276 :   temp=max(t(iz),140._f)
      54             : 
      55             :   ! Convert water vapor concentration to g/cm3:
      56  2914821276 :   gc_cgs = gc(iz, igash2o) / zmet(iz)
      57             : 
      58             :   ! Compute the sulfate composition based on Hanson parameterization
      59             :   ! to temperature and water vapor concentration.
      60  2914821276 :   wtpct(iz) = wtpct_tabaz(carma, temp, gc_cgs, pvapl(iz, igash2o), rc)
      61             : 
      62             :   ! Parameterized fit to Giauque's (1959) enthalpies v. wt %:
      63  2914821276 :   en = 4.184_f * (23624.8_f - 1.14208e8_f / ((wtpct(iz) - 105.318_f)**2 + 4798.69_f))
      64  2914821276 :   en = max(en, 0.0_f)
      65             : 
      66             :   ! Ayers' (1980) fit to sulfuric acid equilibrium vapor pressure:
      67             :   ! (Remember this is the log)
      68             :   ! SULFEQ=-10156/Temp+16.259-En/(8.3143*Temp)
      69             :   !
      70             :   ! Kulmala correction (J. CHEM. PHYS. V.93, No.1, 1 July 1990)
      71  2914821276 :   fk1   = -1._f / temp
      72  2914821276 :   fk4_1 = log(t0_kulm / temp)
      73  2914821276 :   fk4_2 = t0_kulm / temp
      74  2914821276 :   fk4   = 1.0_f + fk4_1 - fk4_2
      75  2914821276 :   factor_kulm = fk1 + fk2 + fk3 * fk4
      76             : 
      77             :   ! This is for pure H2SO4
      78  2914821276 :   sulfeq = fk0 + 10156._f * factor_kulm
      79             : 
      80             :   ! Adjust for WTPCT composition:
      81  2914821276 :   sulfeq = sulfeq - en / (8.3143_f * temp)
      82             : 
      83             :   ! REMEMBER TO TAKE THE EXPONENTIAL!
      84  2914821276 :   sulfeq = exp(sulfeq)
      85             : 
      86             :   ! BUT this is in Atmospheres.  Convert ==> dynes/cm2
      87  2914821276 :   pvap_liq = sulfeq * 1.01325e6_f
      88  2914821276 :   pvap_ice = sulfeq * 1.01325e6_f
      89             : 
      90  2914821276 :   return
      91  2914821276 : end