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 particle production rates due to nucleation <rhompe>:
       6             : !!  binary homogeneous nucleation of sulfuric acid and water only
       7             : !!  Numerical method follows Zhao & Turco, JAS, V.26, No.5, 1995.
       8             : !!
       9             : !!  @author Mike Mills, Chuck Bardeen
      10             : !!  @version Jun-2013
      11   213271518 : subroutine sulfnuc(carma,cstate, iz, rc) 
      12             :   use carma_precision_mod
      13             :   use carma_enums_mod
      14             :   use carma_constants_mod
      15             :   use carma_types_mod
      16             :   use carmastate_mod
      17             :   use carma_mod
      18             :   use sulfate_utils
      19             :   
      20             :   implicit none
      21             :   
      22             :   type(carma_type), intent(in)         :: carma       !! the carma object
      23             :   type(carmastate_type), intent(inout) :: cstate      !! the carma state object
      24             :   integer, intent(in)                  :: iz          !! level index
      25             :   integer, intent(inout)               :: rc          !! return code, negative indicates failure
      26             : 
      27             :   !  Local declarations     
      28             :   integer           :: igroup     ! group index
      29             :   integer           :: ibin       ! bin index
      30             :   integer           :: igas       ! gas index
      31             :   integer           :: iepart     ! concentration element index
      32             :   integer           :: nucbin     ! bin in which nucleation takes place
      33             :   integer           :: ignucto    ! index of target nucleation group
      34             :   integer           :: ienucto    ! index of target nucleation element
      35             :   integer           :: inuc
      36             :   real(kind=f)      :: nucrate    ! nucleation rate (#/x/y/z/s)
      37             :   real(kind=f)      :: h2o        ! H2O concentrations in molec/cm3 
      38             :   real(kind=f)      :: h2so4      ! H2SO4 concentrations in molec/cm3
      39             :   real(kind=f)      :: beta1
      40             :   real(kind=f)      :: beta2
      41             :   real(kind=f)      :: rstar      ! critical radius (cm)
      42             : 
      43             :   ! Cycle through each group, only proceed if BHN
      44   213271518 :   rstar = -1._f
      45             :   
      46   639814554 :   do igroup = 1 , NGROUP
      47             :     
      48   426543036 :     igas = inucgas(igroup)                ! condensing gas
      49             :     
      50   639814554 :     if (igas .ne. 0) then
      51             : 
      52   213271518 :       iepart = ienconc(igroup)              ! particle number density element
      53             : 
      54   213271518 :       if (inucproc(iepart,iepart) .eq. I_HOMNUC) then
      55             :     
      56             :         ! This is where all of the pre calculation needs to go, so that it isn't
      57             :         ! done when the model is not configured for homogeneous nucleation of
      58             :         ! sulfates.
      59   213271518 :         call sulfnucrate(carma, cstate, iz, igroup, h2so4, h2o, beta1, beta2, rstar, nucbin, nucrate, rc)
      60   213271518 :         if (rc /= RC_OK) return
      61             :         
      62             :         ! Do further calculations only if nucleation occurred
      63   213271518 :         if (nucrate .gt. 0._f) then
      64             : 
      65   185917403 :           rhompe(nucbin, iepart) = rhompe(nucbin, iepart) + nucrate
      66             :         
      67             :           ! Since homogeneous nucleation doesn't go through upgxfer or downgxfer, then
      68             :           ! then the effects of latent heat need to be accounted for here.
      69             :   !        rlprod = rlprod + rhompe(nucbin, ielem) * rmass(nucbin,igroup) * rlh_nuc(ielem,ielem) / (CP * rhoa(iz))
      70             :         end if
      71             :       end if
      72             :     end if
      73             :   end do
      74             : 
      75             :   ! Cycle through each group, only proceed if heterogeneous nucleation
      76             :   !
      77             :   ! NOTE: Only do heterogeneous nucleation if an rstar was determined by homogeneous
      78             :   ! nucleation.
      79   213271518 :   if (rstar > 0._f) then
      80   599160330 :     do igroup = 1 , NGROUP
      81             :     
      82   399440220 :       igas = inucgas(igroup)                ! condensing gas
      83             :     
      84   599160330 :       if (igas .ne. 0) then
      85             : 
      86   199720110 :         iepart = ienconc(igroup)              ! particle number density element
      87             : 
      88             :         ! Calculate heterogeneous nucleation loss rates.  Do not allow nucleation into
      89             :         ! an evaporating bin.
      90             :         !
      91             :         ! NOTE: Heterogeneous nucleation assumes that homogeneous nucleation was called
      92             :         ! first to determine the critical cluster size.
      93             :         !
      94             :         ! <ienucto> is index of target nucleation element;
      95             :         ! <ignucto> is index of target nucleation group.
      96   399440220 :         do inuc = 1, nnuc2elem(iepart)
      97             : 
      98   199720110 :           ienucto = inuc2elem(inuc,iepart)
      99             :         
     100   199720110 :           if (ienucto .ne. 0) then
     101   199720110 :             ignucto = igelem(ienucto)
     102             :           else
     103             :             ignucto = 0
     104             :           endif
     105             :         
     106   399440220 :           if (inucproc(iepart,ienucto) .eq. I_HETNUCSULF) then
     107             :     
     108           0 :             do ibin = NBIN, 1, -1
     109             : 
     110             :               ! Bypass calculation if few particles are present
     111           0 :               if (pconmax(iz,igroup) .gt. FEW_PC) then
     112             : 
     113             :                 ! This is where all of the pre calculation needs to go, so that it isn't
     114             :                 ! done when the model is not configured for homogeneous nucleation of
     115             :                 ! sulfates.
     116           0 :                 call sulfhetnucrate(carma, cstate, iz, igroup, ibin, h2so4, h2o, beta1, beta2, nucrate, rc)
     117           0 :                 if (rc /= RC_OK) return
     118             :                     
     119           0 :                 rnuclg(ibin, igroup, ignucto) = rnuclg(ibin, igroup, ignucto) + nucrate
     120             :               end if
     121             :             end do
     122             :           end if
     123             :         end do
     124             :       end if
     125             : 
     126             :     end do
     127             :   end if
     128             :   
     129             :   return
     130   213271518 : end