Suppose that you want to optimize the control parameters in
the Simulink^{®} model `optsim.mdl`

.
(This model can be found in the `optim/optimdemos`

folder.
Note that Simulink must be installed on your system to load this
model.) The model includes a nonlinear process plant modeled as a Simulink block
diagram.

**Plant with Actuator Saturation**

The plant is an under-damped third-order model with actuator
limits. The actuator limits are a saturation limit and a slew rate
limit. The actuator saturation limit cuts off input values greater
than 2 units or less than –2 units. The slew rate limit of
the actuator is 0.8 units/sec. The closed-loop response of the system
to a step input is shown in Closed-Loop Response. You can see this response by opening
the model (type `optsim`

at
the command line or click the model name), and selecting **Run** from the **Simulation** menu.
The response plots to the scope.

**Closed-Loop Response**

The problem is to design a feedback control loop that tracks a unit step input to the system. The closed-loop plant is entered in terms of the blocks where the plant and actuator have been placed in a hierarchical Subsystem block. A Scope block displays output trajectories during the design process.

**Closed-Loop Model**

One way to solve this problem is to minimize the error between the output and the input signal. The variables are the parameters of the Proportional Integral Derivative (PID) controller. If you only need to minimize the error at one time unit, it would be a single objective function. But the goal is to minimize the error for all time steps from 0 to 100, thus producing a multiobjective function (one function for each time step).

The routine `lsqnonlin`

is used to perform a least-squares
fit on the tracking of the output. The tracking is performed via the
function `tracklsq`

, which returns the error signal `yout`

,
the output computed by calling `sim`

,
minus the input signal `1`

. The code for `tracklsq`

is
contained in the file `runtracklsq.m`

, shown below.

The function `runtracklsq`

sets up all the
needed values and then calls `lsqnonlin`

with the
objective function `tracklsq`

, which is nested inside `runtracklsq`

.
The variable `options`

passed to `lsqnonlin`

defines
the criteria and display characteristics. In this case you ask for
output, use the `'levenberg-marquardt'`

algorithm,
and give termination tolerances for the step and objective function
on the order of `0.001`

.

To run the simulation in the model `optsim`

,
the variables `Kp`

, `Ki`

, `Kd`

, `a1`

,
and `a2`

(`a1`

and `a2`

are
variables in the Plant block) must all be defined. `Kp`

, `Ki`

,
and `Kd`

are the variables to be optimized. The function `tracklsq`

is
nested inside `runtracklsq`

so that the variables `a1`

and `a2`

are
shared between the two functions. The variables `a1`

and `a2`

are
initialized in `runtracklsq`

.

The objective function `tracklsq`

runs the
simulation. The simulation can be run either in the base workspace
or the current workspace, that is, the workspace of the function calling `sim`

,
which in this case is the workspace of `tracklsq`

.
In this example, the `SrcWorkspace`

option is set
to `'Current'`

to tell `sim`

to
run the simulation in the current workspace. The simulation is performed
to `100`

seconds.

When the simulation is completed, the `myobj`

object
is created in the current workspace (that is, the workspace of `tracklsq`

).
The Outport block in the block diagram model puts the `yout`

field
of the object into the current workspace at the end of the simulation.

The following is the code for `runtracklsq`

:

function [Kp,Ki,Kd] = runtracklsq % RUNTRACKLSQ demonstrates using LSQNONLIN with Simulink. optsim % Load the model pid0 = [0.63 0.0504 1.9688]; % Set initial values a1 = 3; a2 = 43; % Initialize model plant variables options = optimoptions(@lsqnonlin,'Algorithm','levenberg-marquardt',... 'Display','off','StepTolerance',0.001,'OptimalityTolerance',0.001); pid = lsqnonlin(@tracklsq, pid0, [], [], options); Kp = pid(1); Ki = pid(2); Kd = pid(3); function F = tracklsq(pid) % Track the output of optsim to a signal of 1 % Variables a1 and a2 are needed by the model optsim. % They are shared with RUNTRACKLSQ so do not need to be % redefined here. Kp = pid(1); Ki = pid(2); Kd = pid(3); % Set sim options and compute function value myobj = sim('optsim','SrcWorkspace','Current', ... 'StopTime','100'); F = myobj.get('yout') - 1; end end

Copy the code for `runtracklsq`

to a file named `runtracklsq.m`

,
placed in a folder on your MATLAB^{®} path.

When you run `runtracklsq`

, the optimization
gives the solution for the proportional, integral, and derivative
(`Kp`

, `Ki`

, `Kd`

)
gains of the controller:

[Kp, Ki, Kd] = runtracklsq Kp = 3.1330 Ki = 0.1465 Kd = 14.3918

Here is the resulting closed-loop step response.

**Closed-Loop Response Using lsqnonlin**

The call to `sim`

results
in a call to one of the Simulink ordinary differential equation
(ODE) solvers. A choice must be made about the type of solver to use.
From the optimization point of view, a fixed-step solver is the best choice
if that is sufficient to solve the ODE. However, in the case of a
stiff system, a variable-step method might be required to solve the ODE.

The numerical solution produced by a variable-step solver, however, is not a smooth function of parameters, because of step-size control mechanisms. This lack of smoothness can prevent the optimization routine from converging. The lack of smoothness is not introduced when a fixed-step solver is used. (For a further explanation, see [53].)

Simulink Design Optimization™ software is recommended for solving multiobjective optimization problems in conjunction with Simulink variable-step solvers. It provides a special numeric gradient computation that works with Simulink and avoids introducing a problem of lack of smoothness.