引言

在干燥器、冷水塔、空调和采暖通风等工程设计中，湿球温度是常用的重要气象参数^[1-2]。气象台站以前采用干、湿球温度表进行空气湿度的观测，可以方便地获得湿球温度资料；近年来，由于自动气象站的广泛使用，空气湿度的干、湿球观测方法逐步被自动传感器观测取代，但通过传感器很难直接测量湿球温度^[3]，这给工程设计中湿球温度参数分析和湿球资料的连续应用造成困难，如何获取较为准确的湿球温度资料，成为这类工程设计中亟待解决的问题。荣剑文^[3]研究提出热力学上湿球温度的定义是在定压绝热条件下空气与水直接接触达到稳定热湿平衡时的绝热饱和温度，根据定压绝热条件下湿空气焓值不变的原理，通过水蒸汽饱和压力公式和湿空气焓值计算公式，结合计算机编程来计算工程应用中的湿球温度。宋矗臻等^[4]的研究中指出湿球温度计测得的湿球温度并非是湿空气的绝热饱和温度，吴俊云等^[5]也指出在标准大气压下, 当空气温度超过0 ℃, 干、湿球温度之差小于11℃时, 湿球温度与热力学湿球温度相差小于0.25℃，在工程上可以用等焓线代替等湿球温度线。因此利用湿空气的焓值来计算湿球温度存在一定的局限性，其精度也存在不确定性。赵永胜等^[6]提出了利用大气压强、干球温度、相对湿度和平均风速资料，根据有关湿球温度的计算公式，通过计算机编程来计算无资料地区湿球温度的方法，并对不同风速条件下的湿球温度的计算结果进行了比较，其精度都在1%以内，说明迭代计算法更为精确可靠。关于自动站湿球温度的计算方法，国内发表的研究论文也鲜有介绍。本文根据气象观测用人工观测和自动观测的湿度参量计算公式，利用赤壁市气象站1995—2004年人工观测的资料，研究湿球温度迭代计算的简化方法，并用EXECL电子表格实现自动气象站湿球温度的快速计算。

1 空气湿度测量原理 1.1 干、湿球法测湿原理

干、湿球法测湿的原理是在温度表的水银球体包上脱脂纱布，纱布的下端浸入盛水的容器中，纱布在毛细管作用下经常处于湿润状态，此温度表称为湿球^[7]。湿球纱布中的水分必然向空气中蒸发，即在湿球与通过湿球的空气之间发生湿交换。水的蒸发量与空气中的水汽压平衡并使湿球温度维持在一定的数值。如果再用一支温度表测量当时的气温(此温度表称为干球)，就可以利用干、湿球温度的差值和其他测量条件来计算空气中的水汽压和相对湿度。

1.2 湿敏电容传感器测湿原理

常用的湿敏电容湿度传感器是用有机高分子膜作介质的一种小型电容器^[7]。湿敏电容器上电极是一层多孔金属膜，能透过水汽；下电极为一对刀状或梳状电极，引线由下电极引出。整个感应器是由两个小电容器串联组成传感器置于大气中，当大气中水汽透过上电极进入介电层，介电层吸收水汽后，介电系数发生变化，导致电容器电容量发生变化，而电容量的变化正比于相对湿度，通过这个特性来计算空气中的相对湿度和水汽压。

2 有关湿度计算公式 2.1 干、湿球法

用干、湿球温度求空气中水汽压的计算公式^[7]:

$e = {E_{tw}} - A{P_h}\left( {t - {t_w}} \right)$

(1)

2.2 湿敏电容法

用湿敏电容直接测得相对湿度时，空气中水汽压的计算公式^[7]:

$e = U \times {E_t}/100$

(2)

2.3 饱和水汽压

纯水平液面饱和水汽压采用世界气象组织(WMO)推荐的戈夫-格雷奇(Goff-Gratch)公式^[7]。

$\begin{array}{l} \log {E_t} = 10.79574\left( {1 - {T_1}/T} \right) - 5.028\log \left( {T/{T_1}} \right) + \\ \quad \quad \quad 1.50475 \times {10^{ - 4}}\left[ {1 - {{10}^{ - 8.2969\left( {T/{T_1} - 1} \right)}}} \right] + \\ \quad \quad \quad 0.42873 \times {10^{ - 3}}\left[ {{{10}^{4.76955\left( {1 - {T_1}/T} \right)}} - 1} \right] + \\ \quad \quad \quad 0.78614 \end{array}$

(3)

式(1)~(3) 中，e：水汽压(hPa)；t：干球温度(℃)；t_w：湿球温度(℃)；U：相对湿度(%)；P_h：本站气压(hPa)；E_tw：湿球温度t_w所对应的纯水平液面的饱和水汽压(hPa)，湿球结冰且湿球温度低于0℃时，为纯水平冰面的饱和水汽压；E_t：干球温度t时的纯水平液面饱和水汽压(hPa)；A：干湿表系数(℃^-1)，在百叶箱自然通风情况下A值取0.0007947；T₁=273.16 K(水的三相点温度)，T=273.15+t(绝对温度K)。

3 湿球温度的计算方法

通过自动气象站的气温和相对湿度资料，用式(2) 和式(3) 可以计算出水汽压；通过人工观测的干球温度、湿球温度和气压，用式(1) 和式(3) 也可以计算出水汽压。如果将自动气象站的气温、相对湿度、气压和计算出的水汽压代入式(1) 和式(3)，可以得到等同于人工观测的湿球温度t_w的方程，对方程进行求解，就可以计算出湿球温度值。但式(3) 中同时包含有湿球温度指数和对数的复合函数，方程很难直接一次求解^[6]。在工程计算中，这类可收敛的复杂方程(已证)一般采用迭代法求解^[8-9]，就是不断用变量的旧值递推新值来对方程进行求解的方法，但是利用迭代法求解，需要解决迭代变量、迭代关系式和迭代过程控制三个方面的问题^[10]。

3.1 迭代变量

通过推导得到的方程是关于湿球温度t_w的方程，故迭代变量选用湿球温度t_w。迭代法进行推算时，要先给迭代变量赋予一个初始值，而初始值的选用对迭代过程的计算至关重要。借鉴相对湿度经验计算的研究^[11-13]，利用式(1) 和式(2)，在不考虑气压变化的影响下，湿球温度t_w也可以表示成气温t、相对温度U函数的乘积形式，根据泰勒公式，函数可以近似用泰勒多项展开式表示^[14]，如果将泰勒多项简化成二项式，则湿球温度t_w可以简化成下面的形式:

$\begin{array}{l} {t_w} = {a_1} + {a_2}t + {a_3}{t^2} + \left( {{a_4} + {a_5}t + {a_6}{t^2}} \right)U + \\ \quad \quad \left( {{a_7} + {a_8}t + {a_9}{t^2}} \right){U^2} \end{array}$

(4)

式(4) 中a₁~a₉为参数，可以根据实测的气温和相对湿度资料对参数进行求算。本文利用赤壁市气象站2001—2003年的人工观测资料进行计算，得到了不同值下的多组参数。对相近的参数值求平均，可以得到湿球温度的经验公式一般表达式为：

$\begin{array}{l} {t_w} = - 5.806 + 0.672t - 0.006{t^2} + \left( {0.061 + } \right.\\ \quad \quad \left. {0.004t + 0.000099{t^2}} \right)U + \left( { - 0.000033 - } \right.\quad \\ \quad \quad \left. {0.000005t - 0.0000001{t^2}} \right){U^2} \end{array}$

(5)

通过式(5) 可以直接计算出湿球温度，其湿球温度计算值与真值之间存在一定的误差。经过回代验算，在气温10~20℃，相对湿度为80%~90%时，式(5) 的计算误差稳定在1%左右，而其他情况下，误差很不稳定，误差有时超过10%。但是在迭代关系式计算中，仍可以用式(5) 计算得到的湿球温度作为迭代变量的初始值。

3.2 迭代关系式

迭代关系式是解决迭代问题的关键。湿球温度t_w的方程是含有t_w指数和对数的复合函数，直接推导其迭代关系式非常困难。分析迭代计算过程可以发现，如果将公式(4) 计算的湿球温度初始值t_w0和气温、气压一起代入式(3) 和式(1) 可以计算出水汽压e₀，假设e₀和自动气象站观测得到的水汽压e相同，即Δe₀=(e₀-e)=0，则初始值t_w0就是给定气温、湿度和气压条件下百叶箱内的湿球温度；如果不同，则两者之间必然存在差异Δe₀，如果将差异经过迭代关系式，代入计算湿球温度的新值t_wn，那么当Δe_n的差异趋向于0，则t_wn趋向于湿球温度t_w。利用式(1) 可以得到t_w的方程：

${t_w} = t - 1/A{P_h}\left( {{E_{{t_w}}} - e} \right)$

(6)

根据牛顿迭代法的基本思想，其迭代公式的一般形式可以表示成：X_n+1=X_n-t_nH_ng_n^[8]。其中t_n为最优步长因子，它表示通过从点X_n出发，沿-H_ng_n方向作直线搜索确定X的值，当H=G^-1时(G为g函数的系数矩阵)，又称为阻尼牛顿法的迭代公式。考虑湿球温度的计算精度，将最优步长因子t_n确定为0.5, 根据式(6) 可以构造湿球温度t_w迭代关系式为：

${t_{{w_{n + 1}}}} = {t_{{w_n}}} - \Delta {e_n} \times A{P_h}/2\quad \left( {n = 0,1,2, \cdots } \right)$

(7)

利用迭代关系式(7) 计算得到的t_wn值，经过反复代入式(3) 和式(1) 可以计算出Δe_n=e_n-e (n=0, 1, 2, …)，当Δe_n趋近于0时，则t_wn逼近t_w。在上述迭代计算过程中，可以发现迭代尺度是随Δe_n而变化的，在差异较大时收敛速度较快，在差异较小时收敛速度较慢。

3.3 迭代过程控制

在什么时候结束迭代过程, 是迭代算法必需考虑的问题。气象观测的湿球温度和水汽压精度均为0.1，所以当Δe_n的绝对值小于0.05且Δt_wn+1=t_wn+1-t_wn (n=0, 1, 2, …)的绝对值也小于0.05时，则可以认为t_wn+1为符合气象观测精度要求的湿球温度值。

4 计算分析

为了避免复杂的编程调试困难，根据上面的分析和公式，本文采用EXECL电子表格对计算过程进行实时控制。设计的湿球温度计算电子表格包含有气压、干球、湿球、相对湿度、水汽压、初始湿球温度t_w0、水汽压差Δe_n、湿球温度差Δt_wn+1(n=0，1，2，…)和湿球温度计算值t_w等多个工作表，其中：气压、干球、湿球、相对湿度、水汽压工作表用于读入相应要素的气象资料，而初始的湿球温度t_w0、Δe_n、Δt_wn+1和t_w用于存放迭代计算中间值和结果，同时对Δe_n、Δt_wn+1工作表设计绝对值不小于0.05的红色显示条件控制，当红色显示消失时可以结束计算过程。

利用编制好的湿球温度计算电子表格，本文对赤壁市气象站1995—2004年三次人工观测的资料进行了验算，表 1给出了部分计算实例的过程值和结果。

结果表明：在不考虑湿球结冰情况下，采用上述简单迭代法计算的赤壁市气象站1995—2004年湿球温度值与人工观测值在湿球温度值观测要求精度下基本一致，这说明计算方法完全可行，迭代过程控制精确；对湿球温度实例计算过程分析，从初始到迭代过程结束其迭代步骤只用了5步，说明湿球湿度初始值选用可靠，迭代关系式选用恰当。利用电子表格又对通山县气象站2004—2008年自动气象站24小时资料进行了湿球温度计算，其迭代步骤最多也只用了6步。

5 结语

综上所述，采用迭代法计算湿球温度与传统的湿空气焓值计算法相比，计算结果精度稳定；而采用了一般牛顿迭代法的阻尼公式，计算方法简单；而引入了湿球温度的经验公式计算值作为初始值，使计算量大大减小；利用EXECL进行计算，避免了复杂的计算机编程，计算速度快且可控性好。可广泛应用于自动气象站的湿球温度计算，其他工程设计中不同通风条件下的湿球温度计算也可参考借鉴。