[petsc-dev] Semismooth Newton

[Todd Munson]

On Aug 1, 2012, at 4:13 PM, Jed Brown wrote:

> On Wed, Aug 1, 2012 at 11:10 AM, Todd Munson <tmunson at mcs.anl.gov> wrote:
> 
> Yes, I am all in favor of actually specifying a VI externally and then
> internally formulating the appropriate augmented system.  This is what
> I do in other contexts (i.e. PATH for GAMS/AMPL).  Some things you
> probably want to do with such an interface is differentiate between
> linear and nonlinear constraints -- you do need second derivatives
> for nonlinear constraints and they should be convex.
> 
> Why not treat it as always nonlinear with the trivial second derivative in the linear case?

You can if you want.  I would say that one benefit of splitting is that sometimes
roundoff errors can be an issue and isolating the linear parts (which are just
data) from the nonlinear parts may help.  It depends on your Newton method
though.  Knowing that your constraints are linear also means that you
do not need to worry about the Hessian of the Lagrangian.  Adding and 
dropping linear constraints is would also seem to be easier since
there is no penalty to be paid since the Jacobian does not need
to be modified.

One interesting thing is that for nonlinear constraints, your Jacobian
matrix in the upper right hand block is the summation of a general matrix
(the Jacobian of F(x)) and a symmetric matrix (sum_i lambda_i nabla^2 c_i(x),
where lambda is the Lagrange multiplier).  Is there a good matrix storage
scheme for this combination?

> You can start getting very general (i.e. linear and nonlinear cone
> constraints, which can appear in contact problems with friction),
> but I would advocate for simplicity rather than generality.
> 
> I like simplicity, but I also like not having to change interfaces. I would rather write an interface for a more general case and even temporarily have an implementation that can't do the fully general thing, eventually perhaps doing something more general. We put a rather large amount of effort into making TS work with a convenient (solver-friendly) representation for DAEs, even though lots of methods cannot be used for DAEs. If the general case is not too complicated for the user, I'd rather do it that way and internally reformulate.

The problem with general cones is that they are not easy (you need both a definition
of the cone and the polar cone) and the numerical methods are problematic (second-order
cones, for example, are not differentiable at the origin).  Here you require 
interior-point methods.  So, the cone constraints algorithms themselves are
more of a research issue.

> The addition and deletion of nonlinear constraints in your example is fine.
> For an active set method, the constraints can be implicitly defined and you
> only need to bring in the new active nonlinear constraints and add them to
> the constraint set for the linear algebra.  You can also do this for other
> methods.  So the active nonlinear constraints can change from one Newton
> iteration to the next.  The only place where the inactive constraints
> matter is in the evaluation of the residual for the line search and
> termination tests.
> 
> Cool, I'm glad you're happy with the math. Now we need a suitable interface, which I think means something a fair amount different from what we have now. We need constraint function, its Jacobian, and a Hessian, each of which can change size?

Yep.  I leave it to you to figure this out though.

> For the contact problems with friction, you typically do an explicit
> method though.  In this setting, you know the subset of the possible
> active constraints at each time step and can keep them constant in
> the VI solver.  They only change as you step forward through time.
> So you can optimize for this case and not have to worry about the
> constraints changing within each Newton iteration.
> 
> I'm not sure how you want to design the variable constraint functionality,
> but there is nothing wrong with it technically.
> 
> I can explain the semismooth method with slacks thing in a few lines.
> The simple problem is:
> 
>   0 <= x  perp  f(x) >= 0
> 
> The compact semismooth formulation is:
> 
>   Phi(x,f(x)) = 0
> 
> where Phi is the reformulation used.
> 
> The semismooth formulation with slacks is:
> 
>   f(x) - s = 0
>   Phi(x,s) = 0

> The linear algebra for the semismooth reformulation with slacks
> is very similar to that for the compact formulation, but there
> are subtle differences in the right hand side and the
> subdifferential.
> 
> Do you mean that the slacks are eliminated again in the linear algebra? What are the "subtle differences"?


If you eliminate the slack variables, in the slack formulation, then you get
a system that looks like the semismooth system, but with a different right
hand side.  Apparently, the difference in the right hand side leads to
fewer Newton iterations in general.  Not a big deal.

Cheers, Todd.