From 4982716dce615204474cec9e12b1ac245198f013 Mon Sep 17 00:00:00 2001 From: hartmann Date: Mon, 19 Apr 1999 11:57:45 +0000 Subject: [PATCH] Richardson's iteration method as preconditioned solver. git-svn-id: https://svn.dealii.org/trunk@1177 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/lac/include/lac/solver_richardson.h | 141 ++++++++++++++++++++ 1 file changed, 141 insertions(+) create mode 100644 deal.II/lac/include/lac/solver_richardson.h diff --git a/deal.II/lac/include/lac/solver_richardson.h b/deal.II/lac/include/lac/solver_richardson.h new file mode 100644 index 0000000000..877f3f2d35 --- /dev/null +++ b/deal.II/lac/include/lac/solver_richardson.h @@ -0,0 +1,141 @@ +/*---------------------------- solver_richardson.h ---------------------------*/ +/* $Id$ */ +/* Ralf Hartmann, University of Heidelberg */ +#ifndef __solver_richardson_H +#define __solver_richardson_H +/*---------------------------- solver_richardson.h ---------------------------*/ + + + +#include +#include + + + +/** + * Implementation of the richardson iteration method. The stopping criterion + * is the norm of the residual. + * + * @author Ralf Hartmann + */ +template +class SolverRichardson : public Solver { + public: + /** + * Constructor. + */ + SolverRichardson (SolverControl &cn, VectorMemory &mem) : + Solver (cn,mem), omega(1.) + {}; + + /** + * Solve $Ax=b$ for $x$. + */ + virtual ReturnState solve (const Matrix &A, + Vector &x, + const Vector &b); + + /** + * Set the damping-coefficient. + * Default is 1., i.e. no damping. + */ + void set_omega(double om=1.); + + protected: + /** + * Implementation of the computation of + * the norm of the residual. + */ + virtual double criterion(); + + /** + * Temporary vectors, allocated through + * the #VectorMemory# object at the start + * of the actual solution process and + * deallocated at the end. + */ + Vector *Vr, *Vd; + + /** + * Damping-coefficient. + */ + double omega; ; + + /** + * Within the iteration loop, the + * square of the residual vector is + * stored in this variable. The + * function #criterion# uses this + * variable to compute the convergence + * value, which in this class is the + * norm of the residual vector and thus + * the square root of the #res2# value. + */ + double res2; +}; + + + + +/*----------------- Implementation of the Richardson Method ------------------*/ + + +template +Solver::ReturnState +SolverRichardson::solve (const Matrix &A, + Vector &x, + const Vector &b) { + SolverControl::State conv=SolverControl::iterate; + + // Memory allocation + Vr = memory.alloc(); Vector& r = *Vr; r.reinit(x); + Vd = memory.alloc(); Vector& d = *Vd; d.reinit(x); + + deallog.push("Richardson"); + + // Main loop + for(int iter=0; conv==SolverControl::iterate; iter++) + { + A.residual(r,x,b); + + res2 = r*r; + conv = control().check (iter, criterion()); + if (conv != SolverControl::iterate) + break; + + A.precondition(d,r); + x.add(omega,d); + } + + // Deallocate Memory + memory.free(Vr); + memory.free(Vd); + + deallog.pop(); + // Output + if (conv == SolverControl::failure) + return exceeded; + else + return success; +} + + +template +inline double +SolverRichardson::criterion() +{ + return sqrt(res2); +} + +template +inline void +SolverRichardson::set_omega(double om) +{ + omega=om; +} + + +/*------------------ solver_richardson.h ----------------------*/ +/* end of #ifndef __solver_richardson_H */ +#endif +/*------------------ solver_richardson.h ----------------------*/ -- 2.39.5