From: guido Date: Fri, 7 Jul 2000 22:29:03 +0000 (+0000) Subject: transposed solver X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=0e45a6150c3e5086139b108f7a654fd785d02bbf;p=dealii-svn.git transposed solver git-svn-id: https://svn.dealii.org/trunk@3144 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/lac/include/lac/solver_richardson.h b/deal.II/lac/include/lac/solver_richardson.h index 6aa55bc755..2c000f32a3 100644 --- a/deal.II/lac/include/lac/solver_richardson.h +++ b/deal.II/lac/include/lac/solver_richardson.h @@ -75,9 +75,18 @@ class SolverRichardson : private Solver */ template typename Solver::ReturnState solve (const MATRIX &A, - VECTOR &x, - const VECTOR &b, - const PRECONDITIONER& precondition); + VECTOR &x, + const VECTOR &b, + const PRECONDITIONER& precondition); + + /** + * Solve $A^Tx=b$ for $x$. + */ + template + typename Solver::ReturnState Tsolve (const MATRIX &A, + VECTOR &x, + const VECTOR &b, + const PRECONDITIONER& precondition); /** * Set the damping-coefficient. @@ -188,6 +197,53 @@ SolverRichardson::solve (const MATRIX &A, } +template +template +typename Solver::ReturnState +SolverRichardson::Tsolve (const MATRIX &A, + VECTOR &x, + const VECTOR &b, + const PRECONDITIONER& precondition) +{ + 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++) + { + // Do not use Tresidual, + // but do it in 2 steps + A.Tvmult(r,x); + r.sadd(-1.,1.,b); + res=sqrt(r*r); + + conv = control().check (iter, criterion()); + if (conv != SolverControl::iterate) + break; + + precondition.Tvmult(d,r); + x.add(additional_data.omega,d); + print_vectors(iter,x,r,d); + } + + // Deallocate Memory + memory.free(Vr); + memory.free(Vd); + + deallog.pop(); + // Output + if (conv == SolverControl::failure) + return exceeded; + else + return success; +} + + template void SolverRichardson::print_vectors(const unsigned int,