计算机仿真技术试题含完整答案.docx

一, 数值计算，编程完成以下各题（共20分，每小题5分）1, 脉冲宽度为，周期为的矩形脉冲的傅里叶级数如下式描述：当，绘制出函数的图形。解：syms n t;f=(sin(n*pi/4)/(n*pi/4)*cos(2*pi*n*t);s=symsum(f,n,1,150);y=(1+2*s)/4;x=-0.5:0.01:0.5;Y=subs(y,'t',x);plot(x,Y)2, 画出函数在区间3, 5的图形，求出该函数在区间3, 5中的最小值点与函数的最小值.解：程序如下x=3:0.05:5;y=(sin(5*x).2).*exp(0.05*x.2)-5*(x.5).*cos(1.5*x)+1.5*abs(x+5.5)+x.2.5;mix_where=find(y=min(y);xmin=x(mix_where);hold on;plot(x,y);plot(xmin,min(y),'go','linewidth',5);str=strcat('(',num2str(xmin),',',num2str(min(y),')');text(xmin,min(y),str);Xlabel('x') Ylabel('f(x)')经过运行后得到的图像截图如下：运行后的最小值点=4.6，= -8337.8625 3、 画出函数在1，3区间的图形，并用编程求解该非线性方程的一个根，设初始点为.解： x=1:0.02:3;x0=2;y=(x)(cos(x).2).*exp(-0.3*x)-2.5*abs(x); fplot(y,1,3);Xlabel('x') Ylabel('f(x)') X1=fzero('(cos(x).2).*exp(-0.3*x)-2.5*abs(x)',x0)运行后求得该方程的一个根为z=0.3256。4, 已知非线性方程组如下，编程求方程组的解，设初始点为1 0.5 -1.解：%在新建中建立函数文件fun2_4.mfunction f=fun2_4(x)f=x(1).2+x(1)*sqrt(7)+2;x(1)+5*x(3).2-3;x(2).*x(3)+3;%非线性方程组求解主程序fxxfcz.mx0=1 0.5 -1;fsolve(fun2_4,x0)运行后结果为：ans =-1.3229 3.2264 -0.9298 即是 x=-1.3229 y=3.2264 z=-0.9298 .二, 限制系统仿真（15分）某限制系统的开环传递函数为：，要求：编制一个完整的程序完成以下各小题的要求，所绘制的图形分别定义为四张图。1） 绘制出系统的阶跃信号响应曲线（响应时间为）2） 绘制出系统的脉冲信号响应曲线（响应时间为）3） 绘制出系统的斜坡信号响应曲线（响应时间为）4） 绘制出系统的Bode图（要求频率范围为rad/sec）解：由传递函数知，该传递函数是将其用零极点描述法描述的，将其化为用传递函数表述的形式为：，所以num=0 1.08 9.72 6,den=0.3 6.05 1 0。 %用传递函数编程求解 num=0 1.08 9.72 6; den=0.3 6.05 1 0; sys=tf(num,den); t1=0:0.1:30; figure(1) step(sys) %绘制出系统的阶跃信号响应曲线 t2=0:0.1:20; figure(2) impulse(sys) %绘制出系统的脉冲信号响应曲线 t3=0:0.1:10; figure(3) ramp=t3; lsim(sys,ramp,t3);%绘制出系统的斜坡信号响应曲线 figure(4) w=10(-2):102;bode(sys,w);%绘制出系统的Bode图 fig(1)系统的阶跃信号响应曲线fig(2)系统的脉冲信号响应曲线 fig(3)系统的斜坡信号响应曲线 fig(4)系统的Bode图三, 曲线拟合（15分）已知某型号液力变矩器原始特性参数，要求用多项式拟合的方法编程完成以下各小题：1）用二阶多项式拟合出曲线；用三阶多项式拟合出曲线；用三阶多项式拟合出曲线。2）用不同的颜色与不同的线型，将的原始特性参数数据点与二阶拟合曲线绘制在同一张图形中；将的原始特性参数数据点与三阶拟合曲线绘制在同一张图形中；将的原始特性参数数据点与四阶拟合曲线绘制在同一张图形中。3）运行程序，写出曲线的二阶拟合公式, 曲线的三阶拟合公式与曲线的四阶拟合公式。解：% 曲线拟合（Curve fitting）disp('Input Data-i; Output Data-k(i),eta(i),lambdaB(i):')x=0.065,0.098,0.147,0.187,0.243,0.295,0.344,0.398,0.448,0.499;y1=2.37,2.32,2.23,2.15,2.05,1.96,1.87,1.78,1.69,1.59;y2=0.154,0.227,0.327,0.403,0.497,0.576,0.644,0.707,0.757,0.795;y3=26.775,26.845,27.147,27.549,28.052,28.389,28.645,28.756,28.645,28.243;figure(1)pf1=polyfit(x,y1,2)px1=polyval(pf1,x)plot(x,px1,'k')gridxlabel('转速比i')ylabel('变矩比K')title('二阶多项式拟合k曲线')pause figure(2)pf2=polyfit(x,y2,3)px2=polyval(pf2,x)plot(x,px2,'b')gridxlabel('转速比i')ylabel('效率eta')title('三阶多项式拟合eta 曲线')pause figure(3)pf3=polyfit(x,y3,4)px3=polyval(pf3,x)plot(x,px3,'-r')gridxlabel('转速比i')ylabel('泵轮转矩系数lambdaB')title('四阶多项式拟合lambdaB曲线' )figure(4)pf1=polyfit(x,y1,2)px1=polyval(pf1,x)plot(x,y1,'or',x,px1,'k')gridxlabel('转速比i')ylabel('变矩比K')title('二阶多项式拟合k曲线')Legend('原始数据','拟合曲线')%将的原始特性参数数据点与二阶拟合曲线绘制在同一张图形中pause figure(5)pf2=polyfit(x,y2,3)px2=polyval(pf2,x)plot(x,y2,'*m',x,px2,'b')gridxlabel('转速比i')ylabel('效率eta')title('三阶多项式拟合eta 曲线')Legend('原始数据','拟合曲线',0)%将的原始特性参数数据点与三阶拟合曲线绘制在同一张图形中pause figure(6)pf3=polyfit(x,y3,4)px3=polyval(pf3,x)plot(x,y3,'pk',x,px3,'-r')gridxlabel('转速比i')ylabel('泵轮转矩系数lambdaB')title('四阶多项式拟合lambdaB曲线' )Legend('原始数据','拟合曲线',0)%将的原始特性参数数据点与四阶拟合曲线绘制在同一张图形中y1=poly2str(pf1,'x') %曲线的二阶拟合公式y2=poly2str(pf2,'x') %曲线的三阶拟合公式y3=poly2str(pf3,'x') %曲线的四阶拟合公式运行后的结果如下： 运行后的二阶，三阶，四阶拟合曲线函数为：y1 = 0.01325 x2 - 1.8035 x + 2.491y2 =-0.12713 x3 - 1.6598 x2 + 2.4499 x + 0.0025474y3 =106.7407 x4 - 199.9852 x3 + 95.8404 x2 - 8.7272 x + 26.9754四, 微分方程求解。（25分）自己选择确定一个三阶微分方程，自己设置初始条件，用ode45方法求微分方程的解。要求:（例如：，） 1）仿真时间t=30秒2）结果绘制在一张图中，包括曲线，一阶曲线，二阶曲线，三阶曲线3）用图例吩咐分别说明四条曲线为“”，“”，“” ，“”4）定义横坐标为“时间”，纵坐标为“输出”，图形标题名称为“微分方程的解”解：系统方程为 , 这是一个单变量三阶常微分方程。将上式写成一个一阶方程组的形式，这是函数ode45调用规定的格式。 令： 函数文件程序：function ydot=myfun1(t,y)ydot=y(2);y(3);1-8*y(1)-2*y(3)-4*y(2);主文件程序：t=0 30;y0=0;1;0;tt,yy=ode45(myfun1,t,y0);y=(1-yy(:,3)-2*yy(:,2)-4*yy(:,1)/8;plot(tt,y,'r',tt,yy(:,1),'k',tt,yy(:,2),'-g',tt,yy(:,3),'-.b');legend('y-t','y-t','y-t','y-t') title('微分方程的解')xlabel('时间') ylabel('输出') 运行程序后输出图形如下：五, PID设计（25分）自己选定一个限制系统，（例如：某单位负反馈系统的开环传递函数为），设计一个PID限制器，使系统响应满意较快的上升时间与过渡过程时间, 较小的超调量, 静态误差尽可能小。方法要求：用ZieglerNichols方法对三个参数, , 进行整定，并比较PID限制前后的性能，性能的比较要求编程实现（用未加PID限制的系统闭环传递函数阶跃响应与加PID限制后的闭环传递函数的阶跃响应进行比较）解：1) 分析：用ZieglerNichols方法是一种阅历方法，关键是首先通过根轨迹图找出Km与m，然后利用阅历公式求增益，微分，积分时间常数。程序：ng=400;dg=1 30 200 0;rlocus(ng,dg); %画根轨迹图axis(-30 1 -20 20);gridkm,pole=rlocfind(ng,dg)wm=imag(pole(2)kp=0.6*kmkd=kp*pi/(4*wm)ki=kp*wm/pink=kd kp ki,dk=1 0pausend=conv(nk,ng),dd=conv(dk,dg) n1,d1=feedback(ng,dg,1,1)n2,d2=feedback(nd,dd,1,1);%加PID后的闭环传函figurestep(n1,d1,2)gridhold onpausestep(n2,d2,2)hold off在程序中，首先运用rlocus与rlocfind吩咐求出系统穿越增益Km=12.2961与穿越频率m=13.0220rad/s，然后运用ZN方程求出参数。selected_point =-0.4325 +12.9814ikp =7.3777 kd =0.4450 ki =30.5807为采纳PID限制前后的系统闭环阶跃响应状况比较。图6-1系统的根轨迹图 图6-2 PID限制前后的系统闭环阶跃响应三参数KP，Ki，Kd的整定利用系统的等幅振荡曲线的ZieglerNichols方法限制类型限制器的限制参数KpKiKdP0.5Km0PI0.45Km0.54Km/Tm0PID0.6Km1.2Km/Tm0.072Km/Td2) PID限制系统的开环传函为： 因为式中具有积分项，故假如G(s)是n 型系统，加PID限制后系统变为n+1型，可由下式依据给定的稳态误差指标确定参数Ki。因为是个I型系统，由于系统的开环传递函数中有积分项，故为II型系统，假定单位斜坡输入稳态误差，则可以计算出Ki。即：已知系统性能指标为：系统相角裕量PM=80°，增益穿越频率=4rad/s，故利用这两个参数来求Kp，Kd。程序如下：ng=400;dg=1 30 200 0;ki=5;wgc=4;pm=80;ngv=polyval(ng,j*wgc);dgv=polyval(dg,j*wgc);g=ngv/dgv;thetar=(pm-180)*pi/180;ejtheta=cos(thetar)+j*sin(thetar);eqn=(ejtheta/g)+j*(ki/wgc);x=imag(eqn);r=real(eqn);kp=rkd=x/wgcif ki=0 dk=1 0;nk=kd kp ki;else dk=1;nk=kd kp;endpausend=conv(nk,ng),dd=conv(dk,dg)n1,d1=feedback(ng,dg,1,1)n2,d2=feedback(nd,dd,1,1) %加PID限制后的闭环系统传递函数pauseg1m,p1m,wpc1,wgc1=margin(ng,dg)g2m,p2m,wpc2,wgc2=margin(nd,dd) %幅值裕度，相角裕度，相频曲线穿越-180°时的频率，截止频率w=logspace(-1,2,200);pausefigurebode(ng,dg,w)gridhold onbode(nd,dd,w)hold offfigurestep(n1,d1,5)gridhold onpausestep(n2,d2,5)hold off可以得到:p2m =80.0044 wgc2 =4.0004（即：系统相角裕量PM=80°，增益穿越频率=4rad/s） 图6-3 系统Bode图 图6-4 闭环系统的阶跃响应第 10 页