应用前馈神经网络的大柔性机翼阵风响应分析

2022-08-26刘松林谢长川

计算机测量与控制 2022年8期

常辉，朱靖，安朝，刘松林，谢长川，杨超

(1.北京航空航天大学航空科学与工程学院，北京 100191； 2.中国人民解放军91977部队，北京 100036；3.深圳市百川融创科技有限公司，广东深圳 518110)

0 引言

对于大柔性机翼而言，气动载荷会引起较大的结构变形，这意味着传统线性系统分析方法中的小变形假设将不再适用，结构的力学平衡关系必须建立在变形之后的构型之上，这时气动弹性问题会转化为几何非线性问题[1-2]。从力学本质来看，结构位移与应变呈现非线性关系，而应力应变的本构关系仍满足线性关系[3]。当机翼产生较大的结构变形时，其气动力分布和结构动力学特性与线性情况有明显差异，且会随着结构的变形及载荷变化而改变，因此构建几何非线性气动弹性分析方法，成为大柔性机翼设计发展的主要问题。

结构建模是气动弹性分析中的核心问题，大柔性飞机结构建模需要准确描述结构位移和应变之间的非线性关系，反映大变形带来的静、动刚度变化的非线性效应。针对非线性气动弹性响应，国内外研究人员提出了很多结构建模方法，但在模型计算效率、精确度和模型适用性上无法同时达到最优。因此构建准确的结构大变形预测模型是本文研究的重点。

目前相对成熟的建模方法是以位移为变量的非线性有限元理论，与线性方法不同，位移基非线性有限元理论能合理描述微元变形前后的转动、载荷和应力关系，主要有拉格朗日增量法[4]、基于共旋坐标法(CR)[5]等。拉格朗日增量法计算耗费大、效率较低，CR有限元法虽然效率较高，但是在计算中需要进行复杂的坐标系转换。而大变形梁理论[6-8]计算精度高且易于和二元气动力结合，对于简单模型例如梁等结构具有明显的优势，但对于复杂模型需要进行等效，会造成精度的损失。以应变为变量的有限元理论[9-12]以几何精确梁为基础，可以更方便地描述应变能等物理量，对于梁单元模型计算效率高、收敛速度快，但是对于复杂模型仍需要进行简化，同样会造成精度的损失。除此之外比较常用的方法是非线性结构降阶理论[13-16]，通过给定基函数将结构动力学方程从物理空间转换到低阶空间中。常用结构固有模态作为基函数，将物理坐标方程转换为模态坐标方程进行求解，此方法能够在满足精度求解的基础上提高计算效率，但对于传统的非线性结构使用模态降阶需要在每一次变形之后重新计算结构的模态坐标，且需要进行非线性修正。

实际上，复杂非线性模型构建主要有两种途径：第一个是基于理论的模型架构，即根据物理问题的力学方程，建立理想的数学模型。第二种则是数据驱动方法，即通过分析系统仿真或试验中得到的样本数据，直接构造黑箱或灰箱模型。近年来，随着神经网络等人工智能技术的发展，数据驱动型的机器学习方法也逐渐应用于非线性模型的改进和构建中。机器学习通过算法从数据中建立模型，使之具备一定的判断和预测能力，具有泛化能力强、计算效率高等特点，是目前非线性系统建模的有效工具。机器学习方法在大柔性结构变形预测上的研究还处于探索阶段，现有的研究成果也显示出这一研究方向的发展潜力巨大。

本文以大柔性结构建模为切入点，针对大柔性梁式机翼模型，采用人工神经网络的机器学习算法搭建结构大变形预测模型，结合曲面气动力搭建完整的非线性气动弹性分析框架，计算机翼阵风响应，探究人工神经网络模型应用于大柔性结构建模的可行性，为今后人工智能技术与非线性气动弹性分析结合的研究提供思路。

1 结构大变形预测基础理论

将结构动力学方程进行有限元离散后可表述为：

(1)

其中：M和K分别为物理空间中的线性质量矩阵和刚度矩阵，u(t)为结构位移向量，fnl表示非线性恢复力，f(t)表示结构外力。

在物理空间中结构自由度较高，因此考虑采用模态降阶技术，通过给定基函数将结构运动方程从物理坐标转化为模态坐标，从而降低方程维度。

给定如下假设：

(2)

其中：φi(X)为基函数，qi(t)为待定系数即广义坐标，m为模态阶数。

将方程代入方程中，化简，得到降阶后的结构动力学方程：

(3)

其中：Kq为模态空间中m阶对角刚度矩阵，fq(t)为广义激振力，S(q1,q2,…,qm)是关于所有广义坐标的非线性函数。由于非线性环节S的存在，结构位移的广义坐标、广义速度及广义激振力之间存在非线性关系，很难直接求解动力学方程。

对于结构动力学非线性问题，当结构受到给定静载荷或动载荷时，结构变形的动态响应过程是确定的。因此神经网络模型需要预测的对象是广义坐标的时域响应，大变形结构建模过程变成一个非线性动态系统辨识过程，此时需要辨识的系统的输出是关于输入的非线性泛函，这种泛函的拟合对于处理离散变量的神经网络模型来说是做不到的。此情况下可采用有限个采样值代替时域响应，以差分形式表示系统输入输出，将非线性动态系统离散化。

2 多层前馈神经网络

本文采用多层前馈神经网络进行预测，其拓扑结构如图 1所示。它由一个输入层、一个或多个隐层和一个输出层组成，每一层都包含若干个神经元，输入层、输出层神经元个数由解决的问题本身决定，隐层个数选取没有严格的规定，一般根据输入、输出的维数进行调整，隐层神经元的个数影响着神经网络模型的性能。

图1 多层前馈神经网络结构

网络中各节点之间输入输出关系为[17]：

(4)

oj,k=f(neti,k)

(5)

式中，neti,k为第k层第i个神经元的输入，oi,k为该神经元的输出，bi,k为阈值，wi,j,k为第k层第i个神经元与第k-1层第j个神经元之间的连接权值，f为该神经元的激励函数。对于输入层、输出层神经元，其激励函数一般取线性函数，而对于隐层神经元，激励函数取为Sigmoid等非线性函数，来赋予神经网络非线性拟合能力，Sigmoid函数为S型函数，其表达式为:

(6)

神经网络通过对训练样本的学习，调整神经元之间的连接权值和阈值，使网络输出逐渐逼近样本输出。前馈神经网络的训练可通过多种学习算法进行，其中应用最广泛的是误差反向传播算法(BP，error backpropagation)，其基本思路是根据网络输出值与期待输出值之间误差，沿着误差梯度下降的方向，反向逐层调整个连接权值和阈值，经过反复迭代学习后确定最小误差对应的网格参数。

学习过程由信号的正向传播和误差的反向传播组成。在正向传播过程，输入信号从输入层进入网络，经过隐层的处理后由输出层输出，这个过程只涉及到信号的传播，对神经网络模型没有造成影响，各连接权值和阈值没有发生改变。当网络输出未达到样本预期值时，则转入反向传播过程，将输出值与期待值之间的误差信号沿原通道返回，传递的过程中修正各权值及阈值，使误差最小。

将误差定义为:

(7)

式中，n为神经网络总层数，q为输出神经元个数，dj为给定的预期输出。

网络参数wi,j,k和bi,k的修正沿负梯度方向进行，修正公式为:

(8)

式中，r为迭代步数，η为学习步长。

根据链式求导法则，可得到以下偏导数计算公式：

(9)

式中，神经元的输出oj,k-1由前向通道中计算；输出误差对神经元的偏导数δi,k通过误差的反向传播中逐层递推计算，计算公式如下：

k为输出层时：

δi,k=di-oi,n

(10)

k为其他层时：

δi,k=f′(neti,k)∑j(δj,k+1wj,i,k+1)

(11)

3 神经网络大变形预测模型

构造大变形预测模型的困难在于确定非线性环节的表达式。在非线性结构降阶模型中，非线性函数被定义为广义坐标的二次多项式和三次多项式之和[18]，可通过最小二乘回归分析确定非线性刚度系数，得到结构动力学方程完整表达式。

大变形预测模型的输入变量为结构载荷的广义力，输出变量为结构位移的广义坐标及广义速度的瞬态响应，在模态空间下表示成多维向量，因此该系统为多输入-多输出系统，需要保证同一变量不同模态坐标在时间序列上的阶次保持一致，即广义位移/广义速度/广义力的不同模态阶次相同。另外，由结构动力学方程可得，结构动力学系统在某一时刻的描述不需要考虑结构载荷的动态属性，当前这一时刻结构所受的外载荷可以充分反映系统的瞬态载荷信息，因此该非线性系统的输入阶次为可设为1，因此多输入多输出的模型可写为：

qi(k+1)=f1,i[q1(k),…,ql(k),…,q1(k-n1+1),…,

(12)

结构变形预测模型针对的预测目标主要是当前时刻的结构位移的广义坐标，在动力学问题中由于结构惯性项的加入，在离散时刻下广义位移不再与结构载荷的广义力呈一一映射关系，因此输入特征还需要考虑结构速度、加速度等信息，来增强模型输入与输出的相关性。样本特征选择的思路在于确定系统输出时间序列的阶次n1、n2。分析发现，当输入的广义位移的阶次n1超过3，广义速度的阶次n2超过2时，输入可间接包含加速度信息，由于广义速度可作为样本数据变量，当需要模型输入体现加速度信息时，可将广义速度在时间序列上的阶次设为2。

结合以上分析，根据输入特征包含信息多少，可定义以下三种输入特征组合方案进行对比验证，如表1所示。

表1 输入特征组合方案

在确定神经网络输入和输出特征后，就可以适用神经网络模型对非线性函数f1,i、f2,i进行拟合。对于多输入多输出系统，对输出变量按照数量级大小进行分组来降低神经网络模型使用数量。在本文研究的结构变形预测模型中，输出主要考虑广义位移和广义速度之间的影响，因此将输出变量分为广义位移和广义速度两组分别进行预测，同时增加单个神经网络的模型结构的方案，在仿真中对两种模型结构进行对比验证，分析模型结构对预测精度的影响，变形预测模型结构如图2所示，其中x为模型输入，q、q′分别为多维广义位移和广义速度。

图2 变形预测模型结构

采用梁式平直机翼模型参数如表 2所示，有限元模型如图 3所示，结构模型节点个数为213个。

表2 机翼模型参数

图3 机翼有限元模型

本次仿真验证实验用于神经网络的样本数据由MSC.Nastran有限元软件进行非线性有限元分析得到的，对于给定机翼梁模型，施加不同激励，得到结构在不同工况下的动力学响应，以此作为训练样本集。

首先选择降阶模型所使用的模态。结构线性模态是线性结构降阶模型中常用的基函数，对于非线性结构降阶函数也同样适用。对所需要计算的机翼模型进行结构模态分析后，选择如表3所示的前六阶模态：

表3 前六阶线性模态

接下来选择简单的动载激励方式，即在翼尖加载正弦集中激励，翼尖最大位移约为翼展的27%，满足大变形条件。其激励形式以式表示：

F=Asin(ωt)

(13)

式中，F为施加载荷大小，A为激励幅值大小，ω为激励频率。

样本数据以时间进程排列，记录了载荷作用下机翼5 s内的结构变形动态响应，样本数据形式包括每一时刻各阶模态下结构位移的广义坐标、广义速度及结构载荷的广义力，以2.5 ms为间隔，得到2 000个样本点。

本文神经网络采用Matlab神经网络工具箱，网格参数、激活函数、权值更新算法均为默认值。该网络为三层前馈网络，包括输入层、隐层和输出层。输出层激活函数为线性函数，隐层激活函数为sigmoid函数。

隐层神经元个数选择10个，根据不同样本特征选择方案及不同的模型结构，神经网络输入、输出维数各不同。各方案模型的输入及输出个数由表 4所示，若选择方案一，网络输入维数为12；选择方案二，网络输入维数为18；选择方案三，网络输入维数为24。若选择单个神经网络的结构(a)，网络输出维数为12；若选择多个神经网络的结构(b)，每个网络输出维数为6。

表4 神经网络模型输入、输出量个数

将样本数据随机分为：70%训练集、15%验证集、15%测试集。训练集的作用是用来拟合模型，更新权值和阈值。验证集估计计算误差，当训练集误差降低但验证集误差升高时，停止训练，以防止过拟合现象的出现，并返回具有最小验证机误差的权值和阈值。测试集的作用时通过训练集和验证集得出最优模型后，使用测试集训练网络以衡量该最优模型的性能。

当训练停止后，用均方误差MSE和相关系数R作为衡量模型性能的指标，其表达式分别如下：

(14)

三类样本特征方案训练结果如图 4所示。

图4 三类样本特征选择方案训练结果

根据训练结果可看出，当选用方案一时相关系数最高仅为0.537 47，模型拟合程度较差，均方误差也较大，说明该方案选择的输入特征组合与目标特征的相关性较弱，不适合作为样本输入；选用方案二时，各数据集的相关系数均接近1，均方误差较小，约为0.065 232，模型拟合程度高；当继续增加输入样本维数，选用方案三时模型拟合程度也较高，但均方误差约为0.019 431，相比方案二并没有显著减小，且由于维数增加，模型计算效率降低，因此最佳特征组合方案为方案二。

选定样本特征方案后，对不同神经网络模型结构进行仿真验证。对图 2中(a)单个神经网络组成的模型结构和(b)多个神经网络组成的模型结构进行样本数据集的训练，得到两种模型结构的预测结果与目标值的总均方误差及各阶模态坐标的均方误差，对比结果如表 5所示。

表5 两种模型结构的样本均方误差对比

结果表明与单个神经网络组成的模型结构相比，多个神经网络并联的模型结构对广义位移的预测精度有明显提升，两者的均方误差达到约为3个数量级的差距，广义速度的预测精度虽有所下降，但基本处在同一数量级上，属于可接受范围内，综合考虑选择多个神经网络并联的结构模型作为变形预测模型结构。

4 阵风响应分析

阵风响应是指处于稳定状态的气动弹性系统遭遇阵风扰动时的动态响应。本文采用基于神经网络的大变形预测模型和曲面涡格法[21]气动模型，搭建大柔性机翼气动弹性分析框架。在阵风响应计算中，初始状态为静气弹平衡构型，因此本文使用的静气动弹性平衡构型通过传统的非线性结构降阶模型得到。

如图 5所示，根据结构静平衡构形、阵风强度以及频率等信息，利用曲面涡格法计算该时刻的气动载荷，并进行力插值[22]，将载荷信息传递到神经网络预测模型。神经网络模型可根据当前时刻模态空间下的广义力及上一时刻的广义位移和广义速度，预测当前时刻的结构位移和结构速度，经过位移及速度插值后，更新气动面和流场边界条件，进行下一时间步的非定常气动载荷计算。同时当前时刻预测得到的结构位移和速度作为下一时刻神经网络模型的输入，继续预测下一时刻的结构位移和速度。

图5 大柔性阵风响应分析框架

神经网络模型的预测精度很大程度上取决于训练样本的质量和复杂度。样本质量一般指训练样本分布反应总体分布的程度；样本复杂性是指训练神经网络所需的样本数量[23]，样本数过少会导致模型拟合不充分，样本数过多又会造成模型过度拟合。

要提高神经网络模型对不同工况下的结构变形预测能力，从训练样本选择的角度有两种途径：一是通过增加样本复杂性来覆盖待解决问题对应的样本分布，二是选择能较好反映总体分布的高质量样本。本文建立的神经网络预测模型是一个18维输入、6维输出的系统，其对应的样本空间较为庞大，若通过增加样本复杂性来覆盖待解决问题对应的样本分布，所需的样本数据集会十分复杂和庞大的。因此提高样本质量是训练样本获取和选择的主要途径。

本文从基于结构降阶模型的阵风响应计算中获取机翼的阵风响应数据。给定阵风不同频率和来流速度，计算3°攻角下机翼受到阵风扰动后的时域响应，时间间隔为0.002 5 s，计算2 400步共6 s，记录每一个时间步下的结构位移和速度及结构载荷在模态空间下的广义坐标。各数据集情况如表 6所示，得到6组数据集。由于从开始遭遇突风到达稳定振荡状态之间有一段计算不稳定的过程，在样本集中属于无效样本数据，因此只取2 s后稳定状态的动态响应作为样本数据。，其广义位移的一阶模态坐标的响应如图 6所示。

图6 样本中广义位移的一阶模态坐标响应

表6 训练样本数据集情况

根据前一节得到的神经网络模型，按照方案二整理样本集。将同频率不同风速对应的数据集合并成一个训练样本集，用于该阵风频率下不同来流速度的气动弹性响应的计算，即训练得到两个神经网络模型分别用于计算3 Hz和4 Hz阵风频率下的阵风响应。

以阵风频率为4 Hz、来流速度为12 m/s的阵风响应计算为例，将4 Hz频率下来流速度为10 m/s、12 m/s、14 m/s对应的样本组合成新的数据集，各样本数和比例为1 500∶1 500∶1 500，共4 500组样本数据，训练神经网络模型。使用归一化的均方根误差(nRMSE)对训练集的预测结果进行评估，nRMSE公式为：

(15)

表7 神经网络模型对目标的归一化均方误差值

图7 阵风频率4 Hz、来流速度12 m/s下翼尖垂直位移预测对比

从神经网络模型的训练结果来看，预测值与真实值的归一化均方误差处于较低的水平，说明神经网络模型对样本集的预测精度较高。而从阵风响应计算结果来看，由于神经网络预测模型是从2 s后稳定状态下开始计算的，前2 s翼尖响应于真实值一致。2 s后，翼尖垂直位移计算值逐渐偏离于真实值，偏离值随时间推进不断增大，最后在新的位置稳定振荡。可见即使神经网络模型能对样本进行较高精度的拟合，也不能保证阵风响应计算的准确性。

经过分析，造成神经网络训练结果和阵风响应计算结果不匹配的原因，有两个方面。一是由于误差的积累，在阵风响应计算中，结构模型输出的结构位移和速度的信息，会作为下一时刻的输入，在这个过程中神经网络模型造成的输出误差会传递到下一时间步，并随时间增加不断累积，引起输出值预测值的误差不断增大；二是由于神经网络模型对输入的敏感性，当神经网络拟合的模型与真实模型差别较大时，模型对训练集以外的数据预测能力较差，输入抖动时输出会产生较大偏差，即神经网络模型对输入的敏感性过高。神经网络的敏感性将上一步计算造成的输出误差对模型输出的影响放大，产生更大的输出误差随着计算的迭代进一步累积、增大，使最终的阵风响应结果与真实值有明显的差异。误差的累积特性是由计算模型结构造成的，不可避免，因此要得到准确的阵风响应结果，就需要从优化预测模型性能以改善神经网络的敏感性入手。

在本次实验中，不同风速的训练集在样本空间中的位置较为分散，因此当其他训练集参与神经网络的训练后可能会对拟合效果产生较大影响，导致神经网络敏感性增加，其影响的程度可能与不同训练集之间的比重有关。为了验证上述猜想，可调整不同来流速度的训练集比例训练网络，观察得到的阵风响应结果。代入阵风响应计算模型中求解响应会耗费较大计算资源，考虑到气动力模型求解精度较高且对输入的敏感性较低，可直接给定样本中准确的气动力数据作为神经网络模型输入，而不通过气动力模型，这样可得到近似的广义位移垂直一阶弯曲模态响应曲线，也可以反映模型计算结果准确性。改变总训练集中不同训练集的比例，4 Hz频率下来流速度为10 m/s、12 m/s、14 m/s的训练集比例分别取1 500：1 500：1 500、1 500：1 300：1 300、1 500：1 300：1 000，训练神经网络得到2 s到4.5 s内的近似广义位移垂直一阶弯曲模态响应曲线，与真实值进行对比，如图 8所示。

图8 不同样本比例下近似广义位移垂直一阶弯曲模态响应曲线

图8可反映出不同训练集比例对神经网络性能的影响，当4 Hz频率下来流速度为10 m/s、12 m/s、14 m/s的训练集比例为1 500：1 300：1 000时，近似的阵风响应能稳定振荡且与真实值相近，说明该训练集得到的神经网络模型能够进行准确的阵风响应计算。对3 Hz频率下的训练集进行相同的实验和调整，得到最终调整后10 m/s、12 m/s、14 m/s的训练集比例为1 000：1 000：600。

根据以上训练集比例训练神经网络模型，计算3°攻角机翼受阵风扰动后的时域响应，时间步长为0.002 5 s，计算4 000步共10 s，阵风频率给定3 Hz和4 Hz，风速给定10 m/s、12 m/s、14 m/s。图 9和图 10分别为大柔性机翼模型在3 Hz、10 m/s和4 Hz、10 m/s工况下翼尖位移的动态响应过程与真实值的比较曲线，可以看出目标值与预测值曲线较为一致，说明预测结果较好，其它工况拟合曲线与图 9和图 10类似。表 8统计了上述所有工况下翼尖垂向位移的模型预测值与真实值的归一化均方误差，结果表明神经网络结构模型对训练集拟合程度较高，阵风响应预测结果准确。