To minimize a large-scale quadratic with upper and lower bounds,
you can use the `quadprog`

function
with the `'trust-region-reflective'`

algorithm.

The problem stored in the MAT-file `qpbox1.mat`

is
a positive definite quadratic, and the Hessian matrix `H`

is
tridiagonal, subject to upper (`ub`

) and lower (`lb`

)
bounds.

load qpbox1 % Get H lb = zeros(400,1); lb(400) = -inf; ub = 0.9*ones(400,1); ub(400) = inf; f = zeros(400,1); f([1 400]) = -2;

xstart = 0.5*ones(400,1); options = optimoptions('quadprog','Algorithm','trust-region-reflective'); [x,fval,exitflag,output] = ... quadprog(H,f,[],[],[],[],lb,ub,xstart,options);

Looking at the resulting values of `exitflag`

and `output`

,

exitflag,output exitflag = 3 output = algorithm: 'trust-region-reflective' iterations: 19 constrviolation: 0 firstorderopt: 8.3903e-06 cgiterations: 1673 message: 'Optimization terminated: relative function value changing by le...'

You can see that while convergence occurred in 20 iterations,
the high number of CG iterations indicates that the cost of the linear
system solve is high. In light of this cost, one strategy would be
to limit the number of CG iterations per optimization iteration. The
default number is the dimension of the problem divided by two, 200
for this problem. Suppose you limit it to 50 using the `MaxPCGIter`

flag
in `options`

:

options = optimoptions(options,'MaxPCGIter',50); [x,fval,exitflag,output] = ... quadprog(H,f,[],[],[],[],lb,ub,xstart,options);

This time convergence still occurs and the total number of CG iterations (1547) has dropped:

exitflag,output exitflag = 3 output = algorithm: 'trust-region-reflective' iterations: 36 constrviolation: 0 firstorderopt: 2.3821e-05 cgiterations: 1547 message: 'Optimization terminated: relative function value changing by le...'

A second strategy would be to use a direct solver at each iteration
by setting the `PrecondBandWidth`

option to `inf`

:

options = optimoptions(options,'PrecondBandWidth',inf); [x,fval,exitflag,output] = ... quadprog(H,f,[],[],[],[],lb,ub,xstart,options);

Now the number of iterations has dropped to 10:

exitflag,output exitflag = 3 output = algorithm: 'trust-region-reflective' iterations: 10 constrviolation: 0 firstorderopt: 2.8219e-06 cgiterations: 0 message: 'Optimization terminated: relative function value changing by le...'

Using a direct solver at each iteration usually causes the number
of iterations to decrease, but often takes more time per iteration.
For this problem, the tradeoff is beneficial, as the time for `quadprog`

to solve the problem decreases
by a factor of 10.

You can also use the default `'interior-point-convex'`

algorithm
to solve this convex problem:

options = optimoptions('quadprog','Algorithm','interior-point-convex'); [x,fval,exitflag,output] = ... quadprog(H,f,[],[],[],[],lb,ub,[],options);

Check the exit flag and output structure:

exitflag,output exitflag = 1 output = message: 'Minimum found that satisfies the constraints. Optimization com...' algorithm: 'interior-point-convex' firstorderopt: 1.4120e-06 constrviolation: 0 iterations: 8 cgiterations: []

Was this topic helpful?