CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches

◆ satmedmfvdif_run()

subroutine satmedmfvdif::satmedmfvdif_run ( integer, intent(in) im,
integer, intent(in) km,
integer, intent(in) ntrac,
integer, intent(in) ntcw,
integer, intent(in) ntiw,
integer, intent(in) ntke,
real(kind=kind_phys), intent(in) grav,
real(kind=kind_phys), intent(in) rd,
real(kind=kind_phys), intent(in) cp,
real(kind=kind_phys), intent(in) rv,
real(kind=kind_phys), intent(in) hvap,
real(kind=kind_phys), intent(in) hfus,
real(kind=kind_phys), intent(in) fv,
real(kind=kind_phys), intent(in) eps,
real(kind=kind_phys), intent(in) epsm1,
real(kind=kind_phys), dimension(:,:), intent(inout) dv,
real(kind=kind_phys), dimension(:,:), intent(inout) du,
real(kind=kind_phys), dimension(:,:), intent(inout) tdt,
real(kind=kind_phys), dimension(:,:,:), intent(inout) rtg,
real(kind=kind_phys), dimension(:,:), intent(in) u1,
real(kind=kind_phys), dimension(:,:), intent(in) v1,
real(kind=kind_phys), dimension(:,:), intent(in) t1,
real(kind=kind_phys), dimension(:,:,:), intent(in) q1,
real(kind=kind_phys), dimension(:,:), intent(in) swh,
real(kind=kind_phys), dimension(:,:), intent(in) hlw,
real(kind=kind_phys), dimension(:), intent(in) xmu,
real(kind=kind_phys), dimension(:), intent(in) garea,
real(kind=kind_phys), dimension(:), intent(in) psk,
real(kind=kind_phys), dimension(:), intent(in) rbsoil,
real(kind=kind_phys), dimension(:), intent(in) zorl,
real(kind=kind_phys), dimension(:), intent(in) u10m,
real(kind=kind_phys), dimension(:), intent(in) v10m,
real(kind=kind_phys), dimension(:), intent(in) fm,
real(kind=kind_phys), dimension(:), intent(in) fh,
real(kind=kind_phys), dimension(:), intent(in) tsea,
real(kind=kind_phys), dimension(:), intent(in) heat,
real(kind=kind_phys), dimension(:), intent(in) evap,
real(kind=kind_phys), dimension(:), intent(in) stress,
real(kind=kind_phys), dimension(:), intent(in) spd1,
integer, dimension(:), intent(out) kpbl,
real(kind=kind_phys), dimension(:,:), intent(in) prsi,
real(kind=kind_phys), dimension(:,:), intent(in) del,
real(kind=kind_phys), dimension(:,:), intent(in) prsl,
real(kind=kind_phys), dimension(:,:), intent(in) prslk,
real(kind=kind_phys), dimension(:,:), intent(in) phii,
real(kind=kind_phys), dimension(:,:), intent(in) phil,
real(kind=kind_phys), intent(in) delt,
logical, intent(in) dspheat,
real(kind=kind_phys), dimension(:), intent(out) dusfc,
real(kind=kind_phys), dimension(:), intent(out) dvsfc,
real(kind=kind_phys), dimension(:), intent(out) dtsfc,
real(kind=kind_phys), dimension(:), intent(out) dqsfc,
real(kind=kind_phys), dimension(:), intent(out) hpbl,
integer, dimension(:), intent(in) kinver,
real(kind=kind_phys), intent(in) xkzm_m,
real(kind=kind_phys), intent(in) xkzm_h,
real(kind=kind_phys), intent(in) xkzm_s,
integer, intent(in) index_of_temperature,
integer, intent(in) index_of_x_wind,
integer, intent(in) index_of_y_wind,
integer, intent(in) index_of_process_pbl,
integer, intent(in) ntqv,
integer, intent(in) ntoz,
real(kind=kind_phys), dimension(:,:,:), intent(inout), optional dtend,
integer, dimension(:,:), intent(in) dtidx,
logical, intent(in) gen_tend,
logical, intent(in) ldiag3d,
character(len=*), intent(out) errmsg,
integer, intent(out) errflg )
  1. Compute preliminary variables from input arguments
  • Compute physical height of the layer centers and interfaces from the geopotential height (zi and zl)
  • Compute horizontal grid size (gdx)
  • Initialize tke value at vertical layer centers and interfaces from tracer (tke and tkeh)
  • Compute reciprocal of \( \Delta z \) (rdzt)
  • Compute reciprocal of pressure (tx1, tx2)
    • Compute background vertical diffusivities for scalars and momentum (xkzo and xkzmo)
      • set background diffusivities as a function of horizontal grid size with xkzm_h & xkzm_m for gdx >= 25km and 0.01 for gdx=5m, i.e.,
        xkzm_hx = 0.01 + (xkzm_h - 0.01)/(xkgdx-5.) * (gdx-5.)
        xkzm_mx = 0.01 + (xkzm_h - 0.01)/(xkgdx-5.) * (gdx-5.)
  • Some output variables and logical flags are initialized
  • Compute \(\theta\)(theta), and \(q_l\)(qlx), \(\theta_e\)(thetae), \(\theta_v\)(thvx), \(\theta_{l,v}\) (thlvx) including ice water
  • The background vertical diffusivities in the inversion layers are limited to be less than or equal to xkzinv
  • Compute an empirical cloud fraction based on Xu and Randall (1996) [201]
  • Compute buoyancy modified by clouds
  • Initialize diffusion coefficients to 0 and calculate the total radiative heating rate (dku, dkt, radx)
  • Compute stable/unstable PBL flag (pblflg) based on the total surface energy flux (false if the total surface energy flux is into the surface)
  • Compute critical bulk Richardson number ( \(Rb_{cr}\)) (crb)
    • For the unstable PBL, crb is a constant (0.25)
    • For the stable boundary layer (SBL), \(Rb_{cr}\) varies with the surface Rossby number, \(R_{0}\), as given by Vickers and Mahrt (2004) [196]

      \[ Rb_{cr}=0.16(10^{-7}R_{0})^{-0.18}, R_{0}=\frac{U_{10}}{f_{0}z_{0}} \]

      where \(U_{10}\) is the wind speed at 10m above the ground surface, \(f_0\) the Coriolis parameter, and \(z_{0}\) the surface roughness length. To avoid too much variation, we restrict \(Rb_{cr}\) to vary within the range of 0.15~0.35
  • Compute \(\frac{\Delta t}{\Delta z}\) , \(u_*\)
  • Compute buoyancy \(\frac{\partial \theta_v}{\partial z}\) (bf) and the wind shear squared (shr2)
  • Given the thermal's properties and the critical Richardson number, a loop is executed to find the first level above the surface (kpblx) where the modified Richardson number is greater than the critical Richardson number, using equation 10a from Troen and Mahrt (1996) [193] (also equation 8 from Hong and Pan (1996) [92]):
  • Once the level is found, some linear interpolation is performed to find the exact height of the boundary layer top (where \(R_{i} > Rb_{cr}\)) and the PBL height (hpbl and kpbl) and the PBL top index are saved.
  1. Compute Monin-Obukhov similarity parameters
    • Calculate the Monin-Obukhov nondimensional stability paramter, commonly referred to as \(\zeta\) using the following equation from Businger et al.(1971) [28] (eqn 28):

      \[ \zeta = Ri_{sfc}\frac{F_m^2}{F_h} = \frac{z}{L} \]

      where \(F_m\) and \(F_h\) are surface Monin-Obukhov stability functions calculated in sfc_diff.f and \(L\) is the Obukhov length.
  • Calculate the nondimensional gradients of momentum and temperature ( \(\phi_m\) (phim) and \(\phi_h\)(phih)) are calculated using eqns 5 and 6 from Hong and Pan (1996) [92] depending on the surface layer stability:
    • For the unstable and neutral conditions:

      \[ \phi_m=(1-16\frac{0.1h}{L})^{-1/4} \phi_h=(1-16\frac{0.1h}{L})^{-1/2} \]

    • For the stable regime

      \[ \phi_m=\phi_t=(1+5\frac{0.1h}{L}) \]

  • The \(z/L\) (zol) is used as the stability criterion for the PBL.Currently, strong unstable (convective) PBL for \(z/L < -0.02\) and weakly and moderately unstable PBL for \(0>z/L>-0.02\)
    • Compute the velocity scale \(w_s\) (wscale) (eqn 22 of Han et al. 2019). It is represented by the value scaled at the top of the surface layer:

      \[ w_s=(u_*^3+7\alpha\kappa w_*^3)^{1/3} \]

      where \(u_*\) (ustar) is the surface friction velocity, \(\alpha\) is the ratio of the surface layer height to the PBL height (specified as sfcfrac =0.1), \(\kappa =0.4\) is the von Karman constant, and \(w_*\) is the convective velocity scale defined as eqn23 of Han et al.(2019):

      \[ w_{*}=[(g/T)\overline{(w'\theta_v^{'})}_0h]^{1/3} \]

  • Compute a thermal excess.
  • Compute enhance the pbl height by considering the thermal excess. The PBL height calculation follows the same procedure above, except that is uses an updated virtual potential temperature for the thermal.
  1. Determine whether stratocumulus layers exist and compute quantities
    • Starting at the PBL top and going downward, if the level is less than 2.5 km and \(q_l\geq q_{lcr}\) then set kcld = k (find the cloud top index in the PBL. If no cloud water above the threshold is hound, scuflg is set to F.
  • Starting at the PBL top and going downward, if the level is less than the cloud top, find the level of the minimum radiative heating rate wihin the cloud. If the level of the minimum is the lowest model level or the minimum radiative heating rate is positive, then set scuflg to F.
  1. Compute components for mass flux mixing by large thermals
    • If the PBL is convective, the updraft properties are initialized to be the same as the state variables.
  • Call mfpblt(), which is an EDMF parameterization (Siebesma et al.(2007) [180]) to take into account nonlocal transport by large eddies.
  • Call mfscu(), which is a new mass-flux parameterization for stratocumulus-top-induced turbulence mixing.
  1. Compute Prandtl number \(P_r\) (prn) and exchange coefficient varying with height
    • Within PBL, \(P_r\) is given as \(P_r=\phi_h/\phi_m\), where \(\phi_h\) and \(\phi_m\) are the surface layer ( \(z\leq0.1h\), where \(h\) is the PBL height) non-dimensional gradient functions for heat and momentum, respectively.
      • We assume that \(P_r=0.67\) in the stratocumulus and the unstable layers above PBL (prscu).
      • For the stable layers above PBL, \(P_r\) is given as \(P_r=1+2.1R_i\) (prnum) (Kennedy and Shapiro (1980)[105] )
  1. Compute an asymtotic mixing length
  • Following Bougeault and Lacarrere(1989), the characteristic length scale ( \(l_2\)) (eqn 10 in Han et al.(2019) [75]) is given by:

    \[ l_2=min(l_{up},l_{down}) \]

    and dissipation length scale \(l_d\) is given by:

    \[ l_d=(l_{up}l_{down})^{1/2} \]

    where \(l_{up}\) and \(l_{down}\) are the distances that a parcel having an initial TKE can travel upward and downward before being stopped by buoyancy effects.
  • Compute the surface layer length scale ( \(l_1\)) following Nakanishi (2001) [144] (eqn 9 of Han et al.(2019) [75])
  1. Compute eddy diffusivities
  1. Compute buoyancy and shear productions of TKE
  1. First predict tke due to tke production & dissipation(diss)
  1. Compute updraft & downdraft properties for TKE
  1. Compute tridiagonal matrix elements for turbulent kinetic energy
  1. Call tridit() to solve tridiagonal problem for TKE
  1. Compute tridiagonal matrix elements for heat and moisture
  1. Call tridin() to solve tridiagonal problem for heat and moisture
  1. Add TKE dissipative heating to temperature tendency
  1. Compute tridiagonal matrix elements for momentum
  1. Call tridi2() to solve tridiagonal problem for momentum
  1. Save PBL height for diagnostic purpose

Definition at line 62 of file satmedmfvdif.F.

References mfpblt_mod::mfpblt(), mfscu_mod::mfscu(), tridi_mod::tridi2(), tridi_mod::tridin(), and tridi_mod::tridit().

Here is the call graph for this function: