16 #ifndef dealii__solver_richardson_h 17 #define dealii__solver_richardson_h 20 #include <deal.II/base/config.h> 21 #include <deal.II/base/logstream.h> 22 #include <deal.II/lac/solver.h> 23 #include <deal.II/lac/solver_control.h> 24 #include <deal.II/base/subscriptor.h> 25 #include <deal.II/base/signaling_nan.h> 27 DEAL_II_NAMESPACE_OPEN
61 template <
class VectorType = Vector<
double> >
111 template<
typename MatrixType,
typename PreconditionerType>
113 solve (
const MatrixType &A,
116 const PreconditionerType &precondition);
121 template<
typename MatrixType,
typename PreconditionerType>
123 Tsolve (
const MatrixType &A,
126 const PreconditionerType &precondition);
141 const VectorType &d)
const;
147 virtual typename VectorType::value_type
criterion();
172 typename VectorType::value_type
res;
180 template <
class VectorType>
187 use_preconditioned_residual(use_preconditioned_residual)
191 template <
class VectorType>
204 template <
class VectorType>
212 res(numbers::signaling_nan<double>())
217 template <
class VectorType>
222 template <
class VectorType>
223 template <
typename MatrixType,
typename PreconditionerType>
228 const PreconditionerType &precondition)
232 double last_criterion = -std::numeric_limits<double>::max();
234 unsigned int iter = 0;
244 deallog.
push(
"Richardson");
255 precondition.vmult(d,r);
293 template <
class VectorType>
294 template <
typename MatrixType,
typename PreconditionerType>
299 const PreconditionerType &precondition)
302 double last_criterion = -std::numeric_limits<double>::max();
304 unsigned int iter = 0;
314 deallog.
push(
"RichardsonT");
325 precondition.Tvmult(d,r);
358 template <
class VectorType>
363 const VectorType &)
const 368 template <
class VectorType>
369 inline typename VectorType::value_type
380 template <
class VectorType>
389 DEAL_II_NAMESPACE_CLOSE
bool use_preconditioned_residual
void set_omega(const double om=1.)
AdditionalData additional_data
#define AssertThrow(cond, exc)
SolverRichardson(SolverControl &cn, VectorMemory< VectorType > &mem, const AdditionalData &data=AdditionalData())
Stop iteration, goal reached.
virtual void print_vectors(const unsigned int step, const VectorType &x, const VectorType &r, const VectorType &d) const
VectorType::value_type res
VectorMemory< VectorType > & memory
boost::signals2::signal< SolverControl::State(const unsigned int iteration, const double check_value, const VectorType ¤t_iterate), StateCombiner > iteration_status
void solve(const MatrixType &A, VectorType &x, const VectorType &b, const PreconditionerType &precondition)
void push(const std::string &text)
void Tsolve(const MatrixType &A, VectorType &x, const VectorType &b, const PreconditionerType &precondition)
virtual ~SolverRichardson()
AdditionalData(const double omega=1, const bool use_preconditioned_residual=false)
virtual VectorType::value_type criterion()