摘要
当前,数值模拟作为研究飞机结冰的主要手段之一,在计算结冰冰形时会引入大量参数不确定性,并影响数值模拟的精度和可信度。发展不确定性量化方法,科学量化这种不确定性对评估数值模拟结果具有重要意义。针对传统参数不确定性量化方法难以解决高维输入到输出的问题,基于本征正交分解和误差反向传播神经网络,提出了一种结冰冰形预测代理模型。以水滴中值粒径和温度为例,验证了代理模型单输入参数和双输入参数情况下的精度和泛化能力。最后,在代理模型计算的冰形基础上,结合蒙特卡洛采样,利用准则确定结冰范围,发现水滴中值粒径不确定性主要影响明冰的冰角生长,而温度和水滴中值粒径不确定性的叠加主要作用于霜冰厚度。该研究为后续多结冰条件的影响分析和多维输入到输出的不确定性量化提供了思路。
飞机穿越云层时,空气中的过冷水滴撞击在迎风表面并冻结成冰,会直接改变飞机的气动外形,引发飞机气动特性的变化,比如升力减小,阻力增加等,结冰严重时,甚至会引发飞行安全事
当前,数值模拟是研究飞机结冰的主要手段之
随着数十年的发展,参数不确定性量化方法不断得到完善和拓展。参数不确定性量化方法可以划分为统计型和随机数学型两大
目前,在飞机结冰不确定性量化研究中多采用多项式混沌方法,其将随机空间以扰动的形式表现出来,具有很好的收敛速度。文献[
结冰条件不确定性对冰形的影响是典型的多维输入输出问题,使用多项式混沌方法容易造成“维数灾难”。神经网络(Neural network, NN)适用于此类多维输入到输出问题的建模,但具有大量冗余特征的样本数据会导致模型的训练效率低,精度不高等问题。已有研究表
本文针对结冰数值模拟高维输入与输出的特点,通过POD获取低维状态下的结冰冰形基模态,并利用BPNN的强学习和自适应能力,构建了结冰冰形计算代理模型,经过与CFD模拟结果对比,验证了代理模型的精度。随后,在所提代理模型的基础上,以统计方法获取结冰范围,用于分析来流条件的不确定性影响,并得到了相应的结论。
本文的研究目的是量化来流条件不确定性对结冰冰形的影响,而传统的参数不确定性量化方法通常是以代理模型的形式建立输入与输出之间的映射关系,并根据输出响应量的统计特性来表征这些不确定性。结冰冰形常被离散为物理空间的网格点,使得结冰冰形的计算具有典型的高维输入到高维输出的特点。本文针对飞机结冰高维输入与输出的特点,通过拉丁超立方采样(Latin hypercube sampling,LHS)获取的来流条件参数,并使用CFD模拟采样数据,得到结冰冰形原始样本数据。结合POD的降维特征和BPNN能模拟任意非线性映射的优势,构建计算结冰冰形的POD⁃BPNN代理模型,最后基于已验证的代理模型,以统计方法执行来流条件参数的不确定性量化分析。模型的构建及不确定性分析具体流程如

图1 翼型结冰分析流程图
Fig.1 Airfoil icing analysis flow chart
采用数值模拟计算结冰冰形,主要包括流场计算、水滴场计算和结冰计算3个步
计算流场时,本文是通过求解低速粘流的时均N⁃S(Navier⁃stokes)控制方程来获取空气流场分布,N⁃S方程可以写为
(1) |
式中:从左至右分别代表非稳态项、对流项、扩散项和源项,而各项的表达式及整个方程的数值离散和求解方法可参考文献[
水滴场的计算中,使用拉格朗日方法建立水滴运动轨迹方程,获取水滴撞击轨迹为
(2) |
式中:为水滴质量,为水滴位置,和分别为水滴密度和空气密度,为水滴体积,为重力加速度,为水滴迎风面积,为阻力系数,和分别为当地气流速度和水滴速度。
由拉格朗日法建立的水滴运动轨迹方程形式较为简单,求解方法也较为成熟,本文使用的结冰计算程序是采用一阶Euler数值方法来求解水滴运动方程,详细的求解过程可参见文献[
结冰计算则是通过Messinger模
(3) |
式中:为撞击水质量总和,为上一控制体流入当期控制体的液态水质量,为蒸发的液态水质量,为当前控制体流入下一控制体的液态水质量,为单位时间留在当前控制体的液态水质量。
由于方程中未知数较多,需和能量平衡方程进行联立求解。本文使用改进后的Messinger模型将结冰控制体表面的能量平衡方程表示为
(4) |
式中:为冻结水的能量,为蒸发水所需要的能量,为流出当前控制体溢流水的能量,为流入当前控制体溢流水的能量,为由水滴撞击控制体表面所带来的的能量,为摩擦产生的能量,为对流换热,为冰与水之间热传导的能量。
POD方法常用于处理具有复杂的高维输出的问题,其核心是寻找一组正交基使样本集可由这一组基函数的线性叠加来近似,最佳正交基的选取要求样本向量在这组基上的投影取到最大。POD基模态可以通过如下推导获得。
假设某样本数据集合为,其中为维列向量,则有阶样本矩阵,定义样本均值为
(5) |
将样本矩阵去中心化,通过样本矩阵的每一维度减去该维度的均值,使每一维度上的均值为0,得到标准化的样本矩阵,并由标准化后的样本数据构造协方差矩阵为
(6) |
由定义可知,协方差矩阵为阶实对称矩阵,此时正交基的求解可转变为求解特征值和特征向量的问题。通过对矩阵进行奇异值分解获得从大到小排序的特征值及其对应的特征向量。定义映射矩阵能使原始的样本矩阵在其上的投影最大且满足
(7) |
使用POD方法降维时,文献[
(8) |
由
神经网络是一种把未知系统看成一个黑箱,模拟人脑利用神经元与神经元之间的连接,通过记忆和学习来获取系统信息的网络结构(如

图2 全连接神经网络结构
Fig.2 Fully connected neural network structure
BPNN的学习过程包括正向传播和反向传播两个过程。在正向传播过程中,输入信息通过输入层经隐藏层并完成逐层处理传入输出层,得到输出信息。如果获得的输出层信息达不到期望,则转入反向传播过程,将输出与期望的误差平方和作为误差信息,由输出层逐层反向传播。误差反向传播过程中,通过逐层求解误差信息对各神经元的权值的偏导数,构成误差信息对权值向量的梯度,以此修改网络权值。通过循环迭代不断修正权值,直到整个网络的输出信息达到或者无限接近期望后,完成整个网络的学习过程。值得注意的是,训练过程中,梯度的方向决定了误差传播的方向,因此在权重更新过程中需要对其取反,从而减小权重引起的误差。由于本文研究工作是验证POD⁃BPNN模型模拟翼型结冰冰形的精度,以及利用该模型量化结冰条件的不确定性,这里不再对具体算法作更细节的阐述,算法的具体推导过程可以参考文献[
POD重构是为了实现对于任意一个训练集以外的新的输入参数,能通过训练好的BPNN对基模态系数向量进行预测,随后对预测的基模态系数向量进行POD重构,得到输出预测值。当所选取的前L阶基模态能量占比足够大时(本文取),可通过POD重构,对新的输出进行近似拟合
(9) |
式中:为的伪逆矩阵,为期望得到预测值;为经BPNN预测得到的系数向量。
经POD重构后,可实现对任意输入的预测,而代理模型高精度的预测能力是进行不确定性量化的前提,为了定量分析代理模型的预测精度,本文使用平均相对误差进行评估。定义模型预测的平均相对误差为
(10) |
式中:为坐标平均相对误差;为网格点数量;为第个网格点的坐标的预测数据;为第个网格点的坐标的CFD计算数据。
(11) |
式中:为坐标平均相对误差;为第个网格点的坐标的预测数据;为第个网格点的坐标的CFD计算数据。平均相对误差越小,表明模型预测数据偏离原始数据程度越小,精度越高。
这里,以二维NACA0012翼型为研究对象,进行相关分析。执行POD⁃BPNN训练需要大量的样本数据,本文假设来流条件中的水滴中值粒径和温度均各自服从正态分布(NASA给出的水滴中值粒径分布可以近似为均值,方差的正态分

图3 单参数典型工况CFD模拟结果
Fig.3 CFD simulation results under single parameter input

图4 双参数典型工况CFD模拟结果
Fig.4 CFD simulation results under two-parameter input
本算例在来流速度为44.2 m/s、迎角为3°、液态水含量为1 g/

图5 单参数BPNN训练结果
Fig.5 BPNN training results under single parameter input

图6 不同粒径下冰形的CFD和POD-BPNN模拟
Fig.6 Comparison of numerical simulation and prediction results of ice shape with different
工况 | 数值模拟计算时间/ s | POD⁃BPNN计算时间/ s | x坐标平均相对误差/ % | y坐标平均相对误差/ % |
---|---|---|---|---|
1 | 53.6 | 2.03 | 1.631 | 0.311 |
2 | 53.6 | 2.03 | 0.494 | 0.091 |
3 | 53.6 | 2.03 | 0.135 | 0.042 |
从
在假定来流速度为44.2 m/s、迎角为3°、液态水含量为1 g/

图7 双参数BPNN训练结果
Fig.7 BPNN training results under two-parameter input
通过选取3个不同的特征工况:工况4(MVD=10,t=-6 ℃)、工况5(MVD=15,t=-10 ℃)、工况6(MVD=20,t=-15 ℃),进行CFD和POD⁃BPNN模拟,用于验证POD⁃BPNN模型在双参数输入条件下的精度。

图8 双变量结冰冰形模拟结果对比
Fig.8 Comparison of ice shape numerical simulation and prediction results under two-parameter input
工况 | 数值模拟 计算时间/ s | POD⁃BPNN计算时间/ s | x坐标平均相对误差/ % | y坐标平均相对误差/ % |
---|---|---|---|---|
4 | 53.6 | 2.07 | 0.974 | 0.775 |
5 | 53.6 | 2.07 | 1.389 | 0.211 |
6 | 53.6 | 2.07 | 0.770 | 0.580 |
3个工况的预测精度中,工况4和工况6的x坐标预测平均相对误差均小于1%,对工况5的x坐标预测误差较大,但平均相对误差仅为1.389%;对工况4、工况5和工况6的坐标预测平均相对误差均小于1%,表明在双参数输入条件下,POD⁃BPNN模型仍然具有足够的精度。值得注意的是,不同于单参数冰形算例验证过程,双参数结冰冰形算例验证使用的冰形样本数据同时考虑了温度和水滴中值粒径变化对冰形的影响,这保证了所提模型在双参数输入条件下的精度。
前文通过对比单参数和双参数输入条件下结冰冰形CFD模拟与POD⁃BPNN模拟结果,验证了本文所采用的POD⁃BPNN模型预测结冰冰形可靠性,本节将在POD⁃BPNN模型的基础上,进一步研究随机参数水滴中值粒径和温度的不确定性影响。
参数不确定性影响往往是通过输出响应量的统计特性进行表征,样本数量足够大时,中心极限定理指出样本均值的抽样分布近似为正态分

图9 结冰冰形不确定性分析
Fig.9 Uncertainty analysis of ice formation
将表示结冰冰形的397个离散特征经POD方法降维后,原始样本数据被投影到POD模态构成的完备空间中,并被用于构建POD⁃BPNN结冰冰形代理模型。通过设置不同的工况条件,对比不同工况下的CFD模拟结果,并结合蒙特卡洛采样,在POD⁃BPNN计算的结冰冰形基础上执行结冰条件不确定性量化,得到如下结论:
(1)使用POD⁃BPNN预测模型计算结冰冰形时,所提模型在单参数和双参数输入情况下,都表现出较高的精度和较强的泛化能力。但是,由于神经网络训练过程学习到的是经POD降维后最终保留的冰形共同特征,因此,如果冰形较复杂,并且样本数据中大部分样本都共同存在这一特征,则神经网络会捕捉到这种复杂冰形;而当这种复杂冰形不是共同存在的特征时,POD过程会将这一信息当成噪声,从而被忽略。为进一步提高模型对复杂冰形的预测能力,可以通过提高POD过程中能量占比的取值,从而保留冰形的更多信息,来提升模型预测复杂冰形的能力。
(2)在POD⁃BPNN计算的结冰冰形基础上,对结冰条件执行不确定性量化时发现,水滴中值粒径不确定性主要影响明冰冰角的生长,而在同时考虑水滴中值粒径和温度的不确定性影响时,水滴中值粒径的不确定性影响被削弱,而两者相互叠加的不确定性则主要作用于霜冰厚度。
(3)飞机结冰过程是一个多参数耦合且高度非线性的复杂物理变化过程,仅考虑一维和二维输入变量的不确定性影响仍然不足以精确表示出结冰冰形受不确定性的影响程度。虽然本文采取的POD⁃BPNN模型在量化单参数和双参数输入的不确定性影响时具有良好的效果,并且获得了可靠的结论,但是为了保证模型在更高维度输入参数情况下的预测精度,获取的冰形样本数据同样需要考虑相应个数的结冰条件参数,以此来保证预测模型的精度和可行性。
参考文献
苏媛,徐忠达,吴祯龙. 飞机结冰后若干飞行力学问题综述[J].航空动力学报,2014,29(8):1878-1893. [百度学术]
Su Yuan, Xu Zhongda, Wu Zhenlong. Overview of several ice accretion effects on aircraft flight dynamics[J]. Journal of Aerospace Power, 2014, 29(8): 1878-1893. [百度学术]
李伟斌, 宋超, 赵凡, 等. 基于Kriging模型的冰风洞试验冰形参数化方法[J]. 航空动力学报, 2021,36(2): 369-376. [百度学术]
Li Weibin, Song Chao, Zhao Fan, et al. Parameterization method based on Kriging model for ice shape formed in icing wind tunnel[J]. Journal of Aerospace Power, 2021,36(2): 369-376. [百度学术]
郝云权, 赵大志, 李伟斌, 等. 飞机结冰的不确定性量化研究进展[J]. 航空动力学报, 2022,37(9): 1855-1871. [百度学术]
Hao Yunquan, Zhao Dazhi, Li Weibin, et al. Recent advances in the uncertainty quantification of aircraft icing[J]. Journal of Aerospace Power, 2022, 37(9): 1855-1871. [百度学术]
王鹏, 修东滨. 不确定性量化导论[M]. 北京: 科学出版社, 2019. [百度学术]
Wang Peng, Xiu Dongbin. Introduction to uncertainty quantification[M]. Beijing: Science Press, 2019. [百度学术]
WASSERSTEIN, Ronald L. Monte Carlo: Concepts, algorithms, and applications[J]. Technometrics, 1996, 39(3): 338-338. [百度学术]
Fox B L. Strategies for quasi-Monte Carlo[M]. Dordrecht: Kluwer Academic, 1999. [百度学术]
Loh W L. On Latin hypercube sampling[J]. The Annals of Statistics, 1996,24(5): 2058-2080. [百度学术]
Liu W, Belytschko T, Mani A. Probabilistic finite elements for nonlinear structural dynamics[J]. Computer Methods in Applied Mechanics & Engineering, 1986, 56: 61–81. [百度学术]
Deodatis G, Shinozuka M. Weighted integral method I: Stochasticstiffness matrix[J]. Journal of Engineering Mechanics, 1991,117(8): 1851-1864. [百度学术]
Zhang D X. Stochastic methods for flow in porous media[J]. New York: Academic Press, 2002. [百度学术]
Wang P, Tartakovsky A M. Uncertainty quantification in kinematic-wave models[J]. Journal of Computational Physics, 2012, 231(23): 7868-7880. [百度学术]
Sullivan T J, Introduction to uncertainty quantification[M]. UK: Mathematics Institute University of Warwick Coventry, 2015. [百度学术]
XIU D, Karniadakis G E. The Wiener-Askey polynomial chaos for stochastic differential equations [J]. SIAM Journal on Scientific Computing,2002,24(2): 619-644. [百度学术]
Degennaro A M, Rowley C W, Martinelli L. Uncertainty quantification for airfoil icing using polynomial chaos expansion[J]. Journal of Aircraft, 2015,52(5): 1404-1411. [百度学术]
Degennaro A M. Uncertainty quantification for airfoil icing[D]. Princeton: Princeton University, 2016. [百度学术]
Tabatabaei N, Raisee M, Cervantes M J. Uncertainty quantification of aerodynamic icing losses in wind turbine with polynomial chaos expansion[J]. Journal of Energy Resources Technology, 2019,141(5): 1-17. [百度学术]
王晓东,于佳鑫,房代宝,等. 随机风况下风力机翼型结冰对气动特性的影响研究[J]. 风机技术, 2020, 62(2): 59-66. [百度学术]
Wang Xiaodong, Yu Jiaxin, Fang Daibao, et al. Investigations on the influence of icing on aerodynamics of wind turbine airfoil under stochastic wind conditions[J]. Chinese Journal of Turbomachinery, 2020, 62(2): 59-66. [百度学术]
申晓斌, 郁嘉, 林贵平, 等. 基于特征正交分解法的翼型结冰冰形快速预测[J]. 航空动力学报, 2013, 28(4): 807-812. [百度学术]
Shen Xiaobin, Yu Jia, Lin Guiping, et al. Fast prediction of ice shape based on proper orthogonal decomposition method[J]. Journal of Aerospace Power, 2013, 28(4): 807-812. [百度学术]
柴聪聪, 易贤, 郭磊, 等. 基于BP神经网络的冰形特征参数预测[J]. 实验流体力学, 2021, 35(3): 16-21. [百度学术]
Chai Congcong, Yi Xian, Guo Lei, et al. Prediction of ice shape characteristic parameters based on BP neural network[J]. Journal of Experiments in Fluid Mechanics, 2021, 35(3): 16-21. [百度学术]
贾续毅, 龚春林, 李春娜. 基于POD和BPNN的流场快速计算方法[J]. 西北工业大学学报, 2021, 39(6): 1212-1221. [百度学术]
Jia Xuyi, Gong Chunlin, Li Chunna. Fast flow simulation method based on POD and BPNN[J]. Journal of Northwestern Polytechnical University, 2021, 39(6): 1212-1221. [百度学术]
张珂, 张玮, 阎卫增, 等, 圆度误差的神经网络评定及测量不确定度研究[J]. 机械科学与技术, 2019, 38(3): 428-432. [百度学术]
Zhang Ke, Zhang Wei, Yan Weizeng, et al. Research on evaluation and uncertainty of measurement of circularity errors via neural network algorithm[J]. Mechanical Science and Technology for Aerospace Engineering. 2019, 38(3): 428-432. [百度学术]
程曦, 张志勇. 基于人工神经网络的复杂介质中波的传播不确定性分析方法[J]. 电子与信息学报, 2021, 43(12): 3663-3670. [百度学术]
Chen Xi, Zhang Zhiyong. An uncertainty analysis method of wave propagation in complex media based on artificial neural network[J]. Journal of Electronics & Information Technology, 2021, 43(12): 3663-3670. [百度学术]
易贤, 桂业伟, 朱国林, 等. 运输机翼型结冰的计算和实验[J]. 航空动力学报, 2011, 26(4): 808-813. [百度学术]
Yi Xian, Gui Yewei, Zhu Guolin, et al. Experimental and computational investigation into ice accretion on airfoil of a transport aircraft[J]. Journal of Aerospace Power, 2011, 26(4): 808-813. [百度学术]
易贤. 飞机积冰的数值计算与积冰试验相似准则研究[D]. 绵阳: 中国空气动力研究与发展中心, 2007. [百度学术]
Yi Xian. Numerical computation of aircraft icing and study on icing test scaling law[D]. Mianyang: China Aerodynamics Research and Development Center, 2007. [百度学术]
Messinger B. Equilibrium temperature of an unheated icing surface as a function of air speed[J]. Journal of the Aeronautical Sciences, 1953, 20(1): 29-42. [百度学术]
陈希. 旋翼结冰的高精度数值模拟与防/除冰方法[D]. 南京: 南京航空航天大学, 2018. [百度学术]
Chen Xi. Investigations on high accuracy numerical simulation and anti/de-icing methods for rotor icing[D]. Nanjing: Nanjing University of Aeronautics and Astronautics, 2018. [百度学术]
Ostrowski Z, Bialecki R A. Estimation of constant thermal conductivity by use of proper orthogonal decomposition[J]. Computational Mechanics, 2005, 37(1): 52-59. [百度学术]
周志华. 机器学习[M]. 北京: 清华大学出版社, 2016. [百度学术]
Zhou Zhihua. Machine learing[M]. Beijing: Tsinghua University Press, 2016. [百度学术]
黄丽. BP神经网络算法改进及应用研究[D]. 重庆: 重庆师范大学, 2008. [百度学术]
Huang Li. BP neural network algorithm improvement and application research[D]. Chongqing: Chongqing Normal University, 2008. [百度学术]
Gary A, RUFF B M. Users manual for NASA Lewis ice accretion prediction code(LEWICE)[R]. NASA-CR-185129, 1990. [百度学术]