#ifndef __deal2__mg_base_h
#define __deal2__mg_base_h
+/*
+ * This file contains MGLevelObject and some abstract base classes
+ * used by Multigrid.
+ */
#include <base/config.h>
#include <base/subscriptor.h>
* @author Guido Kanschat, 2002
*/
template <class VECTOR>
-class MGCoarseGrid : public Subscriptor
+class MGCoarseGridBase : public Subscriptor
{
public:
/**
* Virtual destructor.
*/
- virtual ~MGCoarseGrid ();
+ virtual ~MGCoarseGridBase ();
/**
* Solver method implemented by
* @author Wolfgang Bangerth, Guido Kanschat, 1999, 2002
*/
template <class VECTOR>
- class MGTransfer : public Subscriptor
+ class MGTransferBase : public Subscriptor
{
public:
/**
* Destructor. Does nothing here, but
* needs to be declared virtual anyway.
*/
- virtual ~MGTransfer();
+ virtual ~MGTransferBase();
/**
* Prolongate a vector from level
* @author Guido Kanschat, 2002
*/
template <class VECTOR>
-class MGSmoother : public Subscriptor
+class MGSmootherBase : public Subscriptor
{
public:
/**
* Virtual destructor.
*/
- virtual ~MGSmoother();
+ virtual ~MGSmootherBase();
/**
* Smoothing function. This is the
* @author Guido Kanschat, 1999, Ralf Hartmann, 2002.
*/
template<class SOLVER, class MATRIX, class PRECOND, class VECTOR = Vector<double> >
-class MGCoarseGridLACIteration : public MGCoarseGrid<VECTOR>
+class MGCoarseGridLACIteration : public MGCoarseGridBase<VECTOR>
{
public:
/**
template <int dim> class MGDoFHandler;
+/*
+ * MGSmootherBase is defined in mg_base.h
+ */
+
/**
* Smoother doing nothing. This class is not useful for many applications other
* than for testing some multigrid procedures. Also some applications might
* @author Guido Kanschat, 1999, 2002
*/
template <class VECTOR>
-class MGSmootherIdentity : public MGSmoother<VECTOR>
+class MGSmootherIdentity : public MGSmootherBase<VECTOR>
{
public:
/**
* @author Wolfgang Bangerth, Guido Kanschat, 1999, 2002
*/
template <class VECTOR>
-class MGSmootherContinuous : public MGSmoother<VECTOR>
+class MGSmootherContinuous : public MGSmootherBase<VECTOR>
{
private:
/**
* @author Guido Kanschat, 2003
*/
template<class MATRIX, class RELAX, class VECTOR>
-class MGSmootherRelaxation : public MGSmoother<VECTOR>
+class MGSmootherRelaxation : public MGSmootherBase<VECTOR>
{
public:
/**
template <int dim> class MGDoFHandler;
+/*
+ * MGTransferBase is defined in mg_base.h
+ */
/**
* Implementation of the @p{MGTransferBase} interface for which the transfer
* @author Wolfgang Bangerth, Guido Kanschat, 1999, 2000, 2001, 2002
*/
template<typename number>
-class MGTransferPrebuilt : public MGTransfer<Vector<number> >
+class MGTransferPrebuilt : public MGTransferBase<Vector<number> >
{
public:
/**
* @author Guido Kanschat, 2001, 2002
*/
template <typename number>
-class MGTransferBlock : public MGTransfer<BlockVector<number> >,
+class MGTransferBlock : public MGTransferBase<BlockVector<number> >,
private MGTransferBlockBase
{
public:
* @author Guido Kanschat, 2001, 2002
*/
template <typename number>
-class MGTransferSelect : public MGTransfer<Vector<number> >,
+class MGTransferSelect : public MGTransferBase<Vector<number> >,
private MGTransferBlockBase
{
public:
template <int dim>
Multigrid(const MGDoFHandler<dim>& mg_dof_handler,
const MGMatrixBase<VECTOR>& matrix,
- const MGCoarseGrid<VECTOR>& coarse,
- const MGTransfer<VECTOR>& transfer,
- const MGSmoother<VECTOR>& pre_smooth,
- const MGSmoother<VECTOR>& post_smooth,
+ const MGCoarseGridBase<VECTOR>& coarse,
+ const MGTransferBase<VECTOR>& transfer,
+ const MGSmootherBase<VECTOR>& pre_smooth,
+ const MGSmootherBase<VECTOR>& post_smooth,
const unsigned int minlevel = 0,
const unsigned int maxlevel = 1000000);
/**
* The matrix for each level.
*/
- SmartPointer<const MGCoarseGrid<VECTOR> > coarse;
+ SmartPointer<const MGCoarseGridBase<VECTOR> > coarse;
/**
* Object for grid tranfer.
*/
- SmartPointer<const MGTransfer<VECTOR> > transfer;
+ SmartPointer<const MGTransferBase<VECTOR> > transfer;
/**
* The pre-smoothing object.
*/
- SmartPointer<const MGSmoother<VECTOR> > pre_smooth;
+ SmartPointer<const MGSmootherBase<VECTOR> > pre_smooth;
/**
* The post-smoothing object.
*/
- SmartPointer<const MGSmoother<VECTOR> > post_smooth;
+ SmartPointer<const MGSmootherBase<VECTOR> > post_smooth;
/**
* Edge matrix from fine to coarse.
template <int dim>
Multigrid<VECTOR>::Multigrid (const MGDoFHandler<dim>& mg_dof_handler,
const MGMatrixBase<VECTOR>& matrix,
- const MGCoarseGrid<VECTOR>& coarse,
- const MGTransfer<VECTOR>& transfer,
- const MGSmoother<VECTOR>& pre_smooth,
- const MGSmoother<VECTOR>& post_smooth,
+ const MGCoarseGridBase<VECTOR>& coarse,
+ const MGTransferBase<VECTOR>& transfer,
+ const MGSmootherBase<VECTOR>& pre_smooth,
+ const MGSmootherBase<VECTOR>& post_smooth,
const unsigned int min_level,
const unsigned int max_level)
:
template <class VECTOR>
-MGTransfer<VECTOR>::~MGTransfer()
+MGTransferBase<VECTOR>::~MGTransferBase()
{}
template <class VECTOR>
-MGCoarseGrid<VECTOR>::~MGCoarseGrid()
+MGSmootherBase<VECTOR>::~MGSmootherBase()
+{}
+
+
+template <class VECTOR>
+MGCoarseGridBase<VECTOR>::~MGCoarseGridBase()
{}
// Explicit instantiations
-template class MGTransfer<Vector<double> >;
-template class MGTransfer<Vector<float> >;
-template class MGTransfer<BlockVector<double> >;
-template class MGTransfer<BlockVector<float> >;
+template class MGTransferBase<Vector<double> >;
+template class MGTransferBase<Vector<float> >;
+template class MGTransferBase<BlockVector<double> >;
+template class MGTransferBase<BlockVector<float> >;
template class MGMatrixBase<Vector<double> >;
template class MGMatrixBase<Vector<float> >;
template class MGMatrixBase<BlockVector<double> >;
template class MGMatrixBase<BlockVector<float> >;
-template class MGCoarseGrid<Vector<double> >;
-template class MGCoarseGrid<Vector<float> >;
-template class MGCoarseGrid<BlockVector<double> >;
-template class MGCoarseGrid<BlockVector<float> >;
+template MGSmootherBase<Vector<float> >;
+template MGSmootherBase<Vector<double> >;
+template MGSmootherBase<BlockVector<float> >;
+template MGSmootherBase<BlockVector<double> >;
+
+template class MGCoarseGridBase<Vector<double> >;
+template class MGCoarseGridBase<Vector<float> >;
+template class MGCoarseGridBase<BlockVector<double> >;
+template class MGCoarseGridBase<BlockVector<float> >;
#include <algorithm>
-template <class VECTOR>
-MGSmoother<VECTOR>::~MGSmoother()
-{}
-
-
//////////////////////////////////////////////////////////////////////
// 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