兰州理工大学《控制系统计算机仿真》上机报告院系:电信学院班级:基地2班姓名:学号:时间:2011年05月16日电气工程与信息工程学院《控制系统计算机仿真》上机报告兰州理工大学共页第页2-2用MATLAB语言求下列系统的状态方程、传递函数、零极点增益和部分分式形式的模型参数,并分别写出其相应的数学模型表达式:(1)2450351024247)(23423ssssssssG(2)uXX022475.025.075.125.1125.15.025.025.025.125.425.25.025.1525.2Xy2020解:(1):1)num=[172424];den=[110355024];[ABCD]=tf2ss(num,den)A=-10-35-50-24100001000010B=1000C=172424D=0即状态方程为:uXX000101000010000124503510Xy4224712)num=[172424];den=[110355024];[ZPK]=tf2zp(num,den)Z=-2.7306+2.8531i-2.7306-2.8531i-1.5388P=-4.0000《控制系统计算机仿真》上机报告兰州理工大学共页第页-3.0000-2.0000-1.0000K=1即表达式)1)(2)(3)(4()5388.1)(8531.27306.2)(8531.27306.2()(sssssisissG3)[RPH]=residue(num,den)R=4.0000-6.00002.00001.0000P=-4.0000-3.0000-2.0000-1.0000H=[]即表达式G(s)=4/(s+4)-6/(s+3)+2/(s+2)+1/(s+1)(2):1)A=[2.25,-5,-1.25,-0.5;2.25,-4.25,-1.25,-0.250.25,-0.5,-1.25,-1;1.25,-1.75,-0.25,-0.75];B=[4;2;2;0];C=[0,2,0,2];D=0;[num,den]=ss2tf(A,B,C,D)num=04.000014.000022.000015.0000den=1.00004.00006.25005.25002.2500即表达式25.225.525.641522144)(23423ssssssssG2)[Z,P,K]=ss2zp(A,B,C,D)Z=-1.0000+1.2247i-1.0000-1.2247i-1.5000P=-0.5000+0.8660i-0.5000-0.8660i-1.5000《控制系统计算机仿真》上机报告兰州理工大学共页第页-1.5000K=4.0000即表达式)5.1)(866.05.0)(866.05.0()2247.11)(2247.11(4)(sisisisissG3)[R,P,H]=residue(num,den)R=4.0000-0.00000.0000-2.3094i0.0000+2.3094iP=-1.5000-1.5000-0.5000+0.8660i-0.5000-0.8660iH=[]即表达式G(s)=4/(s+1.5)-2.3094i/(s+0.5000-0.8660i)+2.3094i/(s+0.5000+0.8660i)2-3用殴拉法求下列系统的输出响应)(ty在10t上,1.0h时的数值解。yy2',1)0(y要求保留4位小数,并将结果以图形的方式与真解tety2)(比较。t=0:0.1:1t=Columns1through800.10000.20000.30000.40000.50000.60000.7000Columns9through110.80000.90001.0000h=0.1;y(1)=1;fori=1:10y(i+1)=y(i)+h*(-2*y(i));endplot(t,y,'g')holdonm=exp(-2*t)m=Columns1through81.00000.81870.67030.54880.44930.3679《控制系统计算机仿真》上机报告兰州理工大学共页第页0.30120.2466Columns9through110.20190.16530.1353plot(t,m,'r')2-5用四阶龙格-库塔梯形法求解2-3的数值解,并通过与真值及殴拉法的比较,分析其精度。t=0:0.1:1;h=0.1;y(1)=1;fori=1:10k1=-2*y(i);k2=-2*(y(i)+h*k1/2);k3=-2*(y(i)+h*k2/2);k4=-2*(y(i)+h*k3);y(i+1)=y(i)+h*(k1+2*k2+2*k3+k4)/6;endplot(t,y,'g-')holdonm=exp(-2*t)m=Columns1through81.00000.81870.67030.54880.44930.36790.30120.2466Columns9through11《控制系统计算机仿真》上机报告兰州理工大学共页第页0.20190.16530.1353plot(t,m,'r--')00.10.20.30.40.50.60.70.80.910.10.20.30.40.50.60.70.80.914-2设典型闭环结构控制系统如下图所示,当阶跃输入幅值20R时,用sp4_1.m求取输出y(t)的响应。a=[154204.2213.863.5];b=[18750001562000];X0=[0000];V=0.002;n=4;T0=0;Tf=10;h=0.01;R=20;sp4_15.638.2132.2045410562.110875.123466sssssy(t)r(t)_0.002《控制系统计算机仿真》上机报告兰州理工大学共页第页附文件sp4_1.mb=b/a(1);a=a/a(1);A=a(2:n+1);A=[rot90(rot90(eye(n-1,n)));-fliplr(A)];B=[zeros(1,n-1),1]';m1=length(b);C=[fliplr(b),zeros(1,n-m1)];Ab=A-B*C*V;X=X0';y=0;t=T0;N=round(Tf-T0)/h;fori=1:NK1=Ab*X+B*RK2=Ab*(X+h*K1/2)+B*RK3=Ab*(X+h*K2/2)+B*RK4=Ab*(X+h*K3)+B*RX=X+h*(K1+2*K2+2*K3+K4)/6;y=[y,C*X];t=[t,t(i)+h];end[t',y']plot(t,y)4-5下图中,若各环节传递函数已知为:ssG01.011)(1sssG085.017.01)(2,,sssG051.015.01)(3,ssG15.0121.0)(4,ssG01.011.0)(5,ssG01.010044.0)(6,但《控制系统计算机仿真》上机报告兰州理工大学共页第页212.0)(7sG;试列写链接矩阵W、W0和非零元素阵WIJ,将程序sp4_2完善后,应用此程序求输出4y的响应曲线。P=[10.0110;00.08510.17;00.05110.15;10.150.210;10.010.10;10.010.00440];WIJ=[101;211;26-1;321;35-1;431;44-0.212;541;641];n=6;Y0=1;Yt0=[000000];h=0.01;L1=10;T0=0;Tf=50;nout=6;sp4_2;010203040506000.20.40.60.811.21.4附文件sp4_2.mA=diag(P(:,1));B=diag(P(:,2));G2(s)G7(s)G3(s)y0_y4_G1(s)G4(s)G5(s)G6(s)《控制系统计算机仿真》上机报告兰州理工大学共页第页C=diag(P(:,3));D=diag(P(:,4));m=length(WIJ(:,1));W0=zeros(n,1);W=zeros(n,n);fork=1:mif(WIJ(k,2)==0);W0(WIJ(k,1))=WIJ(k,3);elseW(WIJ(k,1),WIJ(k,2))=WIJ(k,3);end;end;Q=B-D*W;Qn=inv(Q);R=C*W-A;V1=C*W0;Ab=Qn*R;b1=Qn*V1;Y=Yt0';y=Y(nout);t=T0;N=round(Tf-T0)/(h*L1);fori=1:Nforj=1:L1;K1=Ab*Y+b1*Y0K2=Ab*(Y+h*K1/2)+b1*Y0K3=Ab*(Y+h*K2/2)+b1*Y0K4=Ab*(Y+h*K3)+b1*Y0Y=Y+h*(K1+2*K2+2*K3+K4)/6;endy=[y,Y(nout)];t=[t,t(i)+h*L1];end[t',y']plot(t,y)4-7用离散相似法仿真程序sp4_3.m重求上题输出4y的数据与曲线,并与四阶龙格-库塔法比较精度。P=[10.0110;00.08510.17;00.05110.15;10.150.210;10.010.10;10.010.00440];WIJ=[101;211;26-1;321;35-1;431;44-0.212;541;641];n=6;Y0=1;Yt0=[000000];h=0.01;L1=10;T0=0;Tf=50;nout=6;sp4_3;plot(t,y)《控制系统计算机仿真》上机报告兰州理工大学共页第页0510152025303540455000.20.40.60.811.21.4附文件sp4_3.mA=diag(P(:,1));B=diag(P(:,2));C=diag(P(:,3));D=diag(P(:,4));m=length(WIJ(:,1));W0=zeros(n,1);W=zeros(n,n);fork=1:mif(WIJ(k,2)==0);W0(WIJ(k,1))=WIJ(k,3);elseW(WIJ(k,1),WIJ(k,2))=WIJ(k,3);end;end;fori=1:nif(A(i,i)==0);FI(i,i)=1;FIM(i,i)=h*C(i,i)/B(i,i);FIJ(i,i)=h*h*C(i,i)/B(i,i)/2;FIC(i,i)=1;FID(i,i)=0;if(D(i,i)~=0);FID(i,i)=D(i,i)/B(i,i);elseendelseFI(i,i)=exp(-h*A(i,i)/B(i,i));FIM(i,i)=(1-FI(i,i))*C(i,i)/A(i,i);FIJ(i,i)=h*C(i,i)/A(i,i)-FIM(i,i)*B(i,i)/A(i,i);《控制系统计算机仿真》上机报告兰州理工大学共页第页FIC(i,i)=1;FID(i,i)=0;if(D(i,i)~=0);FIC(i,i)=C(i,i)/D(i,i)-A(i,i)/B(i,i);FID(i,i)=D(i,i)/B(i,i);elseendendendY=zeros(n,1);X=Y;y=0;Uk=zeros(n,1);Ub=Uk;t=T0:h*L1:Tf;N=length(t);fork=1:N-1fori=1:L1Ub=Uk;Uk=W*Y+W0*Y0;Udot=(Uk-Ub)/h;Uf=2*Uk-Ub;X=FI'*X+FIM'*Uk+FIJ'*Udot;Y=FIC'*X+FID'*Uf;en