4.机电系统动态性能的计算机仿真 4.1 概述 机电系统计算机仿真是目前对复杂机电系统进行分析的重要手段与方法。在进行机电系统分析综合与设计工作过程中,除了需要进行理论分析外,还要对系统的特性进行实验研究。系统性能指标与参数是否达到预期的要求?它的经济性能如何?这些都需要在系统设计中给出明确的结论。对于那些在实际调试过程中存在很大风险或实验费用昂贵的系统,一般不允许对设计好的系统直接进行实验,然而没有经过实验研究是不能将设计好的系统直接放到生产实际中去的,因此就必须对其进行模拟实验研究。当然在有些情况下可以构造一套物理模拟装置来进行实验,但这种方法十分费时而且费用又高,而在有的情况下物理模拟几乎是不可能的。近年来随着计算机的迅速发展,采用计算机对机电系统进行数学仿真的方法已被人们采纳。所谓机电系统计算机仿真就是以机电系统的数学模型为基础,借助计算机对机电系统的动静态过程进行实验研究。这里讲的机电系统计算机仿真是指借助数字计算机实现对机电系统的仿真分析。这种实验研究的特点是:将实际系统的运动规律用数学表达式加以描述,它通常是一组常微分方程或差分方程,然后利用计算机来求解这一数学模型,以达到对系统进行分析研究的目的。 对机电系统进行计算机仿真的基本过程包括:首先建立系统的数学模型,因为数学模型是系统仿真的基本依据,所以数学模型极为重要。然后根据系统的数学模型建立相应的仿真模型,一般需要通过一定的算法或数值积分方法对原系统的数学模型进行离散化处理,从而建立起相应的仿真模型,这是进行机电系统仿真分析的关键步骤;最后根据系统的仿真模型编制相应的仿真程序,在计算机上进行仿真实验研究并对仿真结果加以分析。 机电系统计算机仿真的应用与发展已经过了近40年的历程,进入20世纪80 年代以来,随着微型计算机技术以及软件技术的飞速发展与广泛应用,使得机电系统计算机仿真获得了实质性的发展,并使其走进广大的机电系统生产、设计、研究的第一线。目前在该领域应用最为广泛的软件包就是美国Mathworks公司开发的MATLAB语言软件。虽然它最初并非是为机电系统仿真与设计开发的,但是它的强大的矩阵运算功能与图形处理及绘制能力,以及在MATLAB平台上开发出来的各种工具箱,和面向结构图的SIMULINK 系统分析环境,为机电系统计算机仿真提供了强有力的软件工具,从而为机电系统计算机仿真与辅助设计开辟新的方法与途径提供了充分的条件。 考虑到部分读者不具备MATLAB基础知识,因此,在本书的附录一中对MATLAB语言基础知识作了简单的介绍。 4.2面向微分方程的连续系统仿真与MATLAB实现 对于控制系统的最基本数学描述就是微分方程,因此有必要在介绍机电系统计算机仿真时对这种方法加以介绍。实际上在前面章节已经介绍了几种典型数值积分方法以及它们的MATLAB实现。数值积分的数学基础就是微分方程(状态方程),如果一个系统以微分方程的形式加以描述,就可以利用前面的知识对该系统进行仿真研究,有关这部分内容本章就不再加以讨论了。这里将介绍如何利用MATLAB中为用户所提供的常微分方程解函数ode和对微分方程进行求解计算,以及在面向微分方程的系统仿真研究中如何应用M函数。 4.2.1 基于ode函数的面向微分方程的系统仿真 1. 常微分方程解函数ode的基本格式 (l)ode函数的基本分类。在 MATLAB 5.x中 ode函数可以分为两类五种形式。一类是应用于非刚性微分方程求解,如ode45;另一类属于求解刚性微分方程的ode解函数,如ode23s.他们的具体形式如下: 1 ode45函数用于求解非刚性微分方程.它属于中阶龙格一库塔法; 2 ode23函数用于求解非刚性微分方程,它属于底阶龙格-库塔法; 3 ode113函数用于求解非刚性微分方程的变阶法; 4 ode15s函数用于求解刚性微分方程的变阶法; 5 ode23s函数用于求解刚性微分方程的低阶法。 (2)ode函数的基本格式。上还这五种 ode函数的格式基本上一样.这里主要以最常用的ode45函数为例加以介绍,其结果多数也可用于其他ode函数。ode45的基本格式主要有以下几种。 ① [t,y]=ode45('F',tspan,y0) 其中:'F'是ode文件名字符串; tspan= [T0 TFINAL]为一个向量,它指定了积分的启始时刻和结束时刻; y0为积分运算的初始条件。 该函数对微分方程系统从T0到TFINAL时间段进行积分。函数返回一个列向量。在解向量Y中的每一行与列向量T的返回时间相对应。为了在指定的时刻点T0,T1,…,TFINAL得到解,则使用tspan=[T0,T1,…,TFINAL]。 ② [t , y]=ode45('F', tspan, y0,options) 该函数的解如同前一格式,只是在这里用选项options中的参数来代替前一格式中的缺省值,options中的参数变量由函数odeset予以建立。最常使用的选项参数是标量相对容差'RelTol'(缺省值为1e-3)与绝对容差向量'AbsTol'(缺省的全部向量元素均为1e-6)。 (3)odeset函数的基本格式。odeset函数用来建立或改变options选项结构,其基本格式为 options=odeset('namel',valuel,'name2',value2,…) 通过该odeset函数建立积分器选项options的结构,在options中所命名的属性具有指定值,对于任何未指定的属性均为缺省值。 options=odeset(oldopts,'namel', valuel,…) 通过该odeset函数改变现存选项结构oldopts。 options=odeset(oldopt, newopts) 通过该odeset函数将现有选项结构oldopts与新选项结构newopts结合在一起,新选项中的属性将覆盖对应原有属性。 2.机电系统仿真的ode函数实现 由ode函数可知,当已知系统由微分方程加以描述时,并建立以'F'作为ODE文件名,则可应用 [t,y]=ode45('F',tspan,y0) 实现对系统的求解分析。这里所讲的ODE文件实际就是M函数文件,因此对于以高阶微分方程加以描述的系统,首先应将其写成一组一阶微分方程的形式,当然如果系统是以状态方程描述的,即可直接对其编写M函数文件。函数文件的特征就是文件第一行一定是以function开始,而且在函数文件中的变量均为局部变量,同时由函数文件构造的功能函数可以像MAT-LAB的其他函数一样进行调用。实际上在ode函数[t,y]=ode45('F',tspan,y0)中,对以'F'作为文件名的ODE文件就是这类函数文件。下面通过实例说明对系统仿真过程中的具体应用。 例4-1 之已知一三阶系统,其状态方程表达式为  其中 , 输入为u= 100,试应用ode函数求解该系统。 解: 由于该系统是以状态方程形式给出的,因此即可直接对其编写M函数文件。 function xd= rr(t,x) u=100; xd=[-14*x(1)+9*x(2)+10*x(3)+2.8*u;12*x(1)-9*x(2) +10*x(3)+4*u;24*x(1)-24*x(2)-18*x(3)+12*u]; 该函数文件对所研究系统的状态方程加以描述并形成rr.m文件以便在ode解函数中使用后即可利用ode解函数对系统进行求解分析。考虑状态的零初值条件,即 图4-1 系统状态曲线 clear x0=[0;0;0]; [t,x]=ode45('rr',[0 3],x0); plot(t,x) 对其同样可存一M文件,如li32.m文件,这样即可在命令窗口执行该文件。在运行过程中它会自动调用相关 rr函数文件,实现预定的工作,图 4-1给出了相应的运行结果。 根据ode函数用途可知,该函数应用于求解常微分方程,因此对于含有输出方程的系统状态空间表达式就不能直接利用ode函数求解系统输出响应。考虑到ode函数返回变量的向量格式,可以根据输出方程的具体要求,由状态变量的线性组合直接得到所要求的系统输出。利用ode函数可以方便地对以微分方程形式描述的系统进行仿真分析,然而ode解函数并非是万能的,在应用ode函数时对其进行干预一般是不方便的,因此在很多情况下,对以微分方程形式给出的系统,还必须采用前面所介绍的各种数值积分方法进行仿真分析。 4.2.2 基于 M函数的面向微分方程的系统仿真 M函数是MATLAB中所特有的一类函数,前面介绍ode函数时已经使用了M函数。这里将介绍采用M函数建立描述系统的数学模型并在此基础上编制相应的程序。根据需要调用该M文件,实现所要求的运算。下面通过实例加以说明。 例4-2已知系统状态方程为  应用数值积分方法(预估校正梯形法)实现对系统状态的仿真研究。 解:首先建立描述该系统数学模型的M函数文件,相应的程序如下。 function df=rrl(X,u) A= [-0.5572 -0.7814; 0.7814 0]; B=[1;0] U= 20; df=A*X+B*u 将该函数文件存盘取作rrl.m文件。然后根据预估校正梯形公式编制相应的程序。 h=0.0l;t0=0;tf=15; X=[0 0]'; X1=[0 0]; T=(tf-t0)/h for i=l:T K1=rrl(X,1); X= X+h*K1 K2=rrl(X,1); X= X+h*(K1+K2)/2 X1=[X1;X']; end plot(X1) 在该程序段中两次调用了rr1.m函数文件,在第一次调用rrl.m后,对各状态的计算将产生新值,因此K2就是经过预估计算所得到的相关值,从而即可应用梯形法实现对系统状态的仿真分析。图4-2给出了系统状态响应仿真曲线。 图4-2 系统状态响应 4.3 面向传递函数的连续系统仿真 一个控制系统的动态模型经常是一传递函数的形式给出,特别是当给出开环传递函数, 在进行系统仿真时就要将开环传递函数转换为开环状态方程,再由开环状态方程求出对应的闭环状态方程,然后,再应用数值积分方法进行分解和分析仿真。下面将就这类问题加以介绍并给出MATLAB的实现。 4.3.1 面向传递函数的系统仿真 考虑一个连续系统,它的开环传递函数的一般描述为  其相应的状态空间表达式为  其中A、B、C为相应的系统系数阵,这些系数阵可以写作能控标准型形式或能观标准型形式,当然也可以写成其他的形式,这些都不会影响系统分析的结果。现如果假定系统结构图如图4-3所示,对于SISO系统,其v、u、y 、f 均为标量。由结构图可知u = v-fy,将其代入状态空间表达式得到   图 4-3 系统结构图 式中 (4-1) 输出方程为:  (4-2) (4-1)式即为系统闭环状态方程。因此可以方便地利用各种数值积分方法对其进行求解,实现对以传递函数形式给出的系统进行仿真分析研究。显然,如果系统的数学模型是以微分方程或状态方程形式给出的,一般就可以直接利用前述数值积分方法或ode函数方法对系统进行仿真研究。 4.3.2 面向传递函数的系统仿真的MATLAB实现 下面通过一个例子说明如何利用MATLAB语言进行以传递函数形式给出的系统的仿真分析。 例4-3 设单位反馈系统的开环传递函数为  试给出系统闭环的阶跃响应。 解 根据开环传递函数求出系统开环状态方程,应用MATLAB控制工具箱中的函数,可以方便地求出系统开环状态方程。 num =2*[1,2]; den = conv(conv([1,0],[1,3]),[1,2,2]) G1 = tf(num,den); G2 = ss(G1); [a,b,c,d] = ssdata (G); 这样即可以得到系统系数阵分别为 a = -5.0000 -2.0000 -0.7500 0 4.0000 0 0 0 0 2.0000 0 0 0 0 1.0000 0 b = 1 0 0 0 c = 0 0 0.2500 0.5000 d = 0 进而就可以方便地求出系统闭环状态方程 Ab = a-b*c; Bb = b; Cb = c Db = 0 Ab = -5.0000 -2.0000 -1.0000 -0.5000 4.0000 0 0 0 0 2.0000 0 0 0 0 1.0000 0 Bb = 1 0 0 0 Cb = 0 0 0.2500 0.5000 Db = 0 从而可以应用数值积分方法求取系统的阶跃响应。下面给出了应用MATLAB语言编写的四阶龙格-库塔数值积分方法的仿真程序,相应的仿真曲线由图4-4 给出。 X = [ 0;0;0;0]; Y = 0; t = 0; U = 25; Tf = 15; h=0.01; N = Tf/h; Ab = [ -5 -2 -1 -0.5; 4 0 0 0; 0 2 0 0; 0 0 1 0];  图4-4 四阶龙格-库塔数值积分方法的仿真曲线 Bb = [1; 0; 0;0]; Cb = [0 0 0.25 0.5]; Db = 0; for i = 1:n K1 = Ab * X + Bb*U; K2 = Ab * (X + h*K1/2)+Bb*U; K3 = Ab * (X + h*K2/2)+Bb*U; K4 = Ab * (X + h*K3)+Bb*U; X = X+ h*(K1 + 2*K2+ 2*K3 +K4)/6; Y = [ Y, Cb*X]; t = [t, t(i)+h]; end plot (t,Y) Grid 当由已知系统开环传递函数求取系统开环状态方程时,除了可以利用ss函数直接求取外,还可以使用状态方程相似变换函数ss2ss进行求取,下面给出上例中应用相似变换函数ss2ss求取系统开环状态方程的基本用法。这里主要利用矩阵水平翻转函数fliplr构造出反对角阵作为变换阵,从而得到相应的状态方程。 num = 2 *[1,2]; den = conv(conv([1,0],[1,3]),[1,2,2]); G1= tf (num, den); G2 = ss(G1); T = fliplr (eye(4)); G = ss2ss(G2,T) 运行结果: a = x1 x2 x3 x4 x1 0 1.00000 0 0 x2 0 0 2.00000 0 x3 0 0 0 4.00000 x4 0 -0.75000 -2.00000 -5.00000 b = u1 x1 0 x2 0 x3 0 x4 1.00000 c = x1 x2 x3 x4 y1 0.50000 0.25000 0 0 d = u1 y1 0 从计算结果可以看出,该状态方程十分类似通常定义的系统可控标准型形式。 当然如果在研究系统时,需要状态方程以可控标准型形式给出,则可应用MATLAB编程得到。但不论哪种形式的状态方程,最终所要研究的系统阶跃响应都是一样的。 4.4面向结构图的连续系统仿真 4.4.1 基于典型环节的系统仿真 由于习惯上往往采用结构图形式对系统进行描述,因此在对系统进行仿真时,有必要研究如何以结构图为基础进行。在此将介绍以典型环节为基本环节的面向结构图的系统仿真方法,这种方法的主要优点是: ① 便于研究各环节参数对系统的影响 ② 可以得到每个环节的动态响应; ③ 对于存在系统非线性环节较容易处理; ④ 不但适合单人单出系统,而且对于多入多出系统也适用。 该方法只要建立起各个典型环节的以及各环节之间的连接关系,就能方便地求取任意环节和系统的动态响应。由于直接输入各环节的参数,因此研究参数变化对系统的影响也十分方便。 1.典型环节的确定 控制系统常见的环节有:比例环节 K、积分环节、惯性环节、一阶超前(滞后)环节、二阶振荡环节等。为了仿真的方便,必须选择一种能够代表这些环节的典型环节。经过分析可以看出,只要采用下述形式的环节作为典型环节,就可以满足对系统研究的要求,即典型模型可表示为  (4-3) 显然,当典型环节中时即为惯性环节;时即是积分环节;至于振荡环节则只要用两个典型环节串联并加上一个负反馈即可得到。其中典型环节中的、、、均为对角阵,其对角线元素值由各个环节中的相关参数确定。根据典型环节可以得到第i个环节的输出为  因此得到  (i=1,2,…,n) 将其写为矩阵方程形式,得到 (4-4)  其中 , B= , 2.系统连接矩阵 前面给出了从微观角度对系统环节的描述,然而为了研究一个系统,还必须从宏观角度对系统加以描述,这就需要将各个环节之间的连接关系用一种数学形式加以表示。考察图4-5所示系统,该系统由4个典型环节构成,根据系统结构图可以得到各环节的输入输出关系方程组。 图4-5系统结构图 将其写作矩阵形式  (4-5) 其中 ,  ', ',  一般称与为系统连接矩阵。为一方阵,元素表示结构图中第j 个环节对第i个环节的作用;若第j个环节对第i个环节没有作用,则。表示外部输入信号对系统的连接情况。 3.周环系统的状态方程 将(4-5)式代入(4-4)式,消去U得  对上式展开,整理并取拉氏反变换得  考虑在外输入为阶跃信号时,避免对应的初始导数无穷大,则要求。如果令  (4-6)  (4-7) 可以得到  如果的逆存在,则得到  若令 , (4-8) 则上式可写作  (4-7) 显然上式即为一典型的状态方程,从而可利用前面介绍的数值积分方法对其进行求解,实现对系统的仿真研究。下面以一个结构图形式的系统为例,介绍基于MATLAB的面向结构图的系统仿真实现。 例4-4 已知系统结构图如图4-5所示,采用MATLAB实现面向结构图的系统仿真 。 解:根据系统结构图首先利用MATLAB建立各个环节参数构成的系数阵(此处的系数阵均为对角阵)与系统连接阵,在此基础上再依(4-6)式、(4-7)式、(4-8)式建立闭环系统的相应系数阵。下面给出在MATLAB环境下的相应程序。 AA=diag([1 6.5 1 1]); BB=diag([4 1 3 2 ]); CC=diag([4 1 2.4 5]); DD=diag([0 0 0 1]); W=[0 0 -4 0;1 0 0 -1;0 1 0 0 ;0 0 1 0]; W0=[1 0 0 0 ]’; Q1=BB-DD*W; Q=inv(Q1); R=CC*W-AA; A=Q*R; B=Q*CC*W0; 一旦得到了闭环系统的系数阵,即可应用数值积分方法实现对系统的分析仿真。应注意,在该状态方程中各个状态变量的实际含义。这里应用四阶龙格一库塔法对其进行仿真分析,下面给出了相应的仿真程序。图4-6给出了采用面向结构图系统阶跃响应曲线。 图4-6面向结图仿真(阶跃响应) 图4-7环节参数改变后面向结构图的仿真 X=[0;0;0;0];Y=0;t=0;U=25; Tf=15;h=0.01; n=Tf/h; Cb=[0 01 0]; for i= 1:n K1=A*X+B*U; K2=A*(X+h*K1/2)+ B*U; K3=A*(X+ h*K2/2)+ B*U; K4=A*(X+ h*K3)+ B*U; X=X+h*(K1+2*K2+2*K3+K4)/6; Y=[Y,Cb*X]; t=[t,t(i)+h]; end plot(t,y) grid 在这里由于采用的是面向结构图形式的系统仿真,因此研究因各环节参数发生变化对系统的影响就十分方便,仅对环节系数阵加以修改即可。如果将第三个环节修改为2.4/( s +1)时,只需修改B阵的相关元素即可,其他部分无需改动,图4-7给出了改变环节三参数后的仿真结果。由响应曲线看出,系统的响应快速性明显提高了,同时超调量也有所降低。另外采用面向结构图仿真,还可十分方便地研究各个环节输出的响应情况,只要根据需要修改相应系数阵即可。图4-8(a)、(b)和(c) 图4-8 环节 1(a)、环节2(b)、环节4(c)的响应 分别给出了环节1、环节2和环节4的输出响应。 4.4.2 基于Connect连接函数的系统仿真 除了前面介绍的面向结构图的系统仿真方法外,还可利用MATLAB与控制工具箱所提供的相关函数与命令,实现以结构图形式给出的系统仿真。下面通过实例说明如何以这种方法进行仿真工作。 例4-5 已知系统如图4-9所示,研究该系统的阶跃响应。 解: 显然该系统是以结构图形式给出的较复杂系统;利用MATLAB和控制工具箱可以方便地实现系统仿真。就基于Connect连接函数的仿真而言,一般可分为四步进行。 1.建立系统中间状态模型 首先应对系统各环节进行编号(图4-9),然后即可据此与给定系统参数建立有关中间模型,以便运行时在MATLAB的工作空间建立相应的系统环节模型。其基本语句和格式为 nblocks=10; nl=1;dl= [0.01 1]; n2=[0.17 1];d2=[0.085 0]; n3=1;d3=[0.01 1」; n4=[0.15 1];d4=[0.051 0]; n5=70;d5=[0.0067 1]; 图4-9系统结构图 n6= 0.21;d6= [0.15 1」; n7= 130;d7=[1 0]; n8=-0.212;d8=1; n9= - 0.1 ;d9=[0.011] ; n10=- 0.0044;d10= [0 .01 1」; blkbuild; 其中ni、di表示第i个环节传递函数的分子、分母多项式的系数向量,如果所研究的系统是以状态方程形式给出的各个环节,则应采用ai、bi、ci、di来表示环节状态方程的系数阵。对于每个环节只能取二者之一,MATLAB会自动根据所输入的阵列元素数量加以判定并进行相应的运算。nblocks需根据系统环节个数对其进行赋值, blkbuild语句将对nblocks语句后的各环节输入情况进行校验,并自动生成系统的隐含中间状态方程模型(a、b、c、d系数阵)。 2.建立系统连接矩阵 在MATLAB工作空间建立了系统各环节的模型后,还应该建立反映各环节之间相互连接关系的Q阵,以及输入、输出向量(inputs、outputs),下面结合实例将Q阵与inputs、outputs列写规则加以介绍。根据已知系统可以得到系统连接阵Q和输入、输出向量: Q=[1 0 0; 2 1 10; 3 2 0; 4 3 9; 5 4 0; 6 5 8; 7 6 0; 8 7 0; 9 6 0; 10 7 0]; inputs=1; outputs=7; Q阵的行数与系统环节数相同,Q阵的第一列应为相应环节编号,第i行的第二列及以后的所有元素为进入第i个环节的所有环节编号。输入变量inputs为输入信号进入的环节编号;输出变量ouputs为系统信号输出所对应的环节。如果是多输入多输出系统,输入变量和输出变量应为一向量。 3.建立整体系统模型 完成了上述工作之后,就可以应用connect函数建立整体系统的数学模型。其基本格式为 [ A, B, C,D]=connect(a,b,c,d,Q,inputs,outputs) 或 G= connect(a, b, c, d, Q, inputs, outputs) 在connect函数中的a、b 、c、d系数阵是由blkbuild函数自动生成的,其返回变量既可是系数阵,也可是系统阵G,不论哪种情况都将生成完整的系统模型,从而就可利用它们对该复杂系统进行分析研究。 4.系统仿真 由于已经求出了系统的系数阵,为了获得系统阶跃响应,可以利用Step(A,B,C,D)求取系统的阶跃相应,图4-10即为系统的阶跃响应曲线。当然也可应用数值积分方法对其进行研究,可参考前面内容。  图4-10系统阶跃响应 4.5 Matlab中的系统仿真工具SIMULINK应用 4.5.1 Simulink应用简介 Simulink是MATLAB下的面向结构图方式的仿真环境,对于结构复杂的控制系统,要快速地建立系统模型是较为困难的,Mathworks公司所提供的这种采用图形化方式的对系统建模和仿真的工具为用户带来了极大的方便。Simulink与用户交互接口是基于Windows的图形编程方法,因此非常易于为用户所接受,使用十分灵活方便。特别是它的S-Function使其更增强了功能扩展能力。目前Simulink作为MATLAB的重要组成部分,已经成为仿真研究的重要工具。下面将简要介绍有关Simulink的基本应用。 在进入MATLAB后,键人Simulink命令,就可以打开一个名为Simulink的模型库窗口,如图4-11所示。 图4-12 信号源Source 库 图4-13 输出方式(Sinks)库 Simulink模型库包含功能齐全的子模型库:Sources (信号源库)、Sinks输出方式库)、Discrete(离散模型库)、Linear(线性环节库)、Nonlinear(非线性环节库)、Connections(连接及接口库)、Blocksets and toolboxes(模块建立和工具箱库)和 Demos(示例库)。每个子模型库包含很多的功能模块,双击这些模型库即可打开它们,图4-12、图4-13、图4-14、图4-15分别给出了双击打开的信号源库、输出方式库、线性环节库以及连接模块库等。 每个子模块都可以从相应的库被鼠标拖(即复制)到用户建立的模型窗口中。在用户窗口内,可以用鼠标对各个模块进行连接,并相应设置子模块参数、仿真算法参数,以便进行仿真。设计好的系统模型应加以存盘,在 MATLAB 5.x中 Simulink的文件名以.mdl作为后缀,这一点是与 MATLAB 4.x不同的。为了说明 Simulink的使用,下面以简单示例介绍在 Simulink环境中建立系统模型以及对其进行仿真分析的方法。 例4-6已知系统开环传递函数为 试分析闭环系统的阶跃响应。 解 首先在Simulink的Library窗口中的File中选择New,打开一个新工作平台,然后分别在信号源库、输出方式库以及线性环节库中找到 Step、Scope、Transfer Fun以及 Sum各功能块,将其拖至工作平台并按要求加以连接。在对相关参数进行正确设置后,即可启动仿真命令,如图4-16所示,从而自动得到由图4-17给出的系统仿真结果。 有关Simulink的详细应用请参阅Mathworks的simulink使用手册。 图4-14 线性环节模块(Linear)库 图4-15 连接及接口模块(Connection)库  图4-16 在Simulink 中建立仿真系统 图4-17 Simulink仿真结果 4.5.2 S-Function (S函数)应用 前面已经介绍过Simulink是一类面向结构图的仿真环境,然而由于Simulink提供了S-Function(S函数)的功能,从而极大地提高和丰富了Simulink对系统进行分析与研究的能力。S-Function通过MATLAB语言进行编程与Simulink中的S-Function功能模块有机地结合在一起,实现了用户所需的各种运算和分析功能,因此说 S-Function是Simulink中极为重要的一部分内容。S-Function的使用必须严格按照所规定的格式加以进行,S-Function的运行走向是由Flag标志加以控制的,其实现的功能主要由sys参数及其相关的功能语句块所决定。在高版本(5.x)的MATLAN中,用户可以利用模板程序进行编程,保留所需功能,删掉多余的内容,从而简化了编程工作量,保证了程序编制质量,S-Function的模版程序可以在Blocks目录中找到。 在S-Function中根据Flag提供的标志,执行6种不同的任务(case0,case1,case2,case3,case4与case9),其中最为重要的是下述三种: Case 0, [sys,x0,str, ts]=mdllnitializeSizes; ease l, sys=mdlDerivatives(t,x,u); Case 3, sys=mdloutputs(t,x,u); 它们分别为S-Function提供进行结构初始化设置、状态微分和状态输出功能。因此,欲灵活地使用S-Function,就必须首先掌握它们的应用。在mdllnitializeSizes的返回变量中sys是根据Flag的标志运算得到的解,x0为初始状态,str对于M类文件设置为空阵(S-Function也可由其他类文件构成),ts为列向量,定义了取样时间与偏移量。 下面通过一个示例对S-Function的具体应用加以说明 例4-6 已知系统状态方程的系数阵为 试应用S-Function求取系统的阶跃响应。 解:应用S-Function一般可分为三步进行。①在Simulink窗口构造所需结构图,其中必须包括S-Function功能块。②根据用户的需要编写S-Function程序,这是最为重要的一步,当然可以借助模版程序加以进行。③编写相关的主运行程序。 图4-18 系统结构图 首先建立结构图,如图4-18所示。在结构图中S-Function模块是关键,该模块名应与所编制的S-function程序名相同;To Workspace模块将计算数据写入工作空间;Fcn模块为自定义函数模块,用户可根据需要自行定义。然后可利用模版根据要求编制如下相应程序。 S-function程序:ssq.m function[sys,x0,str,ts]=ssq(t, x,u,flag) A=[-18 19 -24;19 -21 20;44 -24 -36]; B=[0 l l」’; C=[1 0 2]; switch flag, case 0, [sys, x0, str, ts]= mdlInitializeSizes(A,B,C); case 1, sys= mdlDerivatives(t ,x、u, A, B, C); case 2, sys=[ ]; case 3, sys=mdlOutputs(t,x,u,A,B,C); sys=simsizes(sizes); x0=[0 0 0]; str=[ ]; ts=[ 0 0 ]; function sys=mdlDerivatives(t,x,u, A,B,C) ys= A* x + B*u; function sys= mdlUpdate(t,x,u) sys=[ ]; function case 4, sys=[ ]; case 9, sys=[ ]; otherwise error(['Unhandled flag=',num2str(flag)]); end function[sys,x0,str, ts] =mdlInitializeSizes(A,B,C) sizes=simsizes; sizes.NumContstates =3; sizes.NumDiscStates =0; sizes.NumOutputs =1; sizes.NumInputs =1; sizes.DirFeedthrough =1; sizes.NumSampleTimes =1;% at least one sample time is needed sys=mdlOuputs(t,x,u,A,B,C) sys=C* x; function sys=mdlGetTimeOfNextVarHit(t,x,u) sampleTime= 1; sys= t+ sampleTime; function sys=mdlTerminate(t, x, u) sys=[ ]; 无关函数一律赋空阵,函数是根据用户的具体情况加以设置的,因为本系统有三个状态变量,所以对 NumContStates设置 3 ;又因为是单输入单输出系统,所以 NumOutPuts=1, NumInputs=1。 最后编制主程序,在这里主程序的功能就是运行该S-Function并将结果绘制出来。 主程序为 sim('sspaceq') plot(t,simout) 运行主程序后,图4-17示出了系统运行结果。 这里仅给出了一个简单的示例说明S-Function的应用,显然这种功能可以推广到各种复杂的系统,请读者自行加以尝试。 表4-1至表4-7分别给出了Simulink的基本模块库及其相应的功能说明。 无关函数一律赋空阵,函数是根据用户的具体情况加以设置的,因为本系统有三个状态变量,所以对 NumContStates设置 3 ;又因为是单输入单输出系统,所以 NumOutPuts=1, NumInputs=1。 最后编制主程序,在这里主程序的功能就是运行该S-Function并将结果绘制出来。 主程序为 sim('sspaceq') plot(t,simout) 运行主程序后,图4-17示出了系统运行结果。 这里仅给出了一个简单的示例说明S-Function的应用,显然这种功能可以推广到各种复杂的系统,请读者自行加以尝试。 表4-1至表4-7分别给出了Simulink的基本模块库及其相应的功能说明。  图4-19 利用S-Function的仿真结果 表4-1 信号源库 模块名  功能说明  模块名  功能说明   Constant  常量  Chirp Signal  脉冲信号   Signal Genarator  信号发生器  Clock  时钟信号   Step  阶跃信号  Digital Clock  数字时钟信号   Ramp  斜坡信号  From File  来源于文件数据   Sine Wave 正弦信号  From Workspace  来源于工作空间   Repeating Sequence  频率信号  Random Number  随机数  Discrete Pluse Generator  离散脉冲信号源 Uniform Random Number  均匀分布随机数   Pluse Generator  脉冲信号源 Band-Limited White Noise  有界白噪音   表4-2 输出方式库 模块名  功能说明  模块名  功能说明   Scope  示波器  To File  存到文件   XY Graph  XY图形窗口  To Workspace  记录到工作空间   Display  数字显示  Stop Simulation  停止仿真   表4-3 离散模型库 模块名  功能说明  模块名  功能说明   Unit Delay  单位延迟 Discrete State-Space  离散状态方程  Discrete Time Integrator 离散时间积分器  Discrete Filter  滤波器   Zero-order Hold  零阶保持器 Discrete Transter Fcn  离散传递函数   First-order Hold  一阶保持器  Discrete Zero-Pole  离散零极点方程   表4-4 线性环节模型库 模块名  功能说明  模块名  功能说明   Gain  增益  Derivative  微分器   Sum  加法器  Dot Product  点乘   Integrator  积分器  Mrtrix Gain  增益矩阵   State-Space  状态方程  Slider Gain  可变增益   Zero-Pole 零极点方程     表4-5 连接模型库 模块名  功能说明  模块名  功能说明   In  输入  Data Store Memory  数据存储单元   Out  输出  Data Store Write  数据存储写   Mux  综合  Enable  使能   Denux  综合  Trigger  触发   From  从变量读取  Subsystem  子系统   Goto  写到变量  Selector  选择器   Data Store Read  数据存储器     表4-6 模块建立和工具箱模块库 模块名  功能说明  模块名  功能说明   Controls Toolbox  控制工具箱  NNET Blocks  神经网络模块   DSP Bloxkset  信号处理模块  Stateflow  状态流程   Fuzzy  模糊模块  System ID Block  系统辨别模块   NCD Blocks  非线性控制模块     表4-7 非线性环节模型库 模块名  功能说明  模块名  功能说明   Abs  绝对值  Backlash  回滞环节  Trigonometric Function  加法器  Dead Zone  死区环节   Math Function  积分器  Memory  存储环节   Min-Max  状态方程  Transport Delay  传递延迟   Zero-Pole  零极点方程 Variable Transport Delay  变量延迟   Product  微分器  Fcn  函数   Logic Operator  点乘  S-Function  S函数  Relational Operator  增益矩阵  Switch  开关   Sign  可变增益  Manual Switch  自定义开关   Rate Limiter  速度限制环节  Relay  继点环节   Saturation  饱和环节  Multiport Switch  综合开关   Quantizer  编码器  MATLAB Fcn  MATLAB函数