基于S函數的simulink仿真
S函數可以用計算機語言來描述動態系統。在控制系統設計中，S函數可以用來描述控制算法、自適應算法和模型動力學方程。
S函數中使用文本方式輸入公式和方程，適合復雜動態系統的數學描述，并且在仿真過程中可以對仿真參數進行更精確的描述、
1.1 S函數簡介
S函數是系統函數（system function）的簡稱。可以用MATLAB代碼、C、C++等語言來編寫S函數。
1.2 S函數的使用步驟
步驟如下：
- 創建S函數源文件
- 在動態系統的simulink模型框圖中添加S-function模塊，并且進行正確設置
- 在simulink模型框圖中按照定義好的功能連接輸入輸出端口
1.3 S函數的基本功能及重要參數設定
S函數的基本功能及重要參數設定如下：
- S函數功能模塊：各種功能模塊完成不同的任務，這些功能模塊（函數）稱為仿真例程或回調函數（call - back functions），包括初始化（initialization）、導數（mdlDerivative）、輸出（mdlOutput）等。
- NumContStates表示S - 函數描述的模塊中連續狀態的個數。
- NumDiscStates表示離散狀態的個數。
- NumOutputs和NumInputs分別表示模塊輸出和輸入的個數。
- 直接饋通（dirFeedthrough）為輸入信號是否在輸出端出現的標識，取值為0或1。例如，形如 y = k × u y = k×u y=k×u的系統需要輸入（即直接反饋），其中， u u u是輸入， k k k是增益， y y y是輸出，形如等式 y = x , x ˙ = u y = x,\dot{x}=u y=x,x˙=u的系統不需要輸入（即不存在直接反饋），其中， x x x是狀態， u u u是輸入， y y y為輸出。
- NumSampleTimes為模塊采樣周期的個數，S函數支持多采樣周期的系統。 除了sys外，還應設置系統的初始狀態變量 x 0 x_0 x0?、說明變量str和采樣周期變量 t s t_s ts?。 t s t_s ts?變量為雙列矩陣，其中每一行對應一個采樣周期。對連續系統和單個采樣周期的系統來說，該變量為 [ t 1 , t 2 ] [t_1,t_2] [t1?,t2?]， t 1 t_1 t1?為采樣周期， t 1 = ? 1 t_1 = - 1 t1?=?1表示繼承輸入信號的采樣周期， t 2 t_2 t2?為偏移量，一般取為0。對連續系統來說， t s t_s ts?取為 [ ? 1 , 0 ] [-1,0] [?1,0]。
1.4 S函數描述實例
在控制系統設計中，S函數可以用于控制器、自適應律和模型描述。
以模型 J θ ¨ = u + d ( t ) J\ddot{\theta}=u+d(t) Jθ¨=u+d(t)為例，其中， u u u為控制輸入， d ( t ) d(t) d(t)為加在控制輸入端的擾動，模型輸出為 θ 和 θ ˙ \theta和\dot{\theta} θ和θ˙，即轉動角度和角速度， J J J為轉動慣量，該模型可以描述如下：
x ˙ 1 = x 2 x ˙ 2 = 1 J ( u + d ( t ) ) \begin{align*} \dot{x}_1&=x_2\\ \dot{x}_2&=\frac{1}{J}(u + d(t)) \end{align*} x˙1?x˙2??=x2?=J1?(u+d(t))?
其中： x 1 = θ , x 2 = θ ˙ x_1=\theta ,x_2=\dot{\theta} x1?=θ,x2?=θ˙
1 首先，初始化Initialization函數
采用S函數來描述動力學方程，可選取1輸人2輸出系統，如果角度和角速度的初始值取零，則模型初始化參數寫為[0，0]，模型初始化S函數描述如下：（見模板)
2 微分方程描述的mdlDerivative函數
該函數可用于描述微分方程并實現數值求解。在控制系統中，可以采樣該函數來描述被控對象和自適應律等，并通過Simulink環境下選擇數值分析方法來實現對模型的數值求解
取 J = 2 ， d ( t ) = s i n t J=2，d(t)=sint J=2，d(t)=sint，則采用S函數可以實現模型角度和角速度的求解，描述如下：
function sys=mdlDerivatives(t,x,u)J=2;
dt=sin(t);
ut=u(1);
sys(1)=x(2);
sys(2)=1/J*(ut+dt);sys = [dx1;dx2];
3 用于輸出的mdlOutput函數
S函數的mdlOutput函數通常用于描述控制器或模型的輸出。采用S函數的mdlOutput模塊來描述模型角度和角速度的輸出：
function sys=mdlOutputs(t,x,u)sys(1) = x(1);
sys(2) = x(2);
最后，給出S函數模板
function [sys,x0,str,ts,simStateCompliance] = plant(t,x,u,flag,pa)
%SFUNTMPL General MATLAB S-Function Template
% With MATLAB 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 MATLAB S-function syntax is:
% [SYS,X0,STR,TS,SIMSTATECOMPLIANCE] = 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.
% 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.
%
%
% 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 more 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].
%
% SIMSTATECOMPLIANCE = Specifices how to handle this block when saving and
% restoring the complete simulation state of the
% model. The allowed values are: 'DefaultSimState',
% 'HasNoSimState' or 'DisallowSimState'. If this value
% is not speficified, then the block's compliance with
% simState feature is set to 'UknownSimState'.% Copyright 1990-2010 The MathWorks, Inc.%
% The following outlines the general structure of an S-function.
%
switch flag,%%%%%%%%%%%%%%%%%%% Initialization %%%%%%%%%%%%%%%%%%%case 0,[sys,x0,str,ts,simStateCompliance]=mdlInitializeSizes;%%%%%%%%%%%%%%%% Derivatives %%%%%%%%%%%%%%%%case 1,sys=mdlDerivatives(t,x,u,pa);%%%%%%%%%%% 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 %%%%%%%%%%%%%%%%%%%%%otherwiseDAStudio.error('Simulink:blocks:unhandledFlag', num2str(flag));end% end sfuntmpl%
%=============================================================================
% mdlInitializeSizes
% Return the sizes, initial conditions, and sample times for the S-function.
%=============================================================================
%
function [sys,x0,str,ts,simStateCompliance]=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.
%
sizes = simsizes;sizes.NumContStates = 2;
sizes.NumDiscStates = 0;
sizes.NumOutputs = 2;
sizes.NumInputs = 1;
sizes.DirFeedthrough = 0;
sizes.NumSampleTimes = 1; % at least one sample time is neededsys = simsizes(sizes);%
% initialize the initial conditions
%
x0 = [0,0];%
% str is always an empty matrix
%
str = [];%
% initialize the array of sample times
%
ts = [0 0];% Specify the block simStateCompliance. The allowed values are:
% 'UnknownSimState', < The default setting; warn and assume DefaultSimState
% 'DefaultSimState', < Same sim state as a built-in block
% 'HasNoSimState', < No sim state
% 'DisallowSimState' < Error out when saving or restoring the model sim state
simStateCompliance = 'UnknownSimState';% end mdlInitializeSizes%
%=============================================================================
% mdlDerivatives
% Return the derivatives for the continuous states.
%=============================================================================
%
function sys=mdlDerivatives(t,x,u,pa)
k=pa.k;
m=pa.m;x1=x(1);
x2=x(2);dx1=x2;
dx2=-k/m*x1^3+u/m;sys = [dx1;dx2];% 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 = x;% 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
)
