CN112149234B - 一种基于俯仰角速率输入的飞行器质点运动模型设计方法 - Google Patents

一种基于俯仰角速率输入的飞行器质点运动模型设计方法 Download PDF

Info

Publication number
CN112149234B
CN112149234B CN202011072038.1A CN202011072038A CN112149234B CN 112149234 B CN112149234 B CN 112149234B CN 202011072038 A CN202011072038 A CN 202011072038A CN 112149234 B CN112149234 B CN 112149234B
Authority
CN
China
Prior art keywords
pitch
axis
elevator
angle
aircraft
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
CN202011072038.1A
Other languages
English (en)
Other versions
CN112149234A (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.)
Nanjing University of Aeronautics and Astronautics
Original Assignee
Nanjing University of Aeronautics and Astronautics
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 Nanjing University of Aeronautics and Astronautics filed Critical Nanjing University of Aeronautics and Astronautics
Priority to CN202011072038.1A priority Critical patent/CN112149234B/zh
Publication of CN112149234A publication Critical patent/CN112149234A/zh
Application granted granted Critical
Publication of CN112149234B publication Critical patent/CN112149234B/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/10—Geometric CAD
    • G06F30/15—Vehicle, aircraft or watercraft design
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/10—Geometric CAD
    • G06F30/17—Mechanical parametric or variational design
    • 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
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/04—Constraint-based CAD
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Geometry (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Evolutionary Computation (AREA)
  • Computer Hardware Design (AREA)
  • Pure & Applied Mathematics (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Analysis (AREA)
  • Computational Mathematics (AREA)
  • Aviation & Aerospace Engineering (AREA)
  • Automation & Control Theory (AREA)
  • Control Of Position, Course, Altitude, Or Attitude Of Moving Bodies (AREA)

Abstract

本发明提供了一种基于俯仰角速率输入的飞行器质点运动模型设计方法，利用现有的飞行器六自由度刚体运动方程，仅考虑俯仰运动，不考虑滚转与偏航运动，将原方程的输入由原来的升降舵变为俯仰角速率，将原微分方程中状态量俯仰角速率的特性用典型二阶***来描述，并将飞行器的控制能力限制在角速率的微分方程中；本发明针对以迎角为输入的质点动力学方程的不足，建立一种基于俯仰角速率输入的质点运动模型,将飞行器控制能力的限制以及姿态运动的影响直接融入到质点动力学方程中，为轨迹规划和制导设计提供模型基础。

Description

一种基于俯仰角速率输入的飞行器质点运动模型设计方法
技术领域
本发明涉及飞行器建模技术领域，主要涉及一种基于俯仰角速率输入的飞行器质点运动模型设计方法。
背景技术
随着航空航天技术的快速发展，飞行器的飞行速度越来越快、飞行高度越来越高、气动布局越来越先进、气动特性越来越复杂、飞行任务越来越多样，这都对控制***的设计提出了严重的挑战。飞行器飞行过程中马赫数、迎角、高度、动压变化范围很大，不同状态下的气动特性差异大，飞行状态与动力***耦合严重，姿态运动与质点运动耦合严重。研究姿态运动必须考虑质点运动的影响，研究质点运动也必须考虑姿态运动的影响，而姿态运动又取决于飞行器自身的控制能力。另一方面，受总体、结构、防热***的限制，飞行器的控制能力有限，不同的飞行状态控制能力不同。在大迎角、大马赫数飞行阶段操纵能力明显不足，其质点运动受其控制能力的严格限制，同时姿态的变化对动力***影响严重，轨迹规划及制导设计时，必须充分考虑其控制能力的限制，并将姿态变化的影响严格控制在允许的范围内。因此，需要建立考虑控制约束的质点运动模型，以适应姿态运动与质点运动之间的严重耦合。
目前，国内外轨迹规划及制导设计时普遍采用以迎角为输入的质点动力学方程，通过规划迎角规划轨迹剖面，并进行制导律的设计与仿真。文献“高超声速飞行器多约束再入轨迹快速优化”(2019，Vol.40(No.7)：758～767)给出了以迎角和倾侧角为输入的高超声速飞行器无量纲三自由度运动方程，文献“Integration methods for aircraftscheduling and trajectory optimization at a busy terminal manoeuvring area”(2019，Vol.41(No.3)：641～681)给出了以迎角和倾侧角为输入的无量纲三自由度运动方程。
以迎角为输入的质点动力学方程，虽然可以直接体现迎角的变化及变化率，对于无动力飞行器也可以间接反映对角速率的约束，但是在轨迹设计过程中，还存在一些问题。首先，对于带动力的飞行器，质点动力学方程中迎角的变化并不能直接体现姿态的变化，但是姿态的变化对飞行器的质点运动耦合严重；其次，质点动力学方程中没有体现对控制能力的约束，无法直接判断轨迹设计的合理性。因此，需要构建新的动力学方程，直接融入对姿态变化和控制能力的约束。
发明内容
发明目的：本发明提供了一种基于俯仰角速率输入的飞行器质点运动模型设计方法，针对以迎角为输入的质点动力学方程的不足，建立一种基于俯仰角速率输入的质点运动模型,适应面对称飞行器质点运动与姿态运动之间严重的耦合影响，将飞行器控制能力的限制以及姿态运动的影响直接融入到质点动力学方程中，为轨迹规划和制导设计提供模型基础。
技术方案：为实现上述目的，本发明采用的技术方案为：
一种基于俯仰角速率输入的飞行器质点运动模型设计方法，其特征在于，包括以下步骤：
步骤S1、获取飞行器气动数据；所述气动数据包括机体轴x轴的力系数的基本项Cx0和升降舵产生的增量Cxc、机体轴z轴的力系数的基本项Cz0和升降舵产生的增量Czc、俯仰通道的稳定力矩系数Cm0和控制力矩系数Cmc、俯仰通道的俯仰阻尼导数Cmq和平尾下洗流时差阻尼导数
Figure BDA0002715388650000021
步骤S2、分别建立飞行器气动力系数和气动力矩系数的数学模型；
所述气动力系数基本项是马赫数Ma、迎角α和高度H的非线性函数，升降舵产生的气动力系数增量项是马赫数Ma、迎角α、高度H和升降舵δe的非线性函数，气动力系数表示如下：
Cx0＝Cx0(Ma,α,H)
Cxc＝Cxc(Ma,α,H,δe)
Cz0＝Cz0(Ma,α,H)
Czc＝Czc(Ma,α,H,δe)
所述气动力矩系数包括稳定力矩系数和控制力矩系数两部分，所述稳定力矩系数是马赫数Ma、迎角α和高度H的非线性函数；控制力矩系数是马赫数Ma、迎角α、高度H和升降舵δe的非线性函数；气动力矩系数表示如下：
Cm＝Cm0(Ma,α,H)+Cmc(Ma,α,H,δe)
所述俯仰通道的阻尼导数是马赫数Ma和迎角α的非线性函数，表示如下：
Figure BDA0002715388650000031
步骤S3、获取飞行器推力数据，所述推力是时间t、马赫数Ma和高度H的非线性函数，具体表示如下：
T＝T(t,Ma,H)；
步骤S4、根据当前高度H，计算密度ρ和音速VS，并计算出迎角α、速度V、马赫数Ma和动压
Figure BDA0002715388650000037
密度ρ和音速VS表示如下：
Figure BDA0002715388650000032
Figure BDA0002715388650000033
其中g为重力加速度，e为自然常数；
由当前机体轴x轴和z轴的速度U和W，计算当前的迎角α、速度V、马赫数Ma和动压
Figure BDA0002715388650000034
如下：/>
Figure BDA0002715388650000035
步骤S5、根据发动机推力方向与机体轴x轴的夹角η，计算推力在机体轴x轴和z轴上的分量Tx和Tz：
Figure BDA0002715388650000036
步骤S6、计算满足力矩平衡的升降舵配平舵面δe0；
在当前马赫数Ma、迎角α、高度H下，计算满足力矩平衡的升降舵配平舵面δe0，如下：
Cmc(Ma,α,H,δe0)＝-Cm0(Ma,α,H)；
步骤S7、计算机体系x轴和z轴方向的合力Fx和Fz如下：
Figure BDA0002715388650000041
Figure BDA0002715388650000042
其中，S为机翼参考面积；
步骤S8、计算俯仰角加速度最大值
Figure BDA0002715388650000043
和最小值/>
Figure BDA0002715388650000044
根据升降舵的最大值δemax和最小值δemin，计算俯仰力矩最大值Mmax和最小值Mmin如下：
Figure BDA0002715388650000045
Figure BDA0002715388650000046
其中bA为机翼平均气动弦长；
计算俯仰角加速度的最大值
Figure BDA0002715388650000047
和最小值/>
Figure BDA0002715388650000048
如下：
Figure BDA0002715388650000049
Figure BDA00027153886500000410
其中Iyy为绕机体轴y轴的转动惯量；
步骤S9、计算俯仰力矩对俯仰角速率Q、迎角变化率
Figure BDA00027153886500000411
迎角α和升降舵δe的偏导数：
Figure BDA00027153886500000412
Figure BDA00027153886500000413
Figure BDA00027153886500000414
/>
Figure BDA00027153886500000415
其中，Δα为平衡状态下迎角的扰动量，Δδe为平衡状态下升降舵的增量；
步骤S10、计算俯仰角速率控制等效模型的频率ωn和阻尼ξ如下：
Figure BDA0002715388650000051
Figure BDA0002715388650000052
其中Kp为俯仰角速率反馈到升降舵的增益，Kα为迎角反馈到升降舵的增益；
步骤S11、建立俯仰角速率控制等效模型,计算俯仰角速率变化率
Figure BDA0002715388650000053
和俯仰角加速度变化率/>
Figure BDA0002715388650000054
如下：
Figure BDA0002715388650000055
Figure BDA0002715388650000056
Figure BDA0002715388650000057
其中Qc为等效模型的输入；
步骤S12、计算俯仰角变化率
Figure BDA0002715388650000058
如下：
Figure BDA0002715388650000059
步骤S13、计算机体轴x轴和z轴方向的加速度
Figure BDA00027153886500000510
和/>
Figure BDA00027153886500000511
如下：
Figure BDA00027153886500000512
Figure BDA00027153886500000513
步骤S14、计算高度变化率
Figure BDA00027153886500000514
经度变化率/>
Figure BDA00027153886500000515
和纬度变化率/>
Figure BDA00027153886500000516
如下：
Figure BDA00027153886500000517
Figure BDA00027153886500000518
Figure BDA00027153886500000519
其中，R0为地球半径，ψ为偏航角，取固定值；
步骤S15、根据步骤S11-S14的计算结果，构造基于俯仰角速率输入的飞行器质点运动方程如下：
Figure BDA00027153886500000520
其中状态量x和输入量u分别为：
Figure BDA0002715388650000061
/>
有益效果：本发明具备以下优点：
(1)以俯仰角速率为输入建立飞行器质点运动模型，将对角速率的约束直接融入到质点运动模型中，且通过对俯仰角速率的约束实现对迎角的约束。
(2)建立了俯仰角速率数学模型，将俯仰角速率的模型描述为典型的二阶***，将控制能力的约束直接融入到质点运动模型中。
附图说明
图1是本发明提供的俯仰角速率控制结构图；
图2是本发明提供的飞行器质点运动模型设计流程图。
具体实施方式
下面结合附图对本发明作更进一步的说明。
飞行器的刚体运动可以用微分方程描述如下：
Figure BDA0002715388650000062
其中,x,u分别表示***的状态和输入。飞行器刚体运动方程分为运动学方程和动力学方程，运动学方程描述位置和速度的关系，动力学描述加速度和力/力矩的关系。
飞行器的纵向运动主要指俯仰运动和沿速度方向的线运动，在机体坐标系和地理坐标系下描述飞行器纵向运动的微分方程，飞行状态x包括俯仰角速率Q、俯仰角θ、偏航角ψ、机体系x轴的速度U、机体系z轴的速度W、高度H、经度l、纬度λ，输入u为升降舵δe。
Figure BDA0002715388650000063
下面通过线运动学方程、角运动学方程、线动力学方程和角动力学方程描述飞行器的纵向运动微分方程。
线运动学方程：
Figure BDA0002715388650000071
其中，R0为地球半径。
角运动学方程：
Figure BDA0002715388650000072
线动力学方程：
Figure BDA0002715388650000073
/>
其中，Fx为机体系x轴方向的力，Fz为机体系z轴方向的力，m为飞行器质量，g为重力加速度。
角动力学方程：
Figure BDA0002715388650000074
其中，M为俯仰力矩，Iyy为绕机体轴y轴的转动惯量。
飞行器的合力是气动力和推力的总和，合力矩则为气动力矩和推力矩的总和：
Figure BDA0002715388650000075
其中，飞行器的气动力/力矩与当前飞行状态下的迎角α、马赫数Ma、高度H和升降舵δe有关。Tx为推力在机体轴x轴上的分量，Tz为推力在机体轴z轴上的分量。
仅考虑质点运动时，无法直接获取控制产生的升降舵，式(7)中俯仰力矩无法直接计算。因此本发明采用将输入量由升降舵改为俯仰角速率，俯仰角速率的动态特性用等效二阶***描述，并将舵面带来的控制能力的约束体现在二阶***的描述中。因为小扰动情况下，气动舵面变化对气动力的影响较小，直接用升降舵配平舵面计算气动力。
在当前飞行状态下进行配平和小扰动线性化，平衡状态满足：
Figure BDA0002715388650000081
在平衡状态下，进行小扰动线性化，得到线性化的动力学方程：
Figure BDA0002715388650000082
根据线性化方程可以得到传递函数：
Figure BDA0002715388650000083
考虑俯仰角速率控制律：
Δδe＝KpΔQ+KαΔα (11)
图1给出了俯仰角速率控制结构图，其闭环***传递函数可以描述为
Figure BDA0002715388650000084
其中，ξ和ωn分别为二阶环节的阻尼和频率。Mq为俯仰力矩对俯仰角速率Q的偏导数，
Figure BDA0002715388650000085
为俯仰力矩对迎角变化率/>
Figure BDA0002715388650000086
的偏导数，Mα为俯仰力矩对迎角α的偏导数，/>
Figure BDA0002715388650000087
为俯仰力矩对升降舵δe的偏导数，Kp为俯仰角速率反馈到升降舵的增益，Kα为迎角反馈到升降舵的增益。
将闭环***传递函数描述为微分方程的形式：
Figure BDA0002715388650000088
其中，
Figure BDA0002715388650000089
为俯仰角加速度。
小扰动线性方程和配平状态叠加形成全量微分方程：
Figure BDA00027153886500000810
由于升降舵舵偏受限，因此升降舵产生俯仰角加速度的能力受限，需要约束俯仰角加速度
Figure BDA00027153886500000811
俯仰角加速度最大值/>
Figure BDA00027153886500000812
对应最大舵偏产生的俯仰角加速度，俯仰角加速度最小值/>
Figure BDA00027153886500000813
对应最小舵偏产生的俯仰角加速度：
Figure BDA0002715388650000091
其中，Mmax为升降舵能产生的最大俯仰力矩，
Figure BDA0002715388650000092
为最大俯仰角加速度，Mmin为升降舵能产生的最小俯仰力矩，/>
Figure BDA0002715388650000093
为最小俯仰角加速度。
升降舵产生俯仰角加速度
Figure BDA0002715388650000094
需要满足：
Figure BDA0002715388650000095
用式(14)替换式(6)，式(3)、(4)、(5)、(14)和(16)共同构成基于俯仰角速率输入的飞行器质点运动方程。其中状态量和输入量分别为：
Figure BDA0002715388650000096
下面结合图2，以典型的面对称飞行器为例，阐述基于以俯仰角速率为输入的质点运动模型的实施方式。本发明中所提到的“升降舵”是指俯仰通道的所有控制舵面的统称，不是飞行器机体结构上单个的物理舵面，如左升降、右升降、左V尾、右V尾，左升降和右升降可以一起定义为“升降舵”，左V尾和右V尾也可以一起定义为“升降舵”。
步骤S1、获取飞行器气动数据；所述气动数据包括机体轴x轴的力系数的基本项Cx0和升降舵产生的增量Cxc、机体轴z轴的力系数的基本项Cz0和升降舵产生的增量Czc、俯仰通道的稳定力矩系数Cm0和控制力矩系数Cmc、俯仰通道的俯仰阻尼导数Cmq和平尾下洗流时差阻尼导数
Figure BDA0002715388650000097
步骤S2、分别建立飞行器气动力系数和气动力矩系数的数学模型；
所述气动力系数基本项是马赫数Ma、迎角α和高度H的非线性函数，升降舵产生的气动力系数增量项是马赫数Ma、迎角α、高度H和升降舵δe的非线性函数，气动力系数表示如下：
Cx0＝Cx0(Ma,α,H)
Cxc＝Cxc(Ma,α,H,δe)
Cz0＝Cz0(Ma,α,H)
Czc＝Czc(Ma,α,H,δe)
所述气动力矩系数包括稳定力矩系数和控制力矩系数两部分，所述稳定力矩系数是马赫数Ma、迎角α和高度H的非线性函数；控制力矩系数是马赫数Ma、迎角α、高度H和升降舵δe的非线性函数；气动力矩系数表示如下：
Cm＝Cm0(Ma,α,H)+Cmc(Ma,α,H,δe)
所述俯仰通道的阻尼导数是马赫数Ma和迎角α的非线性函数，表示如下：
Figure BDA0002715388650000101
步骤S3、获取飞行器推力数据，所述推力是时间t、马赫数Ma和高度H的非线性函数，具体表示如下：
T＝T(t,Ma,H)；
步骤S4、根据当前高度H，计算密度ρ和音速VS，并计算出迎角α、速度V、马赫数Ma和动压
Figure BDA0002715388650000102
密度ρ和音速VS表示如下：
Figure BDA0002715388650000103
Figure BDA0002715388650000104
其中g为重力加速度，e为自然常数；
由当前机体轴x轴和z轴的速度U和W，计算当前的迎角α、速度V、马赫数Ma和动压
Figure BDA0002715388650000105
如下：
Figure BDA0002715388650000111
步骤S5、根据发动机推力方向与机体轴x轴的夹角η，计算推力在机体轴x轴和z轴上的分量Tx和Tz：
Figure BDA0002715388650000112
步骤S6、计算满足力矩平衡的升降舵配平舵面δe0；
在当前马赫数Ma、迎角α、高度H下，计算满足力矩平衡的升降舵配平舵面δe0，如下：
Cmc(Ma,α,H,δe0)＝-Cm0(Ma,α,H)；
步骤S7、计算机体系x轴和z轴方向的合力Fx和Fz如下：
Figure BDA0002715388650000113
Figure BDA0002715388650000114
其中，S为机翼参考面积；
步骤S8、计算俯仰角加速度最大值
Figure BDA0002715388650000115
和最小值/>
Figure BDA0002715388650000116
根据升降舵的最大值δemax和最小值δemin，计算俯仰力矩最大值Mmax和最小值Mmin如下：
Figure BDA0002715388650000117
Figure BDA0002715388650000118
其中bA为机翼平均气动弦长；
计算俯仰角加速度的最大值
Figure BDA0002715388650000119
和最小值/>
Figure BDA00027153886500001110
如下：
Figure BDA0002715388650000121
Figure BDA0002715388650000122
其中Iyy为绕机体轴y轴的转动惯量；
步骤S9、计算俯仰力矩对俯仰角速率Q、迎角变化率
Figure BDA0002715388650000123
迎角α和升降舵δe的偏导数：
Figure BDA0002715388650000124
Figure BDA0002715388650000125
Figure BDA0002715388650000126
Figure BDA0002715388650000127
其中，Δα为平衡状态下迎角的扰动量，Δδe为平衡状态下升降舵的增量；
步骤S10、计算俯仰角速率控制等效模型的频率ωn和阻尼ξ如下：
Figure BDA0002715388650000128
Figure BDA0002715388650000129
其中Kp为俯仰角速率反馈到升降舵的增益，Kα为迎角反馈到升降舵的增益；
步骤S11、建立俯仰角速率控制等效模型,计算俯仰角速率变化率
Figure BDA00027153886500001210
和俯仰角加速度变化率/>
Figure BDA00027153886500001211
如下：
Figure BDA00027153886500001212
Figure BDA00027153886500001213
Figure BDA00027153886500001214
其中Qc为等效模型的输入；
步骤S12、计算俯仰角变化率
Figure BDA00027153886500001215
如下：
Figure BDA00027153886500001216
步骤S13、计算机体轴x轴和z轴方向的加速度
Figure BDA00027153886500001217
和/>
Figure BDA00027153886500001218
如下：
Figure BDA0002715388650000131
Figure BDA0002715388650000132
步骤S14、计算高度变化率
Figure BDA0002715388650000133
经度变化率/>
Figure BDA0002715388650000134
和纬度变化率/>
Figure BDA0002715388650000135
如下：
Figure BDA0002715388650000136
Figure BDA0002715388650000137
Figure BDA0002715388650000138
其中，R0为地球半径，ψ为偏航角，取固定值；
步骤S15、根据步骤S11-S14的计算结果，构造基于俯仰角速率输入的飞行器质点运动方程如下：
Figure BDA0002715388650000139
其中状态量x和输入量u分别为：
Figure BDA00027153886500001310
以上所述仅是本发明的优选实施方式，应当指出：对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。

Claims (1)

1.一种基于俯仰角速率输入的飞行器质点运动模型设计方法，其特征在于，包括以下步骤：
步骤S1、获取飞行器气动数据；所述气动数据包括机体轴x轴的力系数的基本项Cx0和升降舵产生的增量Cxc、机体轴z轴的力系数的基本项Cz0和升降舵产生的增量Czc、俯仰通道的稳定力矩系数Cm0和控制力矩系数Cmc、俯仰通道的俯仰阻尼导数Cmq和平尾下洗流时差阻尼导数
Figure FDA0002715388640000011
步骤S2、分别建立飞行器气动力系数和气动力矩系数的数学模型；
所述气动力系数基本项是马赫数Ma、迎角α和高度H的非线性函数，升降舵产生的气动力系数增量项是马赫数Ma、迎角α、高度H和升降舵δe的非线性函数，气动力系数表示如下：
Cx0＝Cx0(Ma,α,H)
Cxc＝Cxc(Ma,α,H,δe)
Cz0＝Cz0(Ma,α,H)
Czc＝Czc(Ma,α,H,δe)
所述气动力矩系数包括稳定力矩系数和控制力矩系数两部分，所述稳定力矩系数是马赫数Ma、迎角α和高度H的非线性函数；控制力矩系数是马赫数Ma、迎角α、高度H和升降舵δe的非线性函数；气动力矩系数表示如下：
Cm＝Cm0(Ma,α,H)+Cmc(Ma,α,H,δe)
所述俯仰通道的阻尼导数是马赫数Ma和迎角α的非线性函数，表示如下：
Figure FDA0002715388640000012
步骤S3、获取飞行器推力数据，所述推力是时间t、马赫数Ma和高度H的非线性函数，具体表示如下：
T＝T(t,Ma,H)；
步骤S4、根据当前高度H，计算密度ρ和音速VS，并计算出迎角α、速度V、马赫数Ma和动压
Figure FDA0002715388640000021
密度ρ和音速VS表示如下：
Figure FDA0002715388640000022
Figure FDA0002715388640000023
其中g为重力加速度，e为自然常数；
由当前机体轴x轴和z轴的速度U和W，计算当前的迎角α、速度V、马赫数Ma和动压
Figure FDA0002715388640000024
如下：/>
Figure FDA0002715388640000025
步骤S5、根据发动机推力方向与机体轴x轴的夹角η，计算推力在机体轴x轴和z轴上的分量Tx和Tz：
Figure FDA0002715388640000026
步骤S6、计算满足力矩平衡的升降舵配平舵面δe0；
在当前马赫数Ma、迎角α、高度H下，计算满足力矩平衡的升降舵配平舵面δe0，如下：
Cmc(Ma,α,H,δe0)＝-Cm0(Ma,α,H)；
步骤S7、计算机体系x轴和z轴方向的合力Fx和Fz如下：
Figure FDA0002715388640000031
Figure FDA0002715388640000032
其中，S为机翼参考面积；
步骤S8、计算俯仰角加速度最大值
Figure FDA0002715388640000033
和最小值/>
Figure FDA0002715388640000034
根据升降舵的最大值δemax和最小值δemin，计算俯仰力矩最大值Mmax和最小值Mmin如下：
Figure FDA0002715388640000035
Figure FDA0002715388640000036
其中bA为机翼平均气动弦长；
计算俯仰角加速度的最大值
Figure FDA0002715388640000037
和最小值/>
Figure FDA0002715388640000038
如下：
Figure FDA0002715388640000039
Figure FDA00027153886400000310
其中Iyy为绕机体轴y轴的转动惯量；
步骤S9、计算俯仰力矩对俯仰角速率Q、迎角变化率
Figure FDA00027153886400000311
迎角α和升降舵δe的偏导数：
Figure FDA00027153886400000312
Figure FDA00027153886400000313
Figure FDA00027153886400000314
Figure FDA00027153886400000315
其中，Δα为平衡状态下迎角的扰动量，Δδe为平衡状态下升降舵的增量；
步骤S10、计算俯仰角速率控制等效模型的频率ωn和阻尼ξ如下：
Figure FDA00027153886400000316
Figure FDA00027153886400000317
其中Kp为俯仰角速率反馈到升降舵的增益，Kα为迎角反馈到升降舵的增益；
步骤S11、建立俯仰角速率控制等效模型,计算俯仰角速率变化率
Figure FDA0002715388640000041
和俯仰角加速度变化率/>
Figure FDA0002715388640000042
如下：
Figure FDA0002715388640000043
Figure FDA0002715388640000044
Figure FDA0002715388640000045
其中Qc为等效模型的输入；
步骤S12、计算俯仰角变化率
Figure FDA0002715388640000046
如下：
Figure FDA0002715388640000047
步骤S13、计算机体轴x轴和z轴方向的加速度
Figure FDA0002715388640000048
和/>
Figure FDA0002715388640000049
如下：
Figure FDA00027153886400000410
Figure FDA00027153886400000411
步骤S14、计算高度变化率
Figure FDA00027153886400000412
经度变化率/>
Figure FDA00027153886400000413
和纬度变化率/>
Figure FDA00027153886400000414
如下：
Figure FDA00027153886400000415
Figure FDA00027153886400000416
Figure FDA00027153886400000417
其中，R0为地球半径，ψ为偏航角，取固定值；
步骤S15、根据步骤S11-S14的计算结果，构造基于俯仰角速率输入的飞行器质点运动方程如下：
Figure FDA00027153886400000418
其中状态量x和输入量u分别为：
Figure FDA00027153886400000419
/>
CN202011072038.1A 2020-10-09 2020-10-09 一种基于俯仰角速率输入的飞行器质点运动模型设计方法 Active CN112149234B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011072038.1A CN112149234B (zh) 2020-10-09 2020-10-09 一种基于俯仰角速率输入的飞行器质点运动模型设计方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011072038.1A CN112149234B (zh) 2020-10-09 2020-10-09 一种基于俯仰角速率输入的飞行器质点运动模型设计方法

Publications (2)

Publication Number Publication Date
CN112149234A CN112149234A (zh) 2020-12-29
CN112149234B true CN112149234B (zh) 2023-05-23

Family

ID=73952672

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011072038.1A Active CN112149234B (zh) 2020-10-09 2020-10-09 一种基于俯仰角速率输入的飞行器质点运动模型设计方法

Country Status (1)

Country Link
CN (1) CN112149234B (zh)

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112965014B (zh) * 2021-03-04 2022-09-02 哈尔滨工业大学 一种飞机机械结构变化引起磁干扰的补偿方法及装置
CN113607378B (zh) * 2021-08-02 2022-04-08 厦门大学 一种绳系支撑飞行器模型强迫自由角运动模拟与抑制方法

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
GB1036721A (en) * 1961-06-16 1966-07-20 Smiths Industries Ltd Improvements in or relating to automatic control systems for aircraft
ITTO20080923A1 (it) * 2008-12-11 2010-06-12 Alenia Aeronautica Spa Procedimento di stima dell'angolo di incidenza e dell'angolo di derapata di un aeromobile
CN106707790A (zh) * 2015-11-13 2017-05-24 成都飞机工业（集团）有限责任公司 一种无人机非线性数学模型建立方法
CN110263497A (zh) * 2019-07-19 2019-09-20 南京航空航天大学 一种基于相对增益的气动耦合影响分析方法
CN110316358A (zh) * 2019-03-29 2019-10-11 南京航空航天大学 基于动态逆的战斗机大迎角控制方法
CN111708377A (zh) * 2020-06-21 2020-09-25 西北工业大学 基于惯导/飞控***信息融合的飞行控制方法

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
GB1036721A (en) * 1961-06-16 1966-07-20 Smiths Industries Ltd Improvements in or relating to automatic control systems for aircraft
ITTO20080923A1 (it) * 2008-12-11 2010-06-12 Alenia Aeronautica Spa Procedimento di stima dell'angolo di incidenza e dell'angolo di derapata di un aeromobile
CN106707790A (zh) * 2015-11-13 2017-05-24 成都飞机工业（集团）有限责任公司 一种无人机非线性数学模型建立方法
CN110316358A (zh) * 2019-03-29 2019-10-11 南京航空航天大学 基于动态逆的战斗机大迎角控制方法
CN110263497A (zh) * 2019-07-19 2019-09-20 南京航空航天大学 一种基于相对增益的气动耦合影响分析方法
CN111708377A (zh) * 2020-06-21 2020-09-25 西北工业大学 基于惯导/飞控***信息融合的飞行控制方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
王建华 ; 刘鲁华 ; 汤国建 ; .高超声速飞行器俯冲段制导与姿控***设计.宇航学报.2016,(第08期),全文. *
高正红,焦天峰.飞行器快速俯仰产生大迎角非定常气动力数学模型研究.西北工业大学学报.2001,(第04期),全文. *

Also Published As

Publication number Publication date
CN112149234A (zh) 2020-12-29

Similar Documents

Publication Publication Date Title
CN111766899B (zh) 一种基于干扰观测器的四旋翼无人机集群抗干扰编队控制方法
Liu et al. Control techniques of tilt rotor unmanned aerial vehicle systems: A review
Crowther et al. Kinematic analysis and control design for a nonplanar multirotor vehicle
CN112241125B (zh) 一种基于微分平坦特性的无人机轨迹跟踪方法
KR102294829B1 (ko) 가변 하중을 가지는 무인비행체 비행제어 시스템 및 방법
Lu et al. Real-time simulation system for UAV based on Matlab/Simulink
CN105425812B (zh) 一种基于双模型下的无人机自动着舰轨迹控制方法
CN112149234B (zh) 一种基于俯仰角速率输入的飞行器质点运动模型设计方法
CN109703768B (zh) 一种基于姿态/轨迹复合控制的软式空中加油对接方法
CN110316400B (zh) 一种鸭翼布局固定翼无人机直接升力控制方法
CN109507890A (zh) 一种基于eso的无人机动态逆广义预测控制器
CN109703769B (zh) 一种基于预瞄策略的空中加油对接控制方法
Pfeifle et al. Cascaded incremental nonlinear dynamic inversion for three-dimensional spline-tracking with wind compensation
CN112327926B (zh) 一种无人机编队的自适应滑模控制方法
CN116991170B (zh) 短距起飞和垂直降落飞机的降落阶段自适应控制设计方法
CN113778129A (zh) 一种干扰补偿的高超声速变后掠翼飞行器跟踪控制方法
CN114721266B (zh) 飞机舵面结构性缺失故障情况下的自适应重构控制方法
CN114815888B (zh) 一种仿射形式的制导控制一体化控制方法
Zhu et al. Trajectory linearization control for a miniature unmanned helicopter
Dogan et al. Control and simulation of relative motion for aerial refueling in racetrack maneuvers
Liu et al. Trajectory stabilization control for aerial recovery of cable-drogue-UAV assembly
CN116088549B (zh) 一种尾座式垂直起降无人机姿态控制方法
Raghavan et al. Flight control for flexible, high-aspect-ratio flying wings
CN114265420B (zh) 适于高动态、控制慢响应的制导控制一体化设计方法
CN107703967B (zh) 一种控制受限飞艇航迹控制方法

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