Parameter Sensitivities with OpenModelica

This section describes the use of OpenModelica to compute parameter sensitivities using forward sensitivity analysis together with the Sundials/IDA solver.

Note: this is a very short preliminary description which soon will be considerably improved, since this a rather new feature and will continuous improved.

Note: OpenModelica version 1.10 or newer is required.

Background

Parameter sensitivity analysis aims at analyzing the behavior of the corresponding model states w.r.t. model parameters.

Formally, consider a Modelica model as a DAE system:

F(x, \dot x, y, p, t) = 0 \; x(t_0) = x_0(p)

where x(t) \in \mathbf{R}^n represent state variables, \dot x(t) \in \mathbf{R}^n represent state derivatives, y(t) \in \mathbf{R}^k represent algebraic variables, p \in \mathbf{R}^m model parameters.

For parameter sensitivity analysis the derivatives

\frac{\partial x}{ \partial p}

are required which quantify, according to their mathematical definition, the impact of parameters p on states x. In the Sundials/IDA implementation the derivatives are used to evolve the solution over the time by:

\dot s_i = \frac{\partial x}{ \partial p_i}

An Example

This section demonstrates the usage of the sensitivities analysis in OpenModelica on an example. This module is enabled by the following OpenModelica compiler flag:

>>> setCommandLineOptions("--calculateSensitivities");
Listing 4 LotkaVolterra.mo
model LotkaVolterra
  Real x(start=5, fixed=true),y(start=3, fixed=true);
  parameter Real mu1=5,mu2=2;
  parameter Real lambda1=3,lambda2=1;
equation
  0 = x*(mu1-lambda1*y) - der(x);
  0 = -y* (mu2 -lambda2*x) - der(y);
end LotkaVolterra;

Also for the simulation it is needed to set IDA as solver integration method and add a further simulation flag -idaSensitivity to calculate the parameter sensitivities during the normal simulation.

>>> simulate(LotkaVolterra, method="ida", simflags="-idaSensitivity")
record SimulationResult
    resultFile = "«DOCHOME»/LotkaVolterra_res.mat",
    simulationOptions = "startTime = 0.0, stopTime = 1.0, numberOfIntervals = 500, tolerance = 1e-06, method = 'ida', fileNamePrefix = 'LotkaVolterra', options = '', outputFormat = 'mat', variableFilter = '.*', cflags = '', simflags = '-idaSensitivity'",
    messages = "stdout            | info    | The initialization finished successfully without homotopy method.
stdout            | info    | The simulation finished successfully.
",
    timeFrontend = 0.006732210000000001,
    timeBackend = 0.003632196,
    timeSimCode = 0.05869378400000001,
    timeTemplates = 0.03822418900000001,
    timeCompile = 0.41016682,
    timeSimulation = 0.014152228,
    timeTotal = 0.531881285
end SimulationResult;

Now all calculated sensitivities are stored into the results mat file under the $Sensitivities block, where all currently every top-level parameter of the Real type is used to calculate the sensitivities w.r.t. every state.

_images/LotkaVolterraSensitivities.svg

Figure 67 Results of the sensitivities calculated by IDA solver.

_images/LotkaVolterraResults.svg

Figure 68 Results of the LotkaVolterra equations.