LCOV - code coverage report
Current view: | top level - physics/cam - hetfrz_classnuc_cam.F90 (source / functions) |  |  | Hit | Total | Coverage
Test: | coverage.info |  | Lines: | 21 | 290 | 7.2 %
Date: | 2024-12-17 17:57:11 |  | Functions: | 2 | 5 | 40.0 %

              Line data    Source code
           1             : module hetfrz_classnuc_cam
           2             :
           3             : !---------------------------------------------------------------------------------
           4             : !
           5             : !  CAM Interfaces for hetfrz_classnuc module.
           6             : !
           7             : !---------------------------------------------------------------------------------
           8             :
           9             : use shr_kind_mod,   only: r8=>shr_kind_r8
          10             : use spmd_utils,     only: masterproc
          11             : use ppgrid,         only: pcols, pver
          12             : use physconst,      only: rair, cpair, rh2o, rhoh2o, mwh2o, tmelt, pi
          13             : use constituents,   only: cnst_get_ind
          14             : use physics_types,  only: physics_state
          15             : use physics_buffer, only: physics_buffer_desc, pbuf_get_index, pbuf_old_tim_idx, pbuf_get_field
          16             : use phys_control,   only: use_hetfrz_classnuc
          17             : use physics_buffer, only: pbuf_add_field, dtype_r8, pbuf_old_tim_idx, &
          18             :                           pbuf_get_index, pbuf_get_field
          19             : use cam_history,    only: addfld, add_default, outfld, fieldname_len
          20             : use ref_pres,       only: top_lev => trop_cloud_top_lev
          21             : use wv_saturation,  only: svp_water_vect, svp_ice_vect
          22             : use cam_logfile,    only: iulog
          23             : use error_messages, only: handle_errmsg, alloc_err
          24             : use cam_abortutils, only: endrun
          25             : use string_utils,   only: int2str
          26             : use hetfrz_classnuc,only: hetfrz_classnuc_init, hetfrz_classnuc_calc
          27             :
          28             : use aerosol_properties_mod, only: aerosol_properties, aero_name_len
          29             : use aerosol_state_mod, only: aerosol_state
          30             :
          31             : implicit none
          32             : private
          33             : save
          34             :
          35             : public :: &
          36             :    hetfrz_classnuc_cam_readnl,   &
          37             :    hetfrz_classnuc_cam_register, &
          38             :    hetfrz_classnuc_cam_init,     &
          39             :    hetfrz_classnuc_cam_calc
          40             :
          41             : ! Namelist variables
          42             : logical :: hist_hetfrz_classnuc = .false.
          43             : real(r8) :: hetfrz_bc_scalfac = -huge(1._r8) ! scaling factor for BC
          44             : real(r8) :: hetfrz_dust_scalfac = -huge(1._r8) ! scaling factor for dust
          45             :
          46             : ! Vars set via init method.
          47             : real(r8) :: mincld      ! minimum allowed cloud fraction
          48             :
          49             : ! constituent indices
          50             : integer :: &
          51             :    cldliq_idx = -1, &
          52             :    cldice_idx = -1, &
          53             :    numliq_idx = -1, &
          54             :    numice_idx = -1
          55             :
          56             : ! pbuf indices for fields provided by heterogeneous freezing
          57             : integer :: &
          58             :    frzimm_idx = -1, &
          59             :    frzcnt_idx = -1, &
          60             :    frzdep_idx = -1
          61             :
          62             : ! pbuf indices for fields needed by heterogeneous freezing
          63             : integer :: &
          64             :    ast_idx = -1
          65             :
          66             : type index_t
          67             :    integer :: bin_ndx
          68             :    integer :: spc_ndx
          69             : end type index_t
          70             :
          71             : type(index_t),allocatable :: indices(:)
          72             : character(len=16),allocatable :: types(:)
          73             : character(len=fieldname_len),allocatable :: tot_dens_hnames(:)
          74             : character(len=fieldname_len),allocatable :: cld_dens_hnames(:)
          75             : character(len=fieldname_len),allocatable :: amb_dens_hnames(:)
          76             : character(len=fieldname_len),allocatable :: coated_dens_hnames(:)
          77             : character(len=fieldname_len),allocatable :: uncoated_dens_hnames(:)
          78             : character(len=fieldname_len),allocatable :: cldfn_dens_hnames(:)
          79             : character(len=fieldname_len),allocatable :: coated_frac_hnames(:)
          80             : character(len=fieldname_len),allocatable :: radius_hnames(:)
          81             : character(len=fieldname_len),allocatable :: wactfac_hnames(:)
          82             :
          83             : integer :: tot_num_bins = 0
          84             :
          85             : !===============================================================================
          86             : contains
          87             : !===============================================================================
          88             :
          89        1536 : subroutine hetfrz_classnuc_cam_readnl(nlfile)
          90             :
          91             :   use namelist_utils,  only: find_group_name
          92             :   use units,           only: getunit, freeunit
          93             :   use spmd_utils,      only: mpicom, mstrid=>masterprocid, mpi_logical, mpi_real8, mpi_success
          94             :
          95             :   character(len=*), intent(in) :: nlfile  ! filepath for file containing namelist input
          96             :
          97             :   ! Local variables
          98             :   integer :: unitn, ierr
          99             :   character(len=*), parameter :: subname = 'hetfrz_classnuc_cam_readnl'
         100             :
         101             :   namelist /hetfrz_classnuc_nl/ hist_hetfrz_classnuc, hetfrz_bc_scalfac, hetfrz_dust_scalfac
         102             :
         103             :   !-----------------------------------------------------------------------------
         104             :
         105        1536 :   if (masterproc) then
         106           2 :      unitn = getunit()
         107           2 :      open( unitn, file=trim(nlfile), status='old' )
         108           2 :      call find_group_name(unitn, 'hetfrz_classnuc_nl', status=ierr)
         109           2 :      if (ierr == 0) then
         110           0 :         read(unitn, hetfrz_classnuc_nl, iostat=ierr)
         111           0 :         if (ierr /= 0) then
         112           0 :            call endrun(subname // ':: ERROR reading namelist')
         113             :         end if
         114             :      end if
         115           2 :      close(unitn)
         116           2 :      call freeunit(unitn)
         117             :   end if
         118             :
         119             :   ! Broadcast namelist variables
         120        1536 :   call mpi_bcast(hist_hetfrz_classnuc, 1, mpi_logical, mstrid, mpicom, ierr)
         121        1536 :   if (ierr /= mpi_success) call endrun(subname//" mpi_bcast: hist_hetfrz_classnuc")
         122        1536 :   call mpi_bcast(hetfrz_bc_scalfac, 1, mpi_real8, mstrid, mpicom, ierr)
         123        1536 :   if (ierr /= mpi_success) call endrun(subname//" mpi_bcast: hetfrz_bc_scalfac")
         124        1536 :   call mpi_bcast(hetfrz_dust_scalfac, 1, mpi_real8, mstrid, mpicom, ierr)
         125        1536 :   if (ierr /= mpi_success) call endrun(subname//" mpi_bcast: hetfrz_dust_scalfac")
         126             :
         127        1536 :   if (masterproc) then
         128           2 :      write(iulog,*) subname,': hist_hetfrz_classnuc = ',hist_hetfrz_classnuc
         129           2 :      write(iulog,*) subname,': hetfrz_bc_scalfac = ',hetfrz_bc_scalfac
         130           2 :      write(iulog,*) subname,': hetfrz_dust_scalfac = ',hetfrz_dust_scalfac
         131             :   end if
         132             :
         133        1536 : end subroutine hetfrz_classnuc_cam_readnl
         134             :
         135             : !================================================================================================
         136             :
         137        1536 : subroutine hetfrz_classnuc_cam_register()
         138             :
         139        1536 :    if (.not. use_hetfrz_classnuc) return
         140             :
         141             :    ! pbuf fields provided by hetfrz_classnuc
         142           0 :    call pbuf_add_field('FRZIMM', 'physpkg', dtype_r8, (/pcols,pver/), frzimm_idx)
         143           0 :    call pbuf_add_field('FRZCNT', 'physpkg', dtype_r8, (/pcols,pver/), frzcnt_idx)
         144           0 :    call pbuf_add_field('FRZDEP', 'physpkg', dtype_r8, (/pcols,pver/), frzdep_idx)
         145             :
         146             : end subroutine hetfrz_classnuc_cam_register
         147             :
         148             : !================================================================================================
         149             :
         150           0 : subroutine hetfrz_classnuc_cam_init(mincld_in, aero_props)
         151             :
         152             :    real(r8), intent(in) :: mincld_in
         153             :    class(aerosol_properties), intent(in) :: aero_props
         154             :
         155             :    ! local variables
         156             :    integer :: istat, cnt, ibin, ispc
         157             :    character(len=42) :: tmpstr
         158             :    character(len=aero_name_len) :: species_type
         159             :    character(len=*), parameter :: routine = 'hetfrz_classnuc_cam_init'
         160             :
         161             :    !--------------------------------------------------------------------------------------------
         162             :
         163           0 :    if (.not. use_hetfrz_classnuc) return
         164             :
         165           0 :    cnt = 0
         166           0 :    do ibin = 1, aero_props%nbins()
         167           0 :       do ispc = 1, aero_props%nspecies(ibin)
         168           0 :          if (aero_props%hetfrz_species(ibin,ispc)) then
         169           0 :             cnt = cnt+1
         170             :          end if
         171             :       end do
         172             :    end do
         173             :
         174           0 :    tot_num_bins = cnt
         175             :
         176           0 :    allocate(indices(tot_num_bins), stat=istat)
         177           0 :    call alloc_err(istat, routine, 'indices', tot_num_bins)
         178           0 :    allocate(types(tot_num_bins), stat=istat)
         179           0 :    call alloc_err(istat, routine, 'types', tot_num_bins)
         180             :
         181           0 :    allocate(tot_dens_hnames(tot_num_bins), stat=istat)
         182           0 :    call alloc_err(istat, routine, 'tot_dens_hnames', tot_num_bins)
         183             :
         184           0 :    allocate(cld_dens_hnames(tot_num_bins), stat=istat)
         185           0 :    call alloc_err(istat, routine, 'cld_dens_hnames', tot_num_bins)
         186             :
         187           0 :    allocate(cldfn_dens_hnames(tot_num_bins), stat=istat)
         188           0 :    call alloc_err(istat, routine, 'cldfn_dens_hnames', tot_num_bins)
         189             :
         190           0 :    allocate(amb_dens_hnames(tot_num_bins), stat=istat)
         191           0 :    call alloc_err(istat, routine, 'amb_dens_hnames', tot_num_bins)
         192             :
         193           0 :    allocate(coated_dens_hnames(tot_num_bins), stat=istat)
         194           0 :    call alloc_err(istat, routine, 'coated_dens_hnames', tot_num_bins)
         195             :
         196           0 :    allocate(uncoated_dens_hnames(tot_num_bins), stat=istat)
         197           0 :    call alloc_err(istat, routine, 'uncoated_dens_hnames', tot_num_bins)
         198             :
         199           0 :    allocate(coated_frac_hnames(tot_num_bins), stat=istat)
         200           0 :    call alloc_err(istat, routine, 'coated_frac_hnames', tot_num_bins)
         201             :
         202           0 :    allocate(radius_hnames(tot_num_bins), stat=istat)
         203           0 :    call alloc_err(istat, routine, 'radius_hnames', tot_num_bins)
         204             :
         205           0 :    allocate(wactfac_hnames(tot_num_bins), stat=istat)
         206           0 :    call alloc_err(istat, routine, 'wactfac_hnames', tot_num_bins)
         207             :
         208           0 :    cnt = 0
         209           0 :    do ibin = 1, aero_props%nbins()
         210             :
         211           0 :       do ispc = 1, aero_props%nspecies(ibin)
         212           0 :          if (aero_props%hetfrz_species(ibin,ispc)) then
         213           0 :             call aero_props%species_type(ibin, ispc, species_type)
         214           0 :             cnt = cnt+1
         215           0 :             indices(cnt)%bin_ndx = ibin
         216           0 :             indices(cnt)%spc_ndx = ispc
         217           0 :             types(cnt) = trim(species_type)
         218           0 :             tmpstr = trim(species_type)//trim(int2str(ibin))
         219             :
         220           0 :             cldfn_dens_hnames(cnt) = trim(tmpstr)//'_cld_fn'
         221           0 :             tot_dens_hnames(cnt) = trim(tmpstr)//'_tot_num'
         222           0 :             cld_dens_hnames(cnt) = trim(tmpstr)//'_cld_num'
         223           0 :             amb_dens_hnames(cnt) = trim(tmpstr)//'_amb_num'
         224           0 :             coated_dens_hnames(cnt) = trim(tmpstr)//'_coated'
         225           0 :             uncoated_dens_hnames(cnt) = trim(tmpstr)//'_uncoated'
         226           0 :             coated_frac_hnames(cnt) = trim(tmpstr)//'_coated_frac'
         227           0 :             radius_hnames(cnt) = trim(tmpstr)//'_radius'
         228           0 :             wactfac_hnames(cnt) = trim(tmpstr)//'_wactfac'
         229             :
         230           0 :             call addfld(tot_dens_hnames(cnt),(/ 'lev' /), 'A', '#/cm3', &
         231           0 :                  'total '//trim(tmpstr)//' number density' )
         232           0 :             call addfld(cld_dens_hnames(cnt),(/ 'lev' /), 'A', '#/cm3', &
         233           0 :                  'cloud borne '//trim(tmpstr)//' number density' )
         234           0 :             call addfld(cldfn_dens_hnames(cnt),(/ 'lev' /), 'A', '#/cm3', &
         235           0 :                  'cloud borne '//trim(tmpstr)//' number density derived from fn' )
         236           0 :             call addfld(amb_dens_hnames(cnt),(/ 'lev' /), 'A', '#/cm3', &
         237           0 :                  'ambient '//trim(tmpstr)//' number density' )
         238           0 :             call addfld(coated_dens_hnames(cnt),(/ 'lev' /), 'A', '#/cm3', &
         239           0 :                  'coated '//trim(tmpstr)//' number density' )
         240           0 :             call addfld(uncoated_dens_hnames(cnt),(/ 'lev' /), 'A', '#/cm3', &
         241           0 :                  'uncoated '//trim(tmpstr)//' number density' )
         242           0 :             call addfld(coated_frac_hnames(cnt),(/ 'lev' /), 'A', '#/cm3', &
         243           0 :                  'coated '//trim(tmpstr)//' fraction' )
         244           0 :             call addfld(radius_hnames(cnt),(/ 'lev' /), 'A', 'm', &
         245           0 :                  'ambient '//trim(tmpstr)//' radius' )
         246           0 :             call addfld(wactfac_hnames(cnt),(/ 'lev' /), 'A', ' ', &
         247           0 :                  trim(tmpstr)//' water activity mass factor' )
         248             :
         249             :          end if
         250             :       end do
         251             :
         252             :    end do
         253             :
         254           0 :    mincld = mincld_in
         255             :
         256           0 :    call cnst_get_ind('CLDLIQ', cldliq_idx)
         257           0 :    call cnst_get_ind('CLDICE', cldice_idx)
         258           0 :    call cnst_get_ind('NUMLIQ', numliq_idx)
         259           0 :    call cnst_get_ind('NUMICE', numice_idx)
         260             :
         261             :    ! pbuf fields used by hetfrz_classnuc
         262           0 :    ast_idx      = pbuf_get_index('AST')
         263             :
         264           0 :    call addfld('FRZIMM', (/ 'lev' /), 'A', ' ', 'immersion  freezing')
         265           0 :    call addfld('FRZCNT', (/ 'lev' /), 'A', ' ', 'contact    freezing')
         266           0 :    call addfld('FRZDEP', (/ 'lev' /), 'A', ' ', 'deposition freezing')
         267           0 :    call addfld('FREQIMM', (/ 'lev' /), 'A', 'fraction', 'Fractional occurance of immersion  freezing')
         268           0 :    call addfld('FREQCNT', (/ 'lev' /), 'A', 'fraction', 'Fractional occurance of contact    freezing')
         269           0 :    call addfld('FREQDEP', (/ 'lev' /), 'A', 'fraction', 'Fractional occurance of deposition freezing')
         270           0 :    call addfld('FREQMIX', (/ 'lev' /), 'A', 'fraction', 'Fractional occurance of mixed-phase clouds' )
         271             :
         272           0 :    call addfld('DSTFREZIMM', (/ 'lev' /), 'A', 'm-3s-1', 'dust immersion  freezing rate')
         273           0 :    call addfld('DSTFREZCNT', (/ 'lev' /), 'A', 'm-3s-1', 'dust contact    freezing rate')
         274           0 :    call addfld('DSTFREZDEP', (/ 'lev' /), 'A', 'm-3s-1', 'dust deposition freezing rate')
         275             :
         276           0 :    call addfld('BCFREZIMM', (/ 'lev' /), 'A', 'm-3s-1', 'bc immersion  freezing rate')
         277           0 :    call addfld('BCFREZCNT', (/ 'lev' /), 'A', 'm-3s-1', 'bc contact    freezing rate')
         278           0 :    call addfld('BCFREZDEP', (/ 'lev' /), 'A', 'm-3s-1', 'bc deposition freezing rate')
         279             :
         280             :    call addfld('NIMIX_IMM', (/ 'lev' /), 'A', '#/m3', &
         281           0 :                'Activated Ice Number Concentration due to het immersion freezing in Mixed Clouds')
         282             :    call addfld('NIMIX_CNT', (/ 'lev' /), 'A', '#/m3', &
         283           0 :                'Activated Ice Number Concentration due to het contact freezing in Mixed Clouds')
         284             :    call addfld('NIMIX_DEP', (/ 'lev' /), 'A', '#/m3', &
         285           0 :                'Activated Ice Number Concentration due to het deposition freezing in Mixed Clouds')
         286             :
         287             :    call addfld('DSTNIDEP', (/ 'lev' /), 'A', '#/m3', &
         288           0 :                'Activated Ice Number Concentration due to dst dep freezing in Mixed Clouds')
         289             :    call addfld('DSTNICNT', (/ 'lev' /), 'A', '#/m3', &
         290           0 :                'Activated Ice Number Concentration due to dst cnt freezing in Mixed Clouds')
         291             :    call addfld('DSTNIIMM', (/ 'lev' /), 'A', '#/m3', &
         292           0 :                'Activated Ice Number Concentration due to dst imm freezing in Mixed Clouds')
         293             :
         294             :    call addfld('BCNIDEP', (/ 'lev' /), 'A', '#/m3', &
         295           0 :                'Activated Ice Number Concentration due to bc dep freezing in Mixed Clouds')
         296             :    call addfld('BCNICNT', (/ 'lev' /), 'A', '#/m3', &
         297           0 :                'Activated Ice Number Concentration due to bc cnt freezing in Mixed Clouds')
         298             :    call addfld('BCNIIMM', (/ 'lev' /), 'A', '#/m3', &
         299           0 :                'Activated Ice Number Concentration due to bc imm freezing in Mixed Clouds')
         300             :
         301             :    call addfld('NUMICE10s', (/ 'lev' /), 'A', '#/m3', &
         302           0 :                'Ice Number Concentration due to het freezing in Mixed Clouds during 10-s period')
         303             :    call addfld('NUMIMM10sDST', (/ 'lev' /), 'A', '#/m3', &
         304           0 :                'Ice Number Concentration due to imm freezing by dst in Mixed Clouds during 10-s period')
         305             :    call addfld('NUMIMM10sBC', (/ 'lev' /), 'A', '#/m3', &
         306           0 :                'Ice Number Concentration due to imm freezing by bc in Mixed Clouds during 10-s period')
         307             :
         308           0 :    if (hist_hetfrz_classnuc) then
         309             :
         310           0 :       call add_default('FREQIMM', 1, ' ')
         311           0 :       call add_default('FREQCNT', 1, ' ')
         312           0 :       call add_default('FREQDEP', 1, ' ')
         313           0 :       call add_default('FREQMIX', 1, ' ')
         314             :
         315           0 :       call add_default('DSTFREZIMM', 1, ' ')
         316           0 :       call add_default('DSTFREZCNT', 1, ' ')
         317           0 :       call add_default('DSTFREZDEP', 1, ' ')
         318             :
         319           0 :       call add_default('BCFREZIMM', 1, ' ')
         320           0 :       call add_default('BCFREZCNT', 1, ' ')
         321           0 :       call add_default('BCFREZDEP', 1, ' ')
         322             :
         323           0 :       call add_default('NIMIX_IMM', 1, ' ')
         324           0 :       call add_default('NIMIX_CNT', 1, ' ')
         325           0 :       call add_default('NIMIX_DEP', 1, ' ')
         326             :
         327           0 :       call add_default('DSTNIDEP', 1, ' ')
         328           0 :       call add_default('DSTNICNT', 1, ' ')
         329           0 :       call add_default('DSTNIIMM', 1, ' ')
         330             :
         331           0 :       call add_default('BCNIDEP', 1, ' ')
         332           0 :       call add_default('BCNICNT', 1, ' ')
         333           0 :       call add_default('BCNIIMM', 1, ' ')
         334             :
         335           0 :       call add_default('NUMICE10s', 1, ' ')
         336           0 :       call add_default('NUMIMM10sDST', 1, ' ')
         337           0 :       call add_default('NUMIMM10sBC', 1, ' ')
         338             :
         339             :    end if
         340             :
         341             :    call hetfrz_classnuc_init(rair, cpair, rh2o, rhoh2o, mwh2o, tmelt, iulog, &
         342           0 :         hetfrz_bc_scalfac, hetfrz_dust_scalfac )
         343             :
         344             : end subroutine hetfrz_classnuc_cam_init
         345             :
         346             : !================================================================================================
         347             :
         348           0 : subroutine hetfrz_classnuc_cam_calc(aero_props, aero_state, state, deltatin, factnum, pbuf)
         349             :
         350             :    ! arguments
         351             :    class(aerosol_properties),   intent(in) :: aero_props
         352             :    class(aerosol_state),        intent(in) :: aero_state
         353             :    type(physics_state), target, intent(in) :: state
         354             :    real(r8),                    intent(in) :: deltatin       ! time step (s)
         355             :    real(r8),                    intent(in) :: factnum(:,:,:) ! activation fraction for aerosol number
         356             :    type(physics_buffer_desc),   pointer    :: pbuf(:)
         357             :
         358             :    ! local workspace
         359             :
         360             :    ! outputs shared with the microphysics via the pbuf
         361           0 :    real(r8), pointer :: frzimm(:,:)
         362           0 :    real(r8), pointer :: frzcnt(:,:)
         363           0 :    real(r8), pointer :: frzdep(:,:)
         364             :
         365             :    integer :: itim_old
         366             :    integer :: i, k
         367             :
         368             :    real(r8) :: rho(pcols,pver)          ! air density (kg m-3)
         369             :
         370           0 :    real(r8), pointer :: ast(:,:)
         371             :
         372             :    real(r8) :: lcldm(pcols,pver)
         373             :
         374             :    real(r8) :: esi(pcols), esl(pcols)
         375             :    real(r8) :: con1, r3lx, supersatice
         376             :
         377             :    real(r8) :: qcic
         378             :    real(r8) :: ncic
         379             :
         380             :    real(r8) :: frzbcimm(pcols,pver), frzduimm(pcols,pver)
         381             :    real(r8) :: frzbccnt(pcols,pver), frzducnt(pcols,pver)
         382             :    real(r8) :: frzbcdep(pcols,pver), frzdudep(pcols,pver)
         383             :
         384             :    real(r8) :: freqimm(pcols,pver), freqcnt(pcols,pver), freqdep(pcols,pver), freqmix(pcols,pver)
         385             :
         386             :    real(r8) :: nnuccc_bc(pcols,pver), nnucct_bc(pcols,pver), nnudep_bc(pcols,pver)
         387             :    real(r8) :: nnuccc_dst(pcols,pver), nnucct_dst(pcols,pver), nnudep_dst(pcols,pver)
         388             :    real(r8) :: niimm_bc(pcols,pver), nicnt_bc(pcols,pver), nidep_bc(pcols,pver)
         389             :    real(r8) :: niimm_dst(pcols,pver), nicnt_dst(pcols,pver), nidep_dst(pcols,pver)
         390             :    real(r8) :: numice10s(pcols,pver)
         391             :    real(r8) :: numice10s_imm_dst(pcols,pver)
         392             :    real(r8) :: numice10s_imm_bc(pcols,pver)
         393             :
         394           0 :    real(r8) :: coated(pcols,pver,tot_num_bins)
         395           0 :    real(r8) :: aer_radius(pcols,pver,tot_num_bins)
         396           0 :    real(r8) :: aer_wactfac(pcols,pver,tot_num_bins)
         397             :
         398           0 :    real(r8) :: coated_amb_aer_num(pcols,pver,tot_num_bins)
         399           0 :    real(r8) :: uncoated_amb_aer_num(pcols,pver,tot_num_bins)
         400           0 :    real(r8) :: amb_aer_num(pcols,pver,tot_num_bins)
         401           0 :    real(r8) :: cld_aer_num(pcols,pver,tot_num_bins)
         402           0 :    real(r8) :: tot_aer_num(pcols,pver,tot_num_bins)
         403             :    real(r8) :: fn_cld_aer_num(pcols,pver)
         404           0 :    real(r8) :: fraction_activated(pcols,pver,tot_num_bins)
         405             :
         406             :    character(128) :: errstring   ! Error status
         407             :    !-------------------------------------------------------------------------------
         408             :
         409             :    associate( &
         410             :       lchnk => state%lchnk,             &
         411             :       ncol  => state%ncol,              &
         412             :       t     => state%t,                 &
         413           0 :       qc    => state%q(:pcols,:pver,cldliq_idx), &
         414             :       nc    => state%q(:pcols,:pver,numliq_idx), &
         415             :       pmid  => state%pmid               )
         416             :
         417           0 :    itim_old = pbuf_old_tim_idx()
         418           0 :    call pbuf_get_field(pbuf, ast_idx, ast, start=(/1,1,itim_old/), kount=(/pcols,pver,1/))
         419             :
         420           0 :    rho(:,:) = 0._r8
         421             :
         422           0 :    do k = top_lev, pver
         423           0 :       do i = 1, ncol
         424           0 :          rho(i,k) = pmid(i,k)/(rair*t(i,k))
         425             :       end do
         426             :    end do
         427             :
         428           0 :    do k = top_lev, pver
         429           0 :       do i = 1, ncol
         430           0 :          lcldm(i,k) = max(ast(i,k), mincld)
         431             :       end do
         432             :    end do
         433             :
         434           0 :    do i = 1,tot_num_bins
         435             :
         436           0 :       call aero_state%get_amb_species_numdens( indices(i)%bin_ndx, ncol, pver, types(i), aero_props, rho, amb_aer_num(:,:,i))
         437           0 :       call aero_state%get_cld_species_numdens( indices(i)%bin_ndx, ncol, pver, types(i), aero_props, rho, cld_aer_num(:,:,i))
         438             :
         439           0 :       tot_aer_num(:ncol,:,i) = cld_aer_num(:ncol,:,i) + amb_aer_num(:ncol,:,i)
         440             :
         441           0 :       call outfld(tot_dens_hnames(i), tot_aer_num(:,:,i), pcols, lchnk)
         442           0 :       call outfld(amb_dens_hnames(i), amb_aer_num(:,:,i), pcols, lchnk)
         443           0 :       call outfld(cld_dens_hnames(i), cld_aer_num(:,:,i), pcols, lchnk)
         444             :
         445           0 :       aer_radius(:ncol,:,i) = aero_state%mass_mean_radius( indices(i)%bin_ndx, indices(i)%spc_ndx,  ncol, pver, aero_props, rho )
         446             :
         447           0 :       coated(:ncol,:,i) = aero_state%coated_frac( indices(i)%bin_ndx, types(i), ncol, pver, aero_props, aer_radius(:,:,i) )
         448             :
         449           0 :       call outfld(coated_frac_hnames(i), coated(:,:,i), pcols, lchnk)
         450             :
         451           0 :       coated_amb_aer_num(:ncol,:,i) = amb_aer_num(:ncol,:,i)*coated(:ncol,:,i)
         452           0 :       uncoated_amb_aer_num(:ncol,:,i) = amb_aer_num(:ncol,:,i)*(1._r8-coated(:ncol,:,i))
         453             :
         454           0 :       call outfld(coated_dens_hnames(i), coated_amb_aer_num(:,:,i), pcols, lchnk)
         455           0 :       call outfld(uncoated_dens_hnames(i), uncoated_amb_aer_num(:,:,i), pcols, lchnk)
         456           0 :       call outfld(radius_hnames(i), aer_radius(:ncol,:,i), ncol, lchnk)
         457             :
         458           0 :       call aero_state%watact_mfactor(indices(i)%bin_ndx, types(i), ncol, pver, aero_props, rho, aer_wactfac(:ncol,:,i))
         459           0 :       call outfld(wactfac_hnames(i), aer_wactfac(:,:,i), pcols, lchnk)
         460             :
         461           0 :       fn_cld_aer_num(:ncol,:) = tot_aer_num(:ncol,:,i)*factnum(:ncol,:,indices(i)%bin_ndx)
         462           0 :       call outfld(cldfn_dens_hnames(i), fn_cld_aer_num, pcols, lchnk)
         463             :
         464           0 :       fraction_activated(:ncol,:,i) = factnum(:ncol,:,indices(i)%bin_ndx)
         465             :
         466             :    end do
         467             :
         468             :    ! frzimm, frzcnt, frzdep are the outputs of this parameterization used by the microphysics
         469           0 :    call pbuf_get_field(pbuf, frzimm_idx, frzimm)
         470           0 :    call pbuf_get_field(pbuf, frzcnt_idx, frzcnt)
         471           0 :    call pbuf_get_field(pbuf, frzdep_idx, frzdep)
         472             :
         473           0 :    frzimm(:ncol,:) = 0._r8
         474           0 :    frzcnt(:ncol,:) = 0._r8
         475           0 :    frzdep(:ncol,:) = 0._r8
         476             :
         477           0 :    frzbcimm(:ncol,:) = 0._r8
         478           0 :    frzduimm(:ncol,:) = 0._r8
         479           0 :    frzbccnt(:ncol,:) = 0._r8
         480           0 :    frzducnt(:ncol,:) = 0._r8
         481           0 :    frzbcdep(:ncol,:) = 0._r8
         482           0 :    frzdudep(:ncol,:) = 0._r8
         483             :
         484           0 :    freqimm(:ncol,:) = 0._r8
         485           0 :    freqcnt(:ncol,:) = 0._r8
         486           0 :    freqdep(:ncol,:) = 0._r8
         487           0 :    freqmix(:ncol,:) = 0._r8
         488             :
         489           0 :    numice10s(:ncol,:)         = 0._r8
         490           0 :    numice10s_imm_dst(:ncol,:) = 0._r8
         491           0 :    numice10s_imm_bc(:ncol,:)  = 0._r8
         492             :
         493           0 :    nnuccc_bc(:,:) = 0._r8
         494           0 :    nnucct_bc(:,:) = 0._r8
         495           0 :    nnudep_bc(:,:) = 0._r8
         496             :
         497           0 :    nnuccc_dst(:,:) = 0._r8
         498           0 :    nnucct_dst(:,:) = 0._r8
         499           0 :    nnudep_dst(:,:) = 0._r8
         500             :
         501           0 :    niimm_bc(:,:) = 0._r8
         502           0 :    nicnt_bc(:,:) = 0._r8
         503           0 :    nidep_bc(:,:) = 0._r8
         504             :
         505           0 :    niimm_dst(:,:) = 0._r8
         506           0 :    nicnt_dst(:,:) = 0._r8
         507           0 :    nidep_dst(:,:) = 0._r8
         508             :
         509           0 :    do k = top_lev, pver
         510           0 :       call svp_water_vect(t(1:ncol,k), esl(1:ncol), ncol)
         511           0 :       call svp_ice_vect(t(1:ncol,k), esi(1:ncol), ncol)
         512           0 :       do i = 1, ncol
         513             :
         514           0 :          if (t(i,k) > 235.15_r8 .and. t(i,k) < 269.15_r8) then
         515           0 :             qcic = min(qc(i,k)/lcldm(i,k), 5.e-3_r8)
         516           0 :             ncic = max(nc(i,k)/lcldm(i,k), 0._r8)
         517             :
         518           0 :             con1 = 1._r8/(1.333_r8*pi)**0.333_r8
         519           0 :             r3lx = con1*(rho(i,k)*qcic/(rhoh2o*max(ncic*rho(i,k), 1.0e6_r8)))**0.333_r8 ! in m
         520           0 :             r3lx = max(4.e-6_r8, r3lx)
         521           0 :             supersatice = esl(i)/esi(i)
         522             :
         523             :             call hetfrz_classnuc_calc( tot_num_bins, types, &
         524           0 :                deltatin,  t(i,k),  pmid(i,k),  supersatice,   &
         525             :                fraction_activated(i,k,:),  r3lx,  ncic*rho(i,k)*1.0e-6_r8,  frzbcimm(i,k),  frzduimm(i,k),   &
         526             :                frzbccnt(i,k),  frzducnt(i,k),  frzbcdep(i,k),  frzdudep(i,k),  aer_radius(i,k,:), &
         527             :                aer_wactfac(i,k,:), coated(i,k,:), tot_aer_num(i,k,:),  &
         528             :                uncoated_amb_aer_num(i,k,:), amb_aer_num(i,k,:), &
         529           0 :                cld_aer_num(i,k,:), errstring)
         530             :
         531           0 :             call handle_errmsg(errstring, subname="hetfrz_classnuc_calc")
         532             :
         533           0 :             frzimm(i,k) = frzbcimm(i,k) + frzduimm(i,k)
         534           0 :             frzcnt(i,k) = frzbccnt(i,k) + frzducnt(i,k)
         535           0 :             frzdep(i,k) = frzbcdep(i,k) + frzdudep(i,k)
         536             :
         537           0 :             if (frzimm(i,k) > 0._r8) freqimm(i,k) = 1._r8
         538           0 :             if (frzcnt(i,k) > 0._r8) freqcnt(i,k) = 1._r8
         539           0 :             if (frzdep(i,k) > 0._r8) freqdep(i,k) = 1._r8
         540           0 :             if ((frzimm(i,k) + frzcnt(i,k) + frzdep(i,k)) > 0._r8) freqmix(i,k) = 1._r8
         541             :
         542             :          else
         543           0 :             frzimm(i,k) = 0._r8
         544           0 :             frzcnt(i,k) = 0._r8
         545           0 :             frzdep(i,k) = 0._r8
         546             :          end if
         547             :
         548           0 :          nnuccc_bc(i,k) = frzbcimm(i,k)*1.0e6_r8*ast(i,k)
         549           0 :          nnucct_bc(i,k) = frzbccnt(i,k)*1.0e6_r8*ast(i,k)
         550           0 :          nnudep_bc(i,k) = frzbcdep(i,k)*1.0e6_r8*ast(i,k)
         551             :
         552           0 :          nnuccc_dst(i,k) = frzduimm(i,k)*1.0e6_r8*ast(i,k)
         553           0 :          nnucct_dst(i,k) = frzducnt(i,k)*1.0e6_r8*ast(i,k)
         554           0 :          nnudep_dst(i,k) = frzdudep(i,k)*1.0e6_r8*ast(i,k)
         555             :
         556           0 :          niimm_bc(i,k) = frzbcimm(i,k)*1.0e6_r8*deltatin
         557           0 :          nicnt_bc(i,k) = frzbccnt(i,k)*1.0e6_r8*deltatin
         558           0 :          nidep_bc(i,k) = frzbcdep(i,k)*1.0e6_r8*deltatin
         559             :
         560           0 :          niimm_dst(i,k) = frzduimm(i,k)*1.0e6_r8*deltatin
         561           0 :          nicnt_dst(i,k) = frzducnt(i,k)*1.0e6_r8*deltatin
         562           0 :          nidep_dst(i,k) = frzdudep(i,k)*1.0e6_r8*deltatin
         563             :
         564           0 :          numice10s(i,k) = (frzimm(i,k)+frzcnt(i,k)+frzdep(i,k))*1.0e6_r8*deltatin*(10._r8/deltatin)
         565           0 :          numice10s_imm_dst(i,k) = frzduimm(i,k)*1.0e6_r8*deltatin*(10._r8/deltatin)
         566           0 :          numice10s_imm_bc(i,k) = frzbcimm(i,k)*1.0e6_r8*deltatin*(10._r8/deltatin)
         567             :       end do
         568             :    end do
         569             :
         570           0 :    call outfld('FRZIMM', frzimm, pcols, lchnk)
         571           0 :    call outfld('FRZCNT', frzcnt, pcols, lchnk)
         572           0 :    call outfld('FRZDEP', frzdep, pcols, lchnk)
         573             :
         574           0 :    call outfld('FREQIMM', freqimm, pcols, lchnk)
         575           0 :    call outfld('FREQCNT', freqcnt, pcols, lchnk)
         576           0 :    call outfld('FREQDEP', freqdep, pcols, lchnk)
         577           0 :    call outfld('FREQMIX', freqmix, pcols, lchnk)
         578             :
         579           0 :    call outfld('DSTFREZIMM', nnuccc_dst, pcols, lchnk)
         580           0 :    call outfld('DSTFREZCNT', nnucct_dst, pcols, lchnk)
         581           0 :    call outfld('DSTFREZDEP', nnudep_dst, pcols, lchnk)
         582             :
         583           0 :    call outfld('BCFREZIMM', nnuccc_bc, pcols, lchnk)
         584           0 :    call outfld('BCFREZCNT', nnucct_bc, pcols, lchnk)
         585           0 :    call outfld('BCFREZDEP', nnudep_bc, pcols, lchnk)
         586             :
         587           0 :    call outfld('NIMIX_IMM', niimm_bc+niimm_dst, pcols, lchnk)
         588           0 :    call outfld('NIMIX_CNT', nicnt_bc+nicnt_dst, pcols, lchnk)
         589           0 :    call outfld('NIMIX_DEP', nidep_bc+nidep_dst, pcols, lchnk)
         590             :
         591           0 :    call outfld('DSTNICNT', nicnt_dst, pcols, lchnk)
         592           0 :    call outfld('DSTNIDEP', nidep_dst, pcols, lchnk)
         593           0 :    call outfld('DSTNIIMM', niimm_dst, pcols, lchnk)
         594             :
         595           0 :    call outfld('BCNICNT', nicnt_bc, pcols, lchnk)
         596           0 :    call outfld('BCNIDEP', nidep_bc, pcols, lchnk)
         597           0 :    call outfld('BCNIIMM', niimm_bc, pcols, lchnk)
         598             :
         599           0 :    call outfld('NUMICE10s', numice10s, pcols, lchnk)
         600           0 :    call outfld('NUMIMM10sDST', numice10s_imm_dst, pcols, lchnk)
         601           0 :    call outfld('NUMIMM10sBC', numice10s_imm_bc, pcols, lchnk)
         602             :
         603             :    end associate
         604             :
         605           0 : end subroutine hetfrz_classnuc_cam_calc
         606             :
         607             : !====================================================================================================
         608             :
         609           0 : end module hetfrz_classnuc_cam

Generated by: LCOV version 1.14