0 引言

在飞机的概念和详细设计过程中,气动与结构模型两者之间的耦合作用一直是国内外近些年来考虑弹性效应的载荷设计中重点研究的问题之一 ^[1]。由于气动和结构有限元模型的表面网格单元疏密和大小等各有特点,因此必须找到适合两者网格之间载荷和变形等不同信息的传递方法,才能进行结构和气动之间的耦合计算和分析^[2]。在数据交换过程中要满足几个基本的条件:一是气动载荷要保证插值前后力和力矩的等效作用;二是结构变形须保证插值前后外形的光滑连续;三是尽量分部件进行,在保证整个耦合过程中较高的计算速度的同时防止中间数组畸变引起的插值错误^[3-4]。

气动和结构耦合研究一般包括气动载荷求解和结构变形分析以及相互发生的信息传递等^[5-6]。气动载荷和结构变形等数据插值的质量对整个耦合分析的准确和快慢有决定性的作用。因而考虑静气动弹性影响的载荷设计的主要工作之一就是确定良好的插值方法^[7]。一般来说,表面信息的传递主要有两种方法^[8]。在表面跟踪方法中,主要通过把气动网格和距离最近的结构结点进行关联,根据相应的坐标加权推算两者的传递函数进行计算。而表面拟合方法则是根据模型表面上选取可以表示外形特征的结点构造函数,从而得到相应区域气动网格上的变形量^[9-10]。

本文对飞机的初步和详细设计阶段考虑静气动弹性影响的载荷设计中的问题进行全面分析,同时分别对气动和结构模型之间参数传递和信息交互一般所采用的插值方法进行了原理研究和优劣对比,在此基础上发展和改进了曲面样条函数和形函数面积坐标加权两种插值方法进行了参数定义、公式推导和程序计算分析,分别针对弹性载荷设计中三维气动载荷传递情况和结构变形插值的计算,解决了整个考虑静气动弹性的载荷设计流程中最重要的环节,对于提高整个载荷设计的计算效率、准确性以及载荷计算结果的可视化有特别重要的价值。

1 理论基础

1.1 高阶面元法

飞行器初步设计阶段的弹性飞行载荷分析中使用线性气动力。线性气动力计算方法一般包括平板面元法和高阶面元法。其中高阶面元法由于基于三维真实外形,气动载荷特别是分布载荷的计算相对更加精确。其主要原理如下所示:

小扰动位势流方程在亚音速情况时可写成式(1)的形式^[3]

其中$M_{\mathrm{\infty }}^{}$是自由来流的马赫数,$\varphi$为扰动位函数。平板面元法和高阶面元法均基于小扰动势能理论发展而来。相应的气动模型分别为平面网格和三维曲面网格。前者分别在每个气动网格四分之一和四分之三中心弦线上设定压力点和控制点,而后者无须在相应的地方给出奇点,其压力点和控制点的排列都是根据大量各种翼型和来流工况下的数值模拟得来,如图1所示。

气动模型中,每个网格单元的点源和偶极子的强度分别设定两个不同的边界条件计算相应网格点源和偶极子强度:狄利克雷边界设定${\varphi }_L$=0和诺艾曼边界条件$\frac{\partial \varphi }{\partial n}=-V_{\mathrm{\infty }}\cdot n$,同时附加零力条件$\frac{\partial \varphi }{\partial x}=0$以模拟尾流情况。

1.2 曲面样条函数插值方法

在考虑气动弹性影响的载荷设计中,结构和气动模型之间变形的传递采用曲面样条函数方法进行插值。具体原理如下:

设函数W=W(X)在空间D上连续,那么对任何未知数X _k,相应的有一个确定的值W _k。由此得出,根据已知结构点坐标和相应位置的变形量构建特定的逼近函数,形式如式(2)所示。

式中:

$r_i^2=\sum _{p=1}^N{\left(x_p-x_{pi}\right)}^2$

c₁,c₂,…,c_N+1+n为待定系数;其具体值由下面的方程组来确定:

式中:

$r_{ji}^2=\sum _{p=1}^N{\left(x_{pj}-x_{pi}\right)}^2$

h_j表示第j 个结构点对应的加权系数,在插值前可先行给出。

1.3 形函数面积坐标加权方法

气动和结构模型之间另外一个问题是气动网格变形后载荷向结构模型的传递。本文中的气动载荷^[14]采用了改良版的形函数面积坐标加权方法,以此进行气动力的传递和插值。对气动模型任一网格单元,如图2所示,在积分得到其中心点的气动载荷后,计算比较得到其距离最近的三个结构点。首先获得气动和结构模型的主要网格点作为主要插值点,对每个气动单元中心点,计算所有结构结点和其之间的距离,通过程序获取最近的n个结点,进一步的,在这些点所能构成的所有三角形中,分析得到该气动点其上面的投影点,同时挑选出投影点在其区域内部的结构三角形,假设最终确定的数量为m个。每个气动网格中心点上载荷的分布就在对应的这些结构三角形上进行。通过计算得到这个气动点和这些结构三角形中的任一结点之间的距离,即:

其中,(x_a,y_a,z_a)为气动点(i, j)的坐标,($\mathrm{}{x}_i^k$,$\mathrm{}{y}_i^k$,$\mathrm{}{z}_i^k$)为第k个三角形对应的结点i三个方向的坐标值,通过程序求解能够得出相应的最大距离为:

在得到每个结构三角形中三个结点和该气动点之间的最大距离值后,比较得出其中最小值,从而确定对应的结构三角形。对任一气动网格中心点,相应的气动载荷按这个三角形对应的三个结点的面积坐标来进行分配。

对于结构结点n_i,将所有传递到其上面的所有气动载荷相加可写作:

式(8)中,i_max、 j_max为气动模型分块单元分别在x和y坐标方向的网格中心点序号的最大值,而n_max表示所选的主传力结构结点的序号最大值。

2 耦合数据的准备

2.1 模型的建立

本文中的气动模型为基于算例几何数模的三维表面网格,根据Zonair中的相关要求进行建立,机翼后缘和翼尖分别有cbar和crod单元模拟尾流,同时机翼后缘所在的平面与机身相交的部分须切开以模拟翼身干扰。而结构模型主要基于算例的刚度报告分部件进行建立。算例的气动和结构模型见图3。

2.2 数据的准备

根据已有的气动网格模型和结构有限元模型,编写高阶面元法的计算输入文件。调用Zonair软件,设定相应的计算工况进行分析。挑选典型工况(Ma=0.6,α=1°)导入tecplot,如图4,得到相应情况下的表面cp分布。

3 算例的插值结果

3.1 曲面样条插值方法

在考虑静气动弹性影响的载荷设计流程中,结构和气动模型之间数据的传递需要分部件进行,这样处理有两个好处:一方面分部件之后曲面样条函数的构建可以更加精确,避免不连续曲面构成的插值奇点出现;另外一方面结果的展示和对比也可以更加直观,方便进行插值效果的评估。编写分部件的气动和结构插值点信息文件并结合位移插值程序进行分析,得到不同部件的气动网格单元插值结果。根据算例的结构信息文件和相应的插值结果编写PLT文件并与气动模型在Zonair中计算得出的表面cp分布文件进行对比分析,从而进行位移插值程序效果的检验。算例模型不同部件的表面气动力系数分布插值结果分别如图5~图7所示。

根据以上不同部件表面气动力系数云图的比较可以看出,两种模型的分布差异很小,翼面的前后缘位置和翼身融合部位等典型位置的压力系数对比结果基本一致,从而验证了本文所发展的曲面样条函数方法和程序的有效性。另外,利用本文所开发的程序,数秒即可给出表面气动力系数传递后的全机分布结果和云图,大大提升了弹性载荷设计的精度和效率。

3.2 形函数面积坐标加权法

利用算例模型高阶面元法计算的全机表面气动力系数结果进行积分,得到相应的气动网格分布载荷。在此基础上利用本文所发展的气动载荷传递方法进行插值计算,得到对应的算例结构模型结点上的分布载荷,分别对全机载荷进行积分得到相应的力和力矩,气动和结构模型对应的总载荷对比结果见表1。

从表1结果可看出,气动和结构模型三个方向的总力和力矩基本一致,误差在一般弹性载荷设计载荷插值所允许的范围内,这样插值方法的合理性和程序的准确性就得到了检验。

4 结论

本文在综合分析多种数据传递方法的基本原理和方法的基础上,发展了两种分别适用于结构变形和气动载荷传递的插值方法,同时对其分别进行了公式推导与编程实现,并且利用已有的算例模型进行计算和对比。研究发现:

(1)气动和结构模型表面气动力系数和结构模型变形的传递中,利用曲面样条函数可以高效准确的完成插值计算;

(2)在弹性载荷设计中的气动载荷插值中,利用形函数面积坐标加权方法可以快速得到结构模型的全机分布载荷,并且插值前后和气动模型基本一致;

(3)结构变形和气动载荷插值是整个考虑气动弹性影响的载荷设计过程中的重要步骤,本文发展的两种方法很大的提升了弹性载荷设计的效率和计算精度。

参考文献