扑翼飞行器动力系统建模方法

2021-10-20年鹏宋笔锋宣建林王思琦

航空学报 2021年9期

年鹏，宋笔锋，宣建林,2,*，王思琦

1. 西北工业大学航空学院，西安 710072

2. 西北工业大学长三角研究院，太仓 215400

扑翼飞行器作为一种模仿生物飞行的特殊飞行器，在军民用领域都具有独特的用途。近年来扑翼飞行器都是国内外的研究热点[1-3]，并且随着微电子、微机械技术的发展，一系列扑翼飞行器被成功研制出来[4-7]。扑翼飞行器的优势之一在于低雷诺数下优异的气动性能[8]，因此航时是扑翼飞行器的一项重要的性能指标。但由于扑翼飞行器具有独特的非定常空气动力学机理[9-11]，因此在扑翼飞行器的设计中直接采用数值计算的结果通常不会得到较为准确的飞行性能，如著名的扑翼飞行器Nano Hummingbird[12]和Delfly[13]的设计过程中都大量采用了实验研究的手段，包括气动力的测量实验和自由飞实验。

在改进扑翼飞行器的设计后，为准确评估扑翼飞行器的续航性能，目前多依靠实际外场试飞来采集数据，由于环境的不可控性，试飞结果重复性较差，且外场试飞费时费力，这在一定程度上阻碍了扑翼飞行器的研制。动力系统作为扑翼飞行器的主要功耗源直接影响扑翼飞行器的续航性能，因此有必要建立扑翼飞行器动力系统的仿真模型，利用仿真手段部分替代传统试飞测试的作用。在现有的研究中，Doman等[14]使用拉格朗日分析法建立了包含电机、驱动机构和弹性前缘梁的扑翼飞行器动力系统模型，通过实验验证了模型的精度，该模型可用于载荷与电机的匹配，但模型仅考虑了前缘梁的惯性载荷，而没有对扑动翼的气动载荷进行建模。Kalpathy Venkiteswaran和Sun[15]建立的扑翼动力系统模型包括驱动机构多体动力学模型、简化气动力模型和直流电机模型，该模型可用于扑翼驱动机构功耗的计算和参数优化中，但模型没有考虑电池的影响，且简化气动力模型未经验证。

为较为准确地评估扑动翼改进后对续航性能的影响，便于针对确定的机翼设计进行动力系统的匹配，逐步建立基于模型的扑翼飞行器设计方法，以本团队研制的“信鸽”扑翼飞行器[16]的动力系统建模为例，建立了一套扑翼动力系统各部件的高精度动态建模方法，并通过自由飞测试对仿真系统的精度进行了验证。

1 仿真系统的基本架构

扑翼动力系统仿真需要包含动力能源模型、气动模型、机体动力学模型和控制模型，其结构如图1所示，其中动力能源模型和气动模型是仿真系统的核心模块。

图1 仿真系统基本架构Fig.1 Basic architecture of simulation system

动力能源模型中包含电机动态模型、扑动机构动力学模型、电池动态模型和电调模型，电调模型需要油门控制量来计算电机输入端等效电压系数，气动模型计算得到的瞬时扑动轴负载输入动力能源模型后经扑动机构动力学模型传递至电机动态模型，计算得到瞬时电流作为电池动态模型的负载，电池动态模型的输出电压经过电调输入电机模型，形成动力能源模型的闭环仿真。

气动模型中包含扑动轴负载模型、扑动翼气动力/矩模型和尾翼模型。其中扑动轴负载模型利用动力能源模型输出的机翼瞬时扑动参数和机体动力学模型输出的迎角和空速来计算得到时变的扑动轴负载，扑动翼气动力/矩模型和尾翼模型则用于配合机体动力学模型和控制模型计算仿真所需的配平状态。

2 电机、电调和电池建模

2.1 电机模型

扑翼飞行器目前多采用转速低、扭矩大的直流无刷电机作为主要动力源来驱动扑动机构，直流电机可等效为如图2所示的等效模型[17]。图中：Em为电机输入电压；Im为电机电流；Rm为电机内阻；Lm为电机电感；Ei为感应电动势；Tm为电机输出转矩。

图2 直流电机等效模型[17]Fig.2 Equivalent model of DC motor[17]

由电压平衡有

(1)

由转矩平衡有

(2)

Tm=Km(Im-I0)

(3)

式中：Ke为反向感应电动势常数；Km为转矩常数；I0为空载电流；ωm为电机角速度；Jm.equ为电机轴等效转动惯量；Tload.m为电机轴负载转矩。

关于式(1)中电机电感的影响，电机在扑动翼的交变气动载荷作用下，电机电流存在与载荷同频率的波动(通常小于10 Hz)，通常无刷直流电机的电感为mH量级，这部分电流波动的影响可忽略不计，式(1)可简化为一阶形式,即

Em=RmIm+Keωm

(4)

由式(3)和式(4)知直流电机模型的基本参数包括Ke、Km、Rm和I0，为提高电机模型的精度，本文使用电机实际测试数据进行求解得到电机模型的基本参数。电机测试方法为在给定电压下通过摩擦片连续加载直到全速转动的电机堵转为止，连续采集此过程中电机的负载转矩、转速、电流和电压数据。分别将相邻两个时刻的Em、Im和ωm代入式(4)并联立求解出一组Rm和Ke，遍历整个实验数据，得到不同实验状态下的Rm和Ke。Rm和Ke通常被认为是常数，因此对所有的Rm和Ke取平均值作为电机模型的参数值。使用同样的方法根据式(3)可确定I0和Km的值。表1为使用本文方法获得的“信鸽”扑翼飞行器使用的MN1804电机的模型参数值，图3是给定电压和转矩下电机电流、转速的计算值与实验值的相对误差对比，可知建立的电机模型在实际工作转矩区间(0.01 ～0.03 N·m)相对误差小于10%。

表1 MN1804电机模型的参数值Table 1 Parameter values of MN1804 motor model

图3 MN1804电机模型的相对误差Fig.3 Relative error of MN1804 motor model

2.2 电调模型

将直流无刷电机等效为直流电机进行建模能够获得较高精度的模型，在直流电机模型中通过调节输入电压来控制电机转速，但实际中直流无刷电机是通过电调(ESC)调节电机输入电压的占空比来实现速度调节的，这就需要建立占空比与电机等效输入电压之间的关系。忽略电调的内阻，电调可简化为如图4所示的模型，其中：Eb为电池电压；Ib为电池电流。

图4 简化电调模型Fig.4 Simplified electronic speed controller model

电池的电压Eb经过电调的脉宽调制后的输出给电机的等效电压Em可表示为

Em=σEb

(5)

式中：σ为等效电压系数。

电调的效率通常较高，可认为电调的输入功率PESC.in约等于输出功率PESC.out[18]，得到

(6)

将式(5)代入式(6)可得到电机电流与电池电流的关系

(7)

等效电压系数与电调进行电压脉宽调制的占空比直接相关，调制后的输出电压一般为方波形式，理论上等效电压系数等于占空比信号的有效值。但实际中电调为抑制电机换相过程中的转矩脉动，在换相时刻引入了重叠换相机制[19]，导致不同电调的等效电压系数存在差异。本文实测了恒定电压下不同占空比时电机的转速和电流，再结合式(4)和表1中的电机参数计算得到不同占空比下电机端的等效输入电压，最后利用式(5)计算出不同占空比时的等效电压系数。图5是“信鸽”扑翼飞行器使用的FLYFUN-6A电调的等效电压系数实测曲线，由图可知该电调的占空比与等效电压系数存在指数关系y=aebx+cedx,其中:a=0.830 7;b=0.084 38;c=-0.831 2;d=-6.717;R2=0.999 8。

图5 FLYFUN-6A电调的等效电压系数曲线Fig.5 Equivalent voltage coefficient curve for FLYFUN-6A ESC

2.3 电池模型

锂聚合物电池具有放电倍率大、能量密度高等特点，因此扑翼飞行器常采用锂聚合物电池作为能源。由于锂电池存在倍率容量效应，即随着放电倍率的增大，电池实际能放出的电量呈下降趋势。为模拟这种效应，在无人机航时评估中常使用Peukert经验公式[20]来建立电池模型。由于Peukert公式假设电池电压与放电量存在线性关系，公式在放电初段和放电末段对电压的描述与实际测试结果有较大差异。为更准确地建立整个放电周期内电池电压的变化模型，同时考虑扑翼波动载荷对电池电压的影响，本文使用Tremblay和Dessaint[21]在Shepherd电池模型的基础上改进的动态电池模型。改进的Shepherd电池模型考虑了开环电压与荷电状态的关系，并修正了放电带来的极化电阻，改进的模型能够准确描述放电全过程的电池电压和变化电流下电池的恢复效应。改进的Shepherd电池模型的电压公式为

(8)

图6 电池放电模型等效电路Fig.6 Equivalent circuit of battery discharge model

为获得电池模型的参数，需要对待测电池进行瞬时电压、瞬时电流、放电电荷量等放电数据的实验采集，考虑到扑动翼载荷作用下电池的负载电流变化幅度较大，在采集放电数据时应尽可能地覆盖放电电流的变化范围。基于采集的电池放电数据使用参数辨识算法即可确定模型的参数。以“信鸽”扑翼飞行器使用的定制锂聚合物电池建模为例，首先使用PowerLab电池测试仪采集被测电池在放电电流Ib分别为3、4、5 A时的电压、电流和放电电荷量，选择的放电电流基本涵盖了“信鸽”实际飞行中能达到的电流区间。然后使用最小二乘法进行参数辨识，得到表2中的电池模型参数值，建立的电池模型与实际实验值的相对误差如图7所示。

表2 “信鸽”电池模型的参数值Table 2 Parameter values of Dove’s battery model

由图7可知本文建立的电池模型与实验值相比误差小于6%。此外，电池在负载电流消失后电压迅速回升的恢复效应对于扑翼飞行器这种负载电流动态波动的对象非常重要。本文的电池模型使用如式(9)所示的一阶低通滤波器来模拟电池的恢复效应:

图7 “信鸽”电池模型的相对误差Fig.7 Relative error of Dove’s battery model

(9)

式中：Tbr为电池的动态响应时间。

图8所示为考虑恢复效应的电池模型以4 A放电，直至电池电压低于6 V时将放电电流突然置为零，记录此过程中电池电压变化的仿真数据。与同状态实验数据对比可知模型对恢复效应的仿真与实际结果十分接近。

图8 “信鸽”电池的恢复效应建模结果Fig.8 Modeling result of recovery effect of Dove’s battery

基于上文建立的电池模型，本文还对电池输入为波动电流和对应电流的等效值时放电时间的差异进行了对比。首先根据扑翼飞行器的特点假定波动电流负载为

i(t)=asin(bt)+c

(10)

式中：a为电流的幅值；b为电流的角速度；c为电流的偏移值；对应地式(10)电流的有效值根据功耗相等可计算并化简为

(11)

式中:ieff为电流的有效值。

由式(11)可知电流的有效值与b无关，为较为定量地进行比较，用a/c描述电流的波动程度。使用建立的考虑恢复效应的电池模型，选取不同的a/c值依次使用波动电流和对应电流的有效值进行放电时间仿真(具体取值如表3所示)，结果如图9所示，可知随着电流波动程度的加剧，受电池恢复效应的影响，使用电流有效值进行仿真的误差越来越大。这表明，对于扑翼飞行器这类负载变化较为剧烈的仿真对象，有必要使用实际的动态电流负载进行仿真以提高航时的仿真精度。

表3 波动电流放电时间仿真a/c取值Table 3 Values of a/c for fluctuating current discharge time simulation

图9 波动电流程度与恒定电流放电时间差异的关系Fig.9 Relationship between undulation current and difference of constant current discharge time

3 扑动机构建模

3.1 机构运动学建模

第2节讨论了动态负载对航时仿真精度的重要性，对机构的运动学建模则是动态负载建模的基础。常见的扑动机构可简化为图10，根据前期研究中的建模结果[22]，机构的运动学方程可表示为

(12)

(13)

式中：l1～l4为机构的杆长；G2为杆2的质心位置；lG2为l2杆质心到l1杆铰接点的距离；其余角度在图10中定义；lBD为图10所示铰接点B、D的间距；θ1为机构曲柄的角位移；θ3根据实际情况调整零点位置后即为机构输出的实际扑动角位移。

图10 扑动机构运动学示意图Fig.10 Kinematics diagram of flapping mechanism

机构角速度的矩阵表达形式为

(14)

式中：ω1为机构曲柄l1的角速度，等于电机经过减速器后的输出角速度；ω2为杆2绕杆1旋转的角速度；ω3为机构输出的扑动角速度。

3.2 电机轴负载转矩计算

扑动翼的气动载荷通过机构连杆的传递以转矩的形式施加在减速器上，然后经过减速器最终施加在电机轴上，扑动机构载荷传递示意图如图11 所示。根据力矩平衡有

图11 扑动机构载荷传递示意图Fig.11 Load transfer diagram of flapping mechanism

Twing=ηhingeF23l3sin∠BCD

(15)

Tgear=F21l1sin∠CBA

(16)

式中：Twing为单个扑动翼的翼根转矩；Tgear为减速器的输出转矩；ηhinge为一个铰接副的传动效率；F23为杆2对杆3的力；F21为杆2对杆1的力。

由于杆2两端铰接，忽略杆2的重力，除两端铰接点外其他部位不受力，因此杆2是二力杆，则

(17)

此外，根据减速器两端所受力矩平衡，有

Tload.m=Tgear/(ηgi)

(18)

式中：ηg为减速器的效率；i为减速器的总减速比。

结合式(15)～式(18)，并假设左右扑动翼带来的气动载荷大小相等，可以得到扑动翼的载荷经过扑动机构和减速器传递至电机轴后的转矩为

(19)

式中：

(20)

(21)

3.3 电机轴等效转动惯量计算

为利用式(2)计算电机的转速，将扑动翼、扑动机构、减速器等运动部件的转动惯量统一等效至电机轴。电机轴为整周连续转动，一个扑动周期内运动部件的重力势能没有发生改变，因此等效过程中仅需保证动能守恒，有

(22)

式中：Jm.g1为无刷电机外转子和一级减速器小齿轮组合体绕轴旋转的转动惯量；Jg23为一级减速器大齿轮和二级减速器小齿轮组合体绕轴旋转的转动惯量；Jg4为二级减速器大齿轮绕轴旋转的转动惯量；ml2为杆2的质量；Jl2.G2为杆2绕其质心G2旋转的转动惯量；Jrocker为扑动机构的单侧摇臂绕轴旋转的转动惯量；Jwing为单侧扑动翼绕轴旋转的转动惯量；ηg12为一级减速器齿轮副的传动效率；ηg34为二级减速器齿轮副的传动效率；ωg23为一级减速器大齿轮和二级减速器小齿轮组合体的角速度；VG2为杆2的质心G2点的速度。

此外，按照角速度与减速比的关系，可得到

(23)

式中：i1为一级减速器的减速比。

结合式(22)和式(23)，得到如式(24)所示的电机轴等效转动惯量的表达式为

(24)

式中:定义角速度均以逆时针旋转为正，VG2按照图10中的速度正交分解可表示为

(25)

4 扑动翼气动载荷建模

直接实验采集机翼扑动时产生的动态翼根气动转矩载荷具有一定的难度，本文利用风洞实验中便于采集的数据，提出了一种动态翼根气动转矩的建模方法，分为翼面法向力的建模和其作用点即瞬时展向压力中心位置的建模，将建立的瞬时翼面法向力模型与建立的瞬时展向压力中心位置模型相乘得到翼根转矩模型。

4.1 瞬时翼面法向力建模

影响单侧翼面法向力FN大小的自变量主要有：扑动角位移θflap、扑动角速度ωflap、风速V、迎角α，直接使用数值计算求解瞬时非定常气动力将严重拖慢仿真程序。为获得更为精确、且计算速度更快的气动力模型，基于风洞试验数据提出了瞬时翼面法向力半经验建模方法，即

Fz=Fflap(θflap,ωflap)+Fflow(V,α)

(26)

根据试验数据的规律和实际试飞的经验，扑翼飞行器的纵向气动特性与固定翼类似，体轴系下z方向的气动力Fz可以分解为扑动产生的气动力和来流产生的气动力。其中，扑动产生的气动力主要与扑动参数θflap和ωflap有关；来流产生的气动力主要与V和α有关。实际中由于迎角为零时，周期平均Fz的值接近零，可认为此状态下Fz与来流无关，建立θflap、ωflap与Fz的多项式关系作为扑动产生的气动力。对于其他实验状态下的数据将其与扑动气动力的数据相位对齐后相减，建立V、α与Fz的多项式关系作为来流产生的气动力。

建立Fz的模型后经过式(27)的转换，有

FN=Fz/2cosθflap

(27)

即得到垂直于机翼的单侧翼面法向力FN。本建模方法能够保持实验数据精度的同时，实现实验状态的内插和一定程度的外推。

以信鸽的扑动翼建模为例，分别使用一阶多项式和二阶多项式建立Fflap函数和Fflow函数。图12 为风速10 m/s、扑动频率4 Hz、多个不同迎角状态下，Fz建模结果与对应风洞试验值的单周期瞬时值对比。其中模型的确定系数R2为0.914 8，标准差RMSE为0.453 7，可知建立的Fz模型具有较高的准确性。

图12 “信鸽”扑动翼Fz建模结果Fig.12 Modeling result of Fz of Dove’s flapping wing

4.2 瞬时展向压力中心位置建模

扑动翼的气动载荷在风洞实验中最终体现在实验机构驱动电机的功率消耗上，因此在已知实验电机参数和实验电机输入功率的情况下，可以间接计算得到展向压力中心的位置，具体计算流程见图13。

图13 瞬时展向压力中心位置计算流程图Fig.13 Calculation flow chart of instantaneous position of spanwise center of pressure

首先由实验电机驱动器的输入功率Pin计算得到实验电机的输入功率Pm.in，即

Pm.in=ηESC.expPin

(28)

式中：ηESC.exp为实验电机驱动器的效率。然后根据一阶电机模型有

(29)

Tm.exp=(Im.exp-I0.exp)Km.exp

(30)

式中：Ke.exp为实验电机的反向感应电动势常数；Km.exp为实验电机转矩常数；I0.exp为实验电机空载电流；Rm.exp为实验电机内阻；ωm.exp为实验电机角速度；f为机翼的扑动频率；iexp为实验机构的减速比。

求解式(29)中的一元二次方程并舍去负值解，得到实验电机的电流Im.exp。解出的Im.exp代入式(30)计算得到实验电机的输出转矩Tm.exp。

使用与3.2节相同的计算方法，按照实验机构减速器输入输出能量平衡和减速比关系式(18)，计算得到减速器输出转矩Tg.exp。根据式(19) 解算出单侧翼根转矩Twing的大小，最后利用式(31)得到瞬时展向压力中心的位置lcp.y。

lcp.y=Twing/FN

(31)

扑动翼的瞬时展向压力中心位置同样主要与扑动参数θflap、ωflap和来流参数V、α这4个参数有关，由于参数相互关系较为复杂，本文采用了神经网络拟合的方式来建立以上4个参数与瞬时展向压力中心位置的映射模型。以信鸽的扑动翼建模为例，设置神经元数量为5个，训练神经网络拟合模型。图14为风速10 m/s、扑动频率4 Hz、多个不同迎角状态下，lcp.y建模结果与对应状态解算值的单周期对比。其中模型的确定系数R2为0.896 4，标准差RMSE为0.002 504，可知建立的lcp.y模型具有较高的准确性。

图14 “信鸽”扑动翼lcp.y建模结果Fig.14 Modeling result of lcp.y of Dove’s flapping wing

5 纵向配平与控制系统建模

动力系统仿真时需要首先确定配平状态，动力系统主要影响飞行器的纵向性能，因此本文主要考虑纵向的配平。通过力/矩平衡确定扑翼飞行器纵向配平状态，以提供扑动翼气动载荷模型需要的输入参数(V、α、θflap、ωflap，其中扑动参数θflap和ωflap由配平后的油门控制)。为避免使用瞬时气动力进行配平时配平状态跳动对仿真收敛性产生影响，本文建立了扑动翼周期平均气动模型和平尾气动模型用于配平状态计算。

5.1 扑动翼周期平均气动建模

用于纵向配平计算的扑动翼周期平均气动力模型包括周期平均升力模型、周期平均净推力模型和周期平均俯仰力矩模型。对于柔性扑动翼的气动建模，使用风洞实验数据更为准确，但由于风洞实验状态是离散的，不利于配平状态的求解，本文使用了前期研究[23]中提出的基于风洞试验数据的半经验周期平均扑动翼气动建模方法:

(32)

(33)

并基于相同理念将建模方法推广至周期平均俯仰力矩建模:

(34)

以信鸽的扑动翼气动建模为例，基于风洞试验数据使用式(32)～式(34)建立扑动翼的周期平均气动模型，建立的模型的精度如表4所示，可知建立的模型具有较高的精度。

表4 “信鸽”扑动翼周期平均气动建模精度

5.2 平尾气动建模与纵向平衡

飞行器全机纵向受力示意如图15所示，根据力/矩平衡有

图15 扑翼飞行器全机受力示意图Fig.15 Force diagram of flapping wing vehicle

(35)

扑动翼关于全机重心的俯仰力矩可表示为

(36)

式中：αw为扑动翼的迎角；lcg.x为全机重心距扑动翼前缘水平距离；lcg.z为全机重心距扑动翼平均气动弦的垂直距离；lcp.x为扑动翼弦向压力中心距前缘的距离:

(37)

忽略平尾的阻力和零升力矩，平尾关于全机重心的俯仰力矩可表示为

Mt=-Ltcosαtlht-Ltsinαtlcg.z

(38)

式中：lht为平尾力臂；αt为平尾的迎角，等于机翼迎角αw减下洗角ε。

为便于快速调整平尾的设计参数以适配不同的扑动翼，本文使用估算公式建立平尾的气动模型。平尾的升力模型为

Lt=kqqSt[CLα.t(αw-ε)+CLα.tηeδe]

(39)

式中：kq为速度阻滞系数，视平尾相对机身的位置取0.85～1.0；q为动压；St为平尾面积；δe为升降舵偏角；ηe为升降舵效率系数；CLα.t为平尾升力线斜率。

基于经典升力线理论可以粗略估计低速时平尾处的下洗角为[24]

(40)

式中：CLα.w为扑动翼升力线斜率，可使用风洞试验数据计算得到；AR为扑动翼的展弦比。

低速飞行器的升降舵效率系数为[25]

(41)

忽略平尾处机身直径对平尾升力的影响，直梯形平尾升力线斜率可以表示为[26]

(42)

式中：Kβ为普朗特-葛劳渥压缩性修正因子，即

(43)

其中:AR.t为平尾展弦比；λt为平尾梢根比。

5.3 纵向控制系统建模

为保持配平状态，本文搭建了如图16和图17 所示的简化纵向控制系统模型，按照实际仿真需求纵向控制系统包括高度控制系统和速度控制系统。高度控制系统使用了串级PID控制器，首先PID控制器1根据期望飞行高度和当前飞行高度的差计算出期望迎角，期望迎角与当前迎角的差再经过PID控制器2得到升降舵偏角，升降舵偏角输入平尾模型用于计算配平力矩。速度控制系统使用期望飞行速度和当前飞行速度的差经PID控制器3得到油门大小，油门大小输入电调模型用于计算等效电压系数。控制系统中的姿态和位置反馈数据由刚性飞行器动力学方程[25]解算得到。

图16 高度控制系统示意图Fig.16 Diagram of altitude control system

图17 速度控制系统示意图Fig.17 Diagram of velocity control system

6 仿真模型精度验证

将以上建立的电机、电调、电池、扑动机构、扑动翼、平尾、纵向控制系统等模型集成为扑翼飞行器仿真系统，设定起飞海拔450 m，期望飞行高度50 m，手掷起飞初速度9 m/s，期望巡航速度10 m/s，机载电子设备功耗设定为7 W，进行极限航时仿真。参考试飞数据，将仿真结束条件设定为电池电压小于6.6 V。航时仿真过程中电池电压与“信鸽”扑翼飞行器同状态飞行测试机载数据的对比如图18所示。仿真中的极限航时为1 226.74 s，实际试飞航时为1 196.6 s，相对误差2.52%，证明本文建立的扑翼飞行器动力系统仿真模型具有较高的精度。