Quadratic Minimization with Bound Constraints

This example shows the effects of some option settings on a sparse, bound-constrained, positive definite quadratic problem.

Create the quadratic matrix H as a tridiagonal symmetric matrix of size 400-by-400 with entries +4 on the main diagonal and –2 on the off-diagonals.

Bin = -2*ones(399,1);
H = spdiags(Bin,-1,400,400);
H = H + H';
H = H + 4*speye(400);

Set bounds of [0,0.9] in each component except the 400th. Allow the 400th component to be unbounded.

lb = zeros(400,1);
lb(400) = -inf;
ub = 0.9*ones(400,1);
ub(400) = inf;

Set the linear vector f to zeros, except set f(400) = –2.

f = zeros(400,1);
f(400) = -2;

Trust-Region-Reflective Solution

Solve the quadratic program using the 'trust-region-reflective' algorithm.

options = optimoptions('quadprog','Algorithm',"trust-region-reflective");
tic
[x1,fval1,exitflag1,output1] = ... 
        quadprog(H,f,[],[],[],[],lb,ub,[],options);

Local minimum possible.

quadprog stopped because the relative change in function value is less than the function tolerance.

time1 = toc

time1 = 0.1044

Examine the solution.

fval1,exitflag1,output1.iterations,output1.cgiterations

fval1 = -0.9930

exitflag1 = 3

ans = 18

ans = 1682

The algorithm converges in relatively few iterations, but takes over 1000 CG (conjugate gradient) iterations. To avoid the CG iterations, set options to use a direct solver instead.

options = optimoptions(options,'SubproblemAlgorithm','factorization');
tic
[x2,fval2,exitflag2,output2] = ... 
        quadprog(H,f,[],[],[],[],lb,ub,[],options);

Local minimum possible.

quadprog stopped because the relative change in function value is less than the function tolerance.

time2 = toc

time2 = 0.0185

fval2,exitflag2,output2.iterations,output2.cgiterations

fval2 = -0.9930

exitflag2 = 3

ans = 10

ans = 0

This time, the algorithm takes fewer iterations and no CG iterations. The solution time decreases substantially, despite the relatively time-consuming direct factorization steps, because the solver avoids taking many CG steps.

Interior-Point Solution

The default 'interior-point-convex' algorithm can solve this problem.

tic
[x3,fval3,exitflag3,output3] = ... 
    quadprog(H,f,[],[],[],[],lb,ub); % No options means use the default algorithm

Minimum found that satisfies the constraints.

Optimization completed because the objective function is non-decreasing in 
feasible directions, to within the value of the optimality tolerance,
and constraints are satisfied to within the value of the constraint tolerance.

<stopping criteria details>

time3 = toc

time3 = 0.0402

fval3,exitflag3,output3.iterations

fval3 = -0.9930

exitflag3 = 1

ans = 8

Compare Results

All algorithms give the same objective function value to display precision, –0.9930.

The 'interior-point-convex' algorithm takes the fewest iterations. However, the 'trust-region-reflective' algorithm with the direct subproblem solver reaches the solution fastest.

tt = table([time1;time2;time3],[output1.iterations;output2.iterations;output3.iterations],...
    'VariableNames',["Time" "Iterations"],'RowNames',["TRR" "TRR Direct" "IP"])

tt=3×2 table
                    Time      Iterations
                  ________    __________

    TRR            0.10443        18    
    TRR Direct    0.018544        10    
    IP            0.040204         8