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             : !! This routine evaluates particle loss rates due to nucleation <rnuclg>:
       6             : !! heterogeneous nucleation of glassy aerosols only,.
       7             : !!
       8             : !! The parameterization of glass aerosols is described in Murray et al.
       9             : !! [Nature Geosciences, 2010], and is based upon measurements of the nucleation of
      10             : !! citric acid aerosols at cold temperatures.
      11             : !!
      12             : !! NOTE: This implementation assumes that the aerosol being nucleated is the total
      13             : !! aerosol population and not just the fraction of aerosols that are glassy. To
      14             : !! account for homogenous freezing of the aerosol population, the routine freezaerl
      15             : !! also needs to be called and the overall nucleation rate is the sum of
      16             : !! the rates for homogeneous freezing and for heterogenous nucleation.
      17             : !!
      18             : !! The parameter fglass is the fraction of the total aerosol population that will be
      19             : !! in a glassy state for T <= 212K.
      20             : !! 
      21             : !! The loss rates for all particle elements in a particle group are equal.
      22             : !!
      23             : !! @author Chuck Bardeen, Eric Jensen
      24             : !! @version Apr-2010
      25   213271518 : subroutine freezglaerl_murray2010(carma, cstate, iz, rc)
      26             : 
      27             :   ! types
      28             :   use carma_precision_mod
      29             :   use carma_enums_mod
      30             :   use carma_constants_mod
      31             :   use carma_types_mod
      32             :   use carmastate_mod
      33             :   use carma_mod
      34             : 
      35             :   implicit none
      36             : 
      37             :   type(carma_type), intent(in)         :: carma   !! the carma object
      38             :   type(carmastate_type), intent(inout) :: cstate  !! the carma state object
      39             :   integer, intent(in)                  :: iz      !! z index
      40             :   integer, intent(inout)               :: rc      !! return code, negative indicates failure
      41             : 
      42             :   !  Local declarations
      43             :         ! Define parameters needed for freezing nucleation calculations.
      44             :   real(kind=f), parameter              :: kice1   = 7.7211e-5_f  ! Fit constant from Murray et al.
      45             :   real(kind=f), parameter              :: kice2   = 9.2688e-3_f  ! Fit constant from Murray et al.
      46             :   real(kind=f), parameter              :: ssmin   = 0.21_f       ! Minimum supersaturation for nucleation
      47             :   real(kind=f), parameter              :: ssmax   = 0.7_f        ! Maximum supersaturation for nucleation
      48             :   real(kind=f), parameter              :: tglass  = 212._f       ! Maximum temperature for glassy state
      49             :   real(kind=f), parameter              :: fglass  = 0.5_f        ! Fraction of aerosols that can become glassy
      50             :         
      51             :   integer                              :: igas    ! gas index
      52             :   integer                              :: igroup  ! group index
      53             :   integer                              :: ibin    ! bin index
      54             :   integer                              :: iepart  ! element for condensing group index
      55             :   integer                              :: inuc    ! nucleating element index
      56             :   integer                              :: ienucto ! index of target nucleation element
      57             :   integer                              :: ignucto ! index of target nucleation group
      58             :   integer                              :: inucto  ! index of target nucleation bin
      59             :   real(kind=f)                         :: dfice   ! difference in fraction of aerosol nucleated
      60             :   real(kind=f)                         :: ssi, ssiold
      61             :   
      62             :   ! Assume success.
      63   213271518 :   rc = RC_OK
      64             : 
      65             :   ! Loop over particle groups.
      66   639814554 :   do igroup = 1,NGROUP
      67             : 
      68   426543036 :     igas = inucgas(igroup)                ! condensing gas
      69   426543036 :     iepart = ienconc( igroup )            ! particle number density element
      70             : 
      71   639814554 :     if (igas /= 0) then
      72             : 
      73             :       ! Calculate nucleation loss rates.  Do not allow nucleation into
      74             :       ! an evaporating bin.
      75   426543036 :       do inuc = 1, nnuc2elem(iepart)
      76             : 
      77   213271518 :         ienucto = inuc2elem(inuc,iepart)
      78   426543036 :         if (ienucto /= 0) then
      79   213271518 :           ignucto = igelem(ienucto)
      80             :   
      81             :           ! Only compute nucleation rate for glassy aerosol freezing.
      82   213271518 :               if ((iand(inucproc(iepart,ienucto), I_AF_MURRAY_2010) /= 0)) then
      83             :           
      84             :             ! Is it cold enough for aerosols to be in a glassy state.
      85           0 :             if (t(iz) <= tglass) then
      86             :           
      87             :               !  Loop over particle bins.  Loop from largest to smallest for 
      88             :               !  evaluation of index of smallest bin nucleated during time step <inucstep>.
      89           0 :               do ibin = NBIN, 1, -1
      90             :       
      91             :                 ! Bypass calculation if few particles are present or if it isn't cold enough
      92             :                 ! for the aerosols to be present in a glassy state.
      93           0 :                 if (pconmax(iz,igroup) > FEW_PC) then
      94             :                 
      95             :                   ! Murray et al. [2010] doesn't really give a nucleation rate. Instead it gives
      96             :                   ! a fraction of glassy aerosol particles that have been nucleated as a function
      97             :                   ! of ice supersaturation.
      98             :                   !
      99             :                   ! Since CARMA really wants to work with rates, use the difference in relative
     100             :                   ! humidity and the length of the timestep to come up with an approximation to
     101             :                   ! a nucleation rate.
     102             :                   
     103             :                   ! The supersaturation must be greater than .21 for heterogeneous nucleation to
     104             :                   ! commence. The fraction of glassy aerosol nucleated is:
     105             :                   !
     106             :                   !   fice = 7.7211e-5 * RHi(%) - 9.2688e-3   for 121 % < RHi < 170 %
     107             :                   !
     108             :                   ! To get a pseudo production rate, use
     109             :                   !
     110             :                   !   rnuclg = (fice(RHi) - fice(RHi_old)) / dtime
     111             :                   !
     112           0 :                   ssi    = supsati(iz,igas)
     113           0 :                   ssiold = supsatiold(iz,igas)
     114             :                   
     115           0 :                   if ((ssi >= ssmin) .and. (ssi > ssiold)) then
     116           0 :                     dfice = kice1 * (1._f + min(ssmax, ssi)) * 100._f - kice2
     117             :                     
     118           0 :                     if (ssiold >= ssmin) then
     119           0 :                       dfice = dfice - (kice1 * (1._f + min(ssmax, ssiold)) * 100._f - kice2)
     120             :                     endif
     121             :                     
     122             :                     ! Add the rate of heterogenous freezing to the rate of homogeneous
     123           0 :                     rnuclg(ibin,igroup,ignucto) = rnuclg(ibin,igroup,ignucto) + fglass * dfice / dtime
     124             :                   endif    
     125             :                 endif
     126             :               enddo
     127             :             endif
     128             :           endif
     129             :         endif
     130             :       enddo 
     131             :     endif
     132             :   enddo
     133             : 
     134             :   ! Return to caller with particle loss rates due to nucleation evaluated.
     135   213271518 :   return
     136   213271518 : end