From: David Wells Date: Thu, 12 Nov 2015 20:15:05 +0000 (-0500) Subject: Prefer PreconditionerType to PRECONDITIONER. X-Git-Tag: v8.4.0-rc2~230^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=ef605548bf9168c80949a471604b246e8d9b8888;p=dealii.git Prefer PreconditionerType to PRECONDITIONER. --- diff --git a/include/deal.II/lac/iterative_inverse.h b/include/deal.II/lac/iterative_inverse.h index 53cc554919..8ca447fb97 100644 --- a/include/deal.II/lac/iterative_inverse.h +++ b/include/deal.II/lac/iterative_inverse.h @@ -83,8 +83,8 @@ public: * Initialization function. Provide a matrix and preconditioner for the * solve in vmult(). */ - template - void initialize (const MatrixType &, const PRECONDITION &); + template + void initialize (const MatrixType &, const PreconditionerType &); /** * Delete the pointers to matrix and preconditioner. @@ -124,10 +124,10 @@ private: template -template +template inline void -IterativeInverse::initialize(const MatrixType &m, const PRECONDITION &p) +IterativeInverse::initialize(const MatrixType &m, const PreconditionerType &p) { // dummy variable VectorType *v = 0; diff --git a/include/deal.II/lac/matrix_lib.h b/include/deal.II/lac/matrix_lib.h index adcc672693..3cec3a1edf 100644 --- a/include/deal.II/lac/matrix_lib.h +++ b/include/deal.II/lac/matrix_lib.h @@ -449,9 +449,9 @@ public: * Initialization function. Provide a solver object, a matrix, and another * preconditioner for this. */ - template + template void initialize (const MatrixType &, - const PRECONDITION &); + const PreconditionerType &); /** * Access to the SolverControl object used by the solver. @@ -735,16 +735,17 @@ MeanValueFilter::Tvmult_add(VectorType &, const VectorType &) const //-----------------------------------------------------------------------// template -template +template inline void -InverseMatrixRichardson::initialize (const MatrixType &m, const PRECONDITION &p) +InverseMatrixRichardson::initialize (const MatrixType &m, + const PreconditionerType &p) { if (matrix != 0) delete matrix; matrix = new PointerMatrix(&m); if (precondition != 0) delete precondition; - precondition = new PointerMatrix(&p); + precondition = new PointerMatrix(&p); } diff --git a/include/deal.II/lac/solver_bicgstab.h b/include/deal.II/lac/solver_bicgstab.h index 0ffc4ab648..7b097ebb25 100644 --- a/include/deal.II/lac/solver_bicgstab.h +++ b/include/deal.II/lac/solver_bicgstab.h @@ -128,12 +128,12 @@ public: /** * Solve primal problem only. */ - template + template void - solve (const MatrixType &A, - VectorType &x, - const VectorType &b, - const PRECONDITIONER &precondition); + solve (const MatrixType &A, + VectorType &x, + const VectorType &b, + const PreconditionerType &precondition); protected: /** @@ -253,10 +253,10 @@ private: * The iteration loop itself. The function returns a structure indicating * what happened in this function. */ - template + template IterationResult - iterate(const MatrixType &A, - const PRECONDITIONER &precondition); + iterate(const MatrixType &A, + const PreconditionerType &precondition); }; /*@}*/ @@ -347,10 +347,10 @@ SolverBicgstab::print_vectors(const unsigned int, template -template +template typename SolverBicgstab::IterationResult -SolverBicgstab::iterate(const MatrixType &A, - const PRECONDITIONER &precondition) +SolverBicgstab::iterate(const MatrixType &A, + const PreconditionerType &precondition) { //TODO:[GK] Implement "use the length of the computed orthogonal residual" in the BiCGStab method. SolverControl::State state = SolverControl::iterate; @@ -434,12 +434,12 @@ SolverBicgstab::iterate(const MatrixType &A, template -template +template void -SolverBicgstab::solve(const MatrixType &A, - VectorType &x, - const VectorType &b, - const PRECONDITIONER &precondition) +SolverBicgstab::solve(const MatrixType &A, + VectorType &x, + const VectorType &b, + const PreconditionerType &precondition) { deallog.push("Bicgstab"); Vr = this->memory.alloc(); diff --git a/include/deal.II/lac/solver_cg.h b/include/deal.II/lac/solver_cg.h index 61d1734521..d65703b13f 100644 --- a/include/deal.II/lac/solver_cg.h +++ b/include/deal.II/lac/solver_cg.h @@ -175,12 +175,12 @@ public: /** * Solve the linear system $Ax=b$ for x. */ - template + template void - solve (const MatrixType &A, - VectorType &x, - const VectorType &b, - const PRECONDITIONER &precondition); + solve (const MatrixType &A, + VectorType &x, + const VectorType &b, + const PreconditionerType &precondition); /** * Connect a slot to retrieve the CG coefficients. The slot will be called @@ -453,12 +453,12 @@ SolverCG::compute_eigs_and_cond template -template +template void -SolverCG::solve (const MatrixType &A, - VectorType &x, - const VectorType &b, - const PRECONDITIONER &precondition) +SolverCG::solve (const MatrixType &A, + VectorType &x, + const VectorType &b, + const PreconditionerType &precondition) { SolverControl::State conv=SolverControl::iterate; @@ -521,7 +521,7 @@ SolverCG::solve (const MatrixType &A, return; } - if (types_are_equal::value == false) + if (types_are_equal::value == false) { precondition.vmult(h,g); @@ -553,7 +553,7 @@ SolverCG::solve (const MatrixType &A, if (conv != SolverControl::iterate) break; - if (types_are_equal::value + if (types_are_equal::value == false) { precondition.vmult(h,g); diff --git a/include/deal.II/lac/solver_gmres.h b/include/deal.II/lac/solver_gmres.h index 5a29a11b3c..0e1effa64a 100644 --- a/include/deal.II/lac/solver_gmres.h +++ b/include/deal.II/lac/solver_gmres.h @@ -259,12 +259,12 @@ public: /** * Solve the linear system $Ax=b$ for x. */ - template + template void - solve (const MatrixType &A, - VectorType &x, - const VectorType &b, - const PRECONDITIONER &precondition); + solve (const MatrixType &A, + VectorType &x, + const VectorType &b, + const PreconditionerType &precondition); /** * Connect a slot to retrieve the estimated condition number. @@ -452,12 +452,12 @@ public: /** * Solve the linear system $Ax=b$ for x. */ - template + template void - solve (const MatrixType &A, - VectorType &x, - const VectorType &b, - const PRECONDITIONER &precondition); + solve (const MatrixType &A, + VectorType &x, + const VectorType &b, + const PreconditionerType &precondition); private: @@ -746,12 +746,12 @@ SolverGMRES::compute_eigs_and_cond template -template +template void -SolverGMRES::solve (const MatrixType &A, - VectorType &x, - const VectorType &b, - const PRECONDITIONER &precondition) +SolverGMRES::solve (const MatrixType &A, + VectorType &x, + const VectorType &b, + const PreconditionerType &precondition) { // this code was written a very long time ago by people not associated with // deal.II. we don't make any guarantees to its optimality or that it even @@ -1114,12 +1114,12 @@ SolverFGMRES::SolverFGMRES (SolverControl &cn, template -template +template void -SolverFGMRES::solve (const MatrixType &A, - VectorType &x, - const VectorType &b, - const PRECONDITIONER &precondition) +SolverFGMRES::solve (const MatrixType &A, + VectorType &x, + const VectorType &b, + const PreconditionerType &precondition) { deallog.push("FGMRES"); diff --git a/include/deal.II/lac/solver_minres.h b/include/deal.II/lac/solver_minres.h index 0b37270d3d..23d825839b 100644 --- a/include/deal.II/lac/solver_minres.h +++ b/include/deal.II/lac/solver_minres.h @@ -98,12 +98,12 @@ public: /** * Solve the linear system $Ax=b$ for x. */ - template + template void - solve (const MatrixType &A, - VectorType &x, - const VectorType &b, - const PRECONDITIONER &precondition); + solve (const MatrixType &A, + VectorType &x, + const VectorType &b, + const PreconditionerType &precondition); /** * @addtogroup Exceptions @@ -196,12 +196,12 @@ SolverMinRes::print_vectors(const unsigned int, template -template +template void -SolverMinRes::solve (const MatrixType &A, - VectorType &x, - const VectorType &b, - const PRECONDITIONER &precondition) +SolverMinRes::solve (const MatrixType &A, + VectorType &x, + const VectorType &b, + const PreconditionerType &precondition) { SolverControl::State conv=SolverControl::iterate; diff --git a/include/deal.II/lac/solver_qmrs.h b/include/deal.II/lac/solver_qmrs.h index dc128bb2ef..412d2372bb 100644 --- a/include/deal.II/lac/solver_qmrs.h +++ b/include/deal.II/lac/solver_qmrs.h @@ -125,12 +125,12 @@ public: /** * Solve the linear system $Ax=b$ for x. */ - template + template void - solve (const MatrixType &A, - VectorType &x, - const VectorType &b, - const PRECONDITIONER &precondition); + solve (const MatrixType &A, + VectorType &x, + const VectorType &b, + const PreconditionerType &precondition); /** * Interface for derived class. This function gets the current iteration @@ -198,10 +198,10 @@ private: * The iteration loop itself. The function returns a structure indicating * what happened in this function. */ - template + template IterationResult - iterate (const MatrixType &A, - const PRECONDITIONER &precondition); + iterate (const MatrixType &A, + const PreconditionerType &precondition); /** * Number of the current iteration (accumulated over restarts) @@ -264,12 +264,12 @@ SolverQMRS::print_vectors(const unsigned int, template -template +template void -SolverQMRS::solve (const MatrixType &A, - VectorType &x, - const VectorType &b, - const PRECONDITIONER &precondition) +SolverQMRS::solve (const MatrixType &A, + VectorType &x, + const VectorType &b, + const PreconditionerType &precondition) { deallog.push("QMRS"); @@ -322,10 +322,10 @@ SolverQMRS::solve (const MatrixType &A, template -template +template typename SolverQMRS::IterationResult -SolverQMRS::iterate(const MatrixType &A, - const PRECONDITIONER &precondition) +SolverQMRS::iterate(const MatrixType &A, + const PreconditionerType &precondition) { /* Remark: the matrix A in the article is the preconditioned matrix. * Therefore, we have to precondition x before we compute the first residual. diff --git a/include/deal.II/lac/solver_richardson.h b/include/deal.II/lac/solver_richardson.h index 34d6bb9638..20f07b2bd2 100644 --- a/include/deal.II/lac/solver_richardson.h +++ b/include/deal.II/lac/solver_richardson.h @@ -107,22 +107,22 @@ public: /** * Solve the linear system $Ax=b$ for x. */ - template + template void - solve (const MatrixType &A, - VectorType &x, - const VectorType &b, - const PRECONDITIONER &precondition); + solve (const MatrixType &A, + VectorType &x, + const VectorType &b, + const PreconditionerType &precondition); /** * Solve $A^Tx=b$ for $x$. */ - template + template void - Tsolve (const MatrixType &A, - VectorType &x, - const VectorType &b, - const PRECONDITIONER &precondition); + Tsolve (const MatrixType &A, + VectorType &x, + const VectorType &b, + const PreconditionerType &precondition); /** * Set the damping-coefficient. Default is 1., i.e. no damping. @@ -214,12 +214,12 @@ SolverRichardson::~SolverRichardson() template -template +template void -SolverRichardson::solve (const MatrixType &A, - VectorType &x, - const VectorType &b, - const PRECONDITIONER &precondition) +SolverRichardson::solve (const MatrixType &A, + VectorType &x, + const VectorType &b, + const PreconditionerType &precondition) { SolverControl::State conv=SolverControl::iterate; @@ -285,12 +285,12 @@ SolverRichardson::solve (const MatrixType &A, template -template +template void -SolverRichardson::Tsolve (const MatrixType &A, - VectorType &x, - const VectorType &b, - const PRECONDITIONER &precondition) +SolverRichardson::Tsolve (const MatrixType &A, + VectorType &x, + const VectorType &b, + const PreconditionerType &precondition) { SolverControl::State conv=SolverControl::iterate; double last_criterion = -std::numeric_limits::max(); diff --git a/include/deal.II/multigrid/mg_coarse.h b/include/deal.II/multigrid/mg_coarse.h index 628414b8cc..fa79b357f5 100644 --- a/include/deal.II/multigrid/mg_coarse.h +++ b/include/deal.II/multigrid/mg_coarse.h @@ -50,10 +50,10 @@ public: * Constructor. Store solver, matrix and preconditioning method for later * use. */ - template + template MGCoarseGridLACIteration (SolverType &, const MatrixType &, - const PRECOND &); + const PreconditionerType &); /** * Destructor freeing the pointers. @@ -63,10 +63,10 @@ public: /** * Initialize new data. */ - template + template void initialize (SolverType &, const MatrixType &, - const PRECOND &); + const PreconditionerType &); /** * Clear all pointers. @@ -198,16 +198,16 @@ MGCoarseGridLACIteration template -template +template MGCoarseGridLACIteration -::MGCoarseGridLACIteration (SolverType &s, - const MatrixType &m, - const PRECOND &p) +::MGCoarseGridLACIteration (SolverType &s, + const MatrixType &m, + const PreconditionerType &p) : solver(&s, typeid(*this).name()) { matrix = new PointerMatrix(&m); - precondition = new PointerMatrix(&p); + precondition = new PointerMatrix(&p); } @@ -220,12 +220,12 @@ MGCoarseGridLACIteration template -template +template void MGCoarseGridLACIteration -::initialize (SolverType &s, - const MatrixType &m, - const PRECOND &p) +::initialize (SolverType &s, + const MatrixType &m, + const PreconditionerType &p) { solver = &s; if (matrix) @@ -233,7 +233,7 @@ MGCoarseGridLACIteration matrix = new PointerMatrix(&m); if (precondition) delete precondition; - precondition = new PointerMatrix(&p); + precondition = new PointerMatrix(&p); } diff --git a/include/deal.II/multigrid/mg_smoother.h b/include/deal.II/multigrid/mg_smoother.h index ebc9a66ac6..3dded6693c 100644 --- a/include/deal.II/multigrid/mg_smoother.h +++ b/include/deal.II/multigrid/mg_smoother.h @@ -399,7 +399,7 @@ private: * * @author Guido Kanschat, 2009 */ -template +template class MGSmootherPrecondition : public MGSmoother { public: @@ -416,19 +416,19 @@ public: * matrices and initializes the smoothing operator with the same smoother * for each level. * - * @p additional_data is an object of type @p PRECONDITIONER::AdditionalData + * @p additional_data is an object of type @p PreconditionerType::AdditionalData * and is handed to the initialization function of the relaxation method. */ template void initialize (const MGLevelObject &matrices, - const typename PRECONDITIONER::AdditionalData &additional_data = typename PRECONDITIONER::AdditionalData()); + const typename PreconditionerType::AdditionalData &additional_data = typename PreconditionerType::AdditionalData()); /** * Initialize for matrices. This function stores pointers to the level * matrices and initializes the smoothing operator with the according * smoother for each level. * - * @p additional_data is an object of type @p PRECONDITIONER::AdditionalData + * @p additional_data is an object of type @p PreconditionerType::AdditionalData * and is handed to the initialization function of the relaxation method. */ template @@ -441,7 +441,7 @@ public: * This function stores pointers to the level matrices and initializes the * smoothing operator with the same smoother for each level. * - * @p additional_data is an object of type @p PRECONDITIONER::AdditionalData + * @p additional_data is an object of type @p PreconditionerType::AdditionalData * and is handed to the initialization function of the relaxation method. */ template @@ -456,7 +456,7 @@ public: * This function stores pointers to the level matrices and initializes the * smoothing operator with the according smoother for each level. * - * @p additional_data is an object of type @p PRECONDITIONER::AdditionalData + * @p additional_data is an object of type @p PreconditionerType::AdditionalData * and is handed to the initialization function of the relaxation method. */ template @@ -480,7 +480,7 @@ public: /** * Object containing relaxation methods. */ - MGLevelObject smoothers; + MGLevelObject smoothers; /** * Memory used by this object. @@ -840,9 +840,9 @@ memory_consumption () const //----------------------------------------------------------------------// -template +template inline -MGSmootherPrecondition::MGSmootherPrecondition +MGSmootherPrecondition::MGSmootherPrecondition (const unsigned int steps, const bool variable, const bool symmetric, @@ -853,9 +853,9 @@ MGSmootherPrecondition::MGSmootherPrecon -template +template inline void -MGSmootherPrecondition::clear () +MGSmootherPrecondition::clear () { smoothers.clear(); @@ -867,12 +867,12 @@ MGSmootherPrecondition::clear () -template +template template inline void -MGSmootherPrecondition::initialize +MGSmootherPrecondition::initialize (const MGLevelObject &m, - const typename PRECONDITIONER::AdditionalData &data) + const typename PreconditionerType::AdditionalData &data) { const unsigned int min = m.min_level(); const unsigned int max = m.max_level(); @@ -889,10 +889,10 @@ MGSmootherPrecondition::initialize -template +template template inline void -MGSmootherPrecondition::initialize +MGSmootherPrecondition::initialize (const MGLevelObject &m, const MGLevelObject &data) { @@ -916,10 +916,10 @@ MGSmootherPrecondition::initialize -template +template template inline void -MGSmootherPrecondition::initialize +MGSmootherPrecondition::initialize (const MGLevelObject &m, const DATA &data, const unsigned int row, @@ -940,10 +940,10 @@ MGSmootherPrecondition::initialize -template +template template inline void -MGSmootherPrecondition::initialize +MGSmootherPrecondition::initialize (const MGLevelObject &m, const MGLevelObject &data, const unsigned int row, @@ -969,9 +969,9 @@ MGSmootherPrecondition::initialize -template +template inline void -MGSmootherPrecondition::smooth +MGSmootherPrecondition::smooth (const unsigned int level, VectorType &u, const VectorType &rhs) const @@ -1040,10 +1040,10 @@ MGSmootherPrecondition::smooth -template +template inline std::size_t -MGSmootherPrecondition:: +MGSmootherPrecondition:: memory_consumption () const { return sizeof(*this)