This file contains the CCPP-compliant SATMEDMF scheme (updated version) which computes subgrid vertical turbulence mixing using scale-aware TKE-based moist eddy-diffusion mass-flux (TKE-EDMF) parameterization (by Jongil Han). More...

Functions/Subroutines
subroutine	satmedmfvdifq_can_mod::satmedmfvdifq_can (im, km, ntrac, ntcw, ntrw,

Detailed Description

This file contains the CCPP-compliant SATMEDMF scheme (updated version) which computes subgrid vertical turbulence mixing using scale-aware TKE-based moist eddy-diffusion mass-flux (TKE-EDMF) parameterization (by Jongil Han).

This subroutine contains all of the logic for the scale-aware TKE-based moist eddy-diffusion mass-flux (TKE-EDMF, updated version) scheme. For local turbulence mixing, a TKE closure model is used. Updated version of satmedmfvdif.f to have better low level inversion, to reduce the cold bias in lower troposphere, and to reduce the negative wind speed bias in upper troposphere

Incorporate the LES-based changes for TC simulation (Chen et al.,2022 [37]) with additional improvements on MF working with Cu schemes.

Incorporate the TTE-EDMF; if (tte_edmf=.true.), TKE-EDMF scheme becomes TTE-EDMF scheme and the variable 'te' is read as TTE; if (tte_edmf=.false.), the variable 'te' is read as TKE.

Argument Table

GFS satmedmfvdifq General Algorithm

satmedmfvdifq_can() computes subgrid vertical turbulence mixing using the scale-aware TKE-based moist eddy-diffusion mass-flux (EDMF) parameterization of Han and Bretherton (2019) [77] .

The local turbulent mixing is represented by an eddy-diffusivity scheme which is a function of a prognostic TKE.
For the convective boundary layer, nonlocal transport by large eddies (mfpbltq.f), is represented using a mass flux approach (Siebesma et al.(2007) [179] ).
A mass-flux approach is also used to represent the stratocumulus-top-induced turbulence (mfscuq.f).

GFS satmedmfvdifq Detailed Algorithm

Function/Subroutine Documentation

◆ satmedmfvdifq_can()

subroutine satmedmfvdifq_can_mod::satmedmfvdifq_can	(	integer, intent(in)	im,
		integer, intent(in)	km,
		integer, intent(in)	ntrac,
		integer, intent(in)	ntcw,
		integer, intent(in)	ntrw
	)

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 minimum turbulent mixing length (rlmnz)

Compute background vertical diffusivities for scalars and momentum (xkzo and xkzmo)

set background diffusivities with xkzm_h & xkzm_m for gdx >= xkgdx and as a function of horizontal grid size for gdx < xkgdx
xkzm_hx = xkzm_h * (gdx / xkgdx)
xkzm_mx = xkzm_m * (gdx / xkgdx)

Some output variables and logical flags are initialized

Compute a function for green vegetation fraction and surface roughness. Entrainment rate in updraft is a function of vegetation fraction and surface roughness length

Compute $\theta$(theta), and $q_l$(qlx), $\theta_e$(thetae), $\theta_v$(thvx), $\theta_{l,v}$ (thlvx) including ice water

Compute an empirical cloud fraction based on Xu and Randall (1996) [200]

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)

Calculate the PBL height

The calculation of the boundary layer height follows Troen and Mahrt (1986) [194] section 3. The approach is to find the level in the column where a modified bulk Richardson number exceeds a critical value.

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) [195]
\[ 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) [194] (also equation 8 from Hong and Pan (1996) [95]):

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.

Compute mean tke within pbl for TKE-EDMF

Compute wind shear term as a sink term for updraft and downdraft velocity

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) [29] (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) [95] 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} \]

The counter-gradient terms for temperature and humidity are calculated.

Equation 4 of Hong and Pan (1996) [95] and are used to calculate the "scaled virtual temperature excess near the surface" (equation 9 in Hong and Pan (1996) [95]) for use in the mass-flux algorithm.

Compute mean tke within pbl for TTE-EDMF

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.

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 mfpbltq(), which is an EDMF parameterization (Siebesma et al.(2007) [179]) to take into account nonlocal transport by large eddies. For details of the mfpbltq subroutine, step into its documentation ::mfpbltq

Call mfscuq(), which is a new mass-flux parameterization for stratocumulus-top-induced turbulence mixing. For details of the mfscuq subroutine, step into its documentation ::mfscuq

unify mass fluxes with Cu

taper off MF in high-wind conditions

Compute Prandtl number \f$P_r\f$ (prn) and exchange coefficient varying with height

Compute an asymtotic mixing length

Following Bougeault and Lacarrere(1989), the characteristic length scale ( $l_2$) (eqn 10 in Han et al.(2019) [77]) 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) [145] (eqn 9 of Han et al.(2019) [77])

If sfc_rlm=1, use zk for elm within surface layer

Compute eddy diffusivities

Compute TKE.

Compute a minimum TKE deduced from background diffusivity for momentum.

Compute buoyancy and shear productions of TKE or TTE

First predict te due to te production & dissipation(diss)

Compute updraft & downdraft properties for TKE or TTE

Compute tridiagonal matrix elements for TKE or TTE

Call tridit() to solve tridiagonal problem for TKE

Recover the tendency of TKE or TTE

Compute tridiagonal matrix elements for heat and moisture

Call tridin() to solve tridiagonal problem for heat and moisture

Recover the tendencies of heat and moisture

Add TKE dissipative heating to temperature tendency

Compute tridiagonal matrix elements for momentum

Call tridi2() to solve tridiagonal problem for momentum

Recover the tendencies of momentum

Save PBL height for diagnostic purpose

Definition at line 50 of file satmedmfvdifq_can.F.

References mfpbltq_mod::mfpbltq(), mfscuq_mod::mfscuq(), tridi_mod::tridi2(), tridi_mod::tridin(), and tridi_mod::tridit().

Referenced by canopy_driver::canopy_driver_run().

Here is the call graph for this function:

Here is the caller graph for this function:

Functions/Subroutines