一种基于广义气动力的非定常气动力降阶方法与流程

本发明设计属于飞行器气动弹性试验技术领域，具体为一种基于广义气动力的非定常气动力降阶方法。

背景技术：

颤振是结构在气动力、弹性力和惯性力耦合作用下，发生的一种振幅不衰减的自激振动。目前，研究飞行器结构颤振问题的主要途径为理论计算和颤振试验。

理论计算虽然节省时间和经费，但是建立分析对象的数学模型时，在结构、气动力等方面需要引入一些假设，这些假设可能与真实情况有较大的偏差。

颤振试验分为风洞颤振试验和飞行颤振试验：风洞颤振试验直接考虑气动力的作用，但要对试验模型进行缩比设计，难以完全模拟实物，费用高且有一定风险；飞行颤振试验可以完全考虑试验对象的真实工作环境，但受试验条件的限制，难以获得准确的颤振边界，且费用昂贵、风险大。

地面颤振模拟试验则是一种新兴的颤振试验研究方法，它采用真实的飞行器结构作为试验对象，通过激振器产生的激振力来模拟结构在给定速度下受到的非定常气动力，在地面(风洞外)可以获得真实结构的气动弹性特性。在用激振器进行气动力模拟加载时，由于实际加载到结构上的气动力是通过数量有限的激振器加载的，所以在建立用于地面颤振模拟试验的气动力计算模块时，需要进行两个重要的步骤：

1)将气动面上分布的气动力进行降阶处理，从而利用有限点上的激励力来等效模拟分布的气动力；

2)将计算得到的频域气动力拟合到时域。

从公开的文献中可以看出，地面颤振模拟试验中气动力降阶方法的最新研究进展是：

zeng等在文献(zengj,kingsburyd,ritze,etal.gvt-basedgroundfluttertestwithoutwindtunnel[c]//52ndaiaa/asme/asce/ahs/ascstructures,structuraldynamicsandmaterialsconference19thaiaa/asme/ahsadaptivestructuresconference13t.2011:1942.)中通过面样条插值对气动力进行降阶，在试验气动力近似上，提出了以气动力降阶前后，结构的颤振速度和颤振频率的误差最小为目标的激振点和拾振点位置优化方法，得到的气动力降阶模型的精度很好，但是使用气动力降阶后模型的颤振特性与降阶前一致作为激振点和拾振点的选取目标，与地面颤振模拟试验的初衷，即探索模型的颤振边界，有所矛盾，而且先将气动力模型拟合到时域再降阶，拟合过程会耗费大量时间。

胡巍等在文献(胡巍,杨智春,谷迎松.带操纵面机翼气动弹性地面试验仿真系统中的气动力降阶方法[j].西北工业大学学报,2013,31(5):810-815.)中对带操纵面的机翼气动力降阶问题，提出将气动力模型先进行频域降阶、再时域拟合的方法，分别确定主翼面、操纵面上激振点和拾振点的位置，从而获得气动力降阶模型。先对频域气动力影响系数矩阵a(ω)进行降阶，可以降低矩阵的维数，在此基础上再进行时域拟合，可减少拟合过程耗费的时间，提高计算效率。

吴志刚等在文献(许云涛,吴志刚,杨超.地面颤振模拟试验中的非定常气动力模拟[j].航空学报,2012,33(11):1947-1957.)和文献(wuz,mac,yangc.newapproachtothegroundfluttersimulationtest[j].journalofaircraft,2016,53(5):1578-1580.)中提出在气动力降阶时，以优化颤振关键振型为目标函数对进行激振点拾振点位置优化，但是使用数量较少的激振点和拾振点组合得到的频域降阶气动力模型对应的广义气动力矩阵与降阶前原始广义气动力矩阵可能相差较大，会导致降阶前后结构的颤振边界区别较大。

已有气动力降阶方法都是基于频域气动力计算理论如偶极子格网法(dlm)，zona6，zona7等，在亚音速和超音速范围内具有很高的精度和计算效率，但是对于跨音速颤振问题，由于没有相应的跨音速频域气动力计算方法，上述气动力降阶方法将不再适用。随着计算流体力学技术的快速发展，cfd/csd耦合求解方法成为跨音速颤振求解的一种有效手段，可以考虑跨音速范围内流动的非线性特别是跨音速激波的影响，因此具有较高的准确度。但这种方法计算量大，耗时很长，限制了它在地面颤振模拟试验中的使用。文献(贺顺,杨智春,谷迎松.机翼跨音速颤振特性的频域分析.中国科学:物理学力学天文学,2014,44:285–292。)提出了一种机翼跨音速颤振频域计算方法，首先使机翼按照其某阶模态做微幅简谐运动，采用cfd方法计算出机翼的时域跨音速气动力，再通过模态变换获得时域模态气动力系数，进而计算出给定马赫数下的跨音速广义气动力矩阵，从而可应用频域颤振分析方法求解出机翼的颤振速度及颤振频率；zaero软件中的ztaic方法即等价片条法，输入翼型剖面的定常压力分布，然后使用等价片条理论计算跨声速小扰动方程，得到非定常气动力系数，再用该系数修正含激波效应的非定常压力，反复迭代后也得到考虑跨声速激波效应的广义气动力影响系数矩阵。

技术实现要素：

在地面颤振模拟试验中，通过激振器输出的激振力来模拟分布式的气动力。理论上讲，为了更精确地描述气动面上分布式的气动力，我们都希望使用更多的激振器，但如果激振器数目过多，激振力控制系统的设计难度将大大增加，因此在实际颤振试验中，希望在满足试验精度的前提下，尽量减少气动力的阶数(这样，对应的激振器的数目就会相应减少)。

本发明的目的就是对结构的非定常气动力模型进行降阶，提出一种基于广义气动力的气动力降阶方法，内容包括：坐标变换，插值气动力降阶和气动力时域拟合。首先通过坐标变换将模态坐标下的广义气动力变成物理坐标下有限元模型节点上分布的非定常气动力，然后通过曲面样条插值将分布的非定常气动力进行降阶，等效到有限个加载点处从而获得频域气动力降阶模型，最后使用最小状态法将频域气动力降阶模型拟合到时域。

本发明主要用于地面颤振模拟试验及其仿真中，实现气动力实时计算，具体是通过坐标变换将模态坐标下的广义气动力变成物理坐标下有限元模型节点上分布的非定常气动力，并通过曲面样条插值将分布的非定常气动力进行降阶，等效到有限个加载点处，然后使用最小状态法将降阶气动力模型转换到时域。本发明的优点在于：1)借助cfd跨声速非定常气动力方法，使该降阶方法适用于跨音速颤振分析；2)在尽量减少降阶气动力模型阶数的基础上，提高降阶气动力模型的精度，从而降低激振力控制系统设计的难度。

本发明的技术方案为：

所述一种基于广义气动力的非定常气动力降阶方法，其特征在于：包括以下步骤：

步骤1：针对需要进行地面颤振模拟试验的机翼，建立机翼的有限元模型，进行模态分析，得到机翼有限元模型的质量矩阵m以及机翼的模态振型矩阵φ，并在计算流体力学软件中计算该机翼在给定马赫数下的广义气动力矩阵；

步骤2：坐标变换：

在得到广义气动力矩阵后，根据以下公式

qaa＝mφqhh·φ^tm

得到物理坐标下的气动力影响系数矩阵qaa；其中qhh为步骤1得到的机翼在给定马赫数下的广义气动力矩阵；

步骤3：面样条插值气动力降阶：

根据机翼有限元模型上设定的激振点和拾振点数目，对机翼有限元模型上的激振点及拾振点的位置进行优化，使通过激振点和拾振点表示的气动节点插值振型与气动节点原始振型之间实现最优逼近；

得到机翼有限元模型上的激振点及拾振点位置后，采用插值方法实现从拾振点的位移得到全部结构节点位移的插值变换以及从全部结构节点的气动力到激振点作用力的插值变换；其中从拾振点的位移得到全部结构节点位移的插值变换关系为

x＝[gs]{xs}ns×1

其中x为全部结构节点位移，xs为拾振点的位移，ns为拾振点数目，gs为位移插值矩阵；从全部结构节点的气动力到激振点作用力的插值变换关系为：

{fs}na×1＝[gf]{f}

f为全部结构节点上的气动力，fs为激振点上的作用力，na为激振点数目，gf为力插值矩阵；进而得到降阶后的气动力影响系数矩阵为

[qs]na×ns＝[gf][qaa][gs]

步骤4：将气动力拟合到时域：

采用最小状态法将降阶后的气动力影响系数矩阵qs转换到时域，得到降阶的时域气动力：

其中转换到时域的气动力影响系数矩阵为：

式中，s是拉普拉斯变量，b是机翼的半弦长，v是来流速度，i为单位阵，a0，a1，a2，d，r和e是通过最小状态法求得的系数矩阵；

根据转换到时域的气动力影响系数矩阵，得到降阶的时域气动力为

其中q∞为动压。

进一步的优选方案，所述一种基于广义气动力的非定常气动力降阶方法，其特征在于：步骤3中采用平面薄板样条插值方法实现从拾振点的位移得到全部结构节点位移的插值变换以及从全部结构节点的气动力到激振点作用力的插值变换。

有益效果

本发明提出的基于广义气动力的气动力降阶方法，优点在于在计算广义气动力时，结构有限元模型与气动力模型之间的耦合过程使用了有限元模型的全部节点运动信息，以此广义气动力矩阵为输入，通过坐标变换和样条插值方法实现了气动力降阶过程，避免了已有方法先在物理坐标下对气动力影响系数矩阵降阶导致的降阶气动力模型的广义气动力与降阶前区别较大而引起的系统颤振边界差别较大的问题，而且使得在尽量减少降阶气动力模型阶数的基础上，提高降阶了气动力降阶模型的精度，从而降低了地面颤振模拟试验中激振力控制系统设计的难度。对于地面颤振模拟试验中跨声速气动力降阶问题，通过cfd跨声速气动力计算方法或者zaero的ztaic方法计算得到广义气动力矩阵后，使用本发明提出的气动力降阶方法可以得到用于地面颤振模拟试验的跨音速气动力降阶模型。

本发明的附加方面和优点将在下面的描述中部分给出，部分将从下面的描述中变得明显，或通过本发明的实践了解到。

附图说明

本发明的上述和/或附加的方面和优点从结合下面附图对实施例的描述中将变得明显和容易理解，其中：

图1舵面模型示意图；

图2舵面模型有限元网格；

图34激振点/4拾振点布置图；

图4降阶前后广义气动力矩阵对比图。

具体实施方式

下面详细描述本发明的实施例，所述实施例是示例性的，旨在用于解释本发明，而不能理解为对本发明的限制。

本发明涉及一种基于广义气动力的气动力降阶方法，内容包括：坐标变换，插值气动力降阶和气动力时域拟合。首先通过坐标变换将模态坐标下的广义气动力变成物理坐标下有限元模型节点上分布的非定常气动力，然后通过曲面样条插值将分布的非定常气动力进行降阶，等效到有限个加载点处从而获得频域气动力降阶模型，最后使用最小状态法将频域气动力降阶模型拟合到时域。

步骤1：针对需要进行地面颤振模拟试验的机翼，建立机翼的有限元模型，进行模态分析，得到机翼有限元模型的质量矩阵m以及机翼的模态振型矩阵φ，并在计算流体力学软件中使用文献(贺顺,杨智春,谷迎松.机翼跨音速颤振特性的频域分析.中国科学:物理学力学天文学,2014,44:285–292。)计算该机翼在给定马赫数下的广义气动力矩阵。

步骤2：在得到广义气动力矩阵后，推导物理坐标下的气动力影响系数矩阵：

f＝q∞·qhh·q(1)

其中qhh为广义气动力矩阵，q∞为动压，q为模态坐标，上式以结构的模态位移为输入，模态气动力向量为输出。

当结构的模态振型矩阵φ按质量矩阵m归一化后

φ^tmφ＝i(2)

物理位移和模态位移有如下关系

x＝φq(3)

在(3)式两端左乘φ^tm，可得物理位移和模态位移的关系

q＝φ^tmx(4)

将物理力变换为模态力

f＝φ^tf(5)

进行模态截断后，振型φ不是方阵，φ^t不能直接求逆，由于φ^t为行满秩，求得的广义逆为φ^t的右逆，显然我们需要的是φ^t的左逆。

将(5)式写成矩阵形式

因为振型矩阵按质量矩阵m归一化，由(2)式，可得

如果令将其带入(6)式

所以由模态力得到物理力的变换过程如下

将(4)式和(9)式带入(1)式，可得

f＝q∞·mφqhh·φ^tmx＝q∞·qaax(10)

式中：qaa＝mφqhh·φ^tm为物理坐标下有限元模型全部节点(a-set集)的气动力影响系数矩阵。

即在得到广义气动力矩阵后，可以根据公式

qaa＝mφ·qhh·φ^tm

得到物理坐标下的气动力影响系数矩阵qaa；其中qhh为步骤1得到的机翼在给定马赫数下的广义气动力矩阵。

步骤3：面样条插值气动力降阶：

在确定激振点和拾振点的数目后，选择不同的激振点/拾振点位置时，得到的系统气弹稳定性边界差别很大，且同降阶前的系统气弹稳定性边界相比也有较大的差别，所以需要对激振点/传感点的位置进行优化。

文献(许云涛,吴志刚,杨超.地面颤振模拟试验中的非定常气动力模拟[j].航空学报,2012,33(11):1947-1957.)中提出，降阶气动力与真实气动力的等效逼近，可以转化为使用激振点/拾振点表示的气动节点插值振型与气动节点原始振型之间的最优逼近，定义激振点/拾振点位置优化的目标函数为

式中ηj为j阶模态对颤振的贡献量，为气动节点原始振型φa的元素，为使用激振点/拾振点插值得到的气动节点新振型φ'a的元素，nm为模态截断阶数。obj越小，表明气动力模拟精度越高。

通过对机翼有限元模型上的激振点及拾振点的位置进行优化，使通过激振点和拾振点表示的气动节点插值振型与气动节点原始振型之间实现最优逼近。

得到机翼有限元模型上的激振点及拾振点位置后，通过样条插值方法，将气动力影响系数矩阵qaa插值到少数的激振点和拾振点后处完成气动力降阶：

这里采用平面薄板样条插值方法(thin-platesmoothingspline)来实现上述插值，实现实现从拾振点的位移得到全部结构节点位移的插值变换以及从全部结构节点的气动力到激振点作用力的插值变换。

其中从拾振点的位移得到全部结构节点位移的插值变换关系为

x＝[gs]{xs}ns×1(12)

其中x为全部结构节点位移，xs为拾振点的位移，ns为拾振点数目，gs为位移插值矩阵；从全部结构节点的气动力到激振点作用力的插值变换关系为：

{fs}na×1＝[gf]{f}(13)

f为全部结构节点上的气动力，fs为激振点上的作用力，na为激振点数目，gf为力插值矩阵；进而得到降阶后的气动力影响系数矩阵为

[qs]na×ns＝[gf][qaa][gs](14)

步骤4：将气动力拟合到时域：

采用最小状态法将降阶后的气动力影响系数矩阵qs转换到时域，得到降阶的时域气动力：

其中转换到时域的气动力影响系数矩阵为：

式中，s是拉普拉斯变量，b是机翼的半弦长，v是来流速度，i为单位阵，a0，a1，a2，d，r和e是通过最小状态法求得的系数矩阵；

根据转换到时域的气动力影响系数矩阵，得到降阶的时域气动力为

其中q∞为动压。

下面对图1所示的舵面模型进行气动力降阶处理，气弹分析基于nastran的sol145进行，舵面模型采用根部固支的约束方式，材料为铝，e＝70gpa，ν＝0.3，ρ＝2750kg/m³，对舵面划分结构网格如图2所示。

设定气弹分析的条件为：空气密度ρ空＝1.226kg/m³,马赫数ma＝1.2314，减缩频率k＝1e-7、0.05、0.1、0.15、0.3、0.6、1.0和2.0。半弦长b＝0.316m。气动力采用zona51，计算各减缩频率下的广义气动力矩阵后，计算物理坐标下对应的气动力影响系数矩阵。

使用面样条插值方法，对第一步得到的气动力影响系数矩阵进行降阶，在舵面上各布置4激振点/4拾振点(激振点和拾振点采用对位布置)，通过遗传算法优化激振点/拾振点位置，得到的激振点/拾振点位置如图3所示。

分别将降阶前后气动力矩阵各元素对减缩频率k值绘制散点图，可以定性验证降阶气动力模型的精度。降阶前后，广义气动力矩阵对比图如图4所示，图中各子图的标题表示对应的元素在广义气动力矩阵中的位置，横坐标代表实部，纵坐标代表虚部。其中o代表选择全部结构节点插值得到的广义气动力影响系数矩阵，即降阶前的qhh，x代表使用专利中提出的降阶气动力方法得到的广义气动力影响系数矩阵，☆为使用文献(许云涛,吴志刚,杨超.地面颤振模拟试验中的非定常气动力模拟[j].航空学报,2012,33(11):1947-1957.)气动力降阶方法得到的广义气动力影响系数矩阵。可以看出，对于前两阶广义气动力，本专利中提出的方法得到的降阶气动力模型的广义气动力矩阵更接近原始广义气动力矩阵，这是因为，与方法文献(许云涛,吴志刚,杨超.地面颤振模拟试验中的非定常气动力模拟[j].航空学报,2012,33(11):1947-1957.)中方法相比，新方法在气动力降阶前计算广义气动力矩阵时使用全部结构节点与气动节点耦合插值得到广义气动力矩阵，即使用原始广义气动力矩阵qhh进行气动力降阶，而使用文献(许云涛,吴志刚,杨超.地面颤振模拟试验中的非定常气动力模拟[j].航空学报,2012,33(11):1947-1957.)中方法得到的qhh为使用少数的插值节点(4个激振点/拾振点)与气动节点耦合得到的。

最后将降阶气动力模型拟合到时域，计算模型的颤振特性。降阶前后的颤振速度比较如表1所示：

表1舵面的颤振特性

由表1的结果知，采用本专利中提出的方法，气动力降阶前后舵面的颤振速度误差更小，说明了在保证试验精度的前提下，有效减少激振点数目，满足地面颤振模拟试验的要求。

尽管上面已经示出和描述了本发明的实施例，可以理解的是，上述实施例是示例性的，不能理解为对本发明的限制，本领域的普通技术人员在不脱离本发明的原理和宗旨的情况下在本发明的范围内可以对上述实施例进行变化、修改、替换和变型。