Uzawa iteration

Last updated December 18, 2020

In numerical mathematics, the Uzawa iteration is an algorithm for solving saddle point problems. It is named after Hirofumi Uzawa and was originally introduced in the context of concave programming.^[1]

Basic idea

We consider a saddle point problem of the form

{\begin{pmatrix}A&B\\B^{*}&\end{pmatrix}}{\begin{pmatrix}x_{1}\\x_{2}\end{pmatrix}}={\begin{pmatrix}b_{1}\\b_{2}\end{pmatrix}},

where $A$ is a symmetric positive-definite matrix. Multiplying the first row by $B^{*}A^{-1}$ and subtracting from the second row yields the upper-triangular system

{\begin{pmatrix}A&B\\&-S\end{pmatrix}}{\begin{pmatrix}x_{1}\\x_{2}\end{pmatrix}}={\begin{pmatrix}b_{1}\\b_{2}-B^{*}A^{-1}b_{1}\end{pmatrix}},

where $S:=B^{*}A^{-1}B$ denotes the Schur complement. Since $S$ is symmetric positive-definite, we can apply standard iterative methods like the gradient descent method or the conjugate gradient method to

Sx_{2}=B^{*}A^{-1}b_{1}-b_{2}

in order to compute $x_{2}$ . The vector $x_{1}$ can be reconstructed by solving

Ax_{1}=b_{1}-Bx_{2}.\,

It is possible to update $x_{1}$ alongside $x_{2}$ during the iteration for the Schur complement system and thus obtain an efficient algorithm.

Implementation

We start the conjugate gradient iteration by computing the residual

r_{2}:=B^{*}A^{-1}b_{1}-b_{2}-Sx_{2}=B^{*}A^{-1}(b_{1}-Bx_{2})-b_{2}=B^{*}x_{1}-b_{2},

of the Schur complement system, where

x_{1}:=A^{-1}(b_{1}-Bx_{2})

denotes the upper half of the solution vector matching the initial guess $x_{2}$ for its lower half. We complete the initialization by choosing the first search direction

p_{2}:=r_{2}.\,

In each step, we compute

a_{2}:=Sp_{2}=B^{*}A^{-1}Bp_{2}=B^{*}p_{1}

and keep the intermediate result

p_{1}:=A^{-1}Bp_{2}

for later. The scaling factor is given by

{\displaystyle \alpha

and leads to the updates

x_{2}:=x_{2}+\alpha p_{2},\quad r_{2}:=r_{2}-\alpha a_{2}.

Using the intermediate result $p_{1}$ saved earlier, we can also update the upper part of the solution vector

x_{1}:=x_{1}-\alpha p_{1}.\,

Now we only have to construct the new search direction by the Gram–Schmidt process, i.e.,

{\displaystyle \beta

The iteration terminates if the residual $r_{2}$ has become sufficiently small or if the norm of $p_{2}$ is significantly smaller than $r_{2}$ indicating that the Krylov subspace has been almost exhausted.

Modifications and extensions

If solving the linear system $Ax=b$ exactly is not feasible, inexact solvers can be applied.^[2]^[3]^[4]

If the Schur complement system is ill-conditioned, preconditioners can be employed to improve the speed of convergence of the underlying gradient method.^[2]^[5]

Inequality constraints can be incorporated, e.g., in order to handle obstacle problems.^[5]

Related Research Articles

In numerical analysis, Newton's method, also known as the Newton–Raphson method, named after Isaac Newton and Joseph Raphson, is a root-finding algorithm which produces successively better approximations to the roots of a real-valued function. The most basic version starts with a single-variable function $f$ defined for a real variable $x$ , the function's derivative $f '$ , and an initial guess $x 0$ for a root of $f$ . If the function satisfies sufficient assumptions and the initial guess is close, then

Gradient descent is a first-order iterative optimization algorithm for finding a local minimum of a differentiable function. The idea is to take repeated steps in the opposite direction of the gradient of the function at the current point, because this is the direction of steepest descent. Conversely, stepping in the direction of the gradient will lead to a local maximum of that function; the procedure is then known as gradient ascent.

In linear algebra, a QR decomposition, also known as a QR factorization or QU factorization is a decomposition of a matrix A into a product A = QR of an orthogonal matrix Q and an upper triangular matrix R. QR decomposition is often used to solve the linear least squares problem and is the basis for a particular eigenvalue algorithm, the QR algorithm.

The Gauss–Newton algorithm is used to solve non-linear least squares problems. It is a modification of Newton's method for finding a minimum of a function. Unlike Newton's method, the Gauss–Newton algorithm can only be used to minimize a sum of squared function values, but it has the advantage that second derivatives, which can be challenging to compute, are not required.

In mathematics, the conjugate gradient method is an algorithm for the numerical solution of particular systems of linear equations, namely those whose matrix is symmetric and positive-definite. The conjugate gradient method is often implemented as an iterative algorithm, applicable to sparse systems that are too large to be handled by a direct implementation or other direct methods such as the Cholesky decomposition. Large sparse systems often arise when numerically solving partial differential equations or optimization problems.

In optimization, the line search strategy is one of two basic iterative approaches to find a local minimum $of an objective function . The other approach is trust region.$

Interior-point methods are a certain class of algorithms that solve linear and nonlinear convex optimization problems.

In (unconstrained) minimization, a backtracking line search, a search scheme based on the Armijo–Goldstein condition, is a line search method to determine the amount to move along a given search direction. It involves starting with a relatively large estimate of the step size for movement along the search direction, and iteratively shrinking the step size until a decrease of the objective function is observed that adequately corresponds to the decrease that is expected, based on the local gradient of the objective function.

In the unconstrained minimization problem, the Wolfe conditions are a set of inequalities for performing inexact line search, especially in quasi-Newton methods, first published by Philip Wolfe in 1969.

In numerical optimization, the Broyden–Fletcher–Goldfarb–Shanno (BFGS) algorithm is an iterative method for solving unconstrained nonlinear optimization problems. Like the related Davidon–Fletcher–Powell method, BFGS determines the descent direction by preconditioning the gradient with curvature information. It does so by gradually improving an approximation to the Hessian matrix of the loss function, obtained only from gradient evaluations via a generalized secant method.

In mathematics, preconditioning is the application of a transformation, called the preconditioner, that conditions a given problem into a form that is more suitable for numerical solving methods. Preconditioning is typically related to reducing a condition number of the problem. The preconditioned problem is then usually solved by an iterative method.

Stochastic approximation methods are a family of iterative methods typically used for root-finding problems or for optimization problems. The recursive update rules of stochastic approximation methods can be used, among other things, for solving linear systems when the collected data is corrupted by noise, or for approximating extreme values of functions which cannot be computed directly, but only estimated via noisy observations.

In numerical optimization, the nonlinear conjugate gradient method generalizes the conjugate gradient method to nonlinear optimization. For a quadratic function

Subgradient methods are iterative methods for solving convex minimization problems. Originally developed by Naum Z. Shor and others in the 1960s and 1970s, subgradient methods are convergent when applied even to a non-differentiable objective function. When the objective function is differentiable, sub-gradient methods for unconstrained problems use the same search direction as the method of steepest descent.

In numerical linear algebra, the Alternating Direction Implicit (ADI) method is an iterative method used to solve Sylvester matrix equations. It is a popular method for solving the large matrix equations that arise in systems theory and control, and can be formulated to construct solutions in a memory-efficient, factored form. It is also used to numerically solve parabolic and elliptic partial differential equations, and is a classic method used for modeling heat conduction and solving the diffusion equation in two or more dimensions. It is an example of an operator splitting method.

In numerical analysis, the Schur complement method, named after Issai Schur, is the basic and the earliest version of non-overlapping domain decomposition method, also called iterative substructuring. A finite element problem is split into non-overlapping subdomains, and the unknowns in the interiors of the subdomains are eliminated. The remaining Schur complement system on the unknowns associated with subdomain interfaces is solved by the conjugate gradient method.

A differential equation is a mathematical equation for an unknown function of one or several variables that relates the values of the function itself and of its derivatives of various orders. A matrix differential equation contains more than one function stacked into vector form with a matrix relating the functions to their derivatives.

In statistics, the backfitting algorithm is a simple iterative procedure used to fit a generalized additive model. It was introduced in 1985 by Leo Breiman and Jerome Friedman along with generalized additive models. In most cases, the backfitting algorithm is equivalent to the Gauss–Seidel method algorithm for solving a certain linear system of equations.

Locally Optimal Block Preconditioned Conjugate Gradient (LOBPCG) is a matrix-free method for finding the largest eigenvalues and the corresponding eigenvectors of a symmetric positive definite generalized eigenvalue problem

The conjugate residual method is an iterative numeric method used for solving systems of linear equations. It's a Krylov subspace method very similar to the much more popular conjugate gradient method, with similar construction and convergence properties.

References

↑ Uzawa, H. (1958). "Iterative methods for concave programming". In Arrow, K. J.; Hurwicz, L.; Uzawa, H. (eds.). Studies in linear and nonlinear programming . Stanford University Press.
1 2 Elman, H. C.; Golub, G. H. (1994). "Inexact and preconditioned Uzawa algorithms for saddle point problems". SIAM J. Numer. Anal. 31 (6): 1645–1661. CiteSeerX 10.1.1.307.8178 . doi:10.1137/0731085.
↑ Bramble, J. H.; Pasciak, J. E.; Vassilev, A. T. (1997). "Analysis of the inexact Uzawa algorithm for saddle point problems". SIAM J. Numer. Anal. 34 (3): 1072–1982. CiteSeerX 10.1.1.52.9559 . doi:10.1137/S0036142994273343.
↑ Zulehner, W. (1998). "Analysis of iterative methods for saddle point problems. A unified approach". Math. Comp. 71 (238): 479–505. doi: 10.1090/S0025-5718-01-01324-2 .
1 2 Gräser, C.; Kornhuber, R. (2007). "On Preconditioned Uzawa-type Iterations for a Saddle Point Problem with Inequality Constraints". Domain Decomposition Methods in Science and Engineering XVI. Lec. Not. Comp. Sci. Eng. 55. pp. 91–102. CiteSeerX 10.1.1.72.9238 . doi:10.1007/978-3-540-34469-8_8. ISBN 978-3-540-34468-1.