LCOV - code coverage report
Current view: top level - physics/cam - nucleate_ice_cam.F90 (source / functions) Hit Total Coverage
Test: coverage.info Lines: 286 402 71.1 %
Date: 2025-03-14 01:21:06 Functions: 4 4 100.0 %

          Line data    Source code
       1             : module nucleate_ice_cam
       2             : 
       3             : !---------------------------------------------------------------------------------
       4             : !
       5             : !  CAM Interfaces for nucleate_ice module.
       6             : !
       7             : !  B. Eaton - Sept 2014
       8             : !---------------------------------------------------------------------------------
       9             : 
      10             : use shr_kind_mod,   only: r8=>shr_kind_r8
      11             : use spmd_utils,     only: masterproc
      12             : use ppgrid,         only: pcols, pver
      13             : use physconst,      only: pi, rair, tmelt
      14             : use constituents,   only: pcnst, cnst_get_ind
      15             : use physics_types,  only: physics_state, physics_ptend, physics_ptend_init
      16             : use physics_buffer, only: physics_buffer_desc
      17             : use phys_control,   only: use_hetfrz_classnuc
      18             : use rad_constituents, only: rad_cnst_get_info, rad_cnst_get_aer_mmr, rad_cnst_get_aer_props
      19             : 
      20             : use physics_buffer, only: pbuf_add_field, dtype_r8, pbuf_old_tim_idx, &
      21             :                           pbuf_get_index, pbuf_get_field, &
      22             :                           pbuf_set_field
      23             : use cam_history,    only: addfld, add_default, outfld
      24             : 
      25             : use ref_pres,       only: top_lev => trop_cloud_top_lev
      26             : use wv_saturation,  only: qsat_water
      27             : 
      28             : use cam_logfile,    only: iulog
      29             : use cam_abortutils, only: endrun
      30             : 
      31             : use nucleate_ice,   only: nucleati_init, nucleati
      32             : 
      33             : use aerosol_properties_mod, only: aerosol_properties
      34             : use aerosol_state_mod, only: aerosol_state
      35             : 
      36             : use phys_control,   only: cam_physpkg_is
      37             : 
      38             : implicit none
      39             : private
      40             : save
      41             : 
      42             : public :: &
      43             :    nucleate_ice_cam_readnl,   &
      44             :    nucleate_ice_cam_register, &
      45             :    nucleate_ice_cam_init,     &
      46             :    nucleate_ice_cam_calc
      47             : 
      48             : ! Namelist variables
      49             : logical, public, protected :: use_preexisting_ice = .false.
      50             : logical                    :: hist_preexisting_ice = .false.
      51             : logical                    :: nucleate_ice_incloud = .false.
      52             : logical                    :: nucleate_ice_use_troplev = .false.
      53             : real(r8)                   :: nucleate_ice_subgrid = -1._r8
      54             : real(r8)                   :: nucleate_ice_subgrid_strat = -1._r8
      55             : real(r8)                   :: nucleate_ice_strat = 0.0_r8
      56             : 
      57             : ! Vars set via init method.
      58             : real(r8) :: mincld      ! minimum allowed cloud fraction
      59             : real(r8) :: bulk_scale  ! prescribed aerosol bulk sulfur scale factor
      60             : 
      61             : ! constituent indices
      62             : integer :: &
      63             :    cldliq_idx = -1, &
      64             :    cldice_idx = -1, &
      65             :    numice_idx = -1
      66             : 
      67             : integer :: &
      68             :    naai_idx = -1,     &
      69             :    naai_hom_idx = -1
      70             : 
      71             : integer :: &
      72             :    aist_idx = -1
      73             : 
      74             : integer :: &
      75             :     qsatfac_idx = -1
      76             : 
      77             : ! Bulk aerosols
      78             : character(len=20), allocatable :: aername(:)
      79             : real(r8), allocatable :: num_to_mass_aer(:)
      80             : 
      81             : integer :: naer_all = -1 ! number of aerosols affecting climate
      82             : integer :: idxsul   = -1 ! index in aerosol list for sulfate
      83             : integer :: idxdst1  = -1 ! index in aerosol list for dust1
      84             : integer :: idxdst2  = -1 ! index in aerosol list for dust2
      85             : integer :: idxdst3  = -1 ! index in aerosol list for dust3
      86             : integer :: idxdst4  = -1 ! index in aerosol list for dust4
      87             : integer :: idxbcphi = -1 ! index in aerosol list for Soot (BCPHIL)
      88             : 
      89             : ! MODAL or CARMA aerosols
      90             : logical :: clim_modal_carma = .false.
      91             : logical :: prog_modal_aero = .false.
      92             : 
      93             : logical :: lq(pcnst) = .false. ! set flags true for constituents with non-zero tendencies
      94             : 
      95             : integer, allocatable :: aer_cnst_idx(:,:)
      96             : 
      97             : !===============================================================================
      98             : contains
      99             : !===============================================================================
     100             : 
     101      243456 : subroutine nucleate_ice_cam_readnl(nlfile)
     102             : 
     103             :   use namelist_utils,  only: find_group_name
     104             :   use units,           only: getunit, freeunit
     105             :   use spmd_utils,      only: mpicom, masterprocid, mpi_logical, mpi_real8
     106             : 
     107             :   character(len=*), intent(in) :: nlfile  ! filepath for file containing namelist input
     108             : 
     109             :   ! Local variables
     110             :   integer :: unitn, ierr
     111             :   character(len=*), parameter :: subname = 'nucleate_ice_cam_readnl'
     112             : 
     113             :   namelist /nucleate_ice_nl/ use_preexisting_ice, hist_preexisting_ice, &
     114             :        nucleate_ice_subgrid, nucleate_ice_subgrid_strat, nucleate_ice_strat, &
     115             :        nucleate_ice_incloud, nucleate_ice_use_troplev
     116             : 
     117             :   !-----------------------------------------------------------------------------
     118             : 
     119        1536 :   if (masterproc) then
     120           2 :      unitn = getunit()
     121           2 :      open( unitn, file=trim(nlfile), status='old' )
     122           2 :      call find_group_name(unitn, 'nucleate_ice_nl', status=ierr)
     123           2 :      if (ierr == 0) then
     124           2 :         read(unitn, nucleate_ice_nl, iostat=ierr)
     125           2 :         if (ierr /= 0) then
     126           0 :            call endrun(subname // ':: ERROR reading namelist')
     127             :         end if
     128             :      end if
     129           2 :      close(unitn)
     130           2 :      call freeunit(unitn)
     131             : 
     132             :   end if
     133             : 
     134             :   ! Broadcast namelist variables
     135        1536 :   call mpi_bcast(use_preexisting_ice,  1, mpi_logical,masterprocid, mpicom, ierr)
     136        1536 :   call mpi_bcast(hist_preexisting_ice, 1, mpi_logical,masterprocid, mpicom, ierr)
     137        1536 :   call mpi_bcast(nucleate_ice_subgrid, 1, mpi_real8,  masterprocid, mpicom, ierr)
     138        1536 :   call mpi_bcast(nucleate_ice_subgrid_strat, 1, mpi_real8,  masterprocid, mpicom, ierr)
     139        1536 :   call mpi_bcast(nucleate_ice_strat,   1, mpi_real8,  masterprocid, mpicom, ierr)
     140        1536 :   call mpi_bcast(nucleate_ice_incloud, 1, mpi_logical,masterprocid, mpicom, ierr)
     141        1536 :   call mpi_bcast(nucleate_ice_use_troplev, 1, mpi_logical,masterprocid, mpicom, ierr)
     142             : 
     143        1536 : end subroutine nucleate_ice_cam_readnl
     144             : 
     145             : !================================================================================================
     146             : 
     147        1536 : subroutine nucleate_ice_cam_register()
     148             : 
     149             :    ! global scope for NAAI needed when clubb_do_icesuper=.true.
     150        1536 :    call pbuf_add_field('NAAI',     'global', dtype_r8, (/pcols,pver/), naai_idx)
     151        1536 :    call pbuf_add_field('NAAI_HOM', 'physpkg', dtype_r8, (/pcols,pver/), naai_hom_idx)
     152             : 
     153        1536 : end subroutine nucleate_ice_cam_register
     154             : 
     155             : !================================================================================================
     156             : 
     157        1536 : subroutine nucleate_ice_cam_init(mincld_in, bulk_scale_in, pbuf2d, aero_props)
     158             :    use phys_control, only: phys_getopts
     159             :    use time_manager, only: is_first_step
     160             : 
     161             :    real(r8), intent(in) :: mincld_in
     162             :    real(r8), intent(in) :: bulk_scale_in
     163             :    class(aerosol_properties), optional, intent(in) :: aero_props
     164             : 
     165             :    type(physics_buffer_desc), pointer :: pbuf2d(:,:)
     166             : 
     167             :    ! local variables
     168             :    integer :: iaer
     169             :    integer :: ierr
     170             :    integer :: ispc, ibin
     171             :    integer :: idxtmp
     172             :    integer :: nmodes, nbins
     173             : 
     174             :    character(len=*), parameter :: routine = 'nucleate_ice_cam_init'
     175             :    logical :: history_cesm_forcing
     176             : 
     177             :    character(len=32) :: tmpname
     178             : 
     179             :    !--------------------------------------------------------------------------------------------
     180        1536 :    call phys_getopts(prog_modal_aero_out = prog_modal_aero, history_cesm_forcing_out = history_cesm_forcing)
     181             : 
     182             :    ! clim_modal_aero determines whether modal or carma aerosols are used in the climate calculation.
     183             :    ! The modal aerosols can be either prognostic or prescribed.
     184        1536 :    call rad_cnst_get_info(0, nmodes=nmodes, nbins=nbins)
     185             : 
     186        1536 :    clim_modal_carma = (nmodes > 0) .or. (nbins > 0)
     187             : 
     188        1536 :    mincld     = mincld_in
     189        1536 :    bulk_scale = bulk_scale_in
     190             : 
     191        1536 :    lq(:) = .false.
     192             : 
     193        1536 :    if (clim_modal_carma.and.use_preexisting_ice) then
     194             : 
     195        1536 :       if (.not. present(aero_props)) then
     196           0 :          call endrun(routine//' :  aero_props must be present')
     197             :       end if
     198             : 
     199             :       ! constituent tendencies are calculated only if use_preexisting_ice is TRUE
     200             :       ! set lq for constituent tendencies --
     201             : 
     202       15360 :       allocate(aer_cnst_idx(aero_props%nbins(),0:maxval(aero_props%nspecies())), stat=ierr)
     203        1536 :       if( ierr /= 0 ) then
     204           0 :          call endrun(routine//': aer_cnst_idx allocation failed')
     205             :       end if
     206      204288 :       aer_cnst_idx = -1
     207             : 
     208       10752 :       do ibin = 1, aero_props%nbins()
     209        9216 :          if (aero_props%icenuc_updates_num(ibin)) then
     210             : 
     211             :             ! constituents of this bin will need to be updated
     212             : 
     213        1536 :             if (aero_props%icenuc_updates_mmr(ibin,0)) then ! species 0 indicates bin MMR
     214           0 :                call aero_props%amb_mmr_name( ibin, 0, tmpname)
     215             :             else
     216        1536 :                call aero_props%amb_num_name( ibin, tmpname)
     217             :             end if
     218             : 
     219        1536 :             call cnst_get_ind(tmpname, idxtmp, abort=.false.)
     220        1536 :             aer_cnst_idx(ibin,0) = idxtmp
     221        1536 :             if (idxtmp>0) then
     222        1536 :                lq(idxtmp) = .true.
     223             :             end if
     224             : 
     225             :             ! iterate over the species within the bin
     226        6144 :             do ispc = 1, aero_props%nspecies(ibin)
     227        6144 :                if (aero_props%icenuc_updates_mmr(ibin,ispc)) then
     228             :                   ! this aerosol constituent will be updated
     229        1536 :                   call aero_props%amb_mmr_name( ibin, ispc, tmpname)
     230        1536 :                   call cnst_get_ind(tmpname, idxtmp, abort=.false.)
     231        1536 :                   aer_cnst_idx(ibin,ispc) = idxtmp
     232        1536 :                   if (idxtmp>0) then
     233        1536 :                      lq(idxtmp) = .true.
     234             :                   end if
     235             :                end if
     236             :             end do
     237             : 
     238             :          end if
     239             :       end do
     240             : 
     241             :    end if
     242             : 
     243             :    ! Initialize naai.
     244        1536 :    if (is_first_step()) then
     245         768 :       call pbuf_set_field(pbuf2d, naai_idx, 0.0_r8)
     246             :    end if
     247             : 
     248        1536 :    if( masterproc ) then
     249           2 :       write(iulog,*) 'nucleate_ice parameters:'
     250           2 :       write(iulog,*) '  mincld                     = ', mincld_in
     251           2 :       write(iulog,*) '  bulk_scale                 = ', bulk_scale_in
     252           2 :       write(iulog,*) '  use_preexisiting_ice       = ', use_preexisting_ice
     253           2 :       write(iulog,*) '  hist_preexisiting_ice      = ', hist_preexisting_ice
     254           2 :       write(iulog,*) '  nucleate_ice_subgrid       = ', nucleate_ice_subgrid
     255           2 :       write(iulog,*) '  nucleate_ice_subgrid_strat = ', nucleate_ice_subgrid_strat
     256           2 :       write(iulog,*) '  nucleate_ice_strat         = ', nucleate_ice_strat
     257           2 :       write(iulog,*) '  nucleate_ice_incloud       = ', nucleate_ice_incloud
     258           2 :       write(iulog,*) '  nucleate_ice_use_troplev   = ', nucleate_ice_use_troplev
     259             :    end if
     260             : 
     261        1536 :    call cnst_get_ind('CLDLIQ', cldliq_idx)
     262        1536 :    call cnst_get_ind('CLDICE', cldice_idx)
     263        1536 :    call cnst_get_ind('NUMICE', numice_idx)
     264        1536 :    qsatfac_idx  = pbuf_get_index('QSATFAC', ierr)
     265             : 
     266        1536 :    if (((nucleate_ice_subgrid .eq. -1._r8) .or. (nucleate_ice_subgrid_strat .eq. -1._r8)) .and. (qsatfac_idx .eq. -1)) then
     267           0 :      call endrun(routine//': ERROR qsatfac is required when subgrid = -1 or subgrid_strat = -1')
     268             :    end if
     269             : 
     270        1536 :    if (cam_physpkg_is("cam7")) then
     271             :       ! Updates for PUMAS v1.21+
     272           0 :       call addfld('NIHFTEN',  (/ 'lev' /), 'A', '1/m3/s', 'Activated Ice Number Concentration tendency due to homogenous freezing', sampled_on_subcycle=.true.)
     273           0 :       call addfld('NIDEPTEN', (/ 'lev' /), 'A', '1/m3/s', 'Activated Ice Number Concentration tendency due to deposition nucleation', sampled_on_subcycle=.true.)
     274           0 :       call addfld('NIIMMTEN', (/ 'lev' /), 'A', '1/m3/s', 'Activated Ice Number Concentration tendency due to immersion freezing', sampled_on_subcycle=.true.)
     275           0 :       call addfld('NIMEYTEN', (/ 'lev' /), 'A', '1/m3/s', 'Activated Ice Number Concentration tendency due to meyers deposition', sampled_on_subcycle=.true.)
     276             :    else
     277        3072 :       call addfld('NIHF',  (/ 'lev' /), 'A', '1/m3', 'Activated Ice Number Concentration due to homogenous freezing', sampled_on_subcycle=.true.)
     278        3072 :       call addfld('NIDEP', (/ 'lev' /), 'A', '1/m3', 'Activated Ice Number Concentration due to deposition nucleation', sampled_on_subcycle=.true.)
     279        3072 :       call addfld('NIIMM', (/ 'lev' /), 'A', '1/m3', 'Activated Ice Number Concentration due to immersion freezing', sampled_on_subcycle=.true.)
     280        3072 :       call addfld('NIMEY', (/ 'lev' /), 'A', '1/m3', 'Activated Ice Number Concentration due to meyers deposition', sampled_on_subcycle=.true.)
     281             :    endif
     282             : 
     283        3072 :    call addfld('NIREGM',(/ 'lev' /), 'A', 'C', 'Ice Nucleation Temperature Threshold for Regime', sampled_on_subcycle=.true.)
     284        3072 :    call addfld('NISUBGRID',(/ 'lev' /), 'A', '', 'Ice Nucleation subgrid saturation factor', sampled_on_subcycle=.true.)
     285        3072 :    call addfld('NITROP_PD',(/ 'lev' /), 'A', '', 'Chemical Tropopause probability', sampled_on_subcycle=.true.)
     286        1536 :    if ( history_cesm_forcing ) then
     287           0 :       call add_default('NITROP_PD',8,' ')
     288             :    endif
     289             : 
     290        1536 :    if (use_preexisting_ice) then
     291        3072 :       call addfld('fhom',      (/ 'lev' /), 'A','fraction', 'Fraction of cirrus where homogeneous freezing occur', sampled_on_subcycle=.true.)
     292        3072 :       call addfld ('WICE',     (/ 'lev' /), 'A','m/s','Vertical velocity Reduction caused by preexisting ice', sampled_on_subcycle=.true.)
     293        3072 :       call addfld ('WEFF',     (/ 'lev' /), 'A','m/s','Effective Vertical velocity for ice nucleation', sampled_on_subcycle=.true.)
     294             : 
     295        1536 :       if (cam_physpkg_is("cam7")) then
     296             :          ! Updates for PUMAS v1.21+
     297           0 :          call addfld ('INnso4TEN',   (/ 'lev' /), 'A','1/m3/s','Number Concentration tendency so4 (in) to ice_nucleation', sampled_on_subcycle=.true.)
     298           0 :          call addfld ('INnbcTEN',    (/ 'lev' /), 'A','1/m3/s','Number Concentration tendency bc  (in) to ice_nucleation', sampled_on_subcycle=.true.)
     299           0 :          call addfld ('INndustTEN',  (/ 'lev' /), 'A','1/m3/s','Number Concentration tendency dust (in) ice_nucleation', sampled_on_subcycle=.true.)
     300           0 :          call addfld ('INondustTEN',  (/ 'lev' /), 'A','1/m3/s','Number Concentration tendency dust (out) from ice_nucleation', sampled_on_subcycle=.true.)
     301             :          call addfld ('INhetTEN',    (/ 'lev' /), 'A','1/m3/s', &
     302           0 :               'Tendency for contribution for in-cloud ice number density increase by het nucleation in ice cloud', sampled_on_subcycle=.true.)
     303             :          call addfld ('INhomTEN',    (/ 'lev' /), 'A','1/m3/s', &
     304           0 :               'Tendency for contribution for in-cloud ice number density increase by hom nucleation in ice cloud', sampled_on_subcycle=.true.)
     305             :       else
     306        3072 :          call addfld ('INnso4',   (/ 'lev' /), 'A','1/m3','Number Concentration so4 (in) to ice_nucleation', sampled_on_subcycle=.true.)
     307        3072 :          call addfld ('INnbc',    (/ 'lev' /), 'A','1/m3','Number Concentration bc  (in) to ice_nucleation', sampled_on_subcycle=.true.)
     308        3072 :          call addfld ('INndust',  (/ 'lev' /), 'A','1/m3','Number Concentration dust (in) ice_nucleation', sampled_on_subcycle=.true.)
     309        3072 :          call addfld ('INondust',  (/ 'lev' /), 'A','1/m3','Number Concentration dust (out) from ice_nucleation', sampled_on_subcycle=.true.)
     310             :          call addfld ('INhet',    (/ 'lev' /), 'A','1/m3', &
     311        3072 :               'contribution for in-cloud ice number density increase by het nucleation in ice cloud', sampled_on_subcycle=.true.)
     312             :          call addfld ('INhom',    (/ 'lev' /), 'A','1/m3', &
     313        3072 :               'contribution for in-cloud ice number density increase by hom nucleation in ice cloud', sampled_on_subcycle=.true.)
     314             :       endif
     315             : 
     316        3072 :       call addfld ('INFrehom', (/ 'lev' /), 'A','frequency','hom IN frequency ice cloud', sampled_on_subcycle=.true.)
     317        3072 :       call addfld ('INFreIN',  (/ 'lev' /), 'A','frequency','frequency of ice nucleation occur', sampled_on_subcycle=.true.)
     318             : 
     319        1536 :       if (hist_preexisting_ice) then
     320           0 :          call add_default ('WSUBI   ', 1, ' ')  ! addfld/outfld calls are in microp_aero
     321             : 
     322           0 :          call add_default ('fhom    ', 1, ' ')
     323           0 :          call add_default ('WICE    ', 1, ' ')
     324           0 :          call add_default ('WEFF    ', 1, ' ')
     325           0 :          call add_default ('INnso4  ', 1, ' ')
     326           0 :          call add_default ('INnbc   ', 1, ' ')
     327           0 :          call add_default ('INndust ', 1, ' ')
     328           0 :          call add_default ('INhet   ', 1, ' ')
     329           0 :          call add_default ('INhom   ', 1, ' ')
     330           0 :          call add_default ('INFrehom', 1, ' ')
     331           0 :          call add_default ('INFreIN ', 1, ' ')
     332             :       end if
     333             :    end if
     334             : 
     335        1536 :    if (.not. clim_modal_carma) then
     336             : 
     337             :       ! Props needed for BAM number concentration calcs.
     338             : 
     339           0 :       call rad_cnst_get_info(0, naero=naer_all)
     340             :       allocate( &
     341           0 :          aername(naer_all),        &
     342           0 :          num_to_mass_aer(naer_all) )
     343             : 
     344           0 :       do iaer = 1, naer_all
     345             :          call rad_cnst_get_aer_props(0, iaer, &
     346           0 :             aername         = aername(iaer), &
     347           0 :             num_to_mass_aer = num_to_mass_aer(iaer))
     348             :          ! Look for sulfate, dust, and soot in this list (Bulk aerosol only)
     349           0 :          if (trim(aername(iaer)) == 'SULFATE') idxsul = iaer
     350           0 :          if (trim(aername(iaer)) == 'DUST1') idxdst1 = iaer
     351           0 :          if (trim(aername(iaer)) == 'DUST2') idxdst2 = iaer
     352           0 :          if (trim(aername(iaer)) == 'DUST3') idxdst3 = iaer
     353           0 :          if (trim(aername(iaer)) == 'DUST4') idxdst4 = iaer
     354           0 :          if (trim(aername(iaer)) == 'BCPHI') idxbcphi = iaer
     355             :       end do
     356             :    end if
     357             : 
     358             : 
     359             :    call nucleati_init(use_preexisting_ice, use_hetfrz_classnuc, nucleate_ice_incloud, iulog, pi, &
     360        1536 :         mincld)
     361             : 
     362             :    ! get indices for fields in the physics buffer
     363        1536 :    aist_idx = pbuf_get_index('AIST')
     364             : 
     365        1536 : end subroutine nucleate_ice_cam_init
     366             : 
     367             : !================================================================================================
     368             : 
     369    58544640 : subroutine nucleate_ice_cam_calc( &
     370      241920 :    state, wsubi, pbuf, dtime, ptend, aero_props, aero_state )
     371             : 
     372        1536 :    use tropopause,     only: tropopause_findChemTrop
     373             : 
     374             :    ! arguments
     375             :    type(physics_state), target, intent(in)    :: state
     376             :    real(r8),                    intent(in)    :: wsubi(:,:)
     377             :    type(physics_buffer_desc),   pointer       :: pbuf(:)
     378             :    real(r8),                    intent(in)    :: dtime
     379             :    type(physics_ptend),         intent(out)   :: ptend
     380             :    class(aerosol_properties),optional, intent(in) :: aero_props
     381             :    class(aerosol_state),optional, intent(in) :: aero_state
     382             : 
     383             :    ! local workspace
     384             : 
     385             :    ! naai and naai_hom are the outputs shared with the microphysics
     386      241920 :    real(r8), pointer :: naai(:,:)       ! number of activated aerosol for ice nucleation
     387      241920 :    real(r8), pointer :: naai_hom(:,:)   ! number of activated aerosol for ice nucleation (homogeneous freezing only)
     388             : 
     389             :    integer :: lchnk, ncol
     390             :    integer :: itim_old
     391             :    integer :: i, k, l, m
     392             : 
     393             :    character(len=32) :: spectype
     394             : 
     395      241920 :    real(r8), pointer :: t(:,:)          ! input temperature (K)
     396      241920 :    real(r8), pointer :: qn(:,:)         ! input water vapor mixing ratio (kg/kg)
     397      241920 :    real(r8), pointer :: qc(:,:)         ! cloud water mixing ratio (kg/kg)
     398      241920 :    real(r8), pointer :: qi(:,:)         ! cloud ice mixing ratio (kg/kg)
     399      241920 :    real(r8), pointer :: ni(:,:)         ! cloud ice number conc (1/kg)
     400      241920 :    real(r8), pointer :: pmid(:,:)       ! pressure at layer midpoints (pa)
     401             : 
     402      241920 :    real(r8), pointer :: aer_mmr(:,:)    ! aerosol mass mixing ratio
     403      241920 :    real(r8), pointer :: aist(:,:)
     404             :    real(r8) :: icecldf(pcols,pver)  ! ice cloud fraction
     405      241920 :    real(r8), pointer :: qsatfac(:,:)      ! Subgrid cloud water saturation scaling factor.
     406             : 
     407             :    real(r8) :: rho(pcols,pver)      ! air density (kg m-3)
     408             : 
     409      241920 :    real(r8), allocatable :: naer2(:,:,:)    ! bulk aerosol number concentration (1/m3)
     410      241920 :    real(r8), allocatable :: maerosol(:,:,:) ! bulk aerosol mass conc (kg/m3)
     411             : 
     412             :    real(r8) :: qs(pcols)            ! liquid-ice weighted sat mixing rat (kg/kg)
     413             :    real(r8) :: es(pcols)            ! liquid-ice weighted sat vapor press (pa)
     414             :    real(r8) :: gammas(pcols)        ! parameter for cond/evap of cloud water
     415             :    integer  :: troplev(pcols)       ! tropopause level
     416             : 
     417             :    real(r8) :: relhum(pcols,pver)  ! relative humidity
     418             :    real(r8) :: icldm(pcols,pver)   ! ice cloud fraction
     419             : 
     420             :    real(r8) :: dst_num                               ! total dust aerosol number (#/cm^3)
     421             :    real(r8) :: dso4_num                               ! so4 aerosol number (#/cm^3)
     422             :    real(r8) :: so4_num                               ! so4 aerosol number (#/cm^3)
     423             :    real(r8) :: soot_num                              ! soot (hydrophilic) aerosol number (#/cm^3)
     424             :    real(r8) :: wght
     425             :    real(r8) :: oso4_num
     426             :    real(r8) :: odst_num
     427             :    real(r8) :: osoot_num
     428             :    real(r8) :: so4_num_st_cr_tot
     429             :    real(r8) :: ramp
     430             : 
     431             :    real(r8) :: subgrid(pcols,pver)
     432             :    real(r8) :: trop_pd(pcols,pver)
     433             : 
     434             :    ! For pre-existing ice
     435             :    real(r8) :: fhom(pcols,pver)    ! how much fraction of cloud can reach Shom
     436             :    real(r8) :: wice(pcols,pver)    ! diagnosed Vertical velocity Reduction caused by preexisting ice (m/s), at Shom
     437             :    real(r8) :: weff(pcols,pver)    ! effective Vertical velocity for ice nucleation (m/s); weff=wsubi-wice
     438             :    real(r8) :: INnso4(pcols,pver)   ! #/m3, so4 aerosol number used for ice nucleation
     439             :    real(r8) :: INnbc(pcols,pver)    ! #/m3, bc aerosol number used for ice nucleation
     440             :    real(r8) :: INndust(pcols,pver)  ! #/m3, dust aerosol number used for ice nucleation
     441             :    real(r8) :: INondust(pcols,pver)  ! #/m3, dust aerosol number used for ice nucleation
     442             :    real(r8) :: INhet(pcols,pver)    ! #/m3, ice number from het freezing
     443             :    real(r8) :: INhom(pcols,pver)    ! #/m3, ice number from hom freezing
     444             :    real(r8) :: INFrehom(pcols,pver) !  hom freezing occurence frequency.  1 occur, 0 not occur.
     445             :    real(r8) :: INFreIN(pcols,pver)  !  ice nucleation occerence frequency.   1 occur, 0 not occur.
     446             : 
     447             :    ! history output for ice nucleation
     448             :    real(r8) :: nihf(pcols,pver)  !output number conc of ice nuclei due to heterogenous freezing (1/m3)
     449             :    real(r8) :: niimm(pcols,pver) !output number conc of ice nuclei due to immersion freezing (hetero nuc) (1/m3)
     450             :    real(r8) :: nidep(pcols,pver) !output number conc of ice nuclei due to deoposion nucleation (hetero nuc) (1/m3)
     451             :    real(r8) :: nimey(pcols,pver) !output number conc of ice nuclei due to meyers deposition (1/m3)
     452             :    real(r8) :: regm(pcols,pver)  !output temperature thershold for nucleation regime
     453             : 
     454             :    real(r8) :: size_wghts(pcols,pver)
     455             :    real(r8) :: type_wghts(pcols,pver)
     456             :    real(r8), pointer :: num_col(:,:)
     457             :    real(r8) :: dust_num_col(pcols,pver)
     458             :    real(r8) :: sulf_num_col(pcols,pver)
     459             :    real(r8) :: soot_num_col(pcols,pver)
     460             :    real(r8) :: sulf_num_tot_col(pcols,pver)
     461             : 
     462             :    integer :: idxtmp
     463      241920 :    real(r8), pointer :: amb_num(:,:)
     464      241920 :    real(r8), pointer :: amb_mmr(:,:)
     465             :    real(r8), pointer :: cld_num(:,:)
     466      241920 :    real(r8), pointer :: cld_mmr(:,:)
     467             : 
     468             :    real(r8) :: delmmr, delmmr_sum
     469             :    real(r8) :: delnum, delnum_sum
     470             : 
     471             :    real(r8), parameter :: per_cm3 = 1.e-6_r8 ! factor for m-3 to cm-3 conversions
     472             : 
     473             :    integer :: nbins, nmaxspc
     474      241920 :    real(r8), allocatable :: amb_num_bins(:,:,:)
     475      241920 :    real(r8), allocatable :: size_wght(:,:,:,:)
     476             : 
     477             :    !-------------------------------------------------------------------------------
     478             : 
     479      241920 :    lchnk = state%lchnk
     480      241920 :    ncol  = state%ncol
     481      241920 :    t     => state%t
     482      241920 :    qn    => state%q(:,:,1)
     483      241920 :    qc    => state%q(:,:,cldliq_idx)
     484      241920 :    qi    => state%q(:,:,cldice_idx)
     485      241920 :    ni    => state%q(:,:,numice_idx)
     486      241920 :    pmid  => state%pmid
     487             : 
     488      241920 :    if (present(aero_props)) then
     489      241920 :       nbins = aero_props%nbins()
     490     1451520 :       nmaxspc = maxval(aero_props%nspecies())
     491             : 
     492     1209600 :       allocate(size_wght(ncol,pver,nbins,nmaxspc))
     493      967680 :       allocate(amb_num_bins(ncol,pver,nbins))
     494             :    else
     495      241920 :       nbins = 0
     496      241920 :       nmaxspc = 0
     497             :    endif
     498             : 
     499   119460096 :    rho(:ncol,:) = pmid(:ncol,:)/(rair*t(:ncol,:))
     500             : 
     501      241920 :    if (clim_modal_carma) then
     502             : 
     503      241920 :       call physics_ptend_init(ptend, state%psetcols, 'nucleatei', lq=lq)
     504             : 
     505             :    else
     506             :       ! init number/mass arrays for bulk aerosols
     507             :       allocate( &
     508             :            naer2(pcols,pver,naer_all), &
     509           0 :            maerosol(pcols,pver,naer_all))
     510             : 
     511           0 :       do m = 1, naer_all
     512           0 :          call rad_cnst_get_aer_mmr(0, m, state, pbuf, aer_mmr)
     513           0 :          maerosol(:ncol,:,m) = aer_mmr(:ncol,:)*rho(:ncol,:)
     514             : 
     515           0 :          if (m .eq. idxsul) then
     516           0 :             naer2(:ncol,:,m) = maerosol(:ncol,:,m)*num_to_mass_aer(m)*bulk_scale
     517             :          else
     518           0 :             naer2(:ncol,:,m) = maerosol(:ncol,:,m)*num_to_mass_aer(m)
     519             :          end if
     520             :       end do
     521             : 
     522           0 :       call physics_ptend_init(ptend, state%psetcols, 'nucleatei')
     523             :    end if
     524             : 
     525      241920 :    itim_old = pbuf_old_tim_idx()
     526      967680 :    call pbuf_get_field(pbuf, aist_idx, aist, start=(/1,1,itim_old/), kount=(/pcols,pver,1/))
     527   119460096 :    icecldf(:ncol,:pver) = aist(:ncol,:pver)
     528             : 
     529             :    ! naai and naai_hom are the outputs from this parameterization
     530      241920 :    call pbuf_get_field(pbuf, naai_idx, naai)
     531      241920 :    call pbuf_get_field(pbuf, naai_hom_idx, naai_hom)
     532   119460096 :    naai(1:ncol,1:pver)     = 0._r8
     533   119460096 :    naai_hom(1:ncol,1:pver) = 0._r8
     534             : 
     535             :    ! Use the same criteria that is used in chemistry and in CLUBB (for cloud fraction)
     536             :    ! to determine whether to use tropospheric or stratospheric settings. Include the
     537             :    ! tropopause level so that the cold point tropopause will use the stratospheric values.
     538             :    !REMOVECAM - no longer need this when CAM is retired and pcols no longer exists
     539      241920 :    troplev(:) = 0
     540             :    !REMOVECAM_END
     541      241920 :    call tropopause_findChemTrop(state, troplev)
     542             : 
     543      241920 :    if ((nucleate_ice_subgrid .eq. -1._r8) .or. (nucleate_ice_subgrid_strat .eq. -1._r8)) then
     544           0 :       call pbuf_get_field(pbuf, qsatfac_idx, qsatfac)
     545             :    end if
     546             : 
     547      241920 :    trop_pd(:,:) = 0._r8
     548             : 
     549     7983360 :    do k = top_lev, pver
     550   119460096 :       do i = 1, ncol
     551   111476736 :          trop_pd(i, troplev(i)) = 1._r8
     552             : 
     553   119218176 :          if (k <= troplev(i)) then
     554    50802339 :             if (nucleate_ice_subgrid_strat .eq. -1._r8) then
     555           0 :                subgrid(i, k) = 1._r8 / qsatfac(i, k)
     556             :             else
     557    50802339 :                subgrid(i, k) = nucleate_ice_subgrid_strat
     558             :             end if
     559             :          else
     560    60674397 :             if (nucleate_ice_subgrid .eq. -1._r8) then
     561           0 :                subgrid(i, k) = 1._r8 / qsatfac(i, k)
     562             :             else
     563    60674397 :                subgrid(i, k) = nucleate_ice_subgrid
     564             :             end if
     565             :          end if
     566             :       end do
     567             :    end do
     568             : 
     569             : 
     570             :    ! initialize history output fields for ice nucleation
     571   119460096 :    nihf(1:ncol,1:pver)  = 0._r8
     572   119460096 :    niimm(1:ncol,1:pver) = 0._r8
     573   119460096 :    nidep(1:ncol,1:pver) = 0._r8
     574   119460096 :    nimey(1:ncol,1:pver) = 0._r8
     575             : 
     576   119460096 :    regm(1:ncol,1:pver) = 0._r8
     577             : 
     578      241920 :    if (use_preexisting_ice) then
     579      241920 :       fhom(:,:)     = 0.0_r8
     580      241920 :       wice(:,:)     = 0.0_r8
     581      241920 :       weff(:,:)     = 0.0_r8
     582      241920 :       INnso4(:,:)   = 0.0_r8
     583      241920 :       INnbc(:,:)    = 0.0_r8
     584      241920 :       INndust(:,:)  = 0.0_r8
     585      241920 :       INondust(:,:)  = 0.0_r8
     586      241920 :       INhet(:,:)    = 0.0_r8
     587      241920 :       INhom(:,:)    = 0.0_r8
     588      241920 :       INFrehom(:,:) = 0.0_r8
     589      241920 :       INFreIN(:,:)  = 0.0_r8
     590             :    endif
     591             : 
     592     7983360 :    do k = top_lev, pver
     593             :       ! Get humidity and saturation vapor pressures
     594     7741440 :       call qsat_water(t(1:ncol,k), pmid(1:ncol,k), es(1:ncol), qs(1:ncol), ncol, gam=gammas(1:ncol))
     595             : 
     596   119460096 :       do i = 1, ncol
     597             : 
     598   111476736 :          relhum(i,k) = qn(i,k)/qs(i)
     599             : 
     600             :          ! get cloud fraction, check for minimum
     601   119218176 :          icldm(i,k) = max(icecldf(i,k), mincld)
     602             : 
     603             :       end do
     604             :    end do
     605             : 
     606      241920 :    dust_num_col = 0._r8
     607      241920 :    sulf_num_col = 0._r8
     608      241920 :    sulf_num_tot_col = 0._r8
     609      241920 :    soot_num_col = 0._r8
     610             : 
     611      241920 :    if (clim_modal_carma) then
     612             : 
     613      725760 :       if (.not.(present(aero_props).and.present(aero_state))) then
     614           0 :          call endrun('nucleate_ice_cam_calc: aero_props and aero_state must be present')
     615             :       end if
     616             : 
     617             :       ! collect number densities (#/cm^3) for dust, sulfate, and soot
     618             :       call aero_state%nuclice_get_numdens( aero_props, use_preexisting_ice, ncol, pver, rho, &
     619      241920 :                                            dust_num_col, sulf_num_col, soot_num_col, sulf_num_tot_col )
     620             : 
     621     1451520 :       do m = 1, aero_props%nbins()
     622     1209600 :          call aero_state%get_ambient_num(m, amb_num)
     623   597300480 :          amb_num_bins(:ncol,:,m) = amb_num(:ncol,:)
     624    12579840 :          do l = 1, aero_props%nspecies(m)
     625    11128320 :             call aero_props%species_type(m, l, spectype)
     626    12337920 :             call aero_state%icenuc_size_wght( m, ncol, pver, spectype, use_preexisting_ice, size_wght(:,:,m,l))
     627             : 
     628             :             !size_wght(:ncol,:,m,l) = wght(:ncol,:)
     629             :          end do
     630             :       end do
     631             : 
     632             :    else
     633             :       ! for bulk model
     634           0 :       if (idxdst1 > 0 .and. idxdst2 > 0 .and. idxdst3 > 0 .and. idxdst4 > 0) then
     635           0 :          dust_num_col(:ncol,:) = naer2(:ncol,:,idxdst1)/25._r8 * per_cm3 & ! #/cm3
     636           0 :                                  + naer2(:ncol,:,idxdst2)/25._r8 * per_cm3 &
     637           0 :                                  + naer2(:ncol,:,idxdst3)/25._r8 * per_cm3 &
     638           0 :                                  + naer2(:ncol,:,idxdst4)/25._r8 * per_cm3
     639             :       end if
     640           0 :       if (idxsul > 0) then
     641           0 :          sulf_num_col(:ncol,:) = naer2(:ncol,:,idxsul)/25._r8 * per_cm3
     642             :       end if
     643           0 :       if (idxbcphi > 0) then
     644           0 :          soot_num_col(:ncol,:) = naer2(:ncol,:,idxbcphi)/25._r8 * per_cm3
     645             :       end if
     646             :    endif
     647             : 
     648     7983360 :    kloop: do k = top_lev, pver
     649   119460096 :       iloop: do i = 1, ncol
     650             : 
     651   111476736 :          so4_num_st_cr_tot = 0._r8
     652             : 
     653   119218176 :          freezing: if (t(i,k) < tmelt - 5._r8) then
     654             : 
     655             :             ! set aerosol number for so4, soot, and dust with units #/cm^3
     656    87798058 :             so4_num = sulf_num_col(i,k)
     657    87798058 :             dst_num = dust_num_col(i,k)
     658    87798058 :             so4_num_st_cr_tot=sulf_num_tot_col(i,k)
     659             : 
     660             :             ! *** Turn off soot nucleation ***
     661    87798058 :             soot_num = 0.0_r8
     662             : 
     663    87798058 :             if (cam_physpkg_is("cam7")) then
     664             : 
     665             :                call nucleati( &
     666           0 :                     wsubi(i,k), t(i,k), pmid(i,k), relhum(i,k), icldm(i,k),   &
     667           0 :                     qc(i,k), qi(i,k), ni(i,k), rho(i,k),                      &
     668             :                     so4_num, dst_num, soot_num, subgrid(i,k),                 &
     669           0 :                     naai(i,k), nihf(i,k), niimm(i,k), nidep(i,k), nimey(i,k), &
     670             :                     wice(i,k), weff(i,k), fhom(i,k), regm(i,k),               &
     671             :                     oso4_num, odst_num, osoot_num, &
     672           0 :                     call_frm_zm_in = .false., add_preexisting_ice_in = .false.)
     673             : 
     674             :             else
     675             : 
     676             :                call nucleati( &
     677   175596116 :                     wsubi(i,k), t(i,k), pmid(i,k), relhum(i,k), icldm(i,k),   &
     678    87798058 :                     qc(i,k), qi(i,k), ni(i,k), rho(i,k),                      &
     679             :                     so4_num, dst_num, soot_num, subgrid(i,k),                 &
     680           0 :                     naai(i,k), nihf(i,k), niimm(i,k), nidep(i,k), nimey(i,k), &
     681             :                     wice(i,k), weff(i,k), fhom(i,k), regm(i,k),               &
     682   263394174 :                     oso4_num, odst_num, osoot_num)
     683             : 
     684             :             end if
     685             : 
     686             :             ! Move aerosol used for nucleation from interstial to cloudborne,
     687             :             ! otherwise the same coarse mode aerosols will be available again
     688             :             ! in the next timestep and will supress homogeneous freezing.
     689             : 
     690             : 
     691    87798058 :             if (clim_modal_carma .and. use_preexisting_ice) then
     692             : 
     693             :                ! compute tendencies for transported aerosol constituents
     694             :                ! and update not-transported constituents
     695             : 
     696   526788348 :                do m = 1, aero_props%nbins()
     697             : 
     698   526788348 :                   if (aero_props%icenuc_updates_num(m)) then
     699             : 
     700             :                      ! constituents of this bin will need to be updated
     701             : 
     702    87798058 :                      if (amb_num_bins(i,k,m)>0._r8) then
     703    87798058 :                         delmmr_sum = 0._r8
     704    87798058 :                         delnum_sum = 0._r8
     705             : 
     706             :                         ! iterate over the species within the bin
     707   351192232 :                         do l = 1, aero_props%nspecies(m)
     708   351192232 :                            if (aero_props%icenuc_updates_mmr(m,l)) then
     709             : 
     710    87798058 :                               call aero_props%species_type(m, l, spectype)
     711             : 
     712    87798058 :                               wght = size_wght(i,k,m,l)
     713             : 
     714    87798058 :                               if (wght>0._r8) then
     715             : 
     716             :                                  ! this aerosol constituent will be updated
     717             : 
     718    87798058 :                                  idxtmp = aer_cnst_idx(m,l)
     719             : 
     720    87798058 :                                  call aero_state%get_ambient_mmr(l,m,amb_mmr)
     721    87798058 :                                  call aero_state%get_cldbrne_mmr(l,m,cld_mmr)
     722             : 
     723             :                                  ! determine change in aerosol mass
     724    87798058 :                                  delmmr = 0._r8
     725    87798058 :                                  delnum = 0._r8
     726    87798058 :                                  if (trim(spectype)=='dust') then
     727    87798058 :                                     if (dst_num>0._r8) then
     728    87798058 :                                        delmmr = (odst_num / dst_num) * icldm(i,k) * amb_mmr(i,k) * wght
     729    87798058 :                                        delnum = (odst_num * icldm(i,k)) /rho(i,k)/per_cm3
     730             :                                     endif
     731           0 :                                  elseif (trim(spectype)=='sulfate') then
     732           0 :                                     if (so4_num>0._r8) then
     733           0 :                                        delmmr = (oso4_num / so4_num) * icldm(i,k) * amb_mmr(i,k) * wght
     734           0 :                                        delnum = (oso4_num * icldm(i,k)) /rho(i,k)/per_cm3
     735             :                                     endif
     736             :                                  endif
     737             : 
     738    87798058 :                                  if (idxtmp>0) then
     739             :                                     ! constituent tendency (for transported species)
     740    87798058 :                                     ptend%q(i,k,idxtmp) = -delmmr/dtime
     741             :                                  else
     742             :                                     ! apply change of mass to not-transported species
     743           0 :                                     amb_mmr(i,k) = amb_mmr(i,k) - delmmr
     744             :                                  endif
     745    87798058 :                                  cld_mmr(i,k) = cld_mmr(i,k) + delmmr
     746             : 
     747    87798058 :                                  delmmr_sum = delmmr_sum + delmmr
     748    87798058 :                                  delnum_sum = delnum_sum + delnum
     749             :                               end if
     750             :                            end if
     751             :                         end do
     752             : 
     753    87798058 :                         idxtmp = aer_cnst_idx(m,0)
     754             : 
     755             :                         ! update aerosol state bin and tendency for grid box i,k
     756    87798058 :                         call aero_state%update_bin( m,i,k, delmmr_sum, delnum_sum, idxtmp, dtime, ptend%q )
     757             : 
     758             :                      end if
     759             : 
     760             :                   end if
     761             :                end do
     762             : 
     763             :             end if
     764             : 
     765             : 
     766             :             ! Liu&Penner does not generate enough nucleation in the polar winter
     767             :             ! stratosphere, which affects surface area density, dehydration and
     768             :             ! ozone chemistry. Part of this is that there are a larger number of
     769             :             ! particles in the accumulation mode than in the Aitken mode. In volcanic
     770             :             ! periods, the coarse mode may also be important. As a short
     771             :             ! term work around, include the accumulation and coarse mode particles
     772             :             ! and assume a larger fraction of the sulfates nucleate in the polar
     773             :             ! stratosphere.
     774             :             !
     775             :             ! Do not include the tropopause level, as stratospheric aerosols
     776             :             ! only exist above the tropopause level.
     777             :             !
     778             :             ! NOTE: This may still not represent the proper particles that
     779             :             ! participate in nucleation, because it doesn't include STS and NAT
     780             :             ! particles. It may not represent the proper saturation threshold for
     781             :             ! nucleation, and wsubi from CLUBB is probably not representative of
     782             :             ! wave driven varaibility in the polar stratosphere.
     783    87798058 :             if (nucleate_ice_use_troplev .and. clim_modal_carma) then
     784    87798058 :                if ((k < troplev(i)) .and. (nucleate_ice_strat > 0._r8) .and. (oso4_num > 0._r8)) then
     785         490 :                   dso4_num = max(0._r8, (nucleate_ice_strat*so4_num_st_cr_tot - oso4_num) * 1e6_r8 / rho(i,k))
     786         490 :                   naai(i,k) = naai(i,k) + dso4_num
     787         490 :                   nihf(i,k) = nihf(i,k) + dso4_num
     788             :                endif
     789             :             else
     790             :                ! This maintains backwards compatibility with the previous version.
     791           0 :                if (pmid(i,k) <= 12500._r8 .and. pmid(i,k) > 100._r8 .and. abs(state%lat(i)) >= 60._r8 * pi / 180._r8) then
     792           0 :                   ramp = 1._r8 - min(1._r8, max(0._r8, (pmid(i,k) - 10000._r8) / 2500._r8))
     793             : 
     794           0 :                   if (oso4_num > 0._r8) then
     795           0 :                      dso4_num = (max(oso4_num, ramp * nucleate_ice_strat * so4_num) - oso4_num) * 1e6_r8 / rho(i,k)
     796           0 :                      naai(i,k) = naai(i,k) + dso4_num
     797           0 :                      nihf(i,k) = nihf(i,k) + dso4_num
     798             :                   end if
     799             :                end if
     800             :             end if
     801             : 
     802    87798058 :             if (cam_physpkg_is("cam7")) then
     803             :                !Updates for pumas v1.21+
     804             : 
     805           0 :                naai_hom(i,k) = nihf(i,k)/dtime
     806           0 :                naai(i,k)= naai(i,k)/dtime
     807             : 
     808             :                ! output activated ice (convert from #/kg -> #/m3/s)
     809           0 :                nihf(i,k)     = nihf(i,k) *rho(i,k)/dtime
     810           0 :                niimm(i,k)    = niimm(i,k)*rho(i,k)/dtime
     811           0 :                nidep(i,k)    = nidep(i,k)*rho(i,k)/dtime
     812           0 :                nimey(i,k)    = nimey(i,k)*rho(i,k)/dtime
     813             : 
     814           0 :                if (use_preexisting_ice) then
     815           0 :                   INnso4(i,k) =so4_num*1e6_r8/dtime  ! (convert from #/cm3 -> #/m3/s)
     816           0 :                   INnbc(i,k)  =soot_num*1e6_r8/dtime
     817           0 :                   INndust(i,k)=dst_num*1e6_r8/dtime
     818           0 :                   INondust(i,k)=odst_num*1e6_r8/dtime
     819           0 :                   INFreIN(i,k)=1.0_r8          ! 1,ice nucleation occur
     820           0 :                   INhet(i,k) = (niimm(i,k) + nidep(i,k))   ! #/m3/s, nimey not in cirrus
     821           0 :                   INhom(i,k) = nihf(i,k)                 ! #/m3/s
     822           0 :                   if (INhom(i,k).gt.1e3_r8)   then ! > 1/L
     823           0 :                      INFrehom(i,k)=1.0_r8       ! 1, hom freezing occur
     824             :                   endif
     825             : 
     826             :                   ! exclude  no ice nucleaton
     827           0 :                   if ((INFrehom(i,k) < 0.5_r8) .and. (INhet(i,k) < 1.0_r8))   then
     828           0 :                      INnso4(i,k) =0.0_r8
     829           0 :                      INnbc(i,k)  =0.0_r8
     830           0 :                      INndust(i,k)=0.0_r8
     831           0 :                      INondust(i,k)=0.0_r8
     832           0 :                      INFreIN(i,k)=0.0_r8
     833           0 :                      INhet(i,k) = 0.0_r8
     834           0 :                      INhom(i,k) = 0.0_r8
     835           0 :                      INFrehom(i,k)=0.0_r8
     836           0 :                      wice(i,k) = 0.0_r8
     837           0 :                      weff(i,k) = 0.0_r8
     838           0 :                      fhom(i,k) = 0.0_r8
     839             :                   endif
     840             :                endif
     841             : 
     842             :             else ! Not cam7
     843             : 
     844    87798058 :                naai_hom(i,k) = nihf(i,k)
     845             : 
     846             :                ! output activated ice (convert from #/kg -> #/m3/s)
     847    87798058 :                nihf(i,k)     = nihf(i,k) *rho(i,k)
     848    87798058 :                niimm(i,k)    = niimm(i,k)*rho(i,k)
     849    87798058 :                nidep(i,k)    = nidep(i,k)*rho(i,k)
     850    87798058 :                nimey(i,k)    = nimey(i,k)*rho(i,k)
     851             : 
     852    87798058 :                if (use_preexisting_ice) then
     853    87798058 :                   INnso4(i,k) =so4_num*1e6_r8 ! (convert from #/cm3 -> #/m3/s)
     854    87798058 :                   INnbc(i,k)  =soot_num*1e6_r8
     855    87798058 :                   INndust(i,k)=dst_num*1e6_r8
     856    87798058 :                   INondust(i,k)=odst_num*1e6_r8
     857    87798058 :                   INFreIN(i,k)=1.0_r8          ! 1,ice nucleation occur
     858    87798058 :                   INhet(i,k) = (niimm(i,k) + nidep(i,k))   ! #/m3, nimey not in cirrus
     859    87798058 :                   INhom(i,k) = nihf(i,k)                 ! #/m3
     860    87798058 :                   if (INhom(i,k).gt.1e3_r8)   then ! > 1/L
     861       10034 :                      INFrehom(i,k)=1.0_r8       ! 1, hom freezing occur
     862             :                   endif
     863             : 
     864             :                   ! exclude  no ice nucleaton
     865    87798058 :                   if ((INFrehom(i,k) < 0.5_r8) .and. (INhet(i,k) < 1.0_r8))   then
     866    86192534 :                      INnso4(i,k) =0.0_r8
     867    86192534 :                      INnbc(i,k)  =0.0_r8
     868    86192534 :                      INndust(i,k)=0.0_r8
     869    86192534 :                      INondust(i,k)=0.0_r8
     870    86192534 :                      INFreIN(i,k)=0.0_r8
     871    86192534 :                      INhet(i,k) = 0.0_r8
     872    86192534 :                      INhom(i,k) = 0.0_r8
     873    86192534 :                      INFrehom(i,k)=0.0_r8
     874    86192534 :                      wice(i,k) = 0.0_r8
     875    86192534 :                      weff(i,k) = 0.0_r8
     876    86192534 :                      fhom(i,k) = 0.0_r8
     877             :                   endif
     878             :                end if
     879             : 
     880             :             end if ! cam7
     881             :          end if freezing
     882             :       end do iloop
     883             :    end do kloop
     884             : 
     885      241920 :    if (.not. clim_modal_carma) then
     886           0 :       deallocate( &
     887             :            naer2, &
     888           0 :            maerosol)
     889             :    end if
     890             : 
     891      241920 :    if (cam_physpkg_is("cam7")) then
     892             :       ! Updates for PUMAS v1.21+
     893           0 :       call outfld('NIHFTEN',   nihf, pcols, lchnk)
     894           0 :       call outfld('NIIMMTEN', niimm, pcols, lchnk)
     895           0 :       call outfld('NIDEPTEN', nidep, pcols, lchnk)
     896           0 :       call outfld('NIMEYTEN', nimey, pcols, lchnk)
     897             :    else
     898      241920 :       call outfld('NIHF',   nihf, pcols, lchnk)
     899      241920 :       call outfld('NIIMM', niimm, pcols, lchnk)
     900      241920 :       call outfld('NIDEP', nidep, pcols, lchnk)
     901      241920 :       call outfld('NIMEY', nimey, pcols, lchnk)
     902             :    end if
     903      241920 :    call outfld('NIREGM', regm, pcols, lchnk)
     904      241920 :    call outfld('NISUBGRID', subgrid, pcols, lchnk)
     905      241920 :    call outfld('NITROP_PD', trop_pd, pcols, lchnk)
     906             : 
     907      241920 :    if (use_preexisting_ice) then
     908      241920 :       call outfld( 'fhom' , fhom, pcols, lchnk)
     909      241920 :       call outfld( 'WICE' , wice, pcols, lchnk)
     910      241920 :       call outfld( 'WEFF' , weff, pcols, lchnk)
     911      241920 :       if (cam_physpkg_is("cam7")) then
     912             :          ! Updates for PUMAS v1.21+
     913           0 :          call outfld('INnso4TEN',INnso4 , pcols,lchnk)
     914           0 :          call outfld('INnbcTEN',INnbc  , pcols,lchnk)
     915           0 :          call outfld('INndustTEN',INndust, pcols,lchnk)
     916           0 :          call outfld('INondustTEN',INondust, pcols,lchnk)
     917           0 :          call outfld('INhetTEN',INhet  , pcols,lchnk)
     918           0 :          call outfld('INhomTEN',INhom  , pcols,lchnk)
     919             :       else
     920      241920 :          call outfld('INnso4  ',INnso4 , pcols,lchnk)
     921      241920 :          call outfld('INnbc   ',INnbc  , pcols,lchnk)
     922      241920 :          call outfld('INndust ',INndust, pcols,lchnk)
     923      241920 :          call outfld('INondust ',INondust, pcols,lchnk)
     924      241920 :          call outfld('INhet   ',INhet  , pcols,lchnk)
     925      241920 :          call outfld('INhom   ',INhom  , pcols,lchnk)
     926             :       end if
     927      241920 :       call outfld('INFrehom',INFrehom,pcols,lchnk)
     928      241920 :       call outfld('INFreIN ',INFreIN, pcols,lchnk)
     929             :    end if
     930             : 
     931      241920 :    if (allocated(size_wght)) then
     932      241920 :       deallocate(size_wght)
     933             :    end if
     934      241920 :    if (allocated(amb_num_bins)) then
     935      241920 :       deallocate(amb_num_bins)
     936             :    end if
     937             : 
     938      483840 : end subroutine nucleate_ice_cam_calc
     939             : 
     940             : !================================================================================================
     941             : 
     942             : end module nucleate_ice_cam

Generated by: LCOV version 1.14