The Python commands in the following example may be copied and pasted directly into a Python shell, in some cases with minor modifications. Alternatively, they may be copied into a text file, which also is the recommended way.
For more examples on how to simulate an FMU using JModelica.org's high-level features, see Chapter 7.
This example shows how to use the native JModelica.org FMI interface for simulation of an FMU. The FMU that is to be simulated is the bouncing ball example from Qtronics FMU SDK (http://www.qtronic.de/en/fmusdk.html). This example is written similar to the example in the documentation of the 'Functional Mock-up Interface for Model Exchange' version 1.0 (http://www.functional-mockup-interface.org/). The bouncing ball model is to be simulated using the explicit Euler method with event detection.
The example can also be found in the Python examples catalog in the JModelica.org platform.
The bouncing ball consists of two equations,
and one event function (also commonly called root function),
Where the ball bounces and lose some of its energy according to,
Here, h is the height, g the gravity, v the velocity and e a dimensionless parameter. The starting values are, h=1 and v=0 and for the parameters, e=0.7 and g = 9.81.
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
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
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.