I want to solve a non-linear problem with non-linear equality constrains and I'm using a augmented Lagrangian with a penalty regularization term that, as well known, spoils the condition number of my linearized systems (at each Newton iteration I mean). The bigger the penalty term, the worse the condition number is. Would someone know an efficient way to get rid of this bad conditioning in that specific case ?

To be more specific, I'm using the classical augmented lagrangian because I have lots of constraints which may generally be redundant. So blindly incorporating the constraints direclty into the primal variables is very convenient. I tried other more sophisticated approaches based on variable eliminations or efficient preconditioners directly on the KKT system but, because of constraints redundancy, I have some troubles.

The problem with regard to $\mathbf u =[u_1,\cdots,u_n]$ variables is formulated as follow my Lagrangian as the form $$\mathcal L(\mathbf u,\lambda):= \mathcal W(\mathbf u) + \rho \lambda^T \,c(\mathbf u) + \frac{\rho}{2} c^2(\mathbf u) $$

So generally The goal at each Newton iteration is to solve a problem of the form $$A \Delta u = b$$ With (we drop hessian of the constraint) $$A(\mathbf u,\rho): = \nabla_{\mathbf u}^2 \mathcal W(\mathbf u) + \rho C^T(\mathbf u)C(\mathbf u) $$ and $$b(\mathbf u,\rho) :=- \big(\nabla_{\mathbf u}\mathcal W(\mathbf u) +(\rho +\lambda^Tc(\mathbf u)) \nabla_{\mathbf u}(\mathbf u)\big)$$ and the capital $C$ is meant for $C(\mathbf u) := \nabla_{\mathbf u} c(\mathbf u)$.

Thank you.

  • $\begingroup$ Hi Tom. Welcome to Scicomp. In order to help us answer your question, could you write the equations that you're trying to solve? $\endgroup$
    – Paul
    Oct 4 '12 at 13:46
  • $\begingroup$ Do you mean $A\Delta u=b$? $\endgroup$ Oct 4 '12 at 15:13
  • $\begingroup$ oops sorry. Yep, sure. $\endgroup$
    – Tom
    Oct 4 '12 at 15:15

Depending on the problem structure, you can solve the ill-conditioned Augmented Lagrangian system directly. For example, BDDC/FETI-DP can solve almost-incompressible elasticity in primal form with a convergence rate independent of the Poisson ratio (piecewise constant on subdomains, but with arbitrary jumps). Similarly, multigrid methods that exactly reproduce the volumetric mode can have this property. Such methods are problem-specific and in general, large penalties result in systems that are difficult to precondition.

To allow more flexibility in preconditioner choice, I recommend introducing explicit dual variables and writing the larger saddle point system

$$ \begin{pmatrix}A & C^T \\ C & -\rho^{-1} \end{pmatrix} \begin{pmatrix} x \\ y \end{pmatrix} = \begin{pmatrix} b \\ 0 \end{pmatrix}$$

as suggested by Arnold Neumaier. This system is much better conditioned and allows you to accurately evaluate a residual. If a preconditioner for some penalized system $A - \tilde\rho C^T C$ exists (where $\tilde \rho \le \rho$), you can use it as a block preconditioner for the saddle point system. For an example of this, see Dohrmann and Lehoucq (2006) which preconditions incompressible elasticity in mixed form using BDDC applied to compressible problems. Another popular class of methods are based on approximating the Schur complement $-\rho^{-1} - C A^{-1} C^T$ using "approximate commutator" arguments. There is an extremely diverse range of methods for solving saddle point problems, see Benzi, Golub, and Liesen, Numerical Solution of Saddle Point Problems (2005) for a review. If you are using PETSc, many of the methods described in the review above can be constructed using run-time options via the PCFIELDSPLIT component.

If you can be more specific about the source of your problem (what are you minimizing and what is the constraint), I may be able to suggest more specific references.

  • $\begingroup$ preconditioners for regularized system open some new ways for me ! However I will need some time to digest all that, I might come back to you after a while if you don't mind.Thank you very much both of you for your answers. $\endgroup$
    – Tom
    Oct 5 '12 at 16:24

Introduce extra variables for the spoiling terms in the KT condition, and you can find a bigger symmetric system that is numerically well-behaved, with only the inverse of the penalty factor entering the matrix.

To solve the ill-conditioned system $(A+\rho C^TC)x=b~$ when $\rho$ is large, introduce $y=\rho Cx$ and recast your problem in the form $Ax+C^Ty=b$, $Cx-\rho^{-1}y=0$, which is generically well-conditioned.

  • $\begingroup$ Generally speaking my constraints are of the form $c(\mathbf u) = 0$ where $\mathbf u$ involves generally just very few degrees of freedom. For instance we may have some projection constraints such as $c(x_s,x_1,x_2)= (x_2-x_1) \mathbf n$ corresponding to the projection of the 3D point $x_s$ on the segment $\[x_1,x_2\]$. $\endgroup$
    – Tom
    Oct 4 '12 at 14:42
  • $\begingroup$ @Tom: I didn't mean the nonlinear problem but the ill-conditioned equations you end up with. Please write down (by editing your question) the form of the linear system you want to solve, and how the penalty parameter enters. $\endgroup$ Oct 4 '12 at 14:44
  • $\begingroup$ I'm trying to figure out how introducing extra-variables would do the trick... Could you send me a reference ? Thank you very much! $\endgroup$
    – Tom
    Oct 4 '12 at 15:11
  • $\begingroup$ @Tom: see the edited answer. $\endgroup$ Oct 4 '12 at 15:16
  • $\begingroup$ But if $\rho$ is large then $\rho^{-1} \approx 0$ and the system is very similar to the original KKT system which is indefinite, isn't it ? $\endgroup$
    – Tom
    Oct 4 '12 at 15:49

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.