]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Use of struct 'AdditionalData' in the constructor of each solver
authorhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 3 May 1999 14:43:03 +0000 (14:43 +0000)
committerhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 3 May 1999 14:43:03 +0000 (14:43 +0000)
git-svn-id: https://svn.dealii.org/trunk@1245 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/solver.h
deal.II/lac/include/lac/solver_bicgstab.h
deal.II/lac/include/lac/solver_cg.h
deal.II/lac/include/lac/solver_gmres.h
deal.II/lac/include/lac/solver_richardson.h

index a57b1567a738625ac23f384edae6134ba3554093..7a3a6028ec4bcfec627adbf679f874329bcc71ad 100644 (file)
@@ -90,6 +90,29 @@ template <class Vector> 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 Matrix, class Vector>
 class Solver
index 1b6b5e1da6043f2f559e55a6ce1c8f26a3696744..eb9df13282da8ae581218cb4cf51c54b45bfa3bf 100644 (file)
 
 /**
  * Bicgstab algorithm by van der Vorst.
+ *
+ * The use of the #AdditionalData# struct is described in the #solver#
+ * base class.
  */
 template<class Matrix, class Vector>
 class SolverBicgstab //<Matrix,Vector>
   : public Solver<Matrix,Vector>
 {
   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<Vector> &mem);
+    SolverBicgstab(SolverControl &cn,
+                  VectorMemory<Vector> &mem,
+                  const AdditionalData &data=AdditionalData());
 
                                     /**
                                      * Solve primal problem only.
@@ -132,10 +145,9 @@ class SolverBicgstab //<Matrix,Vector>
 
 template<class Matrix, class Vector>
 SolverBicgstab<Matrix, Vector>::SolverBicgstab(SolverControl &cn,
-                                              VectorMemory<Vector> &mem)
-               :
-               Solver<Matrix,Vector>(cn,mem)
-{}
+                                              VectorMemory<Vector> &mem,
+                                              const AdditionalData &) :
+               Solver<Matrix,Vector>(cn,mem)  {}
 
 
 template<class Matrix, class Vector>
index bac08f6a059b9d16e23f3cc1dd1c30dc211e317d..8b3c58a7b0f49bbb36371002956e5b5620c233f3 100644 (file)
 /**
  * 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 Matrix, class Vector>
 class SolverCG : public Solver<Matrix,Vector>
 {
   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<Vector> &mem) :
-                   Solver<Matrix,Vector>(cn,mem) {};
+    SolverCG (SolverControl &cn,
+             VectorMemory<Vector> &mem,
+             const AdditionalData &data=AdditionalData());
 
                                     /**
                                      * Solver method.
@@ -70,6 +81,13 @@ class SolverCG : public Solver<Matrix,Vector>
 /*------------------------- Implementation ----------------------------*/
  
 
+template<class Matrix, class Vector>
+SolverCG<Matrix,Vector>::SolverCG(SolverControl &cn,
+                                 VectorMemory<Vector> &mem,
+                                 const AdditionalData &) :
+               Solver<Matrix,Vector>(cn,mem) {};
+
+
 template<class Matrix, class Vector>
 double SolverCG<Matrix,Vector>::criterion()
 {
index b57f16f85196e625e4ef1c598d4c218f8c1e6e84..869c45eb51b8a68ce4b88e261df481a7fc8a4d6e 100644 (file)
  * 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 Matrix, class Vector>
 class SolverGMRES : public Solver<Matrix, Vector>
 {
   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<Vector> &mem,
-                const unsigned int    max_n_tmp_vectors);
-    
+                const AdditionalData &data=AdditionalData());
+
                                     /**
                                      * Solver method.
                                      */
@@ -74,7 +94,11 @@ class SolverGMRES : public Solver<Matrix, Vector>
                    << "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<Matrix, Vector>
 template <class Matrix, class Vector>
 SolverGMRES<Matrix,Vector>::SolverGMRES (SolverControl        &cn,
                                         VectorMemory<Vector> &mem,
-                                        const unsigned int    max_n_tmp_vectors) :
+                                        const AdditionalData &data) :
                Solver<Matrix,Vector> (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<Matrix,Vector>::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
index b4b2e8e41db9a2c3f5ed13992c7e34629f0e01e3..8912a1a8e04a51df3615ba5d95eee345343e09d1 100644 (file)
  * 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 Matrix, class Vector>
 class SolverRichardson : public Solver<Matrix, Vector>
 {
   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<Vector> &mem) :
-                   Solver<Matrix,Vector> (cn,mem), omega(1.)
-      {};
+    SolverRichardson (SolverControl &cn,
+                     VectorMemory<Vector> &mem,
+                     const AdditionalData &data=AdditionalData());
 
                                     /**
                                      * Solve $Ax=b$ for $x$.
@@ -60,9 +79,9 @@ class SolverRichardson : public Solver<Matrix, Vector>
     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<Matrix, Vector>
 /*----------------- Implementation of the Richardson Method ------------------*/
 
 
+template<class Matrix, class Vector>
+SolverRichardson<Matrix,Vector>::SolverRichardson(SolverControl &cn,
+                                                 VectorMemory<Vector> &mem,
+                                                 const AdditionalData &data):
+               Solver<Matrix,Vector> (cn,mem),
+               additional_data(data)  {};
+
+
 template<class Matrix,class Vector>
 template<class Preconditioner>
 Solver<Matrix,Vector>::ReturnState 
@@ -110,7 +137,7 @@ SolverRichardson<Matrix,Vector>::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<Matrix,Vector>::criterion()
   return sqrt(res2);
 }
 
+
 template<class Matrix,class Vector>
 inline void
 SolverRichardson<Matrix,Vector>::set_omega(double om)
 {
-  omega=om;
+  additional_data.omega=om;
 }
 
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.