Homework 8
Assigned Mon Apr 11, due Thurs Apr 21

Recommended seminar this week: Fri Apr 15 at 10 am by Jim Burke, University of Washington.

The Primal Barrier Method

Write a Matlab function to implement the primal barrier method (Algorithm 11.1 of BV) for minimizing a convex function subject to convex inequality constraints, that is (11.1) with no equality constraints. The first parameter should be an anonymous function funs which, when it is called, will provide the name of a routine that computes the objective function, the inequality constraints, and their gradients and Hessians. Use cell arrays to represent the constraint gradients and Hessians: thus, for example, con_grad{k} and con_hess{k} could return the kth constraint gradient and Hessian respectively. Other parameters needed by the barrier method are the initial value for the barrier parameter, the scalar mu (μ), the initial guess x0, the target for the final duality gap, a maximum number of iterations, and several parameters to pass to Newton's method: the line search parameters α and β and a tolerance telling it to terminate when either the gradient or the "Newton decrement" (BV, p.486, which is simply the square root of minus the inner product of Newton direction with the gradient) is sufficiently small. Instead of the barrier method calling your existing newtmeth directly, it may be easier to write a specialized Newton barrier minimization method newtMethBarrier which has parameters that include the barrier parameter t, the anonymous function funs described above, the starting point x0, and the other parameters described above. You should solve the Newton equation using chol (Cholesky factorization), requesting a second output argument that is nonzero if the matrix is not numerically positive definite, so you can print an informative message and terminate instead of incurring an error if this happens - which it may when the barrier parameter is sufficiently large.

Note that the full Newton step (using stepsize one) may be infeasible. If you just go ahead and evaluate the log barrier function at an infeasible point, you will get the log of a negative number which gives a complex number that is totally irrelevant. There is an easy fix: the barrier function should be defined as:

• B = tf0(x) - sum log -fi(x) if x is feasible
• B = inf (∞) otherwise
Then the line search will reject the step of 1 and try 1/2, 1/4 etc until it finds a feasible point.

For a test problem, minimize -eTx (that is, maximize the sum of the variables) subject to the quadratic constraints xTAk x ≤ 1, for k=1,...,m. Use n=10 and m=20. Notice that zero is an initial feasible point. Use the Ak stored in the cell array in quadratics.mat.

Generate analogues of Figures 11.4 and 11.5 (p.572-573) in BV showing how the method performs for different values mu. The vertical axis is the duality gap m/t and the horizontal lines indicate the number of Newton iterations needed for that particular choice of the barrier parameter. Thus, in this case, the best of these 3 choices of mu is 50. You are supposed to generate an analogous figure for the quadratically constrained problem described above, not for the LP solved by the book, so the figure might end up looking quite different. Choose values of mu that are suitable for this problem (so that the best one is neither the smallest nor the largest value tried). For the analogue of Fig. 11.5, you don't have to run so many cases as shown in the book. Also, you don't need to generate the figures using exactly the same plotting marks; anything that conveys the same information is fine. Set α=0.25, β=0.5, the initial barrier parameter to 1, and the target duality gap to 10-6. Does the best choice for mu depend strongly on the termination tolerance you use for Newton's method?

Include Matlab code listings in your submission. Make sure that your codes have plenty of comments explaining what they do.