Learn more about Optimization Toolbox

Multiobjective Optimization Examples

On this page…
Example: Using fminimax with a Simulink Model Example: Signal Processing Using fgoalattain

Example: Using fminimax with a Simulink Model

Another approach to optimizing the control parameters in the Simulink model shown in Plant with Actuator Saturation is to use the fminimax function. In this case, rather than minimizing the error between the output and the input signal, you minimize the maximum value of the output at any time t between 0 and 100.

The code for this example, shown below, is contained in the function runtrackmm, in which the objective function is simply the output yout returned by the sim command. But minimizing the maximum output at all time steps might force the output to be far below unity for some time steps. To keep the output above 0.95 after the first 20 seconds, the constraint function trackmmcon contains the constraint yout >= 0.95 from t=20 to t=100. Because constraints must be in the form g ≤ 0, the constraint in the function is g = -yout(20:100)+.95.

Both trackmmobj and trackmmcon use the result yout from sim, calculated from the current PID values. The nonlinear constraint function is always called immediately after the objective function in fmincon, fminimax, fgoalattain, and fseminf with the same values. This way you can avoid calling the simulation twice by using nested functions so that the value of yout can be shared between the objective and constraint functions as long as it is initialized in runtrackmm.

The following is the code for runtrackmm:

function [Kp, Ki, Kd] = runtrackmm

optsim
pid0 = [0.63 0.0504 1.9688];
% a1, a2, yout are shared with TRACKMMOBJ and TRACKMMCON
a1 = 3; a2 = 43; % Initialize plant variables in model
yout = []; % Give yout an initial value
options = optimset('Display','iter',...
    'TolX',0.001,'TolFun',0.001);
pid = fminimax(@trackmmobj,pid0,[],[],[],[],[],[],...
    @trackmmcon,options);
Kp = pid(1); Ki = pid(2); Kd = pid(3);

    function F = trackmmobj(pid)
        % Track the output of optsim to a signal of 1.
        % Variables a1 and a2 are shared with RUNTRACKMM.
        % Variable yout is shared with RUNTRACKMM and 
        % RUNTRACKMMCON.
        
        Kp = pid(1);
        Ki = pid(2);
        Kd = pid(3);

        % Compute function value
        opt = simset('solver','ode5','SrcWorkspace','Current');
        [tout,xout,yout] = sim('optsim',[0 100],opt);
        F = yout;
    end

    function [c,ceq] = trackmmcon(pid)
        % Track the output of optsim to a signal of 1.
        % Variable yout is shared with RUNTRACKMM and
        % TRACKMMOBJ

        % Compute constraints.
        % Objective TRACKMMOBJ is called before this
        % constraint function, so yout is current.
        c = -yout(20:100)+.95;
        ceq=[];
    end
end

When you run the code, it returns the following results:

[Kp,Ki,Kd] = runtrackmm

Done initializing optsim.

              Objective        Max     Line search   Directional 
 Iter F-count      value    constraint   steplength   derivative  Procedure
    0      5           0       1.11982                                  
    1     11       1.184       0.07978         1        0.482     
    2     17       1.012       0.04285         1       -0.236     
    3     23      0.9996       0.00397         1      -0.0195  Hessian modified twice  
    4     29      0.9996    3.464e-005         1     0.000687  Hessian modified  
    5     35      0.9996    2.272e-009         1      -0.0175  Hessian modified twice  

Local minimum possible. Constraints satisfied.

fminimax stopped because the predicted change in the objective function
is less than the selected value of the function tolerance and constraints 
were satisfied to within the default value of the constraint tolerance.

Active inequalities (to within options.TolCon = 1e-006):
  lower      upper     ineqlin   ineqnonlin
                                     1
                                    14
                                   182

Kp =
    0.5894

Ki =
    0.0605

Kd =
    5.5295

The last value shown in the MAX{F,constraints} column of the output shows that the maximum value for all the time steps is 0.9996. The closed loop response with this result is shown in the figure Closed-Loop Response Using fminimax.

This solution differs from the solution obtained in Example: Using lsqnonlin With a Simulink Model because you are solving different problem formulations.

Closed-Loop Response Using fminimax

Back to Top

Example: Signal Processing Using fgoalattain

Consider designing a linear-phase Finite Impulse Response (FIR) filter. The problem is to design a lowpass filter with magnitude one at all frequencies between 0 and 0.1 Hz and magnitude zero between 0.15 and 0.5 Hz.

The frequency response H(f) for such a filter is defined by

(6-119)

where A(f) is the magnitude of the frequency response. One solution is to apply a goal attainment method to the magnitude of the frequency response. Given a function that computes the magnitude, the function fgoalattain will attempt to vary the magnitude coefficients a(n) until the magnitude response matches the desired response within some tolerance. The function that computes the magnitude response is given in filtmin.m. This function takes a, the magnitude function coefficients, and w, the discretization of the frequency domain we are interested in.

To set up a goal attainment problem, you must specify the goal and weights for the problem. For frequencies between 0 and 0.1, the goal is one. For frequencies between 0.15 and 0.5, the goal is zero. Frequencies between 0.1 and 0.15 are not specified, so no goals or weights are needed in this range.

This information is stored in the variable goal passed to fgoalattain. The length of goal is the same as the length returned by the function filtmin. So that the goals are equally satisfied, usually weight would be set to abs(goal). However, since some of the goals are zero, the effect of using weight=abs(goal) will force the objectives with weight 0 to be satisfied as hard constraints, and the objectives with weight 1 possibly to be underattained (see Goal Attainment Method). Because all the goals are close in magnitude, using a weight of unity for all goals will give them equal priority. (Using abs(goal) for the weights is more important when the magnitude of goal differs more significantly.) Also, setting

options = optimset('GoalsExactAchieve',length(goal)); 

specifies that each objective should be as near as possible to its goal value (neither greater nor less than).

Step 1: Write an M-file filtmin.m

function y = filtmin(a,w)
n = length(a);
y = cos(w'*(0:n-1)*2*pi)*a ;

Step 2: Invoke optimization routine

% Plot with initial coefficients
a0 = ones(15,1);
incr = 50;
w = linspace(0,0.5,incr);

y0 = filtmin(a0,w);
clf, plot(w,y0,'.r');
drawnow;

% Set up the goal attainment problem
w1 = linspace(0,0.1,incr) ;
w2 = linspace(0.15,0.5,incr);
w0 = [w1 w2];
goal = [1.0*ones(1,length(w1)) zeros(1,length(w2))];
weight = ones(size(goal)); 

% Call fgoalattain
options = optimset('GoalsExactAchieve',length(goal));
[a,fval,attainfactor,exitflag]=fgoalattain(@(x)filtmin(x,w0),...
    a0,goal,weight,[],[],[],[],[],[],[],options);

% Plot with the optimized (final) coefficients
y = filtmin(a,w);
hold on, plot(w,y,'r')
axis([0 0.5 -3 3])
xlabel('Frequency (Hz)')
ylabel('Magnitude Response (dB)')
legend('initial', 'final')
grid on

Compare the magnitude response computed with the initial coefficients and the final coefficients (Magnitude Response with Initial and Final Magnitude Coefficients). Note that you could use the firpm function in Signal Processing Toolbox™ software to design this filter.

Magnitude Response with Initial and Final Magnitude Coefficients

Back to Top

Multiobjective Optimization

Equation Solving

© 1984-2009- The MathWorks, Inc.    -   Site Help   -   Patents   -   Trademarks   -   Privacy Policy   -   Preventing Piracy   -   RSS