汽車理論課后習題Matlab程序
1.3
確定一輕型貨車的動力性能（貨車可裝用4擋或5擋變速器，任選
其中的一種進行整車性能計算）：
1）繪制汽車驅(qū)動力與行駛阻力平衡圖。
2）求汽車最高車速，最大爬坡度及克服該坡度時相應的附著率。
3）繪制汽車行駛加速度倒數(shù)曲線，用圖解積分法求汽車用2檔起步加速行駛至70km/h的車速－時間曲線，或者用計算機求汽車用2檔起步加速行駛至70km/h的加速時間。
輕型貨車的有關(guān)數(shù)據(jù)：
汽油發(fā)動機使用外特性的Tq-n曲線的擬合公式為
式中，Tq為發(fā)動機轉(zhuǎn)矩（N?m）;n為發(fā)動機轉(zhuǎn)速（r/min）。
發(fā)動機的最低轉(zhuǎn)速nmin=600r/min,最高轉(zhuǎn)速nmax=4000r/min。
裝載質(zhì)量
2000kg
整車整備質(zhì)量
1800kg
總質(zhì)量
3880kg
車輪半徑
0.367m
傳動系機械效率
ηt=0.85
滾動阻力系數(shù)
f=0.013
空氣阻力系數(shù)×迎風面積
CDA=2.77m2
主減速器傳動比
i0=5.83
飛輪轉(zhuǎn)動慣量
If=0.218kg?m2
二前輪轉(zhuǎn)動慣量
Iw1=1.798kg?m2
四后輪轉(zhuǎn)動慣量
Iw2=3.598kg?m2
變速器傳動比
ig(數(shù)據(jù)如下表)
Ⅰ檔
Ⅱ檔
Ⅲ檔
Ⅳ檔
Ⅴ檔
四檔變速器
6.09
3.09
1.71
1.00
五檔變速器
5.56
2.769
1.644
1.00
0.793
軸距
L=3.2m
質(zhì)心至前軸距離（滿載）
a=1.974m
質(zhì)心高（滿載）
hg=0.9m
解：Matlab程序：
(1)
求汽車驅(qū)動力與行駛阻力平衡圖和汽車最高車速程序：
n=[600:10:4000];
Tq=-19.313+295.27*(n/1000)-165.44*(n/1000).^2+40.874*(n/1000).^3-3.8445*(n/1000).^4;
m=3880;g=9.8;nmin=600;nmax=4000;
G=m*g;
ig=[5.56
2.769
1.644
1.00
0.793];nT=0.85;r=0.367;f=0.013;CDA=2.77;i0=5.83;
L=3.2;a=1.947;hg=0.9;If=0.218;Iw1=1.798;Iw2=3.598;
Ft1=Tq*ig(1)*i0*nT/r;
Ft2=Tq*ig(2)*i0*nT/r;
Ft3=Tq*ig(3)*i0*nT/r;
Ft4=Tq*ig(4)*i0*nT/r;
Ft5=Tq*ig(5)*i0*nT/r;
ua1=0.377*r*n/ig(1)/i0;
ua2=0.377*r*n/ig(2)/i0;
ua3=0.377*r*n/ig(3)/i0;
ua4=0.377*r*n/ig(4)/i0;
ua5=0.377*r*n/ig(5)/i0;
ua=[0:5:120];
Ff=G*f;
Fw=CDA*ua.^2/21.15;
Fz=Ff+Fw;
plot(ua1,Ft1,ua2,Ft2,ua3,Ft3,ua4,Ft4,ua5,Ft5,ua,Fz);
title(＇驅(qū)動力-行駛阻力平衡圖＇);
xlabel(＇ua(km/s)＇);
ylabel(＇Ft(N)＇);
gtext(＇Ft1＇),gtext(＇Ft2＇),gtext(＇Ft3＇),gtext(＇Ft4＇),gtext(＇Ft5＇),gtext(＇Ff+Fw＇);
zoom
on;
[x,y]=ginput(1);
zoom
off;
disp(＇汽車最高車速=＇);disp(x);disp(＇km/h＇);
汽車最高車速=
99.3006
km/h
(2)求汽車最大爬坡度程序：
n=[600:10:4000];
Tq=-19.313+295.27*(n/1000)-165.44*(n/1000).^2+40.874*(n/1000).^3-3.8445*(n/1000).^4;
m=3880;g=9.8;nmin=600;nmax=4000;
G=m*g;
ig=[5.56
2.769
1.644
1.00
0.793];nT=0.85;r=0.367;f=0.013;CDA=2.77;i0=5.83;
L=3.2;a=1.947;hg=0.9;If=0.218;Iw1=1.798;Iw2=3.598;
Ft1=Tq*ig(1)*i0*nT/r;
ua1=0.377*r*n/ig(1)/i0;
Ff=G*f;
Fw1=CDA*ua1.^2/21.15;
Fz1=Ff+Fw1;
Fi1=Ft1-Fz1;
Zoom
on;
imax=100*tan(asin(max(Fi1/G)));
disp(＇汽車最大爬坡度=＇);
disp(imax);
disp(＇%＇);
汽車最大爬坡度=
35.2197%
(3)求最大爬坡度相應的附著率和求汽車行駛加速度倒數(shù)曲線程序：
clear
n=[600:10:4000];
Tq=-19.313+295.27*(n/1000)-165.44*(n/1000).^2+40.874*(n/1000).^3-3.8445*(n/1000).^4;
m=3880;g=9.8;nmin=600;nmax=4000;
G=m*g;
ig=[5.56
2.769
1.644
1.00
0.793];nT=0.85;r=0.367;f=0.013;CDA=2.77;i0=5.83;
L=3.2;a=1.947;hg=0.9;If=0.218;Iw1=1.798;Iw2=3.598;
Ft1=Tq*ig(1)*i0*nT/r;
Ft2=Tq*ig(2)*i0*nT/r;
Ft3=Tq*ig(3)*i0*nT/r;
Ft4=Tq*ig(4)*i0*nT/r;
Ft5=Tq*ig(5)*i0*nT/r;
ua1=0.377*r*n/ig(1)/i0;
ua2=0.377*r*n/ig(2)/i0;
ua3=0.377*r*n/ig(3)/i0;
ua4=0.377*r*n/ig(4)/i0;
ua5=0.377*r*n/ig(5)/i0;
Fw1=CDA*ua1.^2/21.15;
Fw2=CDA*ua2.^2/21.15;
Fw3=CDA*ua3.^2/21.15;
Fw4=CDA*ua4.^2/21.15;
Fw5=CDA*ua5.^2/21.15;
Ff=G*f;
deta1=1+(Iw1+Iw2)/(m*r^2)+(If*ig(1)^2*i0^2*nT)/(m*r^2);
deta2=1+(Iw1+Iw2)/(m*r^2)+(If*ig(2)^2*i0^2*nT)/(m*r^2);
deta3=1+(Iw1+Iw2)/(m*r^2)+(If*ig(3)^2*i0^2*nT)/(m*r^2);
deta4=1+(Iw1+Iw2)/(m*r^2)+(If*ig(4)^2*i0^2*nT)/(m*r^2);
deta5=1+(Iw1+Iw2)/(m*r^2)+(If*ig(5)^2*i0^2*nT)/(m*r^2);
a1=(Ft1-Ff-Fw1)/(deta1*m);ad1=1./a1;
a2=(Ft2-Ff-Fw2)/(deta2*m);ad2=1./a2;
a3=(Ft3-Ff-Fw3)/(deta3*m);ad3=1./a3;
a4=(Ft4-Ff-Fw4)/(deta4*m);ad4=1./a4;
a5=(Ft5-Ff-Fw5)/(deta5*m);ad5=1./a5;
plot(ua1,ad1,ua2,ad2,ua3,ad3,ua4,ad4,ua5,ad5);
axis([0
0
10]);
title(＇汽車的加速度倒數(shù)曲線＇);
xlabel(＇ua(km/h)＇);
ylabel(＇1/a＇);
gtext(＇1/a1＇);gtext(＇1/a2＇);gtext(＇1/a3＇);gtext(＇1/a4＇);gtext(＇1/a5＇);
a=max(a1);
af=asin(max(Ft1-Ff-Fw1)/G);
C=tan(af)/(a/L+hg*tan(af)/L);
disp(＇假設后輪驅(qū)動，最大爬坡度相應的附著率=＇);
disp(C);
假設后輪驅(qū)動，最大爬坡度相應的附著率=
0.4219
(4)
>>clear
nT=0.85;r=0.367;f=0.013;CDA=2.77;i0=5.83;If=0.218;
Iw1=1.798;Iw2=3.598;L=3.2;a=1.947;hg=0.9;m=3880;g=9.8;
G=m*g;
ig=[5.56
2.769
1.644
1.00
0.793];
nmin=600;nmax=4000;
u1=0.377*r*nmin./ig/i0;
u2=0.377*r*nmax./ig/i0;
deta=0*ig;
for
i=1:5
deta(i)=1+(Iw1+Iw2)/(m*r^2)+(If*(ig(i))^2*i0^2*nT)/(m*r^2);
end
ua=[6:0.01:99];N=length(ua);n=0;Tq=0;Ft=0;inv_a=0*ua;delta=0*ua;
Ff=G*f;
Fw=CDA*ua.^2/21.15;
for
i=1:N
k=i;
if
ua(i)<=u2(2)
n=ua(i)*(ig(2)*i0/r)/0.377;
Tq=-19.313+295.27*(n/1000)-165.44*(n/1000)^2+40.874*(n/1000)^3-3.8445*(n/1000)^4;
Ft=Tq*ig(2)*i0*nT/r;
inv_a(i)=(deta(2)*m)/(Ft-Ff-Fw(i));
delta(i)=0.01*inv_a(i)/3.6;
elseif
ua(i)<=u2(3)
n=ua(i)*(ig(3)*i0/r)/0.377;
Tq=-19.313+295.27*(n/1000)-165.44*(n/1000)^2+40.874*(n/1000)^3-3.8445*(n/1000)^4;
Ft=Tq*ig(3)*i0*nT/r;
inv_a(i)=(deta(3)*m)/(Ft-Ff-Fw(i));
delta(i)=0.01*inv_a(i)/3.6;
elseif
ua(i)<=u2(4)
n=ua(i)*(ig(4)*i0/r)/0.377;
Tq=-19.313+295.27*(n/1000)-165.44*(n/1000)^2+40.874*(n/1000)^3-3.8445*(n/1000)^4;
Ft=Tq*ig(4)*i0*nT/r;
inv_a(i)=(deta(4)*m)/(Ft-Ff-Fw(i));
delta(i)=0.01*inv_a(i)/3.6;
else
n=ua(i)*(ig(5)*i0/r)/0.377;
Tq=-19.313+295.27*(n/1000)-165.44*(n/1000)^2+40.874*(n/1000)^3-3.8445*(n/1000)^4;
Ft=Tq*ig(5)*i0*nT/r;
inv_a(i)=(deta(5)*m)/(Ft-Ff-Fw(i));
delta(i)=0.01*inv_a(i)/3.6;
end
a=delta(1:k);
t(i)=sum(a);
end
plot(t,ua);
axis([0
0
100]);
title(＇汽車2檔原地起步換擋加速時間曲線＇);
xlabel(＇時間t（s）＇);
ylabel(＇速度ua（km/h）＇);
>>
ginput
ans
=
25.8223
70.0737
25.7467
70.0737
所以汽車2檔原地起步換擋加速行駛至70km/h的加速時間約為25.8s
2.7已知貨車裝用汽油發(fā)動機的負荷特性與萬有特性。負荷特性曲線的擬合公式為：
其中，b為燃油消耗率[g/(kW?h)]；Pe為發(fā)動機凈功率（kW）；擬合式中的系數(shù)隨轉(zhuǎn)速n變化。怠速油耗（怠速轉(zhuǎn)速400r/min）。
計算與繪制題1.3中貨車的1）汽車功率平衡圖。
2）最高檔與次高檔的等速百公里油耗曲線?；蚶糜嬎銠C求貨車按JB3352-83規(guī)定的六工況循環(huán)行駛的百公里油耗。計算中確定燃油消耗值b時，若發(fā)動機轉(zhuǎn)速與負荷特性中給定的轉(zhuǎn)速不相等，可由相鄰轉(zhuǎn)速的兩根曲線用插值法求得。
解：Matlab程序：
(1)
汽車功率平衡圖程序：
clear
n=[600:10:4000];
Tq=-19.313+295.27*(n/1000)-165.44*(n/1000).^2+40.874*(n/1000).^3-3.8445*(n/1000).^4;
m=3880;g=9.8;
G=m*g;
ig=[5.56
2.769
1.644
1.00
0.793];
nT=0.85;r=0.367;f=0.013;CDA=2.77;i0=5.83;
L=3.2;a=1.947;hg=0.9;If=0.218;Iw1=1.798;Iw2=3.598;
ua1=0.377*r*n/ig(1)/i0;
ua2=0.377*r*n/ig(2)/i0;
ua3=0.377*r*n/ig(3)/i0;
ua4=0.377*r*n/ig(4)/i0;
ua5=0.377*r*n/ig(5)/i0;
Pe1=Tq.*ig(1)*i0.*ua1./(3600*r);
Pe2=Tq.*ig(2)*i0.*ua2./(3600*r);
Pe3=Tq.*ig(3)*i0.*ua3./(3600*r);
Pe4=Tq.*ig(4)*i0.*ua4./(3600*r);
Pe5=Tq.*ig(5)*i0.*ua5./(3600*r);
ua=[0:0.35:119];
Ff=G*f;
Fw=CDA*ua.^2/21.15;
Pf=Ff*ua/3600;
Pw=Fw.*ua/3600;
Pe0=(Pf+Pw)./nT;
Pe=max(Pe1);
plot(ua1,Pe1,ua2,Pe2,ua3,Pe3,ua4,Pe4,ua5,Pe5,ua,Pe0,ua,Pe);
axis([0
119
0
100]);
title(＇汽車功率平衡圖＇);
xlabel(＇ua(km/h)＇);
ylabel(＇Pe(kw)＇);
gtext(＇1＇),gtext(＇2＇),gtext(＇3＇),gtext(＇4＇),gtext(＇5＇),gtext(＇(Pf+Pw)/et＇),gtext(＇Pe＇);
（2）最高檔與次高檔的等速百公里油耗曲線程序：
clear
n=600:1:4000;
m=3880;g=9.8;
G=m*g;
ig=[5.56
2.769
1.644
1.00
0.793];
nT=0.85;r=0.367;f=0.013;CDA=2.77;i0=5.83;
L=3.2;a=1.947;hg=0.9;If=0.218;Iw1=1.798;Iw2=3.598;
n0=[815
1207
1614
2012
2603
3006
3403
3804];
B00=[1326.8
1354.7
1284.4
1122.9
1141.0
1051.2
1233.9
1129.7];
B10=[-416.46
-303.98
-189.75
-121.59
-98.893
-73.714
-84.478
-45.291];
B20=[72.379
36.657
14.524
7.0035
4.4763
2.8593
2.9788
0.71113];
B30=[-5.8629
-2.0553
-0.51184
-0.18517
-0.091077
-0.05138
-0.047449
-0.00075215];
B40=[0.17768
0.043072
0.0068164
0.0018555
0.00068906
0.00035032
0.00028230
-0.000038568];
B0=spline(n0,B00,n);
B1=spline(n0,B10,n);
B2=spline(n0,B20,n);
B3=spline(n0,B30,n);
B4=spline(n0,B40,n);
Ff=G*f;
ua4=0.377*r*n/ig(4)/i0;
ua5=0.377*r*n/ig(5)/i0;
Fz4=Ff+CDA*(ua4.^2)/21.15;
Fz5=Ff+CDA*(ua5.^2)/21.15;
Pe4=Fz4.*ua4./(nT*3.6*1000);
Pe5=Fz5.*ua5./(nT*3.6*1000);
for
i=1:1:3401
b4(i)=B0(i)+B1(i)*Pe4(i)+B2(i)*Pe4(i).^2+B3(i)*Pe4(i).^3+B4(i)*Pe4(i).^4;
b5(i)=B0(i)+B1(i)*Pe5(i)+B2(i)*Pe5(i).^2+B3(i)*Pe5(i).^3+B4(i)*Pe5(i).^4;
end
pg=7.0;
Q4=Pe4.*b4./(1.02.*ua4.*pg);
Q5=Pe5.*b5./(1.02.*ua5.*pg);
plot(ua4,Q4,ua5,Q5);
axis([0
30]);
title(＇最高檔與次高檔等速百公里油耗曲線＇);
xlabel(＇ua(km/h)＇);
ylabel(＇百公里油耗（L/100km）＇);
gtext(＇4＇),gtext(＇5＇);
3.1改變1.3題中輕型貨車的主減速器傳動比，做出為5.17、5.43、5.83、6.17、6.33時的燃油經(jīng)濟性—加速時間曲線，討論不同值對汽車性能的影響。
解：Matlab程序：
主程序：
i0=[5.17,5.43,5.83,6.17,6.33];
%輸入主傳動比的數(shù)據(jù)
for
i=1:1:5
y(i)=jiasushijian(i0(i));
%求加速時間
end
y;
for
i=1:1:5
b(i)=youhao(i0(i));
%求對應i0的六工況百公里油耗
end
b;
plot(b,y,＇+r＇)
hold
on
b1=linspace(b(1),b(5),100);
y1=spline(b,y,b1);
%三次樣條插值
plot(b1,y1);
%繪制燃油經(jīng)濟性-加速時間曲線
title(＇燃油經(jīng)濟性—加速時間曲線＇);
xlabel(＇百公里油耗(L/100km)＇);
ylabel(＇加速時間s＇);
gtext(＇i0=5.17＇),gtext(＇i0=5.43＇),gtext(＇i0=5.83＇),gtext(＇i0=6.17＇),gtext(＇i0=6.33＇);
子程序：
(1)
function
y=jiasushijian(i0)
%求加速時間的處理函數(shù)
n1=linspace(0,5000);
%先求各個檔位的驅(qū)動力
nmax=4000;nmin=600;r=0.367;yita=0.85;CDA=2.77;f=0.013;G=(3880)*9.8;ig=[6.09,3.09,1.71,1.00];%i0=5.83
for
i=1:1:4
%i為檔數(shù)
uamax(i)=chesu(nmax,r,ig(i),i0);
%計算各個檔位的最大速度與最小速度
uamin(i)=chesu(nmin,r,ig(i),i0);
ua(i,:)=linspace(uamin(i),uamax(i),100);
n(i,:)=zhuansu(ua(i,:),r,ig(i),i0);
%計算各個檔位的轉(zhuǎn)速范圍
Ttq(i,:)=zhuanju(n(i,:));
%求出各檔位的轉(zhuǎn)矩范圍
Ft(i,:)=qudongli(Ttq(i,:),ig(i),i0,yita,r);
%求出驅(qū)動力
F(i,:)=f*G+CDA*(ua(i,:).^2)/21.15;
%求出滾動阻力和空氣阻力的和
delta(i,:)=1+(1.798+3.598+0.218*(ig(i)^2)*(i0^2)*yita)/(3880*r^2);
%求轉(zhuǎn)動質(zhì)量換算系數(shù)
a(i,:)=1./(delta(i,:).*3880./(Ft(i,:)-F(i,:)));
%求出加速度
F2(i,:)=Ft(i,:)-F(i,:);
end
%下面分各個檔位進行積分，求出加速時間
temp1(1,:)=ua(2,:)/3.6;
temp1(2,:)=1./a(2,:);
n1=1;
for
j1=1:1:100
if
ua(3,j1)>max(ua(2,:))&&ua(3,j1)<=70
temp2(1,n1)=ua(3,j1)/3.6;
temp2(2,n1)=1./a(3,j1);
n1=n1+1;
end
end
n2=1;
for
j1=1:1:100
if
ua(4,j1)>max(ua(3,:))&&ua(4,j1)<=70;
temp3(1,n2)=ua(4,j1)/3.6;
temp3(2,n2)=1./a(4,j1);
n2=n2+1;
end
end
y=temp1(1,1)*temp1(2,1)+qiuji(temp1(1,:),temp1(2,:))+qiuji(temp2(1,:),temp2(2,:))+qiuji(temp3(1,:),temp3(2,:));
end
(2)
function
ua=chesu(n,r,ig,i0);
%由轉(zhuǎn)速計算車速
ua=0.377*r.*n/(ig*i0);
(3)
function
n=zhuansu(ua,r,ig,i0);
%求轉(zhuǎn)速
n=ig*i0.*ua./(0.377*r);
end
(4)
function
y=zhuanju(n);
%求轉(zhuǎn)矩函數(shù)
y=-19.313+295.27.*(n./1000)-165.44.*(n./1000).^2+40.874.*(n./1000).^3-3.8445.*(n./1000).^4;
(5)
function
y=qudongli(Ttq,ig,i0,yita,r);
%求驅(qū)動力函數(shù)
y=(ig*i0*yita.*Ttq)/r;
end
(6)
function
p=qiuji(x0,y0)
%求積分函數(shù)
n0=size(x0);
n=n0(2);
x=linspace(x0(1),x0(n),200)
;
y=spline(x0,y0,x);
%插值
%
figure;plot(x,y);
p=trapz(x,y)
;
end
(7)
%求不同i0下的六工況油耗
function
b=youhao(i0);
global
f
G
CDA
yita
m
r
If
Iw1
Iw2
pg
B0
B1
B2
B3
B4
n
%聲明全局變量
ig=[6.09,3.09,1.71,1.00];r=0.367;
yita=0.85;CDA=2.77;f=0.013;%i0=5.83;
G=(3880)*9.8;If=0.218;Iw1=1.798;Iw2=3.598;m=3880;
%汽車的基本參數(shù)設定
n0=[815
1207
1614
2012
2603
3006
3403
3804];
B00=[1326.8
1354.7
1284.4
1122.9
1141.0
1051.2
1233.9
1129.7];
B10=[-416.46
-303.98
-189.75
-121.59
-98.893
-73.714
-84.478
-45.291];
B20=[72.379
36.657
14.524
7.0035
4.4763
2.8593
2.9788
0.71113];
B30=[-5.8629
-2.0553
-0.51184
-0.18517
-0.091077
-0.05138
-0.047449
-0.00075215];
B40=[0.17768
0.043072
0.0068164
0.0018555
0.00068906
0.00035032
0.00028230
-0.000038568];
n=600:1:4000;
B0=spline(n0,B00,n);
B1=spline(n0,B10,n);
B2=spline(n0,B20,n);
%使用三次樣條插值，保證曲線的光滑連續(xù)
B3=spline(n0,B30,n);
B4=spline(n0,B40,n);
ua4=0.377*r.*n./(i0*ig(4));
%求出發(fā)動機轉(zhuǎn)速范圍內(nèi)對應的III、IV檔車速
F4=f*G+CDA*(ua4.^2)/21.15;
%求出滾動阻力和空氣阻力的和
P_fw4=F4.*ua4./(yita*3.6*1000);
%求出阻力功率
for
i=1:1:3401
%用擬合公式求出各個燃油消耗率
b4(i)=B0(i)+B1(i)*P_fw4(i)+B2(i)*(P_fw4(i))^2+B3(i)*(P_fw4(i))^3+B4(i)*(P_fw4(i))^4;
end
pg=7.06;
%汽油的重度取7.06N/L
ua4_m=[25,40,50];
%勻速階段的車速
s_m=[50,250,250];
%每段勻速走過的距離
b4_m=spline(ua4,b4,ua4_m);
%插值得出對應速度的燃油消耗率
F4_m=f*G+CDA*(ua4_m.^2)/21.15;
%車速對應的阻力
P_fw4_m=F4_m.*ua4_m./(yita*3.6*1000);
%發(fā)動機功率
Q4_m=P_fw4_m.*b4_m.*s_m./(102.*ua4_m.*pg)
;
Q4_a1=jiasu(40,25,ig(4),0.25,ua4,i0);
Q4_a2=jiasu(50,40,ig(4),0.2,ua4,i0);
Qid=0.299;tid=19.3;s=1075;
Q_i=Qid*tid;
%求出減速階段的燃油消耗量
Q4all=(sum(Q4_m)+Q4_a1+Q4_a2+Q_i)*100/s;
%IV檔六工況百公里燃油消耗量
b=Q4all;
(8)加速階段處理函數(shù)
function
q=jiasu(umax,umin,ig,a,ua0,i0);
global
f
G
CDA
yita
m
r
If
Iw1
Iw2
pg
B0
B1
B2
B3
B4
n;
%i0
;
ua1=umin:1:umax;
%把速度范圍以1km/h為間隔進行劃分
delta=1+(Iw1+Iw2)/(m*r^2)+(If*ig^2*i0^2*yita)/(m*r^2);
P0=(G*f.*ua0./3600+CDA.*ua0.^3/76140+(delta*m.*ua0/3600)*a)/yita;
P=(G*f.*ua1/3600+CDA.*ua1.^3/76140+(delta*m.*ua1/3600)*a)/yita;
dt=1/(3.6*a)
;
%速度每增加1km/h所需要的時間
for
i=1:1:3401
%重新利用擬合公式求出b與ua的關(guān)系
b0(i)=B0(i)+B1(i)*P0(i)+B2(i)*(P0(i))^2+B3(i)*(P0(i))^3+B4(i)*(P0(i))^4;
end
b1=interp1(ua0,b0,ua1);
%插值出各個速度節(jié)點的燃油消耗率
Qt=P.*b1./(367.1.*pg);
%求出各個速度節(jié)點的燃油消耗率
i1=size(Qt);
i=i1(2);
Qt1=Qt(2:i-1);
q=(Qt(1)+Qt(i))*dt./2+sum(Qt1)*dt;
%求該加速階段的燃油消耗量
4.3一中型貨車裝有前后制動器分開的雙管路制動系，其有關(guān)參數(shù)如下：
載荷
質(zhì)量（kg）
質(zhì)心高hg/m
軸距L/m
質(zhì)心至前軸距離a/m
制動力分配系數(shù)β
空載
4080
0.845
3.950
2.100
0.38
滿載
9290
1.170
3.950
2.950
0.38
1）
計算并繪制利用附著系數(shù)曲線和制動效率曲線
2）
求行駛車速Ua＝30km/h，在＝0.80路面上車輪不抱死的制動距離。計算時取制動系反應時間＝0.02s，制動減速度上升時間＝0.02s。
3）
求制動系前部管路損壞時汽車的制動距離s，制動系后部管路損壞時汽車的制動距離。
解：Matlab程序：
(1)
求利用附著系數(shù)曲線和制動效率曲線程序：
clear
k=4080;hgk=0.845;Lk=3.950;ak=2.10;betak=0.38;bk=Lk-ak;%空載時的參數(shù)
mm=9290;hgm=1.170;Lm=3.950;am=2.950;betam=0.38;bm=Lm-am;%滿載時的參數(shù)
z=0:0.01:1.0;
figure(1);
fai=z;
fai_fk=betak*z*Lk./(bk+z*hgk);%空載時前軸的φf
fai_fm=betam*z*Lm./(bm+z*hgm);%滿載時前軸的φf
fai_rk=(1-betak)*z*Lk./(ak-z*hgk);%空載時后軸的φr
fai_rm=(1-betam)*z*Lm./(am-z*hgm);%滿載時后軸的φr
plot(z,fai_fk,＇b--＇,z,fai_fm,＇r＇,z,fai_rk,＇b--＇,z,fai_rm,＇r＇,z,fai,＇k＇);
title(＇利用附著系數(shù)與制動強度的關(guān)系曲線＇);
xlabel(＇制動強度(z/g)＇);
ylabel(＇利用附著系數(shù)φ＇);
gtext(＇φr(空載)＇),gtext(＇φr(滿載)＇),gtext(＇φ=z＇),gtext(＇φf(空載)＇),gtext(＇φf(滿載)＇);
figure(2);
Efk=z./fai_fk*100;%空載時前軸的制動效率
Efm=z./fai_fm*100;
Erk=z./fai_rk*100;
Erm=z./fai_rm*100;
plot(fai_fk,Efk,＇b＇,fai_fm,Efm,＇r＇,fai_rk,Erk,＇b＇,fai_rm,Erm,＇r＇);
axis([0
0
100]);
title(＇前.后制動效率曲線＇);
xlabel(＇附著系數(shù)φ＇);
ylabel(＇制動效率%＇);
gtext(＇Ef＇),gtext(＇Er＇),gtext(＇Er＇),gtext(＇滿載＇),gtext(＇空載＇);
(2)
問和(3)問程序：
clear
mk=4080;hgk=0.845;Lk=3.950;ak=2.10;betak=0.38;bk=Lk-ak;%空載時的參數(shù)
mm=9290;hgm=1.170;Lm=3.950;am=2.950;betam=0.38;bm=Lm-am;%滿載時的參數(shù)
z=0:0.01:1;
fai_fk=betak*z*Lk./(bk+z*hgk);%空載時前軸的φf
fai_fm=betam*z*Lm./(bm+z*hgm);%滿載時前軸的φf
fai_rk=(1-betak)*z*Lk./(ak-z*hgk);%空載時后軸的φr
fai_rm=(1-betam)*z*Lm./(am-z*hgm);%滿載時后軸的φr
Efk=z./fai_fk*100;%空載時前軸的制動效率
Efm=z./fai_fm*100;
Erk=z./fai_rk*100;
Erm=z./fai_rm*100;
t1=0.02;t2=0.02;ua0=30;fai=0.80;g=9.8;
ak1=Erk(81)*g*fai/100;
am1=Erm(81)*g*fai/100;
Sk1=(t1+t2/2)*ua0/3.6+ua0^2/(25.92*ak1);%制動距離
Sm1=(t1+t2/2)*ua0/3.6+ua0^2/(25.92*am1);
disp(＇空載時，汽車制動距離Sk1=＇);
disp(Sk1);
disp(＇滿載時，汽車制動距離Sm1=＇);
disp(Sm1);
ak2=fai*g*ak/(Lk+fai*hgk);
am2=fai*g*am/(Lm+fai*hgm);
ak3=fai*g*bk/(Lk-fai*hgk);
am3=fai*g*bm/(Lk-fai*hgm);
Sk2=(t1+t2/2)*ua0/3.6+ua0^2/(25.92*ak2);%制動距離
Sm2=(t1+t2/2)*ua0/3.6+ua0^2/(25.92*am2);
Sk3=(t1+t2/2)*ua0/3.6+ua0^2/(25.92*ak3);
Sm3=(t1+t2/2)*ua0/3.6+ua0^2/(25.92*am3);
disp(＇空載時，前制動器損壞，汽車制動距離Sk2=＇);
disp(Sk2);
disp(＇滿載時，前制動器損壞，汽車制動距離Sm2=＇);
disp(Sm2);
disp(＇空載時，后制動器損壞，汽車制動距離Sk3=＇);
disp(Sk3);
disp(＇滿載時，后制動器損壞，汽車制動距離Sm3=＇);
disp(Sm3);
空載時，汽車制動距離Sk1=
7.8668
滿載時，汽車制動距離Sm1=
5.6354
空載時，前制動器損壞，汽車制動距離Sk2=
10.0061
滿載時，前制動器損壞，汽車制動距離Sm2=
7.5854
空載時，后制動器損壞，汽車制動距離Sk3=
8.0879
滿載時，后制動器損壞，汽車制動距離Sm3=
13.5986
5.11二自由度轎車模型的有關(guān)參數(shù)如下：
總質(zhì)量
m=1818.2kg
繞Oz軸轉(zhuǎn)動慣量
軸距
L=3.048m
質(zhì)心至前軸距離
a=1.463m
質(zhì)心至后軸距離
b=1.585m
前輪總側(cè)偏剛度
k1=-62618N/rad
后輪總側(cè)偏剛度
k2=-110185N/rad
轉(zhuǎn)向系總傳動比
i=20
試求：
1)
穩(wěn)定性因數(shù)K、特征車速uch。
2)
穩(wěn)態(tài)橫擺角速度增益曲線、車速u=22.35m/s時的轉(zhuǎn)向靈敏度。
3)
靜態(tài)儲備系數(shù)S.M.，側(cè)向加速度為0.4g時的前、后輪側(cè)偏角絕對值之差與轉(zhuǎn)彎半徑的比值R/R0(R0=15m)。
4)
車速u=30.56m/s時，瞬態(tài)響應的橫擺角速度波動的固有（圓）頻率、阻尼比、反應時間與峰值反應時間
解：Matlab程序：
m=1818.2;Iz=3885;L=3.048;a=1.463;b=1.585;k1=-62618;k2=-110185;
i=20;g=9.8;R0=15;u1=30.56;
K=m*(a/k2-b/k1)/L^2;
Uch=(1/K)^(1/2);%特征車速
disp(＇穩(wěn)定性因數(shù)(s^2/m^2)K=＇);
disp(K);
disp(＇特征車速(m/s)Uch=＇);
disp(Uch);
u=0:0.05:30;
S=u./(L*(1+K*u.^2));%穩(wěn)態(tài)橫擺角速度增益
plot(u,S);
title(＇汽車穩(wěn)態(tài)橫擺角速度增益曲線＇);
xlabel(＇車速u(m/s)＇);
ylabel(＇穩(wěn)態(tài)橫擺角速度增益＇);
disp(＇u=22.35m/s時，轉(zhuǎn)向靈敏度為＇);
disp(S(448));
SM=k2/(k1+k2)-a/L;
ay=0.4*g;
A=K*ay*L;
B=L/R0;
R=L/(B-A);
C=R/R0;%轉(zhuǎn)彎半徑比
disp(＇靜態(tài)儲備系數(shù)S.M.=＇);
disp(SM);
disp(＇側(cè)向加速度為0.4g時前、后輪側(cè)偏角絕對值之差(rad)
a1-a2=＇);
disp(A);
disp(＇側(cè)向加速度為0.4g時轉(zhuǎn)彎半徑比值R/R0=＇);
disp(C);
W0=L/u1*(k1*k2/(m*Iz)*(1+K*u1^2))^(1/2);%固有（圓）頻率
D=(-m*(k1*a^2+k2*b^2)-Iz*(k1+k2))/(2*L*(m*Iz*k1*k2*(1+K*u1^2))^(1/2));%阻尼比
t=atan((1-D^2)^(1/2)/(-m*u1*a*W0/(L*k2)-D))/(W0*(1-D^2)^(1/2));%反應時間
E=atan((1-D^2)^(1/2)/D)/(W0*(1-D^2)^(1/2))+t;%峰值反應時間
disp(＇車速u=30.56m/s時的瞬態(tài)響應參數(shù)分別為：＇);
disp(＇橫擺角速度波動的固有(圓)頻率(rad)為
＇);
disp(W0);
disp(＇阻尼比為＇);
disp(D);
disp(＇反應時間(s)為＇);
disp(t);
disp(＇峰值反應時間(s)為＇);
disp(E);
穩(wěn)定性因數(shù)(s^2/m^2)K=
0.0024
特征車速(m/s)Uch=
20.6053
u=22.35m/s時，轉(zhuǎn)向靈敏度為
3.3690
靜態(tài)儲備系數(shù)S.M.=
0.1576
側(cè)向加速度為0.4g時前、后輪側(cè)偏角絕對值之差(rad)
a1-a2=
0.0281
側(cè)向加速度為0.4g時轉(zhuǎn)彎半徑比值R/R0=
1.1608
車速u=30.56m/s時的瞬態(tài)響應參數(shù)分別為：
橫擺角速度波動的固有(圓)頻率(rad)為
5.5758
阻尼比為
0.5892
反應時間(s)為
0.1811
峰值反應時間(s)為
0.3899
6.5車身-車輪雙質(zhì)量系統(tǒng)參數(shù)：。
“人體-座椅”系統(tǒng)參數(shù)：。車速，路面不平度系數(shù)，參考空間頻率n0=0.1m-1。
計算時頻率步長，計算頻率點數(shù)。
1)
計算并畫出幅頻特性、、和均方根值譜、、譜圖。進一步計算值
2)
改變“人體-座椅”系統(tǒng)參數(shù)：。分析值隨的變化。
3)
分別改變車身-車輪雙質(zhì)量系統(tǒng)參數(shù)：。繪制三個響應量均方根值隨以上四個系統(tǒng)參數(shù)變化的曲線。
解：Matlab程序
(1)問
yps=0.25;%阻尼比ζ
gama=9;%剛度比γ
mu=10;%質(zhì)量比μ
fs=3;ypss=0.25;g=9.8;a0=10^(-6);f0=1.5;
ua=20;Gqn0=2.56*10^(-8);n0=0.1;detaf=0.2;N=180;
f=detaf*[0:N];lamta=f/f0;lamtas=f/fs;Wf=0*f;
deta=((1-lamta.^2).*(1+gama-1/mu*lamta.^2)-1).^2+4*yps^2*lamta.^2.*(gama-(1/mu+1)*lamta.^2).^2;
z1_q=gama*sqrt(((1-lamta.^2).^2+4*yps^2*lamta.^2)./deta);
z2_z1=sqrt((1+4*yps^2*lamta.^2)./((1-lamta.^2).^2+4*yps^2*lamta.^2));
p_z2=sqrt((1+(2*ypss*lamtas).^2)./((1-lamtas.^2).^2+(2*ypss*lamtas).^2));
z2_q=gama*sqrt((1+4*yps^2*lamta.^2)./deta);
p_q=p_z2.*z2_q;
jfg_Gqddf=4*pi^2*sqrt(Gqn0*n0^2*ua)*f;
jfg_Gzdd1f=z1_q.*jfg_Gqddf;
jfg_Gzdd2f=z2_q.*jfg_Gqddf;
jfg_Gaf=p_q.*jfg_Gqddf;
sigmaqdd=sqrt(trapz(f,jfg_Gqddf.^2));%路面不平度加速度均方根值
sigmazdd1=sqrt(trapz(f,jfg_Gzdd1f.^2));%車輪加速度均方根值
sigmazdd2=sqrt(trapz(f,jfg_Gzdd2f.^2));%車身加速度均方根值
sigmaa=sqrt(trapz(f,jfg_Gaf.^2));%人體加速度均方根值
for
i=1:(N+1)
if
f(i)<=2
Wf(i)=0.5;
elseif
f(i)<=4
Wf(i)=f(i)/4;
elseif
f(i)<=12.5
Wf(i)=1;
else
Wf(i)=12.5/f(i);
end
end
kk=Wf.^2.*jfg_Gaf.^2;
aw=sqrt(trapz(f,kk));%加權(quán)加速度均方根值
Law=20*log10(aw/a0);%加權(quán)振級
disp(＇路面不平度加速度均方根值為＇);disp(sigmaqdd);
disp(＇車輪加速度均方根值為＇);disp(sigmazdd1);
disp(＇車身加速度均方根值為＇);disp(sigmazdd2);
disp(＇人體加速度均方根值為＇);disp(sigmaa);
disp(＇加權(quán)加速度均方根值為＇);disp(aw);
disp(＇加權(quán)振級＇);disp(Law);
figure(1)
plot(f,z1_q),title(＇幅頻特性|z1/q|,(f=1.5Hz,ζ=0.25,γ=9,μ=10)＇),xlabel(＇激振頻率f/Hz＇),ylabel(＇|z1/q|＇);
figure(2)
plot(f,z2_z1),title(＇幅頻特性|z2/z1|,(f=1.5Hz,ζ=0.25,γ=9,μ=10)＇),xlabel(＇激振頻率f/Hz＇),ylabel(＇|z2/z1|＇);
figure(3)
plot(f,p_z2),title(＇幅頻特性|p/z2|,(fs=1.5Hz,ζs=0.25)＇),xlabel(＇激振頻率f/Hz＇),ylabel(＇|p/z2|＇);
figure(4)
plot(f,jfg_Gzdd1f),title(＇車輪加速度均方根值√Gz1(f)譜圖＇),xlabel(＇激振頻率f/Hz＇),ylabel(＇√Gz1(f)＇);
figure(5)
plot(f,jfg_Gzdd2f),title(＇車身加速度均方根值√Gz2(f)譜圖＇),xlabel(＇激振頻率f/Hz＇),ylabel(＇√Gz2(f)＇);
figure(6)
plot(f,jfg_Gaf),title(＇人體加速度均方根值√Ga(f)譜圖＇),xlabel(＇激振頻率f/Hz＇),ylabel(＇√Ga(f)＇);
路面不平度加速度均方根值為
0.3523
車輪加速度均方根值為
0.2391
車身加速度均方根值為
0.0168
人體加速度均方根值為
0.0161
加權(quán)加速度均方根值為
0.0100
加權(quán)振級
80.0291
(2)問
程序1：
clear
gama=9;%剛度比γ
mu=10;%質(zhì)量比μ
f0=1.5;g=9.8;a0=10^(-6);ua=20;
Gqn0=2.56*10^(-8);n0=0.1;detaf=0.2;N=180;
f=detaf*[0:N];lamta=f/f0;Wf=0*f;
for
i=1:(N+1)
if
f(i)<=2
Wf(i)=0.5;
elseif
f(i)<=4
Wf(i)=f(i)/4;
elseif
f(i)<=12.5
Wf(i)=1;
else
Wf(i)=12.5/f(i);
end
end
fs=3;ypss=0.25;
ypss0=[0.125:0.005:0.5];a=0*ypss0;La=0*ypss0;
M=length(ypss0);
for
i=1:M
yps=ypss0(i);
lamtas=f/fs;
deta=((1-lamta.^2).*(1+gama-1/mu*lamta.^2)-1).^2+4*yps^2*lamta.^2.*(gama-(1/mu+1)*lamta.^2).^2;
p_z2=sqrt((1+(2*ypss*lamtas).^2)./((1-lamtas.^2).^2+(2*ypss*lamtas).^2));
z2_q=gama*sqrt((1+4*yps^2*lamta.^2)./deta);
p_q=p_z2.*z2_q;
jfg_Gqddf=4*pi^2*sqrt(Gqn0*n0^2*ua)*f;
jfg_Gaf=p_q.*jfg_Gqddf;
kk=Wf.^2.*jfg_Gaf.^2;
aw(i)=sqrt(trapz(f,kk));
end
Law=20*log10(aw/a0);
figure(1)
plot(ypss0,aw);title(＇aw隨ζs的變化＇),xlabel(＇“人體—座椅”系統(tǒng)的阻尼比ζs＇),ylabel(＇aw/m*s^-2＇);
figure(2)
plot(ypss0,Law);title(＇Law隨ζs的變化＇),xlabel(＇“人體—座椅”系統(tǒng)的阻尼比ζs＇),ylabel(＇Law/dB＇);
程序2：
clear
yps=0.25;%阻尼比ζ
gama=9;%剛度比γ
mu=10;%質(zhì)量比μ
f0=1.5;g=9.8;a0=10^(-6);ua=20;
Gqn0=2.56*10^(-8);n0=0.1;detaf=0.2;N=180;
f=detaf*[0:N];lamta=f/f0;Wf=0*f;
for
i=1:(N+1)
if
f(i)<=2
Wf(i)=0.5;
elseif
f(i)<=4
Wf(i)=f(i)/4;
elseif
f(i)<=12.5
Wf(i)=1;
else
Wf(i)=12.5/f(i);
end
end
ypss=0.25;
fs=[1.5:0.025:6];
M=length(fs);
for
i=1:M
fs0=fs(i);
lamtas=f/fs0;
deta=((1-lamta.^2).*(1+gama-1/mu*lamta.^2)-1).^2+4*yps^2*lamta.^2.*(gama-(1/mu+1)*lamta.^2).^2;
p_z2=sqrt((1+(2*ypss*lamtas).^2)./((1-lamtas.^2).^2+(2*ypss*lamtas).^2));
z2_q=gama*sqrt((1+4*yps^2*lamta.^2)./deta);
p_q=p_z2.*z2_q;
jfg_Gqddf=4*pi^2*sqrt(Gqn0*n0^2*ua)*f;
jfg_Gaf=p_q.*jfg_Gqddf;
kk=Wf.^2.*jfg_Gaf.^2;
aw(i)=sqrt(trapz(f,kk));
end
Law=20*log10(aw/a0);
figure(3)
plot(fs,aw);title(＇aw隨fs的變化＇),xlabel(＇“人體—座椅”系統(tǒng)的固有頻率fs＇),ylabel(＇aw/m*s^-2＇);
figure(4)
plot(fs,Law);title(＇Law隨fs的變化＇),xlabel(＇“人體—座椅”系統(tǒng)的固有頻率fs＇),ylabel(＇Law/dB＇);
（3）問
程序1：
clear
figure(1)
fs=3;yps_s=0.25;g=9.8;
ua=20;Gqn0=2.56*10^(-8);n0=0.1;detaf=0.2;N=180;
f0=1.5;yps=0.25;gama=9;mu=10;
ff0=[0.25:0.05:3];sigmaz2=0*ff0;sigmafd=0*ff0;sigmaFd_G=0*ff0;
M=length(ff0);
for
i=1:M
f0=ff0(i);
f=detaf*[0:N];lamta=f/f0;lamtas=f/fs;
deta=((1-lamta.^2).*(1+gama-1/mu*lamta.^2)-1).^2+4*yps^2*lamta.^2.*(gama-(1/mu+1)*lamta.^2).^2;
z2_qdot=2*pi*f*gama.*sqrt((1+4*yps^2*lamta.^2)./deta);
fd_qdot=gama*lamta.^2./(2*pi*f+eps)./sqrt(deta);
Fd_Gqdot=2*pi*f*gama/g.*sqrt(((lamta.^2/(mu+1)-1).^2+4*yps^2*lamta.^2)./deta);
Gq_dotf=4*pi^2*Gqn0*n0^2*ua;
Gz2f=(z2_qdot).^2*Gq_dotf;
Gfd_qf=(fd_qdot).^2*Gq_dotf;
GFd_Gf=(Fd_Gqdot).^2*Gq_dotf;
sigmaz2(i)=sqrt(trapz(f,Gz2f));
sigmafd(i)=sqrt(trapz(f,Gfd_qf));
sigmaFd_G(i)=sqrt(trapz(f,GFd_Gf));
if
f0==1.5
sgmz2=sigmaz2(i);
sgmfd=sigmafd(i);
sgmFd_G=sigmaFd_G(i);
end
end
sz2=20*log10(sigmaz2/sgmz2);
sfd=20*log10(sigmafd/sgmfd);
sFd_G=20*log10(sigmaFd_G/sgmFd_G);
plot(ff0,sz2,＇r-＇,ff0,sfd,＇b-.＇,ff0,sFd_G,＇k--＇);
axis([0.25
15]);
title(＇三個響應量均方根值隨f0變化的曲線＇),xlabel(＇車身部分固有頻率f0/Hz＇),ylabel(＇σz2/dB,σfd/dB,σFd/G/dB＇);
程序2：
clear
figure(2)
fs=3;yps_s=0.25;g=9.8;
ua=20;Gqn0=2.56*10^(-8);n0=0.1;detaf=0.2;N=180;
f0=1.5;yps=0.25;gama=9;mu=10;
c=(0.5-0.125)/180;
yps0=[0.125:c:0.5];sigmaz2=0*yps0;sigmafd=0*yps0;sigmaFd_G=0*yps0;
M=length(yps0);
for
i=1:M
yps=yps0(i);
f=detaf*[0:N];lamta=f/f0;lamtas=f/fs;
deta=((1-lamta.^2).*(1+gama-1/mu*lamta.^2)-1).^2+4*yps^2*lamta.^2.*(gama-(1/mu+1)*lamta.^2).^2;
z2_qdot=2*pi*f*gama.*sqrt((1+4*yps^2*lamta.^2)./deta);
fd_qdot=gama*lamta.^2./(2*pi*f+eps)./sqrt(deta);
Fd_Gqdot=2*pi*f*gama/g.*sqrt(((lamta.^2/(mu+1)-1).^2+4*yps^2*lamta.^2)./deta);
Gq_dotf=4*pi^2*Gqn0*n0^2*ua;
Gz2f=(z2_qdot).^2*Gq_dotf;
Gfd_qf=(fd_qdot).^2*Gq_dotf;
GFd_Gf=(Fd_Gqdot).^2*Gq_dotf;
sigmaz2(i)=sqrt(trapz(f,Gz2f));
sigmafd(i)=sqrt(trapz(f,Gfd_qf));
sigmaFd_G(i)=sqrt(trapz(f,GFd_Gf));
if
yps==0.25
sgmz2=sigmaz2(i);
sgmfd=sigmafd(i);
sgmFd_G=sigmaFd_G(i);
end
end
sz2=20*log10(sigmaz2/sgmz2);
sfd=20*log10(sigmafd/sgmfd);
sFd_G=20*log10(sigmaFd_G/sgmFd_G);
plot(yps0,sz2,＇r-＇,yps0,sfd,＇b-.＇,yps0,sFd_G,＇k--＇);
axis([0.125
0.5
4]);
title(＇三個響應量均方根值隨ζ變化的曲線＇),xlabel(＇車身部分阻尼比ζ＇),ylabel(＇σz2/dB,σfd/dB,σFd/G/dB＇);
程序3：
clear
figure(3)
fs=3;yps_s=0.25;g=9.8;
ua=20;Gqn0=2.56*10^(-8);n0=0.1;detaf=0.2;N=180;
f0=1.5;yps=0.25;mu=10;
gama0=[4:0.1:19];sigmaz2=0*gama0;sigmafd=0*gama0;sigmaFd_G=0*gama0;
M=length(gama0);
for
i=1:M
gama=gama0(i);
f=detaf*[0:N];lamta=f/f0;lamtas=f/fs;
deta=((1-lamta.^2).*(1+gama-1/mu*lamta.^2)-1).^2+4*yps^2*lamta.^2.*(gama-(1/mu+1)*lamta.^2).^2;
z2_qdot=2*pi*f*gama.*sqrt((1+4*yps^2*lamta.^2)./deta);
fd_qdot=gama*lamta.^2./(2*pi*f+eps)./sqrt(deta);
Fd_Gqdot=2*pi*f*gama/g.*sqrt(((lamta.^2/(mu+1)-1).^2+4*yps^2*lamta.^2)./deta);
Gq_dotf=4*pi^2*Gqn0*n0^2*ua;
Gz2f=(z2_qdot).^2*Gq_dotf;
Gfd_qf=(fd_qdot).^2*Gq_dotf;
GFd_Gf=(Fd_Gqdot).^2*Gq_dotf;
sigmaz2(i)=sqrt(trapz(f,Gz2f));
sigmafd(i)=sqrt(trapz(f,Gfd_qf));
sigmaFd_G(i)=sqrt(trapz(f,GFd_Gf));
if
gama==9
sgmz2=sigmaz2(i);
sgmfd=sigmafd(i);
sgmFd_G=sigmaFd_G(i);
end
end
sz2=20*log10(sigmaz2/sgmz2);
sfd=20*log10(sigmafd/sgmfd);
sFd_G=20*log10(sigmaFd_G/sgmFd_G);
plot(gama0,sz2,＇r-＇,gama0,sfd,＇b-.＇,gama0,sFd_G,＇k--＇);
axis([4
6]);
title(＇三個響應量均方根值隨γ變化的曲線＇),xlabel(＇懸架與輪胎的剛度比γ＇),ylabel(＇σz2/dB,σfd/dB,σFd/G/dB＇);
程序4：
clear
figure(4)
fs=3;yps_s=0.25;g=9.8;
ua=20;Gqn0=2.56*10^(-8);n0=0.1;detaf=0.2;N=180;
f0=1.5;yps=0.25;gama=9;
mu0=[5:0.1:20];sigmaz2=0*mu0;sigmafd=0*mu0;sigmaFd_G=0*mu0;
M=length(mu0);
for
i=1:M
mu=mu0(i);
f=detaf*[0:N];lamta=f/f0;lamtas=f/fs;
deta=((1-lamta.^2).*(1+gama-1/mu*lamta.^2)-1).^2+4*yps^2*lamta.^2.*(gama-(1/mu+1)*lamta.^2).^2;
z2_qdot=2*pi*f*gama.*sqrt((1+4*yps^2*lamta.^2)./deta);
fd_qdot=gama*lamta.^2./(2*pi*f+eps)./sqrt(deta);
Fd_Gqdot=2*pi*f*gama/g.*sqrt(((lamta.^2/(mu+1)-1).^2+4*yps^2*lamta.^2)./deta);
Gq_dotf=4*pi^2*Gqn0*n0^2*ua;
Gz2f=(z2_qdot).^2*Gq_dotf;
Gfd_qf=(fd_qdot).^2*Gq_dotf;
GFd_Gf=(Fd_Gqdot).^2*Gq_dotf;
sigmaz2(i)=sqrt(trapz(f,Gz2f));
sigmafd(i)=sqrt(trapz(f,Gfd_qf));
sigmaFd_G(i)=sqrt(trapz(f,GFd_Gf));
if
mu==10
sgmz2=sigmaz2(i);
sgmfd=sigmafd(i);
sgmFd_G=sigmaFd_G(i);
end
end
sz2=20*log10(sigmaz2/sgmz2);
sfd=20*log10(sigmafd/sgmfd);
sFd_G=20*log10(sigmaFd_G/sgmFd_G);
plot(mu0,sz2,＇r-＇,mu0,sfd,＇b-.＇,mu0,sFd_G,＇k--＇);
axis([5
2]);
title(＇三個響應量均方根值隨μ變化的曲線＇),xlabel(＇車身與車輪部分質(zhì)量比μ＇),ylabel(＇σz2/dB,σfd/dB,σFd/G/dB＇);
點擊咨詢