LCOV - code coverage report
Current view: | top level - chemistry/modal_aero - sox_cldaero_mod.F90 (source / functions) |  |  | Hit | Total | Coverage
Test: | coverage.info |  | Lines: | 151 | 201 | 75.1 %
Date: | 2024-12-17 22:39:59 |  | Functions: | 4 | 4 | 100.0 %

              Line data    Source code
           1             : !----------------------------------------------------------------------------------
           2             : ! Modal aerosol implementation
           3             : !----------------------------------------------------------------------------------
           4             : module sox_cldaero_mod
           5             :
           6             :   use shr_kind_mod,    only : r8 => shr_kind_r8
           7             :   use cam_abortutils,  only : endrun
           8             :   use ppgrid,          only : pcols, pver
           9             :   use mo_chem_utls,    only : get_spc_ndx
          10             :   use cldaero_mod,     only : cldaero_conc_t, cldaero_allocate, cldaero_deallocate
          11             :   use modal_aero_data, only : ntot_amode, modeptr_accum, lptr_so4_cw_amode, lptr_msa_cw_amode
          12             :   use modal_aero_data, only : numptrcw_amode, lptr_nh4_cw_amode
          13             :   use modal_aero_data, only : cnst_name_cw, specmw_so4_amode
          14             :   use chem_mods,       only : adv_mass
          15             :   use physconst,       only : gravit
          16             :   use phys_control,    only : phys_getopts, cam_chempkg_is
          17             :   use cldaero_mod,     only : cldaero_uptakerate
          18             :   use chem_mods,       only : gas_pcnst
          19             :
          20             :   implicit none
          21             :   private
          22             :
          23             :   public :: sox_cldaero_init
          24             :   public :: sox_cldaero_create_obj
          25             :   public :: sox_cldaero_update
          26             :   public :: sox_cldaero_destroy_obj
          27             :
          28             :   integer :: id_msa, id_h2so4, id_so2, id_h2o2, id_nh3
          29             :
          30             :   real(r8), parameter :: small_value = 1.e-20_r8
          31             :
          32             : contains
          33             :
          34             : !----------------------------------------------------------------------------------
          35             : !----------------------------------------------------------------------------------
          36             :
          37        3072 :   subroutine sox_cldaero_init
          38             :
          39             :     integer :: l, m
          40             :     logical :: history_aerosol      ! Output the MAM aerosol tendencies
          41             :
          42        1536 :     id_msa = get_spc_ndx( 'MSA' )
          43        1536 :     id_h2so4 = get_spc_ndx( 'H2SO4' )
          44        1536 :     id_so2 = get_spc_ndx( 'SO2' )
          45        1536 :     id_h2o2 = get_spc_ndx( 'H2O2' )
          46        1536 :     id_nh3 = get_spc_ndx( 'NH3' )
          47             :
          48        1536 :     if (id_h2so4<1 .or. id_so2<1 .or. id_h2o2<1) then
          49             :       call endrun('sox_cldaero_init:MAM mech does not include necessary species' &
          50           0 :                   //' -- should not invoke sox_cldaero_mod ')
          51             :     endif
          52             :
          53        1536 :     call phys_getopts( history_aerosol_out        = history_aerosol   )
          54             :     !
          55             :     !   add to history
          56             :     !
          57             :
          58        1536 :   end subroutine sox_cldaero_init
          59             :
          60             : !----------------------------------------------------------------------------------
          61             : !----------------------------------------------------------------------------------
          62     1489176 :   function sox_cldaero_create_obj(cldfrc, qcw, lwc, cfact, ncol, loffset) result( conc_obj )
          63             :
          64             :     real(r8), intent(in) :: cldfrc(:,:)
          65             :     real(r8), intent(in) :: qcw(:,:,:)
          66             :     real(r8), intent(in) :: lwc(:,:)
          67             :     real(r8), intent(in) :: cfact(:,:)
          68             :     integer,  intent(in) :: ncol
          69             :     integer,  intent(in) :: loffset
          70             :
          71             :     type(cldaero_conc_t), pointer :: conc_obj
          72             :
          73             :
          74             :     integer :: id_so4_1a, id_so4_2a, id_so4_3a, id_so4_4a, id_so4_5a, id_so4_6a
          75             :     integer :: id_nh4_1a, id_nh4_2a, id_nh4_3a, id_nh4_4a, id_nh4_5a, id_nh4_6a
          76             :     integer :: l,n
          77             :     integer :: i,k
          78             :
          79             :     logical :: mode7
          80             :
          81     1489176 :     mode7 = ntot_amode == 7
          82             :
          83     1489176 :     conc_obj => cldaero_allocate()
          84             :
          85   139982544 :     do k = 1,pver
          86  2314006344 :        do i = 1,ncol
          87  2312517168 :           if( cldfrc(i,k) >0._r8) then
          88   464987798 :              conc_obj%xlwc(i,k) = lwc(i,k) *cfact(i,k) ! cloud water L(water)/L(air)
          89   464987798 :              conc_obj%xlwc(i,k) = conc_obj%xlwc(i,k) / cldfrc(i,k) ! liquid water in the cloudy fraction of cell
          90             :           else
          91  1709036002 :              conc_obj%xlwc(i,k) = 0._r8
          92             :           endif
          93             :        enddo
          94             :     enddo
          95             :
          96  2355876432 :     conc_obj%no3c(:,:) = 0._r8
          97             :
          98     1489176 :     if (mode7) then
          99             : #if ( defined MODAL_AERO_7MODE )
         100             : !put ifdef here so ifort will compile
         101             :        id_so4_1a = lptr_so4_cw_amode(1) - loffset
         102             :        id_so4_2a = lptr_so4_cw_amode(2) - loffset
         103             :        id_so4_3a = lptr_so4_cw_amode(4) - loffset
         104             :        id_so4_4a = lptr_so4_cw_amode(5) - loffset
         105             :        id_so4_5a = lptr_so4_cw_amode(6) - loffset
         106             :        id_so4_6a = lptr_so4_cw_amode(7) - loffset
         107             :
         108             :        id_nh4_1a = lptr_nh4_cw_amode(1) - loffset
         109             :        id_nh4_2a = lptr_nh4_cw_amode(2) - loffset
         110             :        id_nh4_3a = lptr_nh4_cw_amode(4) - loffset
         111             :        id_nh4_4a = lptr_nh4_cw_amode(5) - loffset
         112             :        id_nh4_5a = lptr_nh4_cw_amode(6) - loffset
         113             :        id_nh4_6a = lptr_nh4_cw_amode(7) - loffset
         114             : #endif
         115           0 :        conc_obj%so4c(:ncol,:) &
         116           0 :             = qcw(:ncol,:,id_so4_1a) &
         117           0 :             + qcw(:ncol,:,id_so4_2a) &
         118           0 :             + qcw(:ncol,:,id_so4_3a) &
         119           0 :             + qcw(:ncol,:,id_so4_4a) &
         120           0 :             + qcw(:ncol,:,id_so4_5a) &
         121           0 :             + qcw(:ncol,:,id_so4_6a)
         122             :
         123           0 :        conc_obj%nh4c(:ncol,:) &
         124           0 :             = qcw(:ncol,:,id_nh4_1a) &
         125           0 :             + qcw(:ncol,:,id_nh4_2a) &
         126           0 :             + qcw(:ncol,:,id_nh4_3a) &
         127           0 :             + qcw(:ncol,:,id_nh4_4a) &
         128           0 :             + qcw(:ncol,:,id_nh4_5a) &
         129           0 :             + qcw(:ncol,:,id_nh4_6a)
         130             :     else
         131     1489176 :        id_so4_1a = lptr_so4_cw_amode(1) - loffset
         132     1489176 :        id_so4_2a = lptr_so4_cw_amode(2) - loffset
         133     1489176 :        id_so4_3a = lptr_so4_cw_amode(3) - loffset
         134           0 :        conc_obj%so4c(:ncol,:) &
         135     1489176 :             = qcw(:,:,id_so4_1a) &
         136     1489176 :             + qcw(:,:,id_so4_2a) &
         137  2316984696 :             + qcw(:,:,id_so4_3a)
         138             :
         139             :         ! for 3-mode, so4 is assumed to be nh4hso4
         140             :         ! the partial neutralization of so4 is handled by using a
         141             :         !    -1 charge (instead of -2) in the electro-neutrality equation
         142  2314006344 :        conc_obj%nh4c(:ncol,:) = 0._r8
         143             :
         144             :        ! with 3-mode, assume so4 is nh4hso4, and so half-neutralized
         145     1489176 :        conc_obj%so4_fact = 1._r8
         146             :
         147             :     endif
         148             :
         149     1489176 :   end function sox_cldaero_create_obj
         150             :
         151             : !----------------------------------------------------------------------------------
         152             : ! Update the mixing ratios
         153             : !----------------------------------------------------------------------------------
         154     1489176 :   subroutine sox_cldaero_update( &
         155     4467528 :        ncol, lchnk, loffset, dtime, mbar, pdel, press, tfld, cldnum, cldfrc, cfact, xlwc, &
         156     2978352 :        delso4_hprxn, xh2so4, xso4, xso4_init, nh3g, hno3g, xnh3, xhno3, xnh4c,  xno3c, xmsa, xso2, xh2o2, qcw, qin, &
         157     2978352 :        aqso4, aqh2so4, aqso4_h2o2, aqso4_o3, aqso4_h2o2_3d, aqso4_o3_3d)
         158             :
         159             :     ! args
         160             :
         161             :     integer,  intent(in) :: ncol
         162             :     integer,  intent(in) :: lchnk ! chunk id
         163             :     integer,  intent(in) :: loffset
         164             :
         165             :     real(r8), intent(in) :: dtime ! time step (sec)
         166             :
         167             :     real(r8), intent(in) :: mbar(:,:) ! mean wet atmospheric mass ( amu )
         168             :     real(r8), intent(in) :: pdel(:,:)
         169             :     real(r8), intent(in) :: press(:,:)
         170             :     real(r8), intent(in) :: tfld(:,:)
         171             :
         172             :     real(r8), intent(in) :: cldnum(:,:)
         173             :     real(r8), intent(in) :: cldfrc(:,:)
         174             :     real(r8), intent(in) :: cfact(:,:)
         175             :     real(r8), intent(in) :: xlwc(:,:)
         176             :
         177             :     real(r8), intent(in) :: delso4_hprxn(:,:)
         178             :     real(r8), intent(in) :: xh2so4(:,:)
         179             :     real(r8), intent(in) :: xso4(:,:)
         180             :     real(r8), intent(in) :: xso4_init(:,:)
         181             :     real(r8), intent(in) :: nh3g(:,:)
         182             :     real(r8), intent(in) :: hno3g(:,:)
         183             :     real(r8), intent(in) :: xnh3(:,:)
         184             :     real(r8), intent(in) :: xhno3(:,:)
         185             :     real(r8), intent(in) :: xnh4c(:,:)
         186             :     real(r8), intent(in) :: xmsa(:,:)
         187             :     real(r8), intent(in) :: xso2(:,:)
         188             :     real(r8), intent(in) :: xh2o2(:,:)
         189             :     real(r8), intent(in) :: xno3c(:,:)
         190             :
         191             :     real(r8), intent(inout) :: qcw(:,:,:) ! cloud-borne aerosol (vmr)
         192             :     real(r8), intent(inout) :: qin(:,:,:) ! xported species ( vmr )
         193             :
         194             :     real(r8), intent(out) :: aqso4(:,:)                   ! aqueous phase chemistry
         195             :     real(r8), intent(out) :: aqh2so4(:,:)                 ! aqueous phase chemistry
         196             :     real(r8), intent(out) :: aqso4_h2o2(:)                ! SO4 aqueous phase chemistry due to H2O2 (kg/m2)
         197             :     real(r8), intent(out) :: aqso4_o3(:)                  ! SO4 aqueous phase chemistry due to O3 (kg/m2)
         198             :     real(r8), intent(out), optional :: aqso4_h2o2_3d(:,:)                ! SO4 aqueous phase chemistry due to H2O2 (kg/m2)
         199             :     real(r8), intent(out), optional :: aqso4_o3_3d(:,:)                  ! SO4 aqueous phase chemistry due to O3 (kg/m2)
         200             :
         201             :
         202             :     ! local vars ...
         203             :
         204     2978352 :     real(r8) :: dqdt_aqso4(ncol,pver,gas_pcnst), &
         205     2978352 :          dqdt_aqh2so4(ncol,pver,gas_pcnst), &
         206     2978352 :          dqdt_aqhprxn(ncol,pver), dqdt_aqo3rxn(ncol,pver), &
         207             :          sflx(1:ncol)
         208             :
         209     2978352 :     real(r8) :: faqgain_msa(ntot_amode), faqgain_so4(ntot_amode), qnum_c(ntot_amode)
         210             :
         211             :     real(r8) :: delso4_o3rxn, &
         212             :          dso4dt_aqrxn, dso4dt_hprxn, &
         213             :          dso4dt_gasuptk, dmsadt_gasuptk, &
         214             :          dmsadt_gasuptk_tomsa, dmsadt_gasuptk_toso4, &
         215             :          dqdt_aq, dqdt_wr, dqdt
         216             :
         217             :     real(r8) :: fwetrem, sumf, uptkrate
         218             :     real(r8) :: delnh3, delnh4
         219             :
         220             :     integer :: l, n, m
         221             :     integer :: ntot_msa_c
         222             :
         223             :     integer :: i,k
         224             :     real(r8) :: xl
         225             :
         226             :     ! make sure dqdt is zero initially, for budgets
         227 71735685840 :     dqdt_aqso4(:,:,:) = 0.0_r8
         228 71735685840 :     dqdt_aqh2so4(:,:,:) = 0.0_r8
         229  2314006344 :     dqdt_aqhprxn(:,:) = 0.0_r8
         230  2314006344 :     dqdt_aqo3rxn(:,:) = 0.0_r8
         231             :
         232             :     ! Avoid double counting in-cloud sulfur oxidation when running with
         233             :     ! GEOS-Chem. If running with GEOS-Chem then sulfur oxidation
         234             :     ! is performed internally to GEOS-Chem. Here, we just return to the
         235             :     ! parent routine and thus we do not apply tendencies calculated by MAM.
         236     1489176 :     if ( cam_chempkg_is('geoschem_mam4') ) return
         237             :
         238   139982544 :     lev_loop: do k = 1,pver
         239  2314006344 :        col_loop: do i = 1,ncol
         240  2312517168 :           cloud: if (cldfrc(i,k) >= 1.0e-5_r8) then
         241   454195893 :              xl = xlwc(i,k) ! / cldfrc(i,k)
         242             :
         243   454195893 :              IF (XL .ge. 1.e-8_r8) THEN !! WHEN CLOUD IS PRESENTED
         244             :
         245    83838975 :                 delso4_o3rxn = xso4(i,k) - xso4_init(i,k)
         246             :
         247    83838975 :                 if (id_nh3>0) then
         248           0 :                    delnh3 = nh3g(i,k) - xnh3(i,k)
         249           0 :                    delnh4 = - delnh3
         250             :                 endif
         251             :
         252             :                 !-------------------------------------------------------------------------
         253             :                 ! compute factors for partitioning aerosol mass gains among modes
         254             :                 ! the factors are proportional to the activated particle MR for each
         255             :                 ! mode, which is the MR of cloud drops "associated with" the mode
         256             :                 ! thus we are assuming the cloud drop size is independent of the
         257             :                 ! associated aerosol mode properties (i.e., drops associated with
         258             :                 ! Aitken and coarse sea-salt particles are same size)
         259             :                 !
         260             :                 ! qnum_c(n) = activated particle number MR for mode n (these are just
         261             :                 ! used for partitioning among modes, so don't need to divide by cldfrc)
         262             :
         263   419194875 :                 do n = 1, ntot_amode
         264   335355900 :                    qnum_c(n) = 0.0_r8
         265   335355900 :                    l = numptrcw_amode(n) - loffset
         266   419194875 :                    if (l > 0) qnum_c(n) = max( 0.0_r8, qcw(i,k,l) )
         267             :                 end do
         268             :
         269             :                 ! force qnum_c(n) to be positive for n=modeptr_accum or n=1
         270    83838975 :                 n = modeptr_accum
         271    83838975 :                 if (n <= 0) n = 1
         272    83838975 :                 qnum_c(n) = max( 1.0e-10_r8, qnum_c(n) )
         273             :
         274             :                 ! faqgain_so4(n) = fraction of total so4_c gain going to mode n
         275             :                 ! these are proportional to the activated particle MR for each mode
         276    83838975 :                 sumf = 0.0_r8
         277   419194875 :                 do n = 1, ntot_amode
         278   335355900 :                    faqgain_so4(n) = 0.0_r8
         279   419194875 :                    if (lptr_so4_cw_amode(n) > 0) then
         280   251516925 :                       faqgain_so4(n) = qnum_c(n)
         281   251516925 :                       sumf = sumf + faqgain_so4(n)
         282             :                    end if
         283             :                 end do
         284             :
         285    83838975 :                 if (sumf > 0.0_r8) then
         286   419194875 :                    do n = 1, ntot_amode
         287   419194875 :                       faqgain_so4(n) = faqgain_so4(n) / sumf
         288             :                    end do
         289             :                 end if
         290             :                 ! at this point (sumf <= 0.0) only when all the faqgain_so4 are zero
         291             :
         292             :                 ! faqgain_msa(n) = fraction of total msa_c gain going to mode n
         293             :                 ntot_msa_c = 0
         294             :                 sumf = 0.0_r8
         295   419194875 :                 do n = 1, ntot_amode
         296   335355900 :                    faqgain_msa(n) = 0.0_r8
         297   335355900 :                    if (lptr_msa_cw_amode(n) > 0) then
         298           0 :                       faqgain_msa(n) = qnum_c(n)
         299           0 :                       ntot_msa_c = ntot_msa_c + 1
         300             :                    end if
         301   419194875 :                    sumf = sumf + faqgain_msa(n)
         302             :                 end do
         303             :
         304    83838975 :                 if (sumf > 0.0_r8) then
         305           0 :                    do n = 1, ntot_amode
         306           0 :                       faqgain_msa(n) = faqgain_msa(n) / sumf
         307             :                    end do
         308             :                 end if
         309             :                 ! at this point (sumf <= 0.0) only when all the faqgain_msa are zero
         310             :
         311    83838975 :                 uptkrate = cldaero_uptakerate( xl, cldnum(i,k), cfact(i,k), cldfrc(i,k), tfld(i,k),  press(i,k) )
         312             :                 ! average uptake rate over dtime
         313    83838975 :                 uptkrate = (1.0_r8 - exp(-min(100._r8,dtime*uptkrate))) / dtime
         314             :
         315             :                 ! dso4dt_gasuptk = so4_c tendency from h2so4 gas uptake (mol/mol/s)
         316             :                 ! dmsadt_gasuptk = msa_c tendency from msa gas uptake (mol/mol/s)
         317    83838975 :                 dso4dt_gasuptk = xh2so4(i,k) * uptkrate
         318    83838975 :                 if (id_msa > 0) then
         319           0 :                    dmsadt_gasuptk = xmsa(i,k) * uptkrate
         320             :                 else
         321             :                    dmsadt_gasuptk = 0.0_r8
         322             :                 end if
         323             :
         324             :                 ! if no modes have msa aerosol, then "rename" scavenged msa gas to so4
         325    83838975 :                 dmsadt_gasuptk_toso4 = 0.0_r8
         326    83838975 :                 dmsadt_gasuptk_tomsa = dmsadt_gasuptk
         327    83838975 :                 if (ntot_msa_c == 0) then
         328    83838975 :                    dmsadt_gasuptk_tomsa = 0.0_r8
         329    83838975 :                    dmsadt_gasuptk_toso4 = dmsadt_gasuptk
         330             :                 end if
         331             :
         332             :                 !-----------------------------------------------------------------------
         333             :                 ! now compute TMR tendencies
         334             :                 ! this includes the above aqueous so2 chemistry AND
         335             :                 ! the uptake of highly soluble aerosol precursor gases (h2so4, msa, ...)
         336             :                 ! AND the wetremoval of dissolved, unreacted so2 and h2o2
         337             :
         338    83838975 :                 dso4dt_aqrxn = (delso4_o3rxn + delso4_hprxn(i,k)) / dtime
         339    83838975 :                 dso4dt_hprxn = delso4_hprxn(i,k) / dtime
         340             :
         341             :                 ! fwetrem = fraction of in-cloud-water material that is wet removed
         342             :                 ! fwetrem = max( 0.0_r8, (1.0_r8-exp(-min(100._r8,dtime*clwlrat(i,k)))) )
         343    83838975 :                 fwetrem = 0.0_r8 ! don't have so4 & msa wet removal here
         344             :
         345             :                 ! compute TMR tendencies for so4 and msa aerosol-in-cloud-water
         346   419194875 :                 do n = 1, ntot_amode
         347   335355900 :                    l = lptr_so4_cw_amode(n) - loffset
         348   335355900 :                    if (l > 0) then
         349   251516925 :                       dqdt_aqso4(i,k,l) = faqgain_so4(n)*dso4dt_aqrxn*cldfrc(i,k)
         350   251516925 :                       dqdt_aqh2so4(i,k,l) = faqgain_so4(n)* &
         351   251516925 :                            (dso4dt_gasuptk + dmsadt_gasuptk_toso4)*cldfrc(i,k)
         352   251516925 :                       dqdt_aq = dqdt_aqso4(i,k,l) + dqdt_aqh2so4(i,k,l)
         353   251516925 :                       dqdt_wr = -fwetrem*dqdt_aq
         354   251516925 :                       dqdt= dqdt_aq + dqdt_wr
         355   251516925 :                       qcw(i,k,l) = qcw(i,k,l) + dqdt*dtime
         356             :                    end if
         357             :
         358   335355900 :                    l = lptr_msa_cw_amode(n) - loffset
         359   335355900 :                    if (l > 0) then
         360           0 :                       dqdt_aq = faqgain_msa(n)*dmsadt_gasuptk_tomsa*cldfrc(i,k)
         361           0 :                       dqdt_wr = -fwetrem*dqdt_aq
         362           0 :                       dqdt = dqdt_aq + dqdt_wr
         363           0 :                       qcw(i,k,l) = qcw(i,k,l) + dqdt*dtime
         364             :                    end if
         365             :
         366   335355900 :                    l = lptr_nh4_cw_amode(n) - loffset
         367   419194875 :                    if (l > 0) then
         368           0 :                       if (delnh4 > 0.0_r8) then
         369           0 :                          dqdt_aq = faqgain_so4(n)*delnh4/dtime*cldfrc(i,k)
         370           0 :                          dqdt = dqdt_aq
         371           0 :                          qcw(i,k,l) = qcw(i,k,l) + dqdt*dtime
         372             :                       else
         373           0 :                          dqdt = (qcw(i,k,l)/max(xnh4c(i,k),1.0e-35_r8)) &
         374           0 :                               *delnh4/dtime*cldfrc(i,k)
         375           0 :                          qcw(i,k,l) = qcw(i,k,l) + dqdt*dtime
         376             :                       endif
         377             :                    end if
         378             :                 end do
         379             :
         380             :                 ! For gas species, tendency includes
         381             :                 ! reactive uptake to cloud water that essentially transforms the gas to
         382             :                 ! a different species. Wet removal associated with this is applied
         383             :                 ! to the "new" species (e.g., so4_c) rather than to the gas.
         384             :                 ! wet removal of the unreacted gas that is dissolved in cloud water.
         385             :                 ! Need to multiply both these parts by cldfrc
         386             :
         387             :                 ! h2so4 (g) & msa (g)
         388    83838975 :                 qin(i,k,id_h2so4) = qin(i,k,id_h2so4) - dso4dt_gasuptk * dtime * cldfrc(i,k)
         389    83838975 :                 if (id_msa > 0) qin(i,k,id_msa) = qin(i,k,id_msa) - dmsadt_gasuptk * dtime * cldfrc(i,k)
         390             :
         391             :                 ! so2 -- the first order loss rate for so2 is frso2_c*clwlrat(i,k)
         392             :                 ! fwetrem = max( 0.0_r8, (1.0_r8-exp(-min(100._r8,dtime*frso2_c*clwlrat(i,k)))) )
         393    83838975 :                 fwetrem = 0.0_r8 ! don't include so2 wet removal here
         394             :
         395    83838975 :                 dqdt_wr = -fwetrem*xso2(i,k)/dtime*cldfrc(i,k)
         396    83838975 :                 dqdt_aq = -dso4dt_aqrxn*cldfrc(i,k)
         397    83838975 :                 dqdt = dqdt_aq + dqdt_wr
         398    83838975 :                 qin(i,k,id_so2) = qin(i,k,id_so2) + dqdt * dtime
         399             :
         400             :                 ! h2o2 -- the first order loss rate for h2o2 is frh2o2_c*clwlrat(i,k)
         401             :                 ! fwetrem = max( 0.0_r8, (1.0_r8-exp(-min(100._r8,dtime*frh2o2_c*clwlrat(i,k)))) )
         402    83838975 :                 fwetrem = 0.0_r8 ! don't include h2o2 wet removal here
         403             :
         404    83838975 :                 dqdt_wr = -fwetrem*xh2o2(i,k)/dtime*cldfrc(i,k)
         405    83838975 :                 dqdt_aq = -dso4dt_hprxn*cldfrc(i,k)
         406    83838975 :                 dqdt = dqdt_aq + dqdt_wr
         407    83838975 :                 qin(i,k,id_h2o2) = qin(i,k,id_h2o2) + dqdt * dtime
         408             :
         409             :                 ! NH3
         410    83838975 :                 if (id_nh3>0) then
         411           0 :                    dqdt_aq = delnh3/dtime*cldfrc(i,k)
         412           0 :                    dqdt = dqdt_aq
         413           0 :                    qin(i,k,id_nh3) = qin(i,k,id_nh3) + dqdt * dtime
         414             :                 endif
         415             :
         416             :                 ! for SO4 from H2O2/O3 budgets
         417    83838975 :                 dqdt_aqhprxn(i,k) = dso4dt_hprxn*cldfrc(i,k)
         418    83838975 :                 dqdt_aqo3rxn(i,k) = (dso4dt_aqrxn - dso4dt_hprxn)*cldfrc(i,k)
         419             :
         420             :              ENDIF !! WHEN CLOUD IS PRESENTED
         421             :           endif cloud
         422             :        enddo col_loop
         423             :     enddo lev_loop
         424             :
         425             :     !==============================================================
         426             :     ! ... Update the mixing ratios
         427             :     !==============================================================
         428   139982544 :     do k = 1,pver
         429             :
         430   692466840 :        do n = 1, ntot_amode
         431             :
         432   553973472 :           l = lptr_so4_cw_amode(n) - loffset
         433   553973472 :           if (l > 0) then
         434  6937551504 :              qcw(:,k,l) = MAX(qcw(:,k,l), small_value )
         435             :           end if
         436   553973472 :           l = lptr_msa_cw_amode(n) - loffset
         437   553973472 :           if (l > 0) then
         438           0 :              qcw(:,k,l) = MAX(qcw(:,k,l), small_value )
         439             :           end if
         440   553973472 :           l = lptr_nh4_cw_amode(n) - loffset
         441   692466840 :           if (l > 0) then
         442           0 :              qcw(:,k,l) = MAX(qcw(:,k,l), small_value )
         443             :           end if
         444             :
         445             :        end do
         446             :
         447  2312517168 :        qin(:,k,id_so2) =  MAX( qin(:,k,id_so2),    small_value )
         448             :
         449   139982544 :        if ( id_nh3 > 0 ) then
         450           0 :           qin(:,k,id_nh3) =  MAX( qin(:,k,id_nh3),    small_value )
         451             :        endif
         452             :
         453             :     end do
         454             :
         455             :     ! diagnostics
         456             :
         457     7445880 :     do n = 1, ntot_amode
         458     5956704 :        m = lptr_so4_cw_amode(n)
         459     5956704 :        l = m - loffset
         460     7445880 :        if (l > 0) then
         461    74597328 :           aqso4(:,n)=0._r8
         462   419947632 :           do k=1,pver
         463  6942019032 :              do i=1,ncol
         464 19566214200 :                 aqso4(i,n)=aqso4(i,n)+dqdt_aqso4(i,k,l)*adv_mass(l)/mbar(i,k) &
         465 26503765704 :                      *pdel(i,k)/gravit ! kg/m2/s
         466             :              enddo
         467             :           enddo
         468             :
         469    74597328 :           aqh2so4(:,n)=0._r8
         470   419947632 :           do k=1,pver
         471  6942019032 :              do i=1,ncol
         472 19566214200 :                 aqh2so4(i,n)=aqh2so4(i,n)+dqdt_aqh2so4(i,k,l)*adv_mass(l)/mbar(i,k) &
         473 26503765704 :                      *pdel(i,k)/gravit ! kg/m2/s
         474             :              enddo
         475             :           enddo
         476             :        endif
         477             :     end do
         478             :
         479    24865776 :     aqso4_h2o2(:) = 0._r8
         480   139982544 :     do k=1,pver
         481  2314006344 :        do i=1,ncol
         482  4348047600 :           aqso4_h2o2(i)=aqso4_h2o2(i)+dqdt_aqhprxn(i,k)*specmw_so4_amode/mbar(i,k) &
         483  6660564768 :                   *pdel(i,k)/gravit ! kg SO4 /m2/s
         484             :        enddo
         485             :     enddo
         486             :
         487     1489176 :     if (present(aqso4_h2o2_3d)) then
         488           0 :        aqso4_h2o2_3d(:,:) = 0._r8
         489           0 :        do k=1,pver
         490           0 :           do i=1,ncol
         491           0 :              aqso4_h2o2_3d(i,k)=dqdt_aqhprxn(i,k)*specmw_so4_amode/mbar(i,k) &
         492           0 :                                 *pdel(i,k)/gravit ! kg SO4 /m2/s
         493             :           enddo
         494             :        enddo
         495             :     end if
         496             :
         497    24865776 :     aqso4_o3(:)=0._r8
         498   139982544 :     do k=1,pver
         499  2314006344 :        do i=1,ncol
         500  4348047600 :           aqso4_o3(i)=aqso4_o3(i)+dqdt_aqo3rxn(i,k)*specmw_so4_amode/mbar(i,k) &
         501  6660564768 :                   *pdel(i,k)/gravit ! kg SO4 /m2/s
         502             :        enddo
         503             :     enddo
         504             :
         505     1489176 :     if (present(aqso4_o3_3d)) then
         506           0 :        aqso4_o3_3d(:,:)=0._r8
         507           0 :        do k=1,pver
         508           0 :           do i=1,ncol
         509           0 :              aqso4_o3_3d(i,k)=dqdt_aqo3rxn(i,k)*specmw_so4_amode/mbar(i,k) &
         510           0 :                               *pdel(i,k)/gravit ! kg SO4 /m2/s
         511             :           enddo
         512             :        enddo
         513             :     end if
         514             :
         515             :   end subroutine sox_cldaero_update
         516             :
         517             :   !----------------------------------------------------------------------------------
         518             :   !----------------------------------------------------------------------------------
         519     1489176 :   subroutine sox_cldaero_destroy_obj( conc_obj )
         520             :     type(cldaero_conc_t), pointer :: conc_obj
         521             :
         522     1489176 :     call cldaero_deallocate( conc_obj )
         523             :
         524     1489176 :   end subroutine sox_cldaero_destroy_obj
         525             :
         526             : end module sox_cldaero_mod

Generated by: LCOV version 1.14