2. A first example

In this section, a simple optimal control problem will be solved. Consider the optimal control problem for the Van der Pol oscillator model:

optimization VDP_Opt (objective = cost(finalTime),
                      startTime = 0,
                      finalTime = 20)

  // The states
  Real x1(start=0,fixed=true);
  Real x2(start=1,fixed=true);

  // The control signal
  input Real u;

  Real cost(start=0,fixed=true);

equation
  der(x1) = (1 - x2^2) * x1 - x2 + u;
  der(x2) = x1;
  der(cost) = x1^2 + x2^2 + u^2;
constraint 
   u<=0.75;
end VDP_Opt;

Create a new file named VDP_Opt.mop and save it in you working directory. Notice that this model contains both the dynamic system to be optimized and the optimization specification. This is possible since Optimica is an extension of Modelica and thereby supports also Modelica constructs such as variable declarations and equations. In most cases, however, Modelica models are stored separately from the Optimica specifications.

Next, create a Python script file and a write (or copy paste) the following commands:

# Import the function for compilation of models and the JMUModel class
from pymodelica import compile_jmu
from pyjmi import JMUModel

# Import the plotting library
import matplotlib.pyplot as plt

Next, we compile and load the model:

# Compile model
jmu_name = compile_jmu("VDP_Opt","VDP_Opt.mop")
# Load model
vdp = JMUModel(jmu_name)

The function compile_jmu invokes the Optimica compiler and compiles the model into a DLL, which is then loaded when the vdp object is created. This object represents the compiled model and is used to invoke the optimization algorithm:

res = vdp.optimize()

In this case, we use the default settings for the optimization algorithm. The result object can now be used to access the optimization result:

# Extract variable profiles
x1=res['x1']
x2=res['x2']
u=res['u']
t=res['time']

The variable trajectories are returned as numpy arrays and can be used for further analysis of the optimization result or for visualization:

plt.figure(1)
plt.clf()
plt.subplot(311)
plt.plot(t,x1)
plt.grid()
plt.ylabel('x1')

plt.subplot(312)
plt.plot(t,x2)
plt.grid()
plt.ylabel('x2')

plt.subplot(313)
plt.plot(t,u)
plt.grid()
plt.ylabel('u')
plt.xlabel('time')
plt.show()

You should now see the optimization result as shown in Figure 7.1.

Figure 7.1. Optimal profiles for the VDP oscillator

Optimal profiles for the VDP oscillator

Optimal control and state profiles for the Van Der Pol optimal control problem.