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: