CN106844896B - 一种适用于旋成体外形的来流参数确定方法 - Google Patents

一种适用于旋成体外形的来流参数确定方法 Download PDF

Info

Publication number
CN106844896B
CN106844896B CN201611260976.8A CN201611260976A CN106844896B CN 106844896 B CN106844896 B CN 106844896B CN 201611260976 A CN201611260976 A CN 201611260976A CN 106844896 B CN106844896 B CN 106844896B
Authority
CN
China
Prior art keywords
pressure
polynomial
incoming flow
measuring point
revolution
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
Application number
CN201611260976.8A
Other languages
English (en)
Other versions
CN106844896A (zh
Inventor
豆国辉
陈广强
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
China Academy of Aerospace Aerodynamics CAAA
Original Assignee
China Academy of Aerospace Aerodynamics CAAA
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by China Academy of Aerospace Aerodynamics CAAA filed Critical China Academy of Aerospace Aerodynamics CAAA
Priority to CN201611260976.8A priority Critical patent/CN106844896B/zh
Publication of CN106844896A publication Critical patent/CN106844896A/zh
Application granted granted Critical
Publication of CN106844896B publication Critical patent/CN106844896B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Aerodynamic Tests, Hydrodynamic Tests, Wind Tunnels, And Water Tanks (AREA)

Abstract

一种适用于旋成体外形的来流参数确定方法，(1)建立适用于旋成体外形的表面压力近似模型：首先，确定表面压力计算公式，将表面压力表示为来流动压qc与压力系数Cpi的乘积与来流静压p∞的加和；然后，采用多项式的形式表达所述的压力系数Cpi，多项式因子为飞行攻角α、侧滑角β和来流压力比R；最后，通过多组状态下的旋成体表面测点压力数值拟合或者回归或者最小二乘法，得到上述多项式中的系数；(2)获取飞行试验中的旋成体表面测点压力，根据表面测点压力结合上述近似模型进行反算，得到多项式因子，进而得到来流参数。本发明模型可以用于嵌入式大气数据系统，可以有效提高系统的预测精度。

Description

一种适用于旋成体外形的来流参数确定方法
技术领域
本发明涉及一种适用于旋成体外形的表面压力近似模型，可以应用于火箭弹等嵌入式大气数据系统(FADS)当中，属于飞行器飞行状态测量领域。
背景技术
鉴于高超声速的高温加热效应，嵌入式大气数据系统(FADS)在高超声速飞行器上是必要的，高精度一直是大气数据系统追求的目标，除使用神经网络之间建立从测压点压力到飞行参数之间的函数外，一般使用理论方法建立从飞行参数到表面压力的模型，使用时是通过反求不同测压点压力构成的超定方程组获得对应的飞行参数。
理论模型下的嵌入式大气数据系统的精度直接与表面压力模型相关，压力模型的精度越高，大气数据系统的精度也高。最常用的是基于球面建立的牛顿模型和对应的修正模型，这些模型对大钝头体外形的飞行器有较大的应用价值，但对于旋成体类型的导弹类飞行器误差较大。目前针对旋成体一般采用基于神经网络的大气数据系统，但是该系统很难保证在扰动的测量精度。
发明内容
本发明的技术解决问题：针对旋成体类型飞行器的特点，建立了一种基于高精度的表面压力近似模型的来流参数确定方法，模型可以用于嵌入式大气数据系统，可以有效提高系统的预测精度。
本发明的技术方案：一种适用于旋成体外形的来流参数确定方法，通过下述方式实现：
(1)建立适用于旋成体外形的表面压力近似模型：
首先，确定表面压力计算公式，将表面压力表示为来流动压qc与压力系数Cpi的乘积与来流静压p∞的加和；
然后，采用多项式的形式表达所述的压力系数Cpi，多项式因子为飞行攻角α、侧滑角β和来流压力比R；
最后，通过每个表面测点在多组状态下的压力数值进行拟合或者回归或者最小二乘法，得到上述多项式中的系数，得到每个表面测点对应的表面压力近似模型；
(2)获取飞行试验中的旋成体表面测点压力，根据表面测点压力结合上述近似模型进行反算，得到多项式因子，进而得到来流参数。
进一步的，多项式的阶数至少4阶。
进一步的，当表面测点位于旋成体的顶点、机体纵向和横向平面内时，所述的多项式表达形式简化对应的侧滑角和或攻角的奇次项，只保留侧滑角和或攻角其偶次项。
进一步的，当多项式为5阶时，多项式形式如下：
f(α,β,R)＝a0+a1α+a2β+a3R+a4α2+a5β2+a6R2+a7αβ+a8αR+a9βR+a10α3+a11β3
+a12R3+a13α2β+a14αβ2+a15α2R+a16αR2+a17β2R+a18βR2+a19αβR+a20α4
+a21β4+a22R4+a23α3β+a24α2β2+a25αβ3+a26α3R+a27α2R2+a28αR3+a29β3R
+a30β2R2+a31βR3+a32α2βR+a33αβ2R+a34αβR2+a35α5+a36β5+a36R5+a37α4β
+a38α3β2+a39α2β3+a40αβ4+a41α4R+a42α3R2+a43α2R3+a44αR4+a45β4R
+a46β3R2+a46β2R3+a47βR4+a48α3βR+a49αβ3R+a50αβR3+a51α2β2R+a52α2βR2
+a53αβ2R2
式中，a0～a53为多项式中的系数。
当表面测点位于旋成体顶点时，多项式形式如下：
f1(α,β,R)＝a0+a3R+a4α2+a5β2+a6R2+a12R3+a15α2R+a17β2R+a20α4+a21β4+a22R4
+a24α2β2+a27α2R2+a30β2R2+a36R5+a41α4R+a43α2R3+a45β4R+a46β2R3
+a51α2β2R
当表面测点位于旋成体的纵向平面时，多项式形式如下：
f2(α,β,R)＝a0+a1α+a3R+a4α2+a5β2+a6R2+a8αR+a10α3+a12R3+a14αβ2+a15α2R
+a16αR2+a17β2R+a20α4+a21β4+a22R4+a24α2β2+a26α3R+a27α2R2+a28αR3
+a30β2R2+a33αβ2R+a35α5+a36R5+a38α3β2+a40αβ4+a41α4R+a42α3R2+a43α2R3
+a44αR4+a45β4R+a46β2R3+a51α2β2R+a53αβ2R2
当表面测点位于旋成体的横向平面时，多项式形式如下：
f3(α,β,R)＝a0+a2β+a3R+a4α2+a5β2+a6R2+a9βR+a11β3+a12R3+a13α2β+a15α2R
+a17β2R+a18βR2+a20α4+a21β4+a22R4+a24α2β2+a27α2R2+a29β3R+a30β2R2
+a31βR3+a32α2βR+a36β5+a36R5+a37α4β+a39α2β3+a41α4R+a43α2R3+a45β4R
+a46β3R2+a46β2R3+a47βR4+a51α2β2R+a52α2βR2
进一步的，表面测点包含顶点、沿旋成体轴向至少分布2层测点，纵向和横向平面内应至少包含3个测点。
进一步的，反算过程中至少需要3个表面测点压力值。
进一步的，工程上，反算过程中至少需要7个表面测点压力值。
本发明与现有技术相比的优点如下：
(1)现有的表面压力表达式大多是基于球头外形建立的牛顿模型或其改进形式，应用于旋成体外形时表面测压点压力误差较大，本发明中的多项式形式的压力表达式可以应用于旋成体外形并满足精度要求。
(2)相比于基于神经网络的嵌入式大气数据系统，系统扰动精度难以保证的问题，本发明通过分析来流状态和测点位置对测压点压力的影响，难以将测点位置通过余弦函数或其他函数来建立压力表达式，本发明仅考虑来流状态对各测压点建立不同的表达式形式，在多项式阶次高于4阶时能获得较好的拟合精度，在多形式为5阶时基本能满足需求，建立的测点压力表达式可以较容易描述状态参数对测点压力的影响，分析参数的误差影响，进而保证大气数据系统的测量稳定性。
(3)相比于常用的压力表达式，本发明选择除选择攻角和侧滑角外，还使用压力比作为压力表达式的因子而不是使用马赫数，这样的选择更便于获得高精度的拟合结果。
(4)本发明经过了算例试验，并用于设计嵌入式大气数据系统，系统预测精度较高。
附图说明
图1为本发明方法流程图；
图2为某旋成体外形及表面测压点示意图。
具体实施方式
下面结合附图及实例对本发明作详细说明。本发明涉及一种适用于旋成体外形的来流参数确定方法，如图1所示，通过下述方式实现：
(1)建立适用于旋成体外形的表面压力近似模型：
首先，确定表面压力计算公式，将表面压力表示为来流动压qc与压力系数Cpi的乘积与来流静压p∞的加和，即pi＝qcCpi+p∞；
然后，采用多项式的形式表达所述的压力系数Cpi，Cpi＝f(α,β,R)；
压力系数采用多项式的形式，至少需要采用4阶进行表达；多项式因子为飞行攻角α、侧滑角β和来流压力比R＝p∞/(p∞+qc)；
本例中给出5阶的表达形式，由54项参变量乘积或常数项构成。
f(α,β,R)＝a0+a1α+a2β+a3R+a4α2+a5β2+a6R2+a7αβ+a8αR+a9βR+a10α3+a11β3
+a12R3+a13α2β+a14αβ2+a15α2R+a16αR2+a17β2R+a18βR2+a19αβR+a20α4
+a21β4+a22R4+a23α3β+a24α2β2+a25αβ3+a26α3R+a27α2R2+a28αR3+a29β3R
+a30β2R2+a31βR3+a32α2βR+a33αβ2R+a34αβR2+a35α5+a36β5+a36R5+a37α4β
+a38α3β2+a39α2β3+a40αβ4+a41α4R+a42α3R2+a43α2R3+a44αR4+a45β4R
+a46β3R2+a46β2R3+a47βR4+a48α3βR+a49αβ3R+a50αβR3+a51α2β2R+a52α2βR2
+a53αβ2R2
对于飞行器顶点(旋成体顶点)，考虑攻角和侧滑角的对称性，压力系数的函数可以去掉攻角和侧滑角的奇次项获得简化的函数形式，由20项构成。
f1(α,β,R)＝a0+a3R+a4α2+a5β2+a6R2+a12R3+a15α2R+a17β2R+a20α4+a21β4+a22R4
+a24α2β2+a27α2R2+a30β2R2+a36R5+a41α4R+a43α2R3+a45β4R+a46β2R3
+a51α2β2R
对于飞行器纵向平面内的物面点，考虑侧滑角的对称性，压力系数的函
数可以去掉侧滑角的奇次项获得简化的函数形式，由34项构成。
f2(α,β,R)＝a0+a1α+a3R+a4α2+a5β2+a6R2+a8αR+a10α3+a12R3+a14αβ2+a15α2R
+a16αR2+a17β2R+a20α4+a21β4+a22R4+a24α2β2+a26α3R+a27α2R2+a28αR3
+a30β2R2+a33αβ2R+a35α5+a36R5+a38α3β2+a40αβ4+a41α4R+a42α3R2+a43α2R3
+a44αR4+a45β4R+a46β2R3+a51α2β2R+a53αβ2R2
对于飞行器横向平面内的物面点，考虑攻角的对称性，压力系数的函数可以去掉攻角的奇次项获得简化的函数形式，由34项构成。
f3(α,β,R)＝a0+a2β+a3R+a4α2+a5β2+a6R2+a9βR+a11β3+a12R3+a13α2β+a15α2R
+a17β2R+a18βR2+a20α4+a21β4+a22R4+a24α2β2+a27α2R2+a29β3R+a30β2R2
+a31βR3+a32α2βR+a36β5+a36R5+a37α4β+a39α2β3+a41α4R+a43α2R3+a45β4R
+a46β3R2+a46β2R3+a47βR4+a51α2β2R+a52α2βR2
最后，通过每个表面测点在多组状态下的压力数值进行拟合或者回归或者最小二乘法，得到上述多项式中的系数ai，得到每个表面测点对应的表面压力近似模型；
对于旋成体表面测点的分布理论上没有特殊的要求，应用中确定包含顶点、沿旋成体轴向至少分布2层测点，纵向和横向平面内应至少包含3个测点的布局能够提高测量精度。
在使用中需要通过风洞试验数据或数值计算的表面测点压力数据来获得模型中的多项式系数。通过不同飞行状态下的各测压点压力数据采用线性回归方法获得各测压点的压力多项式的系数，由于方法的限制，飞行状态的个数应大于压力系数多项式的系数的个数。
(2)获取飞行试验中的旋成体表面测点压力，根据表面测点压力结合上述近似模型进行反算，得到多项式因子，进而得到来流参数。来流参数包括攻角、侧滑角、马赫数、来流静压等数据，这些参数可以应用于飞行器的飞行控制系统中(例如作为控制增益表的输入参数)。
理论上，上述反算过程中至少需要3个表面测点压力值。工程上，反算过程中至少需要7个表面测点压力值来减小单点测量误差对最终结果的影响，提高系统的鲁棒性；关于反算可以采用目前常用批处理滤波器方法来实现。
实施例：
如图2所示，某旋成体飞行器，头部圆头处理，半径5mm，紧接着为卡门曲线，之后直线段，其锥角为88.3度，这样的曲线旋转获得旋成体外形。在卡门曲线和直线段构成的曲面上分别取4个点，这些点都位于机体的纵向平面(点3、5、7、9)和横向平面(点2、4、6、8)内，和顶点1一起总共9个点，作为测压点。
本例中，使用数值计算的方法获得1000组飞行状态的测压点压力数据，考虑对称性可以获得另外3000组飞行状态下的压力数据，分析这些测压点数据，分别使用线性回归方法获得这些测压点压力模型的系数。
利用本发明近似模型计算得到的测压点压力值与CFD计算的结果进行对比，按照正态分布统计百分比相对误差，其对应的3σ值如下表所示。
测压点 1 2 3 4 5 6 7 8 9
相对误差(％) 0.313 0.471 0.533 0.471 0.533 2.019 2.076 2.019 2.076
从表中可以看出，算例中各测压点相对误差都很小，近似模型的精度比较高，满足模型的设计需要。
综上所述，本发明建立的旋成体外形的表面压力近似模型有效，模型结构简单，线性回归方法可以获得最优的模型多项式系数，五阶多项式的压力模型的精度满足嵌入式大气数据系统(FADS)中对压力模型的设计需要。
本发明未详细说明部分属本领域技术人员公知常识。

Claims (4)

1.一种适用于旋成体外形的来流参数确定方法，其特征在于通过下述方式实现：
(1)建立适用于旋成体外形的表面压力近似模型：
首先，确定表面压力计算公式，将表面压力表示为来流动压qc与压力系数Cpi的乘积与来流静压p∞的加和；
然后，采用多项式的形式表达所述的压力系数Cpi，多项式因子为飞行攻角α、侧滑角β和来流压力比R；
最后，通过每个表面测点在多组状态下的压力数值进行拟合或者回归或者最小二乘法，得到上述多项式中的系数，得到每个表面测点对应的表面压力近似模型；
(2)获取飞行试验中的旋成体表面测点压力，根据表面测点压力结合上述近似模型进行反算，得到多项式因子，进而得到来流参数；
多项式的阶数至少4阶；
当表面测点位于旋成体的顶点、机体纵向和横向平面内时，所述的多项式表达形式简化对应的侧滑角和或攻角的奇次项，只保留侧滑角和或攻角其偶次项；
当多项式为5阶时，多项式形式如下：
f(α,β,R)＝a0+a1α+a2β+a3R+a4α2+a5β2+a6R2+a7αβ+a8αR+a9βR+a10α3+a11β3+a12R3+a13α2β+a14αβ2+a15α2R+a16αR2+a17β2R+a18βR2+a19αβR+a20α4+a21β4+a22R4+a23α3β+a24α2β2+a25αβ3+a26α3R+a27α2R2+a28αR3+a29β3R+a30β2R2+a31βR3+a32α2βR+a33αβ2R+a34αβR2+a35α5+a36β5+a36R5+a37α4β+a38α3β2+a39α2β3+a40αβ4+a41α4R+a42α3R2+a43α2R3+a44αR4+a45β4R+a46β3R2+a46β2R3+a47βR4+a48α3βR+a49αβ3R+a50αβR3+a51α2β2R+a52α2βR2+a53αβ2R2
式中，a0～a53为多项式中的系数；
当表面测点位于旋成体顶点时，多项式形式如下：
f1(α,β,R)＝a0+a3R+a4α2+a5β2+a6R2+a12R3+a15α2R+a17β2R+a20α4+a21β4+a22R4+a24α2β2+a27α2R2+a30β2R2+a36R5+a41α4R+a43α2R3+a45β4R+a46β2R3+a51α2β2R
当表面测点位于旋成体的纵向平面时，多项式形式如下：
f2(α,β,R)＝a0+a1α+a3R+a4α2+a5β2+a6R2+a8αR+a10α3+a12R3+a14αβ2+a15α2R+a16αR2+a17β2R+a20α4+a21β4+a22R4+a24α2β2+a26α3R+a27α2R2+a28αR3+a30β2R2+a33αβ2R+a35α5+a36R5+a38α3β2+a40αβ4+a41α4R+a42α3R2+a43α2R3+a44αR4+a45β4R+a46β2R3+a51α2β2R+a53αβ2R2
当表面测点位于旋成体的横向平面时，多项式形式如下：
f3(α,β,R)＝a0+a2β+a3R+a4α2+a5β2+a6R2+a9βR+a11β3+a12R3+a13α2β+a15α2R+a17β2R+a18βR2+a20α4+a21β4+a22R4+a24α2β2+a27α2R2+a29β3R+a30β2R2+a31βR3+a32α2βR+a36β5+a36R5+a37α4β+a39α2β3+a41α4R+a43α2R3+a45β4R+a46β3R2+a46β2R3+a47βR4+a51α2β2R+a52α2βR2。
2.根据权利要求1所述的方法，其特征在于：表面测点包含顶点、沿旋成体轴向至少分布2层测点，纵向和横向平面内应至少包含3个测点。
3.根据权利要求1所述的方法，其特征在于：反算过程中至少需要3个表面测点压力值。
4.根据权利要求3所述的方法，其特征在于：工程上，反算过程中至少需要7个表面测点压力值。
CN201611260976.8A 2016-12-30 2016-12-30 一种适用于旋成体外形的来流参数确定方法 Active CN106844896B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201611260976.8A CN106844896B (zh) 2016-12-30 2016-12-30 一种适用于旋成体外形的来流参数确定方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201611260976.8A CN106844896B (zh) 2016-12-30 2016-12-30 一种适用于旋成体外形的来流参数确定方法

Publications (2)

Publication Number Publication Date
CN106844896A CN106844896A (zh) 2017-06-13
CN106844896B true CN106844896B (zh) 2020-07-14

Family

ID=59113706

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201611260976.8A Active CN106844896B (zh) 2016-12-30 2016-12-30 一种适用于旋成体外形的来流参数确定方法

Country Status (1)

Country Link
CN (1) CN106844896B (zh)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109033548B (zh) * 2018-07-03 2020-07-07 中国空气动力研究与发展中心高速空气动力研究所 一种计算槽壁边界条件主要系数的拟合方法
CN110046473B (zh) * 2019-05-27 2023-04-18 中国空气动力研究与发展中心高速空气动力研究所 一种飞行器大气参数解算方法、装置及计算机设备
CN111122106A (zh) * 2019-12-19 2020-05-08 中国航空工业集团公司西安飞机设计研究所 一种全压受感器测量误差的修正方法和计算机设备

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105184109A (zh) * 2015-10-27 2015-12-23 中国人民解放军国防科学技术大学 扰动引力作用下弹道助推段状态偏差解析计算方法
CN105628086A (zh) * 2014-10-29 2016-06-01 北京临近空间飞行器系统工程研究所 一种基于锥面压力分布的超声速飞行来流参数解算方法

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105628086A (zh) * 2014-10-29 2016-06-01 北京临近空间飞行器系统工程研究所 一种基于锥面压力分布的超声速飞行来流参数解算方法
CN105184109A (zh) * 2015-10-27 2015-12-23 中国人民解放军国防科学技术大学 扰动引力作用下弹道助推段状态偏差解析计算方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
端壁抽吸位置对大转角扩压叶栅流场及负荷的影响;郭爽等;《推进技术》;20110615;第32卷(第3期);说明书第327页 *

Also Published As

Publication number Publication date
CN106844896A (zh) 2017-06-13

Similar Documents

Publication Publication Date Title
CN106844896B (zh) 一种适用于旋成体外形的来流参数确定方法
CN104697462B (zh) 一种基于中轴线的航空叶片型面特征参数提取方法
CN113505434B (zh) 基于气动力数学模型的飞行器设计制造方法及其飞行器
CN110750855B (zh) 一种外形定尺寸限制下的蜗壳型线设计方法
CN106005475A (zh) 高超声速内外流一体化全乘波飞行器设计方法
CN112069602B (zh) 一种喷管入口总温的反向重构方法、装置、介质及设备
CN114112286B (zh) 一种高超声速风洞轴对称型面喷管拟合喉道段设计方法
CN115358026A (zh) 一种基于多元线性回归与曲面拟合的五孔探针数据处理方法
CN114323536A (zh) 一种提高五孔探针测量精度的插值方法
CN106844966B (zh) 一种螺旋桨叶面叶背精确建模方法
CN110414168B (zh) 基于与前机身耦合优化的高超声速隔离段设计方法及系统
CN114781078B (zh) 一种基于矩阵变换的隐身蛇形进气道设计方法
CN104330068A (zh) 一种降低叶片型面三坐标测量补偿误差的方法
CN111780949A (zh) 基于cfd分析的高速进气道前体风洞实验总压修正方法
CN109325279B (zh) 一种离散的飞机气动载荷参数化的方法
CN105718619A (zh) 一种基于有限元法的飞行器燃油质量特性确定方法
CN114611420A (zh) 非定常气动力计算精度评估及修正方法
CN105373637A (zh) 一种波瓣混合器流量系数的数值计算方法
CN112163271B (zh) 大气数据传感系统的大气参数解算方法
CN113188799A (zh) 基于速度差极值法的航空发动机推力修正方法
CN105760647B (zh) 一种质量加权平均值计算方法
CN108536932B (zh) 基于互扭约束条件下的航空叶片积叠轴垂直度计算方法
CN111693251B (zh) 一种桨对舵的水动力干扰系数测定方法
CN113281001B (zh) 基于集成式微型大气数据模块的全速域大气数据解算方法
CN114571732A (zh) 一种适用于复杂曲面类零件的3d打印分层方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant