太阳能无人机机翼颤振动力学建模与分析

2023-09-27肖良华

计算机测量与控制 2023年9期

冒森，张斌，肖良华，陈斌，王珏

(成都飞机工业(集团)有限责任公司，成都 610073)

0 引言

太阳能飞机是采用太阳能为能源的飞行器。相对比常规布局飞机，太阳能飞机超长的航时和超高的飞行高度使其具有无与伦比的优势，而且契合绿色环保理念，因此得到国内外研究人员的广泛关注。此类飞机的设计，以能源设计为中心，采用超大展弦比的机翼及极低的翼载设计(翼载较常规动力飞机的低2个量级左右)，例如“太阳神”太阳能无人机结构面密度为 3.2 kg/m2[1]，而现阶段最长航时常规动力无人机“全球鹰” 的结构面密度高达 53 kg/m2[2]，因此太阳能飞机工程实现难度非常高。从全机重量分配来看，结构重量和蓄电池重量分别占全机重量的35%以上，且蓄电池重量随总重的增加而增大，因此结构的轻量化是实现全机减重的有效手段。

轻量化结构和大展弦比构型使得太阳能无人机颤振、突风响应等气动弹性问题非常显著。其中颤振是最重要的一类结构动强度问题。颤振现象是是机体结构在受到弹性力、惯性力、气动力共同作用下可能出现的一种自激振动，是一种飞机的固有属性[3]。太阳能无人机因颤振而发生的事故最早被是2003 年 6 月，美国“太阳神”无人机经过多次成功飞行后，在夏威夷附近海域进行低空飞行时，机体突然发生向上弯曲，进而出现俯仰振荡发散，最终导致坠毁[4]。距现在较近的一次太阳能飞机因颤振而发生事故的是上海奥科赛飞机公司设计的太阳能飞机2号机。此飞机的翼展为15米，在2017年7月的一次飞行过程中，在飞行高度为204.2 m、飞行速度为92 km/h时发生了颤振，最终也导致了飞机坠毁[5]。常规飞行器设计中对于结构颤振问题，一般采取试验与仿真相结合的方式进行解决，在风洞中可直接采用1：1的结构模型来进行试验。而太阳能无人机的翼展远大于现有的风洞尺寸，无法直接开展颤振风洞实验，需要设计并加工缩比模型，而且缩比模型的重量缩比量是尺寸缩比的三次方[6]，由于太阳能飞机的面密度极低，导致缩比模型加工难度很高，如此苛刻的条件让此类飞机基本不能开展全模的颤振试验[7]。因此在方案设计阶段，通过仿真计算的手段准确得到太阳能无人机的颤振速度并且给出防颤振设计规律具有重要意义。

目前在大展弦比飞机的颤振分析方面，国内外学者已经开展了广泛而深入的工作，太阳能无人机在结构布局上可归属于大展弦比常规布局飞机，现有的针对大展弦比飞机的气动弹性分析方法也适用于太阳能飞机。文献[8]和[9]首次使用动力学线化理论，将大变形下的机翼非线性有限元分析进行线性化处理，针对大展弦比机翼，分析了机翼的非线性静气弹和动气弹问题。文献[10]中在计算大展弦比机翼的气动载荷时，考虑曲面效应，应用了推广的三维升力线理论，并且在结构/气动耦合计算中首次使用曲面样条插值理论，从而计算出气动力计算所必需的网格节点位移，在此基础上对一个简单的鱼刺模型进行了相应计算，计算结果证明了该方法能很好的解决大展弦比机翼静气动弹性问题。文献[11]中则是国内第一个针对太阳能无人机开展气动弹性分析。首先建立太阳能无人机的结构有限元模型和气动计算模型，采用MSC.Nastran软件计算了静气动弹性变形和颤振速度。通过对计算结果的分析，总结了弹性变形程度对机翼上升力分布的影响，得到了不同飞行条件下的颤振速度和频率。文献[12]针对国内翼展为15米的太阳能飞机，进行了该太阳能飞机的动力学反向建模，并以所建模型为基础，进行颤振计算分析，最终建立了一套完整的关于太阳能飞机颤振分析的设计方法。文献[13]则在理论层面上更深层次地推导了太阳能飞机中机翼在受力变形后的切线刚度矩阵和质量矩阵，在不同飞行工况静气弹平衡位置使用准模态假设，采用P-K法计算了该飞机的在不同飞行载荷下的非线性颤振问题。文献[14]创建了一套使用大变形插值理论计算考虑几何非线性效应时的气动弹性分析框架，并以某太阳能无人机机翼缩比模型为对象使用此方法进行仿真分析，之后又开展地面模态试验与风洞试验，结果表明使用该方法进行仿真计算不仅可以显著提高计算效率，而且能够有效降低常规方法中载荷选择时产生的计算偏差。

可以看出目前针对太阳能无人机的颤振分析研究很多也非常深入，但是上述大多数的研究是根据确定性的结构参数进行正向结构建模和分析，且采用的模型都是简化的鱼刺模型。实际上在太阳能飞机的方案设计阶段，是无法开展地面GVT(ground vibration test)试验获得机翼的准确的刚度数据和重量数据，难以建立准确的鱼刺模型来进行颤振分析。因此在设计阶段，现有的方法不能很好的解决太阳能无人机的气动弹性建模和颤振分析的问题。

此外，太阳能无人机区别于传统飞机的一大特点在于蒙皮采用薄膜结构形式，体现在结构有限元中即为将蒙皮采用刚度相对较小的板单元[15-16]进行模拟，这会导致初始有限元模型在模态分析时出现大量蒙皮单元局部模态，在计算较多阶模态的情况下仍无法获取足够多的整体模态信息。此外，全机复杂有限元模型的自由度一般较多(通常为数十万自由度)，且包含大量用于部件间或集中质量点连接的MPC单元，因此有必要建立一种能够准确反映初始有限元模型动力学特性的工程化建模方法。

另一方面，由于具有明显弹性轴的太阳能机翼通常呈现出“长且柔”的悬臂梁特点，在机翼上布置吊舱可以起到气动力卸载的作用，但因机翼本身已具有较低的固有频率，系挂吊舱可能会进一步降低机翼的固有频率，使得机翼颤振速度发生大的改变[17]。所以研究太阳能机翼上吊舱不同位置的分布对于颤振特性影响也是稳定性分析中不可缺少的一环的。

针对太阳能飞机的独特之处，本文提出一种面向方案阶段的太阳能无人机机翼颤振分析办法，以结构静力模型为基础，通过工程化处理办法，开展太阳能无人机机翼气动弹性建模和颤振研究，在确保消除局部模态的同时不改变全机模态分布特性的基础上，计算不同机翼吊舱布置下机翼颤振速度，并以此为依据，总结提出太阳能飞机防颤振设计规律。

1 动力学建模

1.1 结构设计

该太阳能无人机机翼结构如图1所示。采用常规布局的单梁式机翼。主梁为变截面梁，翼肋是非主承力件。除连接起落架吊舱的肋为金属件外，其余翼肋均为一体成型复材件。能源/动力/起落架吊舱位于单边机翼展向10 m处。吊舱用于连接推进系统、起落架和安装部分储能电池。

图1 机翼结构设计图

蒙皮采用抗撕裂薄膜。结构材料以复合材料为主，局部使用金属材料，金属材料选择常规的铝合金。机翼主梁和接头选用中温固化环氧树脂单向碳纤维预浸料，超轻质肋和泡沫夹芯结构的复材选用中温固化环氧树脂单向碳纤维预浸料。

1.2 有限元建模

对于该机翼使用商业软件MSC.Patran开展有限元建模，有限元模型中按照实际情况建出，如图2所示。

图2 太阳能无人机机翼有限元模型

机翼肋缘条和腹板均按照实际数模建出，如图3所示。

图3 梁和肋有限元建模

太阳能电池阵铺设在机翼上蒙皮，为模拟真实情况，有限元模型中在机翼上蒙皮单元每一个节点建立质量元，平均分配太阳能电池阵的重量，如图4所示。

图4 太阳能薄膜集中质量单元布置

1.3 建模分析

有限元建模完成后，使用MSC.Nastran的101模块进行静力学分析，检验机翼强度。平飞状态的机翼气动力采用面压的方式加载在机翼蒙皮上，气动力为CFD由计算得到，均以外部载荷的形式加载到结构。机翼结构变形如图5所示。从静力学角度可以看出，蒙皮出现了失稳。

图5 机翼受气动载荷变形图

运用MSC.Nastran 103模块开展有限元模态计算，检验建模准确性。计算结果发现在从4阶模态之后，频率2.0 Hz到2.5 Hz之间出现大量的局部模态，如图6所示。这些局部模态，属于蒙皮失稳模态。局部模态极大的影响颤振计算结果，造成颤振计算结果不准确。

图6 蒙皮局部模态

需要注意的是，在开始进行模态计算中发现扭转频率偏低，排查原因是扭转刚度存在问题，在有限元模型中对于主梁采用杆单元进行模拟，肋采用板单元进行模拟，在梁与肋之间采用单点连接，但由于主梁几何尺寸比较大，仅简化成杆单元和板单元之间采用单点连接的形式会造成模型扭转刚度低于实际结构。

2 气动弹性建模修正

2.1 工程化处理方法

为了体现梁肋连接的抗扭能力，在肋腹板处挖出圆孔，将腹板圆孔周围节点与机翼主梁节点采用MPC固联，如图7所示。

图7 梁和肋有限元建模优化

太阳能机翼的薄膜蒙皮只承受拉力，不能承受压力，而且在安装时是进行预紧安装，蒙皮失稳只存在于仿真计算中，实际情况中不存在蒙皮失稳，因此在颤振计算中需要排除蒙皮局部模态的影响。现阶段国内外用于颤振分析的商业软件有ZAERO和MSC.Nastran，采用ZAERO计算颤振需要MSC.Nastran计算的模态信息，而MSC.Nastran无法针对这种情况直接进行颤振计算。另一种处理局部模态的办法是通过MPC约束蒙皮自由度，但是也会增加机翼刚度，造成颤振计算不准确。

本文提出从工程处理角度，对失稳处的蒙皮，采用增加壳单元面外刚度的方法消除局部模态。在薄壳单元有限元计算中，单元面内刚度矩阵Lm和面外刚度矩阵Lb如下[18]：

(1)

(2)

薄板类壳单元的局部模态通常表现为单个单元的弯曲或鼓包，如果增加壳单元的面外刚度，相当于提高了壳单元抵抗弯曲或鼓包变形的能力，这样即可达到消除局部模态的目的。同时，对于全机而言，带有蒙皮单元的机翼或尾翼部件整体模态(弯曲、扭转等)主要受蒙皮单元面内刚度影响，所以上述增大面外刚度的修改方式不会对全机整体模态产生太大影响。Nastran中可通过修改壳单元的PSHELL属性，增大弯曲刚度的比例，消除局部模态[19-20]。对于常规飞机的复材蒙皮或者金属蒙皮，失稳的处理办法一般情况是也是提高蒙皮刚度进行处理，因此具有一定可信度。

2.2 模态计算

工程方法处理之后，机翼模态结果如图8所示。表 1列出了修改刚度前后全机模型模态分析所得前8阶模态频率和振型。可以看到，修改刚度后模态分析结果与原始模型保持高度一致性，各阶模态频率误差均在1%以内。

表1 原始模型与修正后模型的模态结果对比

图8 机翼模态图

3 颤振计算

颇振分析是用来确定气动弹性系统的动稳定性问题，一般情况下气动弹性颤振方程可以写为：

(3)

(4)

其中：p=ω(γ±i)=2πf(γ±i)为方程复特征值，γ为衰减率，g=2γ为结构阻尼比，QIQR为广义气动力阻尼矩阵和广义气动力刚度矩阵。选择MSC.Nastran的SOL145计算模块开展颤振分析，由于高阶模态局部模态偏多，因此值截取中前八阶模态参与颤振计算。

P-K法理论上可进行任意运动形式下的颤振计算分析，但是由于在计算的过程中使用的Theodorson 气动力仍然是结构在谐振荡下产生的气动力，因此在颤振点出的分析结果是较为准确的，而非临界点处得到的计算结果仅能够作为参考数据，并不一定准确[21]。使用P-K法进行颤振分析虽然计算量较大，但是求得的结果与飞行颤振试验和风洞试验结果比较接近，是一个近似真实阻尼解，而且可以得到不同动压下的颤振速度。P-K法具体的计算流程是在给定动压条件下，以来流速度v作为自变量，在颤振计算时对指定主要结构模态根据提前设定好的一系列减缩频率点反复进行迭代计算，最终得到相应的气动力，特征根和特征值。本文中根据MSC.Nastran的SOL145模块得到的计算结果，绘制对应的V-g和V-f图。

在颤振分析中，还需要指定气动网格的划分情况。在本文计算中气动网格划分情况见图9，将机翼划分4个气动区域，内外翼段分开处理，其中内翼面划分为3个区，副翼单独划分为1个区，外翼段1个区，气动网格总数为640个，由1.3节的有限元分析可得，蒙皮承受气动载荷会发生失稳，所以优化插值点选择，只选择在主梁和翼肋上的结构节点作为插值点。图10所示为颤振分析中前5阶结构模态插值至气动面上的结果，表明插值信息无异常。

图9 气动网格

图10 气动面插值情况

3.1 颤振计算结果

计算得到颤振计算结果如图11、图12所示。

图11 V-G

图12 V-F

从V-G图看出两个模态发生穿越，分别是二阶弯曲和三阶弯曲，二阶弯曲穿越速度为26 m/s，三阶弯曲穿越速度为61 m/s，所以颤振速度为26 m/s。水平一弯模态也发了穿越，但是阻尼小于0.001，属于小阻尼模态不是颤振。

从V-F图可以看出二阶弯曲和三阶弯曲对应的颤振是传统的爆发性弯扭颤振。原理是扭转模态正好处于二阶弯曲模态与三阶弯曲模态之间，随着速度的增加这些模态频率接近发生相互耦合，出现颤振。一阶弯曲模态在54 m/s处频率降为0，发生扭转发散现象，扭转发散速度为54 m/s。一阶弯曲模态与扭转模态频率差距比较大，两者之间没发生颤振耦合。

综上可得机翼的颤振速度为26 m/s，颤振速度计算值偏大，原因在于增加了机翼刚度。但是修改后模态频率相差小于1%，真实的颤振计算结果不会出现显著变化。机翼的颤振速度小于扭转发散速度，所以应该优先考虑颤振问题。该无人机设计最大飞行表速14 m/s，考虑颤振安全余量15%，所以该机翼的颤振计算结果满足颤振设计要求。

3.2 吊舱布置对于颤振速度的影响

现有研究已表明，在机翼吊装外挂在某些程度上会改变机翼的颤振特性，分析机翼外挂对颤振的影响成为飞机尤其是具有大展弦比柔性机翼飞机气动弹性分析的一个必要过程。常见带外挂的战斗机通常都采用小展弦比三角翼，翼展较小且刚度高，因此外挂对于颤振影响不是很突出，但是类似于“捕食者”这种大展弦比无人机在带外挂后对颤振速度影响较大，因为外挂会现在影响大展弦比机翼刚度分布。此外对于小展弦比机翼，面内的弯曲变形较小，对气动弹性的影响可以忽略，但是在大展弦比机翼分析中机翼面内的弯曲变形是不可忽略的，因为展弦比的增加会使机翼的平面内特性更加明显，面内变形的程度也极大增加。且对于太阳能飞机这种结构来说，由于柔性增加且翼载减少，在加挂吊舱后，吊舱的重量和相对于机翼的位置对于机翼的模态会有很大的改变，因此会显著影响到机翼的颤振速度。

本文中模型中吊舱总质量约为40 kg，与机翼结构质量相当，如此大的集中质量，吊挂在机翼的某一位置，会对机翼的颤振速度产生很大的影响。吊舱布置对于颤振的影响可分为两部分：一部分是外伸的发动机对于颤振速度的影响；另一部分为吊舱的展向占位于颤振速度的影响。

首先保持吊舱的展向位置和发动机质量不动，更改发动机连杆长度，得到其与机翼颤振速度的关系如图13。

图13 发动机连杆长度与机翼颤振速度关系

保持吊舱的展向位置不动，发动机连杆长度不变，修改发动机质量得到其与机翼颤振速度的关系见图14。

图14 发动机重量与机翼颤振速度关系

从以上算例可以看出，发动机连杆越长，发动机越重，机翼的颤振速度越高。从模态角度看，这是因为发动机越重距离主梁越远，机翼的扭转模态越低，机翼越不容易发生颤振。从另一方面说明，该机翼重心在主梁后方，可以通过增加前缘配重，使得重心前移，减小扭转重心与重心的距离，从而提高颤振速度。值得注意的是发动机是由总体需求确定的，重量不可以随意更改，可以通过在相同位置处增加配重，提高机翼颤振速度。

最后保持吊发动机连杆长度和发动机重量不变，修改吊舱站位，得到其与机翼颤振速度的关系见图15。

图15 吊舱位置与机翼颤振速度关系

可以看出吊舱越靠近机身，颤振速度越高。这是因为吊舱内移，会使得机翼扭转模态频率显著提高，距离相邻的弯曲模态频率差距增加，不容易发生耦合，因此更不容易发生颤振。虽然吊舱内移可以显著提高机翼颤振速度，但是不是越内移越好。吊舱离机身太近会使得其对机翼的气动力卸载作用减弱，导致机翼翼尖变形过大，对于飞行控制和安全带来不利影响。

从以上分析中可以得出结论，二阶水平弯曲是太阳能机翼需要特别注意的一种模态，这种模态与机翼的一阶扭转模态的耦合可能导致颤振速度降低。这种类型的颤振速度随着机翼弯曲变形的增加而减小。在设计太阳能无人机时，必须特别注意扭转刚度的影响，避免二阶弯曲、一阶扭转和三阶弯曲的频率接近，以确保二阶弯曲模态中涉及的颤振模态具有足够高的颤振速度。