钝头机体用FADS-α的求解算法及精度研究

王 鹏

(中国航天空气动力技术研究院， 北京 100074)

引 言

随着测试技术的进步及先进飞行器对飞行参数测试精度需求的提高, 传统的基于Pitot管及风标传感器的大气数据系统局限性越来越大. 尤其对于高超声速飞行器而言, 由于气动加热严重, 探出的传感装置已经不适用于测量飞行参数; 同时, 随着高性能飞行器对隐身性能的需求, 传统的飞行参数测量系统已经满足不了实际飞控需求[1]. 针对某些型号的再入飞行器而言, 再入环境的严酷性(热流及动压水平)使得传统的探出的传感装置不可用. 呈扇形散开配置的传感器(类似于航天飞机上的系统)可以用于大气数据测量,但如何融入飞行器机身及系统的校准使得基于探针的传感系统非常昂贵. 一种更经济、 技术上可行的方法便是嵌入式大气数据传感(flush air data sensing, FADS)系统.

FADS系统的基本思想是, 通过一组非探出的配置在表面的测压孔测得的压力值推出大气数据, 并不需要探针深入周围的流场中计算大气数据. 该方法可以避免由于小的头部曲率半径引起的高热流的影响, 扩大了大气数据系统的应用范围, 从亚、 跨、 超声速直至高超声速领域. FADS系统可以直接集成到飞行器的头部, 并不需要活动部件. 由于FADS系统并不需要探测周围的流场, 而是根据飞行器头部的压力分布来解算飞行参数, 因此更容易校准. FADS系统本质上是基于飞行器表面压力测量的一类方法. 针对飞机及航天飞行器已经对该技术进行了广泛的研究. NASA研究报告表明设计的机载FADS系统已成功用于航天飞机[2-3]及X-33空天飞行器[4]的返回过程中, 并在有人控制飞行器F-14[5]及KC-135[6]的飞行测试中验证了FADS系统的关键技术. 研究结果证实了FADS系统在飞行器机动飞行中应用的可行性. 国外研究表明, FADS系统已成功应用于多种前身呈钝体的飞行器[7], 并在具有尖楔前体的飞行器中进行过可行性试验验证[8-9].

纵观目前用于钝头体FADS系统的攻角解算方法, 主要有4类算法. 第1类是解析算法, 攻角通过测压孔组合配置直接解算得出, 并不需迭代求解, 这类算法的代表就是Cobleigh等[7]发展的经典三点式算法, 该算法简单易行, 且精度较高, 不存在算法收敛问题, 但存在多解的问题, 且对于测压孔配置约束严格; 第2类是拟合算法, 通过选定配置在不同位置的测压孔组合, 建立攻角的拟合关系式. 典型代表是Weiss[10]提到的Bohn等基于线性理论开发的五孔探针算法; 为了系统比对算法的精度, 本文发展了基于非线性理论的五点拟合算法, 统一归类为拟合算法中. 第3类是基于代理模型的算法, 例如神经网络算法、 模糊逻辑建模方法等. 对于该类算法, 由于并非建立在气动模型基础上, 经验参数较多, 目前研究人员进行过小范围数据的测试, 试验点数不多, 并未形成统一性的结论[11-15]; 第4类是迭代算法, 基于气动模型, 通过最小二乘法等迭代算法得到攻角. 该类算法精度较好, 对测压孔配置约束低, 但存在对初值敏感及收敛性问题, 且实时性不好[16-18]. 就工程实施而言, 不占优势, 但从解算算法角度考虑, 本文一并进行了系统论证.

综上所述, 针对钝头机体用FADS系统攻角解算的几种求解算法, 各有优缺, 但尚未进行过系统比对, 未形成统一性的结论. 特别是国外在钝头机体用FADS系统取得巨大发展及应用的情况下, 国内对于钝头机体用FADS系统的求解算法及精度等考量尚不完善, 部分试验论证结果很差. 基于此, 本文从FADS系统现存及发展的4大类5种攻角解算方法实现精度上进行考量, 对涉及的5种算法进行了对比及验证, 对钝头机体用FADS-α求解精度进行了系统的论证, 包括经典三点式及改进三点式算法、 基于线性理论的五孔探针算法、 基于非线性理论的五点拟合算法、 基于神经网络建模的方法及基于压力模型的加权最小二乘迭代算法. 以期为钝头机体用FADS系统的算法优化及改进提供指导, 推进FADS技术进一步向前发展.

1 钝头机体用FADS系统理论模型

FADS系统原理是通过飞行器表面特定区域的压力分布反推得到飞行参数, 须建立将表面压力与飞行参数关联起来的压力模型; 同时, 该模型须满足的条件为: (1)适用于较大的Mach数范围; (2)足够简单, 通过简化的模型描述复杂的流场问题.

因此, FADS系统气动压力模型把势流模型(适用于亚声速条件)与修正的Newton流模型(适用于超声速条件)通过形压系数结合起来. 形压系数综合考虑了气动外形、 系统因素等影响, 可以将其看作是Mach数、 攻角及侧滑角的函数, 飞行前可以通过风洞试验或CFD计算得到.

FADS系统的理论模型[7]为

pi=q(cos2θi+εsin2θi)+p∞

(1)

式中,pi为第i个测压孔测得的表面压力,q为冲击动压,p∞为静压,ε为形压系数, 其为Mach数、 有效攻角及有效侧滑角的函数

ε=func(Ma,αe,βe)

θi为第i个测压孔处的来流入射角(该点的法线方向与来流方向的夹角), 由式(2)确定

cosθi=cosαecosβecosλi+sinβesinφisinλi+sinαecosβecosφisinλi

(2)

式中,αe,βe分别为有效攻角及侧滑角;φi,λi为第i个测压点的圆周角及圆锥角. 测压点i的圆周角φi及圆锥角λi的定义如图1所示.

图1 测压孔i圆周角及圆锥角定义Fig. 1 Definitions of clock angle and cone angle for pressure port i

2 测压孔配置

压力输入作为FADS系统解算的基础, 其测压孔配置至关重要. 针对钝头体的测压孔配置, 已经形成了一套固有的测压孔选取标准. 针对15°钝头体外形, 在Mach数Ma=2.04, 3.02, 5.01, 攻角α=-5°～30°, 无侧滑角的条件下, 测压孔配置如图2所示. 本文主要基于数值计算的方法得到配置在前体表面各个测压孔的压力分布.

图2 钝头体测压孔配置Fig. 2 Pressure port configuration for blunt fore-bodies

测压孔的具体位置信息如表1所示. 其中, 测压孔1为驻点测压孔, 测压孔2, 3, 6, 7, 10, 11位于纵向对称面中心线(攻角平面)上, 测压孔4, 5, 8, 9, 12, 13位于侧向对称面中心线(侧滑角平面)上. 各个测压孔的用途及不同解算方法中用到的测压孔的详细信息将在下述章节中结合相应算法详细介绍.

表1 测压孔位置信息Table 1 Detailed location information for pressure ports

本文测压孔配置原因分析如下， 对于标准钝头体, 位于钝头前缘上且对来流攻角变化最敏感的测压孔位置应该位于λ=45°的位置, 即测压孔2～5的配置, 对其的敏感性分析已在相关文献[19]系统论证过, 本文只引用相关结论, 不做具体分析. 但是为了系统比较不同类算法的精度, 本文在λ=30°及λ=75°两组站位上也配置了测压孔.

因此, 基于本节配置的测压孔, 分别采用经典的三点式算法、 改进三点式算法、 基于线性理论的五孔探针算法、 基于非线性理论的五点拟合算法、 基于神经网络的拟合算法及加权最小二乘迭代算法, 对钝头机体用FADS-α求解精度进行了系统的对比及论证.

3 算法精度分析

3.1 三点式算法及改进三点式算法

本节将三点式算法与改进三点式算法合并考虑, 主要基于改进三点式算法是在三点式算法基础上发展起来的, 但三点式算法需要解析求解, 改进三点式算法需迭代求解, 两者差别较大.

根据式(1)定义的理论模型建立的经典三点式算法是基于沿纵向中心线对称配置的测压孔组合. 一般包括驻点测压孔, 其余两个测压孔沿纵向对称面中心线与驻点测压孔对称配置, 根据压差消去形压系数, 从而得到当地攻角与压差的数学关系. 具体为式(3)所示

(3)

其中,pi,pj,pk为对应编号为i,j,k的测压孔得到的表面压力. 将式(3)展开即得到经典三点式算法的原始方程(4)为

Γikcos2θj+Γjicos2θk+Γkjcos2θi=0

(4)

其中, 压差Γik,Γji,Γkj定义为

Γik=pi-pk，Γji=pj-pi，Γkj=pk-pj

通过给定驻点测压孔及沿纵向对称面中心线与驻点测压孔对称配置的测压孔位置, 即可以解析求解得到攻角数值.

当αe≤45°时, 取

αe=0.5arctan(A/B)

当αe>45°时, 取

αe=0.5(π-arctan(A/B))

其中

A=Γiksin2λj+Γjisin2λk+Γkjsin2λiB=Γikcosλjsinλjcosφj+Γjicosλksinλkcosφk+Γkjcosλisinλicosφi

对于式(3)做出如下变换, 将cosθi的表述式(2)带入, 同时利用数学中的正切函数半角关系, 得到最终的关于攻角的数学表达式为

(5)

其中, [c0c1c2c3c4]的定义参见文献[20],本文不再详述. 式(5)即为改进的三点式算法的原始方程. 该方程很难通过与经典三点法类似的方法求解, 通过带入特殊配置的测压孔位置信息, 得到攻角的解析表达式, 而只能通过迭代求解, 迭代求解就涉及算法的收敛及合理解的选取问题.

图3所示为改进三点式算法的流程图. 本文选定外形的侧滑角为0°, 所以算法中不涉及侧滑角的求解, 仅仅对攻角进行验证. 对于给定的攻角初值, 利用Newton-Raphson方法迭代计算, 仅需要通过8～10步的迭代, 便可使得收敛残差小于10-8. 该算法输出的是当地攻角, 即未进行校准过的攻角, 从而可以客观地比较初始攻角求解精度. 得到当地攻角后, 便可以进行其余飞行参数的求解.

图3 改进三点式算法的FADS-α解算流程Fig. 3 Modified triple algorithm process for FADS-α

采用三点式算法得到的钝头体FADS-α的精度分析, 作者已在相关文献中进行过详细阐述[21-23]. 当地攻角及侧滑角的具体选择方案, 参见文献[7]. 在得到当地攻角及侧滑角后, 必须要进行校准, 得到攻角修正量(上洗角)及侧滑角修正量(侧洗角)的校准曲线, 进而得到真实攻角及侧滑角. 当地攻角及侧滑角的解算是相互独立的, 本文的研究仅涉及攻角的解算. 得到的不同测压孔配置对攻角求解精度的影响如图4所示, 校准曲线如图5所示. 可以看出, 基于三点式算法, 只要测压孔配置合理, 经校准后攻角的求解精度很好.

图4 Ma=2.04, 3.02, 5.01不同测压孔配置得到的攻角偏差曲线Fig. 4 Error distributions for angle of attack based on different pressure port configurations for Ma=2.04, 3.02, 5.01

图5 攻角校准曲线Fig. 5 Calibration curves for angle of attack

针对经典三点式算法及基于不同测压孔选取组合的改进三点式算法的精度进行了分析. 图6所示为基于测压孔对称配置的经典三点式算法原始求解精度与基于不同配置的改进三点式算法的比较. p1,2,3表示基于测压孔(1, 2, 3)组合的经典三点式算法得到的攻角原始求解精度; MTA表示改进三点式算法得到攻角求解精度; 可以看出, 与经典三点式算法相比, 改进算法的求解精度变化不明显, 但是由于采用的测压孔数目明显增加, 对于测压孔的敏感性降低, 如果某一测压孔出现故障, 所解算的攻角精度所受影响不大; 而对于三点式算法, 如果某一测压孔出现故障, 结果是灾难性的.

图6 攻角的算法精度校准曲线Fig. 6 Algorithm accuracy calibration curves for angle of attack

3.2 拟合算法

3.2.1 基于线性理论的五孔探针算法

基于线性理论的五孔探针算法, 基本原理是采用5个测压孔组合, 其中必须包括驻点测压孔, 建立与攻角、 侧滑角及Mach数等飞行参数呈线性变化的函数关系, 从而根据所建立的函数关系来推出飞行参数. 该方法简单易行. 采用的各个参数基本是相互独立的, 相互之间没有相关性, 一方面简化了求解过程, 同时带来的缺陷就是攻角增大时精度变差. 当攻角大于某一临界值时, 线性变差, 从而使得攻角的求解精度较差.

本文建立的3个参数为kα,kMa,kβ, 其定义为

式中,kα仅仅与攻角相关, 与Mach数及侧滑角无关;kMa仅仅与Mach数相关, 与攻角及侧滑角无关,kβ仅仅与侧滑角相关, 与攻角及Mach数无关.α～kα关系如图7所示.

图7 攻角与kα关系曲线Fig. 7 Curves for angle of attack with kα

可以看出,α～kα仅仅与攻角相关, 与Mach数等其他飞行参数无关. 因此, 只要得到kα, 便可以通过线性拟合的方法得到攻角, 简单易行. 但是, 可以明显看出当攻角增大时, 线性变差, 从而使得攻角的求解精度变差. 本文得到的α～kα的线性关系式为

α=0.06kα

图8所示为α～kMa关系曲线,kMa与攻角无关, 仅与Mach数的变化相关. 因此, 只要得到足够多的kMa曲线, 便可以通过线性差值方式得到Mach数等其他飞行参数.

图8 Ma=2.04, 3.02, 5.01时攻角与kMa关系曲线Fig. 8 Curves for kMa with angle of attack for Ma=2.04, 3.02, 5.01

3.2.2 基于非线性理论的五点式算法

本节中的五点式算法原理与3.2.1节中五孔探针算法一致, 都属于拟合算法. 但这里选择测压孔建立的函数关系与相关的飞行参数不是线性关系. 因此, 采用加权最小二乘拟合的算法, 将攻角、 Mach数等飞行参数拟合为所建立的各个参数的函数. 本节仅考虑与攻角相关的函数关系的建立. 定义

pα=(p3-p2)/p1

其中，pα～α的关系如图9所示. 可见, 不考虑侧滑角的情况下,pα仅与Mach数及攻角相关, 且随攻角的变化非常明显, 规律性很好.

图9 pα～α变化曲线Fig. 9 Curves for pα～α

不考虑侧滑角的情况下建立的拟合关系式为

本文假定Mach数已知, 所以上式可以通过最小二乘拟合的方法进一步简化为pα的多项式, 采用下述形式

具体为

(6)

图变化曲线

Mach数已知的情况下, 对于攻角的拟合精度进行了分析. 攻角原始曲线如图11所示. 根据此方法得到的攻角拟合曲线如图12所示. 可以看出, 攻角拟合精度较好.

图11 攻角原始误差曲线Fig. 11 Initial error distribution for angle of attack

图12 攻角校准曲线Fig. 12 Calibration curves for angle of attack

3.3 神经网络算法

神经网络技术作为一种非线性建模预测手段, 由于在自适应性、 模糊推理能力及自学习能力等方面的优势, 非常适用于解决非线性问题. 针对钝头机体用FADS系统的特点, 本文建立了BP神经网络, 以代替FADS系统空气动力学模型, 使得FADS系统的解算不再依赖传统的空气动力学模型(输入量依靠气动模型确定, 解算方法完全替代气动模型), 免去了相关校准参数的标定. 同时鉴于神经网络本身的容错功能,各层中个别单元出现错误也不会出现灾难性的后果.

BP神经网络的普通算法流程如图13所示. 网络模型的拓扑结构分为输入层、 隐含层和输出层, 其信号前向传递, 误差反向传递. 信号经输入层输入, 经隐含层逐层递推, 最后转入输出层. 且每一层的神经元状态只影响下一层神经元状态. 输出层与期望值相差大于设定的误差, 则转入反向传递, 根据设定误差调整网络权值和阈值, 最终使得神经网络模型的输出不断逼近预期输出.

图13 BP神经网络的算法流程图Fig. 13 BP neural network algorithm process

本文采用测压孔得到的表面测压数据作为输入量, 建立表面压力与飞行参数的函数关系. 作为输入量的压力数据, 应与输出量有必然内在变化关系. 尽管神经网络可以训练任意的数据集, 得到与输出的关系, 但如果没有理论上的合理支撑, 那训练结果没有意义. 基于本文建立的测压孔配置及敏感性分析的结果, 如图14所示, 建立的用于解算攻角的神经网络模型如图15所示. 图14仅为测压孔2, 7的压力系数随攻角的变化趋势, 其他测压点的变化趋势图类似.

图14 测压孔2和7压力系数随攻角变化Fig. 14 Pressure coefficient with angle of attack for Port 2 and 7

图15 解算攻角的神经网络算法结构Fig. 15 Neural network architecture for angle of attack

基于测试状态, 生成网络模型训练及验证所需的数据集——训练用数据及测试用数据. 训练用数据用于神经网络中的权值和阈值的训练, 监控计算过程中的训练误差; 测试用数据用于测试最终得到的网络的泛化能力, 得到泛化误差. 为了避免过拟合问题, 选用Bayes调整法来改进Levenberg-Marquardt训练算法. 在Mach数Ma=2.04, 3.02, 5.01, 攻角α=-5°～30°, 无侧滑角的条件下范围内选择计算状态(共计约2 500个数据点)进行数值计算, 获得足够多的样本点. 其中, 2 400个数据点用于神经网络的训练, 另外的约100个数据点用于对训练好的网络进行测试及验证.

鉴于FADS系统精度的需求及其对测压孔的敏感性要求, 本文建立的含有双隐含层的4层神经网络模型如图16所示. 通过试错法来确定合适的隐含层的节点个数, 取值范围设为10～20. 训练数据及测试数据依据计算条件生成.

图16 双隐含层神经网络模型Fig. 16 Neural network model with double hidden layers

其中, 输入层I即为选定的测压孔的压力数据, 将驻点压力作为输入量时， 输入层参数为5个测压点数据(1, 2, 3, 6, 7); 隐含层1为

a2=f2(IW2,1+b2)

隐含层2为

a3=f3(IW3,1a2+b3)

输出层为

a4=y=f4(IW4,3a3+b4)

式中,IW2,1为输入层与隐含层1之间的权值,IW4,3为输出层与隐含层2之间的权值,b为阈值,f2,f3为双曲正切函数,f4为斜率为1截距为0的线性函数. 网络隐含层均采用双曲正切函数作为激活函数, 输出层使用线性函数为功能函数. 作者已在相关文章中对详细的算法流程进行过系统的验证[24]. 攻角的训练误差分布与测试误差分布如图17， 18所示.

图17 攻角训练误差分布Fig. 17 Training error distribution for angle of attack

图18 攻角测试误差分布Fig. 18 Testing error distribution for angle of attack

可以看出, 在包含数据包线的范围内的攻角测试误差小于0.1°, 通过神经网络建模的方法得到的攻角求解精度较高.

3.4 基于压力模型的迭代算法

利用第1节中建立的FADS系统的压力模型, 将表面压力分布用势流理论与修正的Newton流理论通过形压系数结合起来, 利用式(2), 将位置关系带入压力模型(1)中, 得到

pi=q(cos2β(c1cos2α+c2sinαcosα+c3)+sinβcosβ(c4cosα+c5sinα)+c6)+p∞

式中,ci(i=1,2…6)为与测压孔圆周角与圆锥角相关的参数, 利用测压孔配置的对称性可以确定

c2=c4=c5=0

从而减少变量个数, 可以得到

pi=q(c1cos2αcos2β+c3cos2β+c6)+p∞

若不考虑侧滑角, 有

pi=q(c0cos2α+c)+p∞

式中,c0,c为与测压孔圆周角及圆锥角相关的参数. 具体解算流程为: 给定迭代变量的初始值, 然后根据建立的不同Mach数下的压力数据库对所用的测压孔的压力进行修正, 进入迭代算法求解

δpj+1i=pi-Fj(B)

定义

Bj+1=Bj+δBj+1

式中

通过递归最小二乘算法

表2 攻角求解精度Table 2 Solving accuracy for angle of attack

4 结论

采用经典的三点式算法及改进三点式算法、 基于线性理论的五孔探针算法、 基于非线性理论的五点拟合算法、 基于神经网络建模方法及攻角迭代算法, 对钝头机体用FADS-α求解精度进行了系统的论证, 得到的结论为：

(1)三点式算法、 改进三点式与加权最小二乘迭代算法精度相当, 只要测压孔配置适当, 通过这3种方法得到的攻角精度都很高, 都可以比较准确地预测攻角. 三点式算法对测压孔配置约束严格, 改进三点式算法及迭代算法对测压孔配置约束较宽松; 且攻角迭代算法对初始值较敏感, 算法存在收敛性问题, 不利于工程实施.

(2)神经网络算法精度较好, 攻角测试误差小于0.1°, 但是该算法是基于无物理模型的算法, 物理背景薄弱, 且需要大批量数据的训练及验证, 算法涉及的参数较多, 且经验性很大, 因此对神经网络算法在FADS系统的应用, 有待进一步论证.

(3)拟合算法优劣明显. 基于线性理论的五孔探针算法精度在小攻角时与上述几种算法精度相当, 简单易行, 但是攻角增大(大于10°)后精度下降显著. 因此, 基于线性理论的五孔探针算法适合于小攻角的情况, 精度比其他算法稍差. 而基于非线性理论的五点拟合算法精度较高, 但是具体拟合过程复杂繁琐.

综上, 对钝头机机体用FADS-α而言, 三点式算法、 改进三点式算法及加权最小二乘迭代算法都是较好的计算方法, 得到的攻角精度都较高. 具体工程实施中算法的选择须结合相关技术指标考量. 本文仅对于攻角的求解精度进行了系统的分析, 针对其他飞行参数的求解及精度的论证将在其他文章中予以发布.