极端热湿地区围护结构热湿耦合传递模型

罗戴维，刘加平，刘大龙

(西安建筑科技大学 建筑学院, 西安 710055)

极端热湿地区常年处在高温高湿状态，这种气候条件会给围护结构隔热隔湿带来极大负担。以南海西沙群岛为例，其日照丰富，太阳辐射强烈，年平均气温26～27 ℃，6月份平均气温为28.9 ℃，一月份平均气温为22.9 ℃，总降雨量约为1 500 mm。南海全年气温高、温差小，雨量丰沛、时空分布不均，强风、台风、暴雨、大雾等灾害性天气频发，是北半球天气气候变化最敏感的地区之一，属于典型的极端热湿气候区[1]。这样的气候条件将造成其传热、传湿的特性发生改变，导致建筑空调耗电问题异常突出。因此，传热和传湿过程必须被同时考虑，这样才能更准确地预测围护结构的热湿迁移以及室内舒适度。

用来描述围护结构热湿耦合传递的数学模型多种多样，Hou等[2]、Kong等[3]使用了水蒸汽压力为湿驱动势的热湿耦合传递模型，Fang等[4]、Tomáš等[5]、Min等[6]、kerget等[7]使用相对湿度为湿驱动势的热湿耦合传递模型，但这些模型的基础理论都相通，具体使用哪种模型需要根据实际需求来定。此外，刘向伟等[8]、张景欣等[9]、王莹莹[10]针对热湿耦合模型进行了多方面的研究，研究对象包括空调房间、计算方法、室内环境质量。然而，传统的热湿耦合模型针对极端热湿地区的应用较为稀少，尤其针对南海地区的研究更是近乎为零。笔者以温度和毛细压力为驱动势来预测围护结构的热量和水分传递，并根据极端热湿地区的高温高辐射多雨的气候条件对模型进行了优化。该模型充分考虑了极端热湿地区计算参数变物性的需求，使其能更加准确地描述极端热湿地区的特点。

1 热湿耦合模型

湿在围护结构中传递的过程呈现出两种不同相态：液相和气相，在传统计算中，液态水的传递一般被忽略。然而，由于极端热湿地区的气候特征，室外相对湿度常年在80%以上，液态水的传输在模型计算中不可忽略或简化。此外，蒸汽传输可以分为两部分：扩散和对流。蒸汽扩散由蒸汽压梯度驱动，其系数为蒸汽渗透系数，单位为kg/(m·Pa·s)；蒸汽对流由空气流动引起，由于极端热湿地区室外常年处于高温状态，室内很容易保持负压状态，空气压差的方向是由室内向室外，这将减轻热负荷。蒸汽对流过程计算较为复杂，因此，在热工设计计算上可以忽略[11]，这样一来，计算结果将偏安全。液态水的传输则经常用毛细压力作为驱动势来计算，其系数为液态水传导系数，单位为kg/(m·Pa·s)。

1.1 传湿模型

极端热湿地区围护结构传湿过程是一个包含水蒸气和液态水的传质过程，通过建筑围护结构中的蒸汽和液体迁移可以表达为[12]

(1)

(2)

式中：ρm为建筑材料的干密度， kg/m3；t为时间，s；u1为蒸汽含湿量，kg/kg；u2为液体含湿量，kg/kg；j1为蒸汽迁移量，kg/(m2·s)；j2为液体迁移量，kg/(m2·s)；I1为液体汽化引起的蒸汽源，kg/(m2·s)；I2为蒸汽凝结引起的液源，kg/(m2·s)。显然，I1=-I2。

叠加方程式(1)和式(2)，可以得到湿迁移的微分方程式

(3)

由Fick定律可得

(4)

式中：δv为材料的水蒸气渗透系数，kg/(m·Pa·s)；Pv为水蒸气分压力，Pa。此外，水蒸气分压力可以转化成相对湿度的关系式

Pv=φPsat

(5)

式中：φ为相对湿度；Psat为饱和水蒸气压力。

因此，式(4)可以转化为

(6)

由Darcy公式可知

(7)

式中：δl为液态水传导系数，kg/(m·Pa·s)；Pc为毛细吸附压，Pa。

由此可得湿传递方程

(8)

式中：w为材料含湿量，kg/m3。

1.2 传热模型

由能量守恒定律可知

(9)

式中：hi为物质的比焓(hi=ciT)，i=0、1、2，分别代表材料、蒸汽、液体的比焓，J/kg；角码′代表了干燥空气；jq为导热方式迁移的热流，jq=-λ*T，其中λ*为物体的有效导热系数，表示了通过毛细多孔体中水分和物体本身导热方式的热迁移。

由迁移方程(1)和(2)变换可得，

(10)

式中：c为湿物体的比热，J/(kg·K)；r12=h1-h2为蒸发潜热，J/kg；cp为湿空气比热，J/(kg·K)；ω为湿空气对流的速度；jk2为液体的对流迁移流，kg/(m2·s)；jm2为液体的分子扩散流，kg/(m2·s)。

在热湿迁移的一般情况下，方程(10)中方括号内的第3项比方程右边前两项小的多，可以忽略不计。因此有

(11)

式中：c为湿物体的比热，J/(kg·K)；r12为蒸发潜热，J/kg；λ*为物体的有效导热系数，表示了通过毛细多孔体中水分和物体本身导热方式的热迁移。

1.3 模型系数的求解

相对湿度可以转化为毛细压力的关系式，由Kelvin方程[13]可知

(12)

式中：ρw为水的密度，kg/m3；Rv为水蒸气气体常数，J/(kg·K)。

因此

(13)

饱和水蒸气分压力是空气温度的单值函数，已经有许多求饱和水蒸气分压力的公式，其中，联合国世界气象组织(WMO)推荐使用戈夫-格雷奇(Goff-Gratch)公式，中国现行《地面气象观测规范》也采用该公式求解饱和水蒸汽分压力，其公式为[14]

lg(101 3.246)

(14)

此公式较复杂，根据文献[15]研究，当空气温度在20～40 ℃范围内时，采用泰登(Tetens)公式的误差也较小。

(15)

式中：t为空气温度，℃。因此，采用泰登(Tetens)公式来计算饱和水蒸气压力。

2 极端热湿地区模型的特殊性

2.1 边界条件

任何数学模型的仿真模型都离不开边界条件的正确设定。极端热湿地区常年高温多雨，相比其他地区，太阳辐射和雨水在边界条件中的计算不可缺少。

在建筑围护结构内表面(x=0)

g=hmi(Wai(t)-W(0,t))

(16)

q=hci(Ti(t)-T(0,t))+hri(Ti(t)-

T(0,t))+r12hmi(Wai(t)-W(0,t))

(17)

在建筑围护结构外表面(x=L)

g=hme(Wae(t)-W(L,t))+grain

(18)

q=hce(Te(t)-T(L,t))+hre(Te(t)-

T(L,t))+r12hme(Wae(t)-W(L,t))+Qrad+qrain

(19)

式中：hmi、hme分别为室内外表面的质交换系数，kg/(m2·s)；hci、hce分别为室内外表面的对流换热系数，W/(m2·K)；hri、hre分别为室内外表面辐射换热系数，W/(m2·K)；Qrad为太阳辐射得热；grain为风吹动雨水淋在围护结构上造成的湿通量，kg/(m2s)；qrain为雨水带来的热通量，W/m2。

2.2 太阳辐射与雨水的计算

围护结构外表面的对流换热系数和辐射换热系数均与外界的气候因素(风速、天空温度、环境温度)以及外表面的温度有关，他们的值可以分别列式计算。但为了计算工作的简便，在热工设计计算中一般采用总换热系数hα来表示，用于统一计算对流和辐射总换热量。

hα=hc+hr

(20)

中国热工规范[16]规定，围护结构内表面的总换热系数在冬季是8.72 W/(m2·K)，夏季为8.75 W/(m2·K)，外表面的总换热系数冬季为23.3 W/(m2·K)，夏季为18.6 W/(m2·K)。该数值是按建筑物周围风速为3 m/s时实测所得。

根据围护结构外表面的热平衡方程式，可以得出由围护结构外表面向室内侧的传热量q0。

q0=hα(Tz-T0)

(21)

式中：Tz为室外综合温度，综合表达了室外空气温度、太阳辐射、大气长波辐射、地面反射和长波辐射对围护结构外表面的综合热作用。

采用室外综合温度计算太阳辐射对围护结构的影响。此外，根据文献[17]，雨水引起的热通量可以用式(22)表达。

qrain=graincl(Train-Ts)

(22)

式中：cl为水的比热容，约等于4 200 J/(kgK)；Train为雨水的温度，℃或K；Ts为围护结构表面温度，℃或K。

2.3 材料的物性参数

极端热湿地区建筑围护材料常年接触潮湿的空气，因此，其材料的物性参数将与平常状态下大不相同。围护结构材料的导热系数、蒸汽渗透系数和液态水渗透系数是随材料含湿量变化的参数，在计算中需要注意。

3 模型的数值计算与验证

采用COMSOL Multiphysics来求解模型。COMSOL Multiphysics会用有限元法自动对控制方程及相应的边界条件进行离散，并用数值求解器对离散后的代数方程组求解，可以大大简化求解的过程，方便研究人员将注意力集中到物理模型的研究。此外，由于热湿耦合传递是一个瞬态过程，模型的相关参数随着时间不断变化，为了反映这一点，用MATLAB编写了系数计算程序。

3.1 HAMSTAD基准案例

为了提高建筑围护结构的性能，改进建筑围护结构热湿耦合传递计算模型，欧盟在2000年发起了一项研究项目HAMSTAD(Heat, Air and Moisture Standards Development)[18]，其目的在于确定多孔材料的液体水转移特性和发展数值评估方法。该研究吸引了许多著名学校和重点实验室参与进来，其结果可靠，并且提供了详细的设置说明，方便模拟验证。

为了准确地反映极端热湿地区的特性，采用HAMSTAD-benchmark 4模型进行模拟验证。HAMSTAD-benchmark 4模型是一个比较特殊的模型，该模型分析了暴露在非稳态室外气候和相对稳定的室内气候条件下外墙内的温度和含湿量分布，在研究热湿耦合传递的基础上，进一步探索了雨水和太阳辐射对围护结构的影响，同时，整个模型是气密性的，蒸汽将以渗透的形式传输。模拟时长为120 h。此外，其物性参数随含湿量的变化如表1所示，满足研究极端热湿地区的需要。

表1 相关材料的物性参数[18]Table 1 Material properties of the materials

3.2 COMSOL软件设置的研究

为了保证模拟的准确性，先对COMSOL软件模拟的设置对计算结果影响进行了分析，首先进行容差设定分析，COMSOL提供了相对容差(Rtol)和绝对容差(Atol)对误差进行限定。因为绝对容差几乎总是必须被修改，因此，没法给出一个最好的通用缺省值，使用一个很小的绝对容差可能是一个很好的选择，这会使相对容差成为控制误差的决定因素，笔者在绝对容差设定时都采用比预期的变化小一个数量级。相对容差软件默认的是0.01，将相对容差分别设定成10-9、10-3和10-2，对结果进行最大误差分析。

相对容差10-3和10-9之间含湿量的最大相对误差为0.7%，10-2和10-9之间含湿量的最大相对误差为12.8%。从结果上看，相对容差设定成10-9或10-3较为合适，按默认取值会造成传湿计算的较大误差，然而，相对容差设定成10-9的计算时间是其他两个值的10倍左右，为了方便，相对容差推荐设定成比默认值小10倍，即0.001。

此外，对网格划分对结果的影响进行对比分析，网格按极细化、较细化和常规进行划分，最大相对误差都在0.1%以内。结果上看，COMSOL软件默认的网格划分都较细致，但由于计算时间差距不大，推荐根据计算的目标来采用极细化或较细化的网格来减少网格划分带来的误差。

此外，值得注意的是，当材料处于湿饱和状态，毛细压力超过饱和毛细压力时，雨水将不再进入围护结构而是直接落下。因此，在计算的时候，需要将其分别计算。笔者通过在全局变量添加插值函数satu(w)来模拟这一过程。

因此，淋雨的湿通量为grain×satu(w)。

3.3 结果对比

图1～图4中TEST曲线为该模型的计算结果，其他曲线是由HAMSTAD参与学校或组织提供，分别为：Chalmers University of Technology, Coordinator of WP’(CTH)、TNO Building and Construction Research (TNO)、University of Leuven (KUL)、Technical University of Dresden (TUD)、Technion Israel Institute of Technology (Technion)、Eindhoven University of Technology (TUE)、National Research Council of Canada (NRC)、Fraunhofer Institute of Building Physics (IBP)。

图1 外表面温度Fig.1 Temperature related to time at the external surface

图2 内表面温度Fig.2 Temperature related to time at the internal surface

图3 外表面含湿量Fig.3 Moisture content related to time at the external surface

图4 内表面含湿量Fig.4 Moisture content related to time at the internal surface

从结果来看，无论内表面还是外表面，该模型所代表的曲线TEST均可以准确地反映出在雨水和太阳辐射影响下多层围护结构的热湿耦合传递过程，模拟结果准确地反映了太阳辐射影响下围护结构表面的波动，以及在雨水作用下围护结构的传湿。模拟结果与基准值基本一致，该模型模拟结果良好。

4 结论

以毛细压力与温度为驱动势，建立了一种适应于极端热湿地区围护结构热湿耦合传递模型，在该模型的基础上，提出了基于计算机软件COMSOL和MATLAB的热湿耦合模型的求解方法，并与HAMSTAD基准进行了验证，验证结果良好。得到以下结论：

1)COMSOL模拟热湿耦合传递时，相对容差推荐设定为0.001，绝对容差设定为0.000 1，网格按极细化划分。

2)在极端热湿地区，为了准确反映高温和多雨对围护结构的影响，提出了适宜于极端热湿地区的热湿耦合模型边界条件，该方法可以准确地模拟出太阳辐射和雨水对围护结构传热传湿的影响。

3)在计算雨水负荷时，需要注意围护结构的吸附能力，当围护结构表面达到饱和时，多余的雨水应该不参与计算。

参考文献：

[1] 刘毅. 提升南海气象观测能力刻不容缓——专访中国工程院院士丁一汇[N]. 人民日报, 2015-06-21(2).

LIU Y. Enhance the ability of the South China Sea meteorological observations without delay [N]. People's Daily, 2015-06-21(2). (in Chinese)

[2] HOU X D, LIU F S, WANG S J. Coupled heat and moisture transfer in hollow concrete block wall filled with compressed straw bricks [J]. Energy and Buildings, 2017, 135(15): 74-84.

[3] KONG F H, ZHANG Q L. Effect of heat and mass coupled transfer combined with freezing process on building exterior envelope [J]. Energy and Buildings, 2013, 62: 486-495.

[4] FANG X D, ANDREAS K A, PAUL P F. Methodologies for shortening test period of coupled heat-moisture transfer in building envelopes [J]. Applied Thermal Engineering, 2009, 29(4): 787-792.

[6] MIN H G, ZHANG W P, ROBERT C. Coupled heat and moisture transport in damaged concrete under an atmospheric environment [J]. Construction and Building Materials, 2017, 143: 607-620.

[8] 刘向伟,陈友明,陈国杰, 等. 围护结构热湿耦合传递模型及简便求解方法[J]. 土木建筑与环境工程,2016, 38(4): 7-12.

LIU X W, CHEN Y M, CHEN G J, et al.Coupled heat and moisture transfer model and simple solution method for building envelopes [J]. Journal of Civil, Architectural & Environmental Engineering, 2016, 38(4): 7-12. (in Chinese)

[9] 张景欣,郭兴国,陈友明, 等. 墙体内热湿耦合过程的时域递归展开算法[J]. 土木建筑与环境工程,2015, 37(6): 147-152.

ZHANG J X, GUO X G, CHEN Y M, et al.A time domain recursive algorithm for solving the model of coupled heat and moisture transfer in building wall [J]. Journal of Civil, Architectural & Environmental Engineering, 2015, 37(6): 147-152. (in Chinese)

[10] 王莹莹. 围护结构湿迁移对室内热环境及空调负荷影响关系研究[D].西安:西安建筑科技大学, 2013.

WANG Y Y. Research on the effect of the palisade structure moisture transfer on the indoor thermal environment and air-conditioning load [D]. Xi'an: Xi'an University of Architecture and Technology, 2013. (in Chinese)

[11] 彦启森. 建筑热过程[M].北京:中国建筑工业出版社, 1986.

YAN Q S.Building heat process [M]. Beijing: China Architecture & Building Press, 1986. (in Chinese)

[12] A.B.雷柯夫. 建筑热物理理论基础[M]. 北京: 科学出版社, 1965.

LUIKOV A B. Theoretical basis of building thermophysics [M]. Beijing: Science Press, 1965. (in Chinese)

[13] LI Q R, JIWU R, FAZIO P. Development of HAM tool for building envelope analysis [J]. Building and Environment, 2009, 44(5): 1065-1073.

[14] Smithsonian Meteorological Tables [Z]. Washington: Smithsonian Institution,1984.

[15] 周西华,梁茵,王小毛. 饱和水蒸汽分压力经验公式的比较[J]. 辽宁工程技术大学学报,2007(3):331-333.

ZHOU X H, LIANG Y, WANG X M. Comparison of saturation vapor pressure formulas [J]. Journal of Liaoning Technical University, 2007(3): 331-333.(in Cinese)

[16] 中华人民共和国建设部. 民用建筑热工设计规范: GB 50176—93 [S]. 1993.

Ministry of Construction of the People's Republic of China. Thermal design code for civil building: GB 50176-93 [S]. 1993.(in Chinese)

[17] KÜNZEL H M. Simultaneous heat and moisture transport in building components: one and two-dimensional calculation using simple parameters [D]. Institute of Building Physics, Germany,1995.

[18] HAGENTOFT C E. HAMSTAD-Final report: methodology of HAM-modeling，Report R-02:8 [R]. Gothenburg, Department of Building Physics, Chalmers University of Technology, 2002: 8-50.