考虑温度的加筋板结构随机振动疲劳寿命预测方法研究

Study on Fatigue Life Prediction Method of Stiffened Plate Structure under Random Vibration Considering Temperature

摘要: 高超音速飞行器在高速飞行时，机体的表面会达到非常高的温度，结构材料在高温条件下其力学性能会发生重大变化。在热应力、外部激励力的联合作用下结构件会发生一定的疲劳破坏。本文对飞行器加筋壁板的动热特性进行了分析，应用有限元方法进行数值模拟，得到飞行器壁板的温度分布，研究了温度对固有振动频率和固有振型的影响以及对疲劳寿命的影响。应用频域方法估算了飞行器壁板的疲劳寿命，并与时域方法估算出的飞行器壁板的疲劳寿命进行了比较。

Abstract: When the hypersonic vehicle flies at high speed, the surface of the airframe will reach a very high temperature, and the mechanical properties of the structure material will change greatly. Under the combined action of thermal stress and external excitation, fatigue failure will occur to some extent. In this paper, the dynamic and thermal characteristics of the stiffened panel of an aircraft are analyzed, and the finite element method is used for numerical simulation. Then the temperature distribution of the panel is obtained, the effect of temperature on natural vibration frequency, natural vibration mode and fatigue life is studied. The fatigue life of aircraft wallboard is estimated by using the frequency domain method and compared with the time-domain method.

1. 引言

高超音速飞行器是飞行速度高于5倍音速，即5马赫(MH)，能够长距离短时间到达目的地执行攻击、运输和侦察等任务的飞行器。高超音速飞行器的加筋壁板属于飞行器的热防护系统，加筋壁板的服役环境非常复杂。飞行器在巡航过程中，结构会承受较大的振动载荷，对结构产生破坏。

最早提出振动疲劳概念的是上个世纪六十年代的学者S. H. Crandall，他提出振动疲劳是在振动载荷作用下使结构产生不可逆且累积的损伤或破坏 [1]。1981年，Valanis在内时理论的基础上，推导出了机械结构疲劳寿命和激励频率之间的内在关系 [2]。1992年，Sanliturk基于动力学理论，在考虑了结构弹性对疲劳寿命、结构惯性和阻尼因子对疲劳寿命的影响的作用下，提出了一种基于频响函数的振动疲劳寿命估算方法 [3]。2003年，Shiryayev等针对非线性振动系统疲劳试验，提出了一种跟踪固有频率的振动疲劳试验技术，此技术为线性振动中固有频率变化提供了更可靠的数据 [4]。在国内，我国的学者姚起杭等人在上世纪末开创了结构振动疲劳的研究。提出了结构在动态载荷作用下的结构振动疲劳分析方法，同时在结构抗振动疲劳设计、振动疲劳寿命分析等方面已取得了很多优秀成果 [5] [6] [7]。2001年，张积亭提出了一种在随机振动环境下疲劳寿命预测的数据处理方法，和一种利用随机响应功率谱密度求出的特征频率作为平均频率去求疲劳损伤量的方法 [8]。2006年，杨智春等人归纳了目前研究高超飞行器壁板颤振中的几个重要结论，并且提出了高速飞行器壁板颤振在今后的研究中需要解决的一些问题 [9]。在振动疲劳其他方面的研究，陈敏亮 [10]，孙伟 [11]，曹明红 [12] 也做出了一系列工作，陈敏亮提出了飞机结构的随机振动疲劳分析方法，孙伟提出了结构振动疲劳寿命估算方法研究，曹明红提出了涉及双模态应力响应谱的振动疲劳寿命估算方法，一并推动了结构振动疲劳理论的发展历程。

本文对飞行器加筋壁板的动热特性进行了分析，应用有限元方法进行数值模拟，基于理论和有限元的方法得到飞行器壁板的温度分布，研究了温度对固有振动频率和固有振型的影响以及对疲劳寿命的影响。应用频域方法估算了飞行器壁板的疲劳寿命，并与时域方法估算出的飞行器壁板的疲劳寿命进行了比较。

2. 理论分析

2.1. 结构热模态理论

结构瞬态温度场及应力场是结构热模态分析的基础，因此在分析结构热模态之前，必须要得到结构瞬态温度场及应力场。结构热模态一般采用数值有限元分析的方法。基于热传导方程分析传热过程的初始条件和边界条件，并结合材料的热物理性能参数得到结构的温度场分布；进而通过热应力方程结合有限元分析得到结构热应力分布。

根据热能守恒定理，可以推导出热传导固体热传导微分方程。依据变分原理，对热传导微分方程和边界条件进行变换，得到有限元计算方程 [13]：

$C\stackrel{\dot{}}{T}+KT=P$ (1)

式中C是热容矩阵，K是热传导矩阵，P是温度载荷列阵，T是结构温度矩阵。

瞬态温度场有限元法是以温度场有限元为基础，结合时间有限元差分格式得到的。采用差分格式得到的有限元方程在t时刻的求解过程如下：

${\frac{\partial T}{\partial t}|}_t=\frac{1}{\Delta \tau }\left(T_t-T_{t-\Delta t}\right)+{\rm O}\left(\Delta \tau \right)$ (2)

将(2)式代入(1)式中，可得到瞬态温度场的有限元方程

$\left(\frac{C_t}{\Delta \tau }+K_t\right)T_t=P_t+\frac{K_t}{\Delta \tau }{T}_{t-\Delta t}$ (3)

式中 $T_{t-\Delta t}$ 是初始温度场或前一时刻已知的温度场。在已知结构边界条件和初始温度场的情况下，可以计算出后续各时刻的结构瞬态温度场。模态是结构固有的特性，在不考虑阻尼效应影响的情况下，系统的自由振动方程为：

$M\ddot{x}+Kx=0\Rightarrow \left(K-{\omega }^{2}M\right)\Phi =0$ (4)

式中K是结构刚度矩阵，M是结构质量矩阵， $\Phi$ 是模态矩阵， $\omega$ 是与向量 $\Phi$ 对应的振动频率。

2.2. 基于功率谱密度的结构随机振动分析

超高音速飞行器在飞行过程中受到来自发动机和空气随机振动载荷的作用，从随机振动随时间的变化规律来看符合平稳随机振动，即随机过程按统计特性随时间平稳变化，工程上一般采用相关函数的方法在时域范围内进行描述。而功率谱密度是在频域范围内对随机过程的描述，能够反映随机振动过程中能量的分布情况。理论上，自相关函数R_x(τ)和功率谱密度S_x(w)符合维纳–辛钦定理，即

$S_x\left(\omega \right)=\frac{1}{2\pi }{\displaystyle {\int }_{-\infty }^{\infty }R\left(\tau \right)}{e}^{-i\omega t}d\tau$，

$R_x\left(\tau \right)={\displaystyle {\int }_{-\infty }^{\infty }{S}_x\left(\omega \right)}{e}^{i\omega t}d\omega$。

功率谱密度分为单边功率谱和双边功率谱，在实际工程中的振动频率大于零，所以定义单边功率谱密度为：

$G\left(\omega \right)=\{\begin{array}{l}2S_x\left(\omega \right),\omega \ge 0\\ S_x\left(0\right),\omega =0\end{array}$

统计矩是描述概率密度的数字特征，在频域特征参数中引入谱矩来描述随机过程功率谱密度的统计特征参数。为了方便工程运算，根据平稳随机过程 $X\left(t\right)$ 的谱矩 $m_i$ 用单边功率谱密度函数定义：

$m_i={\displaystyle \underset{0}{\overset{\infty }{\int }}{f}^i{W}_x\left(f\right)df}\left(i=0,1,2,3,\cdots \right)$ (5)

2.3. 基于频域分析的加筋板疲劳寿命估算

目前广泛采用的频域方法是Dirlik方法 [14]。该方法在大量Monte Carlo 模拟基础上，得到一个针对简单循环和雨流循环经验表达式，即用随机载荷的4阶矩构造了应力幅值概率密度函数。该方法通过对比循环分布，不必作S-N曲线的斜率的假设。这个经验关系经一些研究者验证，认为具有较高的精度。

Dirlik所使用的方法认为雨流法获得的范围分布是三种分布的总和：一个指数函数和两个瑞利(Rayleigh)函数。假设雨流法得到的循环范围的概率密度分布可由以下得到：

$f_{S_a}\left(S_a\right)=\frac{1}{\sqrt{m_0}}\left(\frac{D_1}{Q}\right)e^{\frac{-z}{Q}{S}_a}+\frac{D_2}{R^2}\left(\frac{1}{{\alpha }^2}\right)e^{\frac{-z^2}{2R^2}{S}_a^2}+D_{3}z{e}^{\frac{-z^2}{2}{S}_a^2}$ (6)

其中：

$z=\frac{1}{\sqrt{m_0}}$ (7)

$\gamma =\frac{m_{\text{2}}}{\sqrt{m_0{m}_4}}$ (8)

$D_1=\frac{\text{2}\left(X_m-{\gamma }^2\right)}{1+{\gamma }^2}$ (9)

$D_2=\frac{1-\gamma -D_1+D_1^2}{1-R}$ (10)

$D_3=1-D_1-D_2$ (11)

$X_m=\frac{m_1}{m_0}\sqrt{\frac{m_2}{m_4}}$ (12)

$R=\frac{\gamma -X_m-D_1^2}{1-\gamma -D_1+D_1^2}$ (13)

$Q=\frac{\text{1}\text{.25}\left(\gamma -D_3-D_{2}R\right)}{D_{\text{1}}}$ (14)

在时间段T内，Dirlik疲劳损伤模型如下：

$D=\frac{E\left[P\right]T}{c}{\displaystyle \underset{0}{\overset{\infty }{\int }}{S}_a^m{f}_{S_a}\left(S_a\right)d{S}_a}$ (15)

3. 数值模拟与寿命计算

为了得到高超音速飞行器表面温度载荷以及随机振动载荷耦合作用下的加筋板危险点的应力变化规律，采用ANSYS Workbench有限元分析软件建立加筋板的热力耦合有限元分析模型，通过定义材料参数、添加固定约束，施加表面温度载荷及随机振动载荷，得到加筋板结构在不同材料系数、温度载荷及随机振动幅值条件下加筋板的模态变化规律，应力应变在频域的响应曲线，具体的考虑温度和随机振动载荷的有限元分析流程如图1所示。

Figure 1. Finite element analysis considering temperature and random vibration load

图1. 考虑温度和随机振动载荷的有限元分析

3.1. 有限元模型建立

飞行器加筋板分为底板、过渡板和侧板三大区域，为了提高加筋板的结构强度和使用寿命，在底板上等间距焊接4段加筋杆，整个加筋板通过侧板的8对螺栓孔与飞行器栓接，如图2所示。飞行器壁板结构尺寸为863.6 mm × 508 mm，加筋厚度为1.65 mm，筋条间隔距离为254 mm，底板厚度为1.65 mm。使用三维建模软件Pro/E建立有限元模型输出STP交互文件格式，导入ANSYS Workbench中进行有限元分析，划分网格后得到27497个节点，3660个单元。

Figure 2. Geometric model of stiffened plate of aircraft

图2. 飞行器加筋板几何模型

为了避免温度载荷与随机振动载荷耦合作用过程中的节点畸变，在螺栓孔的内表面采用Sizing的约束方式细化局部网格，将加筋杆与底板之间设置为Bonded连接，并采用Fixed support的方式在螺栓孔位置施加固定约束。本文采用镍基高温合金GH4169作为加筋板研究对象。GH4169在650℃以下具有高的强度、塑形、持久性能和疲劳性能。此外，还具有良好的焊接性能，以及耐腐蚀、抗氧化和耐辐照性能。该合金在−253至700℃温度范围内组织性能稳定。在航空航天、核能等工业中得到了广泛的应用。GH4169的密度为ρ = $8.24\times {10}^3\text{kg}/m^3$，泊松比为v = 0.321，其他物理参数如表1所示。

Table 1. Thermophysical parameters of GH4169

表1. GH4169热物性参数

从表1可以看出，高超音速飞行器加筋板的温度从100℃到1000℃变化过程中，加筋板的弹性模量随温度的升高而降低，材料导热率、比热容和泊松比随温度的升高而增大。

3.2. 热模态分析以及固有频率随温度的变化规律

在ANSYS仿真分析过程中，模态分析是获得随机振动结果应力谱的基础，而飞行器加筋板上底板持续受到高温载荷的影响，底板高温载荷的持续作用将使得整个加筋板在不同温度区域产生不同的热应力分布，将对结构的模态特性产生很大的影响，图3所示为底板100℃温度载荷条件下的温度场分布，此时对应的沿竖直方向的前6阶模态振形如图4所示。

Figure 3. Temperature distribution under temperature loading at 100℃

图3. 100℃温度载荷时的温度分布

Figure 4. First six modes of vibration under temperature load at 100℃

图4. 100℃温度载荷时的前6阶模态振形

图3中的温度外载荷作用在加筋板的底板面上，采用线性加载的方式在900 s的时间区间内将底板面表面温度从22℃变化为100℃，从图中可以看出加筋杆上的温度已经上升到90℃左右，而侧板区域依旧是温度最低的部位。从图4中可以看出当结构受到100摄氏度的温度载荷时，加筋板的前六阶模态集中在100 Hz到200 Hz之间，第一阶固有频率跟第六阶固有频率之间相差70 Hz左右。为了得到温度载荷对固有模态的影响规律，将底板的表面温度从22℃逐渐增大到600℃左右，得到的加筋板前6阶固有频率随温度的变化曲线如图5所示，加筋板的固有频率随温度先减小后增大，表明材料结构的热屈曲温度在70℃左右，然后在150℃左右开始基本保持稳定，增大速率缓慢。

Figure 5. Variation of the first 6 natural frequencies of stiffened plate with temperature

图5. 加筋板前6阶固有频率随温度的变化

3.3. 随机振动分析及PSD曲线的获取

随机振动分析的目的是得到危险点的应力在频域内的变化规律，继而得到功率谱密度曲线，仿真过程中的随机振动分析是以热模态分析的结果为基础，并在约束区域施加输入载荷谱后的耦合结果。加筋板的热模态随机振动分析以100℃下的热温度场为基础，在螺栓孔上施加沿Z轴方向的加速度载荷谱，如图6所示。

Figure 6. Random vibration analysis of stiffened plate

图6. 加筋板随机振动分析

在图6中，加速度载荷谱的施加方式是梯形曲线，在100℃的温度载荷及随机振动条件下，分析得到危险点的位置位于加筋板的螺栓孔内侧，从危险点的PSD响应曲线可以看出两个峰值点对应的频率分别为112.59 Hz和150.98 Hz，与热模态分析得到的第3阶和第4阶频率基本保持一致，所对应的功率谱密度分别为54 [(MPa)²/Hz]和81 + [(MPa)²/Hz]。基于频域分析的疲劳寿命计算过程中，危险点的PSD响应曲线的峰值点对应频率和峰值大小是加筋板结构寿命计算的必要参数，而危险点的实际温度大小决定了其疲劳寿命曲线的变化规律。因此，首先需要通过有限元仿真建立随机振动载荷条件下，危险点的功率谱密度曲线频率、峰值点及实际温度的对应关系，如图7~图9所示。

表1危险点PSD曲线峰值点频率、温度及大小随表面出现温度变化。

Figure 7. The frequency corresponding to the temperature peak point of different base plate

图7. 不同底板温度峰值点对应的频率

Figure 8. PSD value of temperature peak point of different base plate

图8. 不同底板温度峰值点的PSD值

Figure 9. Peak point temperature corresponding to different base plate temperature

图9. 不同底板温度对应的峰值点温度

对称循环下的S-N曲线(非对称可通过Goodman方法转化成对称)，其表达式为：

$S_{a_i}^m{N}_{fi}=c$ (16)

表达式中的应力参数S_a和m将随危险点温度的变化而变化，对于高超音速飞行器的材料GH4169，参数m随温度变化不明显，此处设定为常数m = −14.3，而应力参数S_a随温度的变化关系为：

$S_a=-0.77*T+1973.8$ (17)

式中，T表示危险点的温度。

3.4. 随机振动疲劳寿命估算

结构随机振动分析得到应力响应功率谱密度如图6中的PSD曲线所示，可以看出结构随机振动响应谱为呈现双峰特点。通过提取PSD曲线峰值点对应的频率、峰值大小以及相应的温度值，结合Dirlik雨流幅值分布模型得到不同温度载荷下高超音速加筋板的随机振动寿命曲线变化规律如图10所示。

考虑到加筋板底板温度变动范围过大，因此以不同工况的寿命值的对数作为纵坐标，可以看出在100℃以内的区域，加筋板的寿命随温度的升高而逐渐降低满足非线性关系，且疲劳寿命的降低速度随温度的升高而加快。

Figure 10. Variation of fatigue life of stiffened plate with temperature

图10. 加筋板的疲劳寿命随温度的变化规律

4. 结论

本文以高超音速飞行器加筋板结构为研究对象，考虑飞行过程中温度载荷及随机振动载荷的耦合作用，建立了考虑温度载荷的加筋板随机振动寿命计算模型，通过BP神经网络预测理论得到了特定工况下加筋板结构参数同疲劳寿命之间的关系，为高超音速飞行器的疲劳寿命分析提供了理论基础。

(1) 本文通过构建加筋板的温度场和随机振动载荷模型，结合ANSYS Workbench有限元仿真技术，对高超音速飞行器的疲劳寿命进行了数值仿真分析，找到了加筋板在实际工况下的危险点位置、应力的PSD曲线以及实际使用寿命，对高超音速飞行器的材料选取和寿命研究具有一定的指导意义；

(2) 从加筋板的前6阶热模态变化曲线可以看出，100℃以内加筋板的固有频率随温度先减小后增大，且结构的热屈曲温度在70℃左右，然后固有频率在150℃左右开始基本保持稳定，增速缓慢；

(3) 从结构随机振动分析得到的危险点功率谱密度曲线可以看出结构随机振动响应谱为呈现双峰特点；

(4) 在随机振动载荷一定的条件下，加筋板的疲劳寿命随温度的增大而逐渐减小且满足非线性关系，且疲劳寿命的降低速度随温度的升高而加快。