puzawa(5)
NAME
puzawa -- Uzawa algorithm.
SYNOPSIS
template <class Matrix, class Vector, class Preconditioner, class Real> int puzawa (const Matrix &A, Vector &x, const Vector &b, const Preconditioner &M, int &max_iter, Real &tol, const Real& rho, std::ostream *p_cerr=0);
EXAMPLE
- The simplest call to 'puzawa' has the folling form:
- size_t max_iter = 100;
double tol = 1e-7;
int status = puzawa(A, x, b, EYE, max_iter, tol, 1.0, &cerr);
DESCRIPTION
puzawa solves the linear system A*x=b using the Uzawa method. The Uzawa
method is a descent method in the direction opposite to the gradient,
with a constant step length 'rho'. The convergence is assured when the
step length 'rho' is small enough. If matrix A is symmetric positive
definite, please uses 'pcg' that computes automatically the optimal
descdnt step length at each iteration.
The return value indicates convergence within max_iter (input) iterations (0), or no convergence within max_iter iterations (1). Upon successful return, output arguments have the following values:
x approximate solution to Ax = b
- max_iter
- the number of iterations performed before the tolerance was reached
- tol the residual after the final iteration
IMPLEMENTATION
- template < class Matrix, class Vector, class Preconditioner, class Real, class Size> int puzawa(const Matrix &A, Vector &x, const Vector &Mb, const Preconditioner &M,
- Size &max_iter, Real &tol, const Real& rho,
std::ostream *p_cerr, std::string label) - {
- Vector b = M.solve(Mb);
Real norm2_b = dot(Mb,b);
Real norm2_r = norm2_b;
if (norm2_b == Real(0)) norm2_b = 1;
if (p_cerr) (*p_cerr) << "[" << label << "] #iteration residue" << std::endl;
for (Size n = 0; n <= max_iter; n++) {Vector Mr = A*x - Mb;
Vector r = M.solve(Mr);
norm2_r = dot(Mr, r);
if (p_cerr) (*p_cerr) << "[" << label << "] " << n << " " << sqrt(norm2_r/norm2_b) << std::endl;
if (norm2_r <= sqr(tol)*norm2_b) {tol = sqrt(norm2_r/norm2_b);
max_iter = n;
return 0;}
x -= rho*r; - }
tol = sqrt(norm2_r/norm2_b);
return 1; - }