#define __multigrid_H
/*---------------------------- multigrid.h ---------------------------*/
-#include <base/subscriptor.h>
+#include <base/smartpointer.h>
#include <lac/forward-declarations.h>
#include <basic/forward-declarations.h>
#include <lac/sparsematrix.h>
#include <vector>
class MGTransferBase;
-
+class MGSmoother;
/**
* Prolongation and restriction object.
*/
SmartPointer<MGTransferBase> transfer;
-
- /** Tranfer from dVector to
- * MGVector.
- *
- * This function copies data from a
- * dVector, that is, data on the
- * locally finest level, into the
- * corresponding levels of an
- * MGVector.
+
+ protected:
+ /**
+ * Number of pre-smoothing steps. Is used
+ * as a parameter to #pre_smooth#.
*/
- void copy_to_mg(vector<Vector<float> >& dst, const Vector<double>& src);
-
+ unsigned n_pre_smooth;
+
/**
- * Transfer from MGVector to
- * Vector<double>.
- *
- * Copies data from active portions
- * of an MGVector into the
- * respective positions of a
- * Vector<double>. All other entries of
- * #src# are zero.
+ * Number of post-smoothing steps Is used
+ * as a parameter to #post_smooth#.
*/
- void copy_from_mg(Vector<double>& dst, const vector<Vector<float> >& src);
-
+ unsigned n_post_smooth;
/**
* The actual v-cycle multigrid method.
* This function is called on the
* #coarse_grid_solution# and
* proceeds back up.
*/
- void level_mgstep(unsigned level);
-
- protected:
- /**
- * Number of pre-smoothing steps.
- */
- unsigned n_pre_smooth;
-
- /**
- * Number of post-smoothing steps
- */
- unsigned n_post_smooth;
- /**
- * The (pre-)smoothing algorithm.
- * This function is required to
- * perform #steps# iterations to
- * smoothen the residual $Ax-b$.
- *
- * When overloading this function
- * in derived classes, make sure
- * that smoothing only applies to
- * interior degrees of freedom of
- * the actual level.
- *
- */
- virtual void smooth(unsigned level,
- Vector<float>& x,
- const Vector<float>& b,
- unsigned steps) = 0;
-
- /**
- * The post-smoothing algorithm.
- * Defaults to #smooth#.
- */
- virtual void post_smooth(unsigned level,
- Vector<float>& dst,
- const Vector<float>& src,
- unsigned steps);
+ void level_mgstep(unsigned level, MGSmoother& smoother);
/**
* Apply residual operator on all
* standard implementation doing
* #10 * (n_pre_smooth +
* n_post_smooth)#
- * smoothing steps.
+ * smoothing steps of #pre_smooth#.
*/
virtual void coarse_grid_solution(unsigned l,
Vector<float>& dst,
const Vector<float> &src) const = 0;
};
+/**
+ * Implementation of multigrid with matrices.
+ * While #MultiGridBase# was only an abstract framework for the v-cycle,
+ * here we have the implementation of the pure virtual functions defined there.
+ * Furthermore, level information is obtained from a triangulation object.
+ *
+ * @author Wolfgang Bangerth, Guido Kanschat, 1999
+ */
+template<int dim, class Matrix>
+class MultiGrid
+{
+ public:
+ MultiGrid(const MGDoFHandler<dim>&);
+
+ /** Transfer from dVector to
+ * MGVector.
+ *
+ * This function copies data from a
+ * dVector, that is, data on the
+ * locally finest level, into the
+ * corresponding levels of an
+ * MGVector.
+ */
+ template<typename number>
+ void copy_to_mg(vector<Vector<float> >& dst,
+ const Vector<number>& src) const;
+
+ /**
+ * Transfer from multi-grid vector to
+ * normal vector.
+ *
+ * Copies data from active portions
+ * of an MGVector into the
+ * respective positions of a
+ * Vector<double>. All other entries of
+ * #src# are zero.
+ */
+ template<typename number>
+ void copy_from_mg(Vector<number>& dst,
+ const vector<Vector<float> >& src) const;
+ private:
+ /**
+ * Associated #MGDoFHandler#.
+ */
+ SmartPointer<MGDoFHandler<dim> > dofs;
+
+ /**
+ * Matrix structures for each level.
+ * These are generated by the constructor
+ * of #MultiGrid# and can be accessed
+ * later.
+ */
+ vector<SparseMatrixStruct> matrix_structures;
+
+ /**
+ * Matrices for each level.
+ * The matrices are prepared by
+ * the constructor of #MultiGrid# and can
+ * be accessed for assembling.
+ */
+ vector<SparseMatrix<float> > matrices;
+};
/**
#define __multigrid_H
/*---------------------------- multigrid.h ---------------------------*/
-#include <base/subscriptor.h>
+#include <base/smartpointer.h>
#include <lac/forward-declarations.h>
#include <basic/forward-declarations.h>
#include <lac/sparsematrix.h>
#include <vector>
class MGTransferBase;
-
+class MGSmoother;
/**
* Prolongation and restriction object.
*/
SmartPointer<MGTransferBase> transfer;
-
- /** Tranfer from dVector to
- * MGVector.
- *
- * This function copies data from a
- * dVector, that is, data on the
- * locally finest level, into the
- * corresponding levels of an
- * MGVector.
+
+ protected:
+ /**
+ * Number of pre-smoothing steps. Is used
+ * as a parameter to #pre_smooth#.
*/
- void copy_to_mg(vector<Vector<float> >& dst, const Vector<double>& src);
-
+ unsigned n_pre_smooth;
+
/**
- * Transfer from MGVector to
- * Vector<double>.
- *
- * Copies data from active portions
- * of an MGVector into the
- * respective positions of a
- * Vector<double>. All other entries of
- * #src# are zero.
+ * Number of post-smoothing steps Is used
+ * as a parameter to #post_smooth#.
*/
- void copy_from_mg(Vector<double>& dst, const vector<Vector<float> >& src);
-
+ unsigned n_post_smooth;
/**
* The actual v-cycle multigrid method.
* This function is called on the
* #coarse_grid_solution# and
* proceeds back up.
*/
- void level_mgstep(unsigned level);
-
- protected:
- /**
- * Number of pre-smoothing steps.
- */
- unsigned n_pre_smooth;
-
- /**
- * Number of post-smoothing steps
- */
- unsigned n_post_smooth;
- /**
- * The (pre-)smoothing algorithm.
- * This function is required to
- * perform #steps# iterations to
- * smoothen the residual $Ax-b$.
- *
- * When overloading this function
- * in derived classes, make sure
- * that smoothing only applies to
- * interior degrees of freedom of
- * the actual level.
- *
- */
- virtual void smooth(unsigned level,
- Vector<float>& x,
- const Vector<float>& b,
- unsigned steps) = 0;
-
- /**
- * The post-smoothing algorithm.
- * Defaults to #smooth#.
- */
- virtual void post_smooth(unsigned level,
- Vector<float>& dst,
- const Vector<float>& src,
- unsigned steps);
+ void level_mgstep(unsigned level, MGSmoother& smoother);
/**
* Apply residual operator on all
* standard implementation doing
* #10 * (n_pre_smooth +
* n_post_smooth)#
- * smoothing steps.
+ * smoothing steps of #pre_smooth#.
*/
virtual void coarse_grid_solution(unsigned l,
Vector<float>& dst,
const Vector<float> &src) const = 0;
};
+/**
+ * Implementation of multigrid with matrices.
+ * While #MultiGridBase# was only an abstract framework for the v-cycle,
+ * here we have the implementation of the pure virtual functions defined there.
+ * Furthermore, level information is obtained from a triangulation object.
+ *
+ * @author Wolfgang Bangerth, Guido Kanschat, 1999
+ */
+template<int dim, class Matrix>
+class MultiGrid
+{
+ public:
+ MultiGrid(const MGDoFHandler<dim>&);
+
+ /** Transfer from dVector to
+ * MGVector.
+ *
+ * This function copies data from a
+ * dVector, that is, data on the
+ * locally finest level, into the
+ * corresponding levels of an
+ * MGVector.
+ */
+ template<typename number>
+ void copy_to_mg(vector<Vector<float> >& dst,
+ const Vector<number>& src) const;
+
+ /**
+ * Transfer from multi-grid vector to
+ * normal vector.
+ *
+ * Copies data from active portions
+ * of an MGVector into the
+ * respective positions of a
+ * Vector<double>. All other entries of
+ * #src# are zero.
+ */
+ template<typename number>
+ void copy_from_mg(Vector<number>& dst,
+ const vector<Vector<float> >& src) const;
+ private:
+ /**
+ * Associated #MGDoFHandler#.
+ */
+ SmartPointer<MGDoFHandler<dim> > dofs;
+
+ /**
+ * Matrix structures for each level.
+ * These are generated by the constructor
+ * of #MultiGrid# and can be accessed
+ * later.
+ */
+ vector<SparseMatrixStruct> matrix_structures;
+
+ /**
+ * Matrices for each level.
+ * The matrices are prepared by
+ * the constructor of #MultiGrid# and can
+ * be accessed for assembling.
+ */
+ vector<SparseMatrix<float> > matrices;
+};
/**
unsigned maxlevel, unsigned minlevel,
unsigned pre_smooth, unsigned post_smooth)
:
+ d(maxlevel-minlevel),
+ s(maxlevel-minlevel),
maxlevel(maxlevel), minlevel(minlevel),
transfer(&transfer),
n_pre_smooth(pre_smooth), n_post_smooth(post_smooth)
{}
-
void
-MultiGridBase::level_mgstep(unsigned level)
+MultiGridBase::level_mgstep(unsigned level, MGSmoother& smoother)
{
if (level == minlevel)
{
}
// smoothing of the residual by modifying s
- smooth(level, d[level], s[level], n_pre_smooth);
+// smooth(level, d[level], s[level], n_pre_smooth);
// 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);
+ level_mgstep(level-1, smoother);
// do coarse grid correction
transfer->prolongate(level, s[level], s[level-1]);
// smoothing (modify s again)
- post_smooth(level, d[level], s[level], n_post_smooth);
+// post_smooth(level, d[level], s[level], n_post_smooth);
}
void
-MultiGridBase::post_smooth(unsigned level,
- Vector<float>& dst, const Vector<float>& src,
- unsigned steps)
-{
- smooth(level, dst, src, steps);
-}
-
-
-
-void
-MultiGridBase::coarse_grid_solution(unsigned level,
- Vector<float>& dst,
- const Vector<float>& src)
+MultiGridBase::coarse_grid_solution(unsigned /*level*/,
+ Vector<float>& /*dst*/,
+ const Vector<float>& /*src*/)
{
- smooth(level, dst, src, 10 * (n_pre_smooth + n_post_smooth));
+// smooth(level, dst, src, 10 * (n_pre_smooth + n_post_smooth));
}
unsigned maxlevel, unsigned minlevel,
unsigned pre_smooth, unsigned post_smooth)
:
+ d(maxlevel-minlevel),
+ s(maxlevel-minlevel),
maxlevel(maxlevel), minlevel(minlevel),
transfer(&transfer),
n_pre_smooth(pre_smooth), n_post_smooth(post_smooth)
{}
-
void
-MultiGridBase::level_mgstep(unsigned level)
+MultiGridBase::level_mgstep(unsigned level, MGSmoother& smoother)
{
if (level == minlevel)
{
}
// smoothing of the residual by modifying s
- smooth(level, d[level], s[level], n_pre_smooth);
+// smooth(level, d[level], s[level], n_pre_smooth);
// 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);
+ level_mgstep(level-1, smoother);
// do coarse grid correction
transfer->prolongate(level, s[level], s[level-1]);
// smoothing (modify s again)
- post_smooth(level, d[level], s[level], n_post_smooth);
+// post_smooth(level, d[level], s[level], n_post_smooth);
}
void
-MultiGridBase::post_smooth(unsigned level,
- Vector<float>& dst, const Vector<float>& src,
- unsigned steps)
-{
- smooth(level, dst, src, steps);
-}
-
-
-
-void
-MultiGridBase::coarse_grid_solution(unsigned level,
- Vector<float>& dst,
- const Vector<float>& src)
+MultiGridBase::coarse_grid_solution(unsigned /*level*/,
+ Vector<float>& /*dst*/,
+ const Vector<float>& /*src*/)
{
- smooth(level, dst, src, 10 * (n_pre_smooth + n_post_smooth));
+// smooth(level, dst, src, 10 * (n_pre_smooth + n_post_smooth));
}