Documentation |
The goal in this problem is to minimize the nonlinear function
$$f(x)=1+{\displaystyle \sum _{i=1}^{n}{\left|\left(3-2{x}_{i}\right){x}_{i}-{x}_{i-1}-{x}_{i+1}+1\right|}^{p}}+{\displaystyle \sum _{i=1}^{n/2}{\left|{x}_{i}+{x}_{i+n/2}\right|}^{p}},$$
such that -10.0 ≤ x_{i} ≤ 10.0, where n is 800 (n should be a multiple of 4), p = 7/3, and x_{0} = x_{n + 1} = 0.
The tbroyfg.m file computes the function value and gradient. This file is long and is not included here. You can see the code for this function using the command
type tbroyfg
The sparsity pattern of the Hessian matrix has been predetermined and stored in the file tbroyhstr.mat. The sparsity structure for the Hessian of this problem is banded, as you can see in the following spy plot.
load tbroyhstr spy(Hstr)
In this plot, the center stripe is itself a five-banded matrix. The following plot shows the matrix more clearly:
spy(Hstr(1:20,1:20))
Use optimoptions to set the HessPattern parameter to Hstr. When a problem as large as this has obvious sparsity structure, not setting the HessPattern parameter requires a huge amount of unnecessary memory and computation. This is because fmincon attempts to use finite differencing on a full Hessian matrix of 640,000 nonzero entries.
You must also set the GradObj parameter to 'on' using optimoptions, since the gradient is computed in tbroyfg.m. Then execute fmincon as shown in Step 2.
fun = @tbroyfg; load tbroyhstr % Get Hstr, structure of the Hessian n = 800; xstart = -ones(n,1); xstart(2:2:n) = 1; lb = -10*ones(n,1); ub = -lb; options = optimoptions('fmincon','GradObj','on','HessPattern',Hstr,... 'Algorithm','trust-region-reflective'); [x,fval,exitflag,output] = ... fmincon(fun,xstart,[],[],[],[],lb,ub,[],options);
After seven iterations, the exitflag, fval, and output values are
exitflag = 3 fval = 270.4790 output = iterations: 7 funcCount: 8 cgiterations: 18 firstorderopt: 0.0163 algorithm: 'trust-region-reflective' message: 'Local minimum possible. fmincon stopped because the final chan...' constrviolation: 0
For bound constrained problems, the first-order optimality is the infinity norm of v.*g, where v is defined as in Box Constraints, and g is the gradient.
Because of the five-banded center stripe, you can improve the solution by using a five-banded preconditioner instead of the default diagonal preconditioner. Using the optimoptions function, reset the PrecondBandWidth parameter to 2 and solve the problem again. (The bandwidth is the number of upper (or lower) diagonals, not counting the main diagonal.)
fun = @tbroyfg; load tbroyhstr % Get Hstr, structure of the Hessian n = 800; xstart = -ones(n,1); xstart(2:2:n,1) = 1; lb = -10*ones(n,1); ub = -lb; options = optimoptions('fmincon','GradObj','on','HessPattern',Hstr, ... 'Algorithm','trust-region-reflective','PrecondBandWidth',2); [x,fval,exitflag,output] = ... fmincon(fun,xstart,[],[],[],[],lb,ub,[],options);
The number of iterations actually goes up by two; however the total number of CG iterations drops from 18 to 15. The first-order optimality measure is reduced by a factor of 1e-3:
exitflag = 3 fval = 270.4790 output = iterations: 9 funcCount: 10 cgiterations: 15 firstorderopt: 7.5340e-05 algorithm: 'trust-region-reflective' message: 'Local minimum possible. fmincon stopped because the final chan...' constrviolation: 0