|
Reference documentation for deal.II version 9.2.0
|
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\)
\(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\)
\(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\)
\(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Go to the documentation of this file.
16 #ifndef dealii_solver_qmrs_h
17 #define dealii_solver_qmrs_h
94 template <
typename VectorType = Vector<
double>>
178 template <
typename MatrixType,
typename PreconditionerType>
180 solve(
const MatrixType &
A,
183 const PreconditionerType &preconditioner);
220 template <
typename MatrixType,
typename PreconditionerType>
225 const PreconditionerType &preconditioner,
244 template <
class VectorType>
247 const double last_residual)
249 , last_residual(last_residual)
254 template <
class VectorType>
257 const AdditionalData & data)
259 , additional_data(data)
263 template <
class VectorType>
265 const AdditionalData &data)
267 , additional_data(data)
271 template <
class VectorType>
279 template <
class VectorType>
280 template <
typename MatrixType,
typename PreconditionerType>
285 const PreconditionerType &preconditioner)
315 deallog <<
"Restart step " << step << std::endl;
316 state = iterate(
A, x,
b, preconditioner, *Vr, *Vu, *Vq, *Vt, *Vd);
327 template <
class VectorType>
328 template <
typename MatrixType,
typename PreconditionerType>
333 const PreconditionerType &preconditioner,
344 double tau, rho, theta = 0;
352 if (additional_data.left_preconditioning)
355 preconditioner.vmult(t, r);
362 preconditioner.vmult(q, t);
381 const double sigma = q * t;
384 if (additional_data.breakdown_testing ==
true &&
385 std::fabs(sigma) < additional_data.breakdown_threshold)
388 const double alpha = rho / sigma;
394 const double theta_old = theta;
397 if (additional_data.left_preconditioning)
400 preconditioner.vmult(t, r);
410 const double psi = 1. / (1. + theta);
414 d.sadd(psi * theta_old, psi * alpha, q);
417 print_vectors(step, x, r,
d);
425 if (res < additional_data.solver_tolerance)
431 state = this->iteration_status(step, res, x);
434 return IterationResult(state, res);
439 if (additional_data.breakdown_testing ==
true &&
440 std::fabs(sigma) < additional_data.breakdown_threshold)
443 const double rho_old = rho;
446 if (additional_data.left_preconditioning)
454 preconditioner.vmult(u, t);
459 const double beta = rho / rho_old;
SolverControl::State state
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
AdditionalData(const bool left_preconditioning=false, const double solver_tolerance=1.e-9, const bool breakdown_testing=true, const double breakdown_threshold=1.e-16)
Expression fabs(const Expression &x)
@ iterate
Continue iteration.
IterationResult(const SolverControl::State state, const double last_residual)
virtual void print_vectors(const unsigned int step, const VectorType &x, const VectorType &r, const VectorType &d) const
#define DEAL_II_NAMESPACE_OPEN
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
double breakdown_threshold
IterationResult iterate(const MatrixType &A, VectorType &x, const VectorType &b, const PreconditionerType &preconditioner, VectorType &r, VectorType &u, VectorType &q, VectorType &t, VectorType &d)
AdditionalData additional_data
inline ::VectorizedArray< Number, width > sqrt(const ::VectorizedArray< Number, width > &x)
SolverQMRS(SolverControl &cn, VectorMemory< VectorType > &mem, const AdditionalData &data=AdditionalData())
void solve(const MatrixType &A, VectorType &x, const VectorType &b, const PreconditionerType &preconditioner)
@ failure
Stop iteration, goal not reached.
@ success
Stop iteration, goal reached.
#define DEAL_II_NAMESPACE_CLOSE
#define AssertThrow(cond, exc)
bool left_preconditioning