Hello there. It’s been a while, right ?

The last two weeks weren’t so great: I was haunted by segmentation faults all the time when trying to run the OpenModelica generated model library. Therefore, I couldn’t achieve any new meaningful result worth blogging about.

The idea presented here in my previous posts of issuing different FMI2 calls for each phase (flag) of computational function processing was somehow causing crashes, even if the Scicos block and FMI2 model data pointers seemed ok, and I simply didn’t have a clue about why.

A New Approach

One thing that I’ve noted while trying different fixes for my problem was that FMI2 calls right after initialization (fmi2Intantiate and fmi2SetupExperiment) always went smoothly. So, in a more radical (or conservative, depending on how you see it) attempt to prevent crashes I’ve reestructured my code to perform model creation and termination at each computational function call:

/* ... */

SCICOS_BLOCK_EXPORT void BLOCK_FUNCTION_NAME(scicos_block* block, const int flag)
{    
    /* ... */

    fmi2Component model = fmi2Instantiate( block->label,          // Instance name
                                           fmi2ModelExchange,     // Exported model type
                                           block->uid,            // Model GUID, have to be obtained
                                           "",                    // Optional FMU resources location
                                           &functions,            // User provided callbacks
                                           fmi2False,             // Interactive mode Off
                                           fmi2True );            // Logging On

    /* FMU processing. Set inputs, get outputs, etc. */
    
    fmi2Terminate( model );       // Terminate simulation for this component
    fmi2FreeInstance( model );    // Deallocate memory
    
    switch( flag )
    {
        /* Return different values for different flags/phases */
    }
    
    return;
}

and finally get some results (Yay !):

(Resulting plot for chaos model simulation)

Refinements

As expected, the process of creating and destroying the FMI2 component structure on every processing step is expensive, mostly, I suppose, due to allocation and deallocation of memory, which also violates the rule of just performing it for flags 4 (initialization) and 5 (termination). I’ve verified that when seeing the simulation running in a much slower speed than with the model library generated by modelicac compiler.

Gladly, keeping the model stored in a static variable, calling fmi2Instantiate only on the first computation step, fmi2Reset on subsequential ones and fmi2FreeInstance only in the end, also works. That change, combined with disabling FMU logs and other text output functions, allows things to run almost as fast as with the older compiler.

As a result, the reestructured code for fmi2_wrapper.h looks like this:

#include <string.h>
#include <stdio.h>
#include <stdarg.h>

#include "scicos_block.h"
#include "scicos_malloc.h"
#include "scicos_free.h"
#ifdef _MSC_VER
#define scicos_print printf
#define SCICOS_BLOCK_EXPORT __declspec(dllexport)
#else
#include "scicos_print.h"
#define SCICOS_BLOCK_EXPORT
#endif

// Simple terminal message logging function
static void print_log_message( fmi2ComponentEnvironment componentEnvironment, fmi2String instanceName,
                               fmi2Status status, fmi2String category, fmi2String message, ... )
{
    char messageBuffer[ 1024 ];

    va_list args;
    va_start( args, message );
    vsnprintf( messageBuffer, 1024, (const char*) message, args );
    va_end( args );
}

static void* scicos_calloc( size_t objects_number, size_t object_size )
{
    return scicos_malloc( objects_number * object_size );
}

// User provided functions for FMI2 data management
const fmi2CallbackFunctions functions = { .logger = print_log_message,
                                          .allocateMemory = scicos_calloc,
                                          .freeMemory = scicos_free,
                                          .stepFinished = NULL,
                                          .componentEnvironment = NULL };

// The computational function itself
SCICOS_BLOCK_EXPORT void BLOCK_FUNCTION_NAME(scicos_block* block, const int flag)
{    
    static fmi2Component model;
  
    double **y = block->outptr;
    double **u = block->inptr;
    
    fmi2EventInfo eventInfo;
    
    int i;
    int enterEventMode, terminateSimulation;
    
    // Flag 4: State initialisation and memory allocation
    if( flag == 4 )
    {
        // Store FMI2 component in the work field
        model = fmi2Instantiate( block->label,          // Instance name
                                 fmi2ModelExchange,     // Exported model type
                                 block->uid,            // Model GUID, have to be obtained
                                 "",                    // Optional FMU resources location
                                 &functions,            // User provided callbacks
                                 fmi2False,             // Interactive mode
                                 fmi2False );           // Logging Off
    }
    else
    {
        fmi2Reset( model );
    }
    
    // Define simulation parameters. Internally calls state and event setting functions,
    // which should be called before any model evaluation/event processing ones
    fmi2SetupExperiment( model,                                   // FMI2 component
                         fmi2False,                               // Tolerance undefined
                         0.0,                                     // Tolerance value (not used)
                         0.0,                                     // Start time
                         fmi2False,                               // Stop time undefined
                         1.0 );                                   // Stop time (not used)
    
    // FMI2 component initialization
    fmi2EnterInitializationMode( model );
    
    if( block->nin > 0 )
    {
        for( i = 0; i < block->nin; i++ )
            inputsList[ i ] = (fmi2Real) u[ i ][ 0 ];
        fmi2SetReal( model, INPUT_REFS_LIST, block->nin, inputsList );
    }
    
    fmi2ExitInitializationMode( model );
    
    fmi2EnterContinuousTimeMode( model );
    
    fmi2SetTime( model, (fmi2Real) get_scicos_time() );
    
    for( i = 0; i < block->nx; i++ )
    {
        statesList[ i ] = (fmi2Real) block->x[ i ];
    }
    fmi2SetContinuousStates( model, statesList, block->nx );
    
    fmi2GetEventIndicators( model, eventIndicatorsList, block->ng );
    
    // Inform the model about an accepted step
    // The second parameter informs that the model simulation won't be set to a prior time instant (fmi2False otherwise)
    fmi2CompletedIntegratorStep( model, fmi2True, &enterEventMode, &terminateSimulation );
    if( terminateSimulation ) set_block_error( -3 );
    
    /* Discrete events handling */ 
    
    fmi2GetDerivatives( model, stateDerivativesList, block->nx );
    
    fmi2GetReal( model, OUTPUT_REFS_LIST, block->nout, outputsList );
    
    fmi2Terminate( model );       // Terminate simulation for this component
    
    // Flag 0: Update continuous state
    if( flag == 0 )
    {
        // Output both residuals for implicit (DAE) solvers
        for( i = 0; i < block->nx; i++ )
        {
            block->res[ i ] = (double) stateDerivativesList[ i ] - block->xd[ i ];
        }
    }
    // Flag 1: Update output state
    else if( flag == 1 )
    {
        for( i = 0; i < block->nout; i++ )
        {
            y[ i ][ 0 ] = (double) outputsList[ i ];
        }
    }
    // Flag 5: simulation end and memory deallocation
    else if( flag == 5 )
    {
        fmi2FreeInstance( model );    // Deallocate memory
    }
    // Flag 7: Define property of continuous time states (algebraic or differential states)
    else if( flag == 7 )
    {
        for( i = 0; i < block->nx; i++ )
            block->xprop[ i ] = 1;
    }
    // Flag 9: Zero crossing computation
    else if( flag == 9 )
    {
        // Get event indicators. If sign changes for one of block->g vector values, block->jroot
        // and block->nevptr will be automatically set before StateUpdate job
        for( i = 0; i < block->ng; i++ )
            block->g[ i ] = (double) eventIndicatorsList[ i ];
    }
    // Flag 10: Jacobian computation
    else if( flag == 10 )
    {
        fmi2Real dx = 1.0;      // State variation for discrete derivatives calculation
        int xd_mat_size = block->nx * block->nx;
      
        for( i = 0; i < block->nx; i++ )
        {
            fmi2GetDirectionalDerivative( model,
                                          STATE_REFS_LIST + i, 1,          // array and number of derivated state references
                                          STATE_DER_REFS_LIST, block->nx,  // array and number of derivative references
                                          &dx, stateDerivativesList );     // state deltas and resulting state derivatives
            
            // Fill matrix column top with acquired values
            for( i = 0; i < block->nx; i++ )
                block->res[ i * block->nx ] = (double) stateDerivativesList[ i ];
            
            fmi2GetDirectionalDerivative( model,
                                          STATE_REFS_LIST + i, 1,         // array and number of derivated state references
                                          OUTPUT_REFS_LIST, block->nout,  // array and number of output references
                                          &dx, outputsList );             // state deltas and resulting output derivatives
            
            // Fill matrix column bottom with acquired values
            for( i = 0; i < block->nx; i++ )
                block->res[ xd_mat_size + i * block->nx ] = (double) outputsList[ i ];
        }
        
        set_block_error( 0 );
    }
    
    return;
}

Final Thoughts

Even if now it is possible to get the simulation to output something usable, segfaults still happen rather frequently, and the results are not exactly the same as with modelicac (which is more visible with other Modelica Xcos demos). From now until the end of GSoC Working Period, my main task will involve using debugging tools (like running Scilab with -debug command-line option) to investigate remaining crashes and trying to match outputs from OpenModelica models to modelicac ones as much as possible.

That’s it for now. Thanks one more time for sticking by. See you soon !