LCOV - code coverage report
Current view: | top level - atmos_phys/schemes/zhang_mcfarlane - zm_conv_evap.F90 (source / functions) |  |  | Hit | Total | Coverage
Test: | coverage.info |  | Lines: | 67 | 88 | 76.1 %
Date: | 2024-12-17 22:39:59 |  | Functions: | 1 | 1 | 100.0 %

              Line data    Source code
           1             : module zm_conv_evap
           2             :
           3             :   use ccpp_kinds, only:  kind_phys
           4             :
           5             : ! CACNOTE - Need to ccpp'ize cloud_fraction
           6             :   use cloud_fraction,  only: cldfrc_fice
           7             :
           8             :   implicit none
           9             :
          10             :   save
          11             :   private                         ! Make default type private to the module
          12             : !
          13             : ! PUBLIC: interfaces
          14             : !
          15             :   public zm_conv_evap_run         ! evaporation of precip from ZM schemea
          16             :
          17             : contains
          18             :
          19             :
          20             : !===============================================================================
          21             : !> \section arg_table_zm_conv_evap_run Argument Table
          22             : !! \htmlinclude zm_conv_evap_run.html
          23             : !!
          24     1495368 : subroutine zm_conv_evap_run(ncol, pver, pverp, &
          25             :      gravit, latice, latvap, tmelt, &
          26             :      cpres, ke, ke_lnd, zm_org, &
          27     2990736 :      t,pmid,pdel,q, &
          28     1495368 :      landfrac, &
          29     1495368 :      tend_s, tend_s_snwprd, tend_s_snwevmlt, tend_q, &
          30     1495368 :      prdprec, cldfrc, deltat,  &
          31     5981472 :      prec, snow, ntprprd, ntsnprd, flxprec, flxsnow, prdsnow)
          32             :
          33             :
          34             : !-----------------------------------------------------------------------
          35             : ! Compute tendencies due to evaporation of rain from ZM scheme
          36             : !--
          37             : ! Compute the total precipitation and snow fluxes at the surface.
          38             : ! Add in the latent heat of fusion for snow formation and melt, since it not dealt with
          39             : ! in the Zhang-MacFarlane parameterization.
          40             : ! Evaporate some of the precip directly into the environment using a Sundqvist type algorithm
          41             : !-----------------------------------------------------------------------
          42             :
          43             : !CACNOTE - Not sure what to do about qsat_water
          44             :     use wv_saturation,  only: qsat
          45             :
          46             : !------------------------------Arguments--------------------------------
          47             :     integer,intent(in) :: ncol                               ! number of columns
          48             :     integer,intent(in) :: pver, pverp
          49             :     real(kind_phys),intent(in) :: gravit                     ! gravitational acceleration (m s-2)
          50             :     real(kind_phys),intent(in) :: latice                     ! Latent heat of fusion (J kg-1)
          51             :     real(kind_phys),intent(in) :: latvap                     ! Latent heat of vaporization (J kg-1)
          52             :     real(kind_phys),intent(in) :: tmelt                      ! Freezing point of water (K)
          53             :     real(kind_phys), intent(in) :: cpres      ! specific heat at constant pressure in j/kg-degk.
          54             :     real(kind_phys), intent(in) :: ke           ! Tunable evaporation efficiency set from namelist input zmconv_ke
          55             :     real(kind_phys), intent(in) :: ke_lnd
          56             :     logical, intent(in)         :: zm_org
          57             :     real(kind_phys),intent(in), dimension(:,:) :: t          ! temperature (K)                              (ncol,pver)
          58             :     real(kind_phys),intent(in), dimension(:,:) :: pmid       ! midpoint pressure (Pa)                       (ncol,pver)
          59             :     real(kind_phys),intent(in), dimension(:,:) :: pdel       ! layer thickness (Pa)                         (ncol,pver)
          60             :     real(kind_phys),intent(in), dimension(:,:) :: q          ! water vapor (kg/kg)                          (ncol,pver)
          61             :     real(kind_phys),intent(in), dimension(:) :: landfrac     ! land fraction                                (ncol)
          62             :     real(kind_phys),intent(inout), dimension(:,:) :: tend_s     ! heating rate (J/kg/s)                     (ncol,pver)
          63             :     real(kind_phys),intent(inout), dimension(:,:) :: tend_q     ! water vapor tendency (kg/kg/s)            (ncol,pver)
          64             :     real(kind_phys),intent(out), dimension(:,:) :: tend_s_snwprd ! Heating rate of snow production        (ncol,pver)
          65             :     real(kind_phys),intent(out), dimension(:,:) :: tend_s_snwevmlt ! Heating rate of evap/melting of snow (ncol,pver)
          66             :
          67             :
          68             :
          69             :     real(kind_phys), intent(in   ) :: prdprec(:,:)! precipitation production (kg/ks/s)                      (ncol,pver)
          70             :     real(kind_phys), intent(in   ) :: cldfrc(:,:) ! cloud fraction                                          (ncol,pver)
          71             :     real(kind_phys), intent(in   ) :: deltat             ! time step
          72             :
          73             :     real(kind_phys), intent(inout) :: prec(:)        ! Convective-scale preciptn rate                       (ncol)
          74             :     real(kind_phys), intent(out)   :: snow(:)        ! Convective-scale snowfall rate                       (ncol)
          75             :
          76             :     real(kind_phys), optional, intent(in), allocatable  :: prdsnow(:,:) ! snow production (kg/ks/s)
          77             :
          78             : !
          79             : !---------------------------Local storage-------------------------------
          80             :
          81     2990736 :     real(kind_phys) :: es    (ncol,pver)    ! Saturation vapor pressure
          82     2990736 :     real(kind_phys) :: fice   (ncol,pver)    ! ice fraction in precip production
          83     2990736 :     real(kind_phys) :: fsnow_conv(ncol,pver) ! snow fraction in precip production
          84     2990736 :     real(kind_phys) :: qs   (ncol,pver)    ! saturation specific humidity
          85             :     real(kind_phys),intent(out) :: flxprec(:,:)   ! Convective-scale flux of precip at interfaces (kg/m2/s) ! (ncol,pverp)
          86             :     real(kind_phys),intent(out) :: flxsnow(:,:)   ! Convective-scale flux of snow   at interfaces (kg/m2/s) ! (ncol,pverp)
          87             :     real(kind_phys),intent(out) :: ntprprd(:,:)   ! net precip production in layer                          ! (ncol,pver)
          88             :     real(kind_phys),intent(out) :: ntsnprd(:,:)   ! net snow production in layer                            ! (ncol,pver)
          89             :     real(kind_phys) :: work1                  ! temp variable (pjr)
          90             :     real(kind_phys) :: work2                  ! temp variable (pjr)
          91             :
          92     2990736 :     real(kind_phys) :: evpvint(ncol)         ! vertical integral of evaporation
          93     2990736 :     real(kind_phys) :: evpprec(ncol)         ! evaporation of precipitation (kg/kg/s)
          94     2990736 :     real(kind_phys) :: evpsnow(ncol)         ! evaporation of snowfall (kg/kg/s)
          95     2990736 :     real(kind_phys) :: snowmlt(ncol)         ! snow melt tendency in layer
          96     2990736 :     real(kind_phys) :: flxsntm(ncol)         ! flux of snow into layer, after melting
          97             :
          98             :     real(kind_phys) :: kemask
          99             :     real(kind_phys) :: evplimit               ! temp variable for evaporation limits
         100             :     real(kind_phys) :: rlat(ncol)
         101             :     real(kind_phys) :: dum
         102             :     real(kind_phys) :: omsm
         103             :
         104             :     integer :: i,k                     ! longitude,level indices
         105             :     logical :: old_snow
         106             :
         107             :
         108             : !-----------------------------------------------------------------------
         109             :
         110             :     ! If prdsnow is passed in and allocated, then use it in the calculation, otherwise
         111             :     ! use the old snow calculation
         112     1495368 :     old_snow=.true.
         113     1495368 :     if (present(prdsnow)) then
         114           0 :        if (allocated(prdsnow)) then
         115           0 :           old_snow=.false.
         116             :        end if
         117             :     end if
         118             :
         119             : ! convert input precip to kg/m2/s
         120    24969168 :     prec(:ncol) = prec(:ncol)*1000._kind_phys
         121             :
         122             : ! determine saturation vapor pressure
         123   140564592 :     do k = 1,pver
         124   140564592 :        call qsat(t(1:ncol,k), pmid(1:ncol,k), es(1:ncol,k), qs(1:ncol,k), ncol)
         125             :     end do
         126             : ! determine ice fraction in rain production (use cloud water parameterization fraction at present)
         127             : !REMOVECAM - no longer need these when CAM is retired and pcols no longer exists
         128  2323627992 :     fice(:,:) = 0._kind_phys
         129  2323627992 :     fsnow_conv(:,:) = 0._kind_phys
         130             : !REMOVECAM_END
         131     1495368 :     call cldfrc_fice(ncol, t(1:ncol,:), fice(1:ncol,:), fsnow_conv(1:ncol,:))
         132             :
         133             : ! zero the flux integrals on the top boundary
         134    24969168 :     flxprec(:ncol,1) = 0._kind_phys
         135    24969168 :     flxsnow(:ncol,1) = 0._kind_phys
         136    24969168 :     evpvint(:ncol)   = 0._kind_phys
         137             :     omsm=0.9999_kind_phys
         138             :
         139   140564592 :     do k = 1, pver
         140  2323627992 :        do i = 1, ncol
         141             :
         142             : ! Melt snow falling into layer, if necessary.
         143  2183063400 :         if( old_snow ) then
         144  2183063400 :           if (t(i,k) > tmelt) then
         145   379833079 :              flxsntm(i) = 0._kind_phys
         146   379833079 :              snowmlt(i) = flxsnow(i,k) * gravit/ pdel(i,k)
         147             :           else
         148  1803230321 :              flxsntm(i) = flxsnow(i,k)
         149  1803230321 :              snowmlt(i) = 0._kind_phys
         150             :           end if
         151             :         else
         152             :           ! make sure melting snow doesn't reduce temperature below threshold
         153           0 :           if (t(i,k) > tmelt) then
         154           0 :               dum = -latice/cpres*flxsnow(i,k)*gravit/pdel(i,k)*deltat
         155           0 :               if (t(i,k) + dum .le. tmelt) then
         156           0 :                 dum = (t(i,k)-tmelt)*cpres/latice/deltat
         157           0 :                 dum = dum/(flxsnow(i,k)*gravit/pdel(i,k))
         158           0 :                 dum = max(0._kind_phys,dum)
         159           0 :                 dum = min(1._kind_phys,dum)
         160             :               else
         161             :                 dum = 1._kind_phys
         162             :               end if
         163           0 :               dum = dum*omsm
         164           0 :               flxsntm(i) = flxsnow(i,k)*(1.0_kind_phys-dum)
         165           0 :               snowmlt(i) = dum*flxsnow(i,k)*gravit/ pdel(i,k)
         166             :           else
         167           0 :              flxsntm(i) = flxsnow(i,k)
         168           0 :              snowmlt(i) = 0._kind_phys
         169             :           end if
         170             :         end if
         171             :
         172             : ! relative humidity depression must be > 0 for evaporation
         173  2183063400 :           evplimit = max(1._kind_phys - q(i,k)/qs(i,k), 0._kind_phys)
         174             :
         175  2183063400 :           if (zm_org) then
         176           0 :              kemask = ke * (1._kind_phys - landfrac(i)) + ke_lnd * landfrac(i)
         177             :           else
         178  2183063400 :              kemask = ke
         179             :           endif
         180             :
         181             : ! total evaporation depends on flux in the top of the layer
         182             : ! flux prec is the net production above layer minus evaporation into environmet
         183  2183063400 :           evpprec(i) = kemask * (1._kind_phys - cldfrc(i,k)) * evplimit * sqrt(flxprec(i,k))
         184             :
         185             : ! Don't let evaporation supersaturate layer (approx). Layer may already be saturated.
         186             : ! Currently does not include heating/cooling change to qs
         187  2183063400 :           evplimit   = max(0._kind_phys, (qs(i,k)-q(i,k)) / deltat)
         188             :
         189             : ! Don't evaporate more than is falling into the layer - do not evaporate rain formed
         190             : ! in this layer but if precip production is negative, remove from the available precip
         191             : ! Negative precip production occurs because of evaporation in downdrafts.
         192  2183063400 :           evplimit   = min(evplimit, flxprec(i,k) * gravit / pdel(i,k))
         193             :
         194             : ! Total evaporation cannot exceed input precipitation
         195  2183063400 :           evplimit   = min(evplimit, (prec(i) - evpvint(i)) * gravit / pdel(i,k))
         196             :
         197  2183063400 :           evpprec(i) = min(evplimit, evpprec(i))
         198  2183063400 :           if( .not.old_snow ) then
         199           0 :             evpprec(i) = max(0._kind_phys, evpprec(i))
         200           0 :             evpprec(i) = evpprec(i)*omsm
         201             :           end if
         202             :
         203             :
         204             : ! evaporation of snow depends on snow fraction of total precipitation in the top after melting
         205  2183063400 :           if (flxprec(i,k) > 0._kind_phys) then
         206             : !            prevent roundoff problems
         207    82509498 :              work1 = min(max(0._kind_phys,flxsntm(i)/flxprec(i,k)),1._kind_phys)
         208    82509498 :              evpsnow(i) = evpprec(i) * work1
         209             :           else
         210  2100553902 :              evpsnow(i) = 0._kind_phys
         211             :           end if
         212             :
         213             : ! vertically integrated evaporation
         214  2183063400 :           evpvint(i) = evpvint(i) + evpprec(i) * pdel(i,k)/gravit
         215             :
         216             : ! net precip production is production - evaporation
         217  2183063400 :           ntprprd(i,k) = prdprec(i,k) - evpprec(i)
         218             : ! net snow production is precip production * ice fraction - evaporation - melting
         219             : ! the small amount added to flxprec in the work1 expression has been increased from
         220             : ! 1e-36 to 8.64e-11 (1e-5 mm/day).  This causes the temperature based partitioning
         221             : ! scheme to be used for small flxprec amounts.  This is to address error growth problems.
         222             :
         223  2183063400 :       if( old_snow ) then
         224  2183063400 :           if (flxprec(i,k).gt.0._kind_phys) then
         225    82509498 :              work1 = min(max(0._kind_phys,flxsnow(i,k)/flxprec(i,k)),1._kind_phys)
         226             :           else
         227             :              work1 = 0._kind_phys
         228             :           endif
         229             :
         230  2183063400 :           work2 = max(fsnow_conv(i,k), work1)
         231  2183063400 :           if (snowmlt(i).gt.0._kind_phys) work2 = 0._kind_phys
         232  2183063400 :           ntsnprd(i,k) = prdprec(i,k)*work2 - evpsnow(i) - snowmlt(i)
         233  2183063400 :           tend_s_snwprd  (i,k) = prdprec(i,k)*work2*latice
         234  2183063400 :           tend_s_snwevmlt(i,k) = - ( evpsnow(i) + snowmlt(i) )*latice
         235             :       else
         236           0 :           ntsnprd(i,k) = prdsnow(i,k) - min(flxsnow(i,k)*gravit/pdel(i,k), evpsnow(i)+snowmlt(i))
         237           0 :           tend_s_snwprd  (i,k) = prdsnow(i,k)*latice
         238           0 :           tend_s_snwevmlt(i,k) = -min(flxsnow(i,k)*gravit/pdel(i,k), evpsnow(i)+snowmlt(i) )*latice
         239             :       end if
         240             :
         241             : ! precipitation fluxes
         242  2183063400 :           flxprec(i,k+1) = flxprec(i,k) + ntprprd(i,k) * pdel(i,k)/gravit
         243  2183063400 :           flxsnow(i,k+1) = flxsnow(i,k) + ntsnprd(i,k) * pdel(i,k)/gravit
         244             :
         245             : ! protect against rounding error
         246  2183063400 :           flxprec(i,k+1) = max(flxprec(i,k+1), 0._kind_phys)
         247  2183063400 :           flxsnow(i,k+1) = max(flxsnow(i,k+1), 0._kind_phys)
         248             :
         249             : ! heating (cooling) and moistening due to evaporation
         250             : ! - latent heat of vaporization for precip production has already been accounted for
         251             : ! - snow is contained in prec
         252  2183063400 :           if( old_snow ) then
         253  2183063400 :              tend_s(i,k)   =-evpprec(i)*latvap + ntsnprd(i,k)*latice
         254             :           else
         255           0 :              tend_s(i,k)   =-evpprec(i)*latvap + tend_s_snwevmlt(i,k)
         256             :           end if
         257  2322132624 :           tend_q(i,k) = evpprec(i)
         258             :        end do
         259             :     end do
         260             :
         261             : ! set output precipitation rates (m/s)
         262    24969168 :     prec(:ncol) = flxprec(:ncol,pver+1) / 1000._kind_phys
         263    24969168 :     snow(:ncol) = flxsnow(:ncol,pver+1) / 1000._kind_phys
         264             :
         265     1495368 :   end subroutine zm_conv_evap_run
         266             :
         267             :
         268             : end module zm_conv_evap

Generated by: LCOV version 1.14