基于移动最小二乘法的气动力数据建模方法

周志高，黄 俊，刘志勤，黎茂锋，李光伟

(西南科技大学 计算机科学与技术学院，四川 绵阳 621010)

0 引言

近年来高超声飞行器的优越性能使其成为各科技强国研究的热点技术，而建立准确的气动力模型是飞行器工程中进行飞行控制、稳定设计以及飞行仿真的重要基础和前提[1]。目前投入应用的气动模型主要来自于风洞试验数据，通过将气动力状态量和控制量以气动数据库表格的形式存储，这样的数据库容量非常大，包含范围很广的飞行包线，可以比较精确地模拟飞行器的整体非线性空气动力学。这种方法通过数值计算和风洞试验手段来获得气动数据将会耗费大量的时间、金钱、计算资源和人力，成本大[2]。因此，对于研制成本和研制周期受限的飞行器来说，尽量减少风洞试验状态，提高数据拓展使用有效性，在最少的试验状态下，建立起建立准确的气动模型显得非常必要[3]。

本文采用具有局部拟合效果的移动最小二乘法(Moving least square method,MLSM)来完成非线性气动力数据拟合。1981年，移动最小二乘法由Lancaster和Salkauskas提出[4]，用于数据拟合和曲面构造，而后又有许多学者对其进行研究改进，其中陈美娟等[5-6]在改进的移动最小二乘法中选取带权重的正交函数作为基函数以及提出复变量移动最小二乘法。崔小曼等[7]在改进的移动最小二乘法中引入Tikhonov正则化，对系数矩阵施加约束项从而得到精确解，避免病态方程组的形成。于成龙等[8]利用遗传算法来求取样本点最优支撑域半径，倪洪杰[9]则利用粒子群优化算法来得到最优支撑域半径，冷亚洪[10]利用对支撑域内抽样点数寻优来获取最佳半径。

经实验研究发现，除了支撑域半径，权函数的选取对模型的精度也会产生很大的影响。因此，本文在移动最小二乘模型的基础上，以模型误差为优化目标，以支撑域半径和权函数影响因子为优化设计变量，采用遗传算法对其进行优化，得到最佳支撑域半径值和权函数影响因子，最后代入气动力模型中验证模型结果的准确度。

1 基于MLS的气动力建模方法

1.1 气动力建模简介

气动力数据建模的目的就是利用CFD模拟数值计算的数据、动态风洞试验数据建立空气动力学数学模型和参数，飞行器气动力和力矩系数一般是攻角数、马赫数、飞行高度等状态参数的函数。即Ci=f(α,Ma,h,δz,…),式中，Ci表示飞行器各气动系数，α为攻角数，Ma为马赫数，h为飞行高度，δz为俯仰舵偏角。

目前，飞行器的气动力建模方法主要有两个方面：1)基于传统的物理特性建模，如多项式模型、三角函数模型等。2)基于人工智能等新型技术的建模方法，常见的有BP神经网络模型、Kriging模型等。

1.2 移动最小二乘近似原理

移动最小二乘法相比传统的最小二乘法引入了支撑域的概念。如图1所示，球形区域表示支撑域，评估点只受支撑域内的样本点的影响,并且距离越近的点对其影响程度越高，这种方法有效解决了分段拟合不能局部化处理的缺点。

图1 支撑域散点示意图

具体原理如下，假设已知一系列离散样本点xi∈Ω(i=1,2,3,…,N)的函数值为f(xi)，移动最小二乘法拟合函数为φ(x)可表示为(本文原理性解释都以一维为例)：

(1)

线性基：PT(x) = [1,x],m=2

二次基：PT(x) = [1,x,x2],m=3

三次基：PT(x) = [1,x,x2,x3],m=4

对于多维问题可以依次类推，本文将会依次对二维问题和四维问题进行研究，随着问题个数增加，计算量也会增加。

系数ai(x)的选取是使近似函数φ(x)在计算点x的邻域内是待求函数f(x)的最佳近似,而拟合函数φ(x)在所有节点的误差加权平方和为：

(2)

对J 求最小值，有：

(3)

由矩阵形式表示可得：

A(x)a(x)=B(x)f(x)

(4)

因此系数向量a(x)为：

a(x)=A-1(x)B(x)f(x)

(5)

式中：

A(x)=PT(x)W(x)P(x)

(6)

B(x)=PT(x)W(x)

(7)

f(x)=[f(x1),f(x2),…,f(xN)]T

(8)

将式(5)代入式(1)后，我们可以得出拟合函数为：

φ(x)=PT(x)a(x)=PT(x)A-1(x)B(x)f(x)

(9)

在MLS拟合函数求解中，为了保证矩阵A可逆，要求支撑域内至少有m个样本点，可以设置合适的支撑域半径。式(2)中，w(x-xi)称为权函数，其只与样本点和预测点的距离有关，即只有在预测点周围某邻域内有值。一般常用的权函数有：高斯函数如式(10)所示，五次样条函数如式(11)所示，指数函数如式(12)所示。

(10)

(11)

(12)

权函数的选取对拟合精度的影响很大[11]，如图2所示,曲线越平缓，权函数的全局性越强，而曲线越陡峭，权函数的紧支撑性越强。支撑域内各样本点对预测点的贡献程度也会随权函数的不同而不同，从而会影响模型准确度。

图2 权函数的全局性

2 移动最小二乘模型优化算法

在MLS近似方法中，支撑域半径的选取以及权函数选择直接影响其对样本点的拟合精度[11-12],如果支撑域半径太大，不能体现局部性，而支撑域半径太小，则会使矩阵A不可逆或者病态。目前主要有两种方法来确定支撑域半径，其中一种是根据经验公式来确定一个相对较优的值；另外一种是在经验公式的基础上对支撑域半径进行优化，进而得到最佳支撑域半径值。

由以上分析可知，支撑域半径影响拟合精度的本质是改变了支撑域内所包含的样本点数，而权函数的选择则改变了不同样本点对评估点的影响程度。因此，为了得到较高精度的气动模型，本文选取高斯函数作为权函数，以高斯权函数影响因子β和支撑域半径为优化设计变量，通过留一法计算出模型误差，并将其作为优化目标，再采用遗传算法对优化模型中设计变量进行优化搜索，直至达到结束条件。

具体算法步骤如下：

1)将气动力数据集划分成均匀训练集和均匀测试集，确定训练集中的节点数N、测试集中的节点数M。

2) 选取高斯权函数为权函数，并确定影响因子β的优化范围为[1,9],在这个范围中可以寻求合适的β来保证局部拟合效果。

3)根据经验公式确定支撑域半径d的优化范围[dmin,dmax]，经验公式如下：

d=sacle×c

(13)

式中，scale是大于1的常量，可根据实验效果进行调整，c是均匀样本点间的距离。

4)将训练集中的初始样本点依次代入MLS模型中，通过留一法循环计算出训练集中各节点的预测值φ(xi)，同时确定优化模型如下：

(14)

式中，e为移动最小二乘模型计算值与真实值的残差平方和,其值越小表示MLS模型精度越高。

5)以式(14)为优化模型，设置种群大小、迭代次数、编码方式等参数，采用遗传算法获取最优支撑域半径d和最优权函数影响因子β。

6)基于优化时使用的训练样本集合，将最优支撑域半径d和最优权函数影响因子β代入模型中,得到所求解的气动力模型。

算法流程图如图3。

图3 算法流程图

3 算例及其结果分析

本文忽略了雷诺数、侧滑角以及升降舵的影响，采用数值计算得到的数据进行建模。对应实验数据有：俯仰舵偏角δz= -20°、-15°、-10°、-5°、0°、5°、10°、15°、20°，飞行高度h=0、10、20、30、40、50(单位：km)，Ma=2、4、6、8、10、12、14、16、18、20，α=-15°～15° (共31组)。因而共有9×5×10×31=13 950个状态数据样本点。

为了验证MLS模型经优化后模型精度的提升(记为改进模型)，同时为了减少计算时间，首先考虑在俯仰舵偏角δz= 0°、飞行高度h=0 km情况下，对攻角数α、马赫数Ma和气动系数之间的关系建模。此时包含31×10=310个状态数据点，将样本数据集均匀划分成训练集和测试集，保证节点之间的距离相等，其中训练集样本用于移动最小二乘建模计算。

通过MLS模型预测一个评估点主要经过两次计算，一是依次计算评估点与训练样本点的距离d，通过权函数公式可以得到权值矩阵W；二是基于训练样本计算出基函数PT。最后利用式(6)～(9)可以计算出对应评估点的预测值。其中由经验公式选取支撑域半径值的方式得到的MLS模型称为经验模型，而改进模型则是经过寻优算法求得最佳设计变量值。

将改进模型分别用来预测训练集和测试集中的数据点，并对比经验模型得到的预测结果对比结果如表1所示，定义当前模型预测误差计算公式如下：

(15)

表1 模型优化前后预测误差对比

表1中，a、b分别表示模型在测试集和训练集中的预测误差，从表中可以看出，改进模型在轴向力系数CA、法向力系数CN和俯仰力矩系数Cmz上有着更好的预测效果，并且MLS模型在测试集上也有较好的预测效果，表明了移动最小二乘法建立的模型具有一定的泛化能力。

参考文献[13]中基于偏最小二乘法(PLS)的建模方法。一般来说，多项式模型的项数越多，模型项中信息保留度也就越大，进而更能反应出真实模型。为了形成对比，本文采用与移动最小二乘法中基函数相同次数的多项式进行主成分分析(PCA)，通过留一法交叉验证确定主成分个数，得到了对比结果。

图4和图5给出了两种建模方法关于3种气动系数的建模结果对比。其中“ErrCA”表示轴向力系数在对应样本点的预测误差值，可以看出偏最小二乘法的建模结果出现明显偏差，并且移动最小二乘法的建模结果整体上优于偏最小二乘法。图 6和图7是移动最小二乘法建模结果和偏最小二乘法建模结果对测试集上样本点的预测，同样可以看出偏最小二乘模型中轴向力系数的预测结果有明显误差。

图4 移动最小二乘法建模结果

图5 偏最小二乘法建模结果

图6 移动最小二乘法预测建模结果及误差

图7 偏最小二乘法预测建模结果及误差

从图中可以看出，对训练集和测试集上的样本点进行预测时，移动最小二乘法建立的模型都表现的更准确。如式(16)，移动最小二乘法的系数为节点x的函数，并且由式(6)可知a(x)融入了权函数带来的局部拟合特性，因而在相同次数的多项式拟合中，移动最小二乘法的建模效果要优于偏最小二乘法。偏最小二乘法拟合方程式的系数ai为常数，同时会根据主成分个数删掉一些信息保留度小的项，因此在偏最小二乘法在计算量上要明显小于移动最小二乘法。

(16)

为了得到适用于大空域、宽速域的气动模型，以飞行高度h、俯仰舵偏角δz、攻角数α、马赫数Ma为自变量，将所有实验数据分为训练数据集和测试数据集，由于参与计算的样本点个数以及问题维数增加，使得优化需要的时间过长。因此在有4个输入情况下的气动建模算例中略去了模型优化步骤，直接利用经验公式确定支撑域半径及影响因子β，得到适用于飞行高度为0～50 km、飞行马赫数在2～20 mach的气动力模型。图8给出了4个输入情况下的的气动系数建模预测结果。为了方便展示，以攻角数α为x坐标，以马赫数Ma为坐标，在同样攻角和马赫上会存在不同飞行高度和俯仰舵偏角的响应值，并以平面图展示各个样本点的预测误差。从图中可以看出，建立的气动模型预测误差较小，效果较优，满足现代飞行器对宽速域和大空域的需求。

图8 四变量输入情况下的建模预测结果

4 结束语

本文利用某飞行器多次飞行仿真的气动力数据作为实验数据，将移动最小二乘法应用于气动力数据建模。考虑到移动最小二乘法的拟合精度对支撑域半径以及权函数比较敏感，提出了一种优化算法。首先以攻角数和马赫数为自变量进行建模，经过实验研究分析，确定支撑域半径和影响因子β的优化范围，接着通过遗传算法获取全局最优设计变量。经算例验证，优化后得到的移动最小二乘模型精度有明显提高，且该方法与偏最小二乘方法相比预测数值也更准确。最后以攻角数、马赫数、飞行高度、俯仰舵偏角为自变量进行建模，由于计算量大，跳过优化步骤，其预测误差也较小，得到的模型能满足“大空域、宽速域”的气动模型要求。移动最小二乘法的建模效果整体上优于偏最小二乘法，但移动最小二乘法的计算时间却明显长于偏最小二乘法，尚有待研究。