subroutine mfpbl | ( | integer | im, |
integer | ix, | ||
integer | km, | ||
integer | ntrac, | ||
real(kind=kind_phys) | delt, | ||
logical, dimension(im) | cnvflg, | ||
real(kind=kind_phys), dimension(im,km) | zl, | ||
real(kind=kind_phys), dimension(im,km+1) | zm, | ||
real(kind=kind_phys), dimension(im,km) | thvx, | ||
real(kind=kind_phys), dimension(ix,km,ntrac) | q1, | ||
real(kind=kind_phys), dimension(ix,km) | t1, | ||
real(kind=kind_phys), dimension(ix,km) | u1, | ||
real(kind=kind_phys), dimension(ix,km) | v1, | ||
real(kind=kind_phys), dimension(im) | hpbl, | ||
integer, dimension(im) | kpbl, | ||
real(kind=kind_phys), dimension(im) | sflx, | ||
real(kind=kind_phys), dimension(im) | ustar, | ||
real(kind=kind_phys), dimension(im) | wstar, | ||
real(kind=kind_phys), dimension(im,km) | xmf, | ||
real(kind=kind_phys), dimension(im,km) | tcko, | ||
real(kind=kind_phys), dimension(im,km,ntrac) | qcko, | ||
real(kind=kind_phys), dimension(im,km) | ucko, | ||
real(kind=kind_phys), dimension(im,km) | vcko | ||
) |
The mfpbl routines works as follows: if the PBL is convective, first, the ascending parcel entrainment rate is calculated as a function of height. Next, a surface parcel is initiated according to surface layer properties and the updraft buoyancy is calculated as a function of height. Next, using the buoyancy and entrainment values, the parcel vertical velocity is calculated using a well known steady-state budget equation. With the profile of updraft vertical velocity, the PBL height is recalculated as the height where the updraft vertical velocity returns to 0, and the entrainment profile is updated with the new PBL height. Finally, the mass flux profile is calculated using the updraft vertical velocity and assumed updraft fraction and the updraft properties are calculated using the updated entrainment profile, surface values, and environmental profiles.
[in] | im | integer, number of used points |
[in] | ix | integer, horizontal dimension |
[in] | km | integer, vertical layer dimension |
[in] | ntrac | integer, number of tracers |
[in] | delt | real, physics time step |
[in] | cnvflg | logical, im, flag to denote a strongly unstable (convective) PBL |
[in] | zl | real, (im, km), height of grid centers |
[in] | zm | real, (im, km+1), height of grid interfaces |
[in] | thvx | real, (im, km), virtual potential temperature at grid centers ( \( K \)) |
[in] | q1 | real, (ix, km, ntrac), layer mean tracer concentration (units?) |
[in] | t1 | real, (ix, km), layer mean temperature ( \( K \)) |
[in] | u1 | real, (ix, km), u component of layer wind ( \( m s^{-1} \)) |
[in] | v1 | real, (ix, km), v component of layer wind ( \( m s^{-1} \)) |
[in,out] | hpbl | real, im, PBL top height (m) |
[in,out] | kpbl | integer, im, PBL top index |
[in] | sflx | real, im, total surface heat flux (units?) |
[in] | ustar | real, im, surface friction velocity |
[in] | wstar | real, im, convective velocity scale |
[out] | xmf | real, (im, km), updraft mass flux |
[in,out] | tcko | real, (im, km), updraft temperature ( \( K \)) |
[in,out] | qcko | real, (im, km, ntrac), updraft tracer concentration (units?) |
[in,out] | ucko | real, (im, km), updraft u component of horizontal momentum ( \( m s^{-1} \)) |
[in,out] | vcko | real, (im, km), updraft v component of horizontal momentum ( \( m s^{-1} \)) |
Since the mfpbl subroutine is called regardless of whether the PBL is convective, a check of the convective PBL flag is performed and the subroutine returns back to moninedmf (with the output variables set to the initialized values) if the PBL is not convective.
Calculate the entrainment rate according to equation 16 in [164] for all levels (xlamue) and a default entrainment rate (xlamax) for use above the PBL top.
Using equations 17 and 7 from [164] along with \(u_*\), \(w_*\), and the previously diagnosed PBL height, the initial \(\theta_v\) of the updraft (and its surface buoyancy) is calculated.
From the second level to the middle of the vertical domain, the updraft virtual potential temperature is calculated using the entraining updraft equation as in equation 10 of [164], discretized as
\[ \frac{\theta_{v,u}^k - \theta_{v,u}^{k-1}}{\Delta z}=-\epsilon^{k-1}\left[\frac{1}{2}\left(\theta_{v,u}^k + \theta_{v,u}^{k-1}\right)-\frac{1}{2}\left(\overline{\theta_{v}}^k + \overline{\theta_v}^{k-1}\right)\right] \]
where the superscript \(k\) denotes model level, and subscript \(u\) denotes an updraft property, and the overbar denotes the grid-scale mean value.
Rather than use the vertical velocity equation given as equation 15 in [164] (which parameterizes the pressure term in terms of the updraft vertical velocity itself), this scheme uses the more widely used form of the steady state vertical velocity equation given as equation 6 in [169] discretized as
\[ \frac{w_{u,k}^2 - w_{u,k-1}^2}{\Delta z} = -2b_1\frac{1}{2}\left(\epsilon_k + \epsilon_{k-1}\right)\frac{1}{2}\left(w_{u,k}^2 + w_{u,k-1}^2\right) + 2b_2B \]
The constants used in the scheme are labeled \(bb1 = 2b_1\) and \(bb2 = 2b_2\) and are tuned to be equal to 1.8 and 3.5, respectively, close to the values proposed by [169] .
Find the level where the updraft vertical velocity is less than zero and linearly interpolate to find the height where it would be exactly zero. Set the PBL height to this determined height.
Recalculate the entrainment rate as before except use the updated value of the PBL height.
Calculate the mass flux:
\[ M = a_uw_u \]
where \(a_u\) is the tunable parameter that represents the fractional area of updrafts (currently set to 0.08). Limit the computed mass flux to be less than \(\frac{\Delta z}{\Delta t}\). This is different than what is done in [164] where the mass flux is the product of a tunable constant and the diagnosed standard deviation of \(w\).
The updraft properties are calculated according to the entraining updraft equation
\[ \frac{\partial \phi}{\partial z}=-\epsilon\left(\phi_u - \overline{\phi}\right) \]
where \(\phi\) is \(T\) or \(q\). The equation is discretized according to
\[ \frac{\phi_{u,k} - \phi_{u,k-1}}{\Delta z}=-\epsilon_{k-1}\left[\frac{1}{2}\left(\phi_{u,k} + \phi_{u,k-1}\right)-\frac{1}{2}\left(\overline{\phi}_k + \overline{\phi}_{k-1}\right)\right] \]
The exception is for the horizontal momentum components, which have been modified to account for the updraft-induced pressure gradient force, and use the following equation, following [76]
\[ \frac{\partial v}{\partial z} = -\epsilon\left(v_u - \overline{v}\right)+d_1\frac{\partial \overline{v}}{\partial z} \]
where \(d_1=0.55\) is a tunable constant.
References physcons::con_cp, and physcons::con_g.
Referenced by hedmf::hedmf_run().