2.10 热交换（Heat Exchange）

      与大气的热交换基于以下四个物理过程计算：  

      • 潜热通量（蒸发导致的热损失）  

      • 显热通量（对流导致的热通量）  

      • 净短波辐射  

      • 净长波辐射  

      潜热、显热与长波辐射假定发生在水面；短波辐射的吸收剖面用 Beer 定律近似。光强衰减用修正 Beer 定律表示：

$$ I(d)=\left(1-\beta\right)I_{0}e^{-\lambda d}\tag{2.100} $$

其中：$I(d)$ 为距水面深度 $d$ 处光强；$I_0$ 为水面下方的光强；$\beta$ 表示近表面吸收（主要为红外）所占比例；$\lambda$ 为消光系数。典型取值：$\beta=0.2\sim0.6$，$\lambda=0.5\sim1.4\,\mathrm{m^{-1}}$；默认值为 $\beta=0.3$、$\lambda=1.0\,\mathrm{m^{-1}}$。近表面被吸收的能量为 $\beta I_0$。净短波辐射 $q_{sr,net}$ 按修正 Beer 定律衰减，因此表面净热通量为：

$$ Q_{n}=q_{v}+q_{c}+\beta q_{sr,net}+q_{lr,net}\tag{2.101} $$

      对三维计算，源项 $\hat{H}$ 为：

$$ \hat{H}=\frac{\partial}{\partial z}\Bigg(\frac{q_{sr,net}(1-\beta)e^{-\lambda(\eta-z)}}{\rho_{0}c_{p}}\Bigg)=\frac{q_{sr,net}(1-\beta)\lambda e^{-\lambda(\eta-z)}}{\rho_{0}c_{p}}\tag{2.102} $$

      对二维计算，源项 $\hat{H}$ 为：

$$ \hat{H}=\frac{q_{v}+q_{c}+q_{sr,net}+q_{lr,net}}{\rho_{0}c_{p}}\tag{2.103} $$

      下文分别给出潜热通量、显热通量、净短波辐射与净长波辐射的计算方法。冰覆盖区域不进行热交换计算。

      Dalton 定律给出蒸发热损失（潜热通量）的关系（见 Sahlberg, 1984）：

$$ q_{v}=LC_{e}(a_{1}+b_{1}W_{2m})\big(Q_{water}-Q_{air}\big)\tag{2.104} $$

其中：$L=2.5\times10^{6}\,\mathrm{J/kg}$ 为水的汽化潜热（文献中亦常用 $L=2.5\times10^{6}-2300\,T_{water}$）；$C_e=1.32\times10^{-3}$ 为水汽交换系数（Dalton number）；$W_{2m}$ 为海面以上 2 m 风速；$Q_{water}$ 为近表面水汽密度；$Q_{air}$ 为大气水汽密度；$a_1,b_1$ 为用户参数，默认 $a_1=0.5$、$b_1=0.9$。

      $Q_{water}$ 与 $Q_{air}$ 通常无法直接测得，可通过水汽压换算：

$$ Q_{i}=\frac{0.2167}{T_{i}+T_{K}}e_{i}\tag{2.105} $$

其中下标 $i$ 同时表示 water 与 air；$T_K=273.15\,^\circ\mathrm{K}$。近海面水汽压（假定近表面空气饱和且与水温相同）为：

$$ e_{water}=6.11e^{K}\Bigg(\Big(\frac{1}{T_{K}}-\frac{1}{T_{water}+T_{K}}\Big)\Bigg)\tag{2.106} $$

其中：$K=5418\,^\circ\mathrm{K}$。空气水汽压可由气温与相对湿度 $R$ 给出：

$$ e_{air}=R\cdot6.11e^{K}\Bigg(\Big(\frac{1}{T_{K}}-\frac{1}{T_{air}+T_{K}}\Big)\Bigg)\tag{2.107} $$

      将 $Q_{water}$ 与 $Q_{air}$ 替换为上述表达式，可将潜热通量写为：

$$ \begin{aligned}q_{v}=-P_{v}(a_{1}+b_{1}W_{2m})\left(\frac{\exp\!\left(K\left(\frac{1}{T_{K}}-\frac{1}{T_{water}+T_{K}}\right)\right)}{T_{water}+T_{K}}-\frac{R\cdot\exp\!\left(K\left(\frac{1}{T_{K}}-\frac{1}{T_{air}+T_{K}}\right)\right)}{T_{air}+T_{K}}\right)\end{aligned}\tag{2.108} $$

其中常数合并为 $P_v=4370\,\mathrm{J\cdot^{\circ}K/m^{3}}$。水面冷却过程中潜热损失影响显著，典型可达 $100\,\mathrm{W/m^{2}}$。

      2 m 风速 $W_2$ 由 10 m 风速 $W_{10}$ 推算。假设风速满足对数剖面：

$$ u(z)=\frac{u_{*}}{\kappa}\ln\left(\frac{z}{z_{0}}\right)\tag{2.109} $$

其中：$u_*$ 为风摩擦速度，$z_0$ 为海面粗糙度，$\kappa=0.4$ 为 von Kármán 常数。$u_*$ 与 $z_0$ 可由 Charnock 关系确定：

$$ z_{0}=z_{Charnock}\frac{u_{*}^{2}}{g}\tag{2.110} $$

$$ u_{*}=\frac{\kappa u(z)}{\ln\left(\frac{z}{z_{0}}\right)}\tag{2.111} $$

其中：$z_{Charnock}$ 为 Charnock 参数，默认 $z_{Charnock}=0.014$。实际计算中用 $z=10\,\mathrm{m}$、$u(z)=W_{10}$ 迭代求解 $z_0$，再得到 $W_2$：

$$ \begin{aligned}W_{2}=W_{10}\frac{\ln\left(\frac{2}{z_{0}}\right)}{\ln\left(\frac{10}{z_{0}}\right)}\quad& W_{10}>0.5\,\mathrm{m/s}\\ W_{2}=W_{10}\quad& W_{10}\le0.5\,\mathrm{m/s}\end{aligned}\tag{2.112} $$

      蒸发热损失既包含风致强迫对流，也包含自由对流；自由对流效应通过式 (2.104) 中参数 $a_1$ 体现。另引入临界风速 $W_{critical}$，使式 (2.112) 使用的风速取 $W_{10}=\max(W_{10},W_{critical})$；默认 $W_{critical}=2\,\mathrm{m/s}$。

      显热通量 $q_c$（$\mathrm{W/m^{2}}$）取决于海气边界层类型。通常边界层为湍流，可用下式表示：

$$ q_{c}=\left\{\begin{aligned}&\rho_{air}c_{air}c_{heating}W_{10}(T_{air}-T_{water})&&T_{air}\ge T_{water}\\&\rho_{air}c_{air}c_{cooling}W_{10}(T_{air}-T_{water})&&T_{air}<T_{water}\end{aligned}\right.\tag{2.113} $$

其中：$\rho_{air}=1.225\,\mathrm{kg/m^{3}}$ 为空气密度；$c_{air}=1007\,\mathrm{J/(kg\cdot^{\circ}K)}$ 为空气定压比热；$c_{heating}=0.0011$、$c_{cooling}=0.0011$ 分别为加热与冷却情况下的显热交换系数（Stanton number）（见 Kantha and Clayson, 2000）；$W_{10}$ 为 10 m 风速；$T_{water}$ 为海表温度；$T_{air}$ 为空气温度。显热通量通常在 $0\sim100\,\mathrm{W/m^{2}}$。同样引入临界风速：$W_{10}=\max(W_{10},W_{critical})$，默认 $W_{critical}=2\,\mathrm{m/s}$。

      太阳辐射包含波长约 1000–30000 Å 的电磁波，其中大部分在臭氧层被吸收，抵达地表的谱分布也会因穿过大气而改变。地表短波辐射主要集中在 4000–9000 Å。其强度与日地距离、赤纬角、纬度、大气外辐射、云量及大气水汽含量有关（见 Iqbal, 1983）。

      太阳轨道离心率因子（eccentricity）$E_0$ 为：

$$ \begin{aligned}E_{0}=\left(\frac{r_{0}}{r}\right)^{2}=1.000110+0.034221\cos(\Gamma)+0.001280\sin(\Gamma)+0.000719\cos(2\Gamma)+0.000077\sin(2\Gamma)\end{aligned}\tag{2.114} $$

其中：$r_0$ 为日地平均距离，$r$ 为实际距离；日角 $\Gamma$（rad）定义为：

$$ \Gamma=\frac{2\pi(d_{n}-1)}{365}\tag{2.115} $$

$d_n$ 为一年中的儒略日。

      季节变化主要由赤纬角 $\delta$（rad）控制，可用下式近似：

$$ \begin{aligned}\delta=0.006918-0.399912\cos(\Gamma)+0.07257\sin(\Gamma)-0.006758\cos(2\Gamma)+0.000907\sin(2\Gamma)-0.002697\cos(3\Gamma)+0.00148\sin(3\Gamma)\end{aligned}\tag{2.116} $$

      白昼长度 $n_d$ 随 $\delta$ 变化。给定纬度 $\phi$（北半球为正），有：

$$ n_{d}=\frac{24}{\pi}\arccos\left(-\tan(\phi)\tan(\delta)\right)\tag{2.117} $$

日出角 $\omega_{sr}$ 与日落角 $\omega_{ss}$（rad）为：

$$ \omega_{sr}=\arccos\left(-\tan(\phi)\tan(\delta)\right)\quad\text{and}\quad\omega_{ss}=-\omega_{sr}\tag{2.118} $$

      将一天内的太阳入射角变化积分，可得大气外日总短波辐射 $H_0$（$\mathrm{MJ/m^{2}/day}$）：

$$ H_{0}=\frac{24}{\pi}q_{sc}E_{0}\cos(\phi)\cos(\delta)\Big(\sin(\omega_{sr})-\omega_{sr}\cos(\omega_{sr})\Big)\tag{2.119} $$

其中：$q_{sc}=4.9212\,\mathrm{MJ/(m^{2}\cdot h)}$ 为太阳常数。

      多云条件下日总辐射 $H$（$\mathrm{MJ/m^{2}/day}$）采用：

$$ \frac{H}{H_{0}}=a_{2}+b_{2}\frac{n}{n_{d}}\tag{2.120} $$

其中：$n$ 为日照小时数；$a_2,b_2$ 为用户参数，默认 $a_2=0.295$、$b_2=0.371$。用户给定的晴空系数对应于 $n/n_d$。于是太阳短波辐射通量 $q_s$（$\mathrm{W/m^{2}}$）可写为：

$$ q_{s}=\left(\frac{H}{H_{0}}\right)q_{0}\left(a_{3}+b_{3}\cos(\omega_{i})\right)\frac{10^{6}}{3600}\tag{2.121} $$

其中：

$$ a_{3}=0.4090+0.5016\sin\left(\omega_{sr}-\frac{\pi}{3}\right)\tag{2.122} $$

$$ b_{3}=0.6609+0.4767\sin\left(\omega_{sr}-\frac{\pi}{3}\right)\tag{2.123} $$

大气外小时辐射强度 $q_0$（$\mathrm{MJ/m^{2}/h}$）与时角 $\omega_i$ 为：

$$ q_{0}=q_{sc}E_{0}\Bigg(\sin(\phi)\sin(\delta)+\frac{24}{\pi}\cos(\phi)\cos(\delta)\cos(\omega_{i})\Bigg)\tag{2.124} $$

$$ \omega_{i}=\frac{\pi}{12}\Bigg(12+\Delta t_{displacement}+\frac{4}{60}(L_{S}-L_{E})-\frac{E_{t}}{60}-t_{local}\Bigg)\tag{2.125} $$

其中：$\Delta t_{displacement}$ 为夏令时位移小时数，$L_S$ 为时区标准经度（deg），二者为用户参数（默认 $\Delta t_{displacement}=0\,\mathrm{h}$、$L_S=0\,\mathrm{deg}$）；$L_E$ 为当地经度（deg）；$E_t$（s）为由于地球轨道导致的时间差，计算为：

$$ E_{t}=\left(0.000075+0.001868\cos(\Gamma)-0.032077\sin(\Gamma)-0.014615\cos(2\Gamma)-0.04089\sin(2\Gamma)\right)\cdot229.18\tag{2.126} $$

$t_{local}$ 为当地时间（小时）。

      入射到海面的太阳辐射并非全部进入水体，部分被反射回去，称为反照率（albedo）。对平滑海面，反射系数 $\alpha$ 可表示为：

$$ \alpha=\frac{1}{2}\Bigg(\frac{\sin^{2}(i-r)}{\sin^{2}(i+r)}+\frac{\tan^{2}(i-r)}{\tan^{2}(i+r)}\Bigg)\tag{2.127} $$

其中：$i$ 为入射角，$r$ 为折射角；$\alpha$ 通常在 5%–40% 之间。$\alpha$ 亦可用下式近似：

$$ \alpha=\left\{\begin{aligned}&0.48\frac{altitude}{5}&&altitude<5\\&0.48-\frac{altitude-5}{25}(0.48-0.05)&&5\le altitude\le30\\&0.05&&altitude>30\end{aligned}\right.\tag{2.128} $$

其中太阳高度角（deg）为：

$$ altitude=90-\frac{180}{\pi}\arccos\Big(\sin(\delta)\sin(\phi)+\cos(\delta)\cos(\phi)\cos(\omega_{i})\Big)\tag{2.129} $$

      因此净短波辐射（$\mathrm{W/m^{2}}$）可表示为：

$$ q_{sr,net}=(1-\alpha)q_{s}\tag{2.130} $$

      净短波辐射 $q_{sr,net}$ 可按上述经验方法计算；也可由用户直接指定 $q_s$ 并用式 (2.130) 求得 $q_{sr,net}$，或直接给定 $q_{sr,net}$。

      任意物体/表面都会在全波段发射电磁能量。长波辐射主要位于 9000–25000 Å，为红外辐射，来源包括大气与海面。海面向大气发射的长波辐射减去大气向海面下行的长波辐射称为净长波辐射，受云量、气温、空气水汽压与相对湿度影响。净向外长波辐射 $q_{lr,net}$（$\mathrm{W/m^{2}}$）可由 Brunt 公式给出（见 Lind and Falkenmark, 1972）：

$$ q_{lr,net}=-\sigma_{sb}(T_{air}+T_{K})^{4}\Big(a-b\sqrt{e_{d}}\Big)\Bigg(c+d\frac{n}{n_{d}}\Bigg)\tag{2.131} $$

其中：$e_d$ 为露点温度对应的水汽压（mb）；$n$ 为日照小时数；$n_d$ 为最大日照小时数；$\sigma_{sb}=5.6697\times10^{-8}\,\mathrm{W/(m^{2}\cdot^{\circ}K^{4})}$ 为 Stefan–Boltzmann 常数；$T_{air}(^\circ\mathrm{C})$ 为空气温度。系数为：

$$ a=0.56,\quad b=0.077\,\mathrm{mb^{-1/2}},\quad c=0.10,\quad d=0.90\tag{2.132} $$

      水汽压按下式确定：

$$ e_{d}=10\cdot R\,e_{saturated}\tag{2.133} $$

其中：$R$ 为相对湿度；饱和水汽压 $e_{saturated}$（kPa）在 -51 至 52 °C 区间可近似为：

$$ \begin{aligned}e_{saturated}=3.38639\left(\left(7.38\times10^{-3}T_{air}+0.8072\right)^{8}-1.9\times10^{-5}\left|1.8T_{air}+48\right|+1.316\times10^{-3}\right)\end{aligned}\tag{2.134} $$

      净长波辐射也可不采用上述经验式，而表示为：

$$ q_{lr,net}=q_{ar,net}-q_{br}\tag{2.135} $$

其中：$q_{ar,net}$ 为用户指定的大气入射净辐射；回辐射 $q_{br}$ 为：

$$ q_{br}=(1-r)\varepsilon\sigma_{sb}T_{K}^{4}\tag{2.136} $$

其中：$r=0.03$ 为反射系数，$\varepsilon=0.985$ 为大气发射率（emissivity）。净长波辐射亦可由用户直接指定。