]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
MGBase moved to LAC
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 28 Apr 1999 16:11:12 +0000 (16:11 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 28 Apr 1999 16:11:12 +0000 (16:11 +0000)
git-svn-id: https://svn.dealii.org/trunk@1220 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/multigrid/mg_base.h [new file with mode: 0644]
deal.II/deal.II/source/multigrid/mg_base.cc [new file with mode: 0644]
deal.II/lac/include/lac/mgbase.h [new file with mode: 0644]
deal.II/lac/source/mgbase.cc [new file with mode: 0644]

diff --git a/deal.II/deal.II/include/multigrid/mg_base.h b/deal.II/deal.II/include/multigrid/mg_base.h
new file mode 100644 (file)
index 0000000..e864500
--- /dev/null
@@ -0,0 +1,303 @@
+// $Id$
+
+
+#ifndef __lac_mg_base_h
+#define __lac_mg_base_h
+
+#include <base/subscriptor.h>
+#include <base/smartpointer.h>
+#include <lac/forward-declarations.h>
+#include <lac/vector.h>
+
+#include <vector>
+
+
+/**
+ * Abstract base class for coarse grid solvers.
+ * The interface of a function call operator is defined to execute coarse
+ * grid solution in a derived class.
+ */
+class MGCoarseGridSolver
+{
+  public:
+                                    /**
+                                     * Coarse grid solving method.
+                                     * This is only the interface for
+                                     * a function defined in a derived
+                                     * class like
+                                     * #MGCoarseGridLACIteration#.
+                                     *
+                                     * Remark that the information
+                                     * about the matrix is removed to
+                                     * that class.
+                                     */
+    virtual void operator() (unsigned int level, Vector<float>& dst,
+                            const Vector<float>& src) const = 0;
+};
+
+
+
+/**
+ * Abstract base class for multigrid smoothers.
+ * In fact, this class only provides the interface of the smoothing function.
+ * Using deal.II grid handling, #MGSmoother# is a good point to start.
+ *
+ * @author Wolfgang Bangerth, Guido Kanschat, 1999
+ */
+class MGSmootherBase
+  :
+  public Subscriptor 
+{  
+  public:
+                                    /**
+                                     * Virtual destructor.
+                                     */
+    virtual ~MGSmootherBase();
+    
+                                    /**
+                                     * Smoothen the residual of #u# on the given
+                                     * level. This function should keep the interior
+                                     * level boundary values, so you may want
+                                     * to call #set_zero_interior_boundary#
+                                     * in the derived class #MGSmoother#
+                                     * somewhere in your derived function,
+                                     * or another function doing similar
+                                     * things.
+                                     */
+    virtual void smooth (const unsigned int   level,
+                        Vector<float>       &u,
+                        const Vector<float> &rhs) const = 0;
+
+};
+
+/**
+ * Base class used to declare the operations needed by a concrete class
+ * implementing prolongation and restriction of vectors in the multigrid
+ * context. This class is an abstract one and has no implementations of
+ * possible algorithms for these operations.
+ *
+ * @author Wolfgang Bangerth, Guido Kanschat, 1999
+ */
+class MGTransferBase  :  public Subscriptor
+{
+  public:
+                                    /**
+                                     * Destructor. Does nothing here, but
+                                     * needs to be declared virtual anyway.
+                                     */
+    virtual ~MGTransferBase();
+
+                                    /**
+                                     * Prolongate a vector from level
+                                     * #to_level-1# to level #to_level#.
+                                     *
+                                     * #src# is assumed to be a vector with
+                                     * as many elements as there are degrees
+                                     * of freedom on the coarser level of
+                                     * the two involved levels, while #src#
+                                     * shall have as many elements as there
+                                     * are degrees of freedom on the finer
+                                     * level.
+                                     */
+    virtual void prolongate (const unsigned int   to_level,
+                            Vector<float>       &dst,
+                            const Vector<float> &src) const = 0;
+
+                                    /**
+                                     * Restrict a vector from level
+                                     * #from_level# to level
+                                     * #from_level-1#.
+                                     *
+                                     * #src# is assumed to be a vector with
+                                     * as many elements as there are degrees
+                                     * of freedom on the finer level of
+                                     * the two involved levels, while #src#
+                                     * shall have as many elements as there
+                                     * are degrees of freedom on the coarser
+                                     * level.
+                                     */
+    virtual void restrict (const unsigned int   from_level,
+                          Vector<float>       &dst,
+                          const Vector<float> &src) const = 0;
+};
+
+/**
+ * An array of matrices for each level.
+ * This class provides the functionality of #vector<MATRIX># combined
+ * with a subscriptor for smart pointers.
+ * @author Guido Kanschat, 1999
+ */
+template<class MATRIX>
+class MGMatrix
+  :
+  public Subscriptor,
+  public Vector<MATRIX>
+{
+  public:
+                                    /**
+                                     * Constructor allowing to initialize the number of levels.
+                                     */
+    MGMatrix(unsigned int n_levels);
+};
+
+
+/**
+ * Coarse grid solver using LAC iterative methods.
+ * This is a little wrapper, transforming a triplet of iterative
+ * solver, matrix and preconditioner into a coarse grid solver.
+ * @author Guido Kanschat, 1999
+ */
+template<class SOLVER, class MATRIX, class PRECOND>
+class MGCoarseGridLACIteration
+  :
+  public MGCoarseGridSolver
+{
+  public:
+                                    /**
+                                     * Constructor.
+                                     * Store solver, matrix and
+                                     * preconditioning method for later
+                                     * use.
+                                     */
+    MGCoarseGridLACIteration(const SOLVER&, const MATRIX&, const PRECOND&);
+    
+                                    /**
+                                     * Implementation of the abstract
+                                     * function.
+                                     * Calls the solver method with
+                                     * matrix, vectors and
+                                     * preconditioner.
+                                     */
+    virtual void operator() (unsigned int level, Vector<float>& dst,
+                            const Vector<float>& src) const;
+  private:
+                                    /**
+                                     * Reference to the solver.
+                                     */
+    const SOLVER& solver;
+                                    /**
+                                     * Reference to the matrix.
+                                     */
+    const MATRIX& matrix;
+    
+                                    /**
+                                     * Reference to the preconditioner.
+                                     */
+    const PRECOND& precondition;
+};
+
+/**
+ * Basic class for multigrid preconditioning.
+ *
+ * The functionality of the multigrid method is restricted to defect
+ * correction. It is <B>not</B> iterative and the start solution is
+ * always zero. Since by this <I>u<SUP>E</SUP><SUB>l</SUB></I> and
+ * <I>u<SUP>A</SUP><SUB>l</SUB></I> (see report on multigrid) are
+ * always zero, restriction is simplified a lot and maybe even the
+ * seam condition on grids is oblivious. Still, I am not sure that
+ * these restrictions on the class might cause numerical
+ * inefficiencies.
+ *
+ * The function #precondition# is the actual multigrid method and
+ * makes use of several operations to be implemented in derived
+ * classes. It takes a defect in #src# and the result is the multigrid
+ * preconditioned defect in #dst#.
+ *
+ * @author Guido Kanschat, 1999
+ * */
+class MGBase
+{
+    MGBase(const MGBase&);
+    const MGBase& operator=(const MGBase&);
+    
+                                    /**
+                                     * Auxiliary vector, defect.
+                                     */
+    vector<Vector<float> > d;
+
+                                    /**
+                                     * Auxiliary vector, solution.
+                                     */
+    vector<Vector<float> > s;
+
+                                    /**
+                                     * Auxiliary vector.
+                                     */
+    Vector<float> t;
+    
+                                    /**
+                                     * Prolongation and restriction object.
+                                     */
+    SmartPointer<const MGTransferBase> transfer;
+    
+  protected:
+                                    /**
+                                     * Highest level of cells.
+                                     */
+    unsigned int maxlevel;
+
+                                    /**
+                                     * Level for coarse grid solution.
+                                     */
+    unsigned int minlevel;
+
+                                    /**
+                                     * The actual v-cycle multigrid method.
+                                     * This function is called on the
+                                     * highest level and recursively
+                                     * invokes itself down to the
+                                     * coarsest. There, it calls
+                                     * #coarse_grid_solver# and
+                                     * proceeds back up.
+                                     */
+    void level_mgstep(unsigned int level,
+                     const MGSmootherBase& pre_smooth,
+                     const MGSmootherBase& post_smooth,
+                     const MGCoarseGridSolver& cgs);  
+
+                                    /**
+                                     * Apply residual operator on all
+                                     * cells of a level.
+                                     * This is implemented in a
+                                     * derived class.
+                                     */
+    virtual void level_residual(unsigned int level,
+                            Vector<float>& dst,
+                            const Vector<float>& src,
+                            const Vector<float>& rhs) = 0;  
+
+  
+  public:
+                                    /**
+                                     * Constructor, subject to change.
+                                     */
+    MGBase(const MGTransferBase& transfer,
+                 unsigned int maxlevel, unsigned int minlevel);
+    virtual ~MGBase();
+    
+};
+
+
+template<class SOLVER, class MATRIX, class PRECOND>
+MGCoarseGridLACIteration<SOLVER, MATRIX, PRECOND>
+::MGCoarseGridLACIteration(const SOLVER& s, const MATRIX& m, const PRECOND& p)
+               :
+               solver(s), matrix(m), precondition(p)
+{}
+
+template<class SOLVER, class MATRIX, class PRECOND>
+void
+MGCoarseGridLACIteration<SOLVER, MATRIX, PRECOND>::operator()
+(unsigned int level, Vector<float>& dst, const Vector<float>& src) const
+{
+  solver.solve(matrix, dst, src, precondition);
+}
+
+template<class MATRIX>
+MGMatrix<MATRIX>::MGMatrix(unsigned int nlevels)
+               :
+               vector<MATRIX>(nlevels)
+{}
+
+
+#endif
diff --git a/deal.II/deal.II/source/multigrid/mg_base.cc b/deal.II/deal.II/source/multigrid/mg_base.cc
new file mode 100644 (file)
index 0000000..1ccc525
--- /dev/null
@@ -0,0 +1,47 @@
+// $Id$
+
+#include <lac/mgbase.h>
+
+MGBase::~MGBase () 
+{};
+
+
+
+MGBase::MGBase(const MGTransferBase& transfer,
+                            unsigned maxlevel, unsigned minlevel)
+               :
+               d(maxlevel-minlevel),
+               s(maxlevel-minlevel),
+               transfer(&transfer),
+               maxlevel(maxlevel), minlevel(minlevel)
+{}
+
+void
+MGBase::level_mgstep(unsigned level,
+                    const MGSmootherBase& pre_smooth,
+                    const MGSmootherBase& post_smooth,
+                    const MGCoarseGridSolver& coarse_grid_solver)
+{
+  if (level == minlevel)
+    {
+      coarse_grid_solver(level, s[level], d[level]);
+      return;
+    }
+  
+                                  // smoothing of the residual by modifying s
+  pre_smooth.smooth(level, s[level], d[level]);
+
+                                  // t = d-As
+  level_residual(level, t, s[level], d[level]);
+  
+                                  // make t rhs of lower level
+  transfer->restrict(level, t, d[level-1]);
+                                  // do recursion
+  level_mgstep(level-1, pre_smooth, post_smooth, coarse_grid_solver);
+                                  // do coarse grid correction
+  transfer->prolongate(level, s[level], s[level-1]);
+
+                                  // smoothing (modify s again)
+  post_smooth.smooth(level, s[level], d[level]);
+}
+
diff --git a/deal.II/lac/include/lac/mgbase.h b/deal.II/lac/include/lac/mgbase.h
new file mode 100644 (file)
index 0000000..f614a30
--- /dev/null
@@ -0,0 +1,303 @@
+// $Id$
+
+
+#ifndef __lac_mgbase_h
+#define __lac_mgbase_h
+
+#include <base/subscriptor.h>
+#include <base/smartpointer.h>
+#include <lac/forward-declarations.h>
+#include <lac/vector.h>
+
+#include <vector>
+
+
+/**
+ * Abstract base class for coarse grid solvers.
+ * The interface of a function call operator is defined to execute coarse
+ * grid solution in a derived class.
+ */
+class MGCoarseGridSolver
+{
+  public:
+                                    /**
+                                     * Coarse grid solving method.
+                                     * This is only the interface for
+                                     * a function defined in a derived
+                                     * class like
+                                     * #MGCoarseGridLACIteration#.
+                                     *
+                                     * Remark that the information
+                                     * about the matrix is removed to
+                                     * that class.
+                                     */
+    virtual void operator() (unsigned int level, Vector<float>& dst,
+                            const Vector<float>& src) const = 0;
+};
+
+
+
+/**
+ * Abstract base class for multigrid smoothers.
+ * In fact, this class only provides the interface of the smoothing function.
+ * Using deal.II grid handling, #MGSmoother# is a good point to start.
+ *
+ * @author Wolfgang Bangerth, Guido Kanschat, 1999
+ */
+class MGSmootherBase
+  :
+  public Subscriptor 
+{  
+  public:
+                                    /**
+                                     * Virtual destructor.
+                                     */
+    virtual ~MGSmootherBase();
+    
+                                    /**
+                                     * Smoothen the residual of #u# on the given
+                                     * level. This function should keep the interior
+                                     * level boundary values, so you may want
+                                     * to call #set_zero_interior_boundary#
+                                     * in the derived class #MGSmoother#
+                                     * somewhere in your derived function,
+                                     * or another function doing similar
+                                     * things.
+                                     */
+    virtual void smooth (const unsigned int   level,
+                        Vector<float>       &u,
+                        const Vector<float> &rhs) const = 0;
+
+};
+
+/**
+ * Base class used to declare the operations needed by a concrete class
+ * implementing prolongation and restriction of vectors in the multigrid
+ * context. This class is an abstract one and has no implementations of
+ * possible algorithms for these operations.
+ *
+ * @author Wolfgang Bangerth, Guido Kanschat, 1999
+ */
+class MGTransferBase  :  public Subscriptor
+{
+  public:
+                                    /**
+                                     * Destructor. Does nothing here, but
+                                     * needs to be declared virtual anyway.
+                                     */
+    virtual ~MGTransferBase();
+
+                                    /**
+                                     * Prolongate a vector from level
+                                     * #to_level-1# to level #to_level#.
+                                     *
+                                     * #src# is assumed to be a vector with
+                                     * as many elements as there are degrees
+                                     * of freedom on the coarser level of
+                                     * the two involved levels, while #src#
+                                     * shall have as many elements as there
+                                     * are degrees of freedom on the finer
+                                     * level.
+                                     */
+    virtual void prolongate (const unsigned int   to_level,
+                            Vector<float>       &dst,
+                            const Vector<float> &src) const = 0;
+
+                                    /**
+                                     * Restrict a vector from level
+                                     * #from_level# to level
+                                     * #from_level-1#.
+                                     *
+                                     * #src# is assumed to be a vector with
+                                     * as many elements as there are degrees
+                                     * of freedom on the finer level of
+                                     * the two involved levels, while #src#
+                                     * shall have as many elements as there
+                                     * are degrees of freedom on the coarser
+                                     * level.
+                                     */
+    virtual void restrict (const unsigned int   from_level,
+                          Vector<float>       &dst,
+                          const Vector<float> &src) const = 0;
+};
+
+/**
+ * An array of matrices for each level.
+ * This class provides the functionality of #vector<MATRIX># combined
+ * with a subscriptor for smart pointers.
+ * @author Guido Kanschat, 1999
+ */
+template<class MATRIX>
+class MGMatrix
+  :
+  public Subscriptor,
+  public Vector<MATRIX>
+{
+  public:
+                                    /**
+                                     * Constructor allowing to initialize the number of levels.
+                                     */
+    MGMatrix(unsigned int n_levels);
+};
+
+
+/**
+ * Coarse grid solver using LAC iterative methods.
+ * This is a little wrapper, transforming a triplet of iterative
+ * solver, matrix and preconditioner into a coarse grid solver.
+ * @author Guido Kanschat, 1999
+ */
+template<class SOLVER, class MATRIX, class PRECOND>
+class MGCoarseGridLACIteration
+  :
+  public MGCoarseGridSolver
+{
+  public:
+                                    /**
+                                     * Constructor.
+                                     * Store solver, matrix and
+                                     * preconditioning method for later
+                                     * use.
+                                     */
+    MGCoarseGridLACIteration(const SOLVER&, const MATRIX&, const PRECOND&);
+    
+                                    /**
+                                     * Implementation of the abstract
+                                     * function.
+                                     * Calls the solver method with
+                                     * matrix, vectors and
+                                     * preconditioner.
+                                     */
+    virtual void operator() (unsigned int level, Vector<float>& dst,
+                            const Vector<float>& src) const;
+  private:
+                                    /**
+                                     * Reference to the solver.
+                                     */
+    const SOLVER& solver;
+                                    /**
+                                     * Reference to the matrix.
+                                     */
+    const MATRIX& matrix;
+    
+                                    /**
+                                     * Reference to the preconditioner.
+                                     */
+    const PRECOND& precondition;
+};
+
+/**
+ * Basic class for multigrid preconditioning.
+ *
+ * The functionality of the multigrid method is restricted to defect
+ * correction. It is <B>not</B> iterative and the start solution is
+ * always zero. Since by this <I>u<SUP>E</SUP><SUB>l</SUB></I> and
+ * <I>u<SUP>A</SUP><SUB>l</SUB></I> (see report on multigrid) are
+ * always zero, restriction is simplified a lot and maybe even the
+ * seam condition on grids is oblivious. Still, I am not sure that
+ * these restrictions on the class might cause numerical
+ * inefficiencies.
+ *
+ * The function #precondition# is the actual multigrid method and
+ * makes use of several operations to be implemented in derived
+ * classes. It takes a defect in #src# and the result is the multigrid
+ * preconditioned defect in #dst#.
+ *
+ * @author Guido Kanschat, 1999
+ * */
+class MGBase
+{
+    MGBase(const MGBase&);
+    const MGBase& operator=(const MGBase&);
+    
+                                    /**
+                                     * Auxiliary vector, defect.
+                                     */
+    vector<Vector<float> > d;
+
+                                    /**
+                                     * Auxiliary vector, solution.
+                                     */
+    vector<Vector<float> > s;
+
+                                    /**
+                                     * Auxiliary vector.
+                                     */
+    Vector<float> t;
+    
+                                    /**
+                                     * Prolongation and restriction object.
+                                     */
+    SmartPointer<const MGTransferBase> transfer;
+    
+  protected:
+                                    /**
+                                     * Highest level of cells.
+                                     */
+    unsigned int maxlevel;
+
+                                    /**
+                                     * Level for coarse grid solution.
+                                     */
+    unsigned int minlevel;
+
+                                    /**
+                                     * The actual v-cycle multigrid method.
+                                     * This function is called on the
+                                     * highest level and recursively
+                                     * invokes itself down to the
+                                     * coarsest. There, it calls
+                                     * #coarse_grid_solver# and
+                                     * proceeds back up.
+                                     */
+    void level_mgstep(unsigned int level,
+                     const MGSmootherBase& pre_smooth,
+                     const MGSmootherBase& post_smooth,
+                     const MGCoarseGridSolver& cgs);  
+
+                                    /**
+                                     * Apply residual operator on all
+                                     * cells of a level.
+                                     * This is implemented in a
+                                     * derived class.
+                                     */
+    virtual void level_residual(unsigned int level,
+                            Vector<float>& dst,
+                            const Vector<float>& src,
+                            const Vector<float>& rhs) = 0;  
+
+  
+  public:
+                                    /**
+                                     * Constructor, subject to change.
+                                     */
+    MGBase(const MGTransferBase& transfer,
+                 unsigned int maxlevel, unsigned int minlevel);
+    virtual ~MGBase();
+    
+};
+
+
+template<class SOLVER, class MATRIX, class PRECOND>
+MGCoarseGridLACIteration<SOLVER, MATRIX, PRECOND>
+::MGCoarseGridLACIteration(const SOLVER& s, const MATRIX& m, const PRECOND& p)
+               :
+               solver(s), matrix(m), precondition(p)
+{}
+
+template<class SOLVER, class MATRIX, class PRECOND>
+void
+MGCoarseGridLACIteration<SOLVER, MATRIX, PRECOND>::operator()
+(unsigned int level, Vector<float>& dst, const Vector<float>& src) const
+{
+  solver.solve(matrix, dst, src, precondition);
+}
+
+template<class MATRIX>
+MGMatrix<MATRIX>::MGMatrix(unsigned int nlevels)
+               :
+               vector<MATRIX>(nlevels)
+{}
+
+
+#endif
diff --git a/deal.II/lac/source/mgbase.cc b/deal.II/lac/source/mgbase.cc
new file mode 100644 (file)
index 0000000..1ccc525
--- /dev/null
@@ -0,0 +1,47 @@
+// $Id$
+
+#include <lac/mgbase.h>
+
+MGBase::~MGBase () 
+{};
+
+
+
+MGBase::MGBase(const MGTransferBase& transfer,
+                            unsigned maxlevel, unsigned minlevel)
+               :
+               d(maxlevel-minlevel),
+               s(maxlevel-minlevel),
+               transfer(&transfer),
+               maxlevel(maxlevel), minlevel(minlevel)
+{}
+
+void
+MGBase::level_mgstep(unsigned level,
+                    const MGSmootherBase& pre_smooth,
+                    const MGSmootherBase& post_smooth,
+                    const MGCoarseGridSolver& coarse_grid_solver)
+{
+  if (level == minlevel)
+    {
+      coarse_grid_solver(level, s[level], d[level]);
+      return;
+    }
+  
+                                  // smoothing of the residual by modifying s
+  pre_smooth.smooth(level, s[level], d[level]);
+
+                                  // t = d-As
+  level_residual(level, t, s[level], d[level]);
+  
+                                  // make t rhs of lower level
+  transfer->restrict(level, t, d[level-1]);
+                                  // do recursion
+  level_mgstep(level-1, pre_smooth, post_smooth, coarse_grid_solver);
+                                  // do coarse grid correction
+  transfer->prolongate(level, s[level], s[level-1]);
+
+                                  // smoothing (modify s again)
+  post_smooth.smooth(level, s[level], d[level]);
+}
+

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.