Homework 4
Assigned Mon Mar 30, due Thurs Apr 9 (any time)

The Primal Barrier Method
  1. Write a Matlab function to implement the primal barrier method 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 constraints, their gradients and their Hessians: thus, for example, congrad{k} and conhess{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 your code newtmeth: the line search parameters &alpha and &beta and a tolerance telling it to terminate when the Newton decrement is sufficiently small (you will need to add this capability to newtmeth.)

    New: It was pointed out to me by a couple of students that, given the limitations of Matlab's anonymous functions, your existing newtmeth code can't be conveniently used to minimize a barrier function based on funs. The easiest solution seems to be to create a specialized Newton barrier minimization method newtmethbarrier which has parameters that include the barrier parameter t, the anonymous function funs described above, and the starting point x0. Also, you should change the call to chol to return a second 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.

    New: Even if Cholesky reports that the Hessian is positive definite, the full Newton step (initial step 1) 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 which is totally irrelevant. There is an easy fix: the barrier function should be defined as: By simply using the perfectly good floating point number inf in this way the line search will reject the step of 1 and try 1/2, 1/4 etc until it finds a feasible point.
    New: I fixed a typo, I had the wrong sign in front of the sum previously.
  2. For a test problem, minimize -eTx (that is, maximize the sum of the variables) subject to the quadratic constraints xTAk x &le 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.

  3. Generate analogues of Figures 11.4 and 11.5 (not 11.6) (p.572-573) in the text, showing how the method performs for different values mu. (Use the same values 2, 50 and 150 for the first and third figure, but you don't have to run so many cases for the second figure.) Also, you don't need to generate the figures using exactly the same plotting marks; anything that conveys the same information is fine. Set &alpha=0.01, &beta=0.5, the Newton decrement tolerance to 10-5, the initial barrier parameter to 1, and the target duality gap to 10-5. What seems to be the best choice for mu?

    New: Let's change this so that you solve the Newton minimization problem more accurately -- I am not sure why the book recommends the large number that it does, actually 10-2.5, since quadratic convergence means that only a few more steps should be needed to get a more accurate answer. You can experiment with this. Also, there was a question about the meaning of Figure 11.4. I interpret it as follows: the vertical axis is, I suppose, 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 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. Same for Figure 11.5. The previous reference to Figure 11.6 was a mistake. Finally, Figure 11.4 indicates a final target duality gap of 10-6, not 10-5, so let's use that.

  4. How does the actual total number of Newton steps needed compare with the complexity predicted by the self-concordance theory?

  5. How many of the quadratic constraints are active at the solution (nearly hold with equality, holding with equality in the limit as the barrier parameter increases to &infin)? Run 100 other randomly generated examples using the same n and m, and plot a histogram (type "help hist") showing the number of active constraints at the solution. Just use one value of mu for this. Generate the Ak by first generating B=randn(n) and then setting Ak = B'*B.
Turn in Matlab code listings as well as answers to the questions above. Make sure that your codes have plenty of comments explaining what they do. If anything is not clear, don't hesitate to contact me by email or by dropping by my office.