My dear community,

I am wondering why BFGS methods are not so widely used for simulating mechanical problems which heavily still relies on inverting the hessian matrix. I am essentially interested by constrained problems (contact on the interface for instance), and I know that some authors proposed method to deal with Quasi-Newton approximation for Augmented Lagrangian, but it seems to me that computing and inverting the regularized hessian matrix is still the method of choice. My question is: why ? Stability ?

Thank you very much.

Kind regards, Tom


For an unconstrained linear SPD problem, BFGS convergence matches unpreconditioned conjugate gradients. Differential operators appearing in structural mechanics problems are usually ill-conditioned with many large eigenvalues, making unpreconditioned iterative methods hopelessly slow. To deal with the troublesome spectrum, structural mechanics solvers depends on either sparse direct methods or preconditioned iterative methods (often using algebraic multigrid or domain decomposition).

In contrast, BFGS was developed in an optimization setting where most problems have a spectrum that looks like an integral operator under grid refinement, e.g., $I + K$ where $K$ is compact. For linear versions of these problems, Krylov methods would converge in a constant number of iterations, independent of resolution. The Hessian is often dense and practically inaccessible, so it is applied in unassembled form and only the identity or a diagonal scaling can be used as the "seed" Hessian $H_0$. BFGS is thus a way to improve a crude approximation to the inverse Hessian.

When solving differential problems, we can use BFGS as a lagging method, starting with a good approximation that is merely expensive to construct. This is useful when nonlinearity represents a part of the difficulty in solving a problem, but ill-conditioning is significant enough that preconditioning is needed. This turns out to work pretty well for some problems like elasticity.

Peter Brune and I wrote a short paper on this recently:

See Table II on page 9 for a large-deformation elasticity example. The standard Newton method needs 10 Jacobian evaluations, 11 residual evaluations, and 77 preconditioner applications (a linear V-cycle of algebraic multigrid; note that we can use fewer total preconditioner applications using inexact solves, but this will increase the number of Jacobian evaluations). LBFGS-6 with the same V-cycle of AMG used for $H_0$ requires 4 Jacobian evaluations, 49 residual evaluations, and 24 preconditioner applications. LBFGS lets us trade Jacobian work for more residual evaluations while retaining faster convergence than normal Jacobian lagging.

You can try this out in PETSc using only command-line options, e.g.

-snes_type qn -snes_qn_restart_type periodic -snes_qn_scale_type jacobian

The utility of this method for constrained elasticity will depend on how quickly the constraints are changing, but I expect it to be useful. We are planning to extend PETSc's support for solving constrained systems so that these methods can be used for elastic contact problems. Please let us know if you're interested in this.

  • $\begingroup$ Thank you Jed for your detailed answer ! Yes definetly I'm interested into LBFGS applied to constrained problems, eventhough some preconditioning (probably a simple Separate Displacement Component ala Gustaffsson) is necessary. In my specific case the constraints are changing rapidly and drive the deformation process so inverting many linear system is very heavy. $\endgroup$
    – tom
    Apr 23 '13 at 13:10
  • $\begingroup$ Just a quick-one related to the topic. Would it be possible to define coarse quasi-newton update using the same framework than AMG to have MG approach ? My (fuzzy) idea would be to use an initial jacobian to create MG-style coarsenings and use them to define quasi-newton update in a similar way that Natwon-MG works ? $\endgroup$
    – tom
    Apr 24 '13 at 12:53
  • $\begingroup$ If I understand your question, this is what we do in the paper. An algebraic or geometric multigrid cycle is used to define $H_0$ in BFGS. We rebuild that multigrid hierarchy each time L-BFGS is restarted. $\endgroup$
    – Jed Brown
    Apr 24 '13 at 13:03

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.