Start by importing the necessary modules,

        import numpy as N 
        import pylab as P #Used for plotting
        from pyfmi import FMUModel #The FMI Interface

Next, the FMU is to be loaded and initialized,

        #Load the FMU by specifying the fmu together with the path.
        bouncing_fmu = FMUModel('/path/to/FMU/bouncingBall.fmu')

        Tstart = 0.5 #The start time.
        Tend   = 3.0 #The final simulation time.
        
        bouncing_fmu.time = Tstart #Set the start time before the initialization.
                                   #(Defaults to 0.0)
        
        bouncing_fmu.initialize() #Initialize the model. Also sets all the start 
                                  #attributes defined in the XML file.

The first line loads the FMU and connects the C-functions of the model to Python together with loading the information from the XML-file. The start time also needs to be specified by setting the property time. The model is also initialized, which must be done before the simulation is started.

Note that if the start time is not specified, FMUModel tries to find the starting time in the XML-file structure 'default experiment' and if successful starts the simulation from that time. Also if the XML-file does not contain any information about the default experiment the simulation is started from time zero.

Then information about the first step is retrieved and stored for later use.

        #Get Continuous States
        x = bouncing_fmu.continuous_states
        #Get the Nominal Values
        x_nominal = bouncing_fmu.nominal_continuous_states
        #Get the Event Indicators
        event_ind = bouncing_fmu.get_event_indicators()
        
        #Values for the solution
        vref  = [bouncing_fmu.get_valueref('h')] + \
                [bouncing_fmu.get_valueref('v')] #Retrieve the valureferences for the
                                                 #values 'h' and 'v't_sol = [Tstart]
        sol = [bouncing_fmu.get_real(vref)]

Here the continuous states together with the nominal values and the event indicators are stored to be used in the integration loop. In our case the nominal values are all equal to one. This information is available in the XML-file. We also create lists which are used for storing the result. The final step before the integration is started is to define the step-size.

        time = Tstart
        Tnext = Tend #Used for time events
        dt = 0.01 #Step-size

We are now ready to create our main integration loop where the solution is advanced using the explicit Euler method.

        #Main integration loop.
        while time < Tend and not bouncing_fmu.get_event_info().terminateSimulation:
            #Compute the derivative of the previous step f(x(n), t(n))
            dx = bouncing_fmu.get_derivatives()
            
            #Advance
            h = min(dt, Tnext-time)
            time = time + h
            
            #Set the time
            bouncing_fmu.time = time
            
            #Set the inputs at the current time (if any)
            #bouncing_fmu.set_real,set_integer,set_boolean,set_string (valueref, values)
            
            #Set the states at t = time (Perform the step using x(n+1)=x(n)+hf(x(n), t(n))
            x = x + h*dx 
            bouncing_fmu.continuous_states = x

This is the integration loop for advancing the solution one step. The loop continues until the final time have been reached or if the FMU reported that the simulation is to be terminated. At the start of the loop the derivatives of the continuous states are retrieved and then the simulation time is incremented by the step-size and set to the model. It could also be the case that the model is depended on inputs which can be set using the set_(real/...) methods.

Note that only variables defined in the XML-file to be inputs can be set using the set_(real/...) methods according to the FMI specification.

The step is performed by calculating the new states (x+h*dx) and setting the values into the model. As our model, the bouncing ball also consist of event functions which needs to be monitored during the simulation, we have to check the indicators which is done below.

            #Get the event indicators at t = time
            event_ind_new = bouncing_fmu.get_event_indicators()
            
            #Inform the model about an accepted step and check for step events
            step_event = bouncing_fmu.completed_integrator_step()
            
            #Check for time and state events
            time_event  = abs(time-Tnext) <= 1.e-10
            state_event = True if True in ((event_ind_new>0.0) != (event_ind>0.0))\
                          else False

Events can be, time, state or step events. The time events are checked by continuously monitor the current time and the next time event (Tnext). State events are checked against sign changes of the event functions. Step events are monitored in the FMU, in the method completed_integrator_step() and return True if any event handling is necessary. If an event have occurred, it needs to be handled, see below.

            #Event handling
            if step_event or time_event or state_event:
                
                eInfo = bouncing_fmu.get_event_info()
                eInfo.iterationConverged = False
                
                #Event iteration
                while eInfo.iterationConverged == False:
                    bouncing_fmu.event_update('0') #Stops at each event iteration
                    eInfo = bouncing_fmu.get_event_info()

                    #Retrieve solutions (if needed)
                    if eInfo.iterationConverged == False:
                        #bouncing_fmu.get_real,get_integer,get_boolean,get_string(valueref)
                        pass
                
                #Check if the event affected the state values and if so sets them
                if eInfo.stateValuesChanged:
                    x = bouncing_fmu.continuous_states
            
                #Get new nominal values.
                if eInfo.stateValueReferencesChanged:
                    atol = 0.01*rtol*bouncing_fmu.nominal_continuous_states
                    
                #Check for new time event
                if eInfo.upcomingTimeEvent:
                    Tnext = min(eInfo.nextEventTime, Tend)
                else:
                    Tnext = Tend

If an event occurred, we enter the iteration loop where we loop until the solution of the new states have converged. During this iteration we can also retrieve the intermediate values with the normal get methods. At this point eInfo contains information about the changes made in the iteration. If the state values have changed, they are retrieved. If the state references have changed, meaning that the state variables no longer have the same meaning as before by pointing to another set of continuous variables in the model, for example in the case with dynamic state selection, new absolute tolerances are calculated with the new nominal values. Finally the model is checked for a new time event.

            event_ind = event_ind_new
        
            #Retrieve solutions at t=time for outputs
            #bouncing_fmu.get_real,get_integer,get_boolean,get_string (valueref)
            
            t_sol += [time]
            sol += [bouncing_fmu.get_real(vref)]

In the end of the loop, the solution is stored and the old event indicators are stored for use in the next loop.

After the loop have finished, by reaching the final time, we plot the simulation results

        #Plot the height
        P.figure(1)
        P.plot(t_sol,N.array(sol)[:,0])
        P.title(bouncing_fmu.get_name())
        P.ylabel('Height (m)')
        P.xlabel('Time (s)')
        #Plot the velocity
        P.figure(2)
        P.plot(t_sol,N.array(sol)[:,1])
        P.title(bouncing_fmu.get_name())
        P.ylabel('Velocity (m/s)')
        P.xlabel('Time (s)')
        P.show()

and the figure below shows the results.

Figure 5.1. Simulation result