15#ifndef dealii_solver_richardson_h
16#define dealii_solver_richardson_h
62template <
typename VectorType = Vector<
double>>
110 template <
typename MatrixType,
typename PreconditionerType>
115 const PreconditionerType &preconditioner);
120 template <
typename MatrixType,
typename PreconditionerType>
125 const PreconditionerType &preconditioner);
142 const VectorType &d)
const;
151 virtual typename VectorType::value_type
152 criterion(
const VectorType &r,
const VectorType &d)
const;
165template <
typename VectorType>
168 const bool use_preconditioned_residual)
170 , use_preconditioned_residual(use_preconditioned_residual)
174template <
typename VectorType>
177 const AdditionalData &data)
179 , additional_data(data)
184template <
typename VectorType>
186 const AdditionalData &data)
188 , additional_data(data)
193template <
typename VectorType>
194template <
typename MatrixType,
typename PreconditionerType>
199 const PreconditionerType &preconditioner)
203 double last_criterion = std::numeric_limits<double>::lowest();
205 unsigned int iter = 0;
227 preconditioner.vmult(d, r);
231 last_criterion = criterion(r, d);
232 conv = this->iteration_status(iter, last_criterion, x);
236 x.add(additional_data.omega, d);
237 print_vectors(iter, x, r, d);
250template <
typename VectorType>
251template <
typename MatrixType,
typename PreconditionerType>
256 const PreconditionerType &preconditioner)
259 double last_criterion = std::numeric_limits<double>::lowest();
261 unsigned int iter = 0;
283 preconditioner.Tvmult(d, r);
285 last_criterion = criterion(r, d);
286 conv = this->iteration_status(iter, last_criterion, x);
290 x.add(additional_data.omega, d);
291 print_vectors(iter, x, r, d);
304template <
typename VectorType>
309 const VectorType &)
const
314template <
typename VectorType>
315inline typename VectorType::value_type
317 const VectorType &d)
const
319 if (!additional_data.use_preconditioned_residual)
326template <
typename VectorType>
330 additional_data.omega = om;
@ iterate
Continue iteration.
@ success
Stop iteration, goal reached.
virtual void print_vectors(const unsigned int step, const VectorType &x, const VectorType &r, const VectorType &d) const
void solve(const MatrixType &A, VectorType &x, const VectorType &b, const PreconditionerType &preconditioner)
void set_omega(const double om=1.)
SolverRichardson(SolverControl &cn, VectorMemory< VectorType > &mem, const AdditionalData &data=AdditionalData())
void Tsolve(const MatrixType &A, VectorType &x, const VectorType &b, const PreconditionerType &preconditioner)
virtual ~SolverRichardson() override=default
virtual VectorType::value_type criterion(const VectorType &r, const VectorType &d) const
AdditionalData additional_data
SolverRichardson(SolverControl &cn, const AdditionalData &data=AdditionalData())
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define AssertThrow(cond, exc)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
AdditionalData(const double omega=1, const bool use_preconditioned_residual=false)
bool use_preconditioned_residual