我国是世界上人造板及其制品的生产大国和消费大国，其中地采暖地板是我国出口创汇的重要组成部分之一。据中国林产工业协会地板专业委员会发布的《2017年我国地板行业销量统计报告》数据显示，2017年我国具有一定规模的地板企业木竹地板总销量约41 520万m²，同比增长约4.64%。随着地采暖地板销量的显著上升，其性能指标检测已成为一个亟需解决的问题。然而，目前国内外尚无专门针对地采暖地板蓄能、缓释和节能等性能的研究成果，也没有专门针对地采暖地板性能的检测检验仪器和评价体系，已有研究成果多为针对木材导热系数的理论推导(陈瑞英等，2005a; 2005b；林铭等，2013；杨文斌等，2006；赵景尧等，2015；Gu et al., 2005；Dupleix et al., 2013)和地采暖地板热辐射系统的建模分析(何毅等，2006；张云斌等，2005；张群力等，2006；冉茂宇，2010；Mi et al., 2015)等。为了向相关质量检测检验部门提供地采暖地板蓄能、缓释、节能等性能的检测方法与仪器，本研究提出一种基于密闭绝热检测室的地采暖地板蓄热性能检测方法，并建立相应传热模型，从而为蓄热性能检测仪器的开发提供理论依据，也将为地采暖地板蓄热性能的检测提供一种有效的研究方法和检测手段，这对我国地采暖地板的产品标准化与市场规范化具有重要意义。

1 检测装置介绍

检测室为圆柱腔体，外壁均作抽真空和防热辐射处理，与外界完全隔热；检测室内由上至下均匀安放温度传感器阵列，共N层，每层按照等距同心圆放置共S个温度传感器，用于测量检测室内各空间节点的温度变化。当检测室达到初始温度后，将加热到一定温度、固定体积的地采暖地板样本放置到检测室底部中心托举网上，然后密闭。检测室内壁和空气温度会随着样本散热而上升，样本温度会随着其自身散热而下降，温度传感器阵列实时记录检测室层间各单元节点温度，传感器阵列通过CAN总线方式连接到系统总线控制器，总线控制器再通过Zigbee无线通信模块将数据传输到系统主控制器上，并存储到存储器中。当检测室温度和样本温度达到平衡时，系统进入稳定状态。此时，可根据检测室内壁和空气的平衡温度，计算地采暖地板样本的热释放量；同时，可根据温度传感器阵列的测量值，动态再现地采暖地板材种及热物特性(材料密度、比热、导热系数等)对检测室内温度场的影响，从而揭示地采暖地板的蓄热性能规律，为地采暖地板的推广应用提供质量检测检验方法与仪器。

2 研究方法 2.1 检测室模型的建立

首先，构建密闭绝热圆柱腔检测室的三维物理模型(张敏等，2010；2011；王治云，2009；Man et al., 2000；Dong et al., 2004)，如图 1所示，圆柱腔内外壁和上下盖板均由双层不锈钢板制成，夹层之间为真空，真空度为10^-4Pa。上下盖板与圆柱腔内壁及上下壁之间采用高绝热橡胶材料密封，使之与外界最大程度绝热。

通常情况下，空气在流动过程中密度变化很小，可视为不可压缩流体，检测室内空气的所有热物性均视为常数，但密度随温度的变化遵循Boussinesq假设，采用这种方法相对于设定密度为温度的函数方法来求解问题要收敛得快。检测室内壁及空气初始温度为T₀。检测室底部中心位置托举网放置实木地板样本，样本为弦切面板材，长度为a，宽度为b，高度为h，假设样本内部温度均匀，初始温度为T_h(T_h>T₀)。

2.2 地采暖地板蓄热规律方程的建立

在图 1所示的圆柱腔检测室内，设u、v、w分别为检测室内x、y、z3个方向的热流速度，t为时间，μ为空气运动学黏度，P为空气压力，T₁为空气温度，ρ₁为空气密度，g为重力加速度，方向垂直向下。根据流体力学定律，检测室内空气流动与换热的质量守恒方程、动量方程和能量守恒方程分别为：

$ \frac{{\partial u}}{{\partial x}} + \frac{{\partial v}}{{\partial y}} + \frac{{\partial w}}{{\partial z}} = 0; $

(1)

$ \begin{array}{l} \frac{{\partial \left({{\rho _1}u} \right)}}{{\partial t}} + u\frac{{\partial \left({{\rho _1}u} \right)}}{{\partial x}} + v\frac{{\partial \left({{\rho _1}u} \right)}}{{\partial y}} + w\frac{{\partial \left({{\rho _1}u} \right)}}{{\partial z}} = \\ \frac{\partial }{{\partial x}}\left({\mu \frac{{\partial u}}{{\partial x}}} \right) + \frac{\partial }{{\partial y}}\left({\mu \frac{{\partial u}}{{\partial y}}} \right) + \frac{\partial }{{\partial z}}\left({\mu \frac{{\partial u}}{{\partial z}}} \right) - \frac{{\partial P}}{{\partial x}}; \end{array} $

(2)

$ \begin{array}{l} \frac{{\partial \left({{\rho _1}v} \right)}}{{\partial t}} + u\frac{{\partial \left({{\rho _1}v} \right)}}{{\partial x}} + v\frac{{\partial \left({{\rho _1}v} \right)}}{{\partial y}} + w\frac{{\partial \left({{\rho _1}v} \right)}}{{\partial z}} = \\ \frac{\partial }{{\partial x}}\left({\mu \frac{{\partial v}}{{\partial x}}} \right) + \frac{\partial }{{\partial y}}\left({\mu \frac{{\partial v}}{{\partial y}}} \right) + \frac{\partial }{{\partial z}}\left({\mu \frac{{\partial v}}{{\partial z}}} \right) - \frac{{\partial P}}{{\partial y}}; \end{array} $

(3)

$ \begin{array}{*{20}{l}} {\frac{{\partial \left( {{\rho _1}w} \right)}}{{\partial t}} + u\frac{{\partial \left( {{\rho _1}w} \right)}}{{\partial x}} + v\frac{{\partial \left( {{\rho _1}w} \right)}}{{\partial y}} + w\frac{{\partial \left( {\rho w} \right)}}{{\partial z}} = }\\ {\frac{\partial }{{\partial x}}\left( {\mu \frac{{\partial w}}{{\partial x}}} \right) + \frac{\partial }{{\partial y}}\left( {\mu \frac{{\partial w}}{{\partial y}}} \right) + \frac{\partial }{{\partial z}}\left( {\mu \frac{{\partial w}}{{\partial z}}} \right) - \frac{{\partial P}}{{\partial z}} + {\rho _1}{g};} \end{array} $

(4)

$ \begin{array}{l} \frac{{\partial {\rho _1}{c_1}{T_1}}}{{\partial t}} + \frac{{\partial \left({{\rho _1}{c_1}u{T_1}} \right)}}{{\partial x}} + \frac{{\partial \left({{\rho _1}{c_1}v{T_1}} \right)}}{{\partial y}} + \frac{{\partial \left({{\rho _1}{c_1}w{T_1}} \right)}}{{\partial z}} = \\ \frac{\partial }{{\partial x}}\left({\lambda \frac{{\partial {T_1}}}{{\partial x}}} \right) + \frac{\partial }{{\partial y}}\left({\lambda \frac{{\partial {T_1}}}{{\partial y}}} \right) + \frac{\partial }{{\partial z}}\left({\lambda \frac{{\partial {T_1}}}{{\partial z}}} \right)。\end{array} $

(5)

设T₂、T₃分别为样本温度和检测室内壁温度，λ₂、λ₃为样本导热系数和检测室内壁导热系数，ρ₂、ρ₃为样本密度和检测室内壁密度，C₂、C₃为样本比热容和检测室内壁比热容。根据能量守恒定律，样本和检测室内壁的能量方程分别为：

$ \frac{{\partial {T_2}}}{{\partial t}} = \frac{{{\lambda _2}}}{{{\rho _2}{C_2}}}\left({\frac{{{\partial ^2}{T_2}}}{{\partial {x^2}}} + \frac{{{\partial ^2}{T_2}}}{{\partial {y^2}}} + \frac{{{\partial ^2}{T_2}}}{{\partial {z^2}}}} \right); $

(6)

$ \frac{{\partial {T_3}}}{{\partial t}} = \frac{{{\lambda _3}}}{{{\rho _3}{C_3}}}\left({\frac{{{\partial ^2}{T_3}}}{{\partial {x^2}}} + \frac{{{\partial ^2}{T_3}}}{{\partial {y^2}}} + \frac{{{\partial ^2}{T_3}}}{{\partial {z^2}}}} \right)。$

(7)

2.3 确定蓄热规律方程的边界条件及初值

首先，设检测室所有内壁均与外界完全绝热，内壁自身的温度会随着空气温度上升而上升。设，边界条件为：

空气温度：t=0时刻，T₁=T₀；

内部热源(样本)：t=0时刻，当－a/2 < x < a/2，且－b/2 < y < b/2，且0 < z < h时，T₂=T_h；

上壁面：z=H，t=0时刻，T₃=T₀；

侧壁面：0 < z < H，t=0时刻，在任意x－y平面内，当$\sqrt {{x^2} + {y^2}} = R$时，T₃=T₀；

下壁面：z=0，t=0时刻，当a/2 < x < R，或－R < x < －a/2，或b/2 < y < R，或－R < y < －b/2时，T₃=T₀。

对检测室进行网格划分，网格节点间隔长度为2 mm。其中，样本、空气采用结构化网格体，检测室内壁采用非结构化网格体。检测室半径为100 mm，高度为200 mm，内壁厚为2 mm，材料为不锈钢，密度为8 030 kg·m^-3，比热容为502 J·kg^-1K^-1，导热系数为16.27 W·m^-1K^-1。样本长为100 mm，宽为60 mm，高为15 mm。空气和内壁初始温度为20 ℃，样本初始温度为70 ℃。

检测样本选用白桦(Betula platyphylla)、水曲柳(Fraxinus mandshurica)、西南桦(Betula alnoides)和柞木(Quercus mongolica) 4种材种做成的地采暖地板，其密度、比热容和传热系数等热物特征参数如表 1所示(林铭等，2015)。

2.4 蓄热规律方程数值求解

根据流体力学原理，反映热传递性能的瑞利数需要准确确定：

$ {R_{\rm{a}}} = \frac{{g\beta \mathit{\Delta }{\mathit{T}_1}{h^3}{\rho _1}}}{{v\alpha }}。$

(8)

式中：α、β分别表示热扩散系数和膨胀系数；g为重力加速度；ΔT₁为空气温差；h为检测室高度；ρ₁为空气密度；υ为空气的运动黏度。$\alpha = \frac{{{\lambda _1}}}{{{\rho _1}{c_{\rm{p}}}}}$、$\beta = - \frac{1}{{{\rho _1}}}\frac{{\partial {\rho _1}}}{{\partial {T_1}}}$，λ₁为空气导热系数，c_p为空气比热容。实际计算过程中，g取9.8 m·s^-1，ρ₁取1.205 kg·m^-3，λ₁取0.0 259 W·m^-1K^-1，υ取1 005 J·kg^-1K^-1。

根据设定的边界条件、初值及计算参数，将得到的地采暖地板蓄热规律方程式通过流体力学计算软件FLUENT进行迭代求解，求解流程如图 2所示。最终求得地采暖地板蓄热性能分析模型，其中，能量方程的收敛标准为1×10⁶，连续性方程和动量方程的收敛标准取0.001。

3 蓄热模型结果分析 3.1 第1 000步时温度分布

图 3为第1 000步(步长为10 s)时检测室内X=0处Y-Z截面的温度分布。由图 3可知，紧邻样本周围区域温度上升最快，由于空气浮力作用，样本正上方的狭窄区域相对于其他区域温度上升得更为显著，越接近于样本表面，温度越高，而狭窄区域两侧空间的温度呈对称分布，由于空气重力作用，越接近检测室顶部其温度越高。其中，样本为白桦时，检测室内空气最高温度为25.5 ℃，检测室顶部温度为24.6 ℃，样本最高温度下降到25.7 ℃(图 3a)。样本为水曲柳时，检测室内空气最高温度达到26.8 ℃以上，检测室顶部温度为25.6 ℃，样本最高温度下降到27.2 ℃(图 3b)。样本为西南桦时，检测室内空气最高温度达到26.8 ℃以上，检测室顶部温度为25.6 ℃，且区域明显大于水曲柳，样本最高温度下降到27.2 ℃(图 3c)。样本为柞木时，检测室内空气最高温度达到27.4 ℃，检测室顶部温度为25.8 ℃，温度达到26 ℃的区域明显大于西南桦，样本最高温度下降到27.6 ℃(图 3d)。

根据图 3所示，检测室内在X=0处Y-Z截面等温线变化趋势可得，样本热物理特性对检测室内的温度分布和采暖地板的蓄热性能具有明显影响(图 4)。即当样本分别为水曲柳和西南桦时，检测室内空气最高温度、检测室顶部温度和样本温度完全相同，只是分布区域略有差别。这说明地采暖地板材种的密度和比热相对于导热系数对检测室内的温度影响更大。当样本为柞木时，空气最高温度、检测室顶部温度和样本温度相对于西南桦分别增大0.6和0.2 ℃。这是由于虽然柞木的密度和导热系数同西南桦几乎相同，但比热容比西南桦大88 J·kg^-1K^-1，因此，密度和导热系数相同时，样本的比热和检测室内温度值的大小呈正比关系。当样本为白桦时，空气最高温度、检测室顶部温度和样本温度均为4种木材中的最小值，这是由于白桦的密度、比热和导热系数均最小，样本的密度和检测室内温度值的大小呈正比关系。

3.2 平衡时刻

随着样本散热继续进行，当检测室内空气和样本的温度最大值和最小值之差小于等于0.1 ℃时，定义当前时刻为系统平衡时刻，当前温度分布状态为系统平衡状态。计算结果表明，白桦约在2 300步时(步长为10 s)，检测室内温度达到平衡状态，室内平均温度为24.4 ℃；水曲柳约在2 700步时(步长为10 s)，检测室内温度达到平衡状态，室内平均温度为25.32 ℃；西南桦和柞木约在2 600步时(步长为10 s)，检测室内温度达到平衡状态，室内平均温度分别为25.35和25.4 ℃(图 5)。总结得知，当样本为水曲柳和西南桦时，室内平衡温度很接近，但西南桦到达平衡状态的时刻比水曲柳少100步，可见，蓄热性能由样本的密度、比热决定，随着样本导热系数增加，样本自身的散热速度也在加快。当样本为白桦时，平衡状态温度最低，这是由于白桦的密度和比热在4种样本中均为最小值。

4 实测结果对比分析

本研究使用自主研制的地采暖地板蓄能规律分析仪分别对白桦、水曲柳、西南桦和柞木4种材种地采暖地板样本进行蓄热性能测试，每种样本随机选取3块试件，测量结果的平均值作为最终试验数据。分析仪检测室在X=0处Y-Z截面内均匀分布10×10个温度传感器，可实时监测检测室内各空间单元的温度变化，传感器测量精度为0.1 ℃。在不同样本所对应的温度计算值与实测值之间做线性回归，计算平均相对误差及相关系数(表 2)(周玉成等，2007)。

图 6为不同样本对应的各传感器节点的温度计算值与实测值曲线。由表 2和图 6可知，对4种不同样本进行蓄热性能测试，其计算结果的平均相对误差均可控制在2%以内，实测值均高于计算值，这是由不锈钢内壁的辐射换热引起的；相关系数则达到了0.88以上，所建立的模型可以准确分析地采暖地板的蓄热规律。继续测试可得，白桦约在2 200步时(步长为10 s)，检测室内温度达到平衡状态，室内平均温度为24.2 ℃；水曲柳约在2 500步时(步长为10 s)，检测室内温度达到平衡状态，室内平均温度为25.1 ℃；西南桦和柞木约在2 400步时(步长为10 s)，检测室内温度达到平衡状态，室内平均温度分别为25.2 ℃和25.3 ℃。平衡时间和平衡温度均小于模拟值，这是由于不锈钢内壁的辐射加快了传热的速度，而实际测试过程中不可避免地存在热量的损耗，从而降低了最终的平衡温度。

5 结论

针对地采暖地板蓄热性能检测需求，构建了密闭绝热圆柱腔结构的检测室物理模型，在检测室底部中心位置分别放置材种为白桦、水曲柳、西南桦和柞木的地采暖地板样本，建立了地采暖地板蓄热规律方程，对检测室内温度分布进行有限元分析，并与所研制的地采暖地板蓄热性能分析仪的实测结果进行对比，以平均相对误差和相关系数为评价依据，对模型计算结果进行分析验证。结果表明，所构建模型和方法能够准确反映采暖木地板材种及热物特性对检测室内温度场的影响，实现地采暖地板蓄热性能的检测，从而揭示了地采暖地板的蓄热规律。该模型的建立可为地采暖地板蓄势性能的检测提供理论基础。