From a56894899e0e89269758cfc246cf1a7a37e3e133 Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Wed, 29 Dec 2021 19:23:51 +0100 Subject: [PATCH] Refactor PreconditionRelaxation --- include/deal.II/lac/precondition.h | 906 ++++++++++++++++------------- 1 file changed, 502 insertions(+), 404 deletions(-) diff --git a/include/deal.II/lac/precondition.h b/include/deal.II/lac/precondition.h index 2334f5b667..4b0377e8a3 100644 --- a/include/deal.II/lac/precondition.h +++ b/include/deal.II/lac/precondition.h @@ -30,6 +30,7 @@ #include #include +#include #include #include @@ -398,7 +399,8 @@ private: * Jacobi, SOR and SSOR preconditioners are implemented. For preconditioning, * refer to derived classes. */ -template > +template , + typename PreconditionerType = IdentityMatrix> class PreconditionRelaxation : public Subscriptor { public: @@ -422,6 +424,11 @@ public: * Relaxation parameter. */ double relaxation; + + /* + * Preconditioner. + */ + std::shared_ptr preconditioner; }; /** @@ -453,6 +460,35 @@ public: size_type n() const; + /** + * Apply preconditioner. + */ + template + void + vmult(VectorType &, const VectorType &) const; + + /** + * Apply transpose preconditioner. Since this is a symmetric preconditioner, + * this function is the same as vmult(). + */ + template + void + Tvmult(VectorType &, const VectorType &) const; + + /** + * Perform one step of the preconditioned Richardson iteration + */ + template + void + step(VectorType &x, const VectorType &rhs) const; + + /** + * Perform one transposed step of the preconditioned Richardson iteration. + */ + template + void + Tstep(VectorType &x, const VectorType &rhs) const; + protected: /** * Pointer to the matrix object. @@ -463,10 +499,338 @@ protected: * Relaxation parameter. */ double relaxation; + + /* + * Preconditioner. + */ + std::shared_ptr preconditioner; }; +#ifndef DOXYGEN + +namespace internal +{ + namespace PreconditionRelaxation + { + template + class PreconditionJacobiImpl + { + public: + PreconditionJacobiImpl(const MatrixType &A, const double relaxation) + : A(&A) + , relaxation(relaxation) + {} + + template + void + vmult(VectorType &dst, const VectorType &src) const + { + this->A->precondition_Jacobi(dst, src, this->relaxation); + } + + template + void + Tvmult(VectorType &dst, const VectorType &src) const + { + // call vmult, since preconditioner is symmetrical + this->vmult(dst, src); + } + + template + void + step(VectorType &dst, const VectorType &src) const + { + this->A->Jacobi_step(dst, src, this->relaxation); + } + + template + void + Tstep(VectorType &dst, const VectorType &src) const + { + // call step, since preconditioner is symmetrical + this->step(dst, src); + } + + private: + const SmartPointer A; + const double relaxation; + }; + + template + class PreconditionSORImpl + { + public: + PreconditionSORImpl(const MatrixType &A, const double relaxation) + : A(&A) + , relaxation(relaxation) + {} + + template + void + vmult(VectorType &dst, const VectorType &src) const + { + this->A->precondition_SOR(dst, src, this->relaxation); + } + + template + void + Tvmult(VectorType &dst, const VectorType &src) const + { + this->A->precondition_TSOR(dst, src, this->relaxation); + } + + template + void + step(VectorType &dst, const VectorType &src) const + { + this->A->SOR_step(dst, src, this->relaxation); + } + + template + void + Tstep(VectorType &dst, const VectorType &src) const + { + this->A->TSOR_step(dst, src, this->relaxation); + } + + private: + const SmartPointer A; + const double relaxation; + }; + + template + class PreconditionSSORImpl + { + public: + using size_type = typename MatrixType::size_type; + + PreconditionSSORImpl(const MatrixType &A, const double relaxation) + : A(&A) + , relaxation(relaxation) + { + // in case we have a SparseMatrix class, we can extract information + // about the diagonal. + const SparseMatrix *mat = + dynamic_cast *>( + &*this->A); + + // calculate the positions first after the diagonal. + if (mat != nullptr) + { + const size_type n = this->A->n(); + pos_right_of_diagonal.resize(n, static_cast(-1)); + for (size_type row = 0; row < n; ++row) + { + // find the first element in this line which is on the right of + // the diagonal. we need to precondition with the elements on + // the left only. note: the first entry in each line denotes the + // diagonal element, which we need not check. + typename SparseMatrix< + typename MatrixType::value_type>::const_iterator it = + mat->begin(row) + 1; + for (; it < mat->end(row); ++it) + if (it->column() > row) + break; + pos_right_of_diagonal[row] = it - mat->begin(); + } + } + } + + template + void + vmult(VectorType &dst, const VectorType &src) const + { + this->A->precondition_SSOR(dst, + src, + this->relaxation, + pos_right_of_diagonal); + } + + template + void + Tvmult(VectorType &dst, const VectorType &src) const + { + this->A->precondition_SSOR(dst, + src, + this->relaxation, + pos_right_of_diagonal); + } + + template + void + step(VectorType &dst, const VectorType &src) const + { + this->A->SSOR_step(dst, src, this->relaxation); + } + + template + void + Tstep(VectorType &dst, const VectorType &src) const + { + // call step, since preconditioner is symmetrical + this->step(dst, src); + } + + private: + const SmartPointer A; + const double relaxation; + + /** + * An array that stores for each matrix row where the first position after + * the diagonal is located. + */ + std::vector pos_right_of_diagonal; + }; + + template + class PreconditionPSORImpl + { + public: + using size_type = typename MatrixType::size_type; + + PreconditionPSORImpl(const MatrixType & A, + const double relaxation, + const std::vector &permutation, + const std::vector &inverse_permutation) + : A(&A) + , relaxation(relaxation) + , permutation(permutation) + , inverse_permutation(inverse_permutation) + {} + + template + void + vmult(VectorType &dst, const VectorType &src) const + { + dst = src; + this->A->PSOR(dst, permutation, inverse_permutation, this->relaxation); + } + + template + void + Tvmult(VectorType &dst, const VectorType &src) const + { + dst = src; + this->A->TPSOR(dst, permutation, inverse_permutation, this->relaxation); + } + + private: + const SmartPointer A; + const double relaxation; + + const std::vector &permutation; + const std::vector &inverse_permutation; + }; + + template + struct has_step + { + private: + static bool + detect(...); + + template + static decltype( + std::declval().step(std::declval(), + std::declval())) + detect(const U &); + + public: + static const bool value = + !std::is_same()))>::value; + }; + + template + const bool has_step::value; + + template + struct has_Tstep + { + private: + static bool + detect(...); + + template + static decltype( + std::declval().Tstep(std::declval(), + std::declval())) + detect(const U &); + + public: + static const bool value = + !std::is_same()))>::value; + }; + + template + const bool has_Tstep::value; + + template < + typename PreconditionerType, + typename VectorType, + typename std::enable_if::value, + PreconditionerType>::type * = nullptr> + void + step(const PreconditionerType &preconditioner, + VectorType & dst, + const VectorType & src) + { + preconditioner.step(dst, src); + } + + template < + typename PreconditionerType, + typename VectorType, + typename std::enable_if::value, + PreconditionerType>::type * = nullptr> + void + step(const PreconditionerType &preconditioner, + VectorType & dst, + const VectorType & src) + { + Assert(false, ExcNotImplemented()); + (void)preconditioner; + (void)dst; + (void)src; + } + + template < + typename PreconditionerType, + typename VectorType, + typename std::enable_if::value, + PreconditionerType>::type * = nullptr> + void + Tstep(const PreconditionerType &preconditioner, + VectorType & dst, + const VectorType & src) + { + preconditioner.Tstep(dst, src); + } + + template < + typename PreconditionerType, + typename VectorType, + typename std::enable_if::value, + PreconditionerType>::type * = nullptr> + void + Tstep(const PreconditionerType &preconditioner, + VectorType & dst, + const VectorType & src) + { + Assert(false, ExcNotImplemented()); + (void)preconditioner; + (void)dst; + (void)src; + } + + } // namespace PreconditionRelaxation +} // namespace internal + +#endif + + + /** * Jacobi preconditioner using matrix built-in function. The * MatrixType class used is required to have a function @@ -494,48 +858,27 @@ protected: * @endcode */ template > -class PreconditionJacobi : public PreconditionRelaxation +class PreconditionJacobi + : public PreconditionRelaxation< + MatrixType, + internal::PreconditionRelaxation::PreconditionJacobiImpl> { -public: - /** - * Declare type for container size. - */ - using size_type = typename PreconditionRelaxation::size_type; + using PreconditionerType = + internal::PreconditionRelaxation::PreconditionJacobiImpl; + using BaseClass = PreconditionRelaxation; +public: /** * An alias to the base class AdditionalData. */ - using AdditionalData = - typename PreconditionRelaxation::AdditionalData; - - /** - * Apply preconditioner. - */ - template - void - vmult(VectorType &, const VectorType &) const; - - /** - * Apply transpose preconditioner. Since this is a symmetric preconditioner, - * this function is the same as vmult(). - */ - template - void - Tvmult(VectorType &, const VectorType &) const; + using AdditionalData = typename BaseClass::AdditionalData; /** - * Perform one step of the preconditioned Richardson iteration. + * @copydoc PreconditionRelaxation::initialize() */ - template void - step(VectorType &x, const VectorType &rhs) const; - - /** - * Perform one transposed step of the preconditioned Richardson iteration. - */ - template - void - Tstep(VectorType &x, const VectorType &rhs) const; + initialize(const MatrixType & A, + const AdditionalData ¶meters = AdditionalData()); }; @@ -585,47 +928,27 @@ public: * @endcode */ template > -class PreconditionSOR : public PreconditionRelaxation +class PreconditionSOR + : public PreconditionRelaxation< + MatrixType, + internal::PreconditionRelaxation::PreconditionSORImpl> { -public: - /** - * Declare type for container size. - */ - using size_type = typename PreconditionRelaxation::size_type; + using PreconditionerType = + internal::PreconditionRelaxation::PreconditionSORImpl; + using BaseClass = PreconditionRelaxation; +public: /** * An alias to the base class AdditionalData. */ - using AdditionalData = - typename PreconditionRelaxation::AdditionalData; + using AdditionalData = typename BaseClass::AdditionalData; /** - * Apply preconditioner. + * @copydoc PreconditionRelaxation::initialize() */ - template void - vmult(VectorType &, const VectorType &) const; - - /** - * Apply transpose preconditioner. - */ - template - void - Tvmult(VectorType &, const VectorType &) const; - - /** - * Perform one step of the preconditioned Richardson iteration. - */ - template - void - step(VectorType &x, const VectorType &rhs) const; - - /** - * Perform one transposed step of the preconditioned Richardson iteration. - */ - template - void - Tstep(VectorType &x, const VectorType &rhs) const; + initialize(const MatrixType & A, + const AdditionalData ¶meters = AdditionalData()); }; @@ -657,25 +980,20 @@ public: * @endcode */ template > -class PreconditionSSOR : public PreconditionRelaxation +class PreconditionSSOR + : public PreconditionRelaxation< + MatrixType, + internal::PreconditionRelaxation::PreconditionSSORImpl> { -public: - /** - * Declare type for container size. - */ - using size_type = typename PreconditionRelaxation::size_type; + using PreconditionerType = + internal::PreconditionRelaxation::PreconditionSSORImpl; + using BaseClass = PreconditionRelaxation; +public: /** * An alias to the base class AdditionalData. */ - using AdditionalData = - typename PreconditionRelaxation::AdditionalData; - - /** - * An alias to the base class. - */ - using BaseClass = PreconditionRelaxation; - + using AdditionalData = typename BaseClass::AdditionalData; /** * Initialize matrix and relaxation parameter. The matrix is just stored in @@ -683,46 +1001,8 @@ public: * zero and smaller than 2 for numerical reasons. It defaults to 1. */ void - initialize(const MatrixType & A, - const typename BaseClass::AdditionalData ¶meters = - typename BaseClass::AdditionalData()); - - /** - * Apply preconditioner. - */ - template - void - vmult(VectorType &, const VectorType &) const; - - /** - * Apply transpose preconditioner. Since this is a symmetric preconditioner, - * this function is the same as vmult(). - */ - template - void - Tvmult(VectorType &, const VectorType &) const; - - - /** - * Perform one step of the preconditioned Richardson iteration - */ - template - void - step(VectorType &x, const VectorType &rhs) const; - - /** - * Perform one transposed step of the preconditioned Richardson iteration. - */ - template - void - Tstep(VectorType &x, const VectorType &rhs) const; - -private: - /** - * An array that stores for each matrix row where the first position after - * the diagonal is located. - */ - std::vector pos_right_of_diagonal; + initialize(const MatrixType & A, + const AdditionalData ¶meters = AdditionalData()); }; @@ -756,13 +1036,20 @@ private: * @endcode */ template > -class PreconditionPSOR : public PreconditionRelaxation +class PreconditionPSOR + : public PreconditionRelaxation< + MatrixType, + internal::PreconditionRelaxation::PreconditionPSORImpl> { + using PreconditionerType = + internal::PreconditionRelaxation::PreconditionPSORImpl; + using BaseClass = PreconditionRelaxation; + public: /** * Declare type for container size. */ - using size_type = typename PreconditionRelaxation::size_type; + using size_type = typename BaseClass::size_type; /** * Parameters for PreconditionPSOR. @@ -780,12 +1067,10 @@ public: * The relaxation parameter should be larger than zero and smaller than 2 * for numerical reasons. It defaults to 1. */ - AdditionalData( - const std::vector &permutation, - const std::vector &inverse_permutation, - const typename PreconditionRelaxation::AdditionalData - ¶meters = - typename PreconditionRelaxation::AdditionalData()); + AdditionalData(const std::vector &permutation, + const std::vector &inverse_permutation, + const typename BaseClass::AdditionalData ¶meters = + typename BaseClass::AdditionalData()); /** * Storage for the permutation vector. @@ -798,7 +1083,7 @@ public: /** * Relaxation parameters */ - typename PreconditionRelaxation::AdditionalData parameters; + typename BaseClass::AdditionalData parameters; }; /** @@ -813,12 +1098,11 @@ public: * for numerical reasons. It defaults to 1. */ void - initialize(const MatrixType & A, - const std::vector &permutation, - const std::vector &inverse_permutation, - const typename PreconditionRelaxation::AdditionalData - ¶meters = - typename PreconditionRelaxation::AdditionalData()); + initialize(const MatrixType & A, + const std::vector & permutation, + const std::vector & inverse_permutation, + const typename BaseClass::AdditionalData ¶meters = + typename BaseClass::AdditionalData()); /** * Initialize matrix and relaxation parameter. The matrix is just stored in @@ -832,30 +1116,6 @@ public: */ void initialize(const MatrixType &A, const AdditionalData &additional_data); - - /** - * Apply preconditioner. - */ - template - void - vmult(VectorType &, const VectorType &) const; - - /** - * Apply transpose preconditioner. - */ - template - void - Tvmult(VectorType &, const VectorType &) const; - -private: - /** - * Storage for the permutation vector. - */ - const std::vector *permutation; - /** - * Storage for the inverse permutation vector. - */ - const std::vector *inverse_permutation; }; @@ -1486,271 +1746,138 @@ PreconditionRichardson::n() const //--------------------------------------------------------------------------- -template +template inline void -PreconditionRelaxation::initialize(const MatrixType & rA, - const AdditionalData ¶meters) +PreconditionRelaxation::initialize( + const MatrixType & rA, + const AdditionalData ¶meters) { - A = &rA; - relaxation = parameters.relaxation; + A = &rA; + relaxation = parameters.relaxation; + preconditioner = parameters.preconditioner; } -template +template inline void -PreconditionRelaxation::clear() +PreconditionRelaxation::clear() { A = nullptr; } -template -inline typename PreconditionRelaxation::size_type -PreconditionRelaxation::m() const +template +inline + typename PreconditionRelaxation::size_type + PreconditionRelaxation::m() const { Assert(A != nullptr, ExcNotInitialized()); return A->m(); } -template -inline typename PreconditionRelaxation::size_type -PreconditionRelaxation::n() const +template +inline + typename PreconditionRelaxation::size_type + PreconditionRelaxation::n() const { Assert(A != nullptr, ExcNotInitialized()); return A->n(); } -//--------------------------------------------------------------------------- - -template +template template inline void -PreconditionJacobi::vmult(VectorType & dst, - const VectorType &src) const +PreconditionRelaxation::vmult( + VectorType & dst, + const VectorType &src) const { - static_assert( - std::is_same::size_type, - typename VectorType::size_type>::value, - "PreconditionJacobi and VectorType must have the same size_type."); + preconditioner->vmult(dst, src); - Assert(this->A != nullptr, ExcNotInitialized()); - this->A->precondition_Jacobi(dst, src, this->relaxation); + if (this->relaxation != 1.0) + dst *= this->relaxation; } - - -template +template template inline void -PreconditionJacobi::Tvmult(VectorType & dst, - const VectorType &src) const +PreconditionRelaxation::Tvmult( + VectorType & dst, + const VectorType &src) const { - static_assert( - std::is_same::size_type, - typename VectorType::size_type>::value, - "PreconditionJacobi and VectorType must have the same size_type."); + preconditioner->Tvmult(dst, src); - Assert(this->A != nullptr, ExcNotInitialized()); - this->A->precondition_Jacobi(dst, src, this->relaxation); + if (this->relaxation != 1.0) + dst *= this->relaxation; } - - -template +template template inline void -PreconditionJacobi::step(VectorType & dst, - const VectorType &src) const +PreconditionRelaxation::step( + VectorType & dst, + const VectorType &src) const { - static_assert( - std::is_same::size_type, - typename VectorType::size_type>::value, - "PreconditionJacobi and VectorType must have the same size_type."); - - Assert(this->A != nullptr, ExcNotInitialized()); - this->A->Jacobi_step(dst, src, this->relaxation); + internal::PreconditionRelaxation::step(*preconditioner, dst, src); } - - -template +template template inline void -PreconditionJacobi::Tstep(VectorType & dst, - const VectorType &src) const +PreconditionRelaxation::Tstep( + VectorType & dst, + const VectorType &src) const { - static_assert( - std::is_same::size_type, - typename VectorType::size_type>::value, - "PreconditionJacobi and VectorType must have the same size_type."); - - step(dst, src); + internal::PreconditionRelaxation::Tstep(*preconditioner, dst, src); } - - //--------------------------------------------------------------------------- template -template -inline void -PreconditionSOR::vmult(VectorType &dst, const VectorType &src) const -{ - static_assert(std::is_same::size_type, - typename VectorType::size_type>::value, - "PreconditionSOR and VectorType must have the same size_type."); - - Assert(this->A != nullptr, ExcNotInitialized()); - this->A->precondition_SOR(dst, src, this->relaxation); -} - - - -template -template inline void -PreconditionSOR::Tvmult(VectorType & dst, - const VectorType &src) const +PreconditionJacobi::initialize(const MatrixType & A, + const AdditionalData ¶meters_in) { - static_assert(std::is_same::size_type, - typename VectorType::size_type>::value, - "PreconditionSOR and VectorType must have the same size_type."); - - Assert(this->A != nullptr, ExcNotInitialized()); - this->A->precondition_TSOR(dst, src, this->relaxation); -} - + Assert(parameters_in.preconditioner == nullptr, ExcInternalError()); + AdditionalData parameters; + parameters.relaxation = 1.0; + parameters.preconditioner = + std::make_shared(A, parameters_in.relaxation); -template -template -inline void -PreconditionSOR::step(VectorType &dst, const VectorType &src) const -{ - static_assert(std::is_same::size_type, - typename VectorType::size_type>::value, - "PreconditionSOR and VectorType must have the same size_type."); - - Assert(this->A != nullptr, ExcNotInitialized()); - this->A->SOR_step(dst, src, this->relaxation); + this->BaseClass::initialize(A, parameters); } - - -template -template -inline void -PreconditionSOR::Tstep(VectorType &dst, const VectorType &src) const -{ - static_assert(std::is_same::size_type, - typename VectorType::size_type>::value, - "PreconditionSOR and VectorType must have the same size_type."); - - Assert(this->A != nullptr, ExcNotInitialized()); - this->A->TSOR_step(dst, src, this->relaxation); -} - - - //--------------------------------------------------------------------------- template inline void -PreconditionSSOR::initialize( - const MatrixType & rA, - const typename BaseClass::AdditionalData ¶meters) +PreconditionSOR::initialize(const MatrixType & A, + const AdditionalData ¶meters_in) { - this->PreconditionRelaxation::initialize(rA, parameters); + Assert(parameters_in.preconditioner == nullptr, ExcInternalError()); - // in case we have a SparseMatrix class, we can extract information about - // the diagonal. - const SparseMatrix *mat = - dynamic_cast *>( - &*this->A); + AdditionalData parameters; + parameters.relaxation = 1.0; + parameters.preconditioner = + std::make_shared(A, parameters_in.relaxation); - // calculate the positions first after the diagonal. - if (mat != nullptr) - { - const size_type n = this->A->n(); - pos_right_of_diagonal.resize(n, static_cast(-1)); - for (size_type row = 0; row < n; ++row) - { - // find the first element in this line which is on the right of the - // diagonal. we need to precondition with the elements on the left - // only. note: the first entry in each line denotes the diagonal - // element, which we need not check. - typename SparseMatrix::const_iterator - it = mat->begin(row) + 1; - for (; it < mat->end(row); ++it) - if (it->column() > row) - break; - pos_right_of_diagonal[row] = it - mat->begin(); - } - } + this->BaseClass::initialize(A, parameters); } +//--------------------------------------------------------------------------- template -template -inline void -PreconditionSSOR::vmult(VectorType & dst, - const VectorType &src) const -{ - static_assert( - std::is_same::size_type, - typename VectorType::size_type>::value, - "PreconditionSSOR and VectorType must have the same size_type."); - - Assert(this->A != nullptr, ExcNotInitialized()); - this->A->precondition_SSOR(dst, src, this->relaxation, pos_right_of_diagonal); -} - - - -template -template inline void -PreconditionSSOR::Tvmult(VectorType & dst, - const VectorType &src) const +PreconditionSSOR::initialize(const MatrixType & A, + const AdditionalData ¶meters_in) { - static_assert( - std::is_same::size_type, - typename VectorType::size_type>::value, - "PreconditionSSOR and VectorType must have the same size_type."); + Assert(parameters_in.preconditioner == nullptr, ExcInternalError()); - Assert(this->A != nullptr, ExcNotInitialized()); - this->A->precondition_SSOR(dst, src, this->relaxation, pos_right_of_diagonal); -} + AdditionalData parameters; + parameters.relaxation = 1.0; + parameters.preconditioner = + std::make_shared(A, parameters_in.relaxation); - - -template -template -inline void -PreconditionSSOR::step(VectorType &dst, const VectorType &src) const -{ - static_assert( - std::is_same::size_type, - typename VectorType::size_type>::value, - "PreconditionSSOR and VectorType must have the same size_type."); - - Assert(this->A != nullptr, ExcNotInitialized()); - this->A->SSOR_step(dst, src, this->relaxation); -} - - - -template -template -inline void -PreconditionSSOR::Tstep(VectorType & dst, - const VectorType &src) const -{ - static_assert( - std::is_same::size_type, - typename VectorType::size_type>::value, - "PreconditionSSOR and VectorType must have the same size_type."); - - step(dst, src); + this->BaseClass::initialize(A, parameters); } @@ -1760,14 +1887,19 @@ PreconditionSSOR::Tstep(VectorType & dst, template inline void PreconditionPSOR::initialize( - const MatrixType & rA, - const std::vector & p, - const std::vector & ip, - const typename PreconditionRelaxation::AdditionalData ¶meters) + const MatrixType & A, + const std::vector & p, + const std::vector & ip, + const typename BaseClass::AdditionalData ¶meters_in) { - permutation = &p; - inverse_permutation = &ip; - PreconditionRelaxation::initialize(rA, parameters); + Assert(parameters_in.preconditioner == nullptr, ExcInternalError()); + + typename BaseClass::AdditionalData parameters; + parameters.relaxation = 1.0; + parameters.preconditioner = + std::make_shared(A, parameters_in.relaxation, p, ip); + + this->BaseClass::initialize(A, parameters); } @@ -1782,46 +1914,12 @@ PreconditionPSOR::initialize(const MatrixType & A, additional_data.parameters); } - -template -template -inline void -PreconditionPSOR::vmult(VectorType & dst, - const VectorType &src) const -{ - static_assert( - std::is_same::size_type, - typename VectorType::size_type>::value, - "PreconditionPSOR and VectorType must have the same size_type."); - - Assert(this->A != nullptr, ExcNotInitialized()); - dst = src; - this->A->PSOR(dst, *permutation, *inverse_permutation, this->relaxation); -} - - - -template -template -inline void -PreconditionPSOR::Tvmult(VectorType & dst, - const VectorType &src) const -{ - static_assert( - std::is_same::size_type, - typename VectorType::size_type>::value, - "PreconditionPSOR and VectorType must have the same size_type."); - - Assert(this->A != nullptr, ExcNotInitialized()); - dst = src; - this->A->TPSOR(dst, *permutation, *inverse_permutation, this->relaxation); -} - template PreconditionPSOR::AdditionalData::AdditionalData( const std::vector &permutation, const std::vector &inverse_permutation, - const typename PreconditionRelaxation::AdditionalData ¶meters) + const typename PreconditionPSOR::BaseClass::AdditionalData + ¶meters) : permutation(permutation) , inverse_permutation(inverse_permutation) , parameters(parameters) @@ -1852,9 +1950,9 @@ PreconditionUseMatrix::vmult( //--------------------------------------------------------------------------- -template -inline PreconditionRelaxation::AdditionalData::AdditionalData( - const double relaxation) +template +inline PreconditionRelaxation::AdditionalData:: + AdditionalData(const double relaxation) : relaxation(relaxation) {} -- 2.39.5