用 MATLAB 绘制汽车的各种动力特性曲线

汽车理论这门课讲到汽车的动力性，书上一堆曲线图看着挺唬人，但其实都是从几个基础公式推出来的。正好课上有一道大题，给了一辆中型货车的全部参数，要求把几张特性曲线都画出来。我干脆全用 MATLAB 算，把整个流程从参数到图像走一遍，记下来。

题目

一辆中型货车，已知数据如下：

参数数值
总质量 m9300 kg
车轮滚动半径 r0.49 m
传动系机械效率（直接挡）0.9
传动系机械效率（其他挡）0.85
滚动阻力系数 f0.015
空气阻力因素 Cd·A4 m²
主减速器传动比 i₀6.33

各挡传动比和旋转质量换算系数：

挡位12345
传动比7.314.312.451.541
旋转质量换算系数3.171.771.271.121.07

发动机外特性（转速 n 对应转矩 Ttq）：

n (r/min)120014001600180020002200240026002800
Ttq (N·m)326334.5343340.5338332326.5313300

要求：

  1. 计算并绘制驱动力与行驶阻力平衡图，求最大车速和最大爬坡度
  2. 计算并绘制直接挡加速时间曲线
  3. 计算并绘制动力特性图
  4. 计算并绘制功率平衡图

第一步：把参数存进 workspace

m=9300; r=0.49; efficiency_direct=0.9; efficiency_other=0.85;
f=0.015; cd_a=4; i_0=6.33; Ff=m*9.8*f;

gear_ratio(1)=7.31; gear_ratio(2)=4.31; gear_ratio(3)=2.45; gear_ratio(4)=1.54; gear_ratio(5)=1;
gear_i(1)=3.17; gear_i(2)=1.77; gear_i(3)=1.27; gear_i(4)=1.12; gear_i(5)=1.07;

engine_n(1)=1200; engine_n(2)=1400; engine_n(3)=1600; engine_n(4)=1800; engine_n(5)=2000;
engine_n(6)=2200; engine_n(7)=2400; engine_n(8)=2600; engine_n(9)=2800;

engine_T(1)=326; engine_T(2)=334.5; engine_T(3)=343; engine_T(4)=340.5; engine_T(5)=338;
engine_T(6)=332; engine_T(7)=326.5; engine_T(8)=313; engine_T(9)=300;

density_air=1.293;

第二步：把发动机外特性插值加密

发动机外特性只给了 9 个点，直接拿来积分（后面算加速时间要用）太稀疏了，先用样条插值把它加密成连续曲线。插值后的转速记为 nn，转矩记为 T_n：

nn=(1200:1:2800);
T_n=spline(engine_n,engine_T,nn);

第三步：逐项计算各个量

接下来按公式一个一个算：驱动力、车速、空气阻力、加速度、发动机功率、爬坡角度、加速时间、动力因数。

驱动力 Ft，公式：

Ft=Ttq ig i0 ηTrF_t = \frac{T_{tq}\, i_g\, i_0\, \eta_T}{r}
Ft_1_int=T_n*i_0*0.85*gear_ratio(1)/r;
Ft_2_int=T_n*i_0*0.85*gear_ratio(2)/r;
Ft_3_int=T_n*i_0*0.85*gear_ratio(3)/r;
Ft_4_int=T_n*i_0*0.85*gear_ratio(4)/r;
Ft_5_int=T_n*i_0*0.9 *gear_ratio(5)/r;   % 5 挡是直接挡，效率取 0.9

车速 ua，公式：

ua=0.377 n rig i0u_a = 0.377\, \frac{n\, r}{i_g\, i_0}
ua_1_int=nn*r*0.377/i_0/gear_ratio(1);
ua_2_int=nn*r*0.377/i_0/gear_ratio(2);
ua_3_int=nn*r*0.377/i_0/gear_ratio(3);
ua_4_int=nn*r*0.377/i_0/gear_ratio(4);
ua_5_int=nn*r*0.377/i_0/gear_ratio(5);

空气阻力 Fw，公式：

Fw=CDA ua221.15F_w = \frac{C_D A\, u_a^{2}}{21.15}
Fw_1_int=cd_a*ua_1_int.^2/21.15;
Fw_2_int=cd_a*ua_2_int.^2/21.15;
Fw_3_int=cd_a*ua_3_int.^2/21.15;
Fw_4_int=cd_a*ua_4_int.^2/21.15;
Fw_5_int=cd_a*ua_5_int.^2/21.15;

加速度 acc，公式：

a=dudt=1δm[Ft−(Ff+Fw)]a = \frac{du}{dt} = \frac{1}{\delta m}\left[F_t - (F_f + F_w)\right]
acc_1=(Ft_1_int-Ff-Fw_1_int)/m/3.17;
acc_2=(Ft_2_int-Ff-Fw_2_int)/m/1.77;
acc_3=(Ft_3_int-Ff-Fw_3_int)/m/1.27;
acc_4=(Ft_4_int-Ff-Fw_4_int)/m/1.12;
acc_5=(Ft_5_int-Ff-Fw_5_int)/m/1.07;

发动机功率 Pe，公式：

Pe=1ηT(Gfua3600+Gi ua3600+CDA ua376140+δm ua3600dudt)P_e = \frac{1}{\eta_T}\left( \frac{G f u_a}{3600} + \frac{G i\, u_a}{3600} + \frac{C_D A\, u_a^{3}}{76140} + \frac{\delta m\, u_a}{3600}\frac{du}{dt} \right)
for j=1:1601
    Pe_1(j)=(m*9.8*0.015*ua_1_int(j)/3600+cd_a*ua_1_int(j)^3/76140+3.17*m*ua_1_int(j)/3600*acc_1(j))/efficiency_other;
    Pe_2(j)=(m*9.8*0.015*ua_2_int(j)/3600+cd_a*ua_2_int(j)^3/76140+1.77*m*ua_2_int(j)/3600*acc_2(j))/efficiency_other;
    Pe_3(j)=(m*9.8*0.015*ua_3_int(j)/3600+cd_a*ua_3_int(j)^3/76140+1.27*m*ua_3_int(j)/3600*acc_3(j))/efficiency_other;
    Pe_4(j)=(m*9.8*0.015*ua_4_int(j)/3600+cd_a*ua_4_int(j)^3/76140+1.12*m*ua_4_int(j)/3600*acc_4(j))/efficiency_other;
    Pe_5(j)=(m*9.8*0.015*ua_5_int(j)/3600+cd_a*ua_5_int(j)^3/76140+1.07*m*ua_5_int(j)/3600*acc_5(j))/efficiency_direct;
end

爬坡角度和爬坡度，公式：

α=arcsin⁡Ft−(Ff+Fw)G\alpha = \arcsin \frac{F_t - (F_f + F_w)}{G}
angle(1,:)=asin((Ft_1_int-Ff-Fw_1_int)/m/9.8);
angle(2,:)=asin((Ft_2_int-Ff-Fw_2_int)/m/9.8);
angle(3,:)=asin((Ft_3_int-Ff-Fw_3_int)/m/9.8);
angle(4,:)=asin((Ft_4_int-Ff-Fw_4_int)/m/9.8);
angle(5,:)=asin((Ft_5_int-Ff-Fw_5_int)/m/9.8);

加速时间 t，用微分的方法逐点累加，相当于一阶龙格-库塔。公式：

t(i)=t(i−1)+ua(i)−ua(i−1)a(i−1)t(i) = t(i-1) + \frac{u_a(i) - u_a(i-1)}{a(i-1)}
t(1)=0;
for i=2:1601
    t(i)=(ua_5_int(i)-ua_5_int(i-1))/acc_5(i-1)+t(i-1);
end

再算各挡总阻力，以及一条连续的速度-阻力曲线（用来在平衡图上画那条阻力线）：

F_1_resist=Fw_1_int+Ff;
F_2_resist=Fw_2_int+Ff;
F_3_resist=Fw_3_int+Ff;
F_4_resist=Fw_4_int+Ff;
F_5_resist=Fw_5_int+Ff;

ua_continue=[0:1:100];
F_resist_continue=Ff+1/2*density_air*ua_continue.^2*cd_a/3.6/3.6;

动力因数 D，公式：

D=ψ+δ dug dtD = \psi + \frac{\delta\, du}{g\, dt}
D(1,:)=f+3.17/9.8*acc_1;
D(2,:)=f+1.77/9.8*acc_2;
D(3,:)=f+1.27/9.8*acc_3;
D(4,:)=f+1.12/9.8*acc_4;
D(5,:)=f+1.07/9.8*acc_5;

第四步：画图

把上面算出来的量画成四类曲线：

% 驱动力-行驶阻力平衡图
plot(ua_1_int,Ft_1_int,ua_2_int,Ft_2_int,ua_3_int,Ft_3_int,ua_4_int,Ft_4_int,ua_5_int,Ft_5_int,ua_continue,F_resist_continue);
xlabel('ua(km/h)')；ylabel('F/N')；title（'汽车的驱动力-行驶阻力平衡图'）

% 最大爬坡度
plot(ua_1_int,tan(angle(1,:))*100,ua_2_int,tan(angle(2,:))*100,ua_3_int,tan(angle(3,:))*100,ua_4_int,tan(angle(4,:))*100,ua_5_int,tan(angle(5,:))*100);
xlabel('ua(km/h)')；ylabel('i(%)')；title（'汽车的爬坡度图'）

% 直接挡加速时间曲线
plot(t，ua_5_int)；xlabel('t/s')；ylabel('ua km/h')；title（'直接挡时候的加速时间曲线'）

% 动力特性图
f_continue=0.015*ones(1,101);
plot(ua_1_int,D(1,:),ua_2_int,D(2,:),ua_3_int,D(3,:),ua_4_int,D(4,:),ua_5_int,D(5,:),ua_continue,f_continue);
xlabel('ua(km/h)')；ylabel('D')；title（'汽车的动力特性图'）

% 功率平衡图
plot(ua_1_int,Pe_1,ua_2_int,Pe_2,ua_3_int,Pe_3,ua_4_int,Pe_4,ua_5_int,Pe_5);
xlabel('ua km/h')；ylabel('Pe/kW')；title（'汽车功率平衡图'）

结果

从图上读出来：

  • 最高车速 81.71 km/h（驱动力曲线和阻力曲线在最高挡的交点）
  • 最大爬坡度 30%（1 挡爬坡度曲线的峰值）

五张图如下：

驱动力-行驶阻力平衡图（各挡驱动力曲线与那条向上翘的阻力曲线，交点就是最高车速）：

各挡驱动力曲线与行驶阻力曲线的平衡图
图 1：驱动力-行驶阻力平衡图，各挡驱动力曲线与向上翘的阻力曲线的交点即对应最高车速

爬坡度图（1 挡峰值约 30%，就是最大爬坡度）：

各挡爬坡度随车速变化的曲线
图 2：爬坡度图，1 挡曲线峰值约 30%，即整车的最大爬坡度

直接挡加速时间曲线：

直接挡加速时间随车速变化的曲线
图 3：直接挡加速时间曲线，由加速度倒数对车速积分得到

动力特性图（动力因数 D，底下那条水平线是滚动阻力系数 f）：

各挡动力因数曲线与滚动阻力系数水平线
图 4：动力特性图，纵轴为动力因数 D，底部水平线为滚动阻力系数 f

功率平衡图：

发动机功率曲线与阻力消耗功率曲线的平衡图
图 5：功率平衡图，发动机输出功率曲线与各项阻力消耗功率曲线的交点对应最高车速

整道题做下来，最大的体会是：汽车的动力性其实就被几个公式管着，发动机外特性 + 传动比 + 阻力模型，剩下的都是把它们组合起来积分、求交点。MATLAB 在这里就是个趁手的计算器，省去了手算插值和积分的麻烦。

  • MATLAB
  • 汽车理论
  • 车辆工程