技术领域
[0001] 本
发明涉及一种大柔性飞行器实时动
力学模型构建方法，特别是一种大柔性飞行器基于在线辨识的建模方法，属于飞行器建模领域。
背景技术
[0002] 高空长航时飞行器是目前国内外航空领域的研究热点，对多种多样的高空飞行任务具有重要的应用价值。这类飞行器一般都采用大展弦比柔性机翼，在机动飞行或者受到外界扰动时会发生几何非线性大结构
变形，进而改变全机升力及力矩分布，影响飞行器的刚体运动。另外，由于结构振动
频率较低，这类飞行器的飞行动力学和结构动力学存在耦合。因此，在进行大展弦比柔性飞行器飞行动力学建模时需要考虑几何非线性结构变形的影响。同时，在飞行中机翼由于大变形，局部剖面
攻角通常较大甚至接近
失速攻角；大变形也会导致流场的变化，这些都会造成
气动非线性。针对大柔性飞行器飞行动力学的机理建模方法主要是结合考虑几何非线性的结构模型和非定常气动力模型，运用多体动力学理论得到系统的动力学模型。这样得到的模型包含大量的柔性模态，阶数较高，不利于控制系统的设计和实现。
[0003] 目前公开文献中提到的大柔性飞行器
控制器设计模型基本可以分为两类：(1)机理模型；(2)机理模型降阶模型。降阶模型虽然解决了系统维数大的问题，但是，该方法从数学模型的角度出发，未考虑真实的物理过程，可能会忽略掉某些重要的物理现象，从而造成灾难性的事故。另外，由于结构非线性和气动非线性的存在，模型存在未建模和不确定现象。为了实现有效控制，需要可靠性更高的面向控制的模型。
[0004] 大柔性飞行器动力学模型是非线性、时变的。用固定模型设计的控制器很难保证飞行器在各种飞行环境下都能安全、稳定和高效的飞行。
发明内容
[0005] 发明目的：为了克服大柔性飞行器控制器设计中遇到的模型阶数高、非线性特性明显和模型存在不确定且时变等问题，本发明提出一种大柔性飞行器基于在线辨识的建模方法。
[0006] 技术方案：一种大柔性飞行器基于在线辨识的建模方法，包括如下步骤：
[0007] 步骤一：用几何非线性梁理论建立大柔性飞行器的大展弦比机翼的非线性结构模型；将大展弦比机翼的总体结构进行离散化，得到单元结构模型；用非定常气动力理论建立气动模型；将外力模型带入单元结构模型，建立单元气弹机理模型；
[0008] 步骤二：根据步骤一中得到的单元气弹机理模型，将可以测量的
节点变量作为单元气弹机理模型的参数，得到单元参数化模型；
[0009] 步骤三：根据步骤二中的单元参数化模型，运用多体动力学理论计算得到考虑结构变形的整体飞行动力学模型；
[0010] 步骤四：步骤二和步骤三中得到的模型包含作为参数的节点变量，通过
传感器数据进行参数在线辨识，得到大柔性飞行器实时飞行动力学模型。
[0011] 进一步的，所述步骤一中，将大展弦比机翼的梁划分为有限个通过节点连接的单元，离散化方法包括有限元法以及有限差分法。
[0012] 进一步的，所述步骤一中，用非定常气动力理论建立气动模型包括如下具体步骤：
[0013] a)用非线性结构模型中的状态量表示气动力计算所需的变量；
[0014] b)根据二维气动模型的表达式和
翼型参数计算机翼的升力、阻力以及
俯仰力矩；
[0015] c)用非定常气动力理论计算机翼运动对气动力的影响，用有限状态入流理论计算机翼运动
对流场的影响；
[0017] e)在二维气动力模型中增加非定常项以及视在质量引起的气动力；
[0018] f)通过转换矩阵得到气动力、力矩在变形标系中的表示，从而得到非定常气动力模型。
[0019] 进一步的，所述步骤一中，外力模型包括非定常气动力模型、重力模型以及集中力模型；其中，重力模型的建立方法为：通过转换矩阵将惯性系中的重力转换到变形
坐标系中，得到重力模型建；集中力模型的建立方法为：将
发动机推力、吊舱气动力、集中质量的重力看作是作用在单个节点的集中力，集中力使得节点左右两侧的内力和内力矩不连续，对存在集中力的节点设置左右两个力，来描述集中力的作用效果。
[0020] 进一步的，所述步骤二中的单元参数模型记为：
[0021]
[0022] 其中， 为单元状态量，pn为可以用传感器测量的单元参数，un为单元可用的控制量，上标n表示第n个单元。
[0023] 进一步的，所述步骤三中，首先运用多体动力学理论推导平均轴系下的动力学方程，然后根据步骤二中的单元参数化模型计算平均轴系动力学方程中的参数，整体飞行动力学模型用平均轴系下的动力学方程描述，记为：
[0024]
[0025] 其中，X为平均轴系的速度和
角速度， 为单元状态向量， 为单元状态向量的时间导数，p为单元参数向量，u为系统控制输入向量。
[0026] 进一步的，所述步骤四中，传感器测量的参数是节点的应变、速度以及角速度。
[0027] 有益效果：本发明有利于解决大柔性飞行器控制器设计中遇到的模型阶数高、非线性特性明显和模型存在不确定且时变的问题，选择合适的传感器可以测量的量作为参数，将单元机理模型参数化，用参数实现单元之间的连接，实现了单元方程间的解耦，通过建立局部的参数化单元模型降低了控
制模型的阶数，用多体动力学理论建立飞行器整体的动力学模型，选择合适的传感器实现参数在线辨识，得到大柔性飞行器实时动力学模型，通过使用
传感器数据辅助建模，使得模型可以包含未建模动态以及其他不确定，提高了模型的可靠性。
附图说明
[0028] 图1为典型大柔性飞行器的外形图；
[0029] 图2为坐标系示意图，图中的F1,F2,F3分别表示局部坐标系下x,y,z方向的梁内力，M1,M2,M3表示局部坐标系下x,y,z方向的内力矩；XB,YB,ZB表示变形坐标系，XR,YR,ZR表示翼根坐标系；
[0030] 图3为传感器安装示意图；
[0031] 图4为大柔性飞行器基于在线辨识的建模
流程图。
具体实施方式
[0032] 下面结合附图对本发明作更进一步的说明。
[0033] 步骤一：大柔性飞行器的大展弦比机翼结构大变形，出现几何非线性现象，用几何非线性梁理论建立机翼的非线性结构模型；为了分析无法找到解析解的连续振动系统，将大柔性飞行器的大展弦比机翼总体结构进行离散化，得到单元结构模型；机翼的运动影响所产生的气动
载荷，用非定常气动力理论建立非定常气动力模型；将气动力等外力模型带入单元结构模型，建立描述单元气弹动力学行为的单元气弹机理模型。
[0034] 具体步骤为：
[0035] s1)大展弦比机翼的非线性结构建模
[0036] 将大展弦比机翼简化为非线性梁模型，用几何非线性梁对大展弦比机翼进行模拟，根据几何精确本征梁理论建立机翼的数学模型。具有扭转、弯曲变形的几何精确本征梁变形坐标系(B系)下的非线性结构模型表示为：
[0037]
[0038]
[0039]
[0040]
[0041] 其中，F,M,V,Ω分别为梁内力、内力矩、速度、角速度向量在B系中的表示，为方程的一级变量；γ,k,P,H分别为力应变、力矩应变、线动量、
角动量向量在B系中的表示，为方T程的二级变量；k表示梁的初始扭转和弯曲；e1＝[1 0 0] ；f,m分别为外力、外力矩向量在B系中的表示。本
实施例中，()′表示对梁坐标s求导，(·)表示对时间t求导，（～)表示该变量对应的反对称矩阵。
[0042] 方程的一级变量和二级变量之间满足如下本构关系：
[0043]
[0044]
[0045] 其中，R,S,T分别是3×3阶的梁截面
刚度系数矩阵；μ,ξ,I分别为线质量、截面质心与刚心的距离以及单位长度的质量惯性矩；ST表示矩阵S的转置，Δ为3×3阶的单位矩阵。
[0046] 上述关系中，前两个是梁的动力学方程，后两个是根据应变-位移关系和速度-位移关系得到的本征运动学方程。
[0047] s2)结构模型有限单元离散
[0048] 梁等连续振动系统的结构动力学模型一般很难求得解析解，为了降低求解难度，先对大柔性飞行器的大展弦比机翼总体结构进行有限单元离散化，将梁划分为有限个通过节点连接的单元，离散化方法包括但不局限于：有限元法以及有限差分法。选择有限差分的离散方式，考虑变量X：
[0049]
[0050] 其中， 表示离散后变量的节点值，上标表示第n个节点；dl为单元长度，X'表示变量X的空间导数。本实施例中，(∧)表示是节点值，(-)表示单元内部变量的值。
[0051] 离散后的单元结构模型表示为：
[0052]
[0053]
[0054]
[0055]
[0056] 本实施中，上标n表示第n个节点或单元。
[0057] 为了进一步说明单元内部变量的值与节点值的关系，以下给出了一级变量的单元内部变量的值和节点值的关系：
[0058]
[0059]
[0060]
[0061]
[0062] 二级变量单元内部变量的值通过本构关系，用一级变量的单元内部变量的值进行表示。
[0063] s3)非定常气动力模型建立
[0064] 基于片条理论，大展弦比机翼的气动载荷可以用二维气动模型进行计算，同时需要考虑机翼运动对气动载荷的影响。片条理论把机翼考虑为由一系列弦向“片条”元素组成，并假定机翼每个弦向片条的升力系数正比于当地攻角，且各个片条上的升力彼此无影响，即假设机翼每个片条上的升力等于这个片条属于一个二维无限展长机翼上一部分时的升力，翼型截面的性质将被使用。当机翼相对气流存在运动时，气动力和力矩将随时间变化，而且机翼运动产生的尾流变化也会影响气动力和力矩。首先用准定常的二维气动力和力矩表达式，然后增加尾流及视在质量的惯性反作用力对气动力、力矩的影响。具体计算方法如下：
[0065] a)用非线性结构模型中的状态量表示气动力计算所需的变量。
[0066] b)根据二维气动模型的表达式，根据翼型参数计算机翼的升力、阻力以及俯仰力矩。
[0067] c)用非定常气动力理论计算机翼运动对气动力的影响，用有限状态入流理论计算机翼运动对流场的影响。
[0068] d)计算视在质量引起的气动力。
[0069] e)在二维气动力模型中增加非定常项以及视在质量引起的气动力。
[0070] f)通过转换矩阵得到气动力、力矩在变形坐标系中的表示，得到非定常气动力模型。
[0071] s4)重力模型建立
[0072] 大展弦比机翼运动需要考虑重力的影响，由于动力学方程建立在变形坐标系中，而我们只知道重力在惯性系中的表达式，所以通过转换矩阵将惯性系中的重力转换到变形坐标系中。重力和重力产生的力矩表示为：
[0073] fg＝μg
[0074]
[0075] 其中，fg是重力，mg是重力产生的力矩，g是重力
加速度向量在B系下的列向量表示。
[0076] B系中重力向量满足以下方程：
[0077]
[0078] 离散化以后可以写为：
[0079]
[0080] 在质心所在的参考节点处满足时间微分方程，其余节点处的g用空间离散方程计算。
[0081] s5)集中力模型建立
[0082] 发动机推力、吊舱气动力、集中质量的重力都可以看作是作用在单个节点的集中力，这些力使得节点左右两侧的内力和内力矩不连续，对存在集中力的节点设置左右两个力，来描述集中力的作用效果。
[0083] s6)单元气弹机理模型建立
[0084] 将非定常气动力模型、重力模型以及集中力模型带入s2)中得到的单元结构模型中，得到基于气弹理论建立的单元气弹机理模型：
[0085]
[0086]
[0087]
[0088]
[0089]
[0090] 其中， 表示B系下的单元气动力和气动力矩， 为单元
重力加速度向量在B系下的表示。本实施例中，下标l表示单元右节点的左值，下标r表示单元左节点的右值，下标aero表示由空气产生的，下标lr表示节点右值到左值的变换。
[0091] 以下的二级方程用于一级变量和二级变量之间的相互计算：
[0092]
[0093]
[0094] 步骤二：单元气弹机理模型参数化
[0095] 单元参数模型记为：
[0096]
[0097] 其中， 为单元状态量，pn为可以用传感器测量的单元参数，un为单元可用的控制量，上标n表示第n个单元。
[0098] 离散后的单元通过节点相互联系，选择合适的传感器可以测量的节点变量作为参数，得到单元参数化模型。
[0099] 参数化后的单元动力学模型表示为：
[0100]
[0101]
[0102]
[0103]
[0104] 以及
[0105]
[0106]
[0107] 其中， 为单元状态量；Pin为可以用传感器测量的单元参数，如图3所示，本实施例中， 表示第n个单元右节点内力的左值， 表示第n个单元左节点内力的右值， 表示第n个单元右节点内力矩的左值， 表示第n个单元左节点内力矩的右值， 表示第n个单元右节点速度左值， 表示第n个单元左节点速度右值， 表示第n个单元右节点角速度左值， 表示第n个单元左节点角速度右值； 为单元气动力， 为单元气动力矩，un为单元可用的控制量。
[0108] 步骤三：整体飞行动力学模型建立
[0109] 首先运用多体动力学理论推导平均轴系下的动力学方程，然后根据步骤二中的单元参数化模型计算平均轴系动力学方程中的参数，整体飞行动力学模型用平均轴系下的动力学方程描述，记为：
[0110]
[0111] 其中，X为平均轴系的速度和角速度， 为单元状态向量， 为单元状态向量的时间导数，p为单元参数向量，u为系统控制输入向量。
[0112] 大柔性飞行器的质量特性和刚度特性会由于变形发生改变，因此无法建立和刚体飞行器一样的
机体坐标系描述飞行器的飞行动力学行为，借助多体动力学理论，用多刚体系统代替柔性系统，建立飞行器平均轴系下的动力学模型，平均轴系的原点位于飞行器的瞬时质心，并且平均轴系下由于变形引起的线动量和角动量为零，模型中的参数用单元气弹机理模型中的状态转换到翼根坐标系R下计算得到，平均轴系下的动力学方程表示为：
[0113]
[0114]
[0115] 其中，
[0116]
[0117]
[0118] 为平均轴系的平动速度和角速度，在翼根坐标系中进行表示。Mtotal为系统总质量，R为对应于平均轴系原点的质量惯性矩阵，P为R坐标系下系统总的线动量， 为相对于平均轴系原点的角动量， 分别为第n个刚体相对于平均轴系原点的
位置矢量以及变形坐标系(B系)到翼根坐标系(R系)的转换矩阵，μn为刚体单元的质量， 为刚体单元的质量惯性矩， 分别为刚体单元的速度和角速度矢量，以上所有变量都可以通过步骤二中单元参数化模型中的变量计算得到。
[0119] 步骤四：在线辨识实现
[0120] 步骤二中得到的单元参数化模型包含了作为节点变量的参数，用传感器在线测量，经过必要的计算辨识参数值，得到实时的单元模型，根据步骤三进一步计算出整体飞行动力学模型。步骤二中作为参数的变量可以是节点的
应力、应力矩、速度以及角速度，节点的应力和应力矩无法直接测量，通过测量节点的应变，根据材料的本构关系可以计算得到应力以及应力矩；速度以及角速度可以使用
空速传感器和
陀螺仪测量。
[0121] 以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。
