CN103577648A - 运输类飞机货物空投时机翼结构载荷的确定方法 - Google Patents

Abstract

本发明属于飞机飞行载荷设计技术领域，涉及对运输类飞机货物空投时机翼结构载荷确定方法的改进。其特征在于，确定机翼结构载荷的步骤如下：飞机空投货物参数初始化；飞机有限元模态分析；飞机非定常气动力系数矩阵计算；飞机空投货物动响应计算；计算机翼结构载荷动响应；确定空投时机翼结构的设计载荷。本发明提供了一种改进的运输类飞机货物空投时机翼结构载荷的确定方法，提高了机翼结构载荷的确定精度，保证了机翼结构的安全性和飞机的飞行安全。

Description

运输类飞机货物空投时机翼结构载荷的确定方法

技术领域

本发明属于飞机飞行载荷设计技术领域，涉及对运输类飞机货物空投时机翼结构载荷的确定方法的改进。

背景技术

空投是运输类飞机的重要工作模式之一。空投时一般需要确定飞机货舱地板、货桥、连接件、机翼等部件的设计载荷。当前的运输类飞机货物空投时将飞机视为刚体，仅考虑气动力的静气动弹性修正，依托刚体飞行动力学方程进行求解，具体技术见文献：“运输机空投的飞行动力学建模及仿真”,杨妙升、屈香菊,《飞行力学》,2010，28卷(第3期),P.9-12。

目前的运输类飞机货物空投时机翼结构载荷的确定方法的缺点是：第一、所确定的机翼结构载荷误差大。由于运输类飞机具有结构柔性大、频率低的特点，货物空投、尤其是重装货物空投过程的激励载荷具有较宽的频带，会激励起飞机较多弹性模态的振动，而这种动响应可能会被气动弹性效应变得更大，进而引起飞机结构载荷的增大，这些对飞机的结构安全性和飞行安全性会引起不利的影响。第二、在弹性飞机货物空投动响应建模过程中，未考虑货物舱内移动过程、并未在空投过程中叠加离散阵风激励载荷，这些工作的不足既很难满足飞机设计规范的要求，也会造成飞机货物空投时机翼结构设计载荷偏小，进而对飞机的结构安全性和飞行安全性引起不利的影响。

发明内容

本发明的目的是：提供一种改进的运输类飞机货物空投时机翼结构载荷的确定方法，以便提高机翼结构载荷的确定精度，保证机翼结构的安全性和飞机的飞行安全性。

本发明的技术方案是：运输类飞机货物空投时机翼结构载荷的确定方法，其特征在于，确定机翼结构载荷的步骤如下：

1、飞机空投货物参数初始化：依据国军标GJB67.2A-2008中的空投阵风准则进行飞机空投货物参数初始化，需要初始化的空投货物参数为：输入飞机的有限元模型FEM、飞机的质量m_a、飞机飞行高度H、飞机飞行速度V、平均气动弦长c、货物的质量m_b、货物在飞机货舱内的初始装载位置x₀、货物的转动惯量I_by、货物在货舱内运动的加速度a_b、模态数目k、离散阵风强度u_w、离散阵风尺度L_w；

2、飞机有限元模态分析:以上述飞机的有限元模型FEM、模态数目k为控制参数，采用Nastran软件对飞机的有限元模型FEM进行模态分析，得到飞机的模态矩阵Φ、每阶模态的振动圆频率、广义质量矩阵M_qq和广义刚度矩阵K_qq，模态矩阵Φ包含飞机的刚体运动模态和弹性振动模态；

3、飞机非定常气动力系数矩阵计算：以上述模态矩阵Φ、飞机飞行高度H、飞机飞行速度V、平均气动弦长c、离散阵风强度u_w和离散阵风尺度L_w为输入变量进行时域非定常气动力系数矩阵计算，采用亚音速偶极子格网法进行计算得到结构振动引起的非定常气动力系数矩阵Q_a，采用阵风混合建模方法进行计算得到阵风引起的激励力系数矩阵Q_w；

4、飞机空投货物动响应计算：

4.1、重构弹性飞机货物空投动响应微分矩阵方程:在步骤1的12个参数、步骤2的模态矩阵Φ、每阶模态的振动圆频率、广义质量矩阵M_qq和广义刚度矩阵K_qq和步骤3的结构振动引起的非定常气动力系数矩阵Q_a及阵风引起的激励力系数矩阵Q_w的基础上，重构弹性飞机货物空投动响应微分矩阵方程：

$\begin{array}{c}(M_{\mathrm{qq}}+{\mathrm{\Φ}}_{\mathrm{bz}}^T{m}_b{\mathrm{\Φ}}_{\mathrm{bz}})\stackrel{\·\·}{\mathrm{\ξ}}+[{\left(\frac{{\∂\mathrm{\Φ}}_{\mathrm{bz}}}{\∂x}\right)}^T{m}_b{V}_{b2a}{\mathrm{\Φ}}_{\mathrm{bz}}+{\mathrm{\Φ}}_{\mathrm{bz}}^T{m}_b{V}_{b2a}\frac{\∂{\mathrm{\Φ}}_{\mathrm{bz}}}{\∂x}]\stackrel{\·}{\mathrm{\ξ}}+K_{\mathrm{qq}}\mathrm{\ξ}=\\ Q_a\mathrm{\ξ}+Q_w{u}_w-\left[\begin{array}{c}{m}_a\\ 0\end{array}\right]g-{\mathrm{\Φ}}_{\mathrm{bz}}^T{m}_{b}g\end{array}---\left(1\right)$

其中，ξ为对应于模态矩阵Φ的广义坐标位移，包含对应于飞机的刚体运动模态的广义坐标位移和弹性振动模态的广义坐标位移；

为对应于模态矩阵Φ的广义坐标速度，包含对应于飞机的刚体运动模态的广义坐标速度和弹性振动模态的广义坐标速度；

为对应于模态矩阵Φ的广义坐标加速度，包含对应于飞机的刚体运动模态的广义坐标加速度和弹性振动模态的广义坐标加速度；Φ_bz为货物所处位置处的飞机模态矩阵分量；V_b2a为货物在货舱内相对飞机的运动速度，通过对货物在货舱内运动的加速度a_b的一次积分可以得到；x为货物在货舱内的移动距离，通过对货物在货舱内运动的加速度a_b的两次积分可以得到；Q_aξ为飞机的结构振动引起的非定常气动力；Q_wu_w为阵风引起的激励力；g为重力加速度； $\left[\begin{array}{c}{m}_a\\ 0\end{array}\right]$ 为一个长度为k的列向量；

为利用Runge-Kutta法求解微分矩阵方程，将式(1)整理成一阶微分矩阵方程：

${\stackrel{\·}{x}}_{\mathrm{ae}}=A_{\mathrm{ae}}{x}_{\mathrm{ae}}+B_{\mathrm{aw}}{u}_w+B_{\mathrm{ag}}{m}_{b}g---\left(2\right)$

其中，A_ae为弹性飞机本体与非定常气动力组成系统的状态矩阵；B_aw为离散阵风对弹性飞机动响应的扰动矩阵；B_ag为由于货物移动对飞机产生作用力的扰动矩阵；u_w为离散阵风强度；x_ae为状态矢量：

$x_{\mathrm{ae}}=\left[\begin{array}{c}{x}_a\\ \mathrm{\ξ}\\ \stackrel{\·}{\mathrm{\ξ}}\end{array}\right]---\left(3\right)$

其中，x_a为非定常气动力的滞后根项目，模拟从飞机翼面上脱落的涡;

4.2、求解飞机1g平飞时的广义坐标位移ξ^Trim:飞机1g平飞时满足以下条件：

4.2.1、广义坐标加速度

和广义坐标速度

均为零向量；

4.2.2、货物相对飞机运动速度和位移均为0；

4.2.3、无离散阵风激励力；

根据上述条件，飞机1g平飞时的广义坐标位移ξ^Trim的计算如公式(4)所示：

${\mathrm{\ξ}}^{\mathrm{Trim}}=-{(K_{\mathrm{qq}}-Q_a)}^{-1}(\left[\begin{array}{c}{m}_a\\ 0\end{array}\right]+{\mathrm{\Φ}}_{\mathrm{bz}0}^T{m}_b)g---\left(4\right)$

其中，Φ_bz0为货物在货舱内初始装载位置处的飞机的模态矩阵分量；

4.3、求解飞机货物空投动响应微分矩阵方程:以飞机1g平飞时的广义坐标位移ξ^Trim为公式(2)的初始求解条件，采用变步长Runge-Kutta法求解飞机货物空投动响应微分矩阵方程，得到飞机的广义坐标位移ξ的动响应、广义坐标速度

的动响应和广义坐标加速度

的动响应:

$y_{\mathrm{ae}}=\left[\begin{array}{c}\mathrm{\ξ}\\ \stackrel{\·}{\mathrm{\ξ}}\\ \stackrel{\·\·}{\mathrm{\ξ}}\end{array}\right]=C_{\mathrm{ae}}{x}_{\mathrm{ae}}+D_{\mathrm{aw}}{u}_w+D_{\mathrm{ag}}{m}_{b}g---\left(5\right)$

其中，C_ae为弹性飞机本体与非定常气动力组成系统的输出矩阵；D_aw为离散阵风对弹性飞机动响应输出的传递矩阵；D_ag为由于货物移动对飞机动响应输出的传递矩阵；

5、计算机翼结构载荷动响应:根据步骤4的广义坐标位移ξ的动响应，采用节点刚度矩阵方法计算机翼结构载荷动响应F_wing(t)，计算公式如式(6)所示：

F_wing(t)=T_wingξ  (6)

其中，T_wing为从模态广义坐标位移到机翼根部单元梁的载荷变换矩阵；

6、确定空投时机翼结构的设计载荷：根据步骤5中得到的机翼结构载荷动响应F_wing(t)，采用包线筛选方法确定机翼结构的空投设计载荷。

本发明的优点是：提供了一种改进的运输类飞机货物空投时机翼结构载荷的确定方法，提高了机翼结构载荷的确定精度，保证了机翼结构的安全性和飞机的飞行安全性。

附图说明

图1是本发明一个实施例的机翼法向载荷的动响应图，横轴为时间，纵轴为机翼法向载荷增量，带圆圈符号的曲线表示空投过程中有阵风激励的机翼法向载荷增量动响应，带三角符号的曲线表示空投过程中无阵风激励的机翼法向载荷增量动响应，从图中可以看出空投过程中有阵风激励的机翼法向载荷增量比无阵风激励的机翼法向载荷增量大。

图2是本发明一个实施例的机翼弯矩增量的动响应图，横轴为时间，纵轴为机翼弯矩增量，带圆圈符号的曲线表示空投过程中有阵风激励的机翼弯矩增量动响应，带三角符号的曲线表示空投过程中无阵风激励的机翼弯矩增量动响应，从图中可以看出空投过程中有阵风激励的机翼弯矩增量比无阵风激励的机翼弯矩增量大。

图3是本发明一个实施例的机翼扭矩增量的动响应图，横轴为时间，纵轴为机翼扭矩增量，带圆圈符号的曲线表示空投过程中有阵风激励的机翼扭矩增量动响应，带三角符号的曲线表示空投过程中无阵风激励的机翼扭矩增量动响应，从图中可以看出空投过程中有阵风激励的机翼扭矩增量比无阵风激励的机翼扭矩增量大。

具体实施方式

下面结合附图对本发明作进一步详细说明。运输类飞机货物空投时机翼结构载荷的确定方法，其特征在于，确定机翼结构载荷的步骤如下：

1、飞机空投货物参数初始化：依据国军标GJB67.2A-2008中的空投阵风准则进行飞机空投货物参数初始化，需要初始化的空投货物参数为：输入飞机的有限元模型FEM、飞机的质量m_a、飞机飞行高度H、飞机飞行速度V、平均气动弦长c、货物的质量m_b、货物在飞机货舱内的初始装载位置x₀、货物的转动惯量I_by、货物在货舱内运动的加速度a_b、模态数目k、离散阵风强度u_w、离散阵风尺度L_w；

在上述12个参数中飞机的有限元模型FEM是本发明的重要的输入，一般飞机的有限元模型FEM为梁架模型，飞机的有限元模型FEM是对真实飞机通过静力缩聚减少自由度（DoF）获得的，其中相应的点集中质量分配到有限元节点上，这样可以大大减少飞机的分析自由度，提高运输类飞机货物空投时机翼结构载荷确定的速度；模态数目k的选取需要覆盖飞机的0～20Hz的模态，这样仿真计算得到的动响应结果会有足够的精度；其它的10个参数按照国军标GJB67.2A-2008中的空投阵风准则要求进行初始化；

2、飞机有限元模态分析:以上述飞机的有限元模型FEM、模态数目k为控制参数，采用Nastran软件对飞机有限元模型进行模态分析，得到飞机的模态矩阵Φ、每阶模态的振动圆频率、广义质量矩阵M_qq和广义刚度矩阵K_qq，模态矩阵Φ包含飞机的刚体运动模态和弹性振动模态；

本发明采用“稳定模态基底法”思想，对飞机有限元模型只进行一次模态分析，得到的模态矩阵Φ，在飞机进行货物空投动响应计算过程中，通过对模态矩阵Φ的线性组合更新方程中的项目；“稳定模态基底法”技术见文献：“A general approach to modal analysis for time-varyingsystem”,Browder A M,AIAA88-2356,1988;

飞机是一个在天空中自由飞行的飞行器，必须考虑其刚体运动模态；飞机受到外载荷作用还会发生变形运动，也需要考虑其弹性振动模态，基于模态叠加理论，可以通过对飞机的刚体运动模态和弹性振动模态线性组合得到飞机上的任何一点的位移和变形，这些位移和变形是飞机的非定常气动力计算的输入；FEM模态分析原理见文献：“Aeroelastic analysis user’s guide”,MSC.Software Corporation,MSC.Software Corporation,2004；

3、飞机非定常气动力系数矩阵计算：以上述模态矩阵Φ、飞机飞行高度H、飞机飞行速度V、平均气动弦长c、离散阵风强度u_w和离散阵风尺度L_w为输入变量进行时域非定常气动力系数矩阵的计算，采用亚音速偶极子格网法进行计算得到结构振动引起的非定常气动力系数矩阵Q_a，采用阵风混合建模方法进行计算得到阵风引起的激励力系数矩阵Q_w，Q_a和Q_w均通过ZAERO软件计算获得，它们的计算原理见文献：“ZAERO theoretical manual”,ZONATechnology,Inc,ZONA Technology,Inc,2008；

4、飞机空投货物动响应计算：

4.1、重构弹性飞机货物空投动响应微分矩阵方程:飞机进行货物空投过程中需要确定机翼结构的空投设计载荷，这就需要对飞机进行货物空投的整个过程的动响应方程进行重构，并进行动力学仿真；本发明提供了一种运输类飞机货物空投动响应方程的重构方法，如式(7)所示，在步骤1的12个参数、步骤2的模态矩阵Φ、每阶模态的振动圆频率、广义质量矩阵M_qq和广义刚度矩阵K_qq和步骤3的结构振动引起的非定常气动力系数矩阵Q_a及阵风引起的激励力系数矩阵Q_w的基础上，重构弹性飞机货物空投动响应微分矩阵方程：

$\begin{array}{c}(M_{\mathrm{qq}}+{\mathrm{\Φ}}_{\mathrm{bz}}^T{m}_b{\mathrm{\Φ}}_{\mathrm{bz}})\stackrel{\·\·}{\mathrm{\ξ}}+[{\left(\frac{{\∂\mathrm{\Φ}}_{\mathrm{bz}}}{\∂x}\right)}^T{m}_b{V}_{b2a}{\mathrm{\Φ}}_{\mathrm{bz}}+{\mathrm{\Φ}}_{\mathrm{bz}}^T{m}_b{V}_{b2a}\frac{\∂{\mathrm{\Φ}}_{\mathrm{bz}}}{\∂x}]\stackrel{\·}{\mathrm{\ξ}}+K_{\mathrm{qq}}\mathrm{\ξ}=\\ Q_a\mathrm{\ξ}+Q_w\mathrm{\ξ}-\left[\begin{array}{c}{m}_a\\ 0\end{array}\right]g-{\mathrm{\Φ}}_{\mathrm{bz}}^T{m}_{b}g\end{array}---\left(7\right)$

其中，ξ为对应于模态矩阵Φ的广义坐标位移，包含对应于飞机的刚体运动模态的广义坐标位移和弹性振动模态的广义坐标位移；

为对应于模态矩阵Φ的广义坐标速度，包含对应于飞机的刚体运动模态的广义坐标速度和弹性振动模态的广义坐标速度；

为对应于模态矩阵Φ的广义坐标加速度，包含对应于飞机的刚体运动模态的广义坐标加速度和弹性振动模态的广义坐标加速度；Φ_bz为货物所处位置处的模态分量，由于货物在飞机货舱内的移动，Φ_bz为随时间变化的变量，通过对模态矩阵Φ插值可以得到，货物离开飞机后Φ_bz＝0，模态插值方法见文献：“飞机气动弹性力学手册”，管德，航空工业出版社，1994；V_b2a为货物在货舱内相对飞机的运动速度，通过对货物在货舱内运动的加速度a_b的一次积分可以得到；x为货物在货舱内的移动距离，通过对货物在货舱内运动的加速度a_b的两次积分可以得到；Q_aξ为飞机的结构振动引起的非定常气动力；Q_wu_w为阵风引起的激励力；g为重力加速度； $\left[\begin{array}{c}{m}_a\\ 0\end{array}\right]$ 为一个长度为k的列向量；

为了利用Runge-Kutta法求解微分矩阵方程，需要将公式（7）整理为一阶微分矩阵方程：

${\stackrel{\·}{x}}_{\mathrm{ae}}=A_{\mathrm{ae}}{x}_{\mathrm{ae}}+B_{\mathrm{aw}}{u}_w+B_{\mathrm{ag}}{m}_{b}g---\left(8\right)$

其中，A_ae为弹性飞机本体与非定常气动力组成系统的状态矩阵；B_aw为离散阵风对弹性飞机动响应的扰动矩阵；B_ag为由于货物移动对飞机产生作用力的扰动矩阵；u_w为离散阵风强度；x_ae为状态矢量：

$x_{\mathrm{ae}}=\left[\begin{array}{c}{x}_a\\ \mathrm{\ξ}\\ \stackrel{\·}{\mathrm{\ξ}}\end{array}\right]---\left(9\right)$

其中，x_a为非定常气动力的滞后根项目，模拟从飞机翼面上脱落的涡，具体计算见文献：“ZAERO theoretical manual”,ZONA Technology,Inc,ZONATechnology,Inc,2008；

4.2、求解飞机1g平飞时的广义坐标位移ξ^Trim:飞机1g平飞时满足以下条件：

4.2.1、广义坐标加速度和广义坐标速度

均为零；

4.2.2、货物相对飞机运动速度和位移均为0；

4.2.3、无离散阵风扰动；

根据上述条件，飞机1g平飞时的广义坐标位移ξ^Trim的计算如公式(10)所示：

${\mathrm{\ξ}}^{\mathrm{Trim}}=-{(K_{\mathrm{qq}}-Q_a)}^{-1}(\left[\begin{array}{c}{m}_a\\ 0\end{array}\right]+{\mathrm{\Φ}}_{\mathrm{bz}0}^T{m}_b)g---\left(10\right)$

其中，Φ_bz0为货物在货舱内初始装载位置处的飞机的模态矩阵分量；

公式(10)为定解代数方程组，通过代数运算即可得到飞机1g平飞时的广义坐标位移ξ^Trim；

4.3、求解飞机货物空投动响应微分矩阵方程:飞机的任何一个部件载荷的确定均为通过广义坐标位移ξ的动响应、广义坐标速度的动响应和广义坐标加速度

的动响应与模态矩阵Φ的运算实现；机翼结构载荷的确定需要通过模态广义坐标位移ξ的动响应实现；以飞机1g平飞时的广义坐标位移ξ^Trim为公式(8)的初始求解条件，采用变步长Runge-Kutta法求解飞机货物空投动响应微分矩阵方程，得到广义坐标位移ξ的动响应、广义坐标速度的动响应和广义坐标加速度

的动响应：

$y_{\mathrm{ae}}=\left[\begin{array}{c}\mathrm{\ξ}\\ \stackrel{\·}{\mathrm{\ξ}}\\ \stackrel{\·\·}{\mathrm{\ξ}}\end{array}\right]=C_{\mathrm{ae}}{x}_{\mathrm{ae}}+D_{\mathrm{aw}}{u}_w+D_{\mathrm{ag}}{m}_{b}g---\left(11\right)$

其中，C_ae为弹性飞机本体与非定常气动力组成系统的输出矩阵；D_aw为离散阵风对弹性飞机动响应输出的传递矩阵；D_ag为由于货物移动对飞机动响应输出的传递矩阵；

5、计算机翼结构载荷动响应:根据步骤4的广义坐标位移ξ的动响应，采用节点刚度矩阵方法计算机翼结构载荷动响应F_wing(t)，计算公式如公式(12)所示：

F_wing(t)=T_wingξ  (12)

其中，T_wing为从模态广义坐标位移ξ到机翼结构载荷F_wing(t)的变换矩阵，计算原理见文献：“Aeroelastic analysis user’s guide”,MSC.SoftwareCorporation,MSC.Software Corporation,2004；

6、确定空投时机翼结构的设计载荷：步骤5中得到的机翼结构载荷动响应F_wing(t)并不是机翼结构的最终设计载荷，机翼结构载荷动响应F_wing(t)包括剪力、弯矩和扭矩分量，采用包线筛选方法，对F_wing(t)分量组合剪力-弯矩、剪力-扭矩和弯矩-扭矩组合进行载荷筛选，最终可以确定出飞机进行货物空投时的机翼结构的设计载荷，供机翼结构强度设计使用，载荷包线筛选方法见文献：“飞机非线性飞行载荷计算方法研究”，邓立东、李天、薛晓春，《航空学报》,2002，23卷(第4期)，P.317-320。

本发明的工作原理是：

针对当前运输类飞机货物空投时机翼结构载荷的确定方法的两个缺点：第一、所确定的机翼结构载荷误差大；第二、在弹性飞机货物空投动响应建模过程中，未考虑货物舱内移动过程、并未在空投过程中叠加离散阵风激励力。本发明采取以下措施提高运输类飞机货物空投时的机翼结构载荷的确定精度：

（1）以飞机的有限元模型FEM为主要基础，对运输类飞机货物空投过程进行动力学方程的重构；基于模态叠加理论，利用飞机的刚体运动模态和飞机的弹性振动模态通过线性组合表示飞机的运动和变形，这种方法提高了对飞机运动与变形描述的精度；

（2）在建立运输类飞机货物空投动响应微分矩阵方程的过程中，考虑了货物在飞机货舱内的移动过程，建立了能够更为真实反映货物运动过程的动响应方程；

（3）在货物空投过程中，对运输机叠加了符合国军标GJB67.2A-2008要求的离散阵风，即阵风速度剖面形状为1-cos形、当量阵风速度为7.6m/s，通过这样的计算，可以获得货物空投时机翼结构更为精准的设计载荷；

（4）飞机有限元模型FEM的模态求解和非定常气动力系数矩阵的计算均采用成熟、稳定、可靠的商业软件，计算结果精度高，为运输类飞机货物空投时机翼结构载荷的确定提供了可靠的输入数据。

实施例

采用本发明，对确定某型运输机进行货物空投时机翼结构载荷进行了仿真计算。

将某型运输机离散成有限元模型FEM，这样就可以将飞机模拟成为弹性飞机；货物与飞机的质量比为0.13；货物在飞机货舱内的初始装载位置为飞机的质心处；货物在货舱内的移动加速度a_b为2.5m/s²；飞机的飞行高度H为30m；飞机的飞行速度V为150m/s；平均气动弦长c为4.5m；选取飞机的前40阶模态；当量离散阵风强度u_w为7.6m/s；离散阵风尺度L_w为112.5m。以上述参数为初始化参数。

利用商业软件Nastran对飞机的有限元模型FEM进行模态分析，得到飞机的前40阶模态矩阵Φ，以及与模态矩阵Φ对应的每阶模态的振动圆频率、广义质量矩阵M_qq和广义刚度矩阵K_qq，模态矩阵Φ的前6阶为运输机的6个刚体运动模态，剩下的34阶模态为运输机的弹性振动模态。

将上述模态矩阵Φ、飞机飞行高度30m、飞机飞行速度150m/s、平均气动弦长4.5m、当量离散阵风强度7.6m/s和离散阵风尺度112.5m输入到商业软件ZAERO中，进行时域非定常气动力系数矩阵计算；在软件ZAERO中选择亚音速偶极子格网法对结构振动引起的非定常气动力系数矩阵Q_a进行计算；在软件ZAERO中选择阵风混合建模方法对阵风引起的激励力系数矩阵Q_w进行计算。

将模态矩阵Φ、广义质量矩阵M_qq、广义刚度矩阵K_qq、货物和飞机的质量、当量离散阵风强度7.6m/s、结构振动引起的非定常气动力系数矩阵Q_a和阵风引起的激励力系数矩阵Q_w代入弹性飞机货物空投动响应微分矩阵方程、即公式(8)中，完成对弹性飞机货物空投动响应微分矩阵方程的重构。

根据货物在飞机货舱内的初始装载位置—飞机质心处，对模态矩阵Φ进行插值得到货物在货舱内初始装载位置处的模态矩阵分量Φ_bz0，进而利用飞机1g平飞时的广义坐标位移ξ^Trim的计算公式(10)计算得到广义坐标位移ξ^Trim。

以飞机1g平飞时的广义坐标位移ξ^Trim为公式(8)的初始求解条件，采用变步长Runge-Kutta法求解飞机货物空投动响应微分矩阵方程，具体过程为：根据货物在货舱内的移动加速度2.5m/s²，积分计算货物当前时刻的运动速度V_b2a和相对货舱的位移x；对模态矩阵Φ进行模态插值得到货物当前时刻所处位置处的模态矩阵分量Φ_bz；基于当前时刻的V_b2a、x和Φ_bz，更新公式(8)中的项目；采用Runge-Kutta法迭代一步，即完成一步仿真计算。同时，利用公式(11)得到飞机的广义坐标位移ξ的动响应、广义坐标速度

的动响应和广义坐标加速度

的动响应。

根据上述计算得到的广义坐标位移ξ的动响应，利用公式(12)计算机翼结构载荷动响应F_wing(t)。

采用包线筛选方法对F_wing(t)分量组合剪力-弯矩、剪力-扭矩和弯矩-扭矩组合进行载荷筛选，最终可以确定出飞机进行货物空投时的机翼结构设计载荷，供机翼结构强度设计使用。

仿真结果如图1、图2和图3所示，图中对比了有/无阵风激励情况下的机翼载荷动响应结果。图1、图2和图3分别为机翼法向载荷增量、机翼弯矩增量和机翼扭矩增量的动响应图，横轴为时间，从图中可以看出空投过程中叠加阵风激励的机翼结构载荷增量比无阵风激励的机翼结构载荷增量大。

通过对实施例的仿真结果分析，表明本发明对运输类飞机货物空投过程时机翼结构载荷的确定结果规律良好、数值量级合适；在进行运输类飞机货物空投过程时机翼结构载荷的确定过程中，必须严格按照国军标GJB67.2A-2008的空投阵风准则，叠加符合规范要求的离散阵风，以达到机翼结构空投设计载荷的边界值。

Claims (1)

1.运输类飞机货物空投时机翼结构载荷的确定方法，其特征在于，确定机翼结构载荷的步骤如下：

1.1、飞机空投货物参数初始化：依据国军标GJB67.2A-2008中的空投阵风准则进行飞机空投货物参数初始化，需要初始化的空投货物参数为：输入飞机的有限元模型FEM、飞机的质量m_a、飞机飞行高度H、飞机飞行速度V、平均气动弦长c、货物的质量m_b、货物在飞机货舱内的初始装载位置x₀、货物的转动惯量I_by、货物在货舱内运动的加速度a_b、模态数目k、离散阵风强度u_w、离散阵风尺度L_w；

1.2、飞机有限元模态分析:以上述飞机的有限元模型FEM、模态数目k为控制参数，采用Nastran软件对飞机的有限元模型FEM进行模态分析，得到飞机的模态矩阵Φ、每阶模态的振动圆频率、广义质量矩阵M_qq和广义刚度矩阵K_qq，模态矩阵Φ包含飞机的刚体运动模态和弹性振动模态；

1.3、飞机非定常气动力系数矩阵计算：以上述模态矩阵Φ、飞机飞行高度H、飞机飞行速度V、平均气动弦长c、离散阵风强度u_w和离散阵风尺度L_w为输入变量进行时域非定常气动力系数矩阵计算，采用亚音速偶极子格网法进行计算得到结构振动引起的非定常气动力系数矩阵Q_a，采用阵风混合建模方法进行计算得到阵风引起的激励力系数矩阵Q_w；

1.4、飞机空投货物动响应计算：

1.4.1、重构弹性飞机货物空投动响应微分矩阵方程:在步骤1.1的12个参数、步骤1.2的模态矩阵Φ、每阶模态的振动圆频率、广义质量矩阵M_qq和广义刚度矩阵K_qq和步骤1.3的结构振动引起的非定常气动力系数矩阵Q_a及阵风引起的激励力系数矩阵Q_w的基础上，重构弹性飞机货物空投动响应微分矩阵方程：

$\begin{array}{c}(M_{\mathrm{qq}}+{\mathrm{\Φ}}_{\mathrm{bz}}^T{m}_b{\mathrm{\Φ}}_{\mathrm{bz}})\stackrel{\·\·}{\mathrm{\ξ}}+[{\left(\frac{{\∂\mathrm{\Φ}}_{\mathrm{bz}}}{\∂x}\right)}^T{m}_b{V}_{b2a}{\mathrm{\Φ}}_{\mathrm{bz}}+{\mathrm{\Φ}}_{\mathrm{bz}}^T{m}_b{V}_{b2a}\frac{\∂{\mathrm{\Φ}}_{\mathrm{bz}}}{\∂x}]\stackrel{\·}{\mathrm{\ξ}}+K_{\mathrm{qq}}\mathrm{\ξ}=\\ Q_a\mathrm{\ξ}+Q_w{u}_w-\left[\begin{array}{c}{m}_a\\ 0\end{array}\right]g-{\mathrm{\Φ}}_{\mathrm{bz}}^T{m}_{b}g\end{array}---\left(1\right)$

其中，ξ为对应于模态矩阵Φ的广义坐标位移，包含对应于飞机的刚体运动模态的广义坐标位移和弹性振动模态的广义坐标位移；

为对应于模态矩阵Φ的广义坐标速度，包含对应于飞机的刚体运动模态的广义坐标速度和弹性振动模态的广义坐标速度；为对应于模态矩阵Φ的广义坐标加速度，包含对应于飞机的刚体运动模态的广义坐标加速度和弹性振动模态的广义坐标加速度；Φ_bz为货物所处位置处的飞机模态矩阵分量；V_b2a为货物在货舱内相对飞机的运动速度，通过对货物在货舱内运动的加速度a_b的一次积分可以得到；x为货物在货舱内的移动距离，通过对货物在货舱内运动的加速度a_b的两次积分可以得到；Q_aξ为飞机的结构振动引起的非定常气动力；Q_wu_w为阵风引起的激励力；g为重力加速度； $\left[\begin{array}{c}{m}_a\\ 0\end{array}\right]$ 为一个长度为k的列向量；

为利用Runge-Kutta法求解微分矩阵方程，将式(1)整理成一阶微分矩阵方程：

${\stackrel{\·}{x}}_{\mathrm{ae}}=A_{\mathrm{ae}}{x}_{\mathrm{ae}}+B_{\mathrm{aw}}{u}_w+B_{\mathrm{ag}}{m}_{b}g---\left(2\right)$

其中，A_ae为弹性飞机本体与非定常气动力组成系统的状态矩阵；B_aw为离散阵风对弹性飞机动响应的扰动矩阵；B_ag为由于货物移动对飞机产生作用力的扰动矩阵；u_w为离散阵风强度；x_ae为状态矢量：

$x_{\mathrm{ae}}=\left[\begin{array}{c}{x}_a\\ \mathrm{\ξ}\\ \stackrel{\·}{\mathrm{\ξ}}\end{array}\right]---\left(3\right)$

其中，x_a为非定常气动力的滞后根项目，模拟从飞机翼面上脱落的涡;

1.4.2、求解飞机1g平飞时的广义坐标位移ξ^Trim:飞机1g平飞时满足以下条件：

1.4.2.1、广义坐标加速度和广义坐标速度

均为零向量；

1.4.2.2、货物相对飞机运动速度和位移均为0；

1.4.2.3、无离散阵风激励力；

根据上述条件，飞机1g平飞时的广义坐标位移ξ^Trim的计算如公式(4)所示：

${\mathrm{\ξ}}^{\mathrm{Trim}}=-{(K_{\mathrm{qq}}-Q_a)}^{-1}(\left[\begin{array}{c}{m}_a\\ 0\end{array}\right]+{\mathrm{\Φ}}_{\mathrm{bz}0}^T{m}_b)g---\left(4\right)$

其中，Φ_bz0为货物在货舱内初始装载位置处的飞机的模态矩阵分量；

1.4.3、求解飞机货物空投动响应微分矩阵方程:以飞机1g平飞时的广义坐标位移ξ^Trim为公式(2)的初始求解条件，采用变步长Runge-Kutta法求解飞机货物空投动响应微分矩阵方程，得到飞机的广义坐标位移ξ的动响应、广义坐标速度的动响应和广义坐标加速度

的动响应:

$y_{\mathrm{ae}}=\left[\begin{array}{c}\mathrm{\ξ}\\ \stackrel{\·}{\mathrm{\ξ}}\\ \stackrel{\·\·}{\mathrm{\ξ}}\end{array}\right]=C_{\mathrm{ae}}{x}_{\mathrm{ae}}+D_{\mathrm{aw}}{u}_w+D_{\mathrm{ag}}{m}_{b}g---\left(5\right)$

其中，C_ae为弹性飞机本体与非定常气动力组成系统的输出矩阵；D_aw为离散阵风对弹性飞机动响应输出的传递矩阵；D_ag为由于货物移动对飞机动响应输出的传递矩阵；

1.5、计算机翼结构载荷动响应:根据步骤1.4的广义坐标位移ξ的动响应，采用节点刚度矩阵方法计算机翼结构载荷动响应F_wing(t)，计算公式如式(6)所示：

F_wing(t)=T_wingξ  (6)

其中，T_wing为从模态广义坐标位移到机翼根部单元梁的载荷变换矩阵；

1.6、确定空投时机翼结构的设计载荷：根据步骤1.5中得到的机翼结构载荷动响应F_wing(t)，采用包线筛选方法确定机翼结构的空投设计载荷。

Images