CN105183996B - 面元修正与网格预先自适应计算方法 - Google Patents
面元修正与网格预先自适应计算方法 Download PDFInfo
- Publication number
- CN105183996B CN105183996B CN201510582632.8A CN201510582632A CN105183996B CN 105183996 B CN105183996 B CN 105183996B CN 201510582632 A CN201510582632 A CN 201510582632A CN 105183996 B CN105183996 B CN 105183996B
- Authority
- CN
- China
- Prior art keywords
- grid
- wing
- cfd
- angle
- attack
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
Landscapes
- Aerodynamic Tests, Hydrodynamic Tests, Wind Tunnels, And Water Tanks (AREA)
Abstract
一种面元修正与网格预先自适应计算方法,用多组不同迎角下的CFD气动力载荷数据对低阶面元法进行分段线化修正,同时采用自适应模拟退火算法对面元计算网格的分布进行优化。通过对面元网格的自适应优化提高了修正面元法的精度,弥补了传统修正面元法依赖于网格分布的缺点。优化后的网格在应用于修正面元法时不仅保持了面元法计算效率高的优点,还能保证机翼整体受力接近CFD数据的精度,有效提高气动弹性优化迭代设计过程中的精度和效率。本发明得到的气动载荷数据与CFD计算结果精度误差在2%以内,且有效延伸到气动载荷随迎角变化的非线性段,提高了结构刚度参数变化过程中的气动载荷计算效率。
Description
技术领域
本发明涉及一种航空器结构设计中升力面定常气动载荷数值计算问题,提出了一种面元修正与网格预先自适应优化算法技术,使在常规飞行状态范围内所得气动载荷分布可达CFD(Computational Fluid Dynamics)数值精度,且在应用于结构静气弹重分析及优化设计时的计算效率远高于CFD方法。
背景技术
在航空器结构静气动弹性约束下的数值优化设计中,需要大量的气动载荷与结构变形耦合数值计算,高精度高效率的静气动弹性耦合数值算法成为关键技术之一。现代工程分析中常采用两种方法计算气动载荷:一是基于完全空间网格数值计算的CFD方法,另一是基于小扰动势流方程的面元法。CFD方法计算精度高,与实验结果吻合性好,得到了广泛的研究与应用,但该方法计算量大,计算时间长,多用于航空器气动布局的精细分析及气动特性研究。虽通过并行技术以及降阶方法等能使CFD计算时间得到大幅缩减,但仍无法满足气动弹性迭代优化设计的效率需求;同时CFD计算结果难于高效率给出关于结构变量变化的设计敏度数据。在现代航空器结构初步设计阶段,工程中仍大量采用面元法计算气动载荷以提高计算效率。面元法几何外形适应性强,计算网格数量少,且气动力影响系数矩阵与结构参数无关,是气动载荷与弹性耦合设计的理想工程方法。但由于采用了线性速势理论,面元法无法计及厚度、黏性和激波等因素的影响,这些缺陷限制了其应用范围。近年来,针对面元法的技术缺陷提出了多种修正方法,计算精度有了明显改善。
文献1“严德,杨超,万志强.应用气动力修正技术的静气动弹性发散计算[J].北京航空航天大学学报,2007,33(10),1146-1149.”利用高阶面元法求解气动弹性问题,在一定程度上提高了面元法的精确性,但高阶面元法增加了计算网格处理的复杂性,且仍不能解决气体流动分离等非线性问题。文献2“PC Chen and D.D.Liu.Investigation onTransonic Correction Method for Unsteady Aerodynamics and AeroelasticAnalyses[J].Journal of Aircraft,2008,45(6),1890-1903.”和文献3“AlbertoVarello,Erasmo Carrera1,Luciano Demasi.Vortex Lattice Method Coupled withAdvanced one-Dimensional Structural Models[J].ASD Journal,2011,2(2),53-78.”)提出的模态修正面元法等亦可提高气动载荷的计算精度,但若结构参数变化时,需要反复调用CFD计算数据,这仍然难于解决航空器结构静气动弹性约束下数值优化设计的计算效率问题。文献4“JiaHuan,Sun Qin,Liu Jie.A corrected panel method for staticaeroelasticity[J].Chinese J Comput Phys,2014,31(1),21-26.”提出的分段斜率修正面元法进一步弥补了面元法无法考虑非线性因素的缺陷,但无法考虑面元网格对修正精度的影响。
利用CFD进行气动载荷计算时,网格划分的优劣是影响计算结果精度的关键因素之一。划分计算网格时,需要在复杂流场部位的翼尖、前后缘等位置进行加密处理。对于跨音速等可能产生激波时的情况,更要采用动网格技术,捕捉激波位置,并在激波产生位置加密网格,这些技术原理启发着对面元网格分布进行自适应优化的必要性。
发明内容
为克服现代面元修正方法改善精度效果不佳的缺点,本发明提出了一种面元修正与网格预先自适应计算方法。
步骤1,计算不同迎角下机翼部件原刚体气动载荷数据。
将机翼部件作为刚性几何体划分CFD计算的空间网格;所述的空间网格采用非结构网格;计算亚音速时机翼迎角在0°到20°之间的CFD数据,迎角间隔为1°;提取机翼物面网格节点上的压力分布数据,得到各机翼物面网格节点上压力分布,从而确定不同迎角下机翼部件原刚体气动载荷数据。
在确定不同迎角下机翼部件原刚体气动载荷数据时,翼面采用三角形单元划分,最小尺度为1mm。附面层用三棱柱网格划分,共20层,附面层之外的空间网格采用四面体单元。从附面层最内层开始网格尺度以1.2倍的体积比率增长。CFD数据由FLUENT计算,采用耦合压力基求解器。空气假设为理想气体,温度为300K,参考压力为海平面大气压,湍流模型采用S-A模型,其余参数选择默认设置。
步骤2,计算变形刚体气动载荷数据。
在机翼升力曲线的非线性段选择一个迎角αref,将步骤1中计算得到的该迎角时的气动载荷加载到机翼上,利用NASTRAN软件计算得到翼尖的扭转角。用CAD软件重新几何造形以使翼尖扭转相同的角度且使扭转角从翼根到翼尖线性变化。再按照步骤1中的方式重新划分CFD计算网格并实施常规CFD数值计算,提取机翼变形后各翼面网格节点上的压力分布,得到机翼的变形刚体气动载荷数据。
步骤3,对升力面的弦平面进行面元网格参数初始化。
采用通用的三次多项式插值法,通过展向三次曲线和弦向三次曲线分别对机翼展向和弦向进行网格点布置。
对所述展向三次曲线和弦向三次曲线分别进行归一化;以给定的展向网格点数量范围、弦向网格点数量范围和每条三次曲线两端的斜率值作为约束条件,需初始化的变量为弦向和展向的网格点数量。
步骤4,用几何算法对机翼弦平面进行面元网格划分。
根据给定的三次曲线参数及展向网格点数量和弦向网格点数量对机翼弦平面进行面元网格节点布置。
在对展向布置网格节点时,先对展向三次曲线的横轴按照展向节点数量进行等距划分,此时横轴上各点在三次曲线上对应的纵坐标即表示此节点的展向比例位置,且横轴正方向代表从机翼翼根到翼尖的方向;在对弦向布置网格节点时,先对弦向三次曲线的横轴按照弦向节点数量进行等距划分,此时横轴上各点在三次曲线上对应的纵坐标即表示此节点的弦向比例位置,且横轴正方向代表从机翼前缘到后缘的方向。
对于弦向三次曲线,若需要在机翼前缘加密,则应使曲线在(0,0)点的斜率值减小;对于展向三次曲线,若需要在翼尖加密,则应使曲线在(1,1)点的斜率值减小;若要采用均匀划分形式,则只需要将(0,0)、(1,1)处的斜率值均设为1即可。
步骤5,面元法计算气动力矩阵。
通过ZAERO软件生成气动力矩阵。具体是,将步骤4得到的网格分布输入到ZAERO软件,计算给定飞行状态下的气动力影响系数矩阵[AIC]。上述飞行状态与步骤1中计算CFD数据的飞行状态一致。
步骤6,采用分段斜率修正面元法对气动载荷分布进行修正,获得修正后的载荷分布。
利用步骤1得到的定常刚体CFD计算数据对面元法的气动力影响系数矩阵[AIC]实施修正,具体是采用公式(5)进行分段线化高精度面元修正:
式中,n为迎角范围划分的子区间数目,[WT1]i为第i个迎角区间采用斜率修正法得到的修正因子矩阵,[Cor]i为各面元在第i个迎角区间内的下洗速度值占总下洗速度的比例系数对角矩阵。
步骤7,计算CFD数据与修正面元法结果的误差范数,判断是否达到收敛条件。。若是,则得到最优网格分布,迭代过程结束;反之采用自适应模拟退火算法更新网格分布参数,更新后转至步骤4,重复步骤4~步骤6;在重复步骤4~步骤6的过程中不断获得新的修正后的载荷分布。对得到的新的修正后的载荷分布通过本步骤计算CFD数据与修正面元法结果的误差范数,判断是否达到收敛条件,直至满足收敛条件。
本发明得到的气动载荷数据与CFD计算结果精度误差在2%以内,且有效延伸到气动载荷随迎角变化的非线性段,有效提高了结构刚度参数变化过程中的气动载荷计算效率。
本发明的核心部分可描述如下:用多组不同迎角下的CFD气动力载荷数据对低阶面元法进行分段线化修正,同时采用自适应模拟退火算法对面元计算网格的分布进行优化,以有效提高在结构变刚度情况下的气动载荷计算精度和效率。
本发明通过对面元网格的自适应优化大大提高了修正面元法的精度,弥补了传统修正面元法依赖于网格分布的缺点。优化后的网格在应用于修正面元法时不仅保持了面元法计算效率高的优点,还能保证机翼整体受力接近CFD数据的精度,有效提高气动弹性优化迭代设计过程中的精度和效率。
附图说明
图1是本发明的示意图;
图2是三条典型的网格分布三次曲线;
图3是机翼结构的整体有限元模型,其中,(a)代表从机翼翼根到翼尖的展向方向,(b)代表从机翼前缘到后缘的弦向方向;
图4是机翼前梁、后梁和翼肋有限元模型;
图5是刚体机翼的CFD表面网格;
图6是机翼翼尖扭转前后的构型,其中,6a是机翼翼尖扭转前的侧视图,6b是翼尖扭转后的侧视图;
图7是不同方法的升力曲线对比图;
图8是优化后的非均匀面元网格分布图;
图9是机翼最大位移和总升力变化曲线,其中,9a是机翼最大位移变化曲线,9b是机翼总升力变化曲线;
图10是本发明的流程图。图中:
1.迎角;2.升力/N;3.CFD曲线;4.未修正面元法曲线;5.修正面元法曲线;6.迭代次数;7.结构最大位移/mm;8.升力/N;9.修正面元法曲线;10.CFD曲线;11.前梁腹板;12.后梁腹板;13.翼肋。
具体实施方式
本实施例是NACA0012对称翼型的M6机翼的面元修正与网格预先自适应计算方法。
所述NACA0012对称翼型的M6机翼根部完全固支。机翼参数如下:翼根弦长为0.8139m,翼尖弦长为0.4573m,展长为1.1963m,前缘后掠角为30°,后缘后掠角为15.8°,无扭转角,蒙皮厚度为0.002m,梁腹板厚度为0.0015m,翼肋厚度为0.0015m,材料为铝合金,弹性模量E=70GPa,泊松比μ=0.30,密度为ρ=2.7×103kg/m3。图3为机翼的整体有限元模型,其中,(a)代表从机翼翼根到翼尖的展向方向,(b)代表从机翼前缘到后缘的弦向方向。图4为机翼前梁腹板11、后梁腹板12和翼肋13的有限元模型。采用三角形和四边形壳单元进行有限元建模,机翼由翼肋、梁腹板和蒙皮组成。11个翼肋沿展向均匀分布。梁腹板包括前梁腹板和后梁腹板,分别位于距离机翼前缘27.0%和63.5%处。弦向各网格节点位置分别为距离前缘0.0%、5.0%、10.0%、27.0%、45.5%、63.5%、82.0和100%处,展向各网格节点位置与翼肋位置一致。
本实施例的具体实施过程是:
步骤1,计算不同迎角下机翼部件原刚体气动载荷数据。
通过CFD软件计算不同迎角下机翼部件原刚体的气动载荷数据。将机翼部件作为刚性几何体划分CFD计算的空间网格,所述的空间网格采用非结构网格。翼面采用三角形单元划分,最小尺度为1mm。附面层用三棱柱网格划分,共20层,附面层之外的空间网格采用四面体单元。从附面层最内层开始网格尺度以1.2倍的体积比率增长。CFD数据由FLUENT计算,采用耦合压力基求解器。空气假设为理想气体,温度为300K,参考压力为海平面大气压,湍流模型采用S-A模型,其余参数选择默认设置,机翼的CFD表面网格见图5。计算亚音速时机翼迎角在0°到20°之间的CFD数据,迎角间隔为1°,提取机翼物面网格节点上的压力分布数据,得到机翼物面网格节点上压力分布的数据库。本实施例中,马赫数取0.6,机翼迎角范围为0°到10°。
步骤2,计算变形刚体气动载荷数据。
在机翼升力曲线的非线性段选择一个迎角αref,将步骤1中计算得到的该迎角时的气动载荷加载到机翼上,利用NASTRAN软件计算得到翼尖的扭转角。用CAD软件重新几何造形以使翼尖扭转相同的角度且使扭转角从翼根到翼尖线性变化。再按照步骤1中的方式重新划分CFD计算网格并实施常规CFD数值计算,提取机翼变形后各翼面网格节点上的压力分布,得到机翼的变形刚体气动载荷数据。
本实施例中,取αref为8°,计算得到的翼尖扭转角为向下扭转4°,翼尖扭转前后的CFD网格如图6所示,其中,6a是机翼翼尖扭转前的侧视图,6b是翼尖扭转后的侧视图。
步骤3,对升力面的弦平面进行面元网格参数初始化。
上述的面元网格参数即指网格加密方法的参数。通用的代数辅助网格加密方法有三次多项式插值法、等差数列法、等比数列法、拉伸函数法。本实施例采用能够实现多种加密方式的三次多项式插值法,需要两条三次曲线分别对机翼展向和弦向进行网格点布置,分别记为展向三次曲线和弦向三次曲线。
给定展向和弦向的网格点数量范围和每条三次曲线两端的斜率值作为约束条件,需初始化的变量为弦向和展向的网格点数量,即待优化变量。
三次曲线能够实现在需要的位置处进行网格加密。图2给出了(a)、(b)和(c)三条典型的三次曲线,图中对横坐标和纵坐标均进行了归一化处理;所述归一化处理指将横坐标最大值和纵坐标最大值均缩放为1。因确定一条三次曲线需要四个约束,其中的两个约束为对横坐标和纵坐标进行归一化后曲线两端的位置(0,0)和(1,1),另外两个约束为所述给定的三次曲线两端的斜率;通过所述四个约束能够完全确定三次曲线的形状。本实施例中,如图2所示,三条给定的三次曲线(a)、(b)、(c)在(0,0,)点处的斜率分别为2.6、1、0.1,在(1,1)点的斜率分别为0.1、1、2.6。
步骤4,用几何算法对机翼弦平面进行面元网格划分。
根据给定的三次曲线参数及展向和弦向网格点数量对机翼弦平面进行面元网格节点布置,即划分机翼弦平面的面元网格。
在对展向布置网格节点时,先对展向三次曲线的横轴按照展向节点数量进行等距划分,此时横轴上各点在三次曲线上对应的纵坐标即表示此节点的展向比例位置,且横轴正方向代表从机翼翼根到翼尖的方向;在对弦向布置网格节点时,先对弦向三次曲线的横轴按照弦向节点数量进行等距划分,此时横轴上各点在三次曲线上对应的纵坐标即表示此节点的弦向比例位置,且横轴正方向代表从机翼前缘到后缘的方向。
对于弦向三次曲线,若需要在机翼前缘加密,则应使曲线在(0,0)点的斜率值减小;对于展向三次曲线,若需要在翼尖加密,则应使曲线在(1,1)点的斜率值减小;若要采用均匀划分形式,则只需要将(0,0)、(1,1)处的斜率值均设为1即可。
步骤5,面元法计算气动力矩阵。
通过ZAERO软件生成气动力矩阵。具体是,将步骤4得到的网格分布输入到ZAERO软件,计算给定飞行状态下的气动力影响系数矩阵[AIC]。上述飞行状态与步骤1中计算CFD数据的飞行状态一致。本实施例中,飞行状态为0.6马赫数,海平面大气压,温度300K。
步骤6,采用分段斜率修正面元法对气动载荷分布进行修正,获得修正后的载荷分布。
利用步骤1得到的定常刚体CFD计算数据对面元法的气动力影响系数矩阵[AIC]实施精细修正。
修正前的面元法求解气动力公式为
{Cp}=[AIC]{W} (1)
其中,[AIC]为空气动力影响系数矩阵,{W}为无量纲下洗列阵,{Cp}为各面元上的压差系数列阵。
面元修正法基本思路是将面元法计算的气动压差系数与CFD计算数据进行对比修偏,获得气动力影响系数矩阵的修正因子矩阵[WT],并将修正因子矩阵前乘于气动力影响系数矩阵[AIC],此方法称为斜率修正法,如下式。
[AIC]*=[WT1][AIC] (2)
式中,[WT1]为修正因子矩阵。[AIC]*为修正气动影响系数矩阵。
斜率修正法将定常压差系数随迎角的变化率作为计算依据来计算修正因子矩阵[WT1],见公式(3)、(4)。公式(4)中[WT1]ii代表CFD定常气动力数据随迎角的变化率与面元法定常气动力数据随迎角的变化率的比值。对于进入非线性段的气动力系数修正,显然相邻迎角越相近,修正效果越好,即高精度面元修正。
{Cp}=[WT1][AIC]{W} (3)
其中,[Cpα]igiven表示第i个面元网格上CFD压差系数关于迎角的差分数值;[Cpα]ipanel为第i个网格面元法计算得出的压差系数关于迎角的差分数值;[WT1]ii表示割线修正因子矩阵的对角元素第i项,每个面元对应的修正因子构成了对角修正矩阵[WT1];α1、α2表示两相邻的不同飞行迎角。
由式(4),似给定两个不同迎角的CFD气动力系数即能够得到精确的修正结果,实际这没有意义,本实施例的目的是将此修正因子应用于任意飞行迎角下的弹性变形条件,这可理解为各剖面的迎角不同。于是,用尽可能多个不同迎角下的CFD数据进行修正,所得修正因子无疑是高精度的,但修正因子将会成为随迎角变化的非线性关系。为此,本实施例采用公式(5)进行分段线化高精度面元修正:
式中,n为迎角范围划分的子区间数目,[WT1]i为第i个迎角区间采用斜率修正法得到的修正因子矩阵,[Cor]i为各面元在第i个迎角区间内的下洗速度值占总下洗速度的比例系数对角矩阵。
步骤7,根据式(6)计算CFD数据与修正面元法结果的误差范数,判断是否达到收敛条件。若是,则得到最优网格分布,迭代过程结束;反之采用自适应模拟退火算法更新网格分布参数,更新后转至步骤4,重复步骤4~步骤6;在重复步骤4~步骤6的过程中不断获得新的修正后的载荷分布。对得到的新的修正后的载荷分布通过本步骤计算CFD数据与修正面元法结果的误差范数,判断是否达到收敛条件,直至满足收敛条件。
本实施例中,收敛条件为连续20步最优解不变。
面元网格自适应优化的目标函数为式(6),即CFD数据与修正面元数据的误差。
其中,N是计算CFD数据时的迎角数量,本实施例中为11,L、MX、MY分别表示机翼的气动力、俯仰力矩以及滚转力矩。所述L、MX、MY的下标CFD均表示CFD计算结果,由CFD计算得到的压力系数在翼面积分得到;所述L、MX、MY的下标panel表示面元修正法的计算结果,由步骤6中的压差系数Cp在弦平面上积分得到。
本实施例采用自适应模拟退化算法来优化网格分布,使误差范数最小化。自适应模拟退化算法的基本过程为:把每种组合状态s看成某一物质系统的微观状态,而将其对应的目标函数C(s)看成该物质系统在状态s下的内能;用控制参数T类比温度,让T从一个足够高的值慢慢下降,对每个T用Metropolis抽样法在计算机上模拟该体系的热平衡态,即对当前状态s作随机扰动以产生一个新状态s',如果C(s')<C(s)则接受s'为下一状态,否则以概率e-[C(s′)-C(s)]/T接受。经过一定次数(Markov链长)的搜索,认为系统在此温度T下达到平衡,则降低温度T再进行搜索,直到满足结束条件(“RAO Jian.A Statistical Physics Studyon Spatial Inhomogeneous Systems with Small Number of Particles[D]HunanUniversity,2013.”)。在本实施例中,组合状态s就是展向网格数量和弦向网格数量组成的向量,C(s)为式(6)定义的目标函数值。
图7中给出了分别采用面元法、斜率修正面元法及CFD三种方法计算得到的刚性M6机翼在0.6马赫数下不同迎角的机翼总升力。其中CFD数据为步骤1中得到的原刚体气动载荷数据。
由图7可以看出,迎角从0°变化到10°时,M6机翼气动力基本呈线性增长。未修正面元法曲线4和CFD曲线3的误差随着迎角增大逐渐变大,这是由于该方法不能考虑机翼附面层、粘性等影响导致的。由修正面元法曲线5可以看出,经过斜率修正面元修正法计算的气动升力与CFD结果吻合良好,能够准确反映气动力随迎角变化的非线性特征。这验证了分段线化面元修正方法的有效性和精度。
网格优化时,分别取机翼展向与弦向网格点数范围为[3,15]及[3,25]两种情况,并采用均匀网格划分及非均匀网格划分两种方法。对于非均匀的网格,曲线起点和终点斜率分别设为0.4和0.5。优化结果如表1及表2所示。表中优化结果第一列为展向网格点数量,第二列为弦向网格点数量。表中的Δ1表示总升力相对误差,Δ2表示总滚转力矩相对误差,Δ3表示总俯仰力矩的相对误差。
由表中结果可以看出,通过面元网格的自适应优化,目标函数最小值对应的并非最细网格情形。通过均匀网格划分与加密网格的比较可知,适当加密机翼前缘及翼尖网格可以较大地提高俯仰力矩及滚转力矩的计算精度,并且采用非均匀网格划分所需网格节点更少,这可以提高气动力计算效率。对于气动弹性优化等需要反复气动力计算的领域来讲,其效果非常可观。
表1均匀网格划分结果
表2非均匀网格划分结果
为验证本发明应用于静气动弹性分析时的效率和精度,在此对机翼在4o迎角下的几何变形进行分析比较,面元网格分布采用表2的第二个结果,对应图8。分别用CFD和高精度面元修正法与结构有限元程序计算耦合,反复迭代计算结构变形与气动力。机翼的气动力总和以及结构最大变形随迭代步变化曲线如图9所示。由修正面元法曲线9和CFD曲线10可以看出,CFD与高精度面元修正法均迭代4步收敛,且收敛结果一致。弹性气动力较刚性机翼损失5.6%。计算时间方面CFD迭代共耗时4小时,修正面元法仅耗时2分钟。高精度面元修正法准确计算了弹性变形后的机翼气动力,与CFD结果几乎一致,这说明优化后的网格结合高精度面元修正法在计算弹性变形后机翼气动力方面是准确可靠的。
Claims (4)
1.一种面元修正与网格预先自适应计算方法,其特征在于,
步骤1,计算不同迎角下机翼部件原刚体气动载荷数据;
将机翼部件作为刚性几何体划分CFD计算的空间网格;所述的空间网格采用非结构网格;计算亚音速时机翼迎角在0°到20°之间的CFD数据,迎角间隔为1°;
提取机翼物面网格节点上的压力分布数据,得到各机翼物面网格节点上压力分布,从而确定不同迎角下机翼部件原刚体气动载荷数据;
步骤2,计算变形刚体气动载荷数据;
在机翼升力曲线的非线性段选择一个迎角αref,将步骤1中计算得到的该迎角时的气动载荷加载到机翼上,利用NASTRAN软件计算得到翼尖的扭转角;用CAD软件重新几何造形以使翼尖扭转相同的角度且使扭转角从翼根到翼尖线性变化;再按照步骤1中的方式重新划分CFD计算网格并实施常规CFD数值计算,提取机翼变形后各翼面网格节点上的压力分布,得到机翼的变形刚体气动载荷数据;
步骤3,对升力面的弦平面进行面元网格参数初始化;
采用通用的三次多项式插值法,通过展向三次曲线对机翼展向进行网格点布置;通过弦向三次曲线对机翼弦向进行网格点布置;
对所述展向三次曲线和弦向三次曲线分别进行归一化;以给定的展向网格点数量范围、弦向网格点数量范围和每条三次曲线两端的斜率值作为约束条件,需初始化的变量为弦向和展向的网格点数量;
步骤4,用几何算法对机翼弦平面进行面元网格分划;
根据给定的三次曲线参数及展向网格点数量和弦向网格点数量对机翼弦平面进行面元网格节点布置;
步骤5,面元法计算气动力矩阵;
通过ZAERO软件生成气动力矩阵;具体是,将步骤4得到的网格分布输入到ZAERO软件,计算给定飞行状态下的气动力影响系数矩阵[AIC];上述飞行状态与步骤1中计算CFD数据的飞行状态一致;
步骤6,采用分段斜率修正面元法对气动载荷分布进行修正,获得修正后的载荷分布;
利用步骤1得到的定常刚体CFD计算数据对面元法的气动力影响系数矩阵[AIC]实施修正,具体是采用公式(5)进行分段线化高精度面元修正:
式中,n为迎角范围划分的子区间数目,[WT1]i为第i个迎角区间采用斜率修正法得到的修正因子矩阵,[Cor]i为各面元在第i个迎角区间内的下洗速度值占总下洗速度的比例系数对角矩阵;{W}为无量纲下洗列阵;{Cp}为各面元上的压差系数列阵;
步骤7,计算CFD数据与修正面元法结果的误差范数,判断是否达到收敛条件;
若是,则得到最优网格分布,迭代过程结束;反之采用自适应模拟退火算法更新网格分布参数,更新后转至步骤4,重复步骤4~步骤6;在重复步骤4~步骤6的过程中不断获得新的修正后的载荷分布;对得到的新的修正后的载荷分布通过本步骤计算CFD数据与修正面元法结果的误差范数,判断是否达到收敛条件,直至满足收敛条件。
2.如权利要求1所述面元修正与网格预先自适应计算方法,其特征在于,在确定不同迎角下机翼部件原刚体气动载荷数据时,翼面采用三角形单元划分,最小尺度为1mm;附面层用三棱柱网格划分,共20层,附面层之外的空间网格采用四面体单元;从附面层最内层开始网格尺度以1.2倍的体积比率增长;CFD数据由FLUENT计算,采用耦合压力基求解器;空气假设为理想气体,温度为300K,参考压力为海平面大气压,湍流模型采用S-A模型,其余参数选择默认设置。
3.如权利要求1所述面元修正与网格预先自适应计算方法,其特征在于,在对展向布置网格节点时,先对展向三次曲线的横轴按照展向节点数量进行等距划分,此时横轴上各点在三次曲线上对应的纵坐标即表示此节点的展向比例位置,且横轴正方向代表从机翼翼根到翼尖的方向;在对弦向布置网格节点时,先对弦向三次曲线的横轴按照弦向节点数量进行等距划分,此时横轴上各点在三次曲线上对应的纵坐标即表示此节点的弦向比例位置,且横轴正方向代表从机翼前缘到后缘的方向。
4.如权利要求1所述面元修正与网格预先自适应计算方法,其特征在于,对于弦向三次曲线,若需要在机翼前缘加密,则应使曲线在(0,0)点的斜率值减小;对于展向三次曲线,若需要在翼尖加密,则应使曲线在(1,1)点的斜率值减小;若要采用均匀划分形式,则只需要将(0,0)、(1,1)处的斜率值均设为1即可。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201510582632.8A CN105183996B (zh) | 2015-09-14 | 2015-09-14 | 面元修正与网格预先自适应计算方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201510582632.8A CN105183996B (zh) | 2015-09-14 | 2015-09-14 | 面元修正与网格预先自适应计算方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN105183996A CN105183996A (zh) | 2015-12-23 |
CN105183996B true CN105183996B (zh) | 2018-07-27 |
Family
ID=54906075
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201510582632.8A Active CN105183996B (zh) | 2015-09-14 | 2015-09-14 | 面元修正与网格预先自适应计算方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN105183996B (zh) |
Families Citing this family (22)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106168529B (zh) * | 2016-07-27 | 2018-06-26 | 江西洪都航空工业集团有限责任公司 | 一种修正飞行器刚性体升力系数的风洞试验方法 |
CN107038296B (zh) * | 2017-04-06 | 2020-12-18 | 深圳数设科技有限公司 | 一种确定飞机翼面中有限元点载荷分布的方法与装置 |
CN107194079B (zh) * | 2017-05-25 | 2020-12-15 | 中国商用飞机有限责任公司 | 像素细分载荷传递方法和系统 |
CN109117584B (zh) * | 2018-09-05 | 2023-01-13 | 四川腾盾科技有限公司 | 一种低速飞机突风载荷系数计算方法及设备 |
CN109342037B (zh) * | 2018-10-17 | 2021-03-26 | 中国特种飞行器研究所 | 一种框架式翼面结构气动载荷加载系统及加载方法 |
CN109165477A (zh) * | 2018-10-22 | 2019-01-08 | 哈尔滨工程大学 | 一种螺旋桨静强度计算方法 |
CN109446660A (zh) * | 2018-10-31 | 2019-03-08 | 潍柴动力股份有限公司 | 一种发动机进气道自动优化方法及装置 |
CN109598025B (zh) * | 2018-11-02 | 2023-05-26 | 中国航空工业集团公司西安飞机设计研究所 | 一种v形尾翼倾角与面积优化设计方法 |
CN109323841B (zh) * | 2018-11-23 | 2020-03-13 | 中国航空工业集团公司沈阳飞机设计研究所 | 基于网格的机翼总载荷与分布载荷的协调方法 |
CN109299579B (zh) * | 2018-11-23 | 2023-05-23 | 中国航空工业集团公司沈阳飞机设计研究所 | 大展弦比飞机风洞测力试验数据修正方法 |
CN109726437B (zh) * | 2018-12-04 | 2023-05-26 | 中国航空工业集团公司西安飞机设计研究所 | 一种舱门气动载荷等效节点力处理方法 |
CN109684674B (zh) * | 2018-12-04 | 2023-05-05 | 中国航空工业集团公司西安飞机设计研究所 | 一种舱门气动载荷处理方法 |
CN110704953B (zh) * | 2019-09-30 | 2020-08-14 | 西北工业大学 | 一种大展弦比机翼静气弹性能设计敏度的分析方法 |
CN110763424B (zh) * | 2019-10-31 | 2021-08-20 | 中电工业互联网有限公司 | 一种机翼表面压力测量方法、系统及装置 |
CN112446067B (zh) * | 2020-11-03 | 2022-12-02 | 中国空气动力研究与发展中心计算空气动力研究所 | 一种基于弹性变形的区域网格动态重构方法 |
CN112347561B (zh) * | 2020-11-10 | 2024-05-28 | 中国商用飞机有限责任公司北京民用飞机技术研究中心 | 飞行器的静气动弹性分析方法、装置、设备和存储介质 |
CN112464372B (zh) * | 2020-11-25 | 2021-08-27 | 西北工业大学 | 一种弹性机翼副翼舵面效率的设计敏度工程数值方法 |
CN113221406B (zh) * | 2021-04-30 | 2022-04-19 | 中国电建集团华东勘测设计研究院有限公司 | 极大压缩变形条件下模型网格自适应修正方法 |
CN114580196B (zh) * | 2022-03-18 | 2022-08-16 | 中国科学院力学研究所 | 一种切割自由面水翼空泡长度预测方法 |
CN116187133B (zh) * | 2023-02-10 | 2023-10-27 | 西北工业大学深圳研究院 | 一种用于旋压移动网格法的维度分离弹簧比拟法 |
CN116225073B (zh) * | 2023-05-09 | 2023-07-07 | 中国空气动力研究与发展中心计算空气动力研究所 | 一种基于cfd方法的僚机飞行迎角快速确定方法 |
CN117540613B (zh) * | 2024-01-10 | 2024-04-05 | 西安羚控电子科技有限公司 | 一种气动载荷二次插值处理方法及装置 |
Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
EP1661805A1 (de) * | 2004-11-24 | 2006-05-31 | Airbus Deutschland GmbH | Deckhaut für eine formvariable aerodynamische Fläche |
EP2416266A2 (en) * | 2009-03-31 | 2012-02-08 | Airbus Operations, S.L. | Computer-assisted method for the cost-optimised calculation of variable distributions on an aerodynamic surface |
CN102745324A (zh) * | 2012-07-07 | 2012-10-24 | 西北工业大学 | 一种柔性蒙皮 |
CN103646131A (zh) * | 2013-11-26 | 2014-03-19 | 北京航空航天大学 | 一种考虑气动弹性约束的复合材料机翼多目标优化设计方法 |
-
2015
- 2015-09-14 CN CN201510582632.8A patent/CN105183996B/zh active Active
Patent Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
EP1661805A1 (de) * | 2004-11-24 | 2006-05-31 | Airbus Deutschland GmbH | Deckhaut für eine formvariable aerodynamische Fläche |
EP2416266A2 (en) * | 2009-03-31 | 2012-02-08 | Airbus Operations, S.L. | Computer-assisted method for the cost-optimised calculation of variable distributions on an aerodynamic surface |
CN102745324A (zh) * | 2012-07-07 | 2012-10-24 | 西北工业大学 | 一种柔性蒙皮 |
CN103646131A (zh) * | 2013-11-26 | 2014-03-19 | 北京航空航天大学 | 一种考虑气动弹性约束的复合材料机翼多目标优化设计方法 |
Non-Patent Citations (4)
Title |
---|
A Comparison of Two Dynamic Mesh Methods in Fluid–Structure Interaction;JIA Huan 等;《2nd International Conference on Electronic & Mechanical Engineering and Information Technology》;20121231;第2171-2175页 * |
Alberto Varello 等.Vortex Lattice Method coupled with advanced one-dimensional structural models.《Journal of Aeroelasticity & Structural Dynamics》.2011,第2卷(第2期),第53-78页. * |
基于修正面元法的机翼焦点位置变化研究;贾欢 等;《航空工程进展》;20150531;第6卷(第2期);第166-170、188页 * |
贾欢 等.基于修正面元法的机翼静气动弹性计算.《计算物理》.2014,第31卷(第1期),第21-26页. * |
Also Published As
Publication number | Publication date |
---|---|
CN105183996A (zh) | 2015-12-23 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN105183996B (zh) | 面元修正与网格预先自适应计算方法 | |
Potsdam et al. | Rotor airloads prediction using loose aerodynamic/structural coupling | |
Abdessemed et al. | Morphing airfoils analysis using dynamic meshing | |
CN107391891A (zh) | 一种基于模型融合方法的大展弦比机翼优化设计方法 | |
Zhang et al. | Improvement to patched grid technique with high-order conservative remapping method | |
Chaderjian et al. | Detached eddy simulation of the UH-60 rotor wake using adaptive mesh refinement | |
Sleesongsom et al. | Aircraft morphing wing design by using partial topology optimization | |
CN111274648B (zh) | 一种民用飞机前缘襟翼的分布式飞行载荷设计方法 | |
Leifsson et al. | Fast low-fidelity wing aerodynamics model for surrogate-based shape optimization | |
Koning et al. | Using RotCFD to Predict Isolated XV-15 Rotor Performance | |
Lee-Rausch et al. | Application of a full reynolds stress model to high lift flows | |
Eller | Fast, unstructured-mesh finite-element method for nonlinear subsonic flow | |
Nikkhoo et al. | Effect of different aero-structural optimization in the commercial airplane | |
CN110826182B (zh) | 基于顶点法和序贯优化策略的飞行器结构气动弹性设计法 | |
Keye et al. | Development of Deformed CAD Geometries of NASA's Common Research Model for the Sixth AIAA CFD Drag Prediction Workshop | |
Ahmad et al. | High-order accurate CFD/CSD simulation of the UH60 rotor in forward flight | |
Tran et al. | Transition Modeling in CREATETM-AV HELIOS v9 | |
Abdol-Hamid et al. | DPW-6 Results Using FUN3D With Focus on k-kL-MEAH2015 Turbulence Model | |
Juvigny et al. | Multigrid algorithms for the chimera method | |
Ghoreyshi et al. | Grid Quality and Resolution Effects on the Aerodynamic Modeling of Parachute Canopies | |
Liguori et al. | Aircraft fuselage effects on transonic wing pressures via Non-Linear Vortex Lattice Method | |
Ma et al. | VPM/CFD-based research on rotor performance and loads of individual blade control rotor system | |
Stickan et al. | Coupled CFD-CSM simulations of the Aerostabil windtunnel experiments, considering sophisticated structural modeling | |
Eisfeld | Numerical simulation of aerodynamic problems with the SSG/LRR-ω Reynolds stress turbulence model using the unstructured TAU code | |
Hristov et al. | Implementation of a Viscous Decambering Routine to Low-Order Aircraft Modeling for Stability and Control Prediction |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |