16 #ifndef dealii_solver_richardson_h
17 #define dealii_solver_richardson_h
63 template <
typename VectorType = Vector<
double>>
111 template <
typename MatrixType,
typename PreconditionerType>
116 const PreconditionerType &preconditioner);
121 template <
typename MatrixType,
typename PreconditionerType>
126 const PreconditionerType &preconditioner);
143 const VectorType &
d)
const;
152 virtual typename VectorType::value_type
166 template <
typename VectorType>
169 const bool use_preconditioned_residual)
171 , use_preconditioned_residual(use_preconditioned_residual)
175 template <
typename VectorType>
178 const AdditionalData &data)
180 , additional_data(data)
185 template <
typename VectorType>
187 const AdditionalData &data)
189 , additional_data(data)
194 template <
typename VectorType>
195 template <
typename MatrixType,
typename PreconditionerType>
200 const PreconditionerType &preconditioner)
204 double last_criterion = std::numeric_limits<double>::lowest();
206 unsigned int iter = 0;
228 preconditioner.vmult(
d, r);
232 last_criterion = criterion(r,
d);
233 conv = this->iteration_status(iter, last_criterion, x);
237 x.add(additional_data.omega,
d);
238 print_vectors(iter, x, r,
d);
251 template <
typename VectorType>
252 template <
typename MatrixType,
typename PreconditionerType>
257 const PreconditionerType &preconditioner)
260 double last_criterion = std::numeric_limits<double>::lowest();
262 unsigned int iter = 0;
284 preconditioner.Tvmult(
d, r);
286 last_criterion = criterion(r,
d);
287 conv = this->iteration_status(iter, last_criterion, x);
291 x.add(additional_data.omega,
d);
292 print_vectors(iter, x, r,
d);
305 template <
typename VectorType>
310 const VectorType &)
const
315 template <
typename VectorType>
316 inline typename VectorType::value_type
318 const VectorType &
d)
const
320 if (!additional_data.use_preconditioned_residual)
327 template <
typename VectorType>
331 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)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
AdditionalData(const double omega=1, const bool use_preconditioned_residual=false)
bool use_preconditioned_residual