北京大兴国际机场航站楼核心区钢网格结构日照非均匀温度场研究

周 勐，樊健生，刘宇飞，张晋勋，段先军，雷素素

(1.清华大学土木工程系，北京 100084；2.北京城建集团有限责任公司，北京 100088)

大跨空间结构具有结构体量与平面尺寸大等特点，同时在漫长的施工阶段往往要承受复杂的环境作用，经历多样的边界条件变化。由于大跨空间结构在施工过程中一般处于露天环境，太阳辐射与日温度变化会在结构中产生瞬态不均匀的温度场，进而引起无法忽略的温度效应。近年来，随着大跨空间结构数量的增加与规模的扩大，结构温度场的研究也变得越来越重要。

有关太阳辐射的理论研究是温度场模拟的重要基础。国内外学者提出了多个适用于不同区域、气候特点的晴朗天空太阳辐射模型[1－4]，并依据实测数据对模型参数进行了标定，可用于后续的温度场计算。围绕大跨空间结构温度场的研究主要可以分为构件温度分析、结构体系温度场模拟、结构温度场实测研究这几个主要研究方向。构件温度方面，已有研究对象包括圆截面构件[5]、H型截面构件[6]及箱型截面构件[7]，结果表明日照下构件的温度场不均匀，辐射强度、构件倾角、辐射吸收系数等参数对温度场的影响较大，但由于影响的因素众多，关系复杂，对于构件的研究结果往往难以在实际结构中直接应用。对体系温度场的研究大多考虑了由太阳主导的瞬态温度场，结果显示结构中的日温度变化非常明显，不可忽略，但除天线结构[8]外的其他结构均未考虑构件间的遮挡[9－11]。结构温度场及温度效应的实测研究[12－13]表明，太阳辐射引起的温度场在空间与时间上的非均匀性都十分显著，但规律性并不明显；同时日温度变化会使构件应力与变形发生周期性的波动。

随着大跨空间结构施工过程的不断推进，不均匀温度场引起的内力、位移变化很可能会由于结构的拼装嵌补而无法释放，不断累积，且最终的积累量与施工方式、各施工步的时间节点均有关系。可见，对于规模大、形状复杂、重要性高的大跨空间钢结构，准确模拟其瞬态温度场，确定结构任意位置任意时刻的温度值是十分必要的，也是进行精确施工模拟的重要前提。由于大跨空间结构的空间造型一般比较复杂，在太阳辐射下构件间会相互遮挡产生阴影，使结构温度场变得更加难以预测。目前，通用有限元软件没有内置的太阳辐射模块，无法考虑太阳辐射的变化，更无法考虑日照下各单元之间的复杂遮挡关系。

本文编制了构件太阳辐射量计算程序，实现了适用于大跨空间钢结构的日照温度场模拟方法，不仅可以充分考虑地理位置、气温季节变化、太阳辐射等环境因素的影响，同时，还可以完成大量圆截面杆系构件的遮挡效应计算，从而，能够高效准确地模拟结构在日照下的温度场，得到结构中任意位置任意时刻的温度值。以北京大兴国际机场为背景，对日照温度场模拟方法的重要求解条件进行了测试与说明，并对空间曲面钢网格结构在考虑日照作用下的温度场特点进行了研究，为大跨空间结构施工过程的精确模拟提供了重要基础。

1 考虑地理气象信息的温度场模型

随着数值计算的发展与计算机性能的提高，数值模拟逐渐成为了温度场研究的主要手段。为得到结构在任意时刻的温度分布，需要根据诸多热边界模型计算其瞬态温度场。

热传递的主要形式有三种，分别是热传导、热辐射与热对流。具体而言，结构瞬态温度场的影响因素包括结构各相连构件之间的热量传递，构件吸收的太阳辐射、日气温变化、构件与空气的对流换热、构件与环境间的热辐射交换几个部分。

1.1 晴朗天空的太阳辐射模型

一个完整的太阳辐射模型包括两部分：一部分是计算任意日期、时刻下的太阳位置；另一部分是计算大气层对太阳辐射的折减，得到任意平面的实际辐射量。

1.1.1 太阳位置的天文地理计算

地球的公转与自转、观测点的经纬度均会影响太阳的位置。其中，地球公转影响太阳赤纬角δ，自转影响太阳时角ω，二者的计算方法如下[1]：

式中：N为一年中的日序数，取值范围为1~365；hours为观测点的太阳时，即经过误差修正的当地时间，取值范围0 h~24 h。

太阳光线与地面任意表面的相对位置关系可用太阳高度角、太阳方位角与表面入射角表示，如图1所示。

二者的计算方法如下：

式中，φ为观测点所在的纬度，取值范围为-90°~90°，北半球为正。

对于任意倾斜表面，太阳光线入射角θ的计算方法如下：

式中：β为表面的倾角取值范围为 0°~180°，超过90°表示表面向下；γ为表面的方向，取值范围与γs相同。

1.1.2 任意倾斜表面的太阳辐射强度

太阳发出的辐射到达地球大气层后，大部分能够直接到达地面，称为直接辐射；小部分能量经过大气中各种分子的散射与吸收后到达地面，称为散射辐射。同时，到达地面的部分太阳辐射会被地面反射。地表任意结构表面受到的太阳辐射即由这三者组成。

图1 太阳光线与任意倾斜表面的相对位置Fig.1 Relative position of the sunlight to any sloping surface

考虑公转引起的日地距离变化，大气层外的太阳辐射强度为：

其中，Gsc为太阳常数，取1367 W/m2，具有1%的不确定度[1]。

考虑大气层的复杂作用及太阳高度角影响后，水平面的直接辐射强度Gb与散射辐射强度Gd可以按下式进行计算[1]：

式中：TL为Linke浊度系数[14]，表示气溶胶与水蒸气对散射的影响，取值范围为 2~8，大城市一般可取4；m为经过大气压修正的大气光学质量，代表大气分子对散射的影响，计算方法如下[3]：

式中：z为观测点的海拔高度；zh为Rayleigh大气高度，取8434.5 m。

对于地表的任意倾斜表面，直接辐射的实际强度需要通过入射角进行折减；散射辐射与地面反射辐射的实际强度需要通过辐射角系数确定。最终得到的晴朗天空下倾斜表面的太阳辐射强度为：

式中：ρr为地面反射率，一般取0.2；当cosθ为负时，太阳光线无法照射到表面的正面，不计入直接辐射的贡献。

1.2 构件日照阴影及辐射量计算

大跨空间钢结构一般具有较为复杂的造型，构件数量多，且一般上下分层。受到太阳光的照射时，构件之间存在不断变化的遮挡关系，形成大量日照阴影，影响太阳辐射的接收量，进而影响结构中的温度分布。已有的遮挡算法研究主要集中在精度要求很高的天线领域[11]，相关研究表明，是否考虑遮挡对构件温度有较为明显的改变，不可忽略。

计算机领域常见的消隐算法用于消除不可见的对象，从而正确显示对象间的位置关系。本文参考图像空间消隐算法的思路，实现了基于包围盒算法的日照阴影与辐射量计算方法，据此分析圆截面杆系构件受到遮挡后的影响。主要计算步骤如下：

1)坐标系变换。将模型从表示实际位置的原坐标系OXYZ变换到光线坐标系oxyz中，将三维问题解耦，降为二维；

2)单元投影的包围盒检测。根据各单元在光线坐标系xoy平面上的投影，进行初步的位置判断，提高计算效率；

3)单元遮挡关系判断。对每组单元对，先判断平面投影的位置关系，再根据交点的高度判断遮挡关系；

4)单元遮挡区间计算。根据单元对的位置关系与单元截面尺寸，计算被遮挡单元表面的阴影长度及位置；

5)单元日照系数计算。对各单元的表面阴影取并集，并用其未受遮挡部分长度占总长度的比例rs来衡量日照阴影对其太阳直射分量的减弱效果。

记构件表面对太阳辐射的吸收率为α，考虑遮挡关系后，构件单位时间吸收的太阳辐射总量为：

式中：Ib、Id和Ir分别为构件单位时间受到的太阳直射、散射和反射辐射量分项；太阳辐射吸收率受构件表面的材料与颜色影响。

1.3 日气温模型

日气温变化的影响因素复杂，难以进行精确的预测模拟。但总体上具有一定的变化规律，晴朗天气下，气温与太阳辐射的关联度很高：每天日出后气温开始升高，一般在 14:00~15:00左右达到当天的最高值；降温过程中，日落前的降温速度要大于夜间；日出前温度达到最低。

已有的相关研究多采用正弦曲线近似模拟日气温随时间的变化。本文结合日气温的实际变化规律[15]，改用正弦曲线(白天)与直线段(夜间)结合的模型，对日气温的变化进行模拟，如图2所示。

图2 日气温变化模型Fig.2 Daily temperature model

模型中需要给出一天的最高气温Tmax与最低气温Tmin，日出与日落时间可通过前文的太阳辐射模型计算得到。

1.4 热对流模型

由于空气在不断流动，物体表面每时每刻都在与周围的空气进行对流热交换，基本形式遵循牛顿冷却定律：

式中，为物体表面的对流换热系数，一般与物体温度、空气温度、换热表面的尺寸、风速等有关，对不同物体需要采用不同的计算方法。

以大跨空间钢结构中经常采用的圆截面钢管为例，其对流换热属于外掠单管换热问题[16]，表面对流换热系数的计算公式如下：

式中：d为圆管外径；k为空气导热系数；u为空气流速；v为运动粘度系数；系数c与n根据雷诺数分段取值，如表1所示。

表1 系数c和n分段取值Table 1 Segmentation values of coefficient c and n

1.5 长波辐射交换

物体与环境之间不断进行热辐射交换，也称为长波辐射交换。

物体周围的环境可假定为温度与气温相同的黑体，且其与物体的辐射角系数为 1。根据黑体的热辐射定律，可以得到物体与环境间长波辐射交换的热流密度：

2 研究对象及温度场模拟方法

2.1 北京大兴国际机场航站楼屋盖结构特点

北京大兴国际机场位于永定河北岸，北京大兴区礼贤镇、榆垡镇和河北廊坊市广阳区之间，定位为大型国际航空枢纽。

航站楼由主楼和五条指廊组成，建筑面积约70 万平方米，屋面最高点50 m，地上五层，地下二层。航站楼的主体结构为全现浇钢筋混凝土框架结构；屋盖为不规则自由曲面钢网格结构；支撑体系由C型钢柱、支撑筒、钢管柱及幕墙柱组成。结构典型示意图如图3所示。

本文主要对航站楼核心区钢屋盖展开研究，其主要结构特点有：

1)造型复杂，空间位形控制精度要求高。屋盖为不规则自由曲面球节点交叉桁架结构，空间坐标及标高落差大，传力路径复杂，如图3(b)。同时构件密度大、管径粗，相互之间的遮挡效果非常明显；

2)结构面积大，跨度大，支承体系少。钢屋盖覆盖范围横向宽504 m，纵向长462 m，总投影面积达18万平方米。各主要支承构件距离60 m~80 m不等，且屋盖外沿最大悬挑达43 m，结构跨度大，如图3(c)；

3)施工过程复杂，温度效应累积明显。屋盖总体的施工原则为分区施工，分区卸载，总体合拢，结构分为7个分区，每个分区内部进一步分块，如图3(d)。在累计多次的提升、拼装、合拢、卸载等施工过程中，温度荷载会产生显著的影响。

图3 北京大兴国际机场Fig.3 Beijing Daxing international airport

综合以上特点，不均匀的温度变化会在结构中产生明显的温度效应，并会随着施工过程不断积累，因此，对结构中任意时刻温度场的精确模拟就显得尤为重要。考虑到整体结构体量太大，本文选取屋盖结构中的典型分区C2-1进行后续研究。

2.3 温度场模拟方法

本文以通用有限元软件ANSYS的热分析模块为平台，结合MATLAB及APDL语言，计算C2-1区屋盖钢结构考虑日照遮挡的瞬态温度场。整体计算流程如图4所示，图中的每个步骤的计算结果均与时间相关。

前面说过，“本心”是一种能动的精神，因此“发明本心”也可以解释成“挖掘或激发内心的能动的精神”。如果这么解释，那我们就可以发现增能理论和陆九渊心学各自的第二个核心观念“增强权能”和“发明本心”同样具有相似性。

主要计算流程可以分为太阳辐射及遮挡模块与热分析模块两大部分，太阳辐射模块的计算结果需要导入温度场模型中，作为热边界条件之一。太阳辐射模块的主要内容与计算过程在前文已经说明，本文用MATLAB实现；热分析模块在ANSYS平台上结合APDL语言实现，相关细节在下文进行说明。

图4 温度场计算流程图Fig.4 Flow chart of temperature field calculation

3 结构日照温度场求解条件取值

3.1 有限元模型

根据北京大兴国际机场 C2-1施工分区的结构组成与特征，建立热分析有限元模型，如图5所示。

图5 结构热分析有限元模型Fig.5 Finite element model for thermal analysis

该分区结构构件均为钢桁架，因此选用三维二节点单元进行模拟。模型采用Link33单元模拟结构杆件，赋予截面面积及导热系数，实现结构内部的热传导；采用Link34单元将所有结构节点与代表空气的节点相连，设定对流换热系数，模拟结构与空气的热对流；采用Link31单元将所有结构节点与代表环境的节点相连，设定表面辐射率，模拟结构与环境间的长波辐射。此外，构件吸收的太阳辐射通过给单元施加等效生热率实现。需要注意的是，结构杆件具有实际尺寸，其表面温度分布也存在一定的不均匀性。考虑到大跨空间结构中杆件截面尺寸普遍较小，其表面温度的不均匀度并不十分显著，采用截面平均温度作为温度荷载已有不错的代表性。同时，大型结构杆件数量往往以万计，若采用壳单元对每根杆件进行精细建模，遮挡判断及温度场的计算代价会变得无法估计。因此，综合考虑计算精度与效率，本文采用杆单元对结构温度场进行高效的计算模拟。

模型中的材料仅有钢材，基本热工参数如表2所示。其中结构表面涂装了灰白色面漆，太阳辐射系数率α与辐射发射率ε通过查表[1]得到。

表2 材料热工参数Table.2 Thermal characteristic parameters of material

3.2 模型变量及求解条件说明

环境参数方面，结构所在的施工场地经度为116.4°E，纬度为39.5°N，海拔高度55 m。

为使结构温度场的特征及变化规律尽可能明显，取施工过程中经历的某个典型夏至日(2017年6 月23日)作为计算日期，计算结构0 h~24 h的温度场。经查询历史气象信息，取当天晴朗情况下的最低气温为23 ℃，最高气温35 ℃。风速取气象记录中的日均风速5.4 m/s。

需要注意的是，在瞬态温度场的求解过程中，重要求解条件如计算的初始状态、计算时间间隔等均会对结果的精度及计算效率产生不小的影响。因此，选取边界条件的重复天数、计算时间间隔、构件遮挡因素这三个重要求解条件作为变量，进行结构温度场的时程分析，并根据计算结果给出建议的取值。

3.3 边界条件重复天数

计算开始时，结构中的温度场是未知且不均匀的，无法赋予其准确的起始条件。有部分研究假设结构在夜间温度与气温相同，对结构施加均匀温度场作为时程计算的初始条件。但如果结构的散热效率不够高，夜间温度场仍不均匀，则这种计算方式会使结果与实际存在一定偏差。

若以天为单位重复施加一段时间的热边界条件，人为设置初始状态的前序求解过程，可以使初始值引入的偏差在前序求解过程中逐渐消除，后续真正需要的计算结果收敛于真实值，得到正确的初始状态与分析结果。按此思路设置了不同的计算工况，典型构件的结果如图6所示。

图6 边界条件重复天数对结果的影响Fig.6 Effect of repeated days of boundary conditions on the results

可以发现，当不设置前序求解过程，以气温作为初始状态时，前几个小时的计算结果确实存在明显偏差，构件温度以较快的速度偏离气温，构件散热效率并不理想。重复1天与重复2天的计算结果重合，说明以同样的边界条件重复计算1天后，结构中各构件温度已经达到收敛，形成闭环。

根据以上结果，进行后续的瞬态温度场模拟时，只需复制第一天的边界条件形成前序求解过程，即可得到所需时间区段内准确的温度场。

3.4 计算时间间隔

有限元在进行瞬态分析时一般采用如中心差分法、Newmark时间积分法等数值分析方法，计算结果仅在离散的时间点上满足瞬态平衡方程。这就使得时间间隔越小，计算结果总体上越逼近实际曲线。但同时计算成本也会随着时间间隔的减小而迅速提高。为尽可能找到计算精度与效率的平衡点，在相同条件下设置了不同的计算时间间隔对结构进行分析，典型构件的结果如图7所示。

结果显示，当取瞬态分析的时间间隔为2 h时，各离散点之间的连线有明显的转折，结果较为粗糙，特征点丢失的可能性过高。间隔为1 h与0.5 h时，计算结果非常接近，曲线基本重合，且在温度变化迅速的时段有相对平滑的曲线。

图7 计算时间间隔对结果的影响Fig.7 Effect of time interval of calculation on the results

综合考虑精度与计算效率，建议取瞬态计算的时间间隔为1 h。

3.5 构件遮挡因素

如前文所述，造型复杂的大跨空间结构中一般存在较为明显的日照遮挡关系。构件的相互遮挡会减弱太阳直射辐射的作用。若考虑日照阴影带来的影响，计算结果会更加精确。但同时遮挡关系的计算过程较为复杂，会在一定程度上影响计算效率。为验证是否有必要考虑单元间的遮挡关系，在考虑遮挡和无遮挡两种情况下对结构进行了计算分析，典型构件的结果如图8所示。

图8 构件遮挡对结果的影响Fig.8 Effect of member occlusion on the results

结果显示，考虑构件之间的遮挡效果后，结构中各构件的温度明显降低，且不同位置构件的遮挡效果有所不同：位于结构底部的C型柱受到上方多层构件的遮挡，温度降低幅度最大约为 9%；顶层的下弦杆受遮挡效果较弱，降温幅度也有约 5%，不可忽略。同时，考虑日照阴影后构件各时刻吸收的太阳辐射量也有明显降低，从图中可以发现，典型C型柱构件在大部分时间被遮挡的有效长度超过长度的一半。

根据以上结果，构件遮挡对温度计算结果的影响较为显著，不可忽略。

4 结构温度场主要计算结果

对求解过程中的重要求解条件按前文给出的建议取值后，对结构瞬态温度场进行计算，主要计算结果如下。

4.1 构件日照系数

在一天中的不同时刻，太阳所在的位置不同，结构中的构件遮挡情况也并不相同。图9给出了上午8 h结构中各构件的日照系数分布情况，其中：0代表构件被完全遮挡；1代表构件表面没有日照阴影；0~1之间被照射的长度比例越高，表示构件被遮挡的长度比例越高。

从图9可以看出，在上午8 h，结构顶部的上弦杆基本均未被遮挡，而下部的C型柱、支撑筒等主要支承构件的被遮挡程度较高，与预想分布规律相同。值得注意的是，在C型柱中部出现一条明显的遮挡与未遮挡区域的分界线，这是由于这一时刻太阳的高度相对较低，结构顶部的悬挑在倾斜的太阳光线下投下了较长的阴影，遮挡住了下侧支承构件高度较高的部分。这一结果从侧面印证了本文日照阴影计算方法的准确性。

图9 8:00各构件的日照系数分布Fig.9 Distribution of sunshine coefficient of each member at 8:00

4.2 温度的时间变化特性

结构的温度变化主要受太阳辐射与气温波动的影响，而后两者随时间的变化是有一定规律的。因此，结构温度在一天内也会呈现出一定的时间变化特性。

结构中最高温度随时间的变化关系见图10，其变化趋势与气温有相似之处。5 h左右即日出前，结构总体温度达到最低，此时结构最高温度为24.92 ℃，比同一时刻的气温高约1.92 ℃。日出后，构件在太阳辐射的作用下开始升温。虽然太阳辐射在12 h左右达到当天最强，但由于结构的散热效率并不高，因此，总体温度还会继续上升，温度峰值较太阳辐射峰值出现一定的滞后。在当天15:00左右，结构总体温度达到最高，此时，结构中的最高温度为 44.60 ℃，比同一时刻的气温高约9.62 ℃。在此之后，结构的温度开始下降，同时降温速度逐渐减小，直至第二天日出之前再次达到最低，完成一个变化周期。

4.3 温度的空间分布特性

大跨空间结构构件数量多，造型复杂，结构温度场在空间上的分布一般并不均匀。结构总体温度达到最低与最高的两个时刻，结构中的温度场分布情况如图11所示。

图10 结构最高温度随时间变化曲线Fig.10 Maximum temperature in structure versus time

图11 结构温度场分布Fig.11 Distribution of structure temperature field

结果显示，整体温度最低时，结构顶部的构件温度基本接近气温，但下部构件的温度较高；整体温度最高时，结构顶部的构件温度最高，而下部的温度则较低。这是因为结构下部以支承构件为主，C型柱、支撑筒等部分的截面尺寸要明显更大，构件比表面积小，对流换热系数也较小，故温度变化较为缓慢，明显滞后于气温的变化；而结构上部的构件以弦杆为主，构件截面细小，升温降温都非常迅速。白天结构上部吸收了大量的太阳辐射，加速了上部构件在日出后的升温速度；夜间失去热源，结构下部构件散热缓慢，赶不上气温下降的速度，因此，温度比气温略高。

从温度场的总体计算结果中可以发现，结构中温度在空间上的分布规律并不明显，且分布情况会受构件的空间位置、遮挡关系、截面尺寸等因素影响，不均匀性较为显著。对于规模较大、形状复杂、重要性高的结构，建议对温度场进行精细模拟，以确定用于设计校核的温度荷载。对于规模较小的结构，在估计其温度分布时，也要考虑构件尺寸所带来的温度不均匀性。

5 结论

本文实现了适用于大跨空间钢结构的日照温度场模拟方法，并以北京大兴国际机场不规则自由曲面钢网格结构为例，进行了瞬态温度场的数值模拟，为大跨空间结构施工过程的精确模拟提供了重要基础。主要结论如下：

(1)在温度场的瞬态求解过程中，经过24 h的前序求解计算即可消除近似起始条件带来的计算误差，使结果收敛于准确解。综合考虑结果精度与计算效率，可取计算时间间隔为1 h。

(2)考虑日照阴影后，典型构件的温度降低幅度最大约为 9%，构件遮挡对温度计算结果的影响较为显著，不可忽略。

(3)大跨空间钢结构的温度场具有明显的时变特征与空间非均匀性，夏至日的最低温度发生于 5时左右，最高温度发生于 15时左右。不同尺寸构件的温度变化特征有所不同。

(4)对于规模较大、形状复杂、重要性高的结构，建议在前期设计时就对结构的温度场进行精确模拟分析，考虑不同阶段的温度场变化在结构中引起的温度荷载，以提高结构在施工、使用过程中的安全性。