From: hartmann Date: Mon, 3 May 1999 14:43:03 +0000 (+0000) Subject: Use of struct 'AdditionalData' in the constructor of each solver X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=108d0ac67baf209c3f5c5943ecdfdf5ac74b2aeb;p=dealii-svn.git Use of struct 'AdditionalData' in the constructor of each solver git-svn-id: https://svn.dealii.org/trunk@1245 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/lac/include/lac/solver.h b/deal.II/lac/include/lac/solver.h index a57b1567a7..7a3a6028ec 100644 --- a/deal.II/lac/include/lac/solver.h +++ b/deal.II/lac/include/lac/solver.h @@ -90,6 +90,29 @@ template class VectorMemory; * double b, const Vector& z); * }; * \end{verbatim} + * + * \subsection{AdditionalData} + * Several solvers need additional data, like the damping parameter #omega# of the + * #SolverRichardson# class or the maximum number of tmp vectors of the #SolverGMRES#. + * To have a standardized constructor for each solver class the #struct AdditionalData# + * has been introduced to each solver class. Some solver needs no additional data, like + * the SolverCG or SolverBicgstab. For these solvers the struct #AdditionalData# is empty + * and the calling + * of the constructor has not change as a default #AdditionalData# is set by default. + * + * Now the generating of a solver looks like + * \begin{verbatim} + * // GMRES with 50 tmp vectors + * SolverGMRES solver_gmres (solver_control, vector_memory, + * SolverGMRES::AdditionalData(50)); + * + * // Richardson with omega=0.8 + * SolverRichardson solver_richardson (solver_control, vector_memory, + * SolverGMRES::AdditionalData(0.8)); + * + * // CG with default AdditionalData + * SolverCG solver_cg (solver_control, vector_memory); + * \end{verbatim} */ template class Solver diff --git a/deal.II/lac/include/lac/solver_bicgstab.h b/deal.II/lac/include/lac/solver_bicgstab.h index 1b6b5e1da6..eb9df13282 100644 --- a/deal.II/lac/include/lac/solver_bicgstab.h +++ b/deal.II/lac/include/lac/solver_bicgstab.h @@ -11,16 +11,29 @@ /** * Bicgstab algorithm by van der Vorst. + * + * The use of the #AdditionalData# struct is described in the #solver# + * base class. */ template class SolverBicgstab // : public Solver { public: + /** + * Standardized data struct to + * pipe additional data to the + * solver. This solver does not + * need additional data yet. + */ + struct AdditionalData {}; + /** * Constructor. */ - SolverBicgstab(SolverControl &cn, VectorMemory &mem); + SolverBicgstab(SolverControl &cn, + VectorMemory &mem, + const AdditionalData &data=AdditionalData()); /** * Solve primal problem only. @@ -132,10 +145,9 @@ class SolverBicgstab // template SolverBicgstab::SolverBicgstab(SolverControl &cn, - VectorMemory &mem) - : - Solver(cn,mem) -{} + VectorMemory &mem, + const AdditionalData &) : + Solver(cn,mem) {} template diff --git a/deal.II/lac/include/lac/solver_cg.h b/deal.II/lac/include/lac/solver_cg.h index bac08f6a05..8b3c58a7b0 100644 --- a/deal.II/lac/include/lac/solver_cg.h +++ b/deal.II/lac/include/lac/solver_cg.h @@ -15,17 +15,28 @@ /** * Preconditioned cg method. * @author Original implementation by G. Kanschat, R. Becker and F.-T. Suttmeier, reworking and documentation by Wolfgang Bangerth + * + * The use of the #AdditionalData# struct is described in the #solver# + * base class. */ template class SolverCG : public Solver { public: + /** + * Standardized data struct to + * pipe additional data to the + * solver. This solver does not + * need additional data yet. + */ + struct AdditionalData {}; /** * Constructor. */ - SolverCG (SolverControl &cn, VectorMemory &mem) : - Solver(cn,mem) {}; + SolverCG (SolverControl &cn, + VectorMemory &mem, + const AdditionalData &data=AdditionalData()); /** * Solver method. @@ -70,6 +81,13 @@ class SolverCG : public Solver /*------------------------- Implementation ----------------------------*/ +template +SolverCG::SolverCG(SolverControl &cn, + VectorMemory &mem, + const AdditionalData &) : + Solver(cn,mem) {}; + + template double SolverCG::criterion() { diff --git a/deal.II/lac/include/lac/solver_gmres.h b/deal.II/lac/include/lac/solver_gmres.h index b57f16f851..869c45eb51 100644 --- a/deal.II/lac/include/lac/solver_gmres.h +++ b/deal.II/lac/include/lac/solver_gmres.h @@ -45,19 +45,39 @@ * computational speed against memory. One of the few other * possibilities is to use a good preconditioner. * + * The use of the #AdditionalData# struct is described in the #solver# + * base class. + * * @author Wolfgang Bangerth */ template class SolverGMRES : public Solver { public: + /** + * Standardized data struct to + * pipe additional data to the + * solver. + */ + struct AdditionalData + { + AdditionalData(const unsigned int max_n_tmp_vectors=30): + max_n_tmp_vectors(max_n_tmp_vectors) {}; + + /** + * Maximum number of + * tmp vectors. + */ + const unsigned int max_n_tmp_vectors; + }; + /** * Constructor. */ SolverGMRES (SolverControl &cn, VectorMemory &mem, - const unsigned int max_n_tmp_vectors); - + const AdditionalData &data=AdditionalData()); + /** * Solver method. */ @@ -74,7 +94,11 @@ class SolverGMRES : public Solver << "any results, and much more for reasonable ones."); protected: - const unsigned int max_n_tmp_vectors; + /** + * Includes the maximum number of + * tmp vectors. + */ + AdditionalData additional_data; /** * Implementation of the computation of @@ -102,11 +126,12 @@ class SolverGMRES : public Solver template SolverGMRES::SolverGMRES (SolverControl &cn, VectorMemory &mem, - const unsigned int max_n_tmp_vectors) : + const AdditionalData &data) : Solver (cn,mem), - max_n_tmp_vectors (max_n_tmp_vectors) + additional_data(data) { - Assert (max_n_tmp_vectors >= 10, ExcTooFewTmpVectors (max_n_tmp_vectors)); + Assert (additional_data.max_n_tmp_vectors >= 10, + ExcTooFewTmpVectors (additional_data.max_n_tmp_vectors)); }; @@ -162,8 +187,8 @@ SolverGMRES::solve (const Matrix& A, // if the size of the matrix is very small, // then only allocate a number of vectors // which is needed - const unsigned int n_tmp_vectors = (A.m()+3 > max_n_tmp_vectors ? - max_n_tmp_vectors : + const unsigned int n_tmp_vectors = (A.m()+3 > additional_data.max_n_tmp_vectors ? + additional_data.max_n_tmp_vectors : A.m()+3); // allocate an array of n_tmp_vectors diff --git a/deal.II/lac/include/lac/solver_richardson.h b/deal.II/lac/include/lac/solver_richardson.h index b4b2e8e41d..8912a1a8e0 100644 --- a/deal.II/lac/include/lac/solver_richardson.h +++ b/deal.II/lac/include/lac/solver_richardson.h @@ -16,18 +16,37 @@ * Implementation of the richardson iteration method. The stopping criterion * is the norm of the residual. * + * The use of the #AdditionalData# struct is described in the #solver# + * base class. + * * @author Ralf Hartmann */ template class SolverRichardson : public Solver { public: + /** + * Standardized data struct to + * pipe additional data to the + * solver. + */ + struct AdditionalData + { + AdditionalData(double omega=1): + omega(omega) {}; + + /** + * Damping parameter. + */ + double omega; + }; + /** * Constructor. */ - SolverRichardson (SolverControl &cn, VectorMemory &mem) : - Solver (cn,mem), omega(1.) - {}; + SolverRichardson (SolverControl &cn, + VectorMemory &mem, + const AdditionalData &data=AdditionalData()); /** * Solve $Ax=b$ for $x$. @@ -60,9 +79,9 @@ class SolverRichardson : public Solver Vector *Vr, *Vd; /** - * Damping-coefficient. + * Damping parameter. */ - double omega; ; + AdditionalData additional_data; /** * Within the iteration loop, the @@ -83,6 +102,14 @@ class SolverRichardson : public Solver /*----------------- Implementation of the Richardson Method ------------------*/ +template +SolverRichardson::SolverRichardson(SolverControl &cn, + VectorMemory &mem, + const AdditionalData &data): + Solver (cn,mem), + additional_data(data) {}; + + template template Solver::ReturnState @@ -110,7 +137,7 @@ SolverRichardson::solve (const Matrix &A, break; precondition(d,r); - x.add(omega,d); + x.add(additional_data.omega,d); } // Deallocate Memory @@ -133,11 +160,12 @@ SolverRichardson::criterion() return sqrt(res2); } + template inline void SolverRichardson::set_omega(double om) { - omega=om; + additional_data.omega=om; }