机动车理论matlab作业.doc

《机动车理论matlab作业.doc》由会员分享，可在线阅读，更多相关《机动车理论matlab作业.doc（34页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、-_一、确定一轻型货车的动力性能。一、确定一轻型货车的动力性能。 1） 绘制汽车驱动力与行驶阻力平衡图；绘制汽车驱动力与行驶阻力平衡图； 2） 求汽车最高车速与最大爬坡度；求汽车最高车速与最大爬坡度； 3） 绘制汽车行驶加速度倒数曲线；用计算机求汽车用绘制汽车行驶加速度倒数曲线；用计算机求汽车用档起档起步加速行驶至步加速行驶至 70km/h 所需所需的加速时间。的加速时间。 已知数据略。已知数据略。 （参见（参见汽车理论汽车理论习题第一章第习题第一章第 3 题）题） 解题程序如下解题程序如下:用用 Matlab 语言语言 (1)绘制汽车驱动力与行驶阻力平衡图绘制汽车驱动力与行驶阻力平衡图m1=

2、2000; m2=1800; mz=3880; g=9.81; r=0.367; CdA=2.77; f=0.013; nT=0.85; ig=5.56 2.769 1.644 1.00 0.793; i0=5.83; If=0.218; Iw1=1.798; Iw2=3.598; Iw=2*Iw1+4*Iw2; for i=1:69 n(i)=(i+11)*50; Ttq(i)=-19.313+295.27*(n(i)/1000)-165.44*(n(i)/1000)2+40.874*(n(i)/1000)3- 3.8445*(n(i)/1000)4; end for j=1:5 for i

3、=1:69 Ft(i,j)=Ttq(i)*ig(j)*i0*nT/r; ua(i,j)=0.377*r*n(i)/(ig(j)*i0); Fz(i,j)=CdA*ua(i,j)2/21.15+mz*g*f; end end plot(ua,Ft,ua,Ff,ua,Ff+Fw) title(汽车驱动力与行驶阻力平衡图); xlabel(ua(km/h); ylabel(Ft(N); gtext(Ft1)-_gtext(Ft2) gtext(Ft3) gtext(Ft4) gtext(Ft5) gtext(Ff+Fw)(2)求最大速度和最大爬坡度求最大速度和最大爬坡度for k=1:175 n1(

4、k)=3300+k*0.1; Ttq(k)=-19.313+295.27*(n1(k)/1000)-165.44*(n1(k)/1000)2 +40.874*(n1(k)/1000)33.8445*(n1(k)/1000)4; Ft(k)=Ttq(k)*ig(5)*i0*nT/r; ua(k)=0.377*r*n1(k)/(ig(5)*i0); Fz(k)=CdA*ua(k)2/21.15+mz*g*f; E(k)=abs(Ft(k)-Fz(k); end for k=1:175 if(E(k)=min(E) disp(汽车最高车速=); disp(ua(k); disp(km/h); end

5、 end for p=1:150 n2(p)=2000+p*0.5; Ttq(p)=-19.313+295.27*(n2(p)/1000)-165.44*(n2(p)/1000)2+40.874*(n2(p)/1000)-_3-3.8445*(n2(p)/1000)4; Ft(p)=Ttq(p)*ig(1)*i0*nT/r; ua(p)=0.377*r*n2(p)/(ig(1)*i0); Fz(p)=CdA*ua(p)2/21.15+mz*g*f; af(p)=asin(Ft(p)-Fz(p)/(mz*g); end for p=1:150 if(af(p)=max(af) i=tan(af(

6、p); disp(汽车最大爬坡度=); disp(i); end end 汽车最高车速汽车最高车速=99.0679km/h 汽车最大爬坡度汽车最大爬坡度=0.3518 （3） 计算计算 2 档起步加速到档起步加速到 70km/h 所需时间所需时间for i=1:69 n(i)=(i+11)*50; Ttq(i)=-19.313+295.27*(n(i)/1000)-165.44*(n(i)/1000)2+40.874*(n(i)/1000)3- 3.8445*(n(i)/1000)4; end for j=1:5 for i=1:69 deta=1+Iw/(mz*r2)+If*ig(j)2*i

7、02*nT/(mz*r2); ua(i,j)=0.377*r*n(i)/(ig(j)*i0); a(i,j)=(Ttq(i)*ig(j)*i0*nT/r-CdA*ua(i,j)2/21.15 -mz*g*f)/(deta*mz); if(a(i,j)0.05) b1(i,j)=a(i,j); u1(i,j)=ua(i,j); else b1(i,j)=a(i-1,j); u1(i,j)=ua(i-1,j); end -_b(i,j)=1/b1(i,j); end end x1=u1(:,1);y1=b(:,1); x2=u1(:,2);y2=b(:,2); x3=u1(:,3);y3=b(:,

8、3); x4=u1(:,4);y4=b(:,4); x5=u1(:,5);y5=b(:,5); plot(x1,y1,x2,y2,x3,y3,x4,y4,x5,y5); title(加速度倒数时间曲线); axis(0 120 0 30); xlabel(ua(km/h); ylabel(1/aj); gtext(1/a1) gtext(1/a2) gtext(1/a3) gtext(1/a4) gtext(1/a5)for i=1:69 A=ua(i,3)-ua(69,2); if (A0) j=i; end -_B=ua(i,4)-ua(69,3); if(B0) k=i; end if(

9、ua(i,4)0.15);g3(i)=z(i)+0.08;endif(z(i)=0.3);g3(i)=0.38+(z(i)-0.3)/0.74;end endz1=0.15:0.01:0.3;g4=z1-0.08;plot(z,gf1,-.,z,gf2,z,gr1,-.,z,gr2,z,g,z,g3,xk,z1,g4,x) axis(0 1 0 1.2) title(利用附着系数与制动强度的关系曲线) xlabel(制动强度z/g)ylabel(利用附着系数g) gtext(空车前轴) gtext(空车后轴) gtext(满载前轴) gtext(满载后轴) gtext(ECE法规)-_C=0:

10、0.05:1; Er1=(a1/L)./(1-beta)+C*hg1/L)*100; Ef=(L-a2)/L./(beta-C*hg2/L)*100; Er=(a2/L)./(1-beta)+C*hg2/L)*100; plot(C,Er,C,Ef,C,Er1) axis(0 1 0 100) title(前后附着效率曲线)xlabel(附着系数C) ylabel(制动效率（%）) gtext(满载)gtext(Ef) gtext(Er) gtext(空载) gtext(Er)-_C1=0.8 E1=(ak1/L)./(1-beta)+C1*hg1/L); E2=(am2/L)/(1-beta

11、)+C1*hg2/L); a1=E1*C1*9.8; a2=E2*C1*9.8; ua=30;i21=0.02;i22=0.02; s1=(i21+i22/2)*ua/3.6+ua2/(25.92*ak1); s2=(i21+i22/2)*ua/3.6+ua2/(25.92*am2);disp(满载时不抱死的制动距离=) disp(s2)disp(空载时不抱死的制动距离=) disp(s1)满载时不抱死的制动距离满载时不抱死的制动距离=5.3319空载时不抱死的制动距离空载时不抱死的制动距离=6.8406beta3=1 beta4=0 Ekr=(a1/L)/(1-beta4)+C1*hg1/L

12、); Ekf=(L-a1)/(beta3*L-C1*hg1); Emf=(L-a2)/L./(beta3-C1*hg2/L); Emr=(a2/L)./(1-beta4)+C1*hg2/L);-_akr=0.8*9.8*Ekr; akf=0.8*9.8*Ekf; amr=0.8*9.8*Emr; amf=0.8*9.8*Emf; skr=(i21+i22/2)*ua/3.6+ua2/(25.92*akr); skf=(i21+i22/2)*ua/3.6+ua2/(25.92*akf); smf=(i21+i22/2)*ua/3.6+ua2/(25.92*amf); smr=(i21+i22/2

13、)*ua/3.6+ua2/(25.92*amr);disp(空车后管路失效时制动距离) disp(skf)disp(空车前管路失效时制动距离) disp(skr)disp(满载后管路失效时制动距离) disp(smf)disp(满载前管路失效时制动距离) disp(smr) 运行结果为：空车后管路失效时制动距离运行结果为：空车后管路失效时制动距离 8.0879 空车前管路失效时制动距离空车前管路失效时制动距离 10.0061 满载后管路失效时制动距离满载后管路失效时制动距离 13.5986 满载前管路失效时制动距离满载前管路失效时制动距离 7.5854五、二自由度轿车模型的有关参数如下：五、二

14、自由度轿车模型的有关参数如下： 总质量总质量 m=1818.2kgm=1818.2kg绕绕 OzOz 轴转动惯量轴转动惯量 23885mkgIz轴距轴距 L=3.048mL=3.048m质心至前轴距离质心至前轴距离 a=1.463ma=1.463m质心至后轴距离质心至后轴距离 b=1.585mb=1.585m前轮总侧偏刚度前轮总侧偏刚度 k k1 1=-62618N/rad=-62618N/rad后轮总侧偏刚度后轮总侧偏刚度 k k2 2=-110185N/rad=-110185N/rad转向系总传动比转向系总传动比 i=20i=20试求：试求：-_1)1)稳定性因数稳定性因数 K K、特征车

15、速、特征车速 u uchch。2)2)稳态横摆角速度增益曲线稳态横摆角速度增益曲线、车速、车速a sru u=22.35m/su=22.35m/s 时的转向灵敏度时的转向灵敏度。swr 3)3)静态储备系数静态储备系数 S.M.S.M.，侧向加速度为，侧向加速度为 0.4g0.4g 时的时的前、后轮侧偏角绝对值之差前、后轮侧偏角绝对值之差与转弯半径的比值与转弯半径的比值21R/RR/R0 0(R(R0 0=15m)=15m)。4)4)车速车速 u=30.56m/su=30.56m/s 时，瞬态响应的横摆角速度波时，瞬态响应的横摆角速度波动的固有（圆）频率动的固有（圆）频率、阻尼比、阻尼比 、反

16、应时间、反应时间 与峰值反与峰值反0应时间应时间MatlabMatlab 程序：程序：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/L2*(a/k2-b/k1); uch=1/K(1/2); u=0:0.05:25; s=u/L./(1+K*u.2);disp(稳定因素 K=); disp(K);disp(特征车速 uch=); disp(uch); plot(u,s); xlabel(ua/(m/s);ylabel(稳态横摆增益); title(汽车

17、的稳态横态摆角速度增益曲线); disp(ua=22.35m/s 时，转向灵敏度为); disp(s(448); SM=k2/(k1+k2)-a/L; A=K*0.4*g*L; Q=L/R0;-_R=L/(Q-A); r=R/R0;disp(静态储备系数 S.M.=); disp(SM);disp(前后轮侧偏角绝对值之差(a1-a2)=); disp(A);disp(转弯半径的比值=); disp(r); w0=L/u1*(k1*k2*(1+K*u12)/(m*Iz)(1/2); n=(-m*(a2*k1+b2*k2)-Iz*(k1+k2)/(2*L*(m*Iz*k1*k2*(1+K*u12)

18、(1/2); t=atan(1-n2)(1/2)/(-m*u1*a*w0/(L*k2)-n)/(w0*(1-n2)(1/2); e=atan(1-n2)(1/2)/n)/(w0*(1-n2)(1/2)+t;disp(横摆角速度波动时的固有频率为); disp(w0);disp(阻尼比为); disp(n);disp(反应时间为); disp(t);disp(达到第一峰值的时间为); disp(e);运算结果：运算结果：稳定因素稳定因素 K= 0.0024特征车速特征车速 uch=-_20.6053ua=22.35m/s 时，转向灵敏度为时，转向灵敏度为 3.3690静态储备系数静态储备系数 S

19、.M.= 0.1576前后轮侧偏角绝对值之差前后轮侧偏角绝对值之差(a1-a2)= 0.0281转弯半径的比值转弯半径的比值= 1.1608横摆角速度波动时的固有频率为横摆角速度波动时的固有频率为 5.5758阻尼比为阻尼比为 0.5892反应时间为反应时间为 0.1811达到第一峰值的时间为达到第一峰值的时间为 0.3899六、车身六、车身- -车轮双质量系统参数：车轮双质量系统参数：。10, 9,25. 0,5 . 10Hzf“人体人体-座椅座椅”系统参数：系统参数：。车速。车速，路，路25. 0,3ssHzfsmu/20面不平度系数面不平度系数，参考空间频率，参考空间频率 n0=0.1m

20、-1。 38 01056. 2mnGq计算时频率步长计算时频率步长，计算频率点数，计算频率点数。Hzf2 . 0180N1) 计算并画出幅频特性计算并画出幅频特性、和均方根值和均方根值qz /112/ zz2/ zq谱谱、谱图。进一步计算谱图。进一步计算 fGz 1 fGz 2 fGa值值awwazzqLa 、 21 2) 改变改变“人体人体-座椅座椅”系统参数：系统参数：。分析。分析值随值随的变化。的变化。5 . 0125. 0,65 . 1ssHzfawwLa 、ssf、3) 分别改变车身分别改变车身- -车轮双质量系统参数：车轮双质量系统参数：，。绘制。绘制5 . 0125. 0,325

21、. 00Hzf205,185 . 4-_三个响应量均方根值随以上四个系统参数变化三个响应量均方根值随以上四个系统参数变化GFdfdz/2、 的曲线。的曲线。Matlab 程序：程序：f0=1.5;g=0.25;r=9;u=10;fs=3;gs=0.25;v=20;Gqn0=2.56e-008;n0=0.1; df=0.2;N=180; f=0:0.2:36; ff0=f/f0; d=(1-ff0.2).*(1+r-1/u*ff0.2-1).2+4*g2*ff0.2.*r-(1/u+1)*ff0.2.2; z1q=r*sqrt(1-ff0.2).2+4*g2*ff0.2)./d); z2z1=s

22、qrt(1+(2*g*ff0).2)./(1-ff0.2).2+(2*g*ff0).2); ffs=f/fs; pz2=sqrt(1+(2*gs*ffs).2)./(1-ffs.2).2+(2*gs*ffs).2) loglog(f,z1q,*,f,z2z1,.,f,pz2); axis(0 100 0 10) set(gca,xtick,0 0.1 1 10 100) set(gca,ytick,0 0.1 1 10)xlabel(激振频率 f/HZ) -_w=2*pi*f; Gz1=w.2.*z1q.*(2*3.14./w).*sqrt(Gqn0*n02*v); loglog(f,Gz1)

23、 hold on Gz2=w.2.*z1q.*z2z1.*(2*3.14./w).*sqrt(Gqn0*n02*v);title(车轮部分加速度均方根值) xlabel(激振频率 f/Hz) ylabel(|Gz1|) loglog(f,Gz2)title(车身部分加速度均方根值) xlabel(激振频率 f/Hz) ylabel(|Gz2|) Gza=w.2.*z1q.*z2z1.*pz2.*(2*3.14./w).*sqrt(Gqn0*n02*v); loglog(f,Gza)title(人体加速度均方根值) xlabel(激振频率 f/Hz) ylabel(|Gza|)-_format

24、long f0=1.5;g=0.25;r=9;u=10;fs=3;gs=0.25;v=20;Gqn0=2.56*10(-8);n0=0.1;-_df=0.2;N=180; f=0.000001:0.2:36; ff0=f/f0; ffs=f/fs; w=2*pi*f; d=(1-ff0.2).*(1+r-ff0.2/u)-1).2+4.*g2.*ff0.2.*(r-(1/u+1).*ff0.2).2; z1q=r*sqrt(1-ff0.2).2+4*g2*ff0.2)./d); z2z1=sqrt(1+(2*g*ff0).2)./(1-ff0.2).2+(2*g*ff0).2); ffs=f/

25、fs; pz2=sqrt(1+(2*gs*ffs).2)./(1-ffs.2).2+(2*gs*ffs).2); ez=sqrt(sum(16.*(3.14)4.*f.2.*Gqn0.*n02.*v*df);disp(路面不平度加速度均方根值=) disp(ez) ez1=sqrt(sum(z1q.2.*16.*pi4.*Gqn0.*v.*f.2*df*n02);disp(车轮加速度均方根值=) disp(ez1) ez2=sqrt(sum(z2z1.2.*z1q.2.*16.*pi4.*Gqn0.*v.*f.2*n02*df);disp(车身加速度均方根值=) disp(ez2) Gza=w

26、.2.*z1q.*z2z1.*pz2.*(2*3.14./w).*sqrt(Gqn0*n02*v); eza=sqrt(sum(Gza.2*df)disp(人体加速度均方根值=) disp(eza) Gza=w.2.*z1q.*z2z1.*pz2.*(2*3.14./w).*sqrt(Gqn0*n02*v);for i=1:180if(f(i)0.5)f1(i)=f(i); Gza(i)=w(i).2.*z1q(i).*z2z1(i).*pz2(i).*(2*3.14./w(i).*sqrt(Gqn0*n02*v);a1(i)=0.5*Gza(i).2*0.2;endif(f(i)2gs2=l

27、inspace(0.125,0.5,10); for j=1:100 f0=1.5;g=0.25;r=9;u=10;fs=3;gs=0.25;v=20;Gqn0=2.56*10(-8);n0=0.1; df=0.2;N=180; f=0.000001:0.2:36; ff0=f/f0; ffs1=f/fs1(j); w=2*pi*f; d=(1-ff0.2).*(1+r-ff0.2/u)-1).2+4.*g2.*ff0.2.*(r-(1/u+1).*ff0.2).2; z1q=r*sqrt(1-ff0.2).2+4*g2*ff0.2)./d); z2z1=sqrt(1+(2*g*ff0).2)

28、./(1-ff0.2).2+(2*g*ff0).2); pz2=sqrt(1+(2*gs*ffs1).2)./(1-ffs1.2).2+(2*gs*ffs1).2); Gza=w.2.*z1q.*z2z1.*pz2.*(2*3.14./w).*sqrt(Gqn0*n02*v); for i=1:180if(f(i)0.5)f1(i)=f(i); Gza(i)=w(i).2.*z1q(i).*z2z1(i).*pz2(i).*(2*3.14./w(i).*sqrt(Gqn0*n02*v);a1(i)=0.5*Gza(i).2*0.2;endif(f(i)2 Gza(i)=w(i).2.*z1q(

29、i).*z2z1(i).*pz2(i).*(2*3.14./w(i).*sqrt(Gqn0*n02*v);a1(i)=0.5*Gza(i).2*0.2;endif(f(i)2a4(i)=12.5/f(i)*Gza(i).2;f4(i)=f(i); end end aw1=trapz(f1,a1); aw2=trapz(f2,a2); aw3=trapz(f3,a3); aw4=trapz(f4,a4); awb(j)=sqrt(aw1+aw2+aw3+aw4); Law2(j)=20*log10(awb(j)*106); endsubplot(2,2,1) plot(fs1,aw) gtext

30、(aw1)ylabel(人体加权加速度均方根值) xlabel(座椅 fs1/HZ) subplot(2,2,2) plot(fs1,Law1)-_gtext(Law1)ylabel(加权振级) xlabel(座椅 fs1/HZ) subplot(2,2,3) plot(gs2,awb) gtext(aw2)ylabel(人体加权加速度均方根值) xlabel(座椅 gs2)subplot(2,2,4) plot(gs2,Law2)ylabel(加权振级) xlabel(座椅 gs2) gtext(Law2)所以人体加权加速度均方根值和加权振级随所以人体加权加速度均方根值和加权振级随 fs1

31、的增大而增大的增大而增大人体加权加速度均方根值和加权振级随人体加权加速度均方根值和加权振级随 gs1 的增大而减小的增大而减小3）f01=linspace(0.25,3,100)for i=1:100 f0=1.5;g=0.25;r=9;u=10;fs=3;gs=0.25;v=20;Gqn0=2.56*10(-8);n0=0.1;-_df=0.2;N=180; f=0.000001:0.2:36; ff0=f/f01(i); ffs1=f/fs; w=2*pi*f; d=(1-ff0.2).*(1+r-ff0.2/u)-1).2+4.*g2.*ff0.2.*(r-(1/u+1).*ff0.2)

32、.2; z1q=r*sqrt(1-ff0.2).2+4*g2*ff0.2)./d); z2z1=sqrt(1+(2*g*ff0).2)./(1-ff0.2).2+(2*g*ff0).2); pz2=sqrt(1+(2*gs*ffs1).2)./(1-ffs1.2).2+(2*gs*ffs1).2); ez2(i)=sqrt(sum(z2z1.2.*z1q.2.*4.*pi2.*Gqn0.*v*n02*n02*df) FdG=(r.*w/9.8).2.*(ff0.2./(1+u)-1).2+4*g2.*ff0.2)./d; eFdG(i)=sqrt(sum(FdG.*4.*pi2.*Gqn0.*

33、v*n02*df); fd=(r*ff0.2./w).2/d; efd(i)=sqrt(sum(fd.*4.*pi2.*Gqn0.*v*n02*df) ; end for i=1:100 Law1(i)=20*log(ez2(i)/ez2(46); Law2(i)=20*log(eFdG(i)/eFdG(46); Law3(i)=20*log(efd(i)/efd(46); end plot(f01,Law1,r,f01,Law2,b,f01,Law3,g)gtext(ez2) gtext(eFdG) gtext(efd) xlabel(f0/Hz) ylabel(efd/eFdG/ez2 d

34、B)title(fd, dFdG,ez2 随 f0 的变化)-_g=linspace(0.125,0.5,100)for i=1:100 f0=1.5;g2=0.25;r=9;u=10;fs=3;gs=0.25;v=20;Gqn0=2.56*10(-8);n0=0.1; df=0.2;N=180; f=0.000001:0.2:36; ff0=f/f0; ffs1=f/fs; w=2*pi*f; d=(1-ff0.2).*(1+r-ff0.2/u)-1).2+4.*g(i)2.*ff0.2.*(r-(1/u+1).*ff0.2).2; z1q=r*sqrt(1-ff0.2).2+4*g(i)2

35、*ff0.2)./d); z2z1=sqrt(1+(2*g(i)*ff0).2)./(1-ff0.2).2+(2*g(i)*ff0).2); pz2=sqrt(1+(2*gs*ffs1).2)./(1-ffs1.2).2+(2*gs*ffs1).2); ez2(i)=sqrt(sum(z2z1.2.*z1q.2.*4*pi2.*Gqn0.*v.*n02*df); FdG=(r.*w/9.8).2.*(ff0.2./(1+u)-1).2+4*g(i)2.*ff0.2)./d; eFdG(i)=sqrt(sum(FdG.*4*pi2.*Gqn0.*v.*n02*df); fd=(r*ff0.2./

36、w).2/d; efd(i)=sqrt(sum(fd.*4*pi2.*Gqn0.*v.*n02*df) ; end for i=1:100 Law1(i)=20*log(ez2(i)/ez2(34); Law2(i)=20*log(eFdG(i)/eFdG(34); Law3(i)=20*log(efd(i)/efd(34); end plot(g,Law1,-.,g,Law2,r,g,Law3) gtext(ez2) gtext(eFdG) gtext(efd)xlabel(阻尼比 g) ylabel(efd/eFdG/ez2 dB)title(fd, dFdG,ez2 随 g 的变化)-_

37、r=linspace(4.5,18,100) for i=1:100 f0=1.5;g=0.25;r1=9;u=10;fs=3;gs=0.25;v=20;Gqn0=2.56*10(-8);n0=0.1; df=0.2;N=180; f=0.000001:0.2:36; ff0=f/f0; ffs1=f/fs; w=2*pi*f; d=(1-ff0.2).*(1+r(i)-ff0.2/u)-1).2+4.*g2.*ff0.2.*(r(i)-(1/u+1).*ff0.2).2; z1q=r(i)*sqrt(1-ff0.2).2+4*g2*ff0.2)./d); z2z1=sqrt(1+(2*g*f

38、f0).2)./(1-ff0.2).2+(2*g*ff0).2); pz2=sqrt(1+(2*gs*ffs1).2)./(1-ffs1.2).2+(2*gs*ffs1).2); ez2(i)=sqrt(sum(z2z1.2.*z1q.2.*14.*pi4.*Gqn0.*v.*n02*df); FdG=(r(i).*w/9.8).2.*(ff0.2./(1+u)-1).2+4*g2.*ff0.2)./d; eFdG(i)=sqrt(sum(FdG.*4.*pi4.*Gqn0.*v.*n02*df); fd=(r(i)*ff0.2./w).2/d; efd(i)=sqrt(sum(fd.*4.*

39、pi4.*Gqn0.*v.*n02*df) ; end for i=1:100 Law1(i)=20*log(ez2(i)/ez2(34); Law2(i)=20*log(eFdG(i)/eFdG(34); Law3(i)=20*log(efd(i)/efd(34); end-_plot(r,Law1,r,r,Law2,b,r,Law3,g)gtext(ez2) gtext(eFdG) gtext(efd) xlabel(r) ylabel(efd/eFdG/ez2 dB)title(fd,edFdG,ez2 随 r 的变化) u=linspace(5,20,151) for i=1:151

40、f0=1.5;g=0.25;r=9;u1=10;fs=3;gs=0.25;v=20;Gqn0=2.56*10(-8);n0=0.1; df=0.2;N=180; f=0.000001:0.2:36; ff0=f/f0; ffs1=f/fs; w=2*pi*f; d=(1-ff0.2).*(1+r-ff0.2/u(i)-1).2+4.*g2.*ff0.2.*(r-(1/u(i)+1).*ff0.2).2; z1q=r*sqrt(1-ff0.2).2+4*g2*ff0.2)./d); z2z1=sqrt(1+(2*g*ff0).2)./(1-ff0.2).2+(2*g*ff0).2); pz2=s

41、qrt(1+(2*gs*ffs1).2)./(1-ffs1.2).2+(2*gs*ffs1).2); ez2(i)=sqrt(sum(z2z1.2.*z1q.2.*4.*pi2.*Gqn0.*v.*n02*df); FdG=(r.*w/9.8).2.*(ff0.2./(1+u(i)-1).2+4*g2.*ff0.2)./d; eFdG(i)=sqrt(sum(FdG.*4.*pi2.*Gqn0.*v.*n02*df); fd=(r*ff0.2./w).2./d;-_efd(i)=sqrt(sum(fd.*4.*pi2.*Gqn0.*v.*n02*df) ; endfor i=1:151 Law1(i)=20*log(ez2(i)/ez2(34); Law2(i)=20*log(eFdG(i)/eFdG(34); Law3(i)=20*log(efd(i)/efd(34); end plot(u,Law1,r,u,Law2,b,u,Law3,g)gtext(ez2) gtext(eFdG) gtext(efd) xlabel(u) ylabel(efd/eFdG/ez2 dB)title(fd, dFdG,ez2 随 u 的变化)