火箭仪器舱结构模型的热振耦合模拟方法

2024-03-11薛景天董龙雷赵建平丁镇军王潇屹

航天器环境工程 2024年1期

薛景天，董龙雷*，刘振，赵建平，丁镇军，王潇屹

（1.西安交通大学复杂服役环境重大装备结构强度与寿命全国重点实验室，西安 710049；2.北京强度环境研究所，北京 100076； 3.火箭军某部，西安 710407）

0 引言

在火箭飞行过程中，其机体结构表面与高速气流摩擦会产生气动加热现象。气动加热环境对火箭的影响主要体现在两方面：一是可能改变火箭结构本身的动力学特性；二是向内传热使结构内部的仪器或元器件升温甚至受损失效，从而影响火箭的可靠性。目前，针对火箭本体结构的气动加热防护已有相对成熟有效的措施[1]。而热环境对航天器结构的动力学特性影响主要体现在对结构固有频率的影响，Heeg 等[2]对典型的航天飞机真实飞行状态进行研究并发现，航天飞机表面温度可以从常温升至2700 ℃以上，相应地，结构固有频率及阻尼也可剧烈变化高达30%以上。工程中常利用数值仿真与试验相结合的方法对结构进行热振分析[3-10]，但在有关研究中普遍存在的问题是：进行传热分析时将结构材料的热特性参数视为常量，但实际上材料的热特性参数（如热导率和比热容）在传热过程中并非恒定不变，而是会随温度的变化发生改变；即使考虑到此影响而将这些材料参数设为变量，依然无法保证真实反映材料参数的变化，导致热振耦合特性预示出现偏差。

本文提出一种材料热特性参数与振动动态耦合的分析方法，在进行结构传热分析中考虑材料热特性随温度的变化，同时结合振动方程，对某火箭仪器舱气动过程进行热振耦合机理分析。

1 热振耦合机理分析

在火箭的飞行过程中，针对热环境对火箭结构的影响：一方面要关注固有频率对结构动力学特性的影响；另一方面要关注结构传热过程，从而对结构内部元器件和部件的温度变化做出预示。常规的结构动力学分析是从振动方程出发[7]。然而，本文所讨论的是结构的固有特性，与外力无关，并且为简化模型而忽略结构自身的阻尼项，因此振动方程中的阻尼和激励项为0，对结构的模态分析可以围绕公式

展开，其中：K为结构的总刚度矩阵；M为结构的质量矩阵。

加入温度场变化的主要影响体现在改变振动方程中的刚度项：一方面，温度的改变使材料的弹性模量发生改变，进而影响结构的刚度，即

式中：B为几何矩阵；DT为与材料弹性模量ET和泊松比μ有关的弹性矩阵。

另外一方面，当结构的温度变化时会产生热应力并影响结构的刚度，即

式中：G为形函数矩阵；Γ为结构热应力矩阵。

在实际试验过程中，温度场和应力场是同时存在的，故通过模态叠加后结构的刚度可以表示为

或

对于某确定结构而言，其泊松比μ是不变的，而弹性模量ET以及结构热应力Γ均随着结构温度T的变化而发生改变，因此掌握结构温度的变化规律是研究热振耦合机理的关键。

在飞行过程中，飞行器结构表面与气体摩擦产生的热量是通过热传导的方式传递到结构内部的。热传导方程是典型的抛物线型二阶偏微分方程（见式(6)），

式中：λ为材料的热导率；ρ为材料的密度；cp为材料的比热容。式(6)为笛卡儿坐标下的非稳态导热方程，反映的是三维非稳态热传导的传热情况。由式(6)可推导出

式中：T为结构某点的温度；x,y,z为该点的坐标位置；τ为传热时间。式(7)揭示了温度与位置和传热时间的关系，可以得到某一时刻结构的热分布以及在传热过程中结构某点的温度变化。

温度T与弹性模量ET间存在映射关系，而材料的热特性参数（如热导率λ和比热容cp）会随温度的变化而发生改变。根据温度-弹性模量以及温度-热特性参数这2 组对应关系，可以建立材料的弹性模量和热特性参数之间的对应关系，

综上，热振耦合机理的研究是从热环境对结构的刚度影响出发，关注影响结构刚度参数中随温度的变化而发生改变的量，从而对热环境下结构的振动特性做出预示。因此，考虑热特性参数变化的仿真计算更加真实、可靠性更高。

2 闪光瞬态法测量材料热特性

为在热振耦合分析中引入材料热特性的变化，需要得到材料热特性与温度之间的对应关系。目前，主要利用稳态法和瞬态法测量材料的热导率[11]。稳态法根据Fourie 方程直接测量热导率，但该方法适用的温度和热导率范围均较窄，主要是在中等温度下测量中低热导率材料。瞬态法则应用范围较为宽广，尤其适合于高热导率材料在高温下的测试。闪光瞬态法测量原理见图1[12]。

图1 闪光瞬态法测量原理示意[12]Fig.1 Measurement principle schematic of flash transient method[12]

如图1 所示，通过仪器向待测试样底部发出激光束，由试样上方的红外检测器跟踪监测试样的升温过程，得到电信号与升温时间的关系曲线，通过转化得到试样材料热扩散系数与升温时间的关系曲线。

式中：α为试样材料的热扩散系数；d为试样厚度；t0为半升温时间，即试样温度达到设定温度一半所用的时间。

得到材料热扩散系数后，根据热扩散系数与热导率及比热容的关系，利用已知材料参数的标样共同进行测试，可推导得到待测试样材料的热导率及比热容。具体方法为：使用一个与待测试样截面形状相同、厚度相近、热物性相近、表面光滑程度相同且比热容值已知的参比标样，与待测试样同时进行表面涂覆；在理想的绝热条件下，依次进行测量得到2 条测试曲线。由于2 试样加热时间相同，因此吸收的热量可以认为是相同的，故可通过公式(10)得到待测试样的比热容。

计算出材料的热导率λ。

3 数值分析

热振耦合机理的研究通常包括以下步骤：1）设计模型并进行数值分析；2）设计试验装置及夹具，并制定试验条件；3）完成热振试验；4）对试验结果进行数据处理，包括材料的热振耦合特性预示。本文详述开展模型设计并进行数值分析的过程，其余步骤会在之后的文章中进行论述。

3.1 试验模型设计

以某型号火箭仪器舱为研究对象，建立其几何模型，如图2 所示。

图2 火箭仪器舱筒形舱段三维模型Fig.2 Three-dimensional model of the rocket cylindrical instrument cabin

筒形舱段由筒与十字梁2 部分结构组成，其中十字梁用来安装火箭内部部件及其他元器件。对该模型进行数值仿真计算，筒的材料为钢，通过闪光瞬态法测得的结构材料参数如表1 所示。

表1 火箭仪器舱筒形舱段的材料参数Table 1 Material parameters of the rocket cylindrical instrument cabin

3.2 数值计算

3.2.1 传热计算

本文对筒形舱段模型的处理作出以下假设：

1）筒外侧受到均匀温度场的影响，即结构纵向不存在温度梯度，可简化为二维平面；

2）筒的厚度远小于舱段的长度，即可忽略筒的厚度，将舱段简化为薄壁壳结构；

3）舱段结构为对称结构，在传热过程中十字梁结构的温度分布具有对称性，可简化为矩形平板。

经上述简化的筒形舱段传热模型如图3 所示，传热过程被简化为第一类边界条件的一维热传导问题：给定模型外侧固定温度值，分析中间矩形结构部分的温度与传热时间的关系，得到

图3 筒形舱段简化传热模型Fig.3 Simplified heat transfer model of the cylindrical cabin

式中：T为结构某点的温度；x为该点的坐标位置；τ为传热时间。

在实际传热过程中，筒外侧边界条件为恒温400 ℃，结构其他部分的初始温度为20 ℃，即传热的边界条件为

式中u(x,τ)为与位置和时间相关的传热方程，l为筒纵向长度。

采用分离变量法求解该传热模型的偏微分方程得到

通过式(14)可以得到模型上各点的温度变化，即结构上任意位置、任意时刻的温度值。

在热振耦合试验中，通常使用石英灯阵对试件进行温度加载来模拟飞行器在飞行过程中所受到的气动加热。因此需要对所使用的石英灯的个数、功率以及加热时间等进行传热仿真分析。同时，需要通过仿真模态分析对试验结果的可靠性进行验证。本研究中，将由闪光瞬态法得到的材料热特性参数代入有限元计算软件ABAQUS 进行仿真计算得到结构的温度分布。热量以热传导方式从筒外侧传递到结构内部，在十字梁上选取3 个位置的特征点与理论模型计算的结果进行对比，特征点（T1～T3）的位置及结构网格如图4 所示。在考虑及不考虑材料热特性（即材料的热导率和比热容）随温度变化而改变的情况下，计算得到50 000 s 传热时间内3 个位置点的温度变化曲线，如图5～图7 所示。

图4 模型中的特征测点位置Fig.4 Positions of characteristic measuring points on the model

图5 T1 位置的温度曲线Fig.5 Temperature curve for T1

图6 T2 位置的温度曲线Fig.6 Temperature curve for T2

图7 T3 位置的温度曲线Fig.7 Temperature curve for T3

可以看出，是否考虑材料参数随温度的变化对传热过程的仿真计算结果有很大的影响，并且越靠近结构中心（即实际安装仪器位置）的误差越大。

3.2.2 热振耦合仿真计算

热环境对结构模态的影响主要体现在对结构刚度的影响，因此结构刚度在热环境中的变化规律能够反映出结构的热振耦合作用机理。用有限元方法对结构进行模态分析时会对结构的刚度矩阵和质量矩阵进行计算，对比不同温度场下输出的结构刚度矩阵可得到结构在热环境中的刚度变化。由于输出矩阵的单元数过多，选出前50 个单元的刚度矩阵大小绘制曲线，如图8 所示。

图8 不同温度下的结构单元刚度曲线Fig.8 Stiffness curves of structural units at different temperatures

从图8 可以看到，典型结构的前50 个单元的刚度随温度变化发生明显改变：一方面，刚度曲线的峰值发生了改变，总体表现在随着温度的升高，刚度的峰值降低，符合仿真分析计算结果；另一方面，刚度曲线峰值的位置发生了改变，甚至有峰值缺失或者增大现象。由此可见，温度对结构刚度的影响不仅仅体现在数值大小的变化上，还会改变结构刚度的最值位置。因此，在实际结构设计中，通过刚度计算可以预示在不同温度场下的结构刚度薄弱位置，为更深层次的热振耦合机理研究奠定基础。

将计算的温度场结果作为预定义场导入热应力计算中，在温度场和热应力场叠加作用下进行模态计算，得到结构的振型和固有频率。在进行网格无关性检验后得到的网格数计算结果收敛性良好。设定4 种温度工况（20 ℃、200 ℃、300 ℃和400 ℃），得到温度场、应力场以及2 个场共同作用对结构固有频率的影响，结果表明模态振型没有明显变化，不同温度下结构的前6 阶固有频率如表2 所示。