数值天气预报模式的行星边界层方案热力预报变量选择

图1 左图为对流边界层典型热力结构示意图，实线代表边界层热泡，右图为虚位温 $Θ_{v}$ 和垂直湍流浮力通量 $\bar{w^{'} θ_{v}^{'}}$ 的垂直廓线示意图（ $z_{i}$ 和B_S 分别代表边界层高度和地表浮力通量，改自赵昭和周博闻^[11]）

Fig.1 (left) Schematic of the classic convective structures in the convective boundary layer，solid lines outline a boundary layer thermal; (Right) the virtual potential temperature and buoyancy flux profiles ( $z_{i}$ and B_S represent the boundary layer depth and the surface buoyancy flux, respectively （Figure adapted from Zhao and Zhou^[11])

依据主导湍流过程的不同，对流边界层自下而上可分为近地层、混合层和夹卷层，边界层顶以上则统称为（不受边界层摩擦所束缚的）自由大气（图1）.各层中的主导物理过程决定了边界层虚位温 $θ_{v} (z)$ 的廓线形态.近地层直接受陆面加热影响，呈现 $θ_{v}$ 向上递减的超绝热廓线.由于受地面的限制，近地层中的湍流以尺度与离地高度相当的接地湍涡（Attached Eddy）为主^［7］；在近地层之上的混合层中，非局地湍流旺盛发展，在强烈的湍流混合作用下， $θ_{v}$ 的垂直梯度接近于零；在与自由大气相邻的夹卷层中，非局地热对流在上升过程中受自由大气中稳定层结的抑制，在消耗完自身的垂直动量后又下沉至自身的中性层结高度.在下沉过程中夹带了部分自由大气中的高熵空气进入混合层（即垂直夹卷过程），造成了混合层上层（一般为 $z / z_{i} > 0.4$ ）的弱层结稳定性^［8-10］.

图1的右半边展示了经典对流边界层虚位温和浮力通量垂直廓线.这里选取虚位温作为边界层的热力变量，在非饱和的条件下（在对流边界层中，只在晨雾等极少数天气条件下，需要考虑饱和条件，这些情况不在本文的研究范围内），依据状态方程可得：

$P = ρ R T = ρ d R d T + ρ v R v T = ρ R d T 1 + r v / ϵ 1 + r v = ρ R d T v$ ( $1$ )

其中， $ρ = ρ_{d} + ρ_{v}$ 为空气密度； $r_{v} = ρ_{v} / ρ_{d}$ 为水汽密度和干空气密度的混合比； $ϵ = R_{d} / R_{v}$ 为干空气和水汽的状态常数比，等同于水汽和干空气的分子质量比； $T_{v} \equiv T \frac{1 + r_{v} / ϵ}{1 + r_{v}}$ 为虚温，即将水汽对气压的作用转化为等压条件下同等密度的干空气温度.由于水汽的分子质量小于干空气的分子质量，即 $ϵ < 1$ ，因此有 $T_{v} > T$ ，即在等压条件下，通过 $T_{v}$ 将更轻的水汽等效视作更热的干空气.

若将式（1）中的 $T$ 通过无量纲压强 $π \equiv T / θ$ 转换为位温 $θ$ ，可得非饱和条件下的虚位温定义：

$θ v = θ 1 + r v / ϵ 1 + r v ≈ θ 1 + 0.61 r v$ (2)

为了更好地理解热力和水汽对虚位温的作用，可构建虚位温的扰动 $θ_{v}^{'} = θ_{v} - \bar{θ_{v}}$ ，其中，上划线代表雷诺平均，也可认为是日平均.若忽略二阶扰动量，则 $θ_{v}^{'}$ 可写成：

$θ v' ≈ θ' 1 + 0.61 r v ¯ + 0.61 θ ¯ r v'$ (3)

在大气中 $β θ_{v}^{'}$ 代表浮力，其中， $β = g / Θ$ 为热力膨胀系数.在大气边界层中，由于 $θ$ 的变化很小， $β$ 几乎是常数 $9.81 / 300$ m·s^-2·K^-1，因此可以认为 $θ_{v}^{'}$ 就代表浮力.需要指出的是，这里浮力指阿基米德浮力，即静力平衡条件下的浮力.在加速条件下，尤其是深对流过程中，浮力的计算需要额外考虑气压的作用，即有效浮力^［12］，这不在本文的考虑范围内.式（3）右边的两项分别代表了热力扰动 $θ^{'}$ 和水汽扰动 $r_{v}^{'}$ 对于浮力 $θ_{v}^{'}$ 的作用.在大陆边界层中，水汽含量一般较小， $θ_{v}^{'}$ 主要由 $θ^{'}$ 主导.但在高土壤湿度条件下， $r_{v}^{'}$ 对 $θ_{v}^{'}$ 的贡献则会大幅提升，若取 ${\bar{θ}}_{v} = 300 K$ ，那么 $r_{v}^{'} = 10 g \cdot k g^{- 1}$ 比湿扰动就等价于 $θ^{'} = 1.8 K$ 的位温扰动.在海洋边界层中，海水的比热容大，海温的日变化小.此时，水汽对于浮力的贡献则占主导地位，驱动并维持海洋边界层中浅积云对流（简称“浅对流”）.

基于类似式（3）的雷诺分解，可以得到垂直虚位温湍流通量 $\bar{w^{'} θ_{v}^{'}}$ （或垂直浮力通量 $β \bar{w^{'} θ_{v}^{'}}$ ）的数学表达式：

$w' θ v' ¯ ≈ w' θ' ¯ + 0.61 r v ¯ w' θ' ¯ + 0.61 θ ¯ w' r v' ¯$ (4)

式（4）中的约等号是因为忽略了三阶小量 $\bar{w^{'} θ^{'} r_{v}^{'}}$ 的贡献.类似浮力扰动，垂直浮力通量同样可以分解为感热通量 $\bar{w^{'} θ^{'}}$ 和水汽通量 $\bar{w^{'} r_{v}^{'}}$ 的贡献.在干边界层中，浮力通量由感热通量主导；在湿边界层中，潜热通量也会贡献相当部分的浮力通量.这里的“干”和“湿”采用了边界层气象中的常用表达方式.干边界层可以认为是式（4）右边最后一项远小于 $\bar{w^{'} θ_{v}^{'}}$ 的情况，即水汽对于浮力的作用很小，可视为被动标量，而并非指水汽含量为零.“湿”则代表水汽的作用不可忽略，如海洋边界层.图2给出了典型海洋浅积云边界层（Shallow Cumulus⁃Topped Boundary Layer，ScTBL）中的浮力通量及其组成分量的垂直廓线.如图所示，水汽通量项在边界层整层内都对浮力产生显著的正贡献，是浮力通量的重要组成部分.平均水汽对于感热通量的修正（即式（4）右边第二项）对浮力通量的贡献几乎则为零.由此可见，水汽对于浮力通量的贡献主要体现在水汽通量本身（即式（4）右边第三项）.

图2

图2 海洋浅积云边界层中的垂直浮力通量及其组成分量的垂直廓线（数据源自对Barbados Oceanographic and Meteorological Experiment (BOMEX)个例的大涡模拟^[13]）

Fig.2 Vertical profiles of the vertical buoyancy flux and its constituents in a typical shallow cumulus⁃topped boundary layer (Data is obtained from large⁃eddy simulation of the Barbados Oceanographic and Meteorological Experiment (BOMEX)^[13])

图1和图2中可以注意到 $\bar{w^{'} θ_{v}^{'}}$ 在边界层内的线性特征，即：

$w' θ v' ¯ z B s = α θ v - 1 z z i + 1$ (5)

其中， $B_{s} = \bar{w^{'} θ_{v_{s}}^{'}}$ 为地表浮力通量， $α_{θ_{v}} = \bar{w^{'} θ_{v_{e}}^{'}} / B_{s}$ 为夹卷通量与地表通量之比，下标 $s$ 和 $e$ 分别代表地表（Surface）和夹卷（Entrainment）.事实上，不止是 $\bar{w^{'} θ_{v}^{'}}$ ，对流边界层中的绝大多数垂直湍流通量在垂直方向都表现出强烈的线性特征，这主要是由快速的湍流混合所造成的^［14］.以 $θ_{v}$ 为例，在水平均匀条件下， $\bar{w^{'} θ_{v}^{'}}$ 的垂直散度是驱动边界层 $\bar{θ_{v}}$ 变化的主要原因，即：

$∂ θ v ¯ ∂ t = - ∂ w' θ v' ¯ ∂ z$ (6)

在日间加热条件下，边界层内的 $\bar{θ_{v}}$ 虽然不断增高，但由于强烈的湍流混合，致使 $\bar{θ_{v}}$ 的垂直梯度 $\partial_{z} \bar{θ_{v}}$ （即其廓线形态）基本不随时间变化，也就是 $\partial_{t} (\partial_{z} \bar{θ_{v}}) \approx$ 0，这也被称为边界层的准定常态，即 $\partial_{t} \bar{θ_{v}} \neq 0$ ， $\partial_{t} (\partial_{z} \bar{θ_{v}}) = 0$ ，这里采用下标对导数进行简写.在准定常条件下，虽然边界层的热源来自地面和夹卷层的加热，但在湍流的作用下，源自上下两端的局地加热能快速均匀地分散到整个边界层中，实现边界层各个高度的均匀加热，也就是式（6）右边为常数，即式（5）中的线性通量.

虽然都具备线性特征，但是和其他湍流通量相比，浮力通量具有独特的相似性特征，即式（5）中的 $α_{θ_{v}}$ 约为-0.2，是一个在不同地表热力和水汽条件下都相对固定的值^［15-16］.尽管 $α_{θ_{v}}$ 会随着夹卷层风切变的增强（确切地说是理查德森数的增大）而变得更负，但其变化范围要远小于其他湍流通量^［17-18］，理查德森数为湍动能的浮力与机械生成项的比值.图3给出了六个理想对流边界层个例的感热通量 $\bar{w^{'} θ^{'}}$ 和浮力通量 $\bar{w^{'} θ_{v}^{'}}$ 的垂直廓线.这六个个例受不同的地表感热 $\bar{w^{'} θ_{s}^{'}^{}}$ 和水汽通量 $\bar{w^{'} r_{v_{s}}^{'}}$ 驱动，形成了从0.02到 $+ \infty$ 的六个波文比（Bowen ratio，B）， $B \equiv (c_{p} \bar{w^{'} θ_{s}^{'}^{}}) / (L_{v} \bar{w^{'} r_{v_{s}}^{'}})$ ， $c_{p}$ 和 $L_{v}$ 分别为空气的感热容和潜热系数.如图3所示，标准化后的浮力通量 $\bar{w^{'} θ_{v}^{'}} / \bar{w^{'} θ_{v_{s}}^{'}}$ 几乎实现了不同个例间的重叠，而标准化后的感热通量 $\bar{w^{'} θ^{'}} / \bar{w^{'} θ_{s}^{'}}$ 在不同波文比条件下，则展现出显著的不同.关于浮力通量的相似性特征及 $α_{θ_{v}}$ 的鲁棒性，至今尚缺乏合理的物理解释，这有可能与边界层热对流的结构和尺度相关^［19］.

图3

图3 六个具有不同波文比的理想对流边界层中的(a)感热通量 $\bar{w^{'} θ^{'}}$ 和(b)浮力通量 $\bar{w^{'} θ_{v}^{'}}$ 垂直廓线

Fig.3 Vertical profiles of the (a) sensible heat fluxes and (b) buoyancy fluxes from 6 idealized CBL cases with different Bowen ratios

在数值天气预报模式和全球气候模式中，大气边界层的湍流混合过程由行星边界层（Planetary Boundary Layer，PBL）参数化方案（简称边界层方案）承担.确切地说，传统的边界层方案主要为模式提供边界层内的垂直湍流动量通量、热量通量、水汽通量，以及云水通量（主要用于处理边界层雾和层积云）；同时边界层方案也计算自由大气中由“干”湍流引起的垂直通量，即不涉及云的部分，云中湍流则由传统的积云方案（Cumulus Scheme）进行参数化.

对于热通量而言，几乎所有的边界层方案计算的都是感热通量 $\bar{w^{'} θ^{'}}$ ，而非浮力通量 $\bar{w^{'} θ_{v}^{'}}$ 或其他热力学变量的通量，比如温度通量 $\bar{w^{'} T^{'}}$ 、相当位温通量 $\bar{w^{'} θ_{e}^{'}}$ 、液态水位温通量 $\bar{w^{'} θ_{l}^{'}}$ 等.这很可能是因为数值预报模式大都采用包含了绝热加热/冷却效果的位温 $θ$ 作为热力预报变量，自然需要边界层方案提供相应的位温方程湍流倾向项.同时，在经典的干边界层中，位温 $θ$ 也是众多的热力学守恒变量（如 $θ_{v}, θ_{e}, θ_{l}, \dots$ ）中最简单易用的一个，主要受温度控制，不受水汽的直接影响.事实上，绝大多数现有的边界层方案在设计之初，对热通量的参数化都是针对干边界层中的位温湍流混合过程来设计，但不考虑水汽对浮力的作用.那么，针对干边界层位温所设计的热力混合是否适用于湿边界层？对于这一问题，前人尚未开展过相关的理论和评估工作，本文将针对这一问题开展研究.

1 理论

1.1　传统感热通量的修正

因非局地热对流的存在，经典的湍流梯度扩散假设对边界层热通量并不适用，这一点从图1可以看出.梯度扩散假设湍流热通量顺着位温梯度，从高位温向低位温传输，即：

\bar{w^{'} θ^{'}} = - K_{H} \frac{\partial \bar{θ}}{\partial z}

(7)

其中， $K_{H}$ 为湍流热交换系数.因为对流边界层受夹卷加热的作用，混合层中上半部分的位温廓线呈弱稳定态，依据式（7），感热通量在 $~ 0.4 z_{i}$ 以上为负，代表热量下传，这显然与实际不符.不管是观测^［15］、实验^［20］还是数值模拟^［21］都指出，在干边界层中，感热通量在 $~ 0.8 z_{i}$ 以下皆为正，即热量上传，这主要是由非局地热对流引起的.因此，Ertel^［22］，Priestley and Swinbank^［23］和Deardorff^［24］提出了针对非局地热通量的逆梯度修正项 $γ_{c g}$ ，将式（7）修正为：

\bar{w^{'} θ^{'}} = - K_{H} (\frac{\partial \bar{θ}}{\partial z} - γ_{c g})

(8)

其中， $γ_{c g} (z)$ 定义在对流边界层内，且始终为正.对比式（7），式（8）右边新增的 $K_{H} γ_{c g}$ 保证了混合层中的正感热通量，因而被边界层方案广泛采用^［10］.

为了确定 $γ_{c g}$ ，Deardorff^［8］提出基于 $\bar{w^{'} θ^{'}}$ 诊断方程的 $γ_{c g}$ 表达式，并被后人不断改善^［25-26］，Holtslag and Moeng^［27］（简称“HM91”）提出相对完善的推导，以下给出简要介绍.从水平均匀环境下的Bousinessq近似 $\bar{w^{'} θ^{'}}$ 预报方程开始：

\frac{\partial \bar{w^{'} θ^{'}}}{\partial t} = \underset{T}{\underset{︸}{- \frac{\partial \bar{w^{'} w^{'} θ^{'}}}{\partial z}}} \underset{M}{\underset{︸}{- \bar{{w^{'}}^{2}} \frac{\partial \bar{θ}}{\partial z}}} \underset{B}{\underset{︸}{+ β \bar{θ^{'} θ_{v}^{'}}}} \underset{P}{\underset{︸}{- \frac{1}{ρ_{0}} \bar{\frac{θ \partial p}{\partial z}}}}

(9)

其中，右边的四项依次代表垂直湍流混合项（T）、位温梯度生成项（M）、浮力项（B）和气压相关项（P）.需要注意的是，在HM91的原文^［27］中，浮力项的表达式为 $β \bar{{θ^{'}}^{2}}$ ，而非式（9）中的 $β \bar{θ^{'} θ_{v}^{'}}$ ，这是因为HM91仅考虑了干边界层湍流，浮力扰动 $θ_{v}^{'}$ 等价于 $θ^{'}$ .如前所述，干边界层假设也被几乎所有的前人工作所采用.

之后，基于一个大涡模拟的干边界层个例，HM91给出了 $γ_{c g}$ 推导过程中的关键假设，即湍流项和梯度项之间只相差一个常数，即：

T (\frac{z}{z_{i}}) \approx P (\frac{z}{z_{i}}) + b \frac{\bar{w^{'} θ_{s}^{'}}}{z_{i} / w_{*}}, 0 < \frac{z}{z_{i}} < 1

(10)

其中， $b \approx$ 2为常数， $z_{i} / w_{*}$ 为对流特征时间.对流特征速度 $w_{*}$ 沿用了位温的定义：

w_{*}^{d} = {(β \bar{w^{'} θ_{s}^{'}} z_{i})}^{1 / 3}

(11)

而非更为普适的虚位温定义：

w_{*}^{m} = {(β \bar{w^{'} θ_{v_{s}}^{'}} z_{i})}^{1 / 3}

(

12

)

作为特征量，这里采用上标 $d$ 和 $m$ 以区分干湿定义.式（12）正确反映了对流速度受浮力 $β \bar{w^{'} θ_{v_{s}}^{'}}$ 控制的基本物理机制，对干湿边界层都适用.另外，需要特别指出，式（10）中 $b$ 系数项 $\bar{w^{'} θ_{s}^{'}} / (z_{i} / w_{*})$ 在HM91原文中写成了 $w_{*}^{2} θ_{*} / z_{i}$ 的特征量乘积形式，这里将它恢复成更一般的通量 $\bar{w^{'} θ_{s}^{'}^{}}$ 除以热对流时间尺度 $z_{i} / w_{*}$ 的形式，以避免进一步引入干、湿特征温度 $θ_{*}^{d} = \bar{w^{'} θ_{s}^{'}} / w_{*}^{d}$ ， $θ_{*}^{m} = \bar{w^{'} θ_{s}^{'}} / w_{*}^{m}$ .HM91采用式（11）的定义，因此式（10）是否适用于湿边界层尚不明确.具体地，若将位温 $θ$ 看作是虚位温 $θ_{v}$ 在干环境中的特殊情况，倘若要应用于湿边界层，式（10）中的 $w_{*}$ 是否需要采用更具物理意义的式（12）定义？若采用了式（12），是否需要对式（10）中的 $b$ 系数进行相应的修正？

而后，HM91又对式（9）中的气压项进行了如下参数化：

- \frac{1}{ρ_{0}} \frac{\bar{θ \partial p}}{\partial z} = - a β \bar{{θ^{'}}^{2}} - \frac{\bar{w^{'} θ^{'}}}{τ}

(13)

其中，右边第一项为浮力压强， $a \approx$ $1 / 2$ ^［28］；第二项为表征湍流作用的回归各向同性（Return⁃to⁃Isotropy）项，也称Rotta项^［29］， $τ$ 为湍流特征时间尺度，常取 $τ \propto e / ε$ ， $e$ 为湍动能， $ε$ 为湍流耗散率.针对湿边界层湍流，式（13）也存在同样的不确定性，即若将浮力压强项修正为 $- β \bar{θ^{'} θ_{v}^{'}}$ （见式（9）），那么式（13）是否依旧成立？系数 $a$ 是否需要调整？这里需要指出的是，系数 $a$ 的选取对于推导 $γ_{c g}$ 的表达式至关重要，若将式（10）和式（13）代入式（9），假设定常状态（即式（9）左侧为零），经过整理，则有：

\frac{\bar{w^{'} θ^{'}}}{τ} = - \frac{\bar{{w^{'}}^{2}}}{2} \frac{\partial \bar{θ}}{\partial z} + \frac{b}{2} \frac{\bar{w^{'} θ_{s}^{'}}}{z_{i} / w_{*}} + (\frac{1}{2} - a) β \bar{θ^{'} θ_{v}^{'}}

(

14

)

若 $a \neq$ $1 / 2$ ，则无法消去浮力项.若 $a$ 取 $1 / 2$ ，则可将式（14）整理成式（8）的形式，得到 $γ_{c g}$ 的最终表达式：

γ_{c g} = b \frac{\bar{w^{'} θ_{s}^{'}} w_{*}}{\bar{{w^{'}}^{2}} z_{i}}

(

15

)

因避免使用 $θ_{*}$ ，式（15）与HM91原文 $γ_{c g} = b w_{*}^{2} θ_{*} / (\bar{{w^{'}}^{2}} z_{i})$ 在形式上略有不同.依据式（8）， $γ_{c g}$ 的量纲与 $\partial \bar{θ} / \partial z$ 一致，因此可取 $\bar{w^{'} θ_{s}^{'}} / w_{*}$ 为 $\partial \bar{θ}$ 的特征值，取 $z_{i}$ 为 $\partial z$ 的特征值，进一步对式（15）进行标准化，则有：

\frac{w_{*} z_{i}}{\bar{w^{'} θ_{s}^{'}}} γ_{c g} = b \frac{w_{*}^{2}}{\bar{{w^{'}}^{2}}}

(

16

)

如前所述，式（16）中 $b$ ， $w_{*}$ 的选取需要针对湿边界层进行评估.

1.2　浮力通量

这里探索前言中的提议.若采用干、湿环境皆可用且相似性更强的垂直浮力通量 $\bar{w^{'} θ_{v}^{'}}$ 取代感热通量 $\bar{w^{'} θ^{'}}$ 作为边界层的热力学预报变量，那么相应的浮力通量逆梯度输送项又是什么？为此，可以参照HM91的推导，基于浮力通量 $\bar{w^{'} θ_{v}^{'}}$ 的收支方程：

\begin{array}{l} \frac{\partial \bar{w^{'} θ_{v}^{'}}}{\partial t} = \\ \underset{T_{v}}{\underset{︸}{- \frac{\partial \bar{w^{'} w^{'} θ_{v}^{'}}}{\partial z}}} \underset{M_{v}}{\underset{︸}{- \bar{{w^{'}}^{2}} \frac{\partial \bar{θ_{v}}}{\partial z}}} \underset{B_{v}}{\underset{︸}{+ β \bar{{θ_{v}^{'}}^{2}}}} \underset{P_{v}}{\underset{︸}{- \frac{1}{ρ_{0}} \bar{\frac{θ_{v} \partial p}{\partial z}}}} \end{array}

(17)

依次展开类似式（10）~（）的推导.由于这里考虑的是浮力通量，在推导过程中将毫无疑问地摈弃式（11）的地表感热通量，而使用式（12）的浮力通量来定义 $w_{*}$ .若将HM91的工作看作边界层的干过程特例，则理应期待获得类似式（8）的浮力通量参数化，即：

\bar{w^{'} θ_{v}^{'}} = - K_{H} (\frac{\partial {\bar{θ}}_{v}}{\partial z} - γ_{c g, v})

(18)

及类似式（15）的逆梯度浮力通量修正项 $γ_{c g, v}$ ，即：

γ_{c g, v} = b \frac{\bar{w^{'} θ_{v_{s}}^{'}} w_{*}}{\bar{{w^{'}}^{2}} z_{i}}

(

19

)

式（19）中的 $w_{*}$ 也将毫无争议地采用式（12）的湿过程定义.但式（18）和式（19）是否具有足够的干湿普适性，是否存在不依赖于干湿过程的常数 $b$ ，还有待于接下来的验证.

2 个例选择和模式设置

本文选择的湿边界层个例是浅积云经典算例Barbados Oceanographic and Meteorological Experiment （BOMEX）^［30］，模拟典型的定常条件下的海上浅积云边界层.个例设置参照Siebesma et al^［31］，关键实验参数如表1所示，其中， $(U_{g}, V_{g})$ 为地转风矢量， $V_{g}$ 恒为零， $L_{O} = - u_{*}^{3} / (κ β \bar{w^{'} θ_{v_{s}}^{'}})$ 为奥布霍夫长度， $κ$ =0.41为冯卡门常数， $z_{i} / |L_{O}|$ 为边界层整体稳定度（Bulk Stability）^［6］，一般当 $(z_{i} / |L_{O}|) < 25$ 时，边界层内的热对流组织方式为滚轴涡旋.在BOMEX个例中， $z_{i}$ 特指混合层厚度，即浮力通量到达最小值的高度.BOMEX个例的初始位温廓线包含了近地500 m的等熵层，500~1500 m高度为湿对流条件不稳定层，1500 m以上为绝对稳定层，初始风速廓线为地转风.个例设置中还加入了大尺度的沉降，结合平流及辐射冷却造成的位温倾向等环境变量，以实现定常态，具体设置详见Siebesma et al^［31］的表1.

表1 BOMEX个例和干对流边界层个例的关键实验参数设置

Table 1 Key model parameters of the BOMEX case and the dry CBL case

个例	$\begin{array}{l} \bar{w^{'} θ_{s}^{'}} \\ (K \cdot m \cdot s^{- 1}) \end{array}$	$\begin{array}{l} \bar{w^{'} r_{v_{s}}^{'}} \\ (m \cdot s^{- 1}) \end{array}$	$\begin{array}{l} \bar{w^{'} θ_{v_{s}}^{'}} \\ (K \cdot m \cdot s^{- 1}) \end{array}$	$\begin{array}{l} u_{*} \\ (m \cdot s^{- 1}) \end{array}$	$\begin{array}{l} U_{g} \\ (m \cdot s^{- 1}) \end{array}$	$\begin{array}{l} z_{i} \\ (m) \end{array}$	$\begin{array}{l} - L_{O} \\ (m) \end{array}$	$\begin{array}{l} z_{i} / \|L_{O}\| \\ (-) \end{array}$
BOMEX	8.0 $\times$ 10^-3	5.2 $\times$ 10^-5	1.8 $\times$ 10^-2	0.28	$- 10 + 1.8 \times 10^{- 3} z$	545	93	5.9
CBL	5.0 $\times$ 10^-2	0	5.0 $\times$ 10^-2	0.50	10	1007	183	5.5

新窗口打开| 下载CSV

为了对比干湿边界层中的非局地输送，本文还选取了一个具有和BOMEX个例相似整体稳定度 $(z_{i} / |L_{O}|)$ 的干边界层（CBL）个例，关键实验参数见表1.CBL个例的设置参照了Shin and Hong^［32］，个例的初始位温条件为0~925 m的等熵层，与925~1075 m的覆盖逆温层以抑制边界层的增长，其上则为自由大气.初始风速廓线为地转风.

两个个例皆采用了大涡模拟，以尽可能真实地分辨湍流混合过程.BOMEX个例的水平和垂直分辨率皆为10 m，采用 $(864,864,300)$ 个网格点，模拟 $8.64 k m \times 8.64 k m \times 3 k m$ 的三维空间.CBL个例采用水平10 m，垂直4 m的分辨率，采用 $(1008,1008,324)$ 个网格点，模拟 $10.08 k m \times 10.08 k m \times 2 k m$ 的三维空间.模式的水平方向采用双周期边界条件，地面采用预设通量作为边界条件（详见表1），模式顶则采用半滑动边界条件，并将模式顶部的1/4高度设置为海绵层，以吸收上传的重力波.BOMEX和CBL个例分别经过4 h和5 h达到定常状态，而后又运行了3 h以采集模式数据用于分析.

3 结果与讨论

图4展示了BOMEX个例中的标准化垂直感热通量 $\bar{w^{'} θ^{'}}$ （式（9））和垂直浮力通量 $\bar{w^{'} θ_{v}^{'}}$ 收支廓线（式（17）），因为本文的研究对象是边界层，y轴仅标注到 $1.2 z_{i}$ 高度， $\bar{w^{'} θ^{'}}$ 和 $\bar{w^{'} θ_{v}^{'}}$ 分别采用了 $H_{s} \equiv \bar{w^{'} θ_{s}^{'}^{}}$ 和 $B_{s} \equiv \bar{w^{'} θ_{v_{s}}^{'}}$ 进行标准化，特征时间尺度统一选取了更具物理意义的 $z_{i} / w_{*}^{m}$ .由于依据大涡模拟诊断的四个通量强迫项 $(T, M, B, P)$ 之和往往不严格为零，因此在图4中引入一个残差项 $R = - (T + M + B + P)$ .残差项不为零，代表通量收支不闭合，其原因可能是由于模式的时空平均不够充分，导致尚有部分时间倾向项残余，即 $\partial_{t} \neq 0$ .模式的数值方法，主要是平流方案带来的数值混合，也是造成收支不平衡的可能原因^［33］.最后，计算的通量收支是基于模式可分辨场的，没有考虑次网格通量的贡献.在靠近地面时，模式次网格通量增大，是造成残差项在近地层增大的主要原因.总体上，图4中的残差项在模式分辨率相对较高的混合层内相比其他项最小，接近零，说明通量收支基本闭合，结果可信.

图4

图4 BOMEX个例的(a)垂直感热通量和(b)垂直浮力通量收支（图注参照式(9)与式(17)）

Fig.4 Vertical profiles of the (a) vertical sensible heat flux and (b) vertical buoyancy flux budget， legends follows Eqs. 9 and 17

图5展示了干对流边界层CBL个例中的垂直感热通量收支，亦可视作垂直浮力通量收支（依据式（4），当 $r_{v} = 0$ 时，有 $\bar{w^{'} θ^{'}} = \bar{w^{'} θ_{v}^{'}}$ ）.对比图4和图5，可见图5中干边界层 $\bar{w^{'} θ^{'}}$ 收支与图4b中的湿边界层 $\bar{w^{'} θ_{v}^{'}}$ 收支高度一致.与图4a中的湿边界层 $\bar{w^{'} θ^{'}}$ 收支相比，则在廓线形态和量级上都有很大的不同.这一结果为浮力通量在干湿边界层中的普适性提供了支持.

图5

图5 如图4，CBL个例的垂直感热通量收支

Fig.5 As Fig.4, but for the CBL vertical sensible heat flux

接下来，依据图4和图5展示的垂直感热和浮力通量收支，对上一节中提出的问题逐一求证.首先是湍流项和气压项之间的关系，即式（10）的干湿适用性.图6a展示了不同条件下计算的BOMEX湿边界层 $b$ 系数（见式（10）），其中，蓝线和绿线针对感热通量，分别采用干 $w_{*}^{d}$ （式（11））、湿 $w_{*}^{m}$ （式（12））两种 $w_{*}$ 定义；红线针对浮力通量，采用了湿 $w_{*}^{m}$ 的定义.对于 $\bar{w^{'} θ^{'}}$ ，无论在式（10）中采用 $w_{*}^{d}$ 或 $w_{*}^{m}$ ，在混合层中所得的 $b$ 系数皆大于2，尽管采用 $w_{*}^{m}$ 时，所得的 $b = 2.8$ 更加接近HM91在干边界层中所得 $b = 2$ .而对于 $\bar{w^{'} θ_{v}^{'}}$ ，所得的 $b$ 系数约为2，不但如此， $b$ 系数取常值的有效垂直范围也要宽于 $\bar{w^{'} θ^{'}}$ ，可达 $~ 0.9 z_{i}$ 高度，而后者只有 $~ 0.8 z_{i}$ ，图6b则展示了CBL个例中的 $b$ 系数，针对的既是感热通量也是浮力通量.在干边界层中，所得的 $b = 2$ ，证实了HM91的结果.图6再次支持了 $\bar{w^{'} θ_{v}^{'}}$ 收支在干湿边界层内的一致性，即若采用浮力通量，则干湿边界层中 $b$ 系数并不变化，都等于2.

图6

图6 (a) BOMEX个例和(b) CBL个例中的无量纲 $b$ 系数（参见式(10)）

蓝线和绿线针对感热通量 $\bar{w^{'} θ^{'}}$ 分别采用了 $w_{*}^{d}$ 和 $w_{*}^{m}$ 定义，红线代表的是浮力通量 $\bar{w^{'} θ_{v}^{'}}$ ，采用 $w_{*}^{m}$ 定义.

Fig.6 Vertical profiles of the dimensionless constant $b$ in Eq. 10 for (a) the BOMEX case and (b) the CBL case

考察式（13）中气压项的参数化.首先，将式（13）中Rotta项的时间尺度 $τ$ 替换为 $e / ε$ 并引入系数 $C_{T}$ ，将式（13）变为可操作的参数化形式，即：

\begin{array}{l} - \frac{1}{ρ_{0}} \frac{\bar{θ \partial p}}{\partial z} = - a β \bar{θ^{'} θ_{v}^{'}} - C_{T}^{θ} \frac{\bar{w^{'} θ^{'}}}{e / ε} \\ - \frac{1}{ρ_{0}} \frac{\bar{θ_{v} \partial p}}{\partial z} = - a β \bar{{θ_{v}^{'}}^{2}} - C_{T}^{θ_{v}} \frac{\bar{w^{'} θ_{v}^{'}}}{e / ε} \end{array}

(20)

其中， $C_{T}^{θ}$ 和 $C_{T}^{θ_{v}}$ 分别为感热通量和浮力通量的 $R o t t a$ 系数.假设式（13）和式（18）中的 $a$ 系数等于0.5（如前文所述，这对式（14）中浮力项的消去至关重要），而后采用最小二乘法在混合层（即 $0.2 z_{i}$ 至 $0.8 z_{i}$ 的高度区间内）计算最优的 $C_{T}$ 系数，并将其代入式（20）的右边，即获得气压项的参数化.图7a中，实线代表大涡模拟诊断所得的BOMEX个例中感热通量气压项 $P_{H}$ （黑色实线）和浮力通量气压项 $P_{B}$ （红色实线）以及两者相应的参数化项.不难看出，相较于黑色虚实线，红色虚实线更加贴近，即在 $a = 0.5$ 的前提下，针对浮力通量所得最优参数化显然优于针对感热通量所得最优参数化，式（20）左右两边的均方根误差更小.图7b则给出了CBL个例中的感热即浮力通量气压项及其参数化，对于CBL气压项的参数化，同样采用最小二乘法，在固定 $a = 0.5$ 的前提下，获得最优的 $C_{T} = 2.2$ ，与BOMEX个例中浮力通量的气压Rotta项系数 $C_{T}^{θ_{v}}$ 高度一致，再次为 $\bar{w^{'} θ_{v}^{'}}$ 收支在干湿边界层内的普适性提供了有力支撑.

图7

图7 (a) BOMEX个例中的感热通量(黑线)和浮力通量(红线)收支中的气压项和

(b) CBL个例中的感热即浮力通量收支中的气压项

Fig.7 Vertical profiles of the normalized pressure correlation term for (a) the sensible heat flux (black lines) and buoyancy flux (red lines) for the BOMEX case， and (b) the sensible heat/buoyance flux for the CBL case

Solid and dashed lines represent the LES benchmark and the parameterized profiles with Eq. 20，respectively.

最后，考察逆梯度修正项 $γ_{c g}$ 的参数化，图8展示了式（16）标准化后的逆梯度修正项.依据图6，对于感热通量 $\bar{w^{'} θ^{'}}$ ，若采用干 $w_{*}^{d}$ 定义（式（11）），那么在BOMEX湿边界层中，需要将系数 $b$ 修正为4.4（即图8中的蓝线）.若依旧使用 $b = 2$ （图8绿线）则会导致逆梯度输送的低估.在这样的情况下， $\bar{w^{'} θ^{'}} (z)$ 依旧受下垫面通量 $\bar{w^{'} θ_{s}^{'}^{}}$ 控制（类似式（5）），而基本保持不变，那么依据式（8），可预见 $γ_{c g}$ 的低估将导致位温梯度 $\partial \bar{θ} / \partial z$ 的增大.针对BOMEX湿边界层中的感热通量 $\bar{w^{'} θ^{'}}$ ，若采用式（12）中湿 $w_{*}^{m}$ 定义，那么 $b$ 系数将相应的修正为2.8（即图8中的红线）；若针对浮力通量 $\bar{w^{'} θ_{v}^{'}}$ ，其 $b$ 系数依旧为2，与HM91中所得的干边界层感热通量 $b$ 系数一致.黑线代表了标准化后的 $γ_{c g, v}$ ，与采用正确系数的标准化 $γ_{c g}$ （图8中的蓝线）几乎重合.对于干边界层CBL个例，其 $γ_{c g}$ 项与图8中的黑线基本重合，这一点从图5和图6就可以推导得出，因此不再重复绘制.图8再次说明，浮力通量 $\bar{w^{'} θ_{v}^{'}}$ 的逆梯度输送项 $γ_{c g, v}$ 在干湿边界层中具有较好的系数一致性，适于作为边界层方案的热力学预报通量.

图8

图8 BOMEX个例中的标准化逆梯度修正项（式（16））及不同的 $b$ 系数

Fig.8 Vertical profiles of the normalized countergradient correction term (Eq. 16) with different $b$ coefficients

4 结论与展望

依据上一节的结果，发现垂直浮力通量 $\bar{w^{'} θ_{v}^{'}}$ 对于干、湿边界具有较好的一致性，而垂直感热通量 $\bar{w^{'} θ^{'}}$ 则不具备这样的性质.在一般日间对流条件下，即边界层整体稳定度 $z_{i} / |L_{O}| > 𝒪 (1)$ 时，不论干、湿边界层， $\bar{w^{'} θ_{v}^{'}}$ 具备极强自相似特征垂直廓线（参见式（5），图3），而 $\bar{w^{'} θ^{'}}$ 的垂直廓线则随着地表波文比而变化，不具备一致性.其次，浮力通量 ${\bar{w^{'} θ^{'}}}_{v}$ 的定常收支相比感热通量的收支具备更好的干湿一致性.

（1）浮力通量收支的湍流和气压强迫项之差（即式（10））的系数 $b$ 在干湿条件下皆为2（图6和图8），而感热通量的 $b$ 系数在湿边界层中要大于干边界层，采用浮力定义的自由对流速度 $w_{*}^{m}$ 会减小这一差距，但湿边界层中的 $b$ 系数依旧要大于2，在BOMEX个例中为2.8（见图6）.

（2）浮力通量的气压项参数化在干湿边界层中具有相同的Rotta系数 $C_{T}$ ，因而可以采用同样的参数化方案，而感热通量的Rotta系数在湿边界层中要大于干边界层（见图7）.而且基于浮力通量的气压参数化效果要略优于感热通量，具有更小的均方根误差.

由此可见，基于垂直浮力通量的边界层方案具备很好的干湿一致性，即方案关键系数不随湿度而改变，因此要优于传统的基于感热通量的边界层方案.基于上述认识，将虚位温作为边界层热力过程预报量，会比传统的位温预报具备更强的鲁棒性和更广的适用性.在程序方面，从位温到虚位温的修改，并不涉及任何对现有边界层方案的改变，只需在输入部分，利用式（2）将 $θ$ 转化为 $θ_{v}$ ，而后在输出部分，将计算所得 $\bar{w^{'} θ_{v}^{'}}$ ，结合方案预报的水汽通量 $\bar{w^{'} r_{v}^{'}}$ ，利用式（4）转换回 $\bar{w^{'} θ^{'}}$ 即可完成.这正是充分利用了浮力通量的干湿一致性，即完全可以认为传统的边界层方案就是针对浮力通量所设计的，只需输入正确的变量即 $θ_{v}$ 就可以输出正确的参数化通量即 $\bar{w^{'} θ_{v}^{'}}$ .我们已在WRF模式中实现了这一变量转换，以后将展示改进效果.

参考文献

原文顺序

文献年度倒序

文中引用次数倒序

被引期刊影响因子

[1]

Stull

R B

. An introduction to boundary layer meteorology. Boston：Kluwer Academic Publishers，1988.

[2]

Hunt

J C R

， Kaimal

J C

， Gaynor

J E

Eddy structure in the convective boundary layer⁃new measurements and new concepts

Quarterly Journal of the Royal Meteorological Society，1988，114(482)：827-858.

[3]

Salesky

S T

， Chamecki

， Bou⁃Zeid

On the nature of the transition between roll and cellular organization in the convective boundary layer

Boundary⁃Layer Meteorology，2017，163(1)：41-68.

[4]

Young

G S

Turbulence structure of the convective boundary layer

Part Ⅱ. Phonenix 78 aircraft observations of thermals and their environment. Journal of the Atmospheric Sciences，1988，45(4)：727-735.

[5]

Warner

， Telford

J W

Convection below cloud base

Journal of the Atmospheric Sciences，1967，24(4)：374-382.

[6]

Lenschow

D H

， Stephens

P L

The role of thermals in the convective boundary layer

Boundary⁃Layer Meteorology，1980，19(4)：509-532.

[7]

Townsend

A A

. The structure of turbulent shear flow. New York：Cambridge University Press，1976，429.

[8]

Deardorff

J W

Theoretical expression for the countergradient vertical heat flux

Journal of Geophysical Research，1972，77(30)：5900-5904.

[9]

X M

， Xue

， Li

X L

The use of high⁃resolution sounding data to evaluate and optimize nonlocal PBL schemes for simulating the slightly stable upper convective boundary layer

Monthly Weather Review，2019，147(10)：3825-3841.

[10]

Zhou

B W

， Sun

S W

， Yao

，et al.

Reexamining the gradient and countergradient representation of the local and nonlocal heat fluxes in the convective boundary layer

Journal of the Atmospheric Sciences，2018，75(7)：2317-2336.

[11]

赵昭，周博闻.

日间对流边界层中的非局地动量混合

气象科学，2021，41(5)：631-643.

Zhao

， Zhou

B W

Non⁃local mixing of momentum in the daytime convective boundary layer

Journal of the Meteorological Sciences，2021，41(5)：631-643.

[12]

Davies⁃Jones

An expression for effective buoyancy in surroundings with horizontal density gradients

Journal of the Atmospheric Sciences，2003，60(23)：2922-2925.

[13]

Wang

Y H

， Cheng

X P

， Fei

J F

，et al.

Modeling the shallow cumulus⁃topped boundary layer at gray zone resolutions

Journal of the Atmospheric Sciences，2022，79(9)：2435-2451.

[14]

Wyngaard

J C

. Turbulence in the atmosphere. New York：Cambridge University Press，2010，393.

[15]

LeMone

M A

， Pennell

W T

The relationship of trade wind cumulus distribution to subcloud layer fluxes and structure

Monthly Weather Review，1976，104(5)：524-539.

[16]

LeMone

M A

， Angevine

W M

， Bretherton

C S

，et al.

100 years of progress in boundary layer meteorology

Meteorological Monographs，2019，59(1)：9.1-9.85.

[17]

Conzemius

R J

， Fedorovich

. Dynamics of sheared convective boundary layer entrainment. Part I：Methodological background and large⁃eddy simu⁃lations. Journal of the Atmospheric Sciences，2006，63(4)：1151-1178.

[18]

Liu

， Sun

J N

， Shen

L D

. Parameterization of sheared entrainment in a well⁃developed CBL. Part I：Evaluation of the scheme through large⁃eddy simu⁃lations. Advances in Atmospheric Sciences，2016，33(10)：1171-1184.

[19]

De Roode

S R

， Duynkerke

P G

， Jonker

H J J

Large⁃eddy simulation：how large is large enough?

Journal of the Atmospheric Sciences，2004，61(4)：403-421.

[20]

Willis

G E

， Deardorff

J W

A laboratory model of diffusion into the convective planetary boundary layer

Quarterly Journal of the Royal Meteorological Society，1976，102(432)：427-445.

[21]

Deardorff

J W

Three⁃dimensional numerical study of the height and mean structure of a heated planetary boundary layer

Boundary⁃Layer Meteorology，1974，7(1)：81-106.

[22]

Ertel

Der vertikale turbulenz⁃wärmestrom in der atmosphäre

Meteorologische Zeitschrift，1942，59：250-253.

[23]

Priestley

C H B

， Swinbank

W C

Vertical transport of heat by turbulence in the atmosphere

Proceedings of the Royal Society of London. Series A. Mathematical and Physical Sciences，1947，189(1019)：543-561.

[24]

Deardorff

J W

The counter⁃gradient heat flux in the lower atmosphere and in the laboratory

Journal of the Atmospheric Sciences，1966，23(5)：503-506.

[25]

Therry

， Lacarrère

Improving the eddy kinetic energy model for planetary boundary layer description

Boundary⁃Layer Meteorology，1983，25(1)：63-88.

[26]

Troen

I B

， Mahrt

A simple model of the atmospheric boundary layer; sensitivity to surface evaporation

Boundary⁃Layer Meteorology，1986，37(1-2)：129-148.

[27]

Holtslag

A A M

， Moeng

C H

Eddy diffusivity and countergradient transport in the convective atmospheric boundary layer

Journal of the Atmospheric Sciences，1991，48(14)：1690-1698.

[28]

Moeng

C H

， Wyngaard

J C

An analysis of closures for pressure⁃scalar covariances in the convective boundary layer

Journal of the Atmospheric Sciences，1986，43(21)：2499-2513.

[29]

Rotta

Statistische theorie nichthomogener turbulenz

Zeitschrift für Physik，1951，129(6)：547-572.

[30]

Holland

J Z

， Rasmusson

E M

Measurements of the atmospheric mass，energy，and momentum budgets over a 500⁃kilometer square of tropical ocean

Monthly Weather Review，1973，101(1)：44-55.

[31]

Pier Siebesma

， Bretherton

C S

， Brown

，et al.

A large eddy simulation intercomparison study of shallow cumulus convection

Journal of the Atmospheric Sciences，2003，60(10)：1201-1219.

[32]

Shin

H H

， Hong

S Y

Analysis of resolved and parameterized vertical transports in convective boundary layers at gray⁃zone resolutions

Journal of the Atmospheric Sciences，2013，70(10)：3248-3261.

[33]

Heinze

， Mironov

， Raasch

Second⁃moment budgets in cloud topped boundary layers：A large⁃eddy simulation study

Journal of Advances in Modeling Earth Systems，2015，7(2)：510-536.