【分享】S-Function 使用及应用举例

chinamaker

对用simulink做仿真的有用：  
  
S-FUNCTIONS的书写之一  
  
s-function也就是system-function的缩写。说得简单，s-function就是用MATLAB所提供的模型不能完全满足用户，而提供给用户自己编写程序来满足自己要求模型的接口。要了解 s-function，必须了解以下知识：  
（1）direct feedthrough  
（2）dynamically sized inputs  
（3）setting sample times and offsets  
[由于我也不知道怎么把上面三句话精确的翻译成中文，因此在此不加翻译，大家自己理解，也许更好]由于上面三部分的重要性，在此详细进行分析。  
  
一.direct feedthrough  
direct feedthrough意思是说系统的输出或可变采样时间是否受到输入的控制。大家清楚有的系统是受到输入控制如：  
y=k*u (u是输入，k是放大因子，y是输出）  
而有的系统输出是不受到输入影响，如：  
输出：y=x  
dx=u  
x表示状态  
  
二.dynamically sized inputs  
主要是给出：输入连续状态数目（size.NumContStates),离散状态数目（size.NumDiscStates)  
,输出数目（size.NumOutputs),输入数目(size.NumInputs),Direct Feedthrough(size.Dir Feedthrough)。  
三.setting sample times and offsets  
  
setting smaple times and offsets主要设置采样时间 MATLAB为了用户方便，已经书写了S-FUNCTIONS模板函数sfuntmpl.m。为了更好的写S-FUNCTIONS，大家来看一下，该函数sfuntmpl.m内容如下：（我通过在该内容加注释来说明,以$开头）  
  
function [sys,x0,str,ts] = sfuntmpl(t,x,u,flag)  
$输出变量就此四个，大家必须注意它的顺序。$输入变量可以为t,x,u,flag,p1,...,pn等，但是前面的四个变量不能变，特此说明。  
  
%SFUNTMPL General M-file S-function template  
% With M-file S-functions, you can define you own ordinary differential  
% equations (ODEs), discrete system equations, and/or just about  
% any type of algorithm to be used within a Simulink block diagram.  
%  
$上面是其功能  
% The general form of an M-File S-function syntax is:  
% [SYS,X0,STR,TS] = SFUNC(T,X,U,FLAG,P1,...,Pn)  
%  
$调用格式  
% What is returned by SFUNC at a given point in time, T, depends on the  
% value of the FLAG, the current state vector, X, and the current  
% input vector, U.  
%  
% FLAG RESULT DESCRIPTION  
% ----- ------ --------------------------------------------  
% 0 [SIZES,X0,STR,TS] Initialization, return system sizes in SYS,  
% initial state in X0, state ordering strings  
% in STR, and sample times in TS.  
$具体怎样实现，大家参看后面的函数mdlInitializeSizes  
% 1 DX Return continuous state derivatives in SYS.  
% 2 DS Update discrete states SYS = X(n+1)  
% 3 Y Return outputs in SYS.  
% 4 TNEXT Return next time hit for variable step sample  
% time in SYS.  
% 5 Reserved for future (root finding).  
% 9 [] Termination, perform any cleanup SYS=[].  
%  
$参看后面相应函数  
%  
% The state vectors, X and X0 consists of continuous states followed  
% by discrete states.  
%  
% Optional parameters, P1,...,Pn can be provided to the S-function and  
% used during any FLAG operation.  
%  
% When SFUNC is called with FLAG = 0, the following information  
% should be returned:  
%  
% SYS(1) = Number of continuous states.  
% SYS(2) = Number of discrete states.  
% SYS(3) = Number of outputs.  
% SYS(4) = Number of inputs.  
% Any of the first four elements in SYS can be specified  
% as -1 indicating that they are dynamically sized. The  
% actual length for all other flags will be equal to the  
% length of the input, U.  
% SYS(5) = Reserved for root finding. Must be zero.  
% SYS(6) = Direct feedthrough flag (1=yes, 0=no). The s-function  
% has direct feedthrough if U is used during the FLAG=3  
% call. Setting this to 0 is akin to making a promise that  
% U will not be used during FLAG=3. If you break the promise  
% then unpredictable results will occur.  
% SYS(7) = Number of sample times. This is the number of rows in TS.  
%  
$需要说明的是sys的顺序不能乱  
%  
% X0 = Initial state conditions or [] if no states.  
%  
% STR = State ordering strings which is generally specified as [].  
%  
% TS = An m-by-2 matrix containing the sample time  
% (period, offset) information. Where m = number of sample  
% times. The ordering of the sample times must be:  
%  
% TS = [0 0, : Continuous sample time.  
% 0 1, : Continuous, but fixed in minor step  
% sample time.  
% PERIOD OFFSET, : Discrete sample time where  
% PERIOD > 0 & OFFSET 〈 PERIOD.  
% -2 0]; : Variable step discrete sample time  
% where FLAG=4 is used to get time of  
% next hit.  
%  
% There can be more than one sample time providing  
% they are ordered such that they are monotonically  
% increasing. Only the needed sample times should be  
% specified in TS. When specifying than one  
% sample time, you must check for sample hits explicitly by  
% seeing if  
% abs(round((T-OFFSET)/PERIOD) - (T-OFFSET)/PERIOD)  
% is within a specified tolerance, generally 1e-8. This  
% tolerance is dependent upon your model's sampling times  
% and simulation time.  
%  
% You can also specify that the sample time of the S-function  
% is inherited from the driving block. For functions which  
% change during minor steps, this is done by  
% specifying SYS(7) = 1 and TS = [-1 0]. For functions which  
% are held during minor steps, this is done by specifying  
% SYS(7) = 1 and TS = [-1 -1].  
  
% Copyright (c) 1990-1998 by The MathWorks, Inc. All Rights Reserved.  
% $Revision: 1.12 $  
  
%  
% The following outlines the general structure of an S-function.  
%  
switch flag,  
   %%%%%%%%%%%%%%%%%%  
   % Initialization %  
   %%%%%%%%%%%%%%%%%%  
   case 0,[sys,x0,str,ts]=mdlInitializeSizes;  
$大家是不是觉得此函数名太长，当然可以根据自己的爱好加以改变，不过后面的相应改。  
   %%%%%%%%%%%%%%%  
   % Derivatives %  
   %%%%%%%%%%%%%%%  
   case 1,sys=mdlDerivatives(t,x,u);  
   %%%%%%%%%%  
   % Update %  
   %%%%%%%%%%  
   case 2,sys=mdlUpdate(t,x,u);  
   %%%%%%%%%%%  
   % Outputs %  
   %%%%%%%%%%%  
   case 3,sys=mdlOutputs(t,x,u);  
   %%%%%%%%%%%%%%%%%%%%%%%  
   % GetTimeOfNextVarHit %  
   %%%%%%%%%%%%%%%%%%%%%%%  
   case 4,sys=mdlGetTimeOfNextVarHit(t,x,u);  
   %%%%%%%%%%%%%  
   % Terminate %  
   %%%%%%%%%%%%%  
   case 9,sys=mdlTerminate(t,x,u);  
   %%%%%%%%%%%%%%%%%%%%  
   % Unexpected flags %  
   %%%%%%%%%%%%%%%%%%%%  
   otherwise, error(['Unhandled flag = ',num2str(flag)]);  
end  
% end sfuntmpl  
  
%  
%=============================================================================  
% mdlInitializeSizes  
% Return the sizes, initial conditions, and sample times for the S-function.  
%=============================================================================  
%  
function [sys,x0,str,ts]=mdlInitializeSizes  
%  
% call simsizes for a sizes structure, fill it in and convert it to a  
% sizes array.  
%  
% Note that in this example, the values are hard coded. This is not a  
% recommended practice as the characteristics of the block are typically  
% defined by the S-function parameters.  
%  
$关于函数simsizes大家必须遵循，因为把是内部函数，不得随便改变，其作用是返回未初始化的size结构。  
sizes = simsizes;  
$number of continuous states  
sizes.NumContStates = 0;  
$number of discrete states  
sizes.NumDiscStates = 0;  
$number of outputs  
sizes.NumOutputs = 0;  
$ number of inputs  
sizes.NumInputs = 0;  
$Flag for direct feedthrough  
sizes.DirFeedthrough = 1;  
$number of sample times  
sizes.NumSampleTimes = 1;  
% at least one sample time is needed  
$  
sys = simsizes(sizes);  
%  
% initialize the initial conditions  
%  
x0 = [];  
%  
% str is always an empty matrix  
%  
str = [];  
%  
% initialize the array of sample times  
%  
ts = [0 0];  
% end mdlInitializeSizes  
%  
%=============================================================================  
% mdlDerivatives  
% Return the derivatives for the continuous states.  
%=============================================================================  
%  
function sys=mdlDerivatives(t,x,u)  
sys = []; % end mdlDerivatives  
%  
%=============================================================================  
% mdlUpdate  
% Handle discrete state updates, sample time hits, and major time step  
% requirements.  
%=============================================================================  
%  
function sys=mdlUpdate(t,x,u)  
sys = []; % end mdlUpdate  
%  
%=============================================================================  
% mdlOutputs  
% Return the block outputs.  
%=============================================================================  
%  
function sys=mdlOutputs(t,x,u)  
sys = []; % end mdlOutputs  
%  
%=============================================================================  
% mdlGetTimeOfNextVarHit  
% Return the time of the next hit for this block. Note that the result is  
% absolute time. Note that this function is only used when you specify a  
% variable discrete-time sample time [-2 0] in the sample time array in  
% mdlInitializeSizes.  
%=============================================================================  
%  
function sys=mdlGetTimeOfNextVarHit(t,x,u)  
sampleTime = 1; % Example, set the next hit to be one second later.  
sys = t + sampleTime; % end mdlGetTimeOfNextVarHit  
%  
%=============================================================================  
% mdlTerminate  
% Perform any end of simulation tasks.  
%=============================================================================  
%  
function sys=mdlTerminate(t,x,u)  
sys = []; % end mdlTerminate  
  
在下面我将写几个具体的实例。  
  
   之二、连续系统例子：   
  
function [sys,x0,str,ts] = csfunc(t,x,u,flag)  
%CSFUNC An example M-file S-function for defining a continuous system.  
% Example M-file S-function implementing continuous equations:  
% x' = Ax + Bu  
% y = Cx + Du  
%  
% See sfuntmpl.m for a general S-function template.  
%  
% See also SFUNTMPL.  
     
% Copyright (c) 1990-1998 by The MathWorks, Inc. All Rights Reserved.  
% $Revision: 1.5 $  
A=[-0.09 -0.01; 1 0];  
B=[ 1 -7; 0 -2];  
C=[ 0 2; 1 -5];  
D=[-3 0; 1 0];  
switch flag,  
   %%%%%%%%%%%%%%%%%%  
   % Initialization %  
   %%%%%%%%%%%%%%%%%%  
   case 0,[sys,x0,str,ts]=mdlInitializeSizes(A,B,C,D);  
   %%%%%%%%%%%%%%%  
   % Derivatives %  
   %%%%%%%%%%%%%%%  
   case 1,sys=mdlDerivatives(t,x,u,A,B,C,D);  
   %%%%%%%%%%%  
   % Outputs %  
   %%%%%%%%%%%  
   case 3,sys=mdlOutputs(t,x,u,A,B,C,D);  
   %%%%%%%%%%%%%%%%%%%  
   % Unhandled flags %  
   %%%%%%%%%%%%%%%%%%%  
   case { 2, 4, 9 },sys = [];  
   %%%%%%%%%%%%%%%%%%%%  
   % Unexpected flags %  
   %%%%%%%%%%%%%%%%%%%%  
   otherwise, error(['Unhandled flag = ',num2str(flag)]);  
end  
% end csfunc  
%  
%=============================================================================  
% mdlInitializeSizes  
% Return the sizes, initial conditions, and sample times for the S-function.  
%=============================================================================  
%  
function [sys,x0,str,ts]=mdlInitializeSizes(A,B,C,D)  
sizes = simsizes;  
sizes.NumContStates = 2; sizes.NumDiscStates = 0; sizes.NumOutputs = 2;  
sizes.NumInputs = 2; sizes.DirFeedthrough = 1; sizes.NumSampleTimes = 1;  
sys = simsizes(sizes); x0 = zeros(2,1); str = []; ts = [0 0];  
% end mdlInitializeSizes  
%  
%=============================================================================  
% mdlDerivatives  
% Return the derivatives for the continuous states.  
%=============================================================================  
%  
function sys=mdlDerivatives(t,x,u,A,B,C,D)  
sys = A*x + B*u; % end mdlDerivatives  
%  
%=============================================================================  
% mdlOutputs  
% Return the block outputs.  
%=============================================================================  
%  
function sys=mdlOutputs(t,x,u,A,B,C,D)  
sys = C*x + D*u; % end mdlOutputs  
  
   之三、有关离散系统   
  
关于离散系统定义的书写在此举一例：  
  
function [sys,x0,str,ts] = dsfunc(t,x,u,flag)  
%DSFUNC An example M-file S-function for defining a discrete system.  
% Example M-file S-function implementing discrete equations:  
% x(n+1) = Ax(n) + Bu(n)  
% y(n) = Cx(n) + Du(n)  
%  
% See sfuntmpl.m for a general S-function template.  
%  
% See also SFUNTMPL. % Copyright (c) 1990-1998 by The MathWorks, Inc. All Rights Reserved.  
% $Revision: 1.13 $  
% Generate a discrete linear system:  
A=[-1.3839 -0.5097 1.0000 0]; B=[-2.5559 0 0 4.2382];  
C=[ 0 2.0761 0 7.7891]; D=[ -0.8141 -2.9334 1.2426 0];  
switch flag,  
% Initialization  
case 0,[sys,x0,str,ts] = mdlInitializeSizes(A,B,C,D);  
% Update % %%%%%%%%%%  
case 2, sys = mdlUpdate(t,x,u,A,B,C,D); %%%%%%%%%%  
% Output % %%%%%%%%%%  
case 3, sys = mdlOutputs(t,x,u,A,C,D); %%%%%%%%%%%%%  
% Terminate % %%%%%%%%%%%%%  
case 9, sys = []; % do nothing %%%%%%%%%%%%%%%%%%%%  
% Unexpected flags % %%%%%%%%%%%%%%%%%%%%  
otherwise, error(['unhandled flag = ',num2str(flag)]);  
end %end dsfunc  
%  
%=======================================================================  
% mdlInitializeSizes  
% Return the sizes, initial conditions, and sample times for the S-function.  
%=======================================================================  
%  
function [sys,x0,str,ts] = mdlInitializeSizes(A,B,C,D)  
sizes = simsizes;  
sizes.NumContStates = 0; sizes.NumDiscStates = size(A,1); sizes.NumOutputs = size(D,1);  
sizes.NumInputs = size(D,2); sizes.DirFeedthrough = 1; sizes.NumSampleTimes = 1;  
sys = simsizes(sizes); x0 = ones(sizes.NumDiscStates,1); str = []; ts = [1 0];  
% end mdlInitializeSizes  
%  
%=======================================================================  
% mdlUpdate  
% Handle discrete state updates, sample time hits, and major time step  
% requirements.  
%=======================================================================  
% function sys = mdlUpdate(t,x,u,A,B,C,D)  
sys = A*x+B*u; %end mdlUpdate  
%  
%=======================================================================  
% mdlOutputs  
% Return Return the output vector for the S-function  
%=======================================================================  
%  
function sys = mdlOutputs(t,x,u,A,C,D)  
sys = C*x+D*u; %end mdlUpdate  
  
   S-FUNCTIONS的书写之四（离散和连续的混合型）   
  
　  
  
function [sys,x0,str,ts] = mixedm(t,x,u,flag)  
%MIXEDM An example integrator followed by unit delay M-file S-function  
% Example M-file S-function implementing a hybrid system consisting  
% of a continuous integrator (1/s) in series with a unit delay (1/z).  
% Sampling period and offset for unit delay.  
dperiod = 1; doffset = 0;  
switch flag  
%%%%%%%%%%%%%%%%%%  
% Initialization %  
%%%%%%%%%%%%%%%%%%  
case 0, [sys,x0,str,ts]=mdlInitializeSizes(dperiod,doffset);  
%%%%%%%%%%%%%%%  
% Derivatives %  
%%%%%%%%%%%%%%%  
case 1, sys=mdlDerivatives(t,x,u);  
%%%%%%%%%%  
% Update %  
%%%%%%%%%%  
case 2, sys=mdlUpdate(t,x,u,dperiod,doffset);  
%%%%%%%%%%  
% Output %  
%%%%%%%%%%  
case 3, sys=mdlOutputs(t,x,u,doffset,dperiod);  
%%%%%%%%%%%%%  
% Terminate %  
%%%%%%%%%%%%%  
case 9, sys = [];  
% do nothing  
otherwise, error(['unhandled flag = ',num2str(flag)]);  
end % end mixedm %  
%=============================================================================  
% mdlInitializeSizes  
% Return the sizes, initial conditions, and sample times for the S-function.  
%=============================================================================  
%  
function [sys,x0,str,ts]=mdlInitializeSizes(dperiod,doffset)  
sizes = simsizes; sizes.NumContStates = 1; sizes.NumDiscStates = 1;  
sizes.NumOutputs = 1; sizes.NumInputs = 1; sizes.DirFeedthrough = 0;  
sizes.NumSampleTimes = 2; sys = simsizes(sizes); x0 = ones(2,1); str = [];  
ts = [0 0; % sample time  
       dperiod doffset];  
% end mdlInitializeSizes  
%  
%=============================================================================  
% mdlDerivatives % Compute derivatives for continuous states.  
%=============================================================================  
%  
function sys=mdlDerivatives(t,x,u)  
sys = u; % end mdlDerivatives  
%  
%=============================================================================  
% mdlUpdate  
% Handle discrete state updates, sample time hits, and major time step % requirements.  
%=============================================================================  
%  
function sys=mdlUpdate(t,x,u,dperiod,doffset)  
% next discrete state is output of the integrator  
if abs(round((t - doffset)/dperiod) - (t - doffset)/dperiod) < 1e-8, sys = x(1);  
else sys = []; end % end mdlUpdate  
%  
%=============================================================================  
% mdlOutputs % Return the output vector for the S-function  
%=============================================================================  
%  
function sys=mdlOutputs(t,x,u,doffset,dperiod)  
% Return output of the unit delay if we have a  
% sample hit within a tolerance of 1e-8. If we  
% don't have a sample hit then return [] indicating  
% that the output shouldn't change.  
if abs(round((t - doffset)/dperiod) - (t - doffset)/dperiod) < 1e-8, sys = x(2); else sys = []; end % end mdlOutputs 关于这些程序的说明，我将在后面讲解。  
  
   S-Function书写之五   
  
在本帖中我对前面列举的关于系统是连续，离散，连续与离散混合的三个例子加以说明，以至于大家在看下面的例子时能更好的理解。  
  
一.函数mdlInitializeSizes  
  
对于描述连续型的函数csfunc.m  
  
sizes.NumContStates = 2; sizes.NumDiscStates = 0;  
  
$表明本函数是描述连续型的  
  
sizes.NumOutputs = 2; sizes.NumInputs = 2;  
sizes.DirFeedthrough = 1; sizes.NumSampleTimes = 1;  
  
对于描述离散型的函数dsfunc.m  
  
sizes.NumContStates = 0; sizes.NumDiscStates = size(A,1);  
  
$表明本函数是描述离散型的  
  
sizes.NumOutputs = size(D,1); sizes.NumInputs = size(D,2);  
sizes.DirFeedthrough = 1; sizes.NumSampleTimes = 1;  
  
对于描述混合型的函数mixedm.m  
  
sizes.NumContStates = 1; sizes.NumDiscStates = 1;  
  
$表明本函数是描述混合型的  
  
sizes.NumOutputs = 1; sizes.NumInputs = 1;  
sizes.DirFeedthrough = 0; sizes.NumSampleTimes = 2;  
  
关于变量sizes.numoutputs和sizes.numinputs则由描述的系统的输出，输入而定。  
关于变量sizes.DirFeedthrough,则可以查看矩阵D,如果D是非空，则为1，否则为0  
关于变量sizes.NumSampleTimes，即是ts这个矩阵的行数。  
  
二.函数mdlDerivatives  
  
仅仅出现在含有连续型的系统中；这部分应该根据系统状态来决定。  
  
三.函数mdlUpdate 仅仅出现在含有离散型的系统中；该函数应该根据系统状态来决定  
  
四.函数mdlOutputs 该函数在各种类型中都有,应该根据描述系统的输出来决定;  
  
五.关于书写s-functions函数的注意之点  
  
(1)我想大家都会觉得上述三个函数都不利于维护,因此在写s-functions函数时不妨把矩阵 A,B,C,D也考虑成变量,作为输入参数;  
(2)关于在csfuncs.m中的型如case(2,4,9)这种形式,大家最好不要采用,不利于维护,可以分 开来写  
(3)关于结构sizes中的各个变量的值,最好从输入参量中得到,而不是通过人为判断来输入  
(4)因此,其实上述三种形式,可以变化成一种形式.  
(5)可以通过对相应于你所要求的系统对上述三个函数加以相应的改进来到达自己的要求,因此可以把上述函数当成模板函数.  
http://www.dytrol.com

lhuaqing

谢谢，请问如果要得到多个输出应该怎么设置，我知道多输入可以将多个输入值合并为一个向量，但多输出的时候就不知道该如何设置了，难道也是合并为一个向量吗，该如何设置，谢谢！

duoting

先谢谢过楼主，考下来慢慢看

qiyulong

好像一些教程上都有这么一大段

小星星

先谢谢过楼主，考下来慢慢看