CN115659856A - 一种飞行器的气动导数预测方法 - Google Patents

Abstract

本发明涉及一种飞行器的气动导数预测方法，属于飞行器技术领域。本发明基于涡格法建立了飞行器机翼模型；基于建立的飞行器机翼模型给出状态空间方程的定常形式；然后基于状态空间方程的定常形式获得飞行器气动力和气动力矩关于飞行器运动变量的导数，并获得飞行器气动力和气动力矩关于飞行器舵偏的导数。本发明的气动导数预测方法，解决了现有的飞行器的气动导数预测方法预测过程复杂、成本高且周期长，而且预测结果不精确的问题。

Description

一种飞行器的气动导数预测方法

技术领域

本发明属于飞行器技术领域，具体涉及一种飞行器的气动导数预测方法。

背景技术

在设计飞行器时，基于气动导数构建飞行器的气动力模型是目前常用方法。传统预测飞行器的气动导数的方法主要有飞行试验、风洞试验、工程经验方法及CFD数值 计算方法。飞行试验方法以飞行器实际飞行过程中获得的气动数据为基础，通过插值、 拟合或参数辨识等手段求解气动导数，该方法难度大、周期长、风险高、成本高，难 以在飞行器设计初期阶段获得指导性的数据。风洞试验方法以飞行器在风洞试验中获 得的气动数据为基础，通过插值、拟合或参数辨识等手段求解气动导数，该方法也具 有成本高，周期长的缺陷，难以应用于飞行器设计初期阶段。工程经验方法通过简单 公式结合飞行器气动外形参数计算气动导数，计算快捷但是精度受到很大限制。CFD 数值计算方法计算精度优于工程经验方法但计算成本较高，为获取气动导数需要进行 大量仿真限制了设计分析效率。

发明内容

鉴于上述分析，本发明提供了一种飞行器的气动导数预测方法，用以解决现有的飞行器的气动导数预测方法预测过程复杂、成本高且周期长，而且预测结果不精确的 问题。

本发明提供了一种飞行器的气动导数预测方法，包括以下步骤：

S1、基于涡格法建立气动模型：

在飞行器机翼的弧面沿弦向和展向划分多个干四边形气动网格，气动网格包括附着于机翼表面的附着涡网格及沿来流方向拖出的机翼后缘的第一排尾涡网格以及其 他尾涡网格；

S2、对步骤S1得到的气动网格建立状态空间形式的气动力方程：

S2-1、在每个气动网格1/4弦线处布置涡线段，四段等强度涡线段中中的相邻两段涡线段首尾相连，构成一个涡环；选择1/4弦线的中点为气动网格气动力的作用点， 选择3/4弦线中点为气动网格的控制点；

S2-2、将气动模型中的附着涡网格、第一排尾涡网格以及其他尾涡网格对应为：机翼表面的附着涡，机翼后缘的第一排尾涡以及其他尾涡，设机翼表面的附着的涡强 度为Γ_b，机翼后缘的第一排尾涡的强度为Γ_w0，其他尾涡的强度为Γ_wl，气动控制方 程为：

K_bΓ_b+K_w0Γ_w0+K_wlΓ_wl＝-w (29)

其中，w为机翼表面诱导速度，w＝(V_∞+V_g)·n；V_∞为来流速度矢量，V_g为来 流扰动速度矢量，n为控制点处的法向量列阵。K_b,K_w0,K_wl分别为附着涡，第一排 尾涡及其他尾涡的诱导系数矩阵；

S2-3、机翼后缘的第一排尾涡在脱出的过程中保持涡强守恒，关系表达为：

其中，

为机翼后缘第一排尾涡的强度关于时间的导数，Δt为时间步长，C₁为 保证后缘第一排的尾涡与附着涡对应关系正确的系数矩阵，只包含0和1两种元素；

S2-4、规定其他尾涡脱出后保持涡强不变，关系表达为：

其中，

为其他尾涡的强度关于时间的导数，C₂,C₃为表征其他涡与机翼后缘 第一排尾涡位置对应关系的常数提取矩阵，只包含0和1两种元素；

S2-5、综合式(1)、(2)和(3)得到涡格法状态空间方程形式为：

其中，尾涡的强度Γ_w＝[Γ_w0 Γ_wl]^T，A_a,B_a为状态空间系数矩阵，与机翼及尾 涡的几何形状及机翼表面的附着涡、机翼后缘的第一排尾涡以及其他尾涡三者的划分 有关：

其中，I为单位阵，O为只包含0元素的零矩阵。

S2-6、给出气动网格上的气动力表达，气动网格的气动力表达为：

其中，V_∞为来流速度矢量，e_Γ为前缘涡切向矢量，ρ为空气密度，

为来流单 位矢量，

为前缘涡切向单位矢量，A为气动网格面积；Γ_i为第i个气动网格中的有 效涡量，数值为沿飞行器飞行方向的前后位置的气动网格的物面涡的涡强之差；

为 第i个气动网格上有效涡量时间导数，数值为沿飞行器飞行方向的前后位置的气动网 格的物面涡的涡强时间导数的平均值。

S2-7、综合式(2)、(3)、(4)和(7)得到气动网格的气动力的矢量F表达式为：

其中，D_a,

C_a为系数矩阵，

为w关于时间的导数。

S3、给出状态空间形式的气动网格的气动力方程的定常形式：

在步骤S2给出的式(4)及式(8)基础上，忽略涡格法状态空间形式中关于涡 强变化的非定常项，状态空间形式的气动网格的气动力方程的定常形式为：

A_aΓ_w+B_aw＝0 (37)

F＝C_aΓ_w+D_aw (38)

其中，在给定来流工况时，系数矩阵A_a,B_a,C_a,D_a与机翼及尾涡的几何形状有 关。

S4、求飞行器全机的气动力和气动力矩关于飞行器运动变量的导数：

S4-1、仅考虑飞行器的刚体运动，飞行器的运动速度V为：

v＝[V₀+ω×r₀] (39)

其中，V₀为飞行器的平飞速度矢量，ω为飞行器的转动角速度矢量，r₀为飞行 器的距重心位置矢量。

S4-2、第i个气动网格的控制点的表面诱导速度w_i为：

w_i＝v_i·n_0i＝(V_0in_0i+(ω×r_0i)n_0i) (40)

其中，v_i为第i个气动网格的运动速度，V_0i为第i个气动网格平飞运动速度,n_0i为第i个气动网格的法向量，r_0i为第i个气动网格到飞行器重心的位置矢量。

S4-3、由式(11)和(12)可得气动网格的气动力为：

其中，G为气动网格的气动力与飞行器机翼结构的插值矩阵。

S4-4、飞行器全机气动力的合外力F_Ar与合外力矩M_Ar通过全部气动网格的气动 力求和得到：

F_Ar＝Φ_tF (42)

M_Ar＝Φ_rF (43)

具体地，求和矩阵Φ_t和Φ_r为：

其中，N_s为气动网格总数，i＝1,2,3…,N_S，

为第i个气动网格气动力的作用点 到取矩点的张量；x_i为第i个气动网格气动力的作用点到取矩点x方向的距离，y_i为 第i个气动网格气动力的作用点到取矩点y方向的距离，z_i为第i个气动网格气动力 的作用点到取矩点z方向的距离，I_3×3为3阶单位矩阵，I_3×3中的取值3×3与气动力 向量的x方向力、y方向力和z方向力相适应，保证了求和矩阵乘出的力和力矩是总 的合外力和合外力矩。

将V₀＝[-V₀ -V₀β -V₀α]^T，ω＝[p q r]^T代入式(15)中，获得全机气动力 及气动力矩对配平变量的导数如下：

其中，c_p,c_q,c_r为列向量；n₀为所有气动网格法向量列阵，

n_0z为法向量n₀的z方向分量，n_0y为法向量n₀的y方向分量，其元素的形式为：

其中，n_0ix,n_0iy,n_0iz分别为第i个气动网格法向量n_0i的x、y和z方向分量；x_i为 第i个气动网格气动力作用点到取矩点x方向的距离，y_i为第i个气动网格气动力作 用点到取矩点y方向的距离，z_i为第i个气动网格气动力作用点到取矩点z方向的距 离。

值得注意的是，全机的气动力及气动力矩表达式(16)和(17)中，系数矩阵 A_a,B_a,C_a,D_a仅与机翼及尾涡的几何形状相关；与飞行器运动速度相关的变量，即：

攻角α、侧滑角β、滚转角速度p、俯仰角速度q和偏航角速度r，都包含在诱导速 度w表达式中。

S5、计算舵偏气动导数：

求全机气动力对控制面偏角δ的气动导数时，仅通过法向量的变化考虑控制面偏转引起的气动网格的气动力，控制面对应的气动网格不作物理偏转，因此，将法向量 n₀改为n＝n₀+cδ，c为控制面偏转轴矢量，设控制面转轴方向为e，则转动角δ引起 的法向量变化量为：

Δn＝sinδ(e×n₀)+(1-cosδ)[e×(e×n₀)] (52)

当在δ小于5度时：

因此，全机气动力对控制面偏转角δ的导数为：

其中，列向量c_δ为：

c_δi＝V₀·c+(ω×r_0i)·c (56)

S6、通过式(17)、(18)、(19)、(20)、(21)、(26)和(27)获得飞行器的气动 导数。

与现有技术相比，本发明可实现如下有益效果：

(1)本发明的预测方法不需要简谐振动假设就能够直接计算任意运动下的非 定常气动力，预测效率高。

(2)本发明的预测方法基于状态空间形式涡格法实现了飞行器机翼(刚体)的 气动导数的预测，预测过程简单，且提高了预测结果的准确性。

附图说明

附图仅用于示出具体实施例的目的，而并不认为是对本发明的限制。

图1为本发明的飞行器的气动导数预测方法的流程图；

图2为本发明的飞行器的气动导数预测方法的示例设计模型图；

图3为本发明的飞行器的气动导数预测方法的示例模型全机机翼和尾翼的气动网格图；

图4为本发明的飞行器的气动导数预测方法的示例模型展向升力计算结果。

具体实施方式

为了能够更清楚地理解本发明的上述目的、特征和优点，下面结合附图和具体 实施方式对本发明进行进一步的详细描述。需要说明的是，在不冲突的情况下，本发 明的实施例及实施例中的特征可以相互组合。

在下面的描述中阐述了很多具体细节以便于充分理解本发明，但是，本发明还 可以采用其他不同于在此描述的其他方式来实施，因此，本发明的保护范围并不受下 面公开的具体实施例的限制。

本实施例采用有翼无人机模型，全机机翼包括左机翼、右机翼及中部机翼，左机翼、右机翼及中部机翼均为矩形平直机翼，展长相等约为1000mm；在左机翼和右机 翼的后缘均设置有副翼，左副翼及右副翼分别贯穿左机翼及右机翼，控制面弦长占全 机机翼弦长的1/3。全机机翼与尾翼之间通过由机身延伸的尾撑杆连接，尾翼距全机 机翼的垂直距离1020mm。尾翼包括平尾尾翼和垂尾尾翼；平尾尾翼为全动平尾，展 长960mm，弦长120mm，在飞行操纵中起到升降舵作用；垂尾尾翼采用双全动立尾设 计方案，展长240mm，弦长120mm，在飞行操纵中起到方向舵作用。

有翼无人机模型部分参数如表1所示：

表1有翼无人机模型参数

参数名	参数值
		全机机翼展长/mm	3000
副翼展长/mm	1000
		垂尾展长/mm	240
平尾展长/mm	960
		尾翼距全机机翼的垂直距离/mm	1020
全机机翼弦长/mm	270
		副翼弦长/mm	90
垂尾弦长/mm	120
		平尾弦长/mm	120
全机质量/kg	8

S1、建立气动模型：

在有翼无人机模型全机机翼的弧面沿弦向划分8个网格，展向划分90个网格， 拖出的尾涡网格弦向划分160个网格，展向划分90个网格；在尾翼的弧面沿弦向划 分8个网格，展向划分46个网格，拖出的尾涡网格弦向划分160个网格，展向划分 46个网格；全机机翼，尾翼网格划分如图3和4所示。

S2、对步骤S1得到的气动网格建立状态空间形式的气动力方程：

涡格法状态空间方程为：

其中，尾涡强度Γ_w＝[Γ_w0 Γ_wl]^T，A_a,B_a为状态空间系数矩阵，w为包含扰动 在内的机翼表面诱导速度。

气动力矢量F表达式为：

其中，D_a,

C_a为系数矩阵，

为w关于时间的导数

S3、给出状态空间形式的气动力方程的定常形式：

在步骤S2给出的方程(4)及方程(8)基础上，给定状态空间方程的定常形式：

A_aΓ_w+B_aw＝0 (59)

F＝C_aΓ_w+D_aw (60)

S4、求气动力和气动力矩关于飞行器运动变量的导数：

步骤S3得到了状态空间方程的定常形式，求气动力及气动力矩关于飞行器运动变量攻角α、侧滑角β、滚转角速度p、俯仰角速度q和偏航角速度r的导数，其解 析表达式见方程(17)、(18)、(19)、(20)和(21)，计算工况取飞行速度20m/s， 全机配平攻角1度，大气密度1.225kg/m³，全机机翼展向升力分布如图4所示，将飞 行速度、配平攻角和大气密度代入解析表达式计算得到气动导数数值。

S5、计算舵偏气动导数：

步骤S3得到了状态空间方程的定常形式，求气动力及气动力矩关于飞行器舵偏 δ的导数，其解析表达式见方程(26)和(27)，计算工况取飞行速度20m/s，配平攻 角1度，大气密度1.225kg/m³，将飞行速度、配平攻角和大气密度)代入解析表达式 计算得到气动导数数值。

以上所述，仅为本发明较佳的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到的变化或替 换，都应涵盖在本发明的保护范围之内。

Claims (1)

1.一种飞行器的气动导数预测方法，其特征在于，包括以下步骤：

S1、基于涡格法建立气动模型：

在飞行器机翼的弧面沿弦向和展向划分多个干四边形气动网格，气动网格包括附着于机翼表面的附着涡网格及沿来流方向拖出的机翼后缘的第一排尾涡网格以及其他尾涡网格；

S2、对步骤S1得到的气动网格建立状态空间形式的气动力方程：

S2-1、在每个气动网格1/4弦线处布置涡线段；选择1/4弦线的中点为气动网格气动力的作用点，选择3/4弦线中点为气动网格的控制点；

S2-2、将气动模型中的附着涡网格、第一排尾涡网格以及其他尾涡网格对应为：机翼表面的附着涡，机翼后缘的第一排尾涡以及其他尾涡，设机翼表面的附着的涡强度为Γ_b，机翼后缘的第一排尾涡的强度为Γ_w0，其他尾涡的强度为Γ_wl，气动控制方程为：

K_bΓ_b+K_w0Γ_w0+K_wlΓ_wl＝-w (1)

其中，w为机翼表面诱导速度，w＝(V_∞+V_g)·n；V_∞为来流速度矢量，V_g为来流扰动速度矢量，n为控制点处的法向量列阵。K_b,K_w0,K_wl分别为附着涡，第一排尾涡及其他尾涡的诱导系数矩阵；

S2-3、机翼后缘的第一排尾涡在脱出的过程中保持涡强守恒，关系表达为：

其中，

为机翼后缘第一排尾涡的强度关于时间的导数，Δt为时间步长，C₁为保证后缘第一排的尾涡与附着涡对应关系正确的系数矩阵，只包含0和1两种元素；

S2-4、规定其他尾涡脱出后保持涡强不变，关系表达为：

其中，

为其他尾涡的强度关于时间的导数，C₂,C₃为表征其他涡与机翼后缘第一排尾涡位置对应关系的常数提取矩阵，只包含0和1两种元素；

S2-5、综合式(1)、(2)和(3)得到涡格法状态空间方程形式为：

其中，尾涡的强度Γ_w＝[Γ_w0 Γ_wl]^T，A_a,B_a为状态空间系数矩阵，与机翼及尾涡的几何形状及机翼表面的附着涡、机翼后缘的第一排尾涡以及其他尾涡三者的划分有关：

其中，I为单位阵，O为只包含0元素的零矩阵；

S2-6、给出气动网格上的气动力表达，气动网格的气动力表达为：

其中，V_∞为来流速度矢量，e_Γ为前缘涡切向矢量，ρ为空气密度，

为来流单位矢量，

为前缘涡切向单位矢量，A为气动网格面积；Γ_i为第i个气动网格中的有效涡量，数值为沿飞行器飞行方向的前后位置的气动网格的物面涡的涡强之差；

为第i个气动网格上有效涡量时间导数，数值为沿飞行器飞行方向的前后位置的气动网格的物面涡的涡强时间导数的平均值；

S2-7、综合式(2)、(3)、(4)和(7)得到气动网格的气动力的矢量F表达式为：

其中，D_a,

C_a为系数矩阵，

为w关于时间的导数；

S3、给出状态空间形式的气动网格的气动力方程的定常形式：

在步骤S2给出的式(4)及式(8)基础上，忽略涡格法状态空间形式中关于涡强变化的非定常项，状态空间形式的气动网格的气动力方程的定常形式为：

A_aΓ_w+B_aw＝0 (9)

F＝C_aΓ_w+D_aw (10)

其中，在给定来流工况时，系数矩阵A_a,B_a,C_a,D_a与机翼及尾涡的几何形状有关；

S4、求飞行器全机的气动力和气动力矩关于飞行器运动变量的导数：

S4-1、飞行器的运动速度V为：

v＝[V₀+ω×r₀] (11)

其中，V₀为飞行器的平飞速度矢量，ω为飞行器的转动角速度矢量，r₀为飞行器的距重心位置矢量；

S4-2、第i个气动网格的控制点的表面诱导速度w_i为：

w_i＝v_i·n_0i＝(V_0in_0i+(ω×r_0i)n_0i) (12)

其中，v_i为第i个气动网格的运动速度，V_0i为第i个气动网格平飞运动速度,n_0i为第i个气动网格的法向量，r_0i为第i个气动网格到飞行器重心的位置矢量；

S4-3、由式(11)和(12)可得气动网格的气动力为：

其中，G为气动网格的气动力与飞行器机翼结构的插值矩阵；

S4-4、飞行器全机气动力的合外力F_Ar与合外力矩M_Ar通过气动网格的气动力求和得到：

F_Ar＝Φ_tF (14)

M_Ar＝Φ_rF (15)

具体地，求和矩阵Φ_t和Φ_r为：

其中，N_s为气动网格总数，i＝1,2,3…,N_S，

为第i个气动网格气动力的作用点到取矩点的张量；x_i为第i个气动网格气动力的作用点到取矩点x方向的距离，y_i为第i个气动网格气动力的作用点到取矩点y方向的距离，z_i为第i个气动网格气动力的作用点到取矩点z方向的距离，I_3×3为3阶单位矩阵，I_3×3中的取值3×3与气动力向量的x方向力、y方向力和z方向力相适应，保证了求和矩阵乘出的力和力矩是总的合外力和合外力矩；

将V₀＝[-V₀ -V₀β -V₀α]^T，ω＝[p q r]^T代入式(15)中，获得飞行器全机气动力及气动力矩对配平变量的导数如下：

其中，c_p,c_q,c_r为列向量；n₀为所有气动网格法向量列阵，

n_0z为法向量n₀的z方向分量，n_0y为法向量n₀的y方向分量，其元素的形式为：

其中，n_0ix,n_0iy,n_0iz分别为第i个气动网格法向量n_0i的x、y和z方向分量；x_i为第i个气动网格气动力作用点到取矩点x方向的距离，y_i为第i个气动网格气动力作用点到取矩点y方向的距离，z_i为第i个气动网格气动力作用点到取矩点z方向的距离；

S5、计算舵偏气动导数：

将法向量n₀改为n＝n₀+cδ，c为控制面偏转轴矢量，设控制面转轴方向为e，则转动角δ引起的法向量变化量为：

Δn＝sinδ(e×n₀)+(1-cosδ)[e×(e×n₀)] (24)

当在δ小于5度时：

因此，飞行器全机的气动力对控制面偏转角δ的导数为：

其中，列向量c_δ为：

c_δi＝V₀·c+(ω×r_0i)·c (28)

S6、通过式(17)、(18)、(19)、(20)、(21)、(26)和(27)获得飞行器的气动导数。

Images