第4-1章-时域离散相似法Matlab编程课件

第4章时域离散相似法仿真Matlab编程4.1连续模型的离散化线性连续状态空间模型为(4.1-1)其采用零阶保持器的离散相似模型为(4.1-2)其中(4.1-3)1第4章时域离散相似法仿真Matlab编程4.1连续模型的用于离散相似的Matlab函数连续模型转换为离散模型：[G,H]=c2d(A,B,T)T为采样周期可选参数离散相似：[G,H,C,D]=c2dm(A,B,C,D,T,’选项’)[numd,dend]=c2dm(num,den,T,’选项’)2用于离散相似的Matlab函数连续模型转换为离散模型：[G,【例4.1.1】对下面的连续系统，在采样周期T=0.1时进行离散化>>A=[0,1;-2,-3];>>B=[0;1];>>[G,H]=c2d(A,B,0.1);G=0.99090.0861-0.17220.7326H=0.00450.0861经计算，可得当T=0.1时，计算结果与c2d函数得到的结果一致。3【例4.1.1】对下面的连续系统，在采样周期T=0.1时进行4.2针对离散状态空间模型的仿真程序离散状态空间模型为对此模型编写仿真程序，只需要按模型方程迭代即可。function[t,y]=w_DiscreteSimu(tstart,tstop,h,x0,u0,G,H,C,D)%函数功能：对线性离散系统x(k+1)=G*X(k)+H*u(k),y(k)=C*x(k)+D*u(k)进行仿真%输入参数：tstart,tstop,h分别是起始时间、结束时间和仿真步%长，是标量%x0,u0是状态、输入的初值，都是列向量%G,H,C,D是线性离散状态空间模型的系数矩阵%输出参数：t是仿真结果的时间序列%y是仿真结果系统的输出序列44.2针对离散状态空间模型的仿真程序离散状态空间模型为对此function[t,y]=w_DiscreteSimu(tstart,tstop,h,x0,u0,G,H,C,D)t=[tstart:h:tstop];%t数一个行序列cntt=size(t,2);%返回列数cnty=size(C,1);%返回y的维数y=zeros(cnty,cntt);%构造一个空矩阵，用来存储结果y0=C*x0+D*u0;y(:,1)=y0;%将y0作为输出的第1列curx=x0;%当前一步的xcuru=u0;%当前一步的ucury=y0;%当前一步的yfori=1:1:cntt-1curx=G*curx+H*curu;%计算下一步的状态cury=C*curx+D*curu;%计算输出y(:,i+1)=cury;%将输出加入到输出序列里end离散状态空间模型仿真程序5function[t,y]=w_DiscreteSimu(【例4.2.1】对下面的连续系统，采用数值积分法和离散相似法仿真。解：该系统连续状态空间模型矩阵为其离散状态空间模型矩阵为也可以采用c2d函数直接计算离散模型矩阵6【例4.2.1】对下面的连续系统，采用数值积分法和离散相似法仿真程序:用c2d进行模型变换，将数值积分与离散相似2种方法的结果对比。数值积分仿真的结果用plot函数绘图，表现连续的效果；离散模型仿真结果用stairs绘图，各数据点之间是零阶保持。functionsimu4_2tstart=0; tstop=6; A=[0,1;-2,-3]; B=[0;1]; C=[2,0]; D=[0];x0=[0;0]; u0=[1];h=0.5; [G,H]=c2d(A,B,h);%得到离散模型矩阵%数值积分仿真[t,y1]=w_LinearSimu(tstart,tstop,h,x0,u0,A,B,C,D,'RK4');%离散相似法仿真[t,y2]=w_DiscreteSimu(tstart,tstop,h,x0,u0,G,H,C,D);plot(t,y1,'-r*');holdon;stairs(t,y2,':b+');holdoff;axis([0,6,0,1.2]);xlabel('x');ylabel('y');legend('数值积分','离散相似');title('数值积分与离散相似仿真对比')7仿真程序:用c2d进行模型变换，将数值积分与离散相似2种方法在步长h=0.5时，两种方法的仿真结果及其接近。离散模型的每个计算点都在连续模型上当步长h=1时，两种方法的结果出现偏差。离散模型的点不一定在连续模型上。8在步长h=0.5时，两种方法的仿真结果及其接近。离散模型的每4.3面向环节离散化仿真(无非线性环节)【例4.3.1】将下面的系统按环节离散化仿真解：将系统分为2个典型环节，如图所示，开环部分传函为1】环节之间的连接方程94.3面向环节离散化仿真(无非线性环节)【例4.3.1】将2】环节内部离散模型（1）积分环节（2）惯性环节环节内部动态方程环节输出方程102】环节内部离散模型（1）积分环节（2）惯性环节环节内部动态3】系统输出方程面向环节离散化仿真时应该注意的问题应该依次计算各环节的输出，即计算出第1个环节，再利用第1个环节的输出计算第2个环节，依次下去。而不要一次计算出各环节的输入，再算出各环节的输出，这样相当于在每个环节之间存在一个步长的延迟环节，模型的准确性会大大降低。下面，我们将对这两种方式做仿真对比，进而说明此问题。113】系统输出方程面向环节离散化仿真时应该注意的问题11function[t,p]=w_DisNodesSimu2(tstart,tstop,h,x0,y0,r0,W,W0,G,H,C,D,Q)%函数功能：面向环节离散化仿真，%环节输入方程：u(k)=W*x(k)+W0*r(k)%环节迭代方程：x(k+1)=G*x(k)+H*u(k)%环节输出方程y(k)=C*x(k)+D*u(k);%整个系统的输出方程p(k)=Q*y(k);%输入参数：tstart,tstop,h分别是起始时间、结束时间和仿真步长%x0,y0是每个环节状态和输出的初值，%r0为系统外部输入%输出参数：t是仿真结果的时间序列%p是仿真结果系统的输出序列1、整体刷新计算的按环节仿真函数整体刷新即根据环节连接方程计算每个环节的输入，再用计算出来的输入向量根据各环节的离散模型计算各环节的输出。各环节之间相当于并行运算。12function[t,p]=w_DisNodesSimu2function[t,p]=w_DisNodesSimu2(tstart,tstop,h,x0,y0,r0,W,W0,G,H,C,D,Q)t=[tstart:h:tstop];%t数一个行序列cntt=size(t,2);%返回列数cnty=size(Q,1);%返回y的维数p=zeros(cnty,cntt);%构造一个空矩阵，用来存储结果p(:,1)=Q*y0;%fori=1:1:cntt-1u=W*y0+W0*r0;x0=G*x0+H*u;y0=C*x0+D*u;p(:,i+1)=Q*y0;%将y0作为输出的第1列end13function[t,p]=w_DisNodesSimu22、逐个环节刷新的仿真程序该函数所有参数与w_DisNodesSimu2的一样。逐个环节刷新是逐个根据环节的编号顺序，计算环节的输入，再计算输出，再计算下一个环节的输入和输出，依次进行，相当于环节之间串连运算。function[t,p]=w_DisNodesSimu(tstart,tstop,h,x0,y0,r0,W,W0,G,H,C,D,Q)t=[tstart:h:tstop];%t数一个行序列cntnodes=size(W,1);%得到环节的个数cntt=size(t,2);%返回列数cnty=size(Q,1);%返回y的维数p=zeros(cnty,cntt);%构造一个空矩阵，用来存储结果p(:,1)=Q*y0;%fori=1:1:cntt-1forj=1:1:cntnodesu(j)=W(j,:)*y0+W0(j,:)*r0;x0(j)=G(j,j)*x0(j)+H(j,j)*u(j);y0(j)=C(j,j)*x0(j)+D(j,j)*u(j);end;p(:,i+1)=Q*y0;%将y0作为输出的第1列end142、逐个环节刷新的仿真程序function[t,p]=w_functionsimu4_3tstart=0;tstop=10;A=[0,1;-2,-3];B=[0;1];C=[2,0];D=[0];x0=[0;0];u0=[1];h=0.5;[G,H]=c2d(A,B,h);%得到离散模型矩阵[t,y1]=w_LinearSimu(tstart,tstop,h,x0,u0,A,B,C,D,'RK4');[t,y2]=w_DiscreteSimu(tstart,tstop,h,x0,u0,G,H,C,D);%离散相似法仿真%面向环节离散化仿真W=[0,-1;1,0];W0=[1,0]';Q=[0,1];C=diag([1,1]);D=zeros(2);G=diag([1,exp(-3*h)]);H=diag([2*h,(1-exp(-3*h))/3]);y0=x0;[t,y3]=w_DisNodesSimu2(tstart,tstop,h,x0,y0,u0,W,W0,G,H,C,D,Q);plot(t,y1,‘-r’);holdon; stairs(t,y2,':k+'); stairs(t,y3,'-b*'); holdoff;axis([0,10,0,1.2]); xlabel('x');ylabel('y');legend('数值积分','整体离散','按环节离散');title('各环节整体刷新（h=0.5）');仿真调用函数数值积分的结果用plot绘图，离散模型仿真结果用stairs绘图。15functionsimu4_3仿真调用函数数值积分的结果用整体刷新时相当于前后环节之间存在一个步长的延迟，导致模型不精确。环节逐个刷新符合系统实际运算意义，模型比较精确16整体刷新时相当于前后环节之间存在一个步长的延迟，导致模型不精4.4面向环节离散化仿真（有非线性环节）环节划分原则1】线性环节都处理为1阶；2】非线性环节附加在线性环节的前面或后面，环节个数等于线性环节的个数。一个环节的计算顺序1】根据环节之间的连接方程计算该环节的输入；2】根据前置非线性环节计算非线性环节的输出；3】前置非线性环节的输出作为线性环节的输入，计算线性环节的输出；4】线性环节的输出作为后置非线性环节的输入，计算整个环节的输出。174.4面向环节离散化仿真（有非线性环节）环节划分原则一个环模型描述1】环节之间连接方程，写出各环节的输入与其他环节的输出和系统输入之间的连续方程。2】线性环节离散模型3】其中18模型描述1】环节之间连接方程，写出各环节的输入与其他环节的输4】非线性环节模型描述根据非线性环节类型不同，分别用不同的参数描述，包括非线性环节类型和必要的参数。非线性环节作为线性环节的前置或后置环节。5】整个系统的输出方程仿真计算顺序根据环节编号顺序，计算一个环节，再计算下一个环节。不能整体计算每个环节的输入，再计算每个环节的输出，这样每个环节之间将产生滞后。194】非线性环节模型描述5】整个系统的输出方程仿真计算顺序1【例4.4.1】面向环节的非线性系统仿真右图所示系统有一个饱和非线性环节和一个线性环节。饱和非线性不能用微分方程表示，因而不能对整个系统采用数值积分法仿真。采用按环节离散化仿真，首先将线性环节分解为两个基本环节，将非线性环节作为第1个线性环节的前置环节。仿真中，取c=20，K＝1，r=2520【例4.4.1】面向环节的非线性系统仿真右图所示系统有一个1】各线性环节内部离散方程（1）积分环节（3）惯性环节用矩阵表示为211】各线性环节内部离散方程（1）积分环节（3）惯性环节用矩阵2】各环节之间连接方程3】系统输出为4】饱和非线性环节作为第1个环节的前置环节222】各环节之间连接方程3】系统输出为4】饱和非线性环节作为第function[t,p]=w_DisNodesSimuNonL(tstart,tstop,h,x0,y0,r0,W,W0,G,H,C,D,Q,NLBefore,NLBeforeParams,NLAfter,NLAfterParams)%函数功能：带有非线性环节的面向环节离散化仿真，%环节输入方程：u(k)=W*x(k)+W0*r(k)%环节迭代方程：x(k+1)=G*x(k)+H*u(k)%环节输出方程y(k)=C*x(k)+D*u(k);%整个系统的输出方程p(k)=Q*y(k);%输入参数：tstart,tstop,h分别是起始时间、结束时间和仿真步长，是标量%x0,y0是每个环节状态和输出的初值，%r0为系统外部输入%NLBefore是前置非线性环节的名称，如果没有就写'None',NLBeforeParams是前置非线性环节的参数，每个环节%的参数占一行。%NLAfter是后置非线性环节的名称，如果没有就写'None',NLAfterParams是后置非线性环节的参数，每个环节%的参数占一行。%输出参数：t是仿真结果的时间序列%p是仿真结果系统的输出序列1、带非线性环节的按环节离散化仿真程序原型23function[t,p]=w_DisNodesSimuNfunction[t,p]=w_DisNodesSimuNonL(tstart,tstop,h,x0,y0,r0,W,W0,G,H,C,D,Q,NLBefore,NLBeforeParams,NLAfter,NLAfterParams)t=[tstart:h:tstop];%t数一个行序列cntnodes=size(W,1);%得到环节的个数cntt=size(t,2);%返回列数cnty=size(Q,1);%返回y的维数p=zeros(cnty,cntt);%构造一个空矩阵，用来存储结果p(:,1)=Q*y0;fori=1:1:cntt-1forj=1:1:cntnodesu(j)=W(j,:)*y0+W0(j,:)*r0;%根据环节之间的连接计算，u(j)=NonlinearNode(NLBefore(j),u(j),NLBeforeParams(j,:));%计算第j个环节的前置非线性环节的输出x0(j)=G(j,j)*x0(j)+H(j,j)*u(j);y0(j)=C(j,j)*x0(j)+D(j,j)*u(j);%根据线性部分计算线性部分的输出y0(j)=NonlinearNode(NLAfter(j),y0(j),NLAfterParams(j,:));%计算第j个环节的后置非线性环节的输出end;p(:,i+1)=Q*y0;%将y0作为输出的第1列end24function[t,p]=w_DisNodesSimuN2、典型非线性环节函数functionu_out=NonlinearNode(NodeType,u_in,params)%函数功能：典型非线性环节的计算%输入参数：NodeType节点类型，可以取'Saturation','Deadzone','Relay','RelayDeadzone',%'None'表示没有非线性环节%u_in,环节的输入；%params,环节的参数行向量%输出结果：u_out，环节的输出252、典型非线性环节函数functionu_out=NonlswitchNodeTypecase'Saturation'c=params(1);k=params(2);if(k*u_in)<-cu_out=-c;elseif(k*u_in)>cu_out=c;elseu_out=k*u_in;end;case'Deadzone'c=params(1);k=params(2);if(u_in<-c)|(u_in>c)u_out=k*u_in;elseu_out=0;end;case'Relay'c=params(1);ifu_in<0u_out=-c;elseifu_in>0u_out=c;elseu_out=0;end;case'RelayDeadzone'c=params(1);h=params(2);ifu_in<-hu_out=-c;elseifu_in>cu_out=c;elseu_out=0;end;otherwiseu_out=u_in;end;函数NonlinearNode的主体26switchNodeTypecase'Relay'函functionsimu_NonNodeststart=0;tstop=20;h=0.5;K=1;c=10;W=[0,-1;1,0];W0=[1,0]';Q=[0,1];C=diag([1,1]);D=zeros(2);G=diag([1,exp(-h)]);H=diag([K*h,1-exp(-h)]);x0=[0;0];u0=[25];y0=x0;NLBefore=['Saturation','None'];NLBeforeParams=[c,1;0,0];NLAfter=['None','None'];NLAfterParams=[0,0;0,0];[t,y]=w_DisNodesSimuNonL(tstart,tstop,h,x0,y0,u0,W,W0,G,H,C,D,Q,NLBefore,NLBeforeParams,NLAfter,NLAfterParams);stairs(t,y,'-r*');xlabel('Time');ylabel('y');title('含非线性环节的面向环节仿真');3、仿真调用函数27functionsimu_NonNodes3、仿真调用函仿真结果：同样的结果，用plot与stairs绘图的对比。仿真中，取c=20，K＝1，r=2528仿真结果：同样的结果，用plot与stairs绘图的对比。仿4.5利用Simulink进行面向环节仿真【例4.5.1】下图所示线性与非线性混合系统，当回环非线性特性c=1时，系统有自持振荡，振幅约为10.8~9.4，振荡周期约为6s.对该系统采用Simulink建模，进行仿真294.5利用Simulink进行面向环节仿真【例4.5.1】Simulink将此系统看做有6个连续环节。将仿真输出用Scope绘图，同时将结果数据输出到Workspace里，以便用plot绘图。从曲线，可以很明显地看到输出震荡。30Simulink将此系统看做有6个连续环节。将仿真输出用Sc【例4.5.2】用Simulik对下图所示的系统进行仿真仿真中，取c=10，K＝1，r=25，并且在t=1时才开始阶跃输入。31【例4.5.2】用Simulik对下图所示的系统进行仿真仿真仿真结果曲线：很明显，饱和环节的输出受限制。32仿真结果曲线：很明显，饱和环节的输出受限制。32【例4.5.3】对【例4.3.1】所示的系统采用Simulink仿真（1）积分环节（2）惯性环节当T=0.5时33【例4.5.3】对【例4.3.1】所示的系统采用Simul本例采用离散化环节进行建模34本例采用离散化环节进行建模34第4章时域离散相似法仿真Matlab编程4.1连续模型的离散化线性连续状态空间模型为(4.1-1)其采用零阶保持器的离散相似模型为(4.1-2)其中(4.1-3)35第4章时域离散相似法仿真Matlab编程4.1连续模型的用于离散相似的Matlab函数连续模型转换为离散模型：[G,H]=c2d(A,B,T)T为采样周期可选参数离散相似：[G,H,C,D]=c2dm(A,B,C,D,T,’选项’)[numd,dend]=c2dm(num,den,T,’选项’)36用于离散相似的Matlab函数连续模型转换为离散模型：[G,【例4.1.1】对下面的连续系统，在采样周期T=0.1时进行离散化>>A=[0,1;-2,-3];>>B=[0;1];>>[G,H]=c2d(A,B,0.1);G=0.99090.0861-0.17220.7326H=0.00450.0861经计算，可得当T=0.1时，计算结果与c2d函数得到的结果一致。37【例4.1.1】对下面的连续系统，在采样周期T=0.1时进行4.2针对离散状态空间模型的仿真程序离散状态空间模型为对此模型编写仿真程序，只需要按模型方程迭代即可。function[t,y]=w_DiscreteSimu(tstart,tstop,h,x0,u0,G,H,C,D)%函数功能：对线性离散系统x(k+1)=G*X(k)+H*u(k),y(k)=C*x(k)+D*u(k)进行仿真%输入参数：tstart,tstop,h分别是起始时间、结束时间和仿真步%长，是标量%x0,u0是状态、输入的初值，都是列向量%G,H,C,D是线性离散状态空间模型的系数矩阵%输出参数：t是仿真结果的时间序列%y是仿真结果系统的输出序列384.2针对离散状态空间模型的仿真程序离散状态空间模型为对此function[t,y]=w_DiscreteSimu(tstart,tstop,h,x0,u0,G,H,C,D)t=[tstart:h:tstop];%t数一个行序列cntt=size(t,2);%返回列数cnty=size(C,1);%返回y的维数y=zeros(cnty,cntt);%构造一个空矩阵，用来存储结果y0=C*x0+D*u0;y(:,1)=y0;%将y0作为输出的第1列curx=x0;%当前一步的xcuru=u0;%当前一步的ucury=y0;%当前一步的yfori=1:1:cntt-1curx=G*curx+H*curu;%计算下一步的状态cury=C*curx+D*curu;%计算输出y(:,i+1)=cury;%将输出加入到输出序列里end离散状态空间模型仿真程序39function[t,y]=w_DiscreteSimu(【例4.2.1】对下面的连续系统，采用数值积分法和离散相似法仿真。解：该系统连续状态空间模型矩阵为其离散状态空间模型矩阵为也可以采用c2d函数直接计算离散模型矩阵40【例4.2.1】对下面的连续系统，采用数值积分法和离散相似法仿真程序:用c2d进行模型变换，将数值积分与离散相似2种方法的结果对比。数值积分仿真的结果用plot函数绘图，表现连续的效果；离散模型仿真结果用stairs绘图，各数据点之间是零阶保持。functionsimu4_2tstart=0; tstop=6; A=[0,1;-2,-3]; B=[0;1]; C=[2,0]; D=[0];x0=[0;0]; u0=[1];h=0.5; [G,H]=c2d(A,B,h);%得到离散模型矩阵%数值积分仿真[t,y1]=w_LinearSimu(tstart,tstop,h,x0,u0,A,B,C,D,'RK4');%离散相似法仿真[t,y2]=w_DiscreteSimu(tstart,tstop,h,x0,u0,G,H,C,D);plot(t,y1,'-r*');holdon;stairs(t,y2,':b+');holdoff;axis([0,6,0,1.2]);xlabel('x');ylabel('y');legend('数值积分','离散相似');title('数值积分与离散相似仿真对比')41仿真程序:用c2d进行模型变换，将数值积分与离散相似2种方法在步长h=0.5时，两种方法的仿真结果及其接近。离散模型的每个计算点都在连续模型上当步长h=1时，两种方法的结果出现偏差。离散模型的点不一定在连续模型上。42在步长h=0.5时，两种方法的仿真结果及其接近。离散模型的每4.3面向环节离散化仿真(无非线性环节)【例4.3.1】将下面的系统按环节离散化仿真解：将系统分为2个典型环节，如图所示，开环部分传函为1】环节之间的连接方程434.3面向环节离散化仿真(无非线性环节)【例4.3.1】将2】环节内部离散模型（1）积分环节（2）惯性环节环节内部动态方程环节输出方程442】环节内部离散模型（1）积分环节（2）惯性环节环节内部动态3】系统输出方程面向环节离散化仿真时应该注意的问题应该依次计算各环节的输出，即计算出第1个环节，再利用第1个环节的输出计算第2个环节，依次下去。而不要一次计算出各环节的输入，再算出各环节的输出，这样相当于在每个环节之间存在一个步长的延迟环节，模型的准确性会大大降低。下面，我们将对这两种方式做仿真对比，进而说明此问题。453】系统输出方程面向环节离散化仿真时应该注意的问题11function[t,p]=w_DisNodesSimu2(tstart,tstop,h,x0,y0,r0,W,W0,G,H,C,D,Q)%函数功能：面向环节离散化仿真，%环节输入方程：u(k)=W*x(k)+W0*r(k)%环节迭代方程：x(k+1)=G*x(k)+H*u(k)%环节输出方程y(k)=C*x(k)+D*u(k);%整个系统的输出方程p(k)=Q*y(k);%输入参数：tstart,tstop,h分别是起始时间、结束时间和仿真步长%x0,y0是每个环节状态和输出的初值，%r0为系统外部输入%输出参数：t是仿真结果的时间序列%p是仿真结果系统的输出序列1、整体刷新计算的按环节仿真函数整体刷新即根据环节连接方程计算每个环节的输入，再用计算出来的输入向量根据各环节的离散模型计算各环节的输出。各环节之间相当于并行运算。46function[t,p]=w_DisNodesSimu2function[t,p]=w_DisNodesSimu2(tstart,tstop,h,x0,y0,r0,W,W0,G,H,C,D,Q)t=[tstart:h:tstop];%t数一个行序列cntt=size(t,2);%返回列数cnty=size(Q,1);%返回y的维数p=zeros(cnty,cntt);%构造一个空矩阵，用来存储结果p(:,1)=Q*y0;%fori=1:1:cntt-1u=W*y0+W0*r0;x0=G*x0+H*u;y0=C*x0+D*u;p(:,i+1)=Q*y0;%将y0作为输出的第1列end47function[t,p]=w_DisNodesSimu22、逐个环节刷新的仿真程序该函数所有参数与w_DisNodesSimu2的一样。逐个环节刷新是逐个根据环节的编号顺序，计算环节的输入，再计算输出，再计算下一个环节的输入和输出，依次进行，相当于环节之间串连运算。function[t,p]=w_DisNodesSimu(tstart,tstop,h,x0,y0,r0,W,W0,G,H,C,D,Q)t=[tstart:h:tstop];%t数一个行序列cntnodes=size(W,1);%得到环节的个数cntt=size(t,2);%返回列数cnty=size(Q,1);%返回y的维数p=zeros(cnty,cntt);%构造一个空矩阵，用来存储结果p(:,1)=Q*y0;%fori=1:1:cntt-1forj=1:1:cntnodesu(j)=W(j,:)*y0+W0(j,:)*r0;x0(j)=G(j,j)*x0(j)+H(j,j)*u(j);y0(j)=C(j,j)*x0(j)+D(j,j)*u(j);end;p(:,i+1)=Q*y0;%将y0作为输出的第1列end482、逐个环节刷新的仿真程序function[t,p]=w_functionsimu4_3tstart=0;tstop=10;A=[0,1;-2,-3];B=[0;1];C=[2,0];D=[0];x0=[0;0];u0=[1];h=0.5;[G,H]=c2d(A,B,h);%得到离散模型矩阵[t,y1]=w_LinearSimu(tstart,tstop,h,x0,u0,A,B,C,D,'RK4');[t,y2]=w_DiscreteSimu(tstart,tstop,h,x0,u0,G,H,C,D);%离散相似法仿真%面向环节离散化仿真W=[0,-1;1,0];W0=[1,0]';Q=[0,1];C=diag([1,1]);D=zeros(2);G=diag([1,exp(-3*h)]);H=diag([2*h,(1-exp(-3*h))/3]);y0=x0;[t,y3]=w_DisNodesSimu2(tstart,tstop,h,x0,y0,u0,W,W0,G,H,C,D,Q);plot(t,y1,‘-r’);holdon; stairs(t,y2,':k+'); stairs(t,y3,'-b*'); holdoff;axis([0,10,0,1.2]); xlabel('x');ylabel('y');legend('数值积分','整体离散','按环节离散');title('各环节整体刷新（h=0.5）');仿真调用函数数值积分的结果用plot绘图，离散模型仿真结果用stairs绘图。49functionsimu4_3仿真调用函数数值积分的结果用整体刷新时相当于前后环节之间存在一个步长的延迟，导致模型不精确。环节逐个刷新符合系统实际运算意义，模型比较精确50整体刷新时相当于前后环节之间存在一个步长的延迟，导致模型不精4.4面向环节离散化仿真（有非线性环节）环节划分原则1】线性环节都处理为1阶；2】非线性环节附加在线性环节的前面或后面，环节个数等于线性环节的个数。一个环节的计算顺序1】根据环节之间的连接方程计算该环节的输入；2】根据前置非线性环节计算非线性环节的输出；3】前置非线性环节的输出作为线性环节的输入，计算线性环节的输出；4】线性环节的输出作为后置非线性环节的输入，计算整个环节的输出。514.4面向环节离散化仿真（有非线性环节）环节划分原则一个环模型描述1】环节之间连接方程，写出各环节的输入与其他环节的输出和系统输入之间的连续方程。2】线性环节离散模型3】其中52模型描述1】环节之间连接方程，写出各环节的输入与其他环节的输4】非线性环节模型描述根据非线性环节类型不同，分别用不同的参数描述，包括非线性环节类型和必要的参数。非线性环节作为线性环节的前置或后置环节。5】整个系统的输出方程仿真计算顺序根据环节编号顺序，计算一个环节，再计算下一个环节。不能整体计算每个环节的输入，再计算每个环节的输出，这样每个环节之间将产生滞后。534】非线性环节模型描述5】整个系统的输出方程仿真计算顺序1【例4.4.1】面向环节的非线性系统仿真右图所示系统有一个饱和非线性环节和一个线性环节。饱和非线性不能用微分方程表示，因而不能对整个系统采用数值积分法仿真。采用按环节离散化仿真，首先将线性环节分解为两个基本环节，将非线性环节作为第1个线性环节的前置环节。仿真中，取c=20，K＝1，r=2554【例4.4.1】面向环节的非线性系统仿真右图所示系统有一个1】各线性环节内部离散方程（1）积分环节（3）惯性环节用矩阵表示为551】各线性环节内部离散方程（1）积分环节（3）惯性环节用矩阵2】各环节之间连接方程3】系统输出为4】饱和非线性环节作为第1个环节的前置环节562】各环节之间连接方程3】系统输出为4】饱和非线性环节作为第function[t,p]=w_DisNodesSimuNonL(tstart,tstop,h,x0,y0,r0,W,W0,G,H,C,D,Q,NLBefore,NLBeforeParams,NLAfter,NLAfterParams)%函数功能：带有非线性环节的面向环节离散化仿真，%环节输入方程：u(k)=W*x(k)+W0*r(k)%环节迭代方程：x(k+1)=G*x(k)+H*u(k)%环节输出方程y(k)=C*x(k)+D*u(k);%整个系统的输出方程p(k)=Q*y(k);%输入参数：tstart,tstop,h分别是起始时间、结束时间和仿真步长，是标量%x0,y0是每个环节状态和输出的初值，%r0为系统外部输入%NLBefore是前置非线性环节的名称，如果没有就写'None',NLBeforeParams是前置非线性环节的参数，每个环节%的参数占一行。%NLAfter是后置非线性环节的名称，如果没有就写'None',NLAfterParams是后置非线性环节的参数，每个环节%的参数占一行。%输出参数：t是仿真结果的时间序列%p是仿真结果系统的输出序列1、带非线性环节的按环节离散化仿真程序原型57function[t,p]=w_DisNodesSimuNfunction[t,p]=w_DisNodesSimuNonL(tstart,tstop,h,x0,y0,r0,W,W0,G,H,C,D,Q,NLBefore,NLBeforeParams,NLAfter,NLAfterParams)t=[tstart:h:tstop];%t数一个行序列cntnodes=size(W,1);%得到环节的个数cntt=size(t,2);%返回列数cnty=size(Q,1);%返回y的维数p=zeros(cnty,cntt);%构造一个空矩阵，用来存储结果p(:,1)=Q*y0;fori=1:1:cntt-1forj=1:1:cntnodesu(j)=W(j,:)*y0+W0(j,:)*r0;%根据环节之间的连接计算，u(j)=NonlinearNode(NLBefore(j),u(j),NLBeforeParams(j,:));%计算第j个环节的前置非线性环节的输出x0(j)=G(j,j)*x0(j)+H(j,j)*u(j);y0(j)=C(j,j)*x0(j)+D(j,j)*u(j);%根据线性部分计算线性部分的输出y0(j)=NonlinearNode(NLAfter(j),y0(j),NLAfterParams(j,:));%计算第j个环节的后置非线性环节的输出end;p(:,i+1)=Q*y0;%将y0作为输出的第1列end58function[t,p]=w_DisNodesSimuN2、典型非线性环节函数functionu_out=NonlinearNode(NodeType,u_in,params)%函数功能：典型非线性环节的计算%输入参数：NodeType节点类型，可以取'Saturation','Deadzone','Relay','RelayDeadzone',%'None'表示没有非线性环节%u_in,环节的输入；%params,环节的参数行向量%输出结果：u_out，环节的输出592、典型非线性环节函数functionu_out=NonlswitchNodeTypecase'Saturation'c=params(1);k=params(2);if(k*u_in)<-cu_out=-c;elseif(k*u_in)>c