Documentation Center

  • Trials
  • Product Updates

C MEX S-Function Examples

About S-Function Examples

All examples are based on the C MEX S-function templates sfuntmpl_basic.c and sfuntmpl_doc.c. Open sfuntmpl_doc.csfuntmpl_doc.c. for a detailed discussion of the S-function template.

Continuous States

The csfunc.ccsfunc.c example shows how to model a continuous system with states using a C MEX S-function. The following Simulink® model uses this S-function.

sfcndemo_csfuncsfcndemo_csfunc

In continuous state integration, the Simulink solvers integrate a set of continuous states using the following equations.

S-functions that contain continuous states implement a state-space equation. The mdlOutputs method contains the output portion and mdlDerivatives method contains the derivative portion of the state-space equation. To visualize how the integration works, see the flowchart in Simulink Engine Interaction with C S-Functions. The output equation corresponds to the mdlOutputs in the major time step. Next, the example enters the integration section of the flowchart. Here the Simulink engine performs a number of minor time steps during which it calls mdlOutputs and mdlDerivatives. Each of these pairs of calls is referred to as an integration stage. The integration returns with the continuous states updated and the simulation time moved forward. Time is moved forward as far as possible, providing that error tolerances in the state are met. The maximum time step is subject to constraints of discrete events such as the actual simulation stop time and the user-imposed limit.

The csfunc.c example specifies that the input port has direct feedthrough. This is because matrix D is initialized to a nonzero matrix. If D is set equal to a zero matrix in the state-space representation, the input signal is not used in mdlOutputs. In this case, the direct feedthrough can be set to 0, which indicates that csfunc.c does not require the input signal when executing mdlOutputs.

matlabroot/simulink/src/csfunc.c

The S-function begins with #define statements for the S-function name and level, and a #include statement for the simstruc.h header. After these statements, the S-function can include or define any other necessary headers, data, etc. The csfunc.c example defines the variable U as a pointer to the first input port's signal and initializes static variables for the state-space matrices.

/*  File    : csfunc.c
 *  Abstract:
 *
 *      Example C S-function for defining a continuous system.  
 *
 *      x' = Ax + Bu
 *      y  = Cx + Du
 *
 *      For more details about S-functions, see simulink/src/sfuntmpl_doc.c.
 * 
 *  Copyright 1990-2007 The MathWorks, Inc.
 */

#define S_FUNCTION_NAME csfunc
#define S_FUNCTION_LEVEL 2

#include "simstruc.h"

#define U(element) (*uPtrs[element])  /* Pointer to Input Port0 */

static real_T A[2][2]={ { -0.09, -0.01 } ,
                        {  1   ,  0    } 
                      };

static real_T B[2][2]={ {  1   , -7    } ,
                        {  0   , -2    } 
                      };

static real_T C[2][2]={ {  0   , 2     } ,
                        {  1   , -5    } 
                      };

static real_T D[2][2]={ { -3   , 0     } ,
                        {  1   , 0     } 
                      };

The required S-function method mdlInitializeSizes then sets up the following S-function characteristics.

  • ssSetNumSFcnParams sets the number of expected S-function dialog parameters to zero.

  • ssGetSFcnParamsCount determines how many parameters the user actually entered into the S-function dialog. If the number of user-specified parameters does not match the number returned by ssGetNumSFcnParams, the S-function errors out.

  • If the S-function parameter count passes, mdlInitializeSizes sets the number of continuous and discrete states using ssSetNumContStates and ssSetNumDiscStates, respectively. This example has two continuous states and zero discrete states.

  • Next, the method configures the S-function to have a single input and output port, each with a width of two to match the dimensions of the state-space matrices. The method passes a value of 1 to ssSetInputPortDirectFeedThrough to indicate the input port has direct feedthrough.

  • ssSetNumSampleTimes initializes one sample time, which the mdlInitializeSampleTimes function configures later.

  • The S-function indicates that no work vectors are used by passing a value of 0 to ssSetNumRWork, ssSetNumIWork, etc. You can omit these lines because zero is the default value for all of these macros. However, for clarity, the S-function explicitly sets the number of work vectors.

  • Lastly, ssSetOptions sets any applicable options. In this case, the only option is SS_OPTION_EXCEPTION_FREE_CODE, which stipulates that the code is exception free.

The mdlInitializeSizes function for this example is shown below.

/*====================*
 * S-function methods *
 *====================*/

/* Function: mdlInitializeSizes ===============================================
 * Abstract:
 *    Determine the S-function block's characteristics:
 *    number of inputs, outputs, states, etc.
 */
static void mdlInitializeSizes(SimStruct *S)
{
    ssSetNumSFcnParams(S, 0);  /* Number of expected parameters */
    if (ssGetNumSFcnParams(S) != ssGetSFcnParamsCount(S)) {
        return; /* Parameter mismatch reported by the Simulink engine*/
    }

    ssSetNumContStates(S, 2);
    ssSetNumDiscStates(S, 0);

    if (!ssSetNumInputPorts(S, 1)) return;
    ssSetInputPortWidth(S, 0, 2);
    ssSetInputPortDirectFeedThrough(S, 0, 1);

    if (!ssSetNumOutputPorts(S, 1)) return;
    ssSetOutputPortWidth(S, 0, 2);

    ssSetNumSampleTimes(S, 1);
    ssSetNumRWork(S, 0);
    ssSetNumIWork(S, 0);
    ssSetNumPWork(S, 0);
    ssSetNumModes(S, 0);
    ssSetNumNonsampledZCs(S, 0);

    /* Take care when specifying exception free code - see sfuntmpl_doc.c */
    ssSetOptions(S, SS_OPTION_EXCEPTION_FREE_CODE);
}

The required S-function method mdlInitializeSampleTimes specifies the S-function sample rates. The value CONTINUOUS_SAMPLE_TIME passed to the ssSetSampleTime macro specifies that the first S-function sample rate be continuous. ssSetOffsetTime then specifies an offset time of zero for this sample rate. The call to ssSetModelReferenceSampleTimeDefaultInheritance tells the solver to use the default rule to determine if referenced models containing this S-function can inherit their sample times from the parent model.

/* Function: mdlInitializeSampleTimes =========================================
 * Abstract:
 *    Specifiy that we have a continuous sample time.
 */
static void mdlInitializeSampleTimes(SimStruct *S)
{
    ssSetSampleTime(S, 0, CONTINUOUS_SAMPLE_TIME);
    ssSetOffsetTime(S, 0, 0.0);
    ssSetModelReferenceSampleTimeDefaultInheritance(S);
}

The optional S-function method mdlInitializeConditions initializes the continuous state vector. The #define statement before this method is required for the Simulink engine to call this function. In the example below, ssGetContStates obtains a pointer to the continuous state vector. The for loop then initializes each state to zero.

#define MDL_INITIALIZE_CONDITIONS
/* Function: mdlInitializeConditions ========================================
 * Abstract:
 *    Initialize both continuous states to zero.
 */
static void mdlInitializeConditions(SimStruct *S)
{
    real_T *x0 = ssGetContStates(S);
    int_T  lp;

    for (lp=0;lp<2;lp++) { 
        *x0++=0.0; 
    }
}

The required mdlOutputs function computes the output signal of this S-function. The beginning of the function obtains pointers to the first output port, continuous states, and first input port. The S-function uses the data in these arrays to solve the output equation y=Cx+Du.

/* Function: mdlOutputs =======================================================
 * Abstract:
 *      y = Cx + Du 
 */
static void mdlOutputs(SimStruct *S, int_T tid)
{
    real_T            *y    = ssGetOutputPortRealSignal(S,0);
    real_T            *x    = ssGetContStates(S);
    InputRealPtrsType uPtrs = ssGetInputPortRealSignalPtrs(S,0);
 
    UNUSED_ARG(tid); /* not used in single tasking mode */

    /* y=Cx+Du */
    y[0]=C[0][0]*x[0]+C[0][1]*x[1]+D[0][0]*U(0)+D[0][1]*U(1);
    y[1]=C[1][0]*x[0]+C[1][1]*x[1]+D[1][0]*U(0)+D[1][1]*U(1);
}

The mdlDerivatives function calculates the continuous state derivatives. Because this function is an optional method, a #define statement must precede the function. The beginning of the function obtains pointers to the S-function continuous states, state derivatives, and first input port. The S-function uses this data to solve the equation dx=Ax+Bu.

#define MDL_DERIVATIVES
/* Function: mdlDerivatives =================================================
 * Abstract:
 *      xdot = Ax + Bu
 */
static void mdlDerivatives(SimStruct *S)
{
    real_T            *dx   = ssGetdX(S);
    real_T            *x    = ssGetContStates(S);
    InputRealPtrsType uPtrs = ssGetInputPortRealSignalPtrs(S,0);

    /* xdot=Ax+Bu */
    dx[0]=A[0][0]*x[0]+A[0][1]*x[1]+B[0][0]*U(0)+B[0][1]*U(1);
    dx[1]=A[1][0]*x[0]+A[1][1]*x[1]+B[1][0]*U(0)+B[1][1]*U(1);
}

The required mdlTerminate function performs any actions, such as freeing memory, necessary at the end of the simulation. In this example, the function is empty.

/* Function: mdlTerminate =====================================================
 * Abstract:
 *    No termination needed, but we are required to have this routine.
 */
static void mdlTerminate(SimStruct *S)
{
    UNUSED_ARG(S); /* unused input argument */
}

The required S-function trailer includes the files necessary for simulation or code generation, as follows.

#ifdef  MATLAB_MEX_FILE    /* Is this file being compiled as a MEX file? */
#include "simulink.c"      /* MEX file interface mechanism */
#else
#include "cg_sfun.h"       /* Code generation registration function */
#endif

    Note   The mdlOutputs and mdlTerminate functions use the UNUSED_ARG macro to indicate that an input argument the callback requires is not used. This optional macro is defined in simstruc_types.hsimstruc_types.h. If used, you must call this macro once for each input argument that a callback does not use.

Discrete States

The dsfunc.cdsfunc.c example shows how to model a discrete system in a C MEX S-function. The following Simulink model uses this S-function.

sfcndemo_dsfuncsfcndemo_dsfunc

Discrete systems can be modeled by the following set of equations.

The dsfunc.c example implements a discrete state-space equation. The mdlOutputs method contains the output portion and the mdlUpdate method contains the update portion of the discrete state-space equation. To visualize how the simulation works, see the flowchart in Simulink Engine Interaction with C S-Functions. The output equation above corresponds to the mdlOutputs in the major time step. The preceding update equation corresponds to the mdlUpdate in the major time step. If your model does not contain continuous elements, the Simulink engine skips the integration phase and time is moved forward to the next discrete sample hit.

matlabroot/simulink/src/dsfunc.c

The S-function begins with #define statements for the S-function name and level, along with a #include statement for the simstruc.h header. After these statements, the S-function can include or define any other necessary headers, data, etc. The dsfunc.c example defines U as a pointer to the first input port's signal and initializes static variables for the state-space matrices.

/*  File    : dsfunc.c
 *  Abstract:
 *
 *      Example C S-function for defining a discrete system.  
 *
 *      x(n+1) = Ax(n) + Bu(n)
 *      y(n)   = Cx(n) + Du(n)
 *
 *      For more details about S-functions, see simulink/src/sfuntmpl_doc.c.
 * 
 *  Copyright 1990-2007 The MathWorks, Inc.
 */

#define S_FUNCTION_NAME dsfunc
#define S_FUNCTION_LEVEL 2

#include "simstruc.h"

#define U(element) (*uPtrs[element])  /* Pointer to Input Port0 */

static real_T A[2][2]={ { -1.3839, -0.5097 } ,
                        {  1     ,  0      }
                      };
 
static real_T B[2][2]={ { -2.5559,  0      } ,
                        {  0     ,  4.2382 }
                      };
 
static real_T C[2][2]={ {  0     ,  2.0761 } ,
                        {  0     ,  7.7891 }
                      };
 
static real_T D[2][2]={ { -0.8141, -2.9334 } ,
                        {  1.2426,  0      }
                      };

The required S-function method mdlInitializeSizes then sets up the following S-function characteristics.

  • ssSetNumSFcnParams sets the number of expected S-function dialog parameters to zero.

  • ssGetSFcnParamsCount determines how many parameters the user actually entered into the S-function dialog. If the number of user-specified parameters does not match the number returned by ssGetNumSFcnParams, the S-function errors out.

  • If the S-function parameter count passes, mdlInitializeSizes next sets the number of continuous and discrete states using ssSetNumContStates and ssSetNumDiscStates, respectively. This example has zero continuous states and two discrete states.

  • Next, the method configures the S-function to have a single input and output port, each with a width of two to match the dimensions of the state-space matrices. The method passes a value of 1 to ssSetInputPortDirectFeedThrough to indicate the input port has direct feedthrough.

  • ssSetNumSampleTimes initializes one sample time, which the mdlInitializeSampleTimes function configures later.

  • The S-function indicates that no work vectors are used by passing a value of 0 to ssSetNumRWork, ssSetNumIWork, etc. You can omit these lines because zero is the default value for all of these macros. However, for clarity, the S-function explicitly sets the number of work vectors.

  • Lastly, ssSetOptions sets any applicable options. In this case, the only option is SS_OPTION_EXCEPTION_FREE_CODE, which stipulates that the code is exception free.

The mdlInitializeSizes function for this example is shown below.

/*====================*
 * S-function methods *
 *====================*/

/* Function: mdlInitializeSizes ===============================================
 * Abstract:
 *    Determine the S-function block's characteristics:
 *    number of inputs, outputs, states, etc.
 */
static void mdlInitializeSizes(SimStruct *S)
{
    ssSetNumSFcnParams(S, 0);  /* Number of expected parameters */
    if (ssGetNumSFcnParams(S) != ssGetSFcnParamsCount(S)) {
        return; /* Parameter mismatch reported by the Simulink engine*/
    }

    ssSetNumContStates(S, 0);
    ssSetNumDiscStates(S, 2);

    if (!ssSetNumInputPorts(S, 1)) return;
    ssSetInputPortWidth(S, 0, 2);
    ssSetInputPortDirectFeedThrough(S, 0, 1);

    if (!ssSetNumOutputPorts(S, 1)) return;
    ssSetOutputPortWidth(S, 0, 2);

    ssSetNumSampleTimes(S, 1);
    ssSetNumRWork(S, 0);
    ssSetNumIWork(S, 0);
    ssSetNumPWork(S, 0);
    ssSetNumModes(S, 0);
    ssSetNumNonsampledZCs(S, 0);

    /* Take care when specifying exception free code - see sfuntmpl_doc.c */
    ssSetOptions(S, SS_OPTION_EXCEPTION_FREE_CODE);
}

The required S-function method mdlInitializeSampleTimes specifies the S-function sample rates. A call to ssSetSampleTime sets this first S-function sample period to 1.0. ssSetOffsetTime then specifies an offset time of zero for the first sample rate. The call to ssSetModelReferenceSampleTimeDefaultInheritance tells the solver to use the default rule to determine if referenced models containing this S-function can inherit their sample times from the parent model.

/* Function: mdlInitializeSampleTimes =========================================
 * Abstract:
 *    Specifiy a sample time 0f 1.0.
 */
static void mdlInitializeSampleTimes(SimStruct *S)
{
    ssSetSampleTime(S, 0, 1.0);
    ssSetOffsetTime(S, 0, 0.0);
    ssSetModelReferenceSampleTimeDefaultInheritance(S);
}

The optional S-function method mdlInitializeConditions initializes the discrete state vector. The #define statement before this method is required for the Simulink engine to call this function. In the example below, ssGetRealDiscStates obtains a pointer to the discrete state vector. The for loop then initializes each discrete state to one.

#define MDL_INITIALIZE_CONDITIONS
/* Function: mdlInitializeConditions ========================================
 * Abstract:
 *    Initialize both discrete states to one.
 */
static void mdlInitializeConditions(SimStruct *S)
{
    real_T *x0 = ssGetRealDiscStates(S);
    int_T  lp;

    for (lp=0;lp<2;lp++) { 
        *x0++=1.0; 
    }
}

The required mdlOutputs function computes the output signal of this S-function. The beginning of the function obtains pointers to the first output port, discrete states, and first input port. The S-function uses the data in these arrays to solve the output equation y=Cx+Du.

/* Function: mdlOutputs =======================================================
 * Abstract:
 *      y = Cx + Du 
 */
static void mdlOutputs(SimStruct *S, int_T tid)
{
    real_T            *y    = ssGetOutputPortRealSignal(S,0);
    real_T            *x    = ssGetRealDiscStates(S);
    InputRealPtrsType uPtrs = ssGetInputPortRealSignalPtrs(S,0);
 
    UNUSED_ARG(tid); /* not used in single tasking mode */

    /* y=Cx+Du */
    y[0]=C[0][0]*x[0]+C[0][1]*x[1]+D[0][0]*U(0)+D[0][1]*U(1);
    y[1]=C[1][0]*x[0]+C[1][1]*x[1]+D[1][0]*U(0)+D[1][1]*U(1);
}

The Simulink engine calls the mdlUpdate function once every major integration time step to update the discrete states' values. Because this function is an optional method, a #define statement must precede the function. The beginning of the function obtains pointers to the S-function discrete states and first input port. The S-function uses the data in these arrays to solve the equation dx=Ax+Bu, which is stored in the temporary variable tempX before being assigned into the discrete state vector x.

#define MDL_UPDATE
/* Function: mdlUpdate ======================================================
 * Abstract:
 *      xdot = Ax + Bu
 */
static void mdlUpdate(SimStruct *S, int_T tid)
{
    real_T            tempX[2] = {0.0, 0.0};
    real_T            *x       = ssGetRealDiscStates(S);
    InputRealPtrsType uPtrs    = ssGetInputPortRealSignalPtrs(S,0);

    UNUSED_ARG(tid); /* not used in single tasking mode */

    /* xdot=Ax+Bu */
    tempX[0]=A[0][0]*x[0]+A[0][1]*x[1]+B[0][0]*U(0)+B[0][1]*U(1);
    tempX[1]=A[1][0]*x[0]+A[1][1]*x[1]+B[1][0]*U(0)+B[1][1]*U(1);
 
    x[0]=tempX[0];
    x[1]=tempX[1];
}

The required mdlTerminate function performs any actions, such as freeing memory, necessary at the end of the simulation. In this example, the function is empty.

/* Function: mdlTerminate =====================================================
 * Abstract:
 *    No termination needed, but we are required to have this routine.
 */
static void mdlTerminate(SimStruct *S)
{
    UNUSED_ARG(S); /* unused input argument */
}

The required S-function trailer includes the files necessary for simulation or code generation, as follows.

#ifdef  MATLAB_MEX_FILE    /* Is this file being compiled as a MEX file? */
#include "simulink.c"      /* MEX file interface mechanism */
#else
#include "cg_sfun.h"       /* Code generation registration function */
#endif

    Note   The mdlOutputs and mdlTerminate functions use the UNUSED_ARG macro to indicate that an input argument the callback requires is not used. This optional macro is defined in simstruc_types.hsimstruc_types.h. If used, you must call this macro once for each input argument that a callback does not use.

Continuous and Discrete States

The mixedm.cmixedm.c example shows a hybrid (a combination of continuous and discrete states) system. The mixedm.c example combines elements of csfunc.c and dsfunc.c. The following Simulink model uses this S-function.

sfcndemo_mixedmsfcndemo_mixedm

If you have a hybrid system, the mdlDerivatives method calculates the derivatives of the continuous states of the state vector, x, and the mdlUpdate method contains the equations used to update the discrete state vector, xD. The mdlOutputs method computes the S-function outputs after checking for sample hits to determine at what point the S-function is being called.

In Simulink block diagram form, the S-function mixedm.c looks like

which implements a continuous integrator followed by a discrete unit delay.

matlabroot/simulink/src/mixedm.c

The S-function begins with #define statements for the S-function name and level, along with a #include statement for the simstruc.h header. After these statements, the S-function can include or define any other necessary headers, data, etc. The mixedm.c example defines U as a pointer to the first input port's signal.

/*  File    : mixedm.c
 *  Abstract:
 *
 *      An example S-function illustrating multiple sample times by implementing
 *         integrator -> ZOH(Ts=1second) -> UnitDelay(Ts=1second) 
 *      with an initial condition of 1.
 *	(e.g. an integrator followed by unit delay operation).
 *
 *      For more details about S-functions, see simulink/src/sfuntmpl_doc.c
 *
 *  Copyright 1990-2007 The MathWorks, Inc.
 */

#define S_FUNCTION_NAME mixedm
#define S_FUNCTION_LEVEL 2

#include "simstruc.h"

#define U(element) (*uPtrs[element])  /* Pointer to Input Port0 */

The required S-function method mdlInitializeSizes then sets up the following S-function characteristics.

  • ssSetNumSFcnParams sets the number of expected S-function dialog parameters to zero.

  • ssGetSFcnParamsCount determines how many parameters the user actually entered into the S-function dialog. If the number of user-specified parameters does not match the number returned by ssGetNumSFcnParams, the S-function errors out.

  • If the S-function parameter count passes, mdlInitializeSizes next sets the number of continuous and discrete states using ssSetNumContStates and ssSetNumDiscStates, respectively. This example has one continuous state and one discrete state.

  • The S-function initializes one floating-point work vector by passing a value of 1 to ssSetNumRWork. No other work vectors are initialized.

  • Next, the method uses ssSetNumInputPorts and ssSetNumOutputPorts to configure the S-function to have a single input and output port, each with a width of one. The method passes a value of 1 to ssSetInputPortDirectFeedThrough to indicate the input port has direct feedthrough.

  • This S-function assigns sample times using a hybrid block-based and port-based method. The macro ssSetNumSampleTimes initializes two block-based sample times, which the mdlInitializeSampleTimes function configures later. The macros ssSetInputPortSampleTime and ssSetInputPortOffsetTime initialize the input port to have a continuous sample time with an offset of zero. Similarly, ssSetOutputPortSampleTime and ssSetOutputPortOffsetTime initialize the output port sample time to 1 with an offset of zero.

  • Lastly, ssSetOptions sets two S-function options. SS_OPTION_EXCEPTION_FREE_CODE stipulates that the code is exception free and SS_OPTION_PORT_SAMPLE_TIMES_ASSIGNED indicates a combination of block-based and port-based sample times.

The mdlInitializeSizes function for this example is shown below.

*====================*
 * S-function methods *
 *====================*/

/* Function: mdlInitializeSizes ===============================================
 * Abstract:
 *    Determine the S-function block's characteristics:
 *    number of inputs, outputs, states, etc.
 */
static void mdlInitializeSizes(SimStruct *S)
{
    ssSetNumSFcnParams(S, 0);  /* Number of expected parameters */
    if (ssGetNumSFcnParams(S) != ssGetSFcnParamsCount(S)) {
        return; /* Parameter mismatch reported by the Simulink engine*/
    }

    ssSetNumContStates(S, 1);
    ssSetNumDiscStates(S, 1);
    ssSetNumRWork(S, 1);  /* for zoh output feeding the delay operator */

    if (!ssSetNumInputPorts(S, 1)) return;
    ssSetInputPortWidth(S, 0, 1);
    ssSetInputPortDirectFeedThrough(S, 0, 1);
    ssSetInputPortSampleTime(S, 0, CONTINUOUS_SAMPLE_TIME);
    ssSetInputPortOffsetTime(S, 0, 0.0);

    if (!ssSetNumOutputPorts(S, 1)) return;
    ssSetOutputPortWidth(S, 0, 1);
    ssSetOutputPortSampleTime(S, 0, 1.0);
    ssSetOutputPortOffsetTime(S, 0, 0.0);

    ssSetNumSampleTimes(S, 2);

    /* Take care when specifying exception free code - see sfuntmpl_doc.c. */
    ssSetOptions(S, (SS_OPTION_EXCEPTION_FREE_CODE |
                     SS_OPTION_PORT_SAMPLE_TIMES_ASSIGNED));

} /* end mdlInitializeSizes */

The required S-function method mdlInitializeSampleTimes specifies the S-function block-based sample rates. The first call to ssSetSampleTime specifies that the first sample rate is continuous, with the subsequent call to ssSetOffsetTime setting the offset to zero. The second call to this pair of macros sets the second sample time to 1 with an offset of zero. The S-function port-based sample times set in mdlInitializeSizes must all be registered as a block-based sample time. The call to ssSetModelReferenceSampleTimeDefaultInheritance tells the solver to use the default rule to determine if referenced models containing this S-function can inherit their sample times from the parent model.

/* Function: mdlInitializeSampleTimes =========================================
 * Abstract:
 *    Two tasks: One continuous, one with discrete sample time of 1.0.
 */
static void mdlInitializeSampleTimes(SimStruct *S)
{
    ssSetSampleTime(S, 0, CONTINUOUS_SAMPLE_TIME);
    ssSetOffsetTime(S, 0, 0.0);

    ssSetSampleTime(S, 1, 1.0);
    ssSetOffsetTime(S, 1, 0.0);
    ssSetModelReferenceSampleTimeDefaultInheritance(S);
} /* end mdlInitializeSampleTimes */

The optional S-function method mdlInitializeConditions initializes the continuous and discrete state vectors. The #define statement before this method is required for the Simulink engine to call this function. In this example, ssGetContStates obtains a pointer to the continuous state vector and ssGetRealDiscStates obtains a pointer to the discrete state vector. The method then sets all states' initial conditions to one.

#define MDL_INITIALIZE_CONDITIONS
/* Function: mdlInitializeConditions ==========================================
 * Abstract:
 *    Initialize both continuous states to one.
 */
static void mdlInitializeConditions(SimStruct *S)
{
    real_T *xC0 = ssGetContStates(S);
    real_T *xD0 = ssGetRealDiscStates(S);

    xC0[0] = 1.0;
    xD0[0] = 1.0;

} /* end mdlInitializeConditions */

The required mdlOutputs function performs computations based on the current task. The macro ssIsContinuousTask checks if the continuous task is executing. If this macro returns true, ssIsSpecialSampleHit then checks if the discrete sample rate is also executing. If this macro also returns true, the method sets the value of the floating-point work vector to the current value of the continuous state, via pointers obtained using ssGetRWork and ssGetContStates, respectively. The mdlUpdate method later uses the floating-point work vector as the input to the zero-order hold. Updating the work vector in mdlOutputs ensures that the correct values are available during subsequent calls to mdlUpdate. Finally, if the S-function is running at its discrete rate, i.e., the call to ssIsSampleHit returns true, the method sets the output to the value of the discrete state.

/* Function: mdlOutputs =======================================================
 * Abstract:
 *      y = xD, and update the zoh internal output.
 */
static void mdlOutputs(SimStruct *S, int_T tid)
{
    /* update the internal "zoh" output */
    if (ssIsContinuousTask(S, tid)) {
        if (ssIsSpecialSampleHit(S, 1, 0, tid)) {
            real_T *zoh = ssGetRWork(S);
            real_T *xC  = ssGetContStates(S);
            *zoh = *xC;
        }
    }

    /* y=xD */
    if (ssIsSampleHit(S, 1, tid)) {
        real_T *y   = ssGetOutputPortRealSignal(S,0);
        real_T *xD  = ssGetRealDiscStates(S);
        y[0]=xD[0];
    }


} /* end mdlOutputs */

The Simulink engine calls the mdlUpdate function once every major integration time step to update the discrete states' values. Because this function is an optional method, a #define statement must precede the function. The call to ssIsSampleHit ensures the body of the method is executed only when the S-function is operating at its discrete rate. If ssIsSampleHit returns true, the method obtains pointers to the S-function discrete state and floating-point work vector and updates the discrete state's value using the value stored in the work vector.

#define MDL_UPDATE
/* Function: mdlUpdate ======================================================
 * Abstract:
 *      xD = xC
 */
static void mdlUpdate(SimStruct *S, int_T tid)
{
    UNUSED_ARG(tid); /* not used in single tasking mode */

    /* xD=xC */
    if (ssIsSampleHit(S, 1, tid)) {
        real_T *xD = ssGetRealDiscStates(S);
        real_T *zoh = ssGetRWork(S);
        xD[0]=*zoh;
    }
} /* end mdlUpdate */

The mdlDerivatives function calculates the continuous state derivatives. Because this function is an optional method, a #define statement must precede the function. The function obtains pointers to the S-function continuous state derivative and first input port then sets the continuous state derivative equal to the value of the first input.

#define MDL_DERIVATIVES
/* Function: mdlDerivatives =================================================
 * Abstract:
 *      xdot = U
 */
static void mdlDerivatives(SimStruct *S)
{
    real_T            *dx   = ssGetdX(S);
    InputRealPtrsType uPtrs = ssGetInputPortRealSignalPtrs(S,0);

    /* xdot=U */
    dx[0]=U(0);

} /* end mdlDerivatives */

The required mdlTerminate function performs any actions, such as freeing memory, necessary at the end of the simulation. In this example, the function is empty.

/* Function: mdlTerminate =====================================================
 * Abstract:
 *    No termination needed, but we are required to have this routine.
 */
static void mdlTerminate(SimStruct *S)
{
    UNUSED_ARG(S); /* unused input argument */
}

The S-function trailer includes the files necessary for simulation or code generation, as follows.

#ifdef  MATLAB_MEX_FILE    /* Is this file being compiled as a MEX file? */
#include "simulink.c"      /* MEX file interface mechanism */
#else
#include "cg_sfun.h"       /* Code generation registration function */
#endif

    Note   The mdlUpdate and mdlTerminate functions use the UNUSED_ARG macro to indicate that an input argument the callback requires is not used. This optional macro is defined in simstruc_types.hsimstruc_types.h. If used, you must call this macro once for each input argument that a callback does not use.

Variable Sample Time

The example S-function vsfunc.cvsfunc.c uses a variable-step sample time. The following Simulink model uses this S-function.

sfcndemo_vsfuncsfcndemo_vsfunc

Variable step-size functions require a call to mdlGetTimeOfNextVarHit, which is an S-function routine that calculates the time of the next sample hit. S-functions that use the variable-step sample time can be used only with variable-step solvers. The vsfunc.c example is a discrete S-function that delays its first input by an amount of time determined by the second input.

The vsfunc.c example outputs the input u delayed by a variable amount of time. mdlOutputs sets the output y equal to state x. mdlUpdate sets the state vector x equal to u, the input vector. This example calls mdlGetTimeOfNextVarHit to calculate and set the time of the next sample hit, that is, the time when vsfunc.c is next called. In mdlGetTimeOfNextVarHit, the macro ssGetInputPortRealSignalPtrs gets a pointer to the input u. Then this call is made:

ssSetTNext(S, ssGetT(S) + U(1));

The macro ssGetT gets the simulation time t. The second input to the block, U(1), is added to t, and the macro ssSetTNext sets the time of the next hit equal to t+U(1), delaying the output by the amount of time set in (U(1)).

matlabroot/simulink/src/vsfunc.c

The S-function begins with #define statements for the S-function name and level, along with a #include statement for the simstruc.h header. After these statements, the S-function can include or define any other necessary headers, data, etc. The vsfunc.c example defines U as a pointer to the first input port's signal.

/*  File    : vsfunc.c
 *  Abstract:
 *
 *      Variable step S-function example.
 *      This example S-function illustrates how to create a variable step
 *      block.  This block implements a variable step delay
 *      in which the first input is delayed by an amount of time determined
 *      by the second input:
 *
 *      dt      = u(2)
 *      y(t+dt) = u(t)
 *
 *      For more details about S-functions, see simulink/src/sfuntmpl_doc.c.
 *
 *  Copyright 1990-2007 The MathWorks, Inc.
 */

#define S_FUNCTION_NAME vsfunc
#define S_FUNCTION_LEVEL 2

#include "simstruc.h"

#define U(element) (*uPtrs[element])  /* Pointer to Input Port0 */

The required S-function method mdlInitializeSizes then sets up the following S-function characteristics.

  • ssSetNumSFcnParams sets the number of expected S-function dialog parameters to zero.

  • ssGetSFcnParamsCount determines how many parameters the user actually entered into the S-function dialog. If the number of user-specified parameters does not match the number returned by ssGetNumSFcnParams, the S-function errors out.

  • If the S-function parameter count passes, mdlInitializeSizes next sets the number of continuous and discrete states using ssSetNumContStates and ssSetNumDiscStates, respectively. This example has no continuous states and one discrete state.

  • Next, the method uses ssSetNumInputPorts and ssSetNumOutputPorts to configure the S-function to have a single input and output port. Calls to ssSetInputPortWidth and ssSetOutputPortWidth assign widths to these input and output ports. The method passes a value of 1 to ssSetInputPortDirectFeedThrough to indicate the input port has direct feedthrough.

  • ssSetNumSampleTimes then initializes one sample time, which the mdlInitializeSampleTimes function configures later.

  • The S-function indicates that no work vectors are used by passing a value of 0 to ssSetNumRWork, ssSetNumIWork, etc. You can omit these lines because zero is the default value for all of these macros. However, for clarity, the S-function explicitly sets the number of work vectors.

  • Next, ssGetSimMode checks if the S-function is being run in a simulation or by the Simulink Coder™ product. If ssGetSimMode returns SS_SIMMODE_RTWGEN and ssIsVariableStepSolver returns false, indicating use with the Simulink Coder product and a fixed-step solver, then the S-function errors out.

  • Lastly, ssSetOptions sets any applicable options. In this case, the only option is SS_OPTION_EXCEPTION_FREE_CODE, which stipulates that the code is exception free.

The mdlInitializeSizes function for this example is shown below.

/* Function: mdlInitializeSizes ===============================================
 * Abstract:
 *    Determine the S-function block's characteristics:
 *    number of inputs, outputs, states, etc.
 */
static void mdlInitializeSizes(SimStruct *S)
{
    ssSetNumSFcnParams(S, 0);  /* Number of expected parameters */
    if (ssGetNumSFcnParams(S) != ssGetSFcnParamsCount(S)) {
        return; /* Parameter mismatch reported by the Simulink engine*/
    }

    ssSetNumContStates(S, 0);
    ssSetNumDiscStates(S, 1);

    if (!ssSetNumInputPorts(S, 1)) return;
    ssSetInputPortWidth(S, 0, 2);
    ssSetInputPortDirectFeedThrough(S, 0, 1);

    if (!ssSetNumOutputPorts(S, 1)) return;
    ssSetOutputPortWidth(S, 0, 1);

    ssSetNumSampleTimes(S, 1);
    ssSetNumRWork(S, 0);
    ssSetNumIWork(S, 0);
    ssSetNumPWork(S, 0);
    ssSetNumModes(S, 0);
    ssSetNumNonsampledZCs(S, 0);

    if (ssGetSimMode(S) == SS_SIMMODE_RTWGEN && !ssIsVariableStepSolver(S)) {
        ssSetErrorStatus(S, "S-function vsfunc.c cannot be used with RTW "
                         "and Fixed-Step Solvers because it contains variable"
                         " sample time");
    }

    /* Take care when specifying exception free code - see sfuntmpl_doc.c */
    ssSetOptions(S, SS_OPTION_EXCEPTION_FREE_CODE);
}

The required S-function method mdlInitializeSampleTimes specifies the S-function sample rates. The input argument VARIABLE_SAMPLE_TIME passed to ssSetSampleTime specifies that this S-function has a variable-step sample time and ssSetOffsetTime specifies an offset time of zero. The call to ssSetModelReferenceSampleTimeDefaultInheritance tells the solver to use the default rule to determine if referenced models containing this S-function can inherit their sample times from the parent model. Because the S-function has a variable-step sample time, vsfunc.c must calculate the time of the next sample hit in the mdlGetTimeOfNextVarHit method, shown later.

/* Function: mdlInitializeSampleTimes =========================================
 * Abstract:
 *    Variable-Step S-function
 */
static void mdlInitializeSampleTimes(SimStruct *S)
{
    ssSetSampleTime(S, 0, VARIABLE_SAMPLE_TIME);
    ssSetOffsetTime(S, 0, 0.0);
    ssSetModelReferenceSampleTimeDefaultInheritance(S);
}

The optional S-function method mdlInitializeConditions initializes the discrete state vector. The #define statement before this method is required for the Simulink engine to call this function. In the example, the method uses ssGetRealDiscStates to obtain a pointer to the discrete state vector and sets the state's initial value to zero.

#define MDL_INITIALIZE_CONDITIONS
/* Function: mdlInitializeConditions ========================================
 * Abstract:
 *    Initialize discrete state to zero.
 */
static void mdlInitializeConditions(SimStruct *S)
{
    real_T *x0 = ssGetRealDiscStates(S);

    x0[0] = 0.0;
}

The optional mdlGetTimeOfNextVarHit method calculates the time of the next sample hit. Because this method is optional, a #define statement precedes it. First, this method obtains a pointer to the first input port's signal using ssGetInputPortRealSignalPtrs. If the input signal's second element is positive, the macro ssGetT gets the simulation time t. The macro ssSetTNext sets the time of the next hit equal to t+(*U[1]), delaying the output by the amount of time specified by the input's second element (*U[1]).

#define MDL_GET_TIME_OF_NEXT_VAR_HIT
static void mdlGetTimeOfNextVarHit(SimStruct *S)
{
    InputRealPtrsType uPtrs = ssGetInputPortRealSignalPtrs(S,0);
    
    /* Make sure input will increase time */
    if (U(1) <= 0.0) {
        /* If not, abort simulation */
        ssSetErrorStatus(S,"Variable step control input must be "
                         "greater than zero");
        return;
    }
    ssSetTNext(S, ssGetT(S)+U(1));
}

The required mdlOutputs function computes the S-function output signal. The function obtains pointers to the first output port and discrete state and then assigns the state's current value to the output.

/* Function: mdlOutputs =======================================================
 * Abstract:
 *      y = x
 */
static void mdlOutputs(SimStruct *S, int_T tid)
{
    real_T *y = ssGetOutputPortRealSignal(S,0);
    real_T *x = ssGetRealDiscStates(S);

    /* Return the current state as the output */
    y[0] = x[0];
}

The mdlUpdate function updates the discrete state's value. Because this method is optional, a #define statement precedes it. The function first obtains pointers to the S-function discrete state and first input port then assigns the value of the first element of the first input port signal to the state.

#define MDL_UPDATE
/* Function: mdlUpdate ========================================================
 * Abstract:
 *    This function is called once for every major integration time step.
 *    Discrete states are typically updated here, but this function is useful
 *    for performing any tasks that should only take place once per integration
 *    step.
 */
static void mdlUpdate(SimStruct *S, int_T tid)
{
    real_T            *x    = ssGetRealDiscStates(S);
    InputRealPtrsType uPtrs = ssGetInputPortRealSignalPtrs(S,0);

    x[0]=U(0);
}

The required mdlTerminate function performs any actions, such as freeing memory, necessary at the end of the simulation. In this example, the function is empty.

/* Function: mdlTerminate =====================================================
 * Abstract:
 *    No termination needed, but we are required to have this routine.
 */
static void mdlTerminate(SimStruct *S)
{
}

The required S-function trailer includes the files necessary for simulation or code generation, as follows.

#ifdef  MATLAB_MEX_FILE    /* Is this file being compiled as a MEX file? */
#include "simulink.c"      /* MEX file interface mechanism */
#else
#include "cg_sfun.h"       /* Code generation registration function */
#endif

Array Inputs and Outputs

The example S-function sfun_matadd.csfun_matadd.c demonstrates how to implement a matrix addition block. The following Simulink model uses this S-function.

sfcndemo_mataddsfcndemo_matadd

The S-function adds signals of various dimensions to a parameter value entered in the S-function. The S-function accepts and outputs 2-D or n-D signals.

matlabroot/simulink/src/sfun_matadd.c

The S-function begins with #define statements for the S-function name and level, along with a #include statement for the simstruc.h header. After these statements, the S-function includes or defines any other necessary headers, data, etc. This example defines additional variables for the number of S-function parameters, the S-function parameter value, and the flag EDIT_OK that indicates if the parameter value can be edited during simulation.

/* SFUN_MATADD matrix support example. 
 *   C MEX S-function for matrix addition with one input port,
 *   one output port, and one parameter.
 *
 *  Input Signal:  2-D or n-D array
 *  Parameter:     2-D or n-D array
 *  Output Signal: 2-D or n-D array
 *
 *  Input   parameter    output
 *  --------------------------------
 *  scalar   scalar      scalar
 *  scalar   matrix     matrix     (input scalar expansion)
 *  matrix   scalar     matrix     (parameter scalar expansion)
 *  matrix   matrix     matrix
 *
 *  Copyright 1990-2007 The MathWorks, Inc.
 */
#define S_FUNCTION_NAME  sfun_matadd
#define S_FUNCTION_LEVEL 2

#include "simstruc.h"

enum {PARAM = 0, NUM_PARAMS};

#define PARAM_ARG ssGetSFcnParam(S, PARAM)

#define EDIT_OK(S, ARG) \
  (!((ssGetSimMode(S) == SS_SIMMODE_SIZES_CALL_ONLY) \
  && mxIsEmpty(ARG)))

The S-function next implements the mdlCheckParameters method to validate the S-function dialog parameters. The #ifdef statement checks that the S-function is compiled as a MEX file, instead of for use with the Simulink Coder product. Because mdlCheckParameters is optional, the S-function code contains a #define statement to register the method. The body of the function checks that the S-function parameter value is not empty. If the parameter check fails, the S-function errors out with a call to ssSetErrorStatus.

#ifdef MATLAB_MEX_FILE 
#define MDL_CHECK_PARAMETERS 
/* Function: mdlCheckParameters ================================

 * Abstract:
 *    Verify parameter settings.
 */
static void mdlCheckParameters(SimStruct *S)
{
    if(EDIT_OK(S, PARAM_ARG)){
        /* Check that parameter value is not empty*/
        if( mxIsEmpty(PARAM_ARG) ) {
          ssSetErrorStatus(S, "Invalid parameter specified. The"
                              "parameter must be non-empty");
          return;
        }      
    }
} /* end mdlCheckParameters */
#endif 

The required S-function method mdlInitializeSizes then sets up the following S-function characteristics.

  • ssSetNumSFcnParams sets the number of expected S-function dialog parameters to one, as defined by the variable NUM_PARAMS.

  • If this S-function is compiled as a MEX file, ssGetSFcnParamsCount determines how many parameters the user actually entered into the S-function dialog. If the number of user-specified parameters matches the number returned by ssGetNumSFcnParams, the method calls mdlCheckParameters to validate the user-entered data. Otherwise, the S-function errors out.

  • If the parameter check passes, the S-function specifies that all S-function parameters are tunable using ssSetSFcnParamTunable.

  • The S-function then invokes ssAllowSignalsWithMoreThan2D to allow the S-function to accept n-D signals.

  • Next, ssSetNumOutputPorts and ssSetNumInputPorts specify that the S-function has a single output port and a single input port.

  • The S-function uses ssSetInputPortDimensionInfo to specify that the input port is dynamically sized. In this case, the S-function needs to implement an mdlSetInputPortDimensionInfo method to set the actual input dimension.

  • The output dimensions depend on the dimensions of the S-function parameter. If the parameter is a scalar, the call to ssSetOutputPortDimensionInfo specifies that the output port dimensions are dynamically sized. If the parameter is a matrix, the output port dimensions are initialized to the dimensions of the S-function parameter. In this case, the macro DECL_AND_INIT_DIMSINFO initializes a dimsInfo structure. The S-function assigns the width, size, and dimensions of the S-function parameter into the dimsInfo structure and then passes this structure to ssSetOutputPortDimensionInfo in order to set the output port dimensions accordingly.

  • The S-function specifies that the input port has direct feedthrough by passing a value of 1 to ssSetInputPortDirectFeedThrough.

  • ssSetNumSampleTimes initializes one sample time, to be configured later in the mdlInitializeSampleTimes method.

  • Lastly, ssSetOptions sets any applicable options. In this case, SS_OPTION_EXCEPTION_FREE_CODE stipulates that the code is exception free and SS_OPTION_WORKS_WITH_CODE_REUSE signifies that this S-function is compatible with the subsystem code reuse feature of the Simulink Coder product.

/* Function: mdlInitializeSizes ================================
 * Abstract:
 *   Initialize the sizes array
 */
static void mdlInitializeSizes(SimStruct *S)
{
    ssSetNumSFcnParams(S, NUM_PARAMS);

   #if defined(MATLAB_MEX_FILE)
      if (ssGetNumSFcnParams(S) != ssGetSFcnParamsCount(S)) {
         return; }
      mdlCheckParameters(S); 
      if (ssGetErrorStatus(S) != NULL) return; 
   #endif

   {
      int iParam = 0;
      int nParam = ssGetNumSFcnParams(S);

      for ( iParam = 0; iParam < nParam; iParam++ )
        {
          ssSetSFcnParamTunable( S, iParam, SS_PRM_TUNABLE );
        }
    }
 
   /* Allow signal dimensions greater than 2 */
   ssAllowSignalsWithMoreThan2D(S);
    
   /* Set number of input and output ports */
   if (!ssSetNumInputPorts( S,1)) return;
   if (!ssSetNumOutputPorts(S,1)) return;

   /* Set dimensions of input and output ports */
   {
      int_T pWidth = mxGetNumberOfElements(PARAM_ARG);
      /* Input can be a scalar or a matrix signal. */
      if(!ssSetInputPortDimensionInfo(S,0,DYNAMIC_DIMENSION)) {
          return; }

      if( pWidth == 1) { 
       /* Scalar parameter: output dimensions are unknown. */
       if(!ssSetOutputPortDimensionInfo(S,0,DYNAMIC_DIMENSION)){
           return; }
       }
      else{
         /* 
          * Non-scalar parameter: output dimensions are the same
          * as the parameter dimensions. To support n-D signals,
          * must use a dimsInfo structure to specify dimensions.
          */
          DECL_AND_INIT_DIMSINFO(di); /*Initializes structure*/
          int_T      pSize = mxGetNumberOfDimensions(PARAM_ARG);
          const int_T *pDims = mxGetDimensions(PARAM_ARG);
          di.width   = pWidth;
          di.numDims = pSize;
          di.dims    = pDims;
          if(!ssSetOutputPortDimensionInfo(S, 0, &di)) return;
          }
    }
    ssSetInputPortDirectFeedThrough(S, 0, 1);

    ssSetNumSampleTimes(S, 1);
    ssSetOptions(S,
                 SS_OPTION_WORKS_WITH_CODE_REUSE |
                 SS_OPTION_EXCEPTION_FREE_CODE);
} /* end mdlInitializeSizes */

The required S-function method mdlInitializeSampleTimes specifies the S-function sample rates. To specify that this S-function inherits its sample time from its driving block, the S-function calls ssSetSampleTime with the input argument INHERITED_SAMPLE_TIME. The call to ssSetModelReferenceSampleTimeDefaultInheritance tells the solver to use the default rule to determine if referenced models containing this S-function can inherit their sample times from the parent model.

/* Function: mdlInitializeSampleTimes ==========================
 * Abstract:
 *    Initialize the sample times array.
 */
static void mdlInitializeSampleTimes(SimStruct *S)
{
    ssSetSampleTime(S, 0, INHERITED_SAMPLE_TIME);
    ssSetOffsetTime(S, 0, 0.0);
    ssSetModelReferenceSampleTimeDefaultInheritance(S);
} /* end mdlInitializeSampleTimes */

The S-function calls the mdlSetWorkWidths method to register its run-time parameters. Because mdlSetWorkWidths is an optional method, a #define statement precedes it. The method first initializes a name for the run-time parameter and then uses ssRegAllTunableParamsAsRunTimeParams to register the run-time parameter.

/* Function: mdlSetWorkWidths ==================================
 * Abstract:
 *    Set up run-time parameter.
 */
#define MDL_SET_WORK_WIDTHS
static void mdlSetWorkWidths(SimStruct *S)
{
    const char_T    *rtParamNames[] = {"Operand"};
    ssRegAllTunableParamsAsRunTimeParams(S, rtParamNames);
} /* end mdlSetWorkWidths */

The S-function mdlOutputs method uses a for loop to calculate the output as the sum of the input and S-function parameter. The S-function handles n-D arrays of data using a single index into the array.

/* Function: mdlOutputs ========================================
 * Abstract:
 *   Compute the outputs of the S-function.
 */
static void mdlOutputs(SimStruct *S, int_T tid)
{
    InputRealPtrsType uPtr = ssGetInputPortRealSignalPtrs(S,0);
    real_T            *y   = ssGetOutputPortRealSignal(S,0);
    const real_T      *p   = mxGetPr(PARAM_ARG);
    
    int_T             uWidth = ssGetInputPortWidth(S,0);
    int_T             pWidth = mxGetNumberOfElements(PARAM_ARG);
    int_T             yWidth = ssGetOutputPortWidth(S,0);
    int               i;
    
    UNUSED_ARG(tid); /* not used in single tasking mode */
    
    /*
     * Note1: Matrix signals are stored in column major order.
     * Note2: Access each matrix element by one index not two
     *        indices. For example, if the output signal is a 
     *        [2x2] matrix signal,
     *        -          -
     *       | y[0]  y[2] |
     *       | y[1]  y[3] |
     *       -           -
     *       Output elements are stored as follows:
     *           y[0] --> row = 0, col = 0
     *           y[1] --> row = 1, col = 0
     *           y[2] --> row = 0, col = 1
     *           y[3] --> row = 1, col = 1
     */
    
    for (i = 0; i < yWidth; i++) {
        int_T uIdx = (uWidth == 1) ? 0 : i;
        int_T pIdx = (pWidth == 1) ? 0 : i;
        
        y[i] = *uPtr[uIdx] + p[pIdx];
    }
} /* end mdlOutputs */

During signal propagation, the S-function calls the optional mdlSetInputPortDimensionInfo method with the candidate input port dimensions stored in dimsInfo. The #if defined statement checks that the S-function is compiled as a MEX file. Because mdlSetInputPortDimensionInfo is an optional method, a #define statement precedes it. In mdlSetInputPortDimensionInfo, the S-function uses ssSetInputPortDimensionInfo to set the dimensions of the input port to the candidate dimensions. If the call to this macro succeeds, the S-function further checks the candidate dimensions to ensure that the input signal is either a 2-D scalar or a matrix. If this condition is met and the output port dimensions are still dynamically sized, the S-function calls ssSetOutputPortDimensionInfo to set the dimension of the output port to the same candidate dimensions. The ssSetOutputPortDimensionInfo macro cannot modify the output port dimensions if they are already specified.

#if defined(MATLAB_MEX_FILE)
#define MDL_SET_INPUT_PORT_DIMENSION_INFO
/* Function: mdlSetInputPortDimensionInfo ======================
 * Abstract:
 *    This routine is called with the candidate dimensions for
 *    an input port with unknown dimensions. If the proposed 
 *    dimensions are acceptable, the routine should go ahead and
 *    set the actual port dimensions. If they are unacceptable
 *    an error should be generated via ssSetErrorStatus.
 *    Note that any other input or output ports whose dimensions
 *    are implicitly defined by virtue of knowing the dimensions
 *    of the given port can also have their dimensions set.
 */
static void mdlSetInputPortDimensionInfo(SimStruct        *S,
   int_T            port,
   const DimsInfo_T *dimsInfo)
{
    int_T  pWidth          = mxGetNumberOfElements(PARAM_ARG);
    int_T  pSize           = mxGetNumberOfDimensions(PARAM_ARG);
    const int_T  *pDims    = mxGetDimensions(PARAM_ARG);
    
    int_T  uNumDims = dimsInfo->numDims;
    int_T  uWidth   = dimsInfo->width;
    int_T  *uDims   = dimsInfo->dims;
    
    int_T numDims;
    boolean_T  isOk = true;
    int iParam = 0;
    int_T outWidth = ssGetOutputPortWidth(S, 0);
    
    /* Set input port dimension */
    if(!ssSetInputPortDimensionInfo(S, port, dimsInfo)) return;
    
    /*
     * The block only accepts 2-D or higher signals. Check 
     * number of dimensions. If the parameter and the input
     * signal are non-scalar, their dimensions must be the same.
     */
    isOk = (uNumDims >= 2) && (pWidth == 1 || uWidth == 1 || 
       pWidth == uWidth);
    numDims = (pSize != uNumDims) ? numDims : uNumDims;
    
    if(isOk && pWidth > 1 && uWidth > 1){
        for ( iParam = 0; iParam < numDims; iParam++ ) {
            isOk = (pDims[iParam] == uDims[iParam]);
            if(!isOk) break;
        }
    }
    
    if(!isOk){
        ssSetErrorStatus(S,"Invalid input port dimensions. The "
        "input signal must be a 2-D scalar signal, or it must "
        "be a matrix with the same dimensions as the parameter "
        "dimensions.");
        return;
    }
    
    /* Set the output port dimensions */
    if (outWidth == DYNAMICALLY_SIZED){
      if(!ssSetOutputPortDimensionInfo(S,port,dimsInfo)) return;
    }
} /* end mdlSetInputPortDimensionInfo */

During signal propagation, if any output ports have unknown dimensions, the S-function calls the optional mdlSetOutputPortDimensionInfo method. Because this method is optional, a #define statement precedes it. In mdlSetOutputPortDimensionInfo, the S-function uses ssSetOutputPortDimensionInfo to set the dimensions of the output port to the candidate dimensions dimsInfo. If the call to this macro succeeds, the S-function further checks the candidate dimensions to ensure that the input signal is either a 2-D or n-D matrix. If this condition is not met, the S-function errors out with a call to ssSetErrorStatus. Otherwise, the S-function calls ssSetInputPortDimensionInfo to set the dimension of the input port to the same candidate dimensions.

# define MDL_SET_OUTPUT_PORT_DIMENSION_INFO
/* Function: mdlSetOutputPortDimensionInfo =====================
 * Abstract:
 *    This routine is called with the candidate dimensions for
 *    an output port with unknown dimensions. If the proposed
 *    dimensions are acceptable, the routine should go ahead and
 *    set the actual port dimensions. If they are unacceptable
 *    an error should be generated via ssSetErrorStatus.
 *    Note that any other input or output ports whose dimensions
 *    are implicitly defined by virtue of knowing the dimensions
 *    of the given port can also have their dimensions set.
 */
static void mdlSetOutputPortDimensionInfo(SimStruct        *S,
   int_T            port,
   const DimsInfo_T *dimsInfo)
{
    /*
     * If the block has scalar parameter, the output dimensions
     * are unknown. Set the input and output port to have the
     * same dimensions.
     */
    if(!ssSetOutputPortDimensionInfo(S, port, dimsInfo)) return;
    
    /* The block only accepts 2-D or n-D signals.
     * Check number of dimensions.
     */
    if (!(dimsInfo->numDims >= 2)){
        ssSetErrorStatus(S, "Invalid output port dimensions. "
        "The output signal must be a 2-D or n-D array (matrix) "
        "signal.");
        return;
    }else{
       /* Set the input port dimensions */
       if(!ssSetInputPortDimensionInfo(S,port,dimsInfo)) return;
    }
} /* end mdlSetOutputPortDimensionInfo */

Because the S-function has ports that are dynamically sized, it must provide an mdlSetDefaultPortDimensionInfo method. The Simulink engine invokes this method during signal propagation when it cannot determine the dimensionality of the signal connected to the block's input port. This situation can happen, for example, if the input port is unconnected. In this example, the mdlSetDefaultPortDimensionInfo method sets the input and output ports dimensions to a scalar.

# define MDL_SET_DEFAULT_PORT_DIMENSION_INFO
/* Function: mdlSetDefaultPortDimensionInfo ====================
 *    This routine is called when the Simulink engine is not able
 *    to find dimension candidates for ports with unknown dimensions.
 *    This function must set the dimensions of all ports with 
 *    unknown dimensions.
 */
static void mdlSetDefaultPortDimensionInfo(SimStruct *S)
{
    int_T outWidth = ssGetOutputPortWidth(S, 0);
    /* Input port dimension must be unknown. Set it to scalar.*/
    if(!ssSetInputPortMatrixDimensions(S, 0, 1, 1)) return;
    if(outWidth == DYNAMICALLY_SIZED){
        /* Output dimensions are unknown. Set it to scalar. */
        if(!ssSetOutputPortMatrixDimensions(S, 0, 1, 1)) return;
    }
} /* end mdlSetDefaultPortDimensionInfo */
#endif

The required mdlTerminate function performs any actions, such as freeing memory, necessary at the end of the simulation. In this example, the function is empty.

/* Function: mdlTerminate ======================================
 * Abstract:
 *    Called when the simulation is terminated.
 */
static void mdlTerminate(SimStruct *S)
{
    UNUSED_ARG(S); /* unused input argument */
    
} /* end mdlTerminate */

The required S-function trailer includes the files necessary for simulation or code generation.

#ifdef	MATLAB_MEX_FILE
  #include "simulink.c"
#else
  #include "cg_sfun.h"
#endif

/* [EOF] sfun_matadd.c */

    Note   The mdlOutputs and mdlTerminate functions use the UNUSED_ARG macro to indicate that an input argument the callback requires is not used. This optional macro is defined in simstruc_types.hsimstruc_types.h. You must call this macro once for each input argument that a callback does not use.

Zero-Crossing Detection

The example S-function sfun_zc_sat.csfun_zc_sat.c demonstrates how to implement a Saturation block. The following Simulink model uses this S-function.

sfcndemo_sfun_zc_satsfcndemo_sfun_zc_sat

The S-function works with either fixed-step or variable-step solvers. When this S-function inherits a continuous sample time and uses a variable-step solver, it uses a zero-crossings algorithm to locate the exact points at which the saturation occurs.

matlabroot/simulink/src/sfun_zc_sat.c

The S-function begins with #define statements for the S-function name and level, along with a #include statement for the simstruc.h header. After these statements, the S-function includes or defines any other necessary headers, data, etc. This example defines various parameters associated with the upper and lower saturation bounds.

/*  File    : sfun_zc_sat.c
 *  Abstract:
 *
 *      Example of an S-function which has nonsampled zero crossings to
 *      implement a saturation function. This S-function is designed to be
 *      used with a variable or fixed step solver.
 *
 *  A saturation is described by three equations
 *
 *    (1)     y = UpperLimit
 *    (2)     y = u
 *    (3)     y = LowerLimit
 *
 *  and a set of inequalities that specify which equation to use
 *
 *    if                          UpperLimit < u    then   use (1)
 *    if       LowerLimit <= u <= UpperLimit        then   use (2)
 *    if   u < LowerLimit                           then   use (3)
 *
 *  A key fact is that the valid equation 1, 2, or 3, can change at
 *  any instant.  Nonsampled zero crossing support helps the variable step
 *  solvers locate the exact instants when behavior switches from one equation
 *  to another.
 *
 *  Copyright 1990-2007 The MathWorks, Inc.
 */


#define S_FUNCTION_NAME  sfun_zc_sat
#define S_FUNCTION_LEVEL 2

#include "simstruc.h"

/*========================*
 * General Defines/macros *
 *========================*/

/* index to Upper Limit */
#define I_PAR_UPPER_LIMIT 0

/* index to Lower Limit */
#define I_PAR_LOWER_LIMIT 1

/* total number of block parameters */
#define N_PAR             2

/*
 *  Make access to mxArray pointers for parameters more readable.
 */
#define P_PAR_UPPER_LIMIT  ( ssGetSFcnParam(S,I_PAR_UPPER_LIMIT) )
#define P_PAR_LOWER_LIMIT  ( ssGetSFcnParam(S,I_PAR_LOWER_LIMIT) )

This S-function next implements the mdlCheckParameters method to check the validity of the S-function dialog parameters. Because this method is optional, a #define statement precedes it. The #if defined statement checks that this function is compiled as a MEX file, instead of for use with the Simulink Coder product. The body of the function performs basic checks to ensure that the user entered real vectors of equal length for the upper and lower saturation limits. If the parameter checks fail, the S-function errors out.

#define     MDL_CHECK_PARAMETERS
#if defined(MDL_CHECK_PARAMETERS) && defined(MATLAB_MEX_FILE)

  /* Function: mdlCheckParameters =============================================
   * Abstract:
   *   Check that parameter choices are allowable.
   */
  static void mdlCheckParameters(SimStruct *S)
  {
      int_T      i;
      int_T      numUpperLimit;
      int_T      numLowerLimit;
      const char *msg = NULL;

      /*
       * check parameter basics
       */
      for ( i = 0; i < N_PAR; i++ ) {
          if ( mxIsEmpty(    ssGetSFcnParam(S,i) ) ||
               mxIsSparse(   ssGetSFcnParam(S,i) ) ||
               mxIsComplex(  ssGetSFcnParam(S,i) ) ||
               !mxIsNumeric( ssGetSFcnParam(S,i) ) ) {
              msg = "Parameters must be real vectors.";
              goto EXIT_POINT;
          }
      }

      /*
       * Check sizes of parameters.
       */
      numUpperLimit = mxGetNumberOfElements( P_PAR_UPPER_LIMIT );
      numLowerLimit = mxGetNumberOfElements( P_PAR_LOWER_LIMIT );

      if ( ( numUpperLimit != 1             ) &&
           ( numLowerLimit != 1             ) &&
           ( numUpperLimit != numLowerLimit ) ) {
          msg = "Number of input and output values must be equal.";
          goto EXIT_POINT;
      }

      /*
       * Error exit point
       */
  EXIT_POINT:
      if (msg != NULL) {
          ssSetErrorStatus(S, msg);
      }
  }
#endif /* MDL_CHECK_PARAMETERS */

The required S-function method mdlInitializeSizes sets up the following S-function characteristics.

  • ssSetNumSFcnParams sets the number of expected S-function dialog parameters to two, as defined previously in the variable N_PAR.

  • If this method is compiled as a MEX file, ssGetSFcnParamsCount determines how many parameters the user actually entered into the S-function dialog. If the number of user-specified parameters matches the number returned by ssGetNumSFcnParams, the method calls mdlCheckParameters to check the validity of the user-entered data. Otherwise, the S-function errors out.

  • If the parameter check passes, the S-function determines the maximum number of elements entered into either the upper or lower saturation limit parameter. This number is needed later to determine the appropriate output width.

  • Next, the number of continuous and discrete states is set using ssSetNumContStates and ssSetNumDiscStates, respectively. This example has no continuous or discrete states.

  • The method specifies that the S-function has a single output port using ssSetNumOutputPorts and sets the width of this output port using ssSetOutputPortWidth. The output port width is either the maximum number of elements in the upper or lower saturation limit or is dynamically sized. Similar code specifies a single input port and indicates the input port has direct feedthrough by passing a value of 1 to ssSetInputPortDirectFeedThrough.

  • ssSetNumSampleTimes initializes one sample time, which the mdlInitializeSampleTimes function configures later.

  • The S-function indicates that no work vectors are used by passing a value of 0 to ssSetNumRWork, ssSetNumIWork, etc. You can omit these lines because zero is the default value for all of these macros. However, for clarity, the S-function explicitly sets the number of work vectors.

  • The method initializes the zero-crossing detection work vectors using ssSetNumModes and ssSetNumNonsampledZCs. The mdlSetWorkWidths method specifies the length of these dynamically sized vectors later.

  • Lastly, ssSetOptions sets any applicable options. In this case, SS_OPTION_EXCEPTION_FREE_CODE stipulates that the code is exception free and SS_OPTION_ALLOW_INPUT_SCALAR_EXPANSION permits scalar expansion of the input without having to provide an mdlSetInputPortWidth function.

The mdlInitializeSizes function for this example is shown below.

/* Function: mdlInitializeSizes ===============================================
 * Abstract:
 *   Initialize the sizes array.
 */
static void mdlInitializeSizes(SimStruct *S)
{
    int_T numUpperLimit, numLowerLimit, maxNumLimit;

    /*
     * Set and Check parameter count
     */
    ssSetNumSFcnParams(S, N_PAR);

#if defined(MATLAB_MEX_FILE)
    if (ssGetNumSFcnParams(S) == ssGetSFcnParamsCount(S)) {
        mdlCheckParameters(S);
        if (ssGetErrorStatus(S) != NULL) {
            return;
        }
    } else {
        return; /* Parameter mismatch reported by the Simulink engine*/
    }
#endif

    /*
     * Get parameter size info.
     */
    numUpperLimit = mxGetNumberOfElements( P_PAR_UPPER_LIMIT );
    numLowerLimit = mxGetNumberOfElements( P_PAR_LOWER_LIMIT );

    if (numUpperLimit > numLowerLimit) {
        maxNumLimit = numUpperLimit;
    } else {
        maxNumLimit = numLowerLimit;
    }

    /*
     * states
     */
    ssSetNumContStates(S, 0);
    ssSetNumDiscStates(S, 0);

    /*
     * outputs
     *   The upper and lower limits are scalar expanded
     *   so their size determines the size of the output
     *   only if at least one of them is not scalar.
     */
    if (!ssSetNumOutputPorts(S, 1)) return;

    if ( maxNumLimit > 1 ) {
        ssSetOutputPortWidth(S, 0, maxNumLimit);
    } else {
        ssSetOutputPortWidth(S, 0, DYNAMICALLY_SIZED);
    }

    /*
     * inputs
     *   If the upper or lower limits are not scalar then
     *   the input is set to the same size.  However, the
     *   ssSetOptions below allows the actual width to
     *   be reduced to 1 if needed for scalar expansion.
     */
    if (!ssSetNumInputPorts(S, 1)) return;

    ssSetInputPortDirectFeedThrough(S, 0, 1 );

    if ( maxNumLimit > 1 ) {
        ssSetInputPortWidth(S, 0, maxNumLimit);
    } else {
        ssSetInputPortWidth(S, 0, DYNAMICALLY_SIZED);
    }

    /*
     * sample times
     */
    ssSetNumSampleTimes(S, 1);

    /*
     * work
     */
    ssSetNumRWork(S, 0);
    ssSetNumIWork(S, 0);
    ssSetNumPWork(S, 0);


    /*
     * Modes and zero crossings:
     * If we have a variable-step solver and this block has a continuous
     * sample time, then
     *   o One mode element will be needed for each scalar output
     *     in order to specify which equation is valid (1), (2), or (3).
     *   o Two ZC elements will be needed for each scalar output
     *     in order to help the solver find the exact instants
     *     at which either of the two possible "equation switches"
     *     One will be for the switch from eq. (1) to (2);
     *     the other will be for eq. (2) to (3) and vice versa.
     * otherwise
     *   o No modes and nonsampled zero crossings will be used.
     *
     */
    ssSetNumModes(S, DYNAMICALLY_SIZED);
    ssSetNumNonsampledZCs(S, DYNAMICALLY_SIZED);

    /*
     * options
     *   o No mexFunctions and no problematic mxFunctions are called
     *     so the exception free code option safely gives faster simulations.
     *   o Scalar expansion of the inputs is desired.  The option provides
     *     this without the need to  write mdlSetOutputPortWidth and
     *     mdlSetInputPortWidth functions.
     */
    ssSetOptions(S, ( SS_OPTION_EXCEPTION_FREE_CODE |
                      SS_OPTION_ALLOW_INPUT_SCALAR_EXPANSION));

} /* end mdlInitializeSizes */

The required S-function method mdlInitializeSampleTimes specifies the S-function sample rates. The input argument INHERITED_SAMPLE_TIME passed to ssSetSampleTime specifies that this S-function inherits its sample time from its driving block. The call to ssSetModelReferenceSampleTimeDefaultInheritance tells the solver to use the default rule to determine if referenced models containing this S-function can inherit their sample times from the parent model.

/* Function: mdlInitializeSampleTimes =========================================
 * Abstract:
 *    Specify that the block is continuous.
 */
static void mdlInitializeSampleTimes(SimStruct *S)
{
    ssSetSampleTime(S, 0, INHERITED_SAMPLE_TIME);
    ssSetOffsetTime(S, 0, 0);
    ssSetModelReferenceSampleTimeDefaultInheritance(S);
}

The optional method mdlSetWorkWidths initializes the size of the zero-crossing detection work vectors. Because this method is optional, a #define statement precedes it. The #if defined statement checks that the S-function is being compiled as a MEX file. Zero-crossing detection can be done only when the S-function is running at a continuous sample rate using a variable-step solver. The if statement uses ssIsVariableStepSolver, ssGetSampleTime, and ssGetOffsetTime to determine if this condition is met. If so, the method sets the number of modes equal to the width of the first output port and the number of nonsampled zero crossings to twice this amount. Otherwise, the method sets both values to zero.

#define     MDL_SET_WORK_WIDTHS
#if defined(MDL_SET_WORK_WIDTHS) && defined(MATLAB_MEX_FILE)
/* Function: mdlSetWorkWidths ===============================================
 *   The width of the Modes and the ZCs depends on the width of the output.
 *   This width is not always known in mdlInitializeSizes so it is handled
 *   here.
 */
static void mdlSetWorkWidths(SimStruct *S)
{
    int nModes;
    int nNonsampledZCs;

    if (ssIsVariableStepSolver(S) &&
        ssGetSampleTime(S,0) == CONTINUOUS_SAMPLE_TIME &&
        ssGetOffsetTime(S,0) == 0.0) {

        int numOutput = ssGetOutputPortWidth(S, 0);

        /*
         * modes and zero crossings 
         *    o One mode element will be needed for each scalar output
         *      in order to specify which equation is valid (1), (2), or (3).
         *    o Two ZC elements will be needed for each scalar output
         *      in order to help the solver find the exact instants
         *      at which either of the two possible "equation switches"
         *      One will be for the switch from eq. (1) to (2);
         *      the other will be for eq. (2) to (3) and vice versa.
         */
        nModes         = numOutput;
        nNonsampledZCs = 2 * numOutput;
    } else {
        nModes         = 0;
        nNonsampledZCs = 0;
    }
    ssSetNumModes(S,nModes);
    ssSetNumNonsampledZCs(S,nNonsampledZCs);
}
#endif /* MDL_SET_WORK_WIDTHS */

After declaring variables for the input and output signals, the mdlOutputs functions uses an if-else statement to create blocks of code used to calculate the output signal based on whether the S-function uses a fixed-step or variable-step solver. The if statement queries the length of the nonsampled zero-crossing vector. If the length, set in mdlWorkWidths, is zero, then no zero-crossing detection is done and the output signals are calculated directly from the input signals. Otherwise, the function uses the mode work vector to determine how to calculate the output signal. If the simulation is at a major time step, i.e., ssIsMajorTimeStep returns true, mdlOutputs determines which mode the simulation is running in, either saturated at the upper limit, saturated at the lower limit, or not saturated. Then, for both major and minor time steps, the function calculates an output based on this mode. If the mode changed between the previous and current time step, then a zero crossing occurred. The mdlZeroCrossings function, not mdlOutputs, indicates this crossing to the solver.

/* Function: mdlOutputs =======================================================
 * Abstract:
 *
 *  A saturation is described by three equations
 *
 *    (1)     y = UpperLimit
 *    (2)     y = u
 *    (3)     y = LowerLimit
 *
 *  When this block is used with a fixed-step solver or it has a noncontinuous
 *  sample time, the equations are used as it
 *
 *  Now consider the case of this block being used with a variable-step solver
 *  and it has a continusous sample time. Solvers work best on smooth problems.
 *  In order for the solver to work without chattering, limit cycles, or
 *  similar problems, it is absolutely crucial that the same equation be used
 *  throughout the duration of a MajorTimeStep. To visualize this, consider
 *  the case of the Saturation block feeding an Integrator block.
 *
 *  To implement this rule, the mode vector is used to specify the
 *  valid equation based on the following:
 *
 *    if                          UpperLimit < u    then   use (1)
 *    if       LowerLimit <= u <= UpperLimit        then   use (2)
 *    if   u < LowerLimit                           then   use (3)
 *
 *  The mode vector is changed only at the beginning of a MajorTimeStep.
 *
 *  During a minor time step, the equation specified by the mode vector
 *  is used without question.  Most of the time, the value of u will agree
 *  with the equation specified by the mode vector.  However, sometimes u's
 *  value will indicate a different equation.  Nonetheless, the equation
 *  specified by the mode vector must be used.
 * 
 *  When the mode and u indicate different equations, the corresponding
 *  calculations are not correct.  However, this is not a problem.  From
 *  the ZC function, the solver will know that an equation switch occurred
 *  in the middle of the last MajorTimeStep.  The calculations for that
 *  time step will be discarded.  The ZC function will help the solver
 *  find the exact instant at which the switch occurred.  Using this knowledge,
 *  the length of the MajorTimeStep will be reduced so that only one equation
 *  is valid throughout the entire time step.
 */
static void mdlOutputs(SimStruct *S, int_T tid)
{
    InputRealPtrsType uPtrs     = ssGetInputPortRealSignalPtrs(S,0);
    real_T            *y        = ssGetOutputPortRealSignal(S,0);
    int_T             numOutput = ssGetOutputPortWidth(S,0);
    int_T             iOutput;

    /*
     * Set index and increment for input signal, upper limit, and lower limit 
     * parameters so that each gives scalar expansion if needed.
     */
    int_T  uIdx          = 0;
    int_T  uInc          = ( ssGetInputPortWidth(S,0) > 1 );
    const real_T *upperLimit   = mxGetPr( P_PAR_UPPER_LIMIT );
    int_T  upperLimitInc = ( mxGetNumberOfElements( P_PAR_UPPER_LIMIT ) > 1 );
    const real_T *lowerLimit   = mxGetPr( P_PAR_LOWER_LIMIT );
    int_T  lowerLimitInc = ( mxGetNumberOfElements( P_PAR_LOWER_LIMIT ) > 1 );

    UNUSED_ARG(tid); /* not used in single tasking mode */

    if (ssGetNumNonsampledZCs(S) == 0) {
        /*
         * This block is being used with a fixed-step solver or it has
         * a noncontinuous sample time, so we always saturate.
         */
        for (iOutput = 0; iOutput < numOutput; iOutput++) {
            if (*uPtrs[uIdx] >= *upperLimit) {
                *y++ = *upperLimit;
            } else if (*uPtrs[uIdx] > *lowerLimit) {
                *y++ = *uPtrs[uIdx];
            } else {
                *y++ = *lowerLimit;
            }

            upperLimit += upperLimitInc;
            lowerLimit += lowerLimitInc;
            uIdx       += uInc;
        }

    } else {
        /*
         * This block is being used with a variable-step solver.
         */
        int_T *mode = ssGetModeVector(S);

        /*
         * Specify indices for each equation.
         */
        enum { UpperLimitEquation, NonLimitEquation, LowerLimitEquation };

        /*
         * Update the Mode Vector ONLY at the beginning of a MajorTimeStep
         */
        if ( ssIsMajorTimeStep(S) ) {
            /*
             * Specify the mode, ie the valid equation for each output scalar.
             */
            for ( iOutput = 0; iOutput < numOutput; iOutput++ ) {
                if ( *uPtrs[uIdx] > *upperLimit ) {
                    /*
                     * Upper limit eq is valid.
                     */
                    mode[iOutput] = UpperLimitEquation;
                } else if ( *uPtrs[uIdx] < *lowerLimit ) {
                    /*
                     * Lower limit eq is valid.
                     */
                    mode[iOutput] = LowerLimitEquation;
                } else {
                    /*
                     * Nonlimit eq is valid.
                     */
                    mode[iOutput] = NonLimitEquation;
                }
                /*
                 * Adjust indices to give scalar expansion if needed.
                 */
                uIdx       += uInc;
                upperLimit += upperLimitInc;
                lowerLimit += lowerLimitInc;
            }

            /*
             * Reset index to input and limits.
             */
            uIdx       = 0;
            upperLimit = mxGetPr( P_PAR_UPPER_LIMIT );
            lowerLimit = mxGetPr( P_PAR_LOWER_LIMIT );

        } /* end IsMajorTimeStep */

        /*
         * For both MinorTimeSteps and MajorTimeSteps calculate each scalar
         * output using the equation specified by the mode vector.
         */
        for ( iOutput = 0; iOutput < numOutput; iOutput++ ) {
            if ( mode[iOutput] == UpperLimitEquation ) {
                /*
                 * Upper limit eq.
                 */
                *y++ = *upperLimit;
            } else if ( mode[iOutput] == LowerLimitEquation ) {
                /*
                 * Lower limit eq.
                 */
                *y++ = *lowerLimit;
            } else {
                /*
                 * Nonlimit eq.
                 */
                *y++ = *uPtrs[uIdx];
            }

            /*
             * Adjust indices to give scalar expansion if needed.
             */
            uIdx       += uInc;
            upperLimit += upperLimitInc;
            lowerLimit += lowerLimitInc;
        }
    }
} /* end mdlOutputs */

The mdlZeroCrossings method determines if a zero crossing occurred between the previous and current time step. The method obtains a pointer to the input signal using ssGetInputPortRealSignalPtrs. A comparison of this signal's value to the value of the upper and lower saturation limits determines values for the elements of the nonsampled zero-crossing vector. If any element of the nonsampled zero-crossing vector switches from negative to positive, or positive to negative, a zero crossing occurred. In the event of a zero crossing, the Simulink engine modifies the step size and recalculates the outputs to try to locate the exact zero crossing.

#define     MDL_ZERO_CROSSINGS
#if defined(MDL_ZERO_CROSSINGS) && (defined(MATLAB_MEX_FILE) || defined(NRT))

/* Function: mdlZeroCrossings =================================================
 * Abstract:
 *  This will only be called if the number of nonsampled zero crossings is
 *  greater than 0 which means this block has a continuous sample time and the
 *  model is using a variable-step solver.
 *
 *  Calculate zero crossing (ZC) signals that help the solver find the
 *  exact instants at which equation switches occur:
 *
 *    if                          UpperLimit < u    then   use (1)
 *    if       LowerLimit <= u <= UpperLimit        then   use (2)
 *    if   u < LowerLimit                           then   use (3)
 *
 *  The key words are help find. There is no choice of a function that will
 *  direct the solver to the exact instant of the change. The solver will
 *  track the zero crossing signal and do a bisection style search for the
 *  exact instant of equation switch.
 *
 *  There is generally one ZC signal for each pair of signals that can
 *  switch.  The three equations above would break into two pairs (1)&(2)
 *  and (2)&(3).  The  possibility of a "long jump" from (1) to (3) does
 *  not need to be handled as a separate case.  It is implicitly handled.
 *
 *  When ZCs are calculated, the value is normally used twice.  When it is
 *  first calculated, it is used as the end of the current time step.  Later,
 *  it will be used as the beginning of the following step.
 *
 *  The sign of the ZC signal always indicates an equation from the pair.  For
 *  S-functions, which equation is associated with a positive ZC and which is
 *  associated with a negative ZC doesn't really matter.  If the ZC is positive
 *  at the beginning and at the end of the time step, this implies that the
 *  "positive" equation was valid throughout the time step.  Likewise, if the
 *  ZC is negative at the beginning and at the end of the time step, this
 *  implies that the "negative" equation was valid throughout the time step.
 *  Like any other nonlinear solver, this is not foolproof, but it is an
 *  excellent indicator.  If the ZC has a different sign at the beginning and
 *  at the end of the time step, then a equation switch definitely occurred
 *  during the time step.
 *
 *  Ideally, the ZC signal gives an estimate of when an equation switch
 *  occurred.  For example, if the ZC signal is -2 at the beginning and +6 at
 *  the end, then this suggests that the switch occurred 
 *  25% = 100%*(-2)/(-2-(+6)) of the way into the time step.  It will almost
 *  never be true that 25% is perfectly correct.  There is no perfect choice
 *  for a ZC signal, but there are some good rules.  First, choose the ZC
 *  signal to be continuous.  Second, choose the ZC signal to give a monotonic
 *  measure of the "distance" to a signal switch; strictly monotonic is ideal.
 */
static void mdlZeroCrossings(SimStruct *S)
{
    int_T             iOutput;
    int_T             numOutput = ssGetOutputPortWidth(S,0);
    real_T            *zcSignals = ssGetNonsampledZCs(S);
    InputRealPtrsType uPtrs      = ssGetInputPortRealSignalPtrs(S,0);

    /*
     * Set index and increment for the input signal, upper limit, and lower 
     * limit parameters so that each gives scalar expansion if needed.
     */
    int_T  uIdx          = 0;
    int_T  uInc          = ( ssGetInputPortWidth(S,0) > 1 );
    real_T *upperLimit   = mxGetPr( P_PAR_UPPER_LIMIT );
    int_T  upperLimitInc = ( mxGetNumberOfElements( P_PAR_UPPER_LIMIT ) > 1 );
    real_T *lowerLimit   = mxGetPr( P_PAR_LOWER_LIMIT );
    int_T  lowerLimitInc = ( mxGetNumberOfElements( P_PAR_LOWER_LIMIT ) > 1 );

    /*
     * For each output scalar, give the solver a measure of "how close things
     * are" to an equation switch.
     */
    for ( iOutput = 0; iOutput < numOutput; iOutput++ ) {

        /*  The switch from eq (1) to eq (2)
         *
         *    if                          UpperLimit < u    then   use (1)
         *    if       LowerLimit <= u <= UpperLimit        then   use (2)
         *
         *  is related to how close u is to UpperLimit.  A ZC choice
         *  that is continuous, strictly monotonic, and is
         *    u - UpperLimit 
         *  or it is negative.
         */
        zcSignals[2*iOutput] = *uPtrs[uIdx] - *upperLimit;

        /*  The switch from eq (2) to eq (3)
         *
         *    if       LowerLimit <= u <= UpperLimit        then   use (2)
         *    if   u < LowerLimit                           then   use (3)
         *
         *  is related to how close u is to LowerLimit.  A ZC choice
         *  that is continuous, strictly monotonic, and is
         *    u - LowerLimit.
         */
        zcSignals[2*iOutput+1] = *uPtrs[uIdx] - *lowerLimit;

        /*
         * Adjust indices to give scalar expansion if needed.
         */
        uIdx       += uInc;
        upperLimit += upperLimitInc;
        lowerLimit += lowerLimitInc;
    }
}

#endif /* end mdlZeroCrossings */

The S-function concludes with the required mdlTerminate function. In this example, the function is empty.

/* Function: mdlTerminate =====================================================
 * Abstract:
 *    No termination needed, but we are required to have this routine.
 */
static void mdlTerminate(SimStruct *S)
{
    UNUSED_ARG(S); /* unused input argument */
}

The required S-function trailer includes the files necessary for simulation or code generation, as follows.

#ifdef  MATLAB_MEX_FILE    /* Is this file being compiled as a MEX file? */
#include "simulink.c"      /* MEX file interface mechanism */
#else
#include "cg_sfun.h"       /* Code generation registration function */
#endif

    Note   The mdlOutputs and mdlTerminate functions use the UNUSED_ARG macro to indicate that an input argument the callback requires is not used. This optional macro is defined in simstruc_types.hsimstruc_types.h. If used, you must call this macro once for each input argument that a callback does not use.

Discontinuities in Continuous States

The example S-function stvctf.cstvctf.c demonstrates a time-varying continuous transfer function. The following Simulink model uses this S-function.

sfcndemo_stvctfsfcndemo_stvctf

The S-function demonstrates how to work with the solvers so that the simulation maintains consistency, which means that the block maintains smooth and consistent signals for the integrators although the equations that are being integrated are changing.

matlabroot/simulink/src/stvctf.c

The S-function begins with #define statements for the S-function name and level, along with a #include statement for the simstruc.h header. After these statements, the S-function includes or defines any other necessary headers, data, etc. This example defines parameters for the transfer function's numerator and denominator, which are entered into the S-function dialog. The comments at the beginning of this S-function provide additional information on the purpose of the work vectors in this example.

/*
 * File : stvctf.c
 * Abstract:
 *      Time Varying Continuous Transfer Function block
 *
 *      This S-function implements a continuous time transfer function
 *      whose transfer function polynomials are passed in via the input
 *      vector.  This is useful for continuous time adaptive control
 *      applications.
 *
 *      This S-function is also an example of how to use banks to avoid
 *      problems with computing derivatives when a continuous output has 
 *      discontinuities. The consistency checker can be used to verify that 
 *      your S-function is correct with respect to always maintaining smooth 
 *      and consistent signals for the integrators. By consistent we mean that
 *      two mdlOutputs calls at major time t and minor time t are always the 
 *      same. The consistency checker is enabled on the diagnostics page of the
 * 		Configuraion parameters dialog box. The update method of this S-function
 *      modifies the coefficients of the transfer function, which cause the
 *      output to "jump." To have the simulation work properly, we need to let
 *      the solver know of these discontinuities by setting 
 *      ssSetSolverNeedsReset and then we need to use multiple banks of 
 *      coefficients so the coefficients used in the major time step output 
 *      and the minor time step outputs are the same. In the simulation loop 
 *      we have:
 *        Loop:
 *          o Output in major time step at time t
 *          o Update in major time step at time t
 *          o Integrate (minor time step):
 *              o Consistency check: recompute outputs at time t and compare
 *                with current outputs.
 *              o Derivatives at time t
 *              o One or more Output,Derivative evaluations at time t+k
 *                where k <= step_size to be taken.
 *              o Compute state, x
 *              o t = t + step_size
 *            End_Integrate
 *        End_Loop
 *      Another purpose of the consistency checker is to verify that when 
 *      the solver needs to try a smaller step_size, the recomputing of 
 *      the output and derivatives at time t doesn't change. Step size 
 *      reduction occurs when tolerances aren't met for the current step size.
 *      The ideal ordering would be to update after integrate. To achieve 
 *      this we have two banks of coefficients. And the use of the new 
 *      coefficients, which were computed in update, is delayed until after 
 *      the integrate phase is complete.
 *
 *	This block has multiple sample times and will not work correctly
 *	in a multitasking environment. It is designed to be used in
 *	a single tasking (or variable step) simulation environment.
 *	Because this block accesses the input signal in both tasks,
 *	it cannot specify the sample times of the input and output ports
 *	(SS_OPTION_PORT_SAMPLE_TIMES_ASSIGNED).
 *
 * See simulink/src/sfuntmpl_doc.c.
 *
 * Copyright 1990-7 The MathWorks, Inc.
 */

#define S_FUNCTION_NAME  stvctf
#define S_FUNCTION_LEVEL 2

#include "simstruc.h"

/*
 * Defines for easy access to the numerator and denominator polynomials
 * parameters
 */
#define NUM(S)  ssGetSFcnParam(S, 0)
#define DEN(S)  ssGetSFcnParam(S, 1)
#define TS(S)   ssGetSFcnParam(S, 2)
#define NPARAMS 3

This S-function implements the mdlCheckParameters method to check the validity of the S-function dialog parameters. Because this method is optional, a #define statement precedes it. The #if defined statement checks that this function is compiled as a MEX file, instead of for use with the Simulink Coder product. The body of the function performs basic checks to ensure that the user entered real vectors for the numerator and denominator, and that the denominator has a higher order than the numerator. If the parameter check fails, the S-function errors out.

#define MDL_CHECK_PARAMETERS
#if defined(MDL_CHECK_PARAMETERS) && defined(MATLAB_MEX_FILE)
  /* Function: mdlCheckParameters =============================================
   * Abstract:
   *    Validate our parameters to verify:
   *     o The numerator must be of a lower order than the denominator.
   *     o The sample time must be a real positive nonzero value.
   */
  static void mdlCheckParameters(SimStruct *S)
  {
      int_T i;

      for (i = 0; i < NPARAMS; i++) {
          real_T *pr;
          int_T   el;
          int_T   nEls;
          if (mxIsEmpty(    ssGetSFcnParam(S,i)) ||
              mxIsSparse(   ssGetSFcnParam(S,i)) ||
              mxIsComplex(  ssGetSFcnParam(S,i)) ||
              !mxIsNumeric( ssGetSFcnParam(S,i)) ) {
              ssSetErrorStatus(S,"Parameters must be real finite vectors");
              return;
          } 
          pr   = mxGetPr(ssGetSFcnParam(S,i));
          nEls = mxGetNumberOfElements(ssGetSFcnParam(S,i));
          for (el = 0; el < nEls; el++) {
              if (!mxIsFinite(pr[el])) {
                  ssSetErrorStatus(S,"Parameters must be real finite vectors");
                  return;
              }
          }
      }

      if (mxGetNumberOfElements(NUM(S)) > mxGetNumberOfElements(DEN(S)) &&
          mxGetNumberOfElements(DEN(S)) > 0  && *mxGetPr(DEN(S)) != 0.0) {
          ssSetErrorStatus(S,"The denominator must be of higher order than "
                           "the numerator, nonempty and with first "
                           "element nonzero");
          return;
      }

      /* xxx verify finite */
      if (mxGetNumberOfElements(TS(S)) != 1 || mxGetPr(TS(S))[0] <= 0.0) {
          ssSetErrorStatus(S,"Invalid sample time specified");
          return;
      }
  }
#endif /* MDL_CHECK_PARAMETERS */

The required S-function method mdlInitializeSizes then sets up the following S-function characteristics.

  • ssSetNumSFcnParams sets the number of expected S-function dialog parameters to three, as defined previously in the variable NPARAMS.

  • If this method is compiled as a MEX file, ssGetSFcnParamsCount determines how many parameters the user entered into the S-function dialog. If the number of user-specified parameters matches the number returned by ssGetNumSFcnParams, the method calls mdlCheckParameters to check the validity of the user-entered data. Otherwise, the S-function errors out.

  • If the parameter check passes, the S-function specifies the number of continuous and discrete states using ssSetNumContStates and ssSetNumDiscStates, respectively. This example has no discrete states and sets the number of continuous states based on the number of coefficients in the transfer function's denominator.

  • Next, ssSetNumInputPorts specifies that the S-function has a single input port and sets its width to one plus twice the length of the denominator using ssSetInputPortWidth. The method uses the value provided by the third S-function dialog parameter as the input port's sample time. This parameter indicates the rate at which the transfer function is modified during simulation. The S-function specifies that the input port has direct feedthrough by passing a value of 1 to ssSetInputPortDirectFeedThrough.

  • ssSetNumOutputPorts specifies that the S-function has a single output port. The method uses ssSetOutputPortWidth to set the width of this output port, ssSetOutputPortSampleTime to specify that the output port has a continuous sample time, and ssSetOutputPortOffsetTime to set the offset time to zero.

  • ssSetNumSampleTimes then initializes two sample times, which the mdlInitializeSampleTimes function configures later.

  • The method passes a value of four times the number of denominator coefficients to ssSetNumRWork in order to set the length of the floating-point work vector. ssSetNumIWork then sets the length of the integer work vector to two. The RWork vectors store two banks of transfer function coefficients, while the IWork vector indicates which bank in the RWork vector is currently in use. The S-function sets the length of all other work vectors to zero. You can omit these lines because zero is the default value for these macros. However, for clarity, the S-function explicitly sets the number of work vectors.

  • Lastly, ssSetOptions sets any applicable options. In this case, SS_OPTION_EXCEPTION_FREE_CODE stipulates that the code is exception free.

The mdlInitializeSizes function for this example is shown below.

/* Function: mdlInitializeSizes ===============================================
 * Abstract:
 *    Determine the S-function block's characteristics:
 *    number of inputs, outputs, states, etc.
 */
static void mdlInitializeSizes(SimStruct *S)
{
    int_T nContStates;
    int_T nCoeffs;

    /* See sfuntmpl_doc.c for more details on the macros below. */
 
    ssSetNumSFcnParams(S, NPARAMS);  /* Number of expected parameters. */
#if defined(MATLAB_MEX_FILE)
    if (ssGetNumSFcnParams(S) == ssGetSFcnParamsCount(S)) {
        mdlCheckParameters(S);
        if (ssGetErrorStatus(S) != NULL) {
            return;
        }
    } else {
        return; /* Parameter mismatch reported by the Simulink engine*/
    }
#endif


    /*
     * Define the characteristics of the block:
     *
     *   Number of continuous states:     length of denominator - 1
     *   Inputs port width                2 * (NumContStates+1) + 1
     *   Output port width                1
     *   DirectFeedThrough:               0 (Although this should be computed.
     *                                       We'll assume coefficients entered
     *                                       are strictly proper).
     *   Number of sample times:          2 (continuous and discrete)
     *   Number of Real work elements:    4*NumCoeffs
     *                                    (Two banks for num and den coeff's:
     *                                     NumBank0Coeffs
     *                                     DenBank0Coeffs
     *                                     NumBank1Coeffs
     *                                     DenBank1Coeffs)
     *   Number of Integer work elements: 2 (indicator of active bank 0 or 1
     *                                       and flag to indicate when banks
     *                                       have been updated).
     *
     * The number of inputs arises from the following:
     *   o 1 input (u)
     *   o the numerator and denominator polynomials each have NumContStates+1
     *     coefficients
     */
    nCoeffs     = mxGetNumberOfElements(DEN(S));
    nContStates = nCoeffs - 1;

    ssSetNumContStates(S, nContStates);
    ssSetNumDiscStates(S, 0);

    if (!ssSetNumInputPorts(S, 1)) return;
    ssSetInputPortWidth(S, 0, 1 + (2*nCoeffs));
    ssSetInputPortDirectFeedThrough(S, 0, 0);
    ssSetInputPortSampleTime(S, 0, mxGetPr(TS(S))[0]);
    ssSetInputPortOffsetTime(S, 0, 0);

    if (!ssSetNumOutputPorts(S,1)) return;
    ssSetOutputPortWidth(S, 0, 1);
    ssSetOutputPortSampleTime(S, 0, CONTINUOUS_SAMPLE_TIME);
    ssSetOutputPortOffsetTime(S, 0, 0);

    ssSetNumSampleTimes(S, 2);

    ssSetNumRWork(S, 4 * nCoeffs);
    ssSetNumIWork(S, 2);
    ssSetNumPWork(S, 0);

    ssSetNumModes(S, 0);
    ssSetNumNonsampledZCs(S, 0);

    /* Take care when specifying exception free code - see sfuntmpl_doc.c */
    ssSetOptions(S, (SS_OPTION_EXCEPTION_FREE_CODE));

} /* end mdlInitializeSizes */

The required S-function method mdlInitializeSampleTimes specifies the S-function sample rates. The first call to ssSetSampleTime specifies that the first sample rate is continuous and the subsequent call to ssSetOffsetTime sets the offset to zero. The second call to this pair of macros sets the second sample time to the value of the third S-function parameter with an offset of zero. The call to ssSetModelReferenceSampleTimeDefaultInheritance tells the solver to use the default rule to determine if referenced models containing this S-function can inherit their sample times from the parent model.

/* Function: mdlInitializeSampleTimes =========================================
 * Abstract:
 *      This function is used to specify the sample time(s) for the
 *      S-function.  This S-function has two sample times.  The
 *      first, a continuous sample time, is used for the input to the
 *      transfer function, u.  The second, a discrete sample time
 *      provided by the user, defines the rate at which the transfer
 *      function coefficients are updated.
 */
static void mdlInitializeSampleTimes(SimStruct *S)
{
    /*
     * the first sample time, continuous
     */
    ssSetSampleTime(S, 0, CONTINUOUS_SAMPLE_TIME);
    ssSetOffsetTime(S, 0, 0.0);

    /*
     * the second, discrete sample time, is user provided
     */
    ssSetSampleTime(S, 1, mxGetPr(TS(S))[0]);
    ssSetOffsetTime(S, 1, 0.0);
    ssSetModelReferenceSampleTimeDefaultInheritance(S);


} /* end mdlInitializeSampleTimes */

The optional S-function method mdlInitializeConditions initializes the continuous state vector and the initial numerator and denominator vectors. The #define statement before this method is required for the Simulink engine to call this function. The function initializes the continuous states to zero. The numerator and denominator coefficients are initialized from the first two S-function parameters, normalized by the first denominator coefficient. The function sets the value stored in the IWork vector to zero, to indicate that the first bank of numerator and denominator coefficients stored in the RWork vector is currently in use.

#define MDL_INITIALIZE_CONDITIONS
/* Function: mdlInitializeConditions ==========================================
 * Abstract:
 *      Initalize the states, numerator and denominator coefficients.
 */
static void mdlInitializeConditions(SimStruct *S)
{
    int_T  i;
    int_T  nContStates = ssGetNumContStates(S);
    real_T *x0           = ssGetContStates(S);
    int_T  nCoeffs       = nContStates + 1;
    real_T *numBank0     = ssGetRWork(S);
    real_T *denBank0     = numBank0 + nCoeffs;
    int_T *activeBank    = ssGetIWork(S);

    /*
     * The continuous states are all initialized to zero.
     */
    for (i = 0; i < nContStates; i++) {
        x0[i]       = 0.0;
        numBank0[i] = 0.0;
        denBank0[i] = 0.0;
    }
    numBank0[nContStates] = 0.0;
    denBank0[nContStates] = 0.0;

    /*
     * Set up the initial numerator and denominator.
     */
    {
        const real_T *numParam   = mxGetPr(NUM(S));
        int          numParamLen = mxGetNumberOfElements(NUM(S));

        const real_T *denParam   = mxGetPr(DEN(S));
        int          denParamLen = mxGetNumberOfElements(DEN(S));
        real_T       den0        = denParam[0];

        for (i = 0; i < denParamLen; i++) {
            denBank0[i] = denParam[i] / den0;
        }

        for (i = 0; i < numParamLen; i++) {
            numBank0[i] = numParam[i] / den0;
        }
    }

    /*
     * Normalize if this transfer function has direct feedthrough.
     */
    for (i = 1; i < nCoeffs; i++) {
        numBank0[i] -= denBank0[i]*numBank0[0];
    }

    /*
     * Indicate bank0 is active (i.e. bank1 is oldest).
     */
    *activeBank = 0;

} /* end mdlInitializeConditions */

The mdlOutputs function calculates the S-function output signals when the S-function is simulating in a continuous task, i.e., ssIsContinuousTask is true. If the simulation is also at a major time step, mdlOutputs checks if the numerator and denominator coefficients need to be updated, as indicated by a switch in the active bank stored in the IWork vector. At both major and minor time steps, the S-function calculates the output using the numerator coefficients stored in the active bank.

/* Function: mdlOutputs =======================================================
 * Abstract:
 *      The outputs for this block are computed by using a controllable state-
 *      space representation of the transfer function.
 */
static void mdlOutputs(SimStruct *S, int_T tid)
{
    if (ssIsContinuousTask(S,tid)) {
        int               i;
        real_T            *num;
        int               nContStates = ssGetNumContStates(S);
        real_T            *x          = ssGetContStates(S);
        int_T             nCoeffs     = nContStates + 1;
        InputRealPtrsType uPtrs       = ssGetInputPortRealSignalPtrs(S,0);
        real_T            *y          = ssGetOutputPortRealSignal(S,0);
        int_T             *activeBank = ssGetIWork(S);

        /*
         * Switch banks because we've updated them in mdlUpdate and we're no 
         * longer in a minor time step.
         */
        if (ssIsMajorTimeStep(S)) {
            int_T *banksUpdated = ssGetIWork(S) + 1;
            if (*banksUpdated) {
                *activeBank = !(*activeBank);
                *banksUpdated = 0;
                /*
                 * Need to tell the solvers that the derivatives are no
                 * longer valid.
                 */
                ssSetSolverNeedsReset(S);
            }
        }
        num = ssGetRWork(S) + (*activeBank) * (2*nCoeffs);

        /*
         * The continuous system is evaluated using a controllable state space
         * representation of the transfer function.  This implies that the
         * output of the system is equal to:
         *
         *     y(t) = Cx(t) + Du(t)
         *          = [ b1 b2 ... bn]x(t) + b0u(t)
         *
         * where b0, b1, b2, ... are the coefficients of the numerator 
         * polynomial:
         *
         *    B(s) = b0 s^n + b1 s^n-1 + b2 s^n-2 + ... + bn-1 s + bn
         */
        *y = *num++ * (*uPtrs[0]);
        for (i = 0; i < nContStates; i++) {
            *y += *num++ * *x++;
        }
    }

} /* end mdlOutputs */

Although this example has no discrete states, the method still implements the mdlUpdate function to update the transfer function coefficients at every major time step. Because this method is optional, a #define statement precedes it. The method uses ssGetInputPortRealSignalPtrs to obtain a pointer to the input signal. The input signal's values become the new transfer function coefficients, which the S-function stores in the bank of the inactive RWork vector. When the mdlOutputs function is later called at this major time step, it updates the active bank to be this updated bank of coefficients.

#define MDL_UPDATE
/* Function: mdlUpdate ========================================================
 * Abstract:
 *      Every time through the simulation loop, update the
 *      transfer function coefficients. Here we update the oldest bank.
 */
static void mdlUpdate(SimStruct *S, int_T tid)
{
    if (ssIsSampleHit(S, 1, tid)) {
        int_T             i;
        InputRealPtrsType uPtrs        = ssGetInputPortRealSignalPtrs(S,0);
        int_T             uIdx         = 1;/*1st coeff is after signal input*/
        int_T             nContStates  = ssGetNumContStates(S);
        int_T             nCoeffs      = nContStates + 1;
        int_T             bankToUpdate = !ssGetIWork(S)[0];
        real_T            *num         = ssGetRWork(S)+bankToUpdate*2*nCoeffs;
        real_T            *den         = num + nCoeffs;

        real_T            den0;
        int_T             allZero;


        /*
         * Get the first denominator coefficient.  It will be used
         * for normalizing the numerator and denominator coefficients.
         *
         * If all inputs are zero, we probably could have unconnected
         * inputs, so use the parameter as the first denominator coefficient.
         */
        den0 = *uPtrs[uIdx+nCoeffs];
        if (den0 == 0.0) {
            den0 = mxGetPr(DEN(S))[0];
        }

        /*
         * Grab the numerator.
         */
        allZero = 1;
        for (i = 0; (i < nCoeffs) && allZero; i++) {
            allZero &= *uPtrs[uIdx+i] == 0.0;
        }

        if (allZero) { /* if numerator is all zero */
            const real_T *numParam   = mxGetPr(NUM(S));
            int_T        numParamLen = mxGetNumberOfElements(NUM(S));
            /*
             * Move the input to the denominator input and
             * get the denominator from the input parameter.
             */
            uIdx += nCoeffs;
            num += nCoeffs - numParamLen;
            for (i = 0; i < numParamLen; i++) {
                *num++ = *numParam++ / den0;
            }
        } else {
            for (i = 0; i < nCoeffs; i++) {
                *num++ = *uPtrs[uIdx++] / den0;
            }
        }

        /*
         * Grab the denominator.
         */
        allZero = 1;
        for (i = 0; (i < nCoeffs) && allZero; i++) {
            allZero &= *uPtrs[uIdx+i] == 0.0;
        }

        if (allZero) {  /* If denominator is all zero. */
            const real_T *denParam   = mxGetPr(DEN(S));
            int_T        denParamLen = mxGetNumberOfElements(DEN(S));

            den0 = denParam[0];
            for (i = 0; i < denParamLen; i++) {
                *den++ = *denParam++ / den0;
            }
        } else {
            for (i = 0; i < nCoeffs; i++) {
                *den++ = *uPtrs[uIdx++] / den0;
            }
        }

        /*
         * Normalize if this transfer function has direct feedthrough.
         */
        num = ssGetRWork(S) + bankToUpdate*2*nCoeffs;
        den = num + nCoeffs;
        for (i = 1; i < nCoeffs; i++) {
            num[i] -= den[i]*num[0];
        }

        /*
         * Indicate oldest bank has been updated.
         */
        ssGetIWork(S)[1] = 1;
    }

} /* end mdlUpdate */

The mdlDerivatives function calculates the continuous state derivatives. The function uses the coefficients from the active bank to solve a controllable state-space representation of the transfer function.

#define MDL_DERIVATIVES
/* Function: mdlDerivatives ===================================================
 * Abstract:
 *      The derivatives for this block are computed by using a controllable 
 *      state-space representation of the transfer function.
 */
static void mdlDerivatives(SimStruct *S) 
{
    int_T             i;
    int_T             nContStates = ssGetNumContStates(S);
    real_T            *x          = ssGetContStates(S);
    real_T            *dx         = ssGetdX(S);
    int_T             nCoeffs     = nContStates + 1;
    int_T             activeBank  = ssGetIWork(S)[0];
    const real_T      *num        = ssGetRWork(S) + activeBank*(2*nCoeffs);
    const real_T      *den        = num + nCoeffs;
    InputRealPtrsType uPtrs       = ssGetInputPortRealSignalPtrs(S,0);

    /*
     * The continuous system is evaluated using a controllable state-space
     * representation of the transfer function.  This implies that the
     * next continuous states are computed using:
     *
     *     dx = Ax(t) + Bu(t)
     *        = [-a1 -a2 ... -an] [x1(t)] + [u(t)]
     *          [  1  0  ...   0] [x2(t)] + [0]
     *          [  0  1  ...   0] [x3(t)] + [0]
     *          [  .  .  ...   .]   .     +  .
     *          [  .  .  ...   .]   .     +  .
     *          [  .  .  ...   .]   .     +  .
     *          [  0  0  ... 1 0] [xn(t)] + [0]
     *
     * where a1, a2, ... are the coefficients of the numerator polynomial:
     *
     *    A(s) = s^n + a1 s^n-1 + a2 s^n-2 + ... + an-1 s + an
     */
    dx[0] = -den[1] * x[0] + *uPtrs[0];
    for (i = 1; i < nContStates; i++) {
        dx[i] = x[i-1];
        dx[0] -= den[i+1] * x[i];
    }

} /* end mdlDerivatives */

The required mdlTerminate function performs any actions, such as freeing memory, necessary at the end of the simulation. In this example, the function is empty.

/* Function: mdlTerminate =====================================================
 * Abstract:
 *      Called when the simulation is terminated.
 *      For this block, there are no end of simulation tasks.
 */
static void mdlTerminate(SimStruct *S)
{
    UNUSED_ARG(S); /* unused input argument */
} /* end mdlTerminate */

The required S-function trailer includes the files necessary for simulation or code generation, as follows.

#ifdef  MATLAB_MEX_FILE    /* Is this file being compiled as a MEX file? */
#include "simulink.c"      /* MEX file interface mechanism */
#else
#include "cg_sfun.h"       /* Code generation registration function */
#endif

    Note   The mdlTerminate function uses the UNUSED_ARG macro to indicate that an input argument the callback requires is not used. This optional macro is defined in simstruc_types.hsimstruc_types.h. If used, you must call this macro once for each input argument that a callback does not use.

Was this topic helpful?