第5-1章 频域离散相似法仿真Matlab编程

第5章快速数字仿真Matlab编程5.1面向Z传递函数的仿真针对下式Z传递函数描述的离散系统编写程序(5.1-1)以k为当前步，则有将y(k)的系数归一化(5.1-2)1仿真中，应该定义n+1个元素的数组来存储y序列，且对应关系为同样定义m+1个元素存储u序列2function[t,y]=GZSimu(tstart,tstop,h,num,den,iniu,iniy,u)%函数功能：对Z传递函数进行仿真%输入参数：tstart,tstop,h分别是起始时间、结束时间和仿真步长，是标量%num,den是Z传递函数的分子分母多项式系数%iniu,iniy是u,y的初始值,维数与num,den一致，第1个元素存最新的%u是系统输入%输出参数：t是仿真结果的时间序列%y是仿真结果系统的输出序列函数原型定义程序的主要功能是：1】将由num,den描述的Z传递函数(5.1-1)转换为迭代计算式(5.1-2)2】根据iniu,iniy和u进行迭代计算3function[t,y]=GZSimu(tstart,tstop,h,num,den,iniu,iniy,u)t=[tstart:h:tstop];%t数一个行序列cntt=size(t,2);%返回列数a=den(1);%得到分母最高次项的系数by=den(2:(size(den,2)));%得到分母除最高次项外的余项by=-by/a;%将分母余项除以a，并移到差分式的右边bu=num/a;%将分子除以acntu=size(num,2);%u的存储步数cnty=size(den,2);%y的存储步数cury=iniy(2:cnty);%得到用于循环的y序列curu=iniu;%得到用于循环的u序列y=zeros(1,cntt);%构造一个空矩阵，用来存储结果y(:,1)=iniy(1);%作为输出的第1列fori=1:1:cntt-1forj=cntu:-1:2%将存储的u后移curu(j)=curu(j-1);end;curu(1)=u;%curu(1)存最新的uy0=by*(cury')+bu*(curu');%计算最新的yy(:,i+1)=y0;%将输出加入到输出序列里forj=cnty-1:-1:2%将存储的y后移cury(j)=cury(j-1);end;cury(1)=y0;%cury(1)存最新的yend4【例5.1】对于如下由传递函数描述的系统，假设初值为零，输入为单位阶跃作用。采用双线性变换和RK4进行仿真对比。解：将双线性变换带入上式，得其中当步长取0.5时，a=0.16,b=0.6，替换法所得差分方程为：5functionsimu5_1tstart=0;tstop=10;h=0.5;num=[0.16,0,-0.16];den=[1,-1.2,0.36];iniu=[0,0,0];iniy=[0,0,0];u=1;[t,y]=GZSimu(tstart,tstop,h,num,den,iniu,iniy,u);%用差分方程仿真A=[0,1;-1,-2];B=[0,1]';C=[0,1];D=[0];x0=[0,0]';u0=[1];[t,y1]=w_LinearSimu(tstart,tstop,h,x0,u0,A,B,C,D,'RK4');stairs(t,y,'b');holdon;plot(t,y1,'-r*');holdoff;xlabel('x');ylabel('y');仿真程序：对此系统采用差分仿真和RK4两种方法仿真6对连续模型采用RK4仿真；模型离散化采用双线性替换法；两种仿真方法存在一定差别，但动态特性相似，静态特性一致，满足匹配原理。75.2程序化替换法编程程序化替换法的详细原理见教材。编写一个程序，实现程序化双线性替换，由S传递函数求Z传递函数。Function[numz,denz]=w_Gs2Gz(nums,dens,T)%函数功能：采用程序化双线性变换，求S传递函数对应的Z传递函数%输入参数：nums,dens是S传递函数的分子和分母多项式系数%T是离散化步长%输出参数：numz,denz是求得的Z传递函数的分子和分母多项式系数8function[numz,denz]=w_Gs2Gz(nums,dens,T)m=length(nums);%分子阶数n=length(dens);%分母阶数ifm>nerror('s传递函数的分子阶数不能高于分子阶数！');end;ifm<nnums=[zeros(1,n-m),nums];end;XX=ones(n,n);baseitem=[1];fori=1:1:n-1%计算（z+1）^n的系数baseitem=conv(baseitem,[1,1]);end;XX(1,:)=baseitem;fori=2:nforj=2:nXX(i,j)=XX(i-1,j)-XX(i,j-1)-XX(i-1,j-1);end;end;diagT=diag((2/T).^(0:n-1));H=diagT*XX;numz=fliplr(nums)*H;denz=fliplr(dens)*H;双线性变换函数实现9用于模型离散化的Matlab函数S传递函数到Z传递函数转换的Matlab函数为：[numd,dend]=c2dm(num,den,T,’选项’)其中，num,den为S传递函数的分子分母多项式系数；T为离散化步长；Numd,dend为Z传递函数的分子分母多项式系数10【例5.2】采用c2dm和w_Gs2Gz求如下系统的离散模型。>>num=[1,0];>>den=[1,2,1];>>[numz,numz]=c2dm(num,den,0.5,'tustin')numz=0.16000-0.1600numz=1.0000-1.20000.3600>>[numz,numz]=c2dm(num,den,0.5,'matched')numz=00.3284-0.3284numz=1.0000-1.21310.3679>>[numz,denz]=w_Gs2Gz(num,den,0.5)numz=0.16000-0.1600denz=1.0000-1.20000.360011functionsimu5_2tstart=0;tstop=10;h=0.5;num=[1,0];den=[1,2,1];iniu=[0,0,0];iniy=[0,0,0];u=1;[numd,dend]=c2dm(num,den,h,'tustin');[t,y]=GZSimu(tstart,tstop,h,numd,dend,iniu,iniy,u);%用差分方程仿真[numd,dend]=c2dm(num,den,h,'matched');[t,y2]=GZSimu(tstart,tstop,h,numd,dend,iniu,iniy,u);%匹配法模型plot(t,y,'-b*',t,y2,':k+');xlabel('x');ylabel('y');将由c2dm得到的两种模型进行仿真对比125.3可调整积分法仿真【例5.3】对下图所示的系统采用离散型补偿，并仿真。下图中，图a是原系统，图b是等价的离散系统。等效离散系统的图中，引入采样保持器、连续补偿器和延迟一拍环节。反馈回路中引入延迟一拍环节，是因为在计算反馈回路之前，必须先计算正向回路，因此将反馈回路延迟一拍。13原系统闭环传递函数为采用单零点补偿经计算得到系统的补偿离散模型为取K=2,步长0.5时，14替换法差分方程取K=2,步长0.5时，15functionsimu5_4tstart=0;tstop=5;h=0.5;num=[1,1];den=[6,-2];iniu=[0,0];iniy=[0,0];u=1;[t,y]=GZSimu(tstart,tstop,h,num,den,iniu,iniy,u);%用替换法模型num=[0.3161,0];den=[1,-0.3679];[t,y2]=GZSimu(tstart,tstop,h,num,den,iniu,iniy,u);%用补偿模型A=[-2];B=[1];C=[1];D=[0];x0=[0];u0=[1];[t,y1]=