]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
MGSmoother base class introduced as template
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 14 Aug 2002 15:23:23 +0000 (15:23 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 14 Aug 2002 15:23:23 +0000 (15:23 +0000)
git-svn-id: https://svn.dealii.org/trunk@6326 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/multigrid/mg_smoother.h
deal.II/deal.II/include/multigrid/mg_smoother.templates.h
deal.II/deal.II/include/multigrid/multigrid.h
deal.II/deal.II/source/multigrid/mg_smoother.cc

index b4808ebf21fade0854c01d12d41d2d85e1dde542..b01defeb2a418bd7fb27abed28182339d2fd3967 100644 (file)
 
 template <int dim> class MGDoFHandler;
 
+
+/**
+ * Base class for multigrid smoothers. Does nothing but defining the
+ * interface used by multigrid methods.
+ *
+ * @author Guido KAnschat, 2002
+ */
+template <class VECTOR>
+class MGSmoother : public Subscriptor
+{
+  public:
+                                  /**
+                                   * Virtual destructor.
+                                   */
+  virtual ~MGSmoother();
+
+                                  /**
+                                   * Smoothing function. This is the
+                                   * function used in multigrid
+                                   * methods.
+                                   */
+  virtual void smooth (const unsigned int level,
+                      VECTOR&            u,
+                      const VECTOR&      rhs) const = 0;  
+};
+
+
 /**
  * Smoother doing nothing. This class is not useful for many applications other
  * than for testing some multigrid procedures. Also some applications might
  * get convergence without smoothing and then this class brings you the
  * cheapest possible multigrid.
  *
- * @author Guido Kanschat, 1999
+ * @author Guido Kanschat, 1999, 2002
  */
-class MGSmootherIdentity : public Subscriptor
+template <class VECTOR>
+class MGSmootherIdentity : public MGSmoother<VECTOR>
 {
   public:
                                     /**
@@ -42,10 +70,9 @@ class MGSmootherIdentity : public Subscriptor
                                      * operator equals the null
                                      * operator.
                                      */
-  template <class VECTOR>
-  void smooth (const unsigned int level,
-              VECTOR&            u,
-              const VECTOR&      rhs) const;
+  virtual void smooth (const unsigned int level,
+                      VECTOR&            u,
+                      const VECTOR&      rhs) const;
 };
 
 
@@ -56,9 +83,10 @@ class MGSmootherIdentity : public Subscriptor
  * triangulation) to zero, by building a list of these degrees of
  * freedom's indices at construction time.
  *
- * @author Wolfgang Bangerth, Guido Kanschat, 1999
+ * @author Wolfgang Bangerth, Guido Kanschat, 1999, 2002
  */
-class MGSmootherContinuous : public Subscriptor
+template <class VECTOR>
+class MGSmootherContinuous : public MGSmoother<VECTOR>
 {
   private:
                                     /**
@@ -122,7 +150,6 @@ class MGSmootherContinuous : public Subscriptor
                                      * has no interior boundaries, this
                                      * function does nothing in this case.
                                      */
-    template <class VECTOR>
     void set_zero_interior_boundary (const unsigned int level,
                                     VECTOR&            u) const;
     
@@ -173,7 +200,7 @@ class MGSmootherContinuous : public Subscriptor
  * @author Wolfgang Bangerth, Guido Kanschat, 1999
  */
 template<class MATRIX, class VECTOR>
-class MGSmootherRelaxation : public MGSmootherContinuous
+class MGSmootherRelaxation : public MGSmootherContinuous<VECTOR>
 {
   public:
                                     /**
@@ -214,9 +241,9 @@ class MGSmootherRelaxation : public MGSmootherContinuous
                                      * work and find a way to keep
                                      * the boundary values.
                                      */
-    void smooth (const unsigned int level,
-                VECTOR&            u,
-                const VECTOR&      rhs) const;
+    virtual void smooth (const unsigned int level,
+                        VECTOR&            u,
+                        const VECTOR&      rhs) const;
   private:
                                     /**
                                      * Pointer to the matrices.
@@ -236,12 +263,12 @@ class MGSmootherRelaxation : public MGSmootherContinuous
                                     /**
                                      * Auxiliary vector.
                                      */
-    mutable Vector<double> h1;
+    mutable VECTOR h1;
     
                                     /**
                                      * Auxiliary vector.
                                      */
-    mutable Vector<double> h2;
+    mutable VECTOR h2;
 };
 
 
@@ -249,21 +276,25 @@ class MGSmootherRelaxation : public MGSmootherContinuous
 
 template <class VECTOR>
 inline void
-MGSmootherIdentity::smooth (const unsigned int, VECTOR&, const VECTOR&) const
+MGSmootherIdentity<VECTOR>::smooth (
+  const unsigned int, VECTOR&,
+  const VECTOR&) const
 {};
 
 /*----------------------------------------------------------------------*/
 
+template <class VECTOR>
 inline
 void
-MGSmootherContinuous::set_steps(unsigned int i)
+MGSmootherContinuous<VECTOR>::set_steps(unsigned int i)
 {
   steps = i;
 }
 
+template <class VECTOR>
 inline
 unsigned int
-MGSmootherContinuous::get_steps() const
+MGSmootherContinuous<VECTOR>::get_steps() const
 {
   return steps;
 }
index c481c04359b0cf81897c7a4f6e4323f8748ba4f2..3700a525fbf2d10cc60481b8bd8802689e3f0992 100644 (file)
@@ -24,7 +24,7 @@ MGSmootherRelaxation<MATRIX, VECTOR>
                        const unsigned int           steps,
                        const double                 omega)
                :
-               MGSmootherContinuous(mg_dof, steps),
+               MGSmootherContinuous<VECTOR>(mg_dof, steps),
                matrix(&matrix),
                relaxation(relaxation),
                omega(omega)
index 809103a0e013505fd92cbc9d44dfb6212a5726d1..eb890b259a048a4a0e7775361cf5d886e5d12b45 100644 (file)
@@ -87,6 +87,9 @@ class Multigrid : public Subscriptor
                                      */
   template <int dim>
     Multigrid(const MGDoFHandler<dim>& mg_dof_handler,
+             const MGTransfer<VECTOR>& transfer,
+             const MGSmoother<VECTOR>& pre_smooth,
+             const MGSmoother<VECTOR>& post_smooth,
              const unsigned int            minlevel = 0,
              const unsigned int            maxlevel = 1000000);
     
@@ -111,11 +114,8 @@ class Multigrid : public Subscriptor
                                      * function is done in
                                      * @p{level_mgstep}.
                                      */
-    template <class MATRIX, class SMOOTHER, class COARSE>
+    template <class MATRIX, class COARSE>
     void vcycle(const MGLevelObject<MATRIX>& matrix,
-               const MGTransfer<VECTOR>&    transfer,
-               const SMOOTHER&              pre_smooth,
-               const SMOOTHER&              post_smooth,
                const COARSE&                cgs);
 
                                     /**
@@ -143,12 +143,9 @@ class Multigrid : public Subscriptor
                                      * is used to solve the matrix of
                                      * this level.
                                      */
-    template <class MATRIX, class SMOOTHER, class COARSE>
+    template <class MATRIX, class COARSE>
     void level_mgstep (const unsigned int           level,
                       const MGLevelObject<MATRIX>& matrix,
-                      const MGTransfer<VECTOR>& transfer,
-                      const SMOOTHER&              pre_smooth,
-                      const SMOOTHER&              post_smooth,
                       const COARSE&                cgs);
                                     /**
                                      * Level for coarse grid solution.
@@ -181,6 +178,20 @@ class Multigrid : public Subscriptor
 
 
   private:
+                                    /**
+                                     * Object for grid tranfer.
+                                     */
+    SmartPointer<const MGTransfer<VECTOR> > transfer;
+                                    /**
+                                     * The pre-smoothing object.
+                                     */
+    SmartPointer<const MGSmoother<VECTOR> > pre_smooth;
+
+                                    /**
+                                     * The post-smoothing object.
+                                     */
+    SmartPointer<const MGSmoother<VECTOR> > post_smooth;
+    
     
 
                                     /**
@@ -189,7 +200,7 @@ class Multigrid : public Subscriptor
     DeclException2(ExcSwitchedLevels, int, int,
                   << "minlevel and maxlevel switched, should be: "
                   << arg1 << "<=" << arg2);
-    template<int dim, class MATRIX, class VECTOR2, class SMOOTHER, class TRANSFER, class COARSE> friend class PreconditionMG;
+    template<int dim, class MATRIX, class VECTOR2, class TRANSFER, class COARSE> friend class PreconditionMG;
 };
 
 
@@ -205,7 +216,7 @@ class Multigrid : public Subscriptor
  *
  * @author Guido Kanschat, 1999, 2000, 2001, 2002
  */
-template<int dim, class MATRIX, class VECTOR, class TRANSFER, class SMOOTHER, class COARSE>
+template<int dim, class MATRIX, class VECTOR, class TRANSFER, class COARSE>
 class PreconditionMG : public Subscriptor
 {
   public:
@@ -219,8 +230,6 @@ class PreconditionMG : public Subscriptor
                   Multigrid<VECTOR>&           mg,
                   const MGLevelObject<MATRIX>& matrix,
                   const TRANSFER& transfer,
-                  const SMOOTHER&              pre,
-                  const SMOOTHER&              post,
                   const COARSE&                coarse);
     
                                     /**
@@ -287,16 +296,6 @@ class PreconditionMG : public Subscriptor
                                    */
     SmartPointer<const TRANSFER> transfer;
   
-                                    /**
-                                     * The pre-smoothing object.
-                                     */
-    SmartPointer<const SMOOTHER> pre;
-
-                                    /**
-                                     * The post-smoothing object.
-                                     */
-    SmartPointer<const SMOOTHER> post;
-    
                                     /**
                                      * The coarse grid solver.
                                      */
@@ -312,27 +311,30 @@ class PreconditionMG : public Subscriptor
 template <class VECTOR>
 template <int dim>
 Multigrid<VECTOR>::Multigrid (const MGDoFHandler<dim>& mg_dof_handler,
+                             const MGTransfer<VECTOR>& transfer,
+                             const MGSmoother<VECTOR>& pre_smooth,
+                             const MGSmoother<VECTOR>& post_smooth,
                              const unsigned int min_level,
                              const unsigned int max_level)
                :
                minlevel(min_level),
-                                    maxlevel(std::min(mg_dof_handler.get_tria().n_levels()-1,
-                                                      max_level)),
-                                    defect(minlevel,maxlevel),
-                                    solution(minlevel,maxlevel),
-                                    t(minlevel,maxlevel)
+  maxlevel(std::min(mg_dof_handler.get_tria().n_levels()-1,
+                   max_level)),
+  defect(minlevel,maxlevel),
+  solution(minlevel,maxlevel),
+  t(minlevel,maxlevel),
+  transfer(&transfer),
+  pre_smooth(&pre_smooth),
+  post_smooth(&post_smooth)
 {
 };
 
 
 template <class VECTOR>
-template <class MATRIX, class SMOOTHER, class COARSE>
+template <class MATRIX, class COARSE>
 void
 Multigrid<VECTOR>::level_mgstep(const unsigned int        level,
                                const MGLevelObject<MATRIX>& matrix,
-                               const MGTransfer<VECTOR>& transfer,
-                               const SMOOTHER     &pre_smooth,
-                               const SMOOTHER     &post_smooth,
                                const COARSE &coarse_grid_solver)
 {
 #ifdef MG_DEBUG
@@ -354,7 +356,7 @@ Multigrid<VECTOR>::level_mgstep(const unsigned int        level,
     }
 
                           // smoothing of the residual by modifying s
-  pre_smooth.smooth(level, solution[level], defect[level]);
+  pre_smooth->smooth(level, solution[level], defect[level]);
 
 #ifdef MG_DEBUG
   sprintf(name, "MG%d-pre",level);
@@ -373,7 +375,7 @@ Multigrid<VECTOR>::level_mgstep(const unsigned int        level,
   for (unsigned int l = level;l>minlevel;--l)
     {
       t[l-1] = 0.;
-      transfer.restrict_and_add (l, t[l-1], t[l]);
+      transfer->restrict_and_add (l, t[l-1], t[l]);
       defect[l-1] += t[l-1];
     }
 
@@ -381,7 +383,7 @@ Multigrid<VECTOR>::level_mgstep(const unsigned int        level,
 //  edge_vmult(level, defect[level-1], defect[level]);
   
                                   // do recursion
-  level_mgstep(level-1, matrix, transfer, pre_smooth, post_smooth, coarse_grid_solver);
+  level_mgstep(level-1, matrix, coarse_grid_solver);
 
                                   // reset size of the auxiliary
                                   // vector, since it has been
@@ -391,7 +393,7 @@ Multigrid<VECTOR>::level_mgstep(const unsigned int        level,
 
                                   // do coarse grid correction
 
-  transfer.prolongate(level, t[level], solution[level-1]);
+  transfer->prolongate(level, t[level], solution[level-1]);
 
 #ifdef MG_DEBUG
   sprintf(name, "MG%d-cgc",level);
@@ -401,7 +403,7 @@ Multigrid<VECTOR>::level_mgstep(const unsigned int        level,
   solution[level] += t[level];
   
                                   // post-smoothing
-  post_smooth.smooth(level, solution[level], defect[level]);
+  post_smooth->smooth(level, solution[level], defect[level]);
 
 #ifdef MG_DEBUG
   sprintf(name, "MG%d-post",level);
@@ -413,12 +415,9 @@ Multigrid<VECTOR>::level_mgstep(const unsigned int        level,
 
 
 template <class VECTOR>
-template <class MATRIX, class SMOOTHER, class COARSE>
+template <class MATRIX, class COARSE>
 void
 Multigrid<VECTOR>::vcycle(const MGLevelObject<MATRIX>& matrix,
-                         const MGTransfer<VECTOR>& transfer,
-                         const SMOOTHER     &pre_smooth,
-                         const SMOOTHER     &post_smooth,
                          const COARSE &coarse_grid_solver)
 {
                                   // The defect vector has been
@@ -433,7 +432,7 @@ Multigrid<VECTOR>::vcycle(const MGLevelObject<MATRIX>& matrix,
       t[level].reinit(defect[level]);
     }
 
-  level_mgstep (maxlevel, matrix, transfer, pre_smooth, post_smooth, coarse_grid_solver);
+  level_mgstep (maxlevel, matrix, coarse_grid_solver);
 //  abort ();
 }
 
@@ -441,29 +440,25 @@ Multigrid<VECTOR>::vcycle(const MGLevelObject<MATRIX>& matrix,
 /* ------------------------------------------------------------------- */
 
 
-template<int dim, class MATRIX, class VECTOR, class TRANSFER, class SMOOTHER, class COARSE>
-PreconditionMG<dim, MATRIX, VECTOR, TRANSFER, SMOOTHER, COARSE>
+template<int dim, class MATRIX, class VECTOR, class TRANSFER, class COARSE>
+PreconditionMG<dim, MATRIX, VECTOR, TRANSFER, COARSE>
 ::PreconditionMG(const MGDoFHandler<dim>&      mg_dof_handler,
                 Multigrid<VECTOR>& mg,
                 const MGLevelObject<MATRIX>& matrix,
                 const TRANSFER&              transfer,
-                const SMOOTHER&              pre,
-                const SMOOTHER&              post,
                 const COARSE&                coarse)
                :
   mg_dof_handler(&mg_dof_handler),
   multigrid(&mg),
   level_matrices(&matrix),
   transfer(&transfer),
-  pre(&pre),
-  post(&post),
   coarse(&coarse)
 {}
 
 
-template<int dim, class MATRIX, class VECTOR, class TRANSFER, class SMOOTHER, class COARSE>
+template<int dim, class MATRIX, class VECTOR, class TRANSFER, class COARSE>
 void
-PreconditionMG<dim, MATRIX, VECTOR, TRANSFER, SMOOTHER, COARSE>::vmult (
+PreconditionMG<dim, MATRIX, VECTOR, TRANSFER, COARSE>::vmult (
   VECTOR& dst,
   const VECTOR& src) const
 {
@@ -471,9 +466,7 @@ PreconditionMG<dim, MATRIX, VECTOR, TRANSFER, SMOOTHER, COARSE>::vmult (
                       multigrid->defect,
                       src);
   
-  multigrid->vcycle(*level_matrices,
-                   *transfer,
-                   *pre, *post, *coarse);
+  multigrid->vcycle(*level_matrices, *coarse);
 
   transfer->copy_from_mg(*mg_dof_handler,
                         dst,
@@ -481,9 +474,9 @@ PreconditionMG<dim, MATRIX, VECTOR, TRANSFER, SMOOTHER, COARSE>::vmult (
 }
 
 
-template<int dim, class MATRIX, class VECTOR, class TRANSFER, class SMOOTHER, class COARSE>
+template<int dim, class MATRIX, class VECTOR, class TRANSFER, class COARSE>
 void
-PreconditionMG<dim, MATRIX, VECTOR, TRANSFER, SMOOTHER, COARSE>::vmult_add (
+PreconditionMG<dim, MATRIX, VECTOR, TRANSFER, COARSE>::vmult_add (
   VECTOR& dst,
   const VECTOR& src) const
 {
@@ -491,9 +484,7 @@ PreconditionMG<dim, MATRIX, VECTOR, TRANSFER, SMOOTHER, COARSE>::vmult_add (
                       multigrid->defect,
                       src);
   
-  multigrid->vcycle(*level_matrices,
-                   *transfer,
-                   *pre, *post, *coarse);
+  multigrid->vcycle(*level_matrices, *coarse);
 
   transfer->copy_from_mg_add(*mg_dof_handler,
                             dst,
@@ -501,9 +492,9 @@ PreconditionMG<dim, MATRIX, VECTOR, TRANSFER, SMOOTHER, COARSE>::vmult_add (
 }
 
 
-template<int dim, class MATRIX, class VECTOR, class TRANSFER, class SMOOTHER, class COARSE>
+template<int dim, class MATRIX, class VECTOR, class TRANSFER, class COARSE>
 void
-PreconditionMG<dim, MATRIX, VECTOR, TRANSFER, SMOOTHER, COARSE>::Tvmult (
+PreconditionMG<dim, MATRIX, VECTOR, TRANSFER, COARSE>::Tvmult (
   VECTOR&,
   const VECTOR&) const
 {
@@ -511,9 +502,9 @@ PreconditionMG<dim, MATRIX, VECTOR, TRANSFER, SMOOTHER, COARSE>::Tvmult (
 }
 
 
-template<int dim, class MATRIX, class VECTOR, class TRANSFER, class SMOOTHER, class COARSE>
+template<int dim, class MATRIX, class VECTOR, class TRANSFER, class COARSE>
 void
-PreconditionMG<dim, MATRIX, VECTOR, TRANSFER, SMOOTHER, COARSE>::Tvmult_add (
+PreconditionMG<dim, MATRIX, VECTOR, TRANSFER, COARSE>::Tvmult_add (
   VECTOR&,
   const VECTOR&) const
 {
index 827b0b16729d5133504e26e9422c5e455cd3ece0..18e0cef626226fa1398f47c8f46efa7474986099 100644 (file)
@@ -18,7 +18,6 @@
 #include <multigrid/mg_dof_accessor.h>
 #include <grid/tria_iterator.h>
 #include <fe/fe.h>
-//#include <multigrid/multigrid.h>
 #include <multigrid/mg_smoother.h>
 #include <multigrid/mg_smoother.templates.h>
 #include <lac/vector.h>
 #include <algorithm>
 
 
+template <class VECTOR>
+MGSmoother<VECTOR>::~MGSmoother()
+{}
+
+
 //////////////////////////////////////////////////////////////////////
 
 
 
 #if deal_II_dimension == 1
 
-MGSmootherContinuous::MGSmootherContinuous (const MGDoFHandler<1> &/*mg_dof*/,
-                                           const unsigned int steps)
+template <class VECTOR>
+MGSmootherContinuous<VECTOR>::MGSmootherContinuous (
+  const MGDoFHandler<1> &/*mg_dof*/,
+  const unsigned int steps)
                :
                steps(steps)
 {
@@ -48,9 +54,11 @@ MGSmootherContinuous::MGSmootherContinuous (const MGDoFHandler<1> &/*mg_dof*/,
 
 #if deal_II_dimension > 1
 
+template <class VECTOR>
 template <int dim>
-MGSmootherContinuous::MGSmootherContinuous (const MGDoFHandler<dim> &mg_dof,
-                                           const unsigned int steps)
+MGSmootherContinuous<VECTOR>::MGSmootherContinuous (
+  const MGDoFHandler<dim> &mg_dof,
+  const unsigned int steps)
                :
                steps(steps)
 {
@@ -115,7 +123,7 @@ MGSmootherContinuous::MGSmootherContinuous (const MGDoFHandler<dim> &mg_dof,
 
 template <class VECTOR>
 void
-MGSmootherContinuous::set_zero_interior_boundary (const unsigned int level,
+MGSmootherContinuous<VECTOR>::set_zero_interior_boundary (const unsigned int level,
                                                  VECTOR&            u) const
 {
   if (level==0)
@@ -130,28 +138,45 @@ MGSmootherContinuous::set_zero_interior_boundary (const unsigned int level,
 
 
 // explicit instantiations
+
+template MGSmoother<Vector<float> >::~MGSmoother();
+template MGSmoother<Vector<double> >::~MGSmoother();
+template MGSmoother<BlockVector<float> >::~MGSmoother();
+template MGSmoother<BlockVector<double> >::~MGSmoother();
+
 // don't do the following instantiation in 1d, since there is a specialized
 // function there
 #if deal_II_dimension > 1
-template MGSmootherContinuous::MGSmootherContinuous (const MGDoFHandler<deal_II_dimension>&,
-                                                    const unsigned int);
+template MGSmootherContinuous<Vector<float> >::MGSmootherContinuous (
+  const MGDoFHandler<deal_II_dimension>&,
+  const unsigned int);
+template MGSmootherContinuous<Vector<double> >::MGSmootherContinuous (
+  const MGDoFHandler<deal_II_dimension>&,
+  const unsigned int);
+template MGSmootherContinuous<BlockVector<float> >::MGSmootherContinuous (
+  const MGDoFHandler<deal_II_dimension>&,
+  const unsigned int);
+template MGSmootherContinuous<BlockVector<double> >::MGSmootherContinuous (
+  const MGDoFHandler<deal_II_dimension>&,
+  const unsigned int);
 #endif
 
 template
-void MGSmootherContinuous::set_zero_interior_boundary<Vector<double> > (const unsigned int,
-                                                        Vector<double>&) const;
-
+void MGSmootherContinuous<Vector<double> >::set_zero_interior_boundary (
+  const unsigned int,
+  Vector<double>&) const;
 template
-void MGSmootherContinuous::set_zero_interior_boundary<Vector<float> > (const unsigned int,
-                                                        Vector<float>&) const;
-
+void MGSmootherContinuous<Vector<float> >::set_zero_interior_boundary (
+  const unsigned int,
+  Vector<float>&) const;
 template
-void MGSmootherContinuous::set_zero_interior_boundary<BlockVector<double> > (const unsigned int,
-                                                        BlockVector<double>&) const;
-
+void MGSmootherContinuous<BlockVector<double> >::set_zero_interior_boundary (
+  const unsigned int,
+  BlockVector<double>&) const;
 template
-void MGSmootherContinuous::set_zero_interior_boundary<BlockVector<float> > (const unsigned int,
-                                                        BlockVector<float>&) const;
+void MGSmootherContinuous<BlockVector<float> >::set_zero_interior_boundary (
+  const unsigned int,
+  BlockVector<float>&) const;
 
 
   

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.