August 2012 - petsc-dev - mailman3.cels.anl.gov

Re: [petsc-dev] Semismooth Newton
by Dmitry Karpeev 02 Aug '12

02 Aug '12

On Wed, Aug 1, 2012 at 5:01 PM, Todd Munson <tmunson(a)mcs.anl.gov> wrote: > > On Aug 1, 2012, at 4:13 PM, Jed Brown wrote: > > > On Wed, Aug 1, 2012 at 11:10 AM, Todd Munson <tmunson(a)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? > How does one precondition a system like that? Presumably, the user may at best be able to provide a preconditioner for JF. If JF responds well to multigrid, then it may be possible to use AMG-KKT to coarsen the constraint Hessian, but how does MG perform for the sum JF + Hc? There is also the problem of smoothing, although, if JF and Jc are sparse, then DD may be an option with a direct solve on subdomains. Anything more definitive on preconditioning these systems? Dmitry. > > > 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. > >

1 0

Re: [petsc-dev] Semismooth Newton
by Todd Munson 01 Aug '12

01 Aug '12

On Aug 1, 2012, at 4:13 PM, Jed Brown wrote: > On Wed, Aug 1, 2012 at 11:10 AM, Todd Munson <tmunson(a)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.

1 0

Re: [petsc-dev] Semismooth Newton
by Jed Brown 01 Aug '12

01 Aug '12

On Wed, Aug 1, 2012 at 11:10 AM, Todd Munson <tmunson(a)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 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 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? > > 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"?

1 0

Re: [petsc-dev] Semismooth Newton
by Todd Munson 01 Aug '12

01 Aug '12

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. 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. 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. 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. Cheers, Todd.

1 0

Re: [petsc-dev] Semismooth Newton
by Jed Brown 01 Aug '12

01 Aug '12

On Tue, Jul 31, 2012 at 5:38 AM, Todd Munson <tmunson(a)mcs.anl.gov> wrote: > Knowing the constraint structure and forming the augmented system is fine > and can be > helpful information to know when determining the active set. I would > concentrate on > linear constraints and inequalities though; they are nicer to deal with > and you do > not have to be concerned with constraint qualifications. > > Nonlinear constraints are possible, but you need to satisfy a constraint > qualification > and you need second derivatives of the constraints for the augmented > system. (In > optimization terms, the Hessian of the Lagrangian.) > I would be in favor of making the interface support this more general form. I'm not wild about having the user add a bunch of extra variables because then they have to explain enough about that reformulation to decompose it again for an efficient multigrid solver. > > > There are some reformulations for polyhedral constraints, but they are, > in my opinion, > > a bit unwieldy. > > > > Why unwieldy? > > The standard semismooth formulation using the Fischer-Burmeister function > and its > cousins works only for bound constraints. Designing a semismooth > reformulation > for general polyhedral constraints is difficult. It can be done, but its > really not easy. You can define the normal map using projections, but > then you need to be able to do the projections onto the constraint > set. > > Of course, the augmented system results in a box constrained problem and > the > standard formulations work. > I'm not at all opposed to having the library use the augmented system internally, I just don't like making the user reformulate because it's more intrusive to their code and because it locks the library into solving the problem that way. > > One thing we should try is adding slacks in the semismooth algorithm for > box > constraints and then applying the reformulation to the redefined system. > There is some evidence that this will work better. > If you have a reference for a version that you like, I'll implement it when I start messing with that code. > > > Note that changing the size of the constrained size is also important, > and box constraints strike me as even more confusing in that context. > > I do not follow this at all. Are you talking about something like a time > dependent > problem and the constraints in the VI are a function of time? That > becomes a > somewhat odd problem and I'd have to see an example. Normally the > constraint set is fixed for all time and just the activities > change. > Suppose we are doing unstructured elastic contact. Any surface point could end up in contact with any other. If we reformulate with bound contstraints, we end up adding N^2 auxiliary variables. If we can use nonlinear constraints where the number of constraints can change each iteration, then we only need to deal with the O(N) viable contact possibilities. I'm happy to differentiate the constraint functions, but I certainly want the flexibility.

1 0

Re: [petsc-dev] Slides on Barcelona "Mont Blanc" project
by Jed Brown 31 Jul '12

31 Jul '12

On Tue, Jul 31, 2012 at 7:17 AM, Matthew Knepley <knepley(a)gmail.com> wrote: > The problem with the programming model they propose, namely unroll the > entire > task graph and schedule, is that it has been proposed at least 20 times > before > by smart, well-funded people. These slides just tell us how the system is > supposed > to work, which is bullshit. They need to tell us how it is different from > other similar > systems which failed, and why that will make it succeed. > I feel like a lot of those projects don't seem to recognize that the granularity we work with is very small. For strong scaling, we need every grid cell to be processed once every 10-100 microseconds. There is a data dependence and ghost cell cost that needs to be amortized over a large enough patch. The task graph isn't a great model if you have only one task per thread. Also, it tends not to follow data locality well.

1 0