四旋翼数学模型—

一、动力模型

动力模型研究四旋翼机体系 xyzxyzxyz 轴上的力矩与四个电机输入油门的关系。这个过程中，假设电池电压不变，那么油门（即PWM波占空比）越大，则电调输入的等效电压越大，电机与螺旋桨转速越快，产生的升力越大；每个螺旋桨通过升力与反扭力转化为力矩与反扭矩，在 xyzxyzxyz 三个轴上合成新的力矩，以下做定量分析。

1.1 符号说明

变量	符号	单位
每个电机的油门	σ1,σ2,σ3,σ4\sigma_1, \sigma_2, \sigma_3, \sigma_4σ1,σ2,σ3,σ4	归一化 [0,1][0,1][0,1]
每个螺旋桨产生的升力	F1,F2,F3,F4F_1, F_2, F_3, F_4F1,F2,F3,F4	N\rm NN
每个螺旋桨产生的诱导阻力	f1,f2,f3,f4f_1, f_2, f_3, f_4f1,f2,f3,f4	N\rm NN
每个螺旋桨升力产生的力矩	M1,M2,M3,M4M_1, M_2, M_3, M_4M1,M2,M3,M4	N⋅m\rm N\cdot mN⋅m
每个螺旋桨产生的反扭矩	Mr1,Mr3,Mr3,Mr4M_{r1}, M_{r3}, M_{r3}, M_{r4}Mr1,Mr3,Mr3,Mr4	N⋅m\rm N\cdot mN⋅m
xyzxyzxyz 轴上的合力矩	τx,τy,τz\tau_x, \tau_y, \tau_zτx,τy,τz	N⋅m\rm N\cdot mN⋅m
合升力	TTT	N\rm NN

注：为简洁起见，以下省略 1,2,3,41,2,3,41,2,3,4 下标。

常量	符号	取值	含义
电机时间常数	TmT_mTm	260ms260 \rm ms260ms	电机转速达到稳定转速63.2%所需时间
电机与四旋翼重心距离	ddd	0.225m0.225 \rm m0.225m	450四旋翼，机臂长度一半 225mm
电机转速斜率	CmC_mCm	706.1rad/s706.1 \rm rad/s706.1rad/s	油门从增加1，电机转速增加量
电机转速截距	ϖm\varpi_mϖm	170.47rad/s170.47 \rm rad/s170.47rad/s	油门为0时，电机转速
升力系数	cTc_TcT	1.201×10−5N/(rad/s)21.201\times 10^{-5} \rm N/(rad/s)^21.201×10−5N/(rad/s)2	单个螺旋桨转速增加 1rad/s1\rm rad/s1rad/s，升力增加大小
反扭矩矩系数	cMc_McM	1.574×10−7N⋅m/(rad/s)21.574\times 10^{-7}\rm N \cdot m / (rad/s)^21.574×10−7N⋅m/(rad/s)2	单个螺旋桨转速增加 1rad/s1\rm rad/s1rad/s，扭矩增加大小

1.2 电机模型

1.2.1 原理

给电机的油门越大，电机与螺旋桨转速（以下简称电机转速）越快，油门与电机稳态转速之间接近线性关系。但是给电机一个油门之后，电机并不能立即达到对应的转速，可以把电机近似为一阶系统

Tmdϖ(t)dt+ϖ(t)=Cmσ(t)+ϖm(1.1)T_m \frac {\mathrm{d}\varpi(t)}{\mathrm{d}t} + \varpi(t) = C_{m} \sigma(t) + \varpi_{m} \tag{1.1} Tmdtdϖ(t)+ϖ(t)=Cmσ(t)+ϖm(1.1)

其中，CmC_mCm为电机转速斜率，定义为油门从增加1，电机转速增加量；ϖm\varpi_mϖm 为电机转速截距，定义为油门为0时，电机的转速；Cmσ(t)+ϖmC_{m} \sigma(t) + \varpi_{m}Cmσ(t)+ϖm 代表电机的稳定转速，与油门线性相关。

为方便理解，式(1.1) 可转化为

ϖ˙(t)=dϖ(t)dt=1Tm(Cmσ(t)+ϖm−ϖ(t))\dot \varpi(t) = \frac {\mathrm{d}\varpi(t)}{\mathrm{d}t} = \frac{1}{T_m} (C_{m} \sigma(t) + \varpi_{m} - \varpi(t)) ϖ˙(t)=dtdϖ(t)=Tm1(Cmσ(t)+ϖm−ϖ(t))

可见，随着电机转速 ϖ(t)\varpi(t)ϖ(t) 的增加，电机转速增量 ϖ˙(t)\dot \varpi(t)ϖ˙(t) 逐渐减小，电机转速逐渐达到给定转速。这个过程中，如果电机时间常数 TmT_mTm 越大，则每次转速的增量越小，达到稳定转速的时间越长。

1.2.2 仿真

假设四个电机的输入油门分别时 0.7,0.6,0.5,0.40.7, 0.6, 0.5, 0.40.7,0.6,0.5,0.4，观察 2s2\rm s2s 内四个电机转速变化。

%% 模型测试
global dt Tm Cm varpimdt = 1e-3;              % 仿真时间步长
Cm = 706.01;            % 油门增大1，电机转速变化(RPM)
varpim = 170.47;        % 零占空比时电机转速(RPM)
Tm = 0.260;             % 电机时间常数N = 2000;
t = 0:dt:dt*(N-1);
sigma = [0.7; 0.6; 0.5; 0.4];
varpi = zeros(N, 4);k=1;
for tt=0:dt:(N-2)*dtk = k+1;% 动力单元模型varpi(k, 1) = motor(sigma(1), varpi(k-1, 1));       % 电机1转速varpi(k, 2) = motor(sigma(2), varpi(k-1, 2));       % 电机2转速varpi(k, 3) = motor(sigma(3), varpi(k-1, 3));       % 电机3转速varpi(k, 4) = motor(sigma(4), varpi(k-1, 4));       % 电机4转速
endfigure(1);plot(t, varpi(:,1), 'LineWidth', 1.5); hold on
plot(t, varpi(:,2), 'LineWidth', 1.5);
plot(t, varpi(:,3), 'LineWidth', 1.5);
plot(t, varpi(:,4), 'LineWidth', 1.5); hold offlegend(['\sigma_1=' num2str(sigma(1))], ['\sigma_2=' num2str(sigma(2))],['\sigma_3=' num2str(sigma(3))],['\sigma_4=' num2str(sigma(4))]);
xlabel('时间 t (s)');ylabel('转速 \varpi  (rad/s)');title('电机模型测试'); grid on; grid minor%% 电机模型
% 输入：油门大小 sigma（0-1）
%       电机上一时刻的转速（rad/s）
% 输出：此时刻电机转速（rad/s）function varpi = motor(sigma, varpi_)global dt Tm Cm varpim;dvarpi = (Cm * sigma + varpim - varpi_) / Tm * dt;varpi = varpi_ + dvarpi;
end

可见，油门越大，电机稳定转速越快，电机在 260ms260 \rm ms260ms 达到稳定转速的63.2% 。

1.3 动力合成模型

1.3.1 合升力

螺旋桨转动产生升力，该升力和螺旋桨转速的平方成正比

f=cTϖ2(1.2)f = c_T \varpi ^2 \tag{1.2} f=cTϖ2(1.2)

其中，cTc_TcT为升力系数。

四个螺旋桨产生的总升力大小为

T=cT(ϖ12+ϖ22+ϖ32+ϖ42)(1.3)T = c_T(\varpi_1^2 + \varpi_2^2 + \varpi_3^2 + \varpi_4^2) \tag{1.3} T=cT(ϖ12+ϖ22+ϖ32+ϖ42)(1.3)

1.3.2 x,yx,yx,y 轴力矩

螺旋桨产生的升力作用在四旋翼上，会导致机体转动，也即产生了力矩。力矩的定义为

M=d×F(1.4)M = d \times F \tag{1.4} M=d×F(1.4)

其中，ddd 为支点到力作用点的矢量，FFF 为力矢量；MMM为力矩，单位 N⋅m\rm N\cdot mN⋅m，其大小为 ∣d∣⋅∣F∣⋅sin⁡θ|d|\cdot |F| \cdot \sin\theta∣d∣⋅∣F∣⋅sinθ，方向服从右手法则垂直于 ddd 与 FFF 所在的平面。

如下图所示，OOO为四旋翼重心，F1F_1F1 为1号螺旋桨产生的升力，垂直纸面向外，容易得出力矩

M1=d×F1(1.5)M_1 = d\times F_1 \tag{1.5} M1=d×F1(1.5)

M1M_1M1 在机体系 x,yx,yx,y 轴上的分量分别为

M1x=−22dF1=−22dcTϖ12(1.6)M_{1x} = -\frac{\sqrt {2}}{2} dF_1 =-\frac{\sqrt {2}}{2} d c_T \varpi_1^2 \tag{1.6}M1x=−22dF1=−22dcTϖ12(1.6)

M1y=22dF1=22dcTϖ12(1.7)M_{1y} = \frac{\sqrt {2}}{2} dF_1 =\frac{\sqrt {2}}{2} d c_T \varpi_1^2 \tag{1.7}M1y=22dF1=22dcTϖ12(1.7)

规定：刚体绕轴转动，力矩为正，刚体绕轴逆时针转动，力矩为负，刚体绕轴为顺时针转动（将坐标系箭头朝向自己，判断顺逆时针）。

1号螺旋桨产生力矩在 xxx 轴上分量为负，导致机体绕 xxx 轴逆时针转动，也即滚转角减小；在 yyy 轴上分量为正，机体系绕 yyy 逆时针转动，俯仰角增大。事实上，1号螺旋桨拉力很大将导致机体右前方升高，即俯仰角增大，滚转角减小，因此上面的规定是合理的。

同理可得其他三个螺旋桨产生的力矩 M2,M3,M4M_2, M_3, M_4M2,M3,M4，通过力矩的矢量合成，易得 x,yx,yx,y 方向的合力矩为

τx=22dcT(−ϖ12+ϖ22+ϖ32−ϖ42)(1.8)\tau_x = \frac{\sqrt {2}}{2} d c_T (-\varpi_1^2 + \varpi_2^2 + \varpi_3^2 - \varpi_4^2) \tag{1.8} τx=22dcT(−ϖ12+ϖ22+ϖ32−ϖ42)(1.8)

τy=22dcT(ϖ12+ϖ22−ϖ32−ϖ42)(1.9)\tau_y = \frac{\sqrt {2}}{2} d c_T (\varpi_1^2 + \varpi_2^2 - \varpi_3^2 - \varpi_4^2) \tag{1.9} τy=22dcT(ϖ12+ϖ22−ϖ32−ϖ42)(1.9)

1.3.3 zzz轴力矩

螺旋桨在转动的时候，迅速冲撞着螺旋桨平面内的空气。根据牛顿第三定律，空气也会给螺旋桨一个反方向的阻力。如下图所示，俯视逆时针转动的1号螺旋桨，螺旋桨转动方向为 vvv，将受到空气的阻力 f1f_1f1，产生力矩 Mr1M_{r1}Mr1

Mr1=r×f1(1.10)M_{r1} = r \times f_1 \tag{1.10} Mr1=r×f1(1.10)

其中，螺旋桨各点所受空气阻力并不相等，f1f_1f1 为等效空气阻力，rrr 为螺旋桨长度的一半。

由于该力矩使得无人机转动，是反作用力效果，故称该力矩为反扭矩，该空气阻力称为反扭力。根据右手螺旋定则容易得出，叉乘后的方向垂直纸面向内，也即 zzz 轴分量为正，这将导致四旋翼逆时针转动（注意将 zzz 轴箭头朝向自己），俯视而看，四旋翼将顺时针转动。

同理，Mr3M_{r3}Mr3 也将导致四旋翼顺时针转动（俯视），偏航角增大；而Mr2,Mr4M_{r2}, M_{r4}Mr2,Mr4 将导致四旋翼逆时针转动（俯视），偏航角减小。

以上说明反扭矩的来源及对四旋翼偏航角的影响，但根据式(1.10) 很难计算的反扭矩，因为平均的反扭力 f1f_1f1 难以得出，因此，采用实验方式得出反扭矩。实验表明，反扭矩也和螺旋桨转速的平方成正比

Mr=cMϖ2(1.11)M_{r} = c_M \varpi^2 \tag{1.11} Mr=cMϖ2(1.11)

其中，cMc_McM 称为反扭矩系数，代表单个螺旋桨转速增加 1rad/s1\rm rad/s1rad/s，反扭矩增加的大小。

因此，四个螺旋桨产生的反扭矩之和为

τz=cM(ϖ12−ϖ22+ϖ32−ϖ42)(1.12)\tau_z = c_M (\varpi_1^2 - \varpi_2^2 + \varpi_3^2 - \varpi_4^2) \tag{1.12} τz=cM(ϖ12−ϖ22+ϖ32−ϖ42)(1.12)

1.3.4 模型小结与仿真

综上，得到完整的动力合成模型如下

T=cT(ϖ12+ϖ22+ϖ32+ϖ42)τx=22dcT(−ϖ12+ϖ22+ϖ32−ϖ42)τy=22dcT(ϖ12+ϖ22−ϖ32−ϖ42)τz=cM(ϖ12−ϖ22+ϖ32−ϖ42)\begin{aligned} &T = c_T(\varpi_1^2 + \varpi_2^2 + \varpi_3^2 + \varpi_4^2) \\ &\tau_x = \frac{\sqrt {2}}{2} d c_T (-\varpi_1^2 + \varpi_2^2 + \varpi_3^2 - \varpi_4^2) \\ &\tau_y = \frac{\sqrt {2}}{2} d c_T (\varpi_1^2 + \varpi_2^2 - \varpi_3^2 - \varpi_4^2) \\ &\tau_z = c_M (\varpi_1^2 - \varpi_2^2 + \varpi_3^2 - \varpi_4^2) \\ \end{aligned}T=cT(ϖ12+ϖ22+ϖ32+ϖ42)τx=22dcT(−ϖ12+ϖ22+ϖ32−ϖ42)τy=22dcT(ϖ12+ϖ22−ϖ32−ϖ42)τz=cM(ϖ12−ϖ22+ϖ32−ϖ42)

假设四个电机的输入油门分别时 0.7,0.6,0.5,0.40.7, 0.6, 0.5, 0.40.7,0.6,0.5,0.4，计算 2s2\rm s2s 内四个螺旋桨产生的合升力及力矩大小。

%% 模型测试
global dt Tm Cm varpim d cT cMdt = 1e-3;              % 仿真时间步长
Cm = 706.01;            % 油门增大1，电机转速变化(RPM)
varpim = 170.47;        % 零占空比时电机转速(RPM)
Tm = 0.260;             % 电机时间常数
d = 0.225;              % 450mm/2
cT = 1.201e-5;          % 升力系数
cM = 1.574e-7;          % 反扭力系数N = 2000;
t = 0:dt:dt*(N-1);
sigma = [0.7; 0.6; 0.5; 0.4];
varpi = zeros(N, 4);
T = zeros(N, 1);
tau = zeros(N, 3);k=1;
for tt=0:dt:(N-2)*dtk = k+1;% 电机模型varpi(k, 1) = motor(sigma(1), varpi(k-1, 1));       % 电机1转速varpi(k, 2) = motor(sigma(2), varpi(k-1, 2));       % 电机2转速varpi(k, 3) = motor(sigma(3), varpi(k-1, 3));       % 电机3转速varpi(k, 4) = motor(sigma(4), varpi(k-1, 4));       % 电机4转速[T(k), tau(k,:)] = power_mix(varpi(k, :));
endfigure(1);subplot(211); plot(t, T, 'linewidth', 1.5); title('动力合成模型');ylabel('升力 (N)');
subplot(212);plot(t, tau(:,1), 'linewidth', 1.5);hold on
plot(t, tau(:,2),'linewidth', 1.5);plot(t, tau(:,3),'linewidth', 1.5);hold off
ylabel('力矩 (N\cdotm)');xlabel('时间 (t)'); legend('\tau_x', '\tau_y', '\tau_z');%% 电机模型
% 输入：油门大小 sigma（0-1）
%       电机上一时刻的转速（rad/s）
% 输出：此时刻电机转速（rad/s）function varpi = motor(sigma, varpi_)global dt Tm Cm varpim;dvarpi = (Cm * sigma + varpim - varpi_) / Tm * dt;varpi = varpi_ + dvarpi;
end%% 动力合成模型
% 输入：四个电机转速
% 输出：合升力与三轴力矩
function [T, tau] = power_mix(varpi)global cT cM d;T = cT * sum(varpi.^2);tau(1) = sqrt(2)/2 * d * cT * (-varpi(1)^2 + varpi(2)^2 + varpi(3)^2 - varpi(4)^2);tau(2) = sqrt(2)/2 * d * cT * ( varpi(1)^2 + varpi(2)^2 - varpi(3)^2 - varpi(4)^2);tau(3) = cM * (varpi(1)^2 - varpi(2)^2 + varpi(3)^2 - varpi(4)^2);
end

可见，产生了大约1.5kg 的拉力，一般450四旋翼空载质量也大约1.5kg，产生这样的拉力是合理的。绕 yyy 轴的力矩很大，并且是正数，说明这会使得无人机绕 yyy 逆时针运动，俯仰角增大。实际上，由于1号2号电机油门很大，相应的螺旋桨也会产生很大拉力，是的无人机俯仰角增大，与仿真结果一致。同理可分析绕 xxx 轴力矩较小，且会让滚转角稍微变小。同时，由于反扭矩系数远小于升力系数，绕 zzz 的力矩会比绕 x,yx,yx,y 的力矩小一到两个数量级。这说明相比于改变四旋翼俯仰角或滚转角，改变无人机偏航角更为困难，需要电机油门变化量更大，在设计控制系统时应注意。

1.4 动力模型小结

动力模型假设了电池电压不变，建立了四旋翼机体系下升力与力矩与四个电机油门的关系。完整模型如下

Tmϖ˙1(t)+ϖ1(t)=Cmσ1(t)+ϖmTmϖ˙2(t)+ϖ2(t)=Cmσ2(t)+ϖmTmϖ˙3(t)+ϖ3(t)=Cmσ3(t)+ϖmTmϖ˙4(t)+ϖ4(t)=Cmσ4(t)+ϖmT=cT(ϖ12+ϖ22+ϖ32+ϖ42)τx=2/2dcT(−ϖ12+ϖ22+ϖ32−ϖ42)τy=2/2dcT(ϖ12+ϖ22−ϖ32−ϖ42)τz=cM(ϖ12−ϖ22+ϖ32−ϖ42)\begin{aligned} &T_m \dot\varpi_1(t) + \varpi_1(t) = C_{m} \sigma_1(t) + \varpi_{m} \\ &T_m \dot\varpi_2(t) + \varpi_2(t) = C_{m} \sigma_2(t) + \varpi_{m} \\ &T_m \dot\varpi_3(t) + \varpi_3(t) = C_{m} \sigma_3(t) + \varpi_{m} \\ &T_m \dot\varpi_4(t) + \varpi_4(t) = C_{m} \sigma_4(t) + \varpi_{m} \\ &T = c_T(\varpi_1^2 + \varpi_2^2 + \varpi_3^2 + \varpi_4^2) \\ &\tau_x = \sqrt 2/2 d c_T (-\varpi_1^2 + \varpi_2^2 + \varpi_3^2 - \varpi_4^2) \\ &\tau_y = \sqrt 2/2 d c_T (\varpi_1^2 + \varpi_2^2 - \varpi_3^2 - \varpi_4^2) \\ &\tau_z = c_M (\varpi_1^2 - \varpi_2^2 + \varpi_3^2 - \varpi_4^2) \\ \end{aligned}Tmϖ˙1(t)+ϖ1(t)=Cmσ1(t)+ϖmTmϖ˙2(t)+ϖ2(t)=Cmσ2(t)+ϖmTmϖ˙3(t)+ϖ3(t)=Cmσ3(t)+ϖmTmϖ˙4(t)+ϖ4(t)=Cmσ4(t)+ϖmT=cT(ϖ12+ϖ22+ϖ32+ϖ42)τx=2/2dcT(−ϖ12+ϖ22+ϖ32−ϖ42)τy=2/2dcT(ϖ12+ϖ22−ϖ32−ϖ42)τz=cM(ϖ12−ϖ22+ϖ32−ϖ42)

专门针对四旋翼建模与控制出了一期教程哦，欢迎购买支持，课程链接火力教育四旋翼建模与控制