+++ /dev/null
-//---------------------------- mg_base.h ---------------------------
-// $Id$
-// Version: $Name$
-//
-// Copyright (C) 1998, 1999, 2000 by the deal.II authors
-//
-// This file is subject to QPL and may not be distributed
-// without copyright and license information. Please refer
-// to the file deal.II/doc/license.html for the text and
-// further information on this license.
-//
-//---------------------------- mg_base.h ---------------------------
-#ifndef __deal2__mg_base_h
-#define __deal2__mg_base_h
-
-
-/*---------------------------- mgbase.h ---------------------------*/
-
-
-#include <base/subscriptor.h>
-#include <base/smartpointer.h>
-#include <lac/vector.h>
-
-#include <lac/forward_declarations.h>
-#include <grid/forward_declarations.h>
-
-#include <vector>
-
-
-/**
- * An array with an object for each level. The purpose of this class
- * is mostly to allow access by level number, even if the lower levels
- * are not used and therefore have no object at all; this is done by
- * simply shifting the given index by the minimum level we have
- * stored.
- *
- * In most cases, the objects which are stored on each levels, are
- * either matrices or vectors.
- *
- * @author Wolfgang Bangerth, Guido Kanschat, 1999
- */
-template<class Object>
-class MGLevelObject : public Subscriptor
-{
- public:
- /**
- * Constructor allowing to
- * initialize the number of
- * levels. By default, the object
- * is created empty.
- */
- MGLevelObject (const unsigned int minlevel = 0,
- const unsigned int maxlevel = 0);
-
- /**
- * Access object on level #level#.
- */
- Object & operator[] (const unsigned int level);
-
- /**
- * Access object on level
- * #level#. Constant version.
- */
- const Object & operator[] (const unsigned int level) const;
-
- /**
- * Delete all previous contents
- * of this object and reset its
- * size according to the values
- * of #new_minlevel# and
- * #new_maxlevel#.
- */
- void resize (const unsigned int new_minlevel,
- const unsigned int new_maxlevel);
-
- /**
- * Call #clear# on all objects
- * stored by this object. This
- * function is only implemented
- * for some #Object# classes,
- * most notably for vectors and
- * matrices. Note that if
- * #Object==Vector<T>#, #clear#
- * will set all entries to zero,
- * while if
- * #Object==std::vector<T>#,
- * #clear# deletes the elements
- * of the vectors. This class
- * might therefore not be useful
- * for STL vectors.
- */
- void clear();
-
- //////
- unsigned int get_minlevel () const;
- unsigned int get_maxlevel () const;
-
- private:
- /**
- * Level of first component.
- */
- unsigned int minlevel;
-
- /**
- * Array of the objects to be held.
- */
- vector<Object> objects;
-};
-
-
-/**
- * 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.
- *
- * @author Guido Kanschat, 1999
- */
-class MGCoarseGridSolver : public Subscriptor
-{
- public:
- /**
- * Virtual destructor. Does
- * nothing in particular, but
- * needs to be declared anyhow.
- */
- virtual ~MGCoarseGridSolver();
-
- /**
- * Coarse grid solving method.
- * This is only the interface for
- * a function defined in a derived
- * class like
- * #MGCoarseGridLACIteration#.
- *
- * Note that the information
- * about the matrix is removed to
- * that class.
- */
- virtual void operator() (const unsigned int level,
- Vector<double> &dst,
- const Vector<double> &src) const = 0;
-};
-
-
-/**
- * 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.
- *
- * The type of the matrix (i.e. the template parameter #MATRIX#)
- * should be derived from #Subscriptor# to allow for the use of a
- * smart pointer to it.
- *
- * @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 (SOLVER &,
- const MATRIX &,
- const PRECOND &);
-
- /**
- * Implementation of the abstract
- * function.
- * Calls the solver method with
- * matrix, vectors and
- * preconditioner.
- */
- virtual void operator() (const unsigned int level,
- Vector<double> &dst,
- const Vector<double> &src) const;
- private:
- /**
- * Reference to the solver.
- */
- SOLVER& solver;
-
- /**
- * Reference to the matrix.
- */
- const SmartPointer<const MATRIX> matrix;
-
- /**
- * Reference to the preconditioner.
- */
- const PRECOND& precondition;
-};
-
-
-/**
- * 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#. The previous
- * content of #dst# is
- * overwritten.
- *
- * #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<double> &dst,
- const Vector<double> &src) const = 0;
-
- /**
- * Restrict a vector from level
- * #from_level# to level
- * #from_level-1# and add this
- * restriction to
- * #dst#. Obviously, if the
- * refined region on level
- * #from_level# is smaller than
- * that on level #from_level-1#,
- * some degrees of freedom in
- * #dst# are not covered and will
- * not be altered. For the other
- * degress of freedom, the result
- * of the restriction is added.
- *
- * #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_and_add (const unsigned int from_level,
- Vector<double> &dst,
- const Vector<double> &src) const = 0;
-};
-
-
-/**
- * Basic class for preconditioning by multigrid.
- *
- * The functionality of the multigrid method is restricted to defect
- * correction. It is {\bf not} iterative and the start solution is
- * always zero. Since by this $u^E_l$ and $u^A_l$ (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 : public Subscriptor
-{
- private:
- /**
- * Copy constructor. Made private
- * to prevent use.
- */
- MGBase(const MGBase&);
-
- /**
- * Copy operator. Made private
- * to prevent use.
- */
- const MGBase& operator=(const MGBase&);
-
- public:
- /**
- * Constructor, subject to change.
- */
- MGBase (const MGTransferBase &transfer,
- const unsigned int minlevel,
- const unsigned int maxlevel);
-
- /**
- * Virtual destructor.
- */
- virtual ~MGBase();
-
- /**
- * Execute one step of the
- * v-cycle algorithm. This
- * function assumes, that the
- * vector #d# is properly filled
- * with the residual in the outer
- * defect correction
- * scheme. After execution of
- * #vcycle()#, the result is in
- * the vector #s#. We propose to
- * write functions #copy_from_mg#
- * and #copy_to_mg# in derived
- * classes of #MGBase# to access
- * #d# and #s#.
- *
- * The actual work for this
- * function is done in
- * #level_mgstep#.
- */
- void vcycle(const MGSmootherBase &pre_smooth,
- const MGSmootherBase &post_smooth,
- const MGCoarseGridSolver &cgs);
-
- /**
- * Exception.
- */
- DeclException2(ExcSwitchedLevels, int, int,
- << "minlevel and maxlevel switched, should be: "
- << arg1 << "<=" << arg2);
-
- protected:
- /**
- * Highest level of cells.
- */
- unsigned int maxlevel;
-
- /**
- * Level for coarse grid solution.
- */
- unsigned int minlevel;
-
- /**
- * Auxiliary vector.
- */
- MGLevelObject<Vector<double> > defect;
-
- /**
- * Auxiliary vector.
- */
- MGLevelObject<Vector<double> > solution;
-
- /**
- * Prolongation and restriction object.
- */
- SmartPointer<const MGTransferBase> transfer;
-
- /**
- * The actual v-cycle multigrid method.
- * #level# is the level the
- * function should work on. It
- * will usually be called for the
- * highest level from outside,
- * but will then call itself
- * recursively for #level-1#,
- * unless we are on #minlevel#
- * where instead of recursing
- * deeper, the coarse grid solver
- * is used to solve the matrix of
- * this level.
- */
- void level_mgstep (const unsigned int level,
- const MGSmootherBase &pre_smooth,
- const MGSmootherBase &post_smooth,
- const MGCoarseGridSolver &cgs);
-
- /**
- * Apply negative #vmult# on all
- * cells of a level.
- * This is implemented in a
- * derived class.
- */
- virtual void level_vmult (const unsigned int level,
- Vector<double> &dst,
- const Vector<double> &src,
- const Vector<double> &rhs) = 0;
-
- private:
- /**
- * Auxiliary vector.
- */
- Vector<double> t;
-};
-
-
-/**
- * Multi-level preconditioner.
- * Here, we collect all information needed for multi-level preconditioning
- * and provide the standard interface for LAC iterative methods.
- *
- * The template parameter class #MG# is required to inherit #MGBase#.
- * Furthermore, it needs functions #void copy_to_mg(const VECTOR&)#
- * to store #src# in the right hand side of the multi-level method and
- * #void copy_from_mg(VECTOR&)# to store the result of the v-cycle in #dst#.
- *
- * @author Guido Kanschat, 1999
- */
-template<class MG, class VECTOR = Vector<double> >
-class PreconditionMG
-{
- public:
- /**
- * Constructor.
- * Arguments are the multigrid object,
- * pre-smoother, post-smoother and
- * coarse grid solver.
- */
- PreconditionMG(MG &mg,
- const MGSmootherBase &pre,
- const MGSmootherBase &post,
- const MGCoarseGridSolver &coarse);
-
- /**
- * Preconditioning operator,
- * calling the #vcycle# function
- * of the #MG# object passed to
- * the constructor.
- *
- * This is the operator used by
- * LAC iterative solvers.
- */
- void operator() (VECTOR &dst,
- const VECTOR &src) const;
-
- private:
- /**
- * The multigrid object.
- */
- SmartPointer<MG> multigrid;
-
- /**
- * The pre-smoothing object.
- */
- SmartPointer<const MGSmootherBase> pre;
-
- /**
- * The post-smoothing object.
- */
- SmartPointer<const MGSmootherBase> post;
-
- /**
- * The coarse grid solver.
- */
- SmartPointer<const MGCoarseGridSolver> coarse;
-};
-
-
-/* ------------------------------------------------------------------- */
-
-
-template<class Object>
-MGLevelObject<Object>::MGLevelObject(const unsigned int min,
- const unsigned int max)
- :
- minlevel(0)
-{
- resize (min, max);
-};
-
-
-template<class Object>
-Object &
-MGLevelObject<Object>::operator[] (const unsigned int i)
-{
- Assert((i>=minlevel) && (i<minlevel+objects.size()),
- ExcIndexRange (i, minlevel, minlevel+objects.size()));
- return objects[i-minlevel];
-}
-
-
-template<class Object>
-const Object &
-MGLevelObject<Object>::operator[] (const unsigned int i) const
-{
- Assert((i>=minlevel) && (i<minlevel+objects.size()),
- ExcIndexRange (i, minlevel, minlevel+objects.size()));
- return objects[i-minlevel];
-}
-
-
-template<class Object>
-void
-MGLevelObject<Object>::resize (const unsigned int new_minlevel,
- const unsigned int new_maxlevel)
-{
- Assert (new_minlevel <= new_maxlevel, ExcInternalError());
- objects.clear ();
-
- minlevel = new_minlevel;
- objects.resize (new_maxlevel - new_minlevel + 1);
-};
-
-
-template<class Object>
-void
-MGLevelObject<Object>::clear ()
-{
- typename vector<Object>::iterator v;
- for (v = objects.begin(); v != objects.end(); ++v)
- v->clear();
-};
-
-
-template<class Object>
-unsigned int
-MGLevelObject<Object>::get_minlevel () const
-{
- return minlevel;
-};
-
-
-template<class Object>
-unsigned int
-MGLevelObject<Object>::get_maxlevel () const
-{
- return minlevel + objects.size() - 1;
-};
-
-
-/* ------------------ Functions for MGCoarseGridLACIteration ------------ */
-
-
-template<class SOLVER, class MATRIX, class PRECOND>
-MGCoarseGridLACIteration<SOLVER,MATRIX,PRECOND>::MGCoarseGridLACIteration(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() (const unsigned int /* level */,
- Vector<double> &dst,
- const Vector<double> &src) const
-{
- solver.solve(*matrix, dst, src, precondition);
-}
-
-
-/* ------------------------------------------------------------------- */
-
-
-template<class MG, class VECTOR>
-PreconditionMG<MG, VECTOR>::PreconditionMG(MG &mg,
- const MGSmootherBase &pre,
- const MGSmootherBase &post,
- const MGCoarseGridSolver &coarse)
- :
- multigrid(&mg),
- pre(&pre),
- post(&post),
- coarse(&coarse)
-{}
-
-
-template<class MG, class VECTOR>
-void
-PreconditionMG<MG,VECTOR>::operator() (VECTOR &dst,
- const VECTOR &src) const
-{
- multigrid->copy_to_mg(src);
- multigrid->vcycle(*pre, *post, *coarse);
- multigrid->copy_from_mg(dst);
-}
-
-
-/*---------------------------- mgbase.h ---------------------------*/
-
-#endif
-/*---------------------------- mgbase.h ---------------------------*/
+++ /dev/null
-//---------------------------- mg_dof_accessor.h ---------------------------
-// $Id$
-// Version: $Name$
-//
-// Copyright (C) 1998, 1999, 2000 by the deal.II authors
-//
-// This file is subject to QPL and may not be distributed
-// without copyright and license information. Please refer
-// to the file deal.II/doc/license.html for the text and
-// further information on this license.
-//
-//---------------------------- mg_dof_accessor.h ---------------------------
-#ifndef __deal2__mg_dof_accessor_h
-#define __deal2__mg_dof_accessor_h
-
-
-#include <dofs/dof_accessor.h>
-
-
-/**
- * Define the basis for accessors to the degrees of freedom for
- * a multigrid DoF handler object.
- *
- * Note that it is allowed to construct an object of which the
- * #mg_dof_handler# pointer is a Null pointer. Such an object would
- * result in a strange kind of behaviour, though every reasonable
- * operating system should disallow access through that pointer.
- * The reason we do not check for the null pointer in the
- * constructor which gets passed the #DoFHandler# pointer is that
- * if we did we could not make dof iterators member of other classes
- * (like in the #FEValues# class) if we did not know about the
- * #DoFHandler# object to be used upon construction of that object.
- * Through the way this class is implemented here, we allow the
- * creation of a kind of virgin object which only gets useful if
- * assigned to from another object before first usage.
- *
- * Opposite to construction, it is not possible to copy an object
- * which has an invalid dof handler pointer. This is to guarantee
- * that every iterator which is once assigned to is a valid
- * object. However, this assertion only holds in debug mode, when
- * the #Assert# macro is switched on.
- *
- * @author Wolfgang Bangerth, 1998
- */
-template <int dim>
-class MGDoFAccessor {
- public:
- /**
- * Constructor
- */
- MGDoFAccessor () : mg_dof_handler(0) {
- Assert (false, ExcInvalidObject());
- };
-
- /**
- * This should be the default constructor.
- * We cast away the #const#ness of the
- * pointer which clearly is EVIL but
- * we can't help without making all
- * functions which could somehow use
- * iterators (directly or indirectly) make
- * non-const, even if they preserve
- * constness.
- */
- MGDoFAccessor (const MGDoFHandler<dim> *mg_dof_handler) :
- mg_dof_handler(const_cast<MGDoFHandler<dim>*>(mg_dof_handler)) {};
-
- /**
- * Reset the DoF handler pointer.
- */
- void set_mg_dof_handler (MGDoFHandler<dim> *dh) {
- Assert (dh != 0, DoFAccessor<dim>::ExcInvalidObject());
- mg_dof_handler = dh;
- };
-
- /**
- * Exception for child classes
- */
- DeclException0 (ExcInvalidObject);
-
- /**
- * Copy operator.
- */
- MGDoFAccessor<dim> & operator = (const MGDoFAccessor<dim> &da) {
- set_dof_handler (da.mg_dof_handler);
- return *this;
- };
-
- protected:
- /**
- * Store the address of the #MGDoFHandler# object
- * to be accessed.
- */
- MGDoFHandler<dim> *mg_dof_handler;
-};
-
-
-/* -------------------------------------------------------------------------- */
-
-/**
- * This is a switch class which only declares a #typdef#. It is meant to
- * determine which class a #MGDoFObjectAccessor# class is to be derived
- * from. By default, #MGDoFObjectAccessor<celldim,dim># derives from
- * the #typedef# in the general #MGDoFObjectAccessor_Inheritance<celldim,dim>#
- * class, which is #DoFObjectAccessor<celldim,dim>#,
- * but if #celldim==dim#, then the specialization #MGDoFObjectAccessor_Inheritance<dim,dim>#
- * is used which declares its local type to be #DoFCellAccessor<dim>#. Therefore,
- * the inheritance is automatically chosen to be from #DoFCellAccessor# if the
- * object under consideration has full dimension, i.e. constitutes a cell.
- *
- * @author Wolfgang Bangerth, 1999
- */
-template <int celldim, int dim>
-class MGDoFObjectAccessor_Inheritance
-{
- /**
- * Declaration of the #typedef#.
- * See the full documentation for
- * more information.
- */
- typedef DoFObjectAccessor<celldim,dim> BaseClass;
-};
-
-
-/**
- * This is a switch class which only declares a #typdef#. It is meant to
- * determine which class a #DoFObjectAccessor# class is to be derived
- * from. By default, #DoFObjectAccessor<celldim,dim># derives from
- * the #typedef# in the general #DoFObjectAccessor_Inheritance<celldim,dim>#
- * class, which is #TriaObjectAccessor<celldim,dim>#,
- * but if #celldim==dim#, then the specialization #TriaObjectAccessor<dim,dim>#
- * is used which declares its local type to be #CellAccessor<dim>#. Therefore,
- * the inheritance is automatically chosen to be from #CellAccessor# if the
- * object under consideration has full dimension, i.e. constitutes a cell.
- *
- * @author Wolfgang Bangerth, 1999
- */
-template <int dim>
-class MGDoFObjectAccessor_Inheritance<dim,dim>
-{
- /**
- * Declaration of the #typedef#.
- * See the full documentation for
- * more information.
- */
- typedef DoFCellAccessor<dim> BaseClass;
-};
-
-
-/* -------------------------------------------------------------------------- */
-
-
-/**
- * Common template for line, quad, hex.
- *
- * Internal: inheritance is necessary for the general template due to
- * a compiler error.
- * @author Guido Kanschat, 1999
- */
-template <int celldim, int dim>
-class MGDoFObjectAccessor : public MGDoFAccessor<dim>,
- public MGDoFObjectAccessor_Inheritance<celldim, dim>::BaseClass
-{};
-
-
-/**
- * Closure class.
- */
-template<int dim>
-class MGDoFObjectAccessor<0, dim>
-{
- public:
- typedef void AccessorData;
-
- /**
- * Constructor. Should never be called
- * and thus throws an error.
- */
- MGDoFObjectAccessor (Triangulation<dim> *,
- const int,
- const int,
- const AccessorData *)
- {
- Assert (false, ExcInternalError());
- }
-};
-
-
-/**
- * Grant access to the degrees of freedom located on lines.
- * This class follows mainly the route laid out by the accessor library
- * declared in the triangulation library (\Ref{TriaAccessor}). It enables
- * the user to access the degrees of freedom on the lines (there are similar
- * versions for the DoFs on quads, etc), where the dimension of the underlying
- * triangulation does not really matter (i.e. this accessor works with the
- * lines in 1D-, 2D-, etc dimensions).
- *
- *
- * \subsection{Usage}
- *
- * The \Ref{DoFDimensionInfo} classes inherited by the \Ref{DoFHandler} classes
- * declare typedefs to iterators using the accessors declared in this class
- * hierarchy tree. Usage is best to happens through these typedefs, since they
- * are more secure to changes in the class naming and template interface as well
- * as they provide easier typing (much less complicated names!).
- *
- *
- * \subsection{Notes about the class hierarchy structure}
- *
- * Inheritance from #MGDoFObjectAccessor_Inheritance<1,dim>::BaseClass# yields
- * inheritance from #DoFCellAccessor<1># if #dim==1# and from
- * #TriaObjectAccessor<1,dim># for all other #dim# values. Thus, an object
- * of this class shares all features of cells in one dimension, but behaves
- * like an ordinary line in all other cases.
- *
- * @author Wolfgang Bangerth, 1998
- */
-template <int dim>
-class MGDoFObjectAccessor<1, dim> : public MGDoFAccessor<dim>,
- public MGDoFObjectAccessor_Inheritance<1, dim>::BaseClass
-{
- public:
- /**
- * Declare the data type that this accessor
- * class expects to get passed from the
- * iterator classes.
- */
- typedef MGDoFHandler<dim> AccessorData;
-
- /**
- * Default constructor, unused thus
- * not implemented.
- */
- MGDoFObjectAccessor ();
-
- /**
- * Constructor. The #local_data#
- * argument is assumed to be a pointer
- * to a #MGDoFHandler<dim># object.
- */
- MGDoFObjectAccessor (Triangulation<dim> *tria,
- const int level,
- const int index,
- const AccessorData *local_data);
-
- /**
- * Return the index of the #i#th degree
- * of freedom of this line on the level
- * this line lives on.
- */
- unsigned int mg_dof_index (const unsigned int i) const;
-
- /**
- * Set the index of the #i#th degree
- * of freedom of this line on the
- * level this line lives on to #index#.
- */
- void set_mg_dof_index (const unsigned int i,
- const unsigned int index) const;
-
- /**
- * Return the index of the #i#th degree
- * on the #vertex#th vertex for the
- * level this line lives on.
- */
- unsigned int mg_vertex_dof_index (const unsigned int vertex,
- const unsigned int i) const;
-
- /**
- * Set the index of the #i#th degree
- * on the #vertex#th vertex to #index#
- * for the level this line lives on.
- */
- void set_mg_vertex_dof_index (const unsigned int vertex,
- const unsigned int i,
- const unsigned int index) const;
-
- /**
- * Return the indices of the dofs of this
- * line in the standard ordering: dofs
- * on vertex 0, dofs on vertex 1,
- * dofs on line.
- *
- * It is assumed that the vector already
- * has the right size beforehand. The
- * indices refer to the local numbering
- * for the level this line lives on.
- */
- void get_mg_dof_indices (vector<unsigned int> &dof_indices) const;
-
- /**
- * Return the value of the given vector
- * restricted to the dofs of this
- * cell in the standard ordering: dofs
- * on vertex 0, dofs on vertex 1, etc,
- * dofs on line 0, dofs on line 1, etc,
- * dofs on quad 0, etc.
- *
- * It is assumed that the vector already
- * has the right size beforehand. The
- * indices refer to the multilevel
- * numbering local to the present
- * level of this cell. The vector shall
- * therefore have the same number of
- * entries as there are degrees of
- * freedom on this level.
- */
- template <typename number>
- void get_mg_dof_values (const Vector<number> &values,
- Vector<number> &dof_values) const;
-
- /**
- * Return the #i#th child as a MGDoF line
- * iterator. This function is needed since
- * the child function of the base
- * class returns a line accessor without
- * access to the MG data.
- */
- TriaIterator<dim,MGDoFObjectAccessor<1, dim> > child (const unsigned int) const;
-
- /**
- * Implement the copy operator needed
- * for the iterator classes.
- */
- void copy_from (const MGDoFObjectAccessor<1, dim> &a);
-};
-
-
-/**
- * Grant access to the multilevel degrees of freedom located on quads.
- *
- * @see DoFLineAccessor
- */
-template <int dim>
-class MGDoFObjectAccessor<2, dim> : public MGDoFAccessor<dim>,
- public MGDoFObjectAccessor_Inheritance<2, dim>::BaseClass
-{
- public:
- /**
- * Declare the data type that this accessor
- * class expects to get passed from the
- * iterator classes.
- */
- typedef MGDoFHandler<dim> AccessorData;
-
- /**
- * Default constructor, unused thus
- * not implemented.
- */
- MGDoFObjectAccessor ();
-
- /**
- * Constructor. The #local_data#
- * argument is assumed to be a pointer
- * to a #DoFHandler<dim># object.
- */
- MGDoFObjectAccessor (Triangulation<dim> *tria,
- const int level,
- const int index,
- const AccessorData *local_data);
-
- /**
- * Return the index of the #i#th degree
- * of freedom of this quad on the level
- * this quad lives on.
- */
- unsigned int mg_dof_index (const unsigned int i) const;
-
- /**
- * Set the index of the #i#th degree
- * of freedom of this quad on the
- * level this quad lives on to #index#.
- */
- void set_mg_dof_index (const unsigned int i,
- const unsigned int index) const;
-
- /**
- * Return the index of the #i#th degree
- * on the #vertex#th vertex for the level
- * this quad lives on.
- */
- unsigned int mg_vertex_dof_index (const unsigned int vertex,
- const unsigned int i) const;
-
- /**
- * Set the index of the #i#th degree
- * on the #vertex#th vertex for the
- * level this quad lives on to #index#.
- */
- void set_mg_vertex_dof_index (const unsigned int vertex,
- const unsigned int i,
- const unsigned int index) const;
-
- /**
- * Return the indices of the dofs of this
- * quad in the standard ordering: dofs
- * on vertex 0, dofs on vertex 1, etc,
- * dofs on line 0, dofs on line 1, etc,
- * dofs on quad 0, etc.
- *
- * It is assumed that the vector already
- * has the right size beforehand. The
- * indices refer to the local numbering
- * for the level this quad lives on.
- */
- void get_mg_dof_indices (vector<unsigned int> &dof_indices) const;
-
- /**
- * Return the value of the given vector
- * restricted to the dofs of this
- * cell in the standard ordering: dofs
- * on vertex 0, dofs on vertex 1, etc,
- * dofs on line 0, dofs on line 1, etc,
- * dofs on quad 0, etc.
- *
- * It is assumed that the vector already
- * has the right size beforehand. The
- * indices refer to the multilevel
- * numbering local to the present
- * level of this cell. The vector shall
- * therefore have the same number of
- * entries as there are degrees of
- * freedom on this level.
- */
- template <typename number>
- void get_mg_dof_values (const Vector<number> &values,
- Vector<number> &dof_values) const;
-
- /**
- * Return a pointer to the #i#th line
- * bounding this #Quad#.
- */
- TriaIterator<dim,MGDoFObjectAccessor<1, dim> >
- line (const unsigned int i) const;
-
- /**
- * Return the #i#th child as a DoF quad
- * iterator. This function is needed since
- * the child function of the base
- * class returns a quad accessor without
- * access to the DoF data.
- */
- TriaIterator<dim,MGDoFObjectAccessor<2, dim> >
- child (const unsigned int) const;
-
- /**
- * Implement the copy operator needed
- * for the iterator classes.
- */
- void copy_from (const MGDoFObjectAccessor<2, dim> &a);
-};
-
-
-/**
- * Grant access to the multilevel degrees of freedom located on hexhedra.
- *
- * @see DoFLineAccessor
- */
-template <int dim>
-class MGDoFObjectAccessor<3, dim> : public MGDoFAccessor<dim>,
- public MGDoFObjectAccessor_Inheritance<3, dim>::BaseClass
-{
- public:
- /**
- * Declare the data type that this accessor
- * class expects to get passed from the
- * iterator classes.
- */
- typedef MGDoFHandler<dim> AccessorData;
-
- /**
- * Default constructor, unused thus
- * not implemented.
- */
- MGDoFObjectAccessor ();
-
- /**
- * Constructor. The #local_data#
- * argument is assumed to be a pointer
- * to a #DoFHandler<dim># object.
- */
- MGDoFObjectAccessor (Triangulation<dim> *tria,
- const int level,
- const int index,
- const AccessorData *local_data);
-
- /**
- * Return the index of the #i#th degree
- * of freedom of this hex on the level
- * this quad lives on.
- */
- unsigned int mg_dof_index (const unsigned int i) const;
-
- /**
- * Set the index of the #i#th degree
- * of freedom of this hex on the
- * level this hex lives on to #index#.
- */
- void set_mg_dof_index (const unsigned int i,
- const unsigned int index) const;
-
- /**
- * Return the index of the #i#th degree
- * on the #vertex#th vertex for the level
- * this hex lives on.
- */
- unsigned int mg_vertex_dof_index (const unsigned int vertex,
- const unsigned int i) const;
-
- /**
- * Set the index of the #i#th degree
- * on the #vertex#th vertex for the
- * level this hex lives on to #index#.
- */
- void set_mg_vertex_dof_index (const unsigned int vertex,
- const unsigned int i,
- const unsigned int index) const;
-
- /**
- * Return the indices of the dofs of this
- * hex in the standard ordering: dofs
- * on vertex 0, dofs on vertex 1, etc,
- * dofs on line 0, dofs on line 1, etc,
- * dofs on quad 0, etc.
- *
- * It is assumed that the vector already
- * has the right size beforehand. The
- * indices refer to the local numbering
- * for the level this hex lives on.
- */
- void get_mg_dof_indices (vector<unsigned int> &dof_indices) const;
-
- /**
- * Return the value of the given vector
- * restricted to the dofs of this
- * cell in the standard ordering: dofs
- * on vertex 0, dofs on vertex 1, etc,
- * dofs on line 0, dofs on line 1, etc,
- * dofs on quad 0, etc.
- *
- * It is assumed that the vector already
- * has the right size beforehand. The
- * indices refer to the multilevel
- * numbering local to the present
- * level of this cell. The vector shall
- * therefore have the same number of
- * entries as there are degrees of
- * freedom on this level.
- */
- template <typename number>
- void get_mg_dof_values (const Vector<number> &values,
- Vector<number> &dof_values) const;
-
- /**
- * Return a pointer to the #i#th line
- * bounding this #Hex#.
- */
- TriaIterator<dim,MGDoFObjectAccessor<1, dim> >
- line (const unsigned int i) const;
-
- /**
- * Return a pointer to the #i#th quad
- * bounding this #Hex#.
- */
- TriaIterator<dim,MGDoFObjectAccessor<2, dim> >
- quad (const unsigned int i) const;
-
- /**
- * Return the #i#th child as a DoF quad
- * iterator. This function is needed since
- * the child function of the base
- * class returns a hex accessor without
- * access to the DoF data.
- */
- TriaIterator<dim,MGDoFObjectAccessor<3, dim> > child (const unsigned int) const;
-
- /**
- * Implement the copy operator needed
- * for the iterator classes.
- */
- void copy_from (const MGDoFObjectAccessor<3, dim> &a);
-};
-
-
-/**
- * Grant access to the degrees of freedom on a cell. In fact, since all
- * access to the degrees of freedom has been enabled by the classes
- * #DoFObjectAccessor<1, 1># and #DoFObjectAccessor<2, 2># for the space dimension
- * one and two, respectively, this class only collects the pieces
- * together by deriving from the appropriate #DoF*Accessor# and the
- * right #CellAccessor<dim># and finally adding two functions which give
- * access to the neighbors and children as #DoFCellAccessor# objects
- * rather than #CellAccessor# objects (the latter function was inherited
- * from the #CellAccessor<dim># class).
- *
- * @author Wolfgang Bangerth, 1998
- */
-template <int dim>
-class MGDoFCellAccessor : public MGDoFObjectAccessor<dim, dim> {
- public:
- /**
- * Type of faces.
- */
- typedef
- TriaIterator<dim, MGDoFObjectAccessor<dim-1, dim> >
- face_iterator;
- /**
- * Declare the data type that this accessor
- * class expects to get passed from the
- * iterator classes.
- */
- typedef typename MGDoFObjectAccessor<dim, dim>::AccessorData AccessorData;
-
- /**
- * Constructor
- */
- MGDoFCellAccessor (Triangulation<dim> *tria,
- const int level,
- const int index,
- const AccessorData *local_data) :
- MGDoFObjectAccessor<dim, dim> (tria,level,index,local_data) {};
-
- /**
- * Return the #i#th neighbor as a MGDoF cell
- * iterator. This function is needed since
- * the neighbor function of the base
- * class returns a cell accessor without
- * access to the MGDoF data.
- */
- TriaIterator<dim,MGDoFCellAccessor<dim> > neighbor (const unsigned int) const;
-
- /**
- * Return the #i#th child as a MGDoF cell
- * iterator. This function is needed since
- * the child function of the base
- * class returns a cell accessor without
- * access to the DoF data.
- */
- TriaIterator<dim,MGDoFCellAccessor<dim> > child (const unsigned int) const;
-
- /**
- * Return an iterator to the #i#th face
- * of this cell.
- *
- * This function is not implemented in 1D,
- * and maps to MGDoFObjectAccessor<2, dim>::line in 2D.
- */
- face_iterator
- face (const unsigned int i) const;
-
- /**
- * Exception
- */
- DeclException0 (ExcNotUsefulForThisDimension);
-};
-
-
-#endif
+++ /dev/null
-//---------------------------- mg_dof_handler.h ---------------------------
-// $Id$
-// Version: $Name$
-//
-// Copyright (C) 1998, 1999, 2000 by the deal.II authors
-//
-// This file is subject to QPL and may not be distributed
-// without copyright and license information. Please refer
-// to the file deal.II/doc/license.html for the text and
-// further information on this license.
-//
-//---------------------------- mg_dof_handler.h ---------------------------
-#ifndef __deal2__mg_dof_handler_h
-#define __deal2__mg_dof_handler_h
-
-
-/*---------------------------- mg_dof.h ---------------------------*/
-
-
-#include <dofs/dof_handler.h>
-
-
-/**
- * Define some types which differ between the dimensions. This class
- * is analogous to the \Ref{TriaDimensionInfo} class hierarchy.
- *
- * @see MGDoFDimensionInfo<1>
- * @see MGDoFDimensionInfo<2>
- */
-template <int dim>
-class MGDoFDimensionInfo;
-
-
-/**
- * Define some types for the DoF handling in one dimension.
- *
- * The types have the same meaning as those declared in \Ref{TriaDimensionInfo<2>}.
- */
-class MGDoFDimensionInfo<1> {
- public:
- typedef TriaRawIterator<1,MGDoFCellAccessor<1> > raw_line_iterator;
- typedef TriaIterator<1,MGDoFCellAccessor<1> > line_iterator;
- typedef TriaActiveIterator<1,MGDoFCellAccessor<1> > active_line_iterator;
-
- typedef void * raw_quad_iterator;
- typedef void * quad_iterator;
- typedef void * active_quad_iterator;
-
- typedef void * raw_hex_iterator;
- typedef void * hex_iterator;
- typedef void * active_hex_iterator;
-
- typedef raw_line_iterator raw_cell_iterator;
- typedef line_iterator cell_iterator;
- typedef active_line_iterator active_cell_iterator;
-
- typedef void * raw_face_iterator;
- typedef void * face_iterator;
- typedef void * active_face_iterator;
-};
-
-
-/**
- * Define some types for the DoF handling in two dimensions.
- *
- * The types have the same meaning as those declared in \Ref{TriaDimensionInfo<2>}.
- */
-class MGDoFDimensionInfo<2> {
- public:
- typedef TriaRawIterator<2,MGDoFObjectAccessor<1, 2> > raw_line_iterator;
- typedef TriaIterator<2,MGDoFObjectAccessor<1, 2> > line_iterator;
- typedef TriaActiveIterator<2,MGDoFObjectAccessor<1, 2> > active_line_iterator;
-
- typedef TriaRawIterator<2,MGDoFCellAccessor<2> > raw_quad_iterator;
- typedef TriaIterator<2,MGDoFCellAccessor<2> > quad_iterator;
- typedef TriaActiveIterator<2,MGDoFCellAccessor<2> > active_quad_iterator;
-
- typedef void * raw_hex_iterator;
- typedef void * hex_iterator;
- typedef void * active_hex_iterator;
-
- typedef raw_quad_iterator raw_cell_iterator;
- typedef quad_iterator cell_iterator;
- typedef active_quad_iterator active_cell_iterator;
-
- typedef raw_line_iterator raw_face_iterator;
- typedef line_iterator face_iterator;
- typedef active_line_iterator active_face_iterator;
-};
-
-
-/**
- * Define some types for the DoF handling in two dimensions.
- *
- * The types have the same meaning as those declared in \Ref{TriaDimensionInfo<2>}.
- */
-class MGDoFDimensionInfo<3> {
- public:
- typedef TriaRawIterator<3,MGDoFObjectAccessor<1, 3> > raw_line_iterator;
- typedef TriaIterator<3,MGDoFObjectAccessor<1, 3> > line_iterator;
- typedef TriaActiveIterator<3,MGDoFObjectAccessor<1, 3> > active_line_iterator;
-
- typedef TriaRawIterator<3,MGDoFObjectAccessor<2, 3> > raw_quad_iterator;
- typedef TriaIterator<3,MGDoFObjectAccessor<2, 3> > quad_iterator;
- typedef TriaActiveIterator<3,MGDoFObjectAccessor<2, 3> > active_quad_iterator;
-
- typedef TriaRawIterator<3,MGDoFCellAccessor<3> > raw_hex_iterator;
- typedef TriaIterator<3,MGDoFCellAccessor<3> > hex_iterator;
- typedef TriaActiveIterator<3,MGDoFCellAccessor<3> > active_hex_iterator;
-
- typedef raw_hex_iterator raw_cell_iterator;
- typedef hex_iterator cell_iterator;
- typedef active_hex_iterator active_cell_iterator;
-
- typedef raw_quad_iterator raw_face_iterator;
- typedef quad_iterator face_iterator;
- typedef active_quad_iterator active_face_iterator;
-};
-
-
-/**
- * This class manages degrees of freedom for a multilevel hierarchy of
- * grids. It does mostly the same as does the #DoDHandler# class, but
- * it uses a separate enumeration of the degrees of freedom on each
- * level. For example, a vertex has several DoF numbers, one for each
- * level of the triangulation on which it exists.
- *
- * At present, multilevel algorithms are not fully functional, so this
- * documentation is still very brief.
-//TODO: Extend MGDoFHandler doc
- *
- * @author Wolfgang Bangerth, 1998, 1999
- */
-template <int dim>
-class MGDoFHandler : public DoFHandler<dim>
-{
- public:
- typedef typename MGDoFDimensionInfo<dim>::raw_line_iterator raw_line_iterator;
- typedef typename MGDoFDimensionInfo<dim>::line_iterator line_iterator;
- typedef typename MGDoFDimensionInfo<dim>::active_line_iterator active_line_iterator;
-
- typedef typename MGDoFDimensionInfo<dim>::raw_quad_iterator raw_quad_iterator;
- typedef typename MGDoFDimensionInfo<dim>::quad_iterator quad_iterator;
- typedef typename MGDoFDimensionInfo<dim>::active_quad_iterator active_quad_iterator;
-
- typedef typename MGDoFDimensionInfo<dim>::raw_hex_iterator raw_hex_iterator;
- typedef typename MGDoFDimensionInfo<dim>::hex_iterator hex_iterator;
- typedef typename MGDoFDimensionInfo<dim>::active_hex_iterator active_hex_iterator;
-
- typedef typename MGDoFDimensionInfo<dim>::raw_cell_iterator raw_cell_iterator;
- typedef typename MGDoFDimensionInfo<dim>::cell_iterator cell_iterator;
- typedef typename MGDoFDimensionInfo<dim>::active_cell_iterator active_cell_iterator;
-
- typedef typename MGDoFDimensionInfo<dim>::raw_face_iterator raw_face_iterator;
- typedef typename MGDoFDimensionInfo<dim>::face_iterator face_iterator;
- typedef typename MGDoFDimensionInfo<dim>::active_face_iterator active_face_iterator;
-
- /**
- * Constructor. Take #tria# as the
- * triangulation to work on.
- */
- MGDoFHandler (Triangulation<dim> &tria);
-
- /**
- * Destructor
- */
- virtual ~MGDoFHandler ();
-
- /**
- * Go through the triangulation and
- * distribute the degrees of freedoms
- * needed for the given finite element
- * according to the given distribution
- * method. We first call the #DoFHandler#'s
- * function and then distribute the
- * levelwise numbers.
- *
- * A copy of the transferred finite
- * element is stored.
- *
- * This function uses the user flags of the
- * triangulation object, so make sure you
- * don't need them after calling this
- * function, or if so store them.
- */
- virtual void distribute_dofs (const FiniteElement<dim> &, unsigned int offset=0);
-
- /**
- * Actually do the renumbering based on
- * a list of new dof numbers for all the
- * dofs.
- *
- * #new_numbers# is an array of integers
- * with size equal to the number of dofs
- * on the present level. It stores the new
- * indices after renumbering in the
- * order of the old indices.
- */
- void renumber_dofs (const unsigned int level,
- const vector<unsigned int> &new_numbers);
-
- /*--------------------------------------*/
-
- /**
- * @name Cell iterator functions
- */
- /*@{*/
- /**
- * Iterator to the first cell, used
- * or not, on level #level#. If a level
- * has no cells, a past-the-end iterator
- * is returned.
- *
- * This function calls #begin_raw_line#
- * in 1D and #begin_raw_quad# in 2D.
- */
- raw_cell_iterator begin_raw (const unsigned int level = 0) const;
-
- /**
- * Iterator to the first used cell
- * on level #level#.
- *
- * This function calls #begin_line#
- * in 1D and #begin_quad# in 2D.
- */
- cell_iterator begin (const unsigned int level = 0) const;
-
- /**
- * Iterator to the first active
- * cell on level #level#.
- *
- * This function calls #begin_active_line#
- * in 1D and #begin_active_quad# in 2D.
- */
- active_cell_iterator begin_active(const unsigned int level = 0) const;
-
- /**
- * Iterator past the end; this
- * iterator serves for comparisons of
- * iterators with past-the-end or
- * before-the-beginning states.
- *
- * This function calls #end_line#
- * in 1D and #end_quad# in 2D.
- */
- raw_cell_iterator end () const;
-
- /**
- * Return an iterator which is the first
- * iterator not on level. If #level# is
- * the last level, then this returns
- * #end()#.
- */
- cell_iterator end (const unsigned int level) const;
-
- /**
- * Return a raw iterator which is the first
- * iterator not on level. If #level# is
- * the last level, then this returns
- * #end()#.
- */
- raw_cell_iterator end_raw (const unsigned int level) const;
-
- /**
- * Return an active iterator which is the
- * first iterator not on level. If #level#
- * is the last level, then this returns
- * #end()#.
- */
- active_cell_iterator end_active (const unsigned int level) const;
-
- /**
- * Return an iterator pointing to the
- * last cell, used or not.
- *
- * This function calls #last_raw_line#
- * in 1D and #last_raw_quad# in 2D.
- */
- raw_cell_iterator last_raw () const;
-
- /**
- * Return an iterator pointing to the last
- * cell of the level #level#, used or not.
- *
- * This function calls #last_raw_line#
- * in 1D and #last_raw_quad# in 2D.
- */
- raw_cell_iterator last_raw (const unsigned int level) const;
-
- /**
- * Return an iterator pointing to the last
- * used cell.
- *
- * This function calls #last_line#
- * in 1D and #last_quad# in 2D.
- */
- cell_iterator last () const;
-
- /**
- * Return an iterator pointing to the last
- * used cell on level #level#.
- *
- * This function calls #last_line#
- * in 1D and #last_quad# in 2D.
- */
- cell_iterator last (const unsigned int level) const;
-
- /**
- * Return an iterator pointing to the last
- * active cell.
- *
- * This function calls #last_active_line#
- * in 1D and #last_active_quad# in 2D.
- */
- active_cell_iterator last_active () const;
-
- /**
- * Return an iterator pointing to the last
- * active cell on level #level#.
- *
- * This function calls #last_active_line#
- * in 1D and #last_active_quad# in 2D.
- */
- active_cell_iterator last_active (const unsigned int level) const;
- //@}
-
- /*---------------------------------------*/
-
- /**
- * @name Face iterator functions
- */
- /*@{*/
- /**
- * Iterator to the first face, used
- * or not, on level #level#. If a level
- * has no faces, a past-the-end iterator
- * is returned.
- *
- * This function calls #begin_raw_line#
- * in 2D and #begin_raw_quad# in 3D.
- */
- raw_face_iterator begin_raw_face (const unsigned int level = 0) const;
-
- /**
- * Iterator to the first used face
- * on level #level#.
- *
- * This function calls #begin_line#
- * in 2D and #begin_quad# in 3D.
- */
- face_iterator begin_face (const unsigned int level = 0) const;
-
- /**
- * Iterator to the first active
- * face on level #level#.
- *
- * This function calls #begin_active_line#
- * in 2D and #begin_active_quad# in 3D.
- */
- active_face_iterator begin_active_face(const unsigned int level = 0) const;
-
- /**
- * Iterator past the end; this
- * iterator serves for comparisons of
- * iterators with past-the-end or
- * before-the-beginning states.
- *
- * This function calls #end_line#
- * in 2D and #end_quad# in 3D.
- */
- raw_face_iterator end_face () const;
-
- /**
- * Return an iterator which is the first
- * iterator not on level. If #level# is
- * the last level, then this returns
- * #end()#.
- */
- face_iterator end_face (const unsigned int level) const;
-
- /**
- * Return a raw iterator which is the first
- * iterator not on level. If #level# is
- * the last level, then this returns
- * #end()#.
- */
- raw_face_iterator end_raw_face (const unsigned int level) const;
-
- /**
- * Return an active iterator which is the
- * first iterator not on level. If #level#
- * is the last level, then this returns
- * #end()#.
- */
- active_face_iterator end_active_face (const unsigned int level) const;
-
- /**
- * Return an iterator pointing to the
- * last face, used or not.
- *
- * This function calls #last_raw_line#
- * in 2D and #last_raw_quad# in 3D.
- */
- raw_face_iterator last_raw_face () const;
-
- /**
- * Return an iterator pointing to the last
- * face of the level #level#, used or not.
- *
- * This function calls #last_raw_line#
- * in 2D and #last_raw_quad# in 3D.
- */
- raw_face_iterator last_raw_face (const unsigned int level) const;
-
- /**
- * Return an iterator pointing to the last
- * used face.
- *
- * This function calls #last_line#
- * in 2D and #last_quad# in 3D.
- */
- face_iterator last_face () const;
-
- /**
- * Return an iterator pointing to the last
- * used face on level #level#.
- *
- * This function calls #last_line#
- * in 2D and #last_quad# in 3D.
- */
- face_iterator last_face (const unsigned int level) const;
-
- /**
- * Return an iterator pointing to the last
- * active face.
- *
- * This function calls #last_active_line#
- * in 2D and #last_active_quad# in 3D.
- */
- active_face_iterator last_active_face () const;
-
- /**
- * Return an iterator pointing to the last
- * active face on level #level#.
- *
- * This function calls #last_active_line#
- * in 2D and #last_active_quad# in 3D.
- */
- active_face_iterator last_active_face (const unsigned int level) const;
- //@}
-
-
-/*---------------------------------------*/
-
- /**
- * @name Line iterator functions
- */
- /*@{*/
- /**
- * Iterator to the first line, used
- * or not, on level #level#. If a level
- * has no lines, a past-the-end iterator
- * is returned.
- */
- raw_line_iterator begin_raw_line (const unsigned int level = 0) const;
-
- /**
- * Iterator to the first used line
- * on level #level#.
- */
- line_iterator begin_line (const unsigned int level = 0) const;
-
- /**
- * Iterator to the first active
- * line on level #level#.
- */
- active_line_iterator begin_active_line(const unsigned int level = 0) const;
-
- /**
- * Iterator past the end; this
- * iterator serves for comparisons of
- * iterators with past-the-end or
- * before-the-beginning states.
- */
- raw_line_iterator end_line () const;
-
- /**
- * Return an iterator which is the first
- * iterator not on level. If #level# is
- * the last level, then this returns
- * #end()#.
- */
- line_iterator end_line (const unsigned int level) const;
-
- /**
- * Return a raw iterator which is the first
- * iterator not on level. If #level# is
- * the last level, then this returns
- * #end()#.
- */
- raw_line_iterator end_raw_line (const unsigned int level) const;
-
- /**
- * Return an active iterator which is the
- * first iterator not on level. If #level#
- * is the last level, then this returns
- * #end()#.
- */
- active_line_iterator end_active_line (const unsigned int level) const;
-
-
-/**
- * Return an iterator pointing to the
- * last line, used or not.
- */
- raw_line_iterator last_raw_line () const;
-
- /**
- * Return an iterator pointing to the last
- * line of the level #level#, used or not.
-
- */
- raw_line_iterator last_raw_line (const unsigned int level) const;
-
- /**
- * Return an iterator pointing to the last
- * used line.
- */
- line_iterator last_line () const;
-
- /**
- * Return an iterator pointing to the last
- * used line on level #level#.
- */
- line_iterator last_line (const unsigned int level) const;
-
- /**
- * Return an iterator pointing to the last
- * active line.
- */
- active_line_iterator last_active_line () const;
-
- /**
- * Return an iterator pointing to the last
- * active line on level #level#.
- */
- active_line_iterator last_active_line (const unsigned int level) const;
- /*@}*/
-
- /*---------------------------------------*/
-
- /**
- * @name Quad iterator functions*/
- /*@{
- */
- /**
- * Iterator to the first quad, used
- * or not, on level #level#. If a level
- * has no quads, a past-the-end iterator
- * is returned.
- */
- raw_quad_iterator begin_raw_quad (const unsigned int level = 0) const;
-
- /**
- * Iterator to the first used quad
- * on level #level#.
- */
- quad_iterator begin_quad (const unsigned int level = 0) const;
-
- /**
- * Iterator to the first active
- * quad on level #level#.
- */
- active_quad_iterator begin_active_quad(const unsigned int level = 0) const;
-
- /**
- * Iterator past the end; this
- * iterator serves for comparisons of
- * iterators with past-the-end or
- * before-the-beginning states.
- */
- raw_quad_iterator end_quad () const;
-
- /**
- * Return an iterator which is the first
- * iterator not on level. If #level# is
- * the last level, then this returns
- * #end()#.
- */
- quad_iterator end_quad (const unsigned int level) const;
-
- /**
- * Return a raw iterator which is the first
- * iterator not on level. If #level# is
- * the last level, then this returns
- * #end()#.
- */
- raw_quad_iterator end_raw_quad (const unsigned int level) const;
-
- /**
- * Return an active iterator which is the
- * first iterator not on level. If #level#
- * is the last level, then this returns
- * #end()#.
- */
- active_quad_iterator end_active_quad (const unsigned int level) const;
-
-
-/**
- * Return an iterator pointing to the
- * last quad, used or not.
- */
- raw_quad_iterator last_raw_quad () const;
-
- /**
- * Return an iterator pointing to the last
- * quad of the level #level#, used or not.
-
- */
- raw_quad_iterator last_raw_quad (const unsigned int level) const;
-
- /**
- * Return an iterator pointing to the last
- * used quad.
- */
- quad_iterator last_quad () const;
-
- /**
- * Return an iterator pointing to the last
- * used quad on level #level#.
- */
- quad_iterator last_quad (const unsigned int level) const;
-
- /**
- * Return an iterator pointing to the last
- * active quad.
- */
- active_quad_iterator last_active_quad () const;
-
- /**
- * Return an iterator pointing to the last
- * active quad on level #level#.
- */
- active_quad_iterator last_active_quad (const unsigned int level) const;
- /*@}*/
-
- /*---------------------------------------*/
-
- /**
- * @name Hex iterator functions*/
- /*@{
- */
- /**
- * Iterator to the first hex, used
- * or not, on level #level#. If a level
- * has no hexs, a past-the-end iterator
- * is returned.
- */
- raw_hex_iterator begin_raw_hex (const unsigned int level = 0) const;
-
- /**
- * Iterator to the first used hex
- * on level #level#.
- */
- hex_iterator begin_hex (const unsigned int level = 0) const;
-
- /**
- * Iterator to the first active
- * hex on level #level#.
- */
- active_hex_iterator begin_active_hex(const unsigned int level = 0) const;
-
- /**
- * Iterator past the end; this
- * iterator serves for comparisons of
- * iterators with past-the-end or
- * before-the-beginning states.
- */
- raw_hex_iterator end_hex () const;
-
- /**
- * Return an iterator which is the first
- * iterator not on level. If #level# is
- * the last level, then this returns
- * #end()#.
- */
- hex_iterator end_hex (const unsigned int level) const;
-
- /**
- * Return a raw iterator which is the first
- * iterator not on level. If #level# is
- * the last level, then this returns
- * #end()#.
- */
- raw_hex_iterator end_raw_hex (const unsigned int level) const;
-
- /**
- * Return an active iterator which is the
- * first iterator not on level. If #level#
- * is the last level, then this returns
- * #end()#.
- */
- active_hex_iterator end_active_hex (const unsigned int level) const;
-
-
-/**
- * Return an iterator pointing to the
- * last hex, used or not.
- */
- raw_hex_iterator last_raw_hex () const;
-
- /**
- * Return an iterator pointing to the last
- * hex of the level #level#, used or not.
-
- */
- raw_hex_iterator last_raw_hex (const unsigned int level) const;
-
- /**
- * Return an iterator pointing to the last
- * used hex.
- */
- hex_iterator last_hex () const;
-
- /**
- * Return an iterator pointing to the last
- * used hex on level #level#.
- */
- hex_iterator last_hex (const unsigned int level) const;
-
- /**
- * Return an iterator pointing to the last
- * active hex.
- */
- active_hex_iterator last_active_hex () const;
-
- /**
- * Return an iterator pointing to the last
- * active hex on level #level#.
- */
- active_hex_iterator last_active_hex (const unsigned int level) const;
- /*@}*/
-
- /*---------------------------------------*/
-
- /**
- * Return the number of degrees of freedom
- * on the specified level.
- * Included in this number are those
- * DoFs which are constrained by
- * hanging nodes.
- */
- unsigned int n_dofs (const unsigned int level) const;
-
- /**
- * Exception.
- */
- DeclException1 (ExcInvalidLevel,
- int,
- << "The level index " << arg1 << "was not in the valid range.");
-
- private:
-
- /** We need for each vertex a list of the
- * degree of freedom indices on each of the
- * levels this vertex lives on. Since most
- * vertices live only on a few levels, it
- * is not economical to reserve indices for
- * all the levels there are; rather, we
- * create an object which holds the indices
- * on those levels only where the vertex
- * lives. To construct such an array, it
- * is necessary to know beforehand which
- * is the coarsest level the vertex lives
- * on, how many levels it lives on and
- * how many dofs there are on each vertex.
- * If we have this information, we can
- * allocate exactly the amount of memory
- * which is needed and need not handle
- * growing arrays and the like.
- */
- class MGVertexDoFs {
- public:
-
- /**
- * Constructor. This one is empty
- * because it is difficult to make it
- * efficient to use vector<>'s and
- * still construct the object using
- * the constructor. Use the #init#
- * function to really allocate
- * memory.
- */
- MGVertexDoFs ();
-
- /**
- * Allocate memory and
- * set all indices to #-1#.
- */
- void init (const unsigned int coarsest_level,
- const unsigned int finest_level,
- const unsigned int dofs_per_vertex);
-
- /**
- * Destructor
- */
- ~MGVertexDoFs ();
-
- /**
- * Set the index with number
- * #dof_number# of this vertex on
- * #level# to the given index. To
- * compute the position in the array,
- * one has to specify how many
- * dofs per vertex there are. It is
- * not checked that the level number
- * is below the number of the finest
- * level this vertex lives on.
- *
- * The function is inline, so should
- * be reasonably fast.
- */
- void set_index (const unsigned int level,
- const unsigned int dof_number,
- const unsigned int dofs_per_vertex,
- const unsigned int index);
-
- /**
- * Return the index with number
- * #dof_number# of this vertex on
- * #level#. To
- * compute the position in the array,
- * one has to specify how many
- * dofs per vertex there are. It is
- * not checked that the level number
- * is below the number of the finest
- * level this vertex lives on.
- *
- * The function is inline, so should
- * be reasonably fast.
- */
- unsigned int get_index (const unsigned int level,
- const unsigned int dof_number,
- const unsigned int dofs_per_vertex) const;
-
- /**
- * Return the index of the coarsest
- * level this vertex lives on.
- */
- unsigned int get_coarsest_level () const;
-
- /**
- * Return the index of the finest
- * level this vertex lives on.
- */
- unsigned int get_finest_level () const;
-
- /**
- * Exception.
- */
- DeclException0 (ExcNoMemory);
- /**
- * Exception.
- */
- DeclException1 (ExcInvalidLevel,
- int,
- << "The given level number " << arg1 << " is outside "
- << "the range of levels this vertex lives on.");
- /**
- * Exception.
- */
- DeclException0 (ExcInternalError);
-
- private:
- /**
- * Store the coarsest level this
- * vertex lives on. This is used
- * as an offset when accessing the
- * dofs of a specific level.
- */
- unsigned int coarsest_level;
-
- /**
- * Finest level this level lives on.
- * This is mostly used for error
- * checking but can also be accessed
- * by the function #get_finest_level#.
- */
- unsigned int finest_level;
-
- /**
- * Array holding the indices.
- */
- unsigned int *indices;
- };
-
-
-/**
- * Distribute dofs on the given cell,
- * with new dofs starting with index
- * #next_free_dof#. Return the next
- * unused index number. The finite
- * element used is the one given to
- * #distribute_dofs#, which is copied
- * to #selected_fe#.
- *
- * This function is excluded from the
- * #distribute_dofs# function since
- * it can not be implemented dimension
- * independent.
- *
- * Note that unlike for the usual dofs,
- * here all cells and not only active
- * ones are allowed.
- */
- unsigned int distribute_dofs_on_cell (cell_iterator &cell,
- unsigned int next_free_dof);
-
- /**
- * Reserve enough space for the MG dof
- * indices for a given triangulation.
- */
- void reserve_space ();
-
- /**
- * Space to store the DoF numbers for the
- * different levels. Unlike the #levels#
- * object in the #DoFHandler#, these are
- * not global numbers but rather are
- * numbers which start from zero on each
- * level.
- */
- vector<DoFLevel<dim>*> mg_levels;
-
- /**
- * For each vertex there is a list of
- * indices of the degrees of freedom indices
- * on the different levels it lives on and
- * which are these levels.
- */
- vector<MGVertexDoFs> mg_vertex_dofs;
-
- /**
- * Vectors storing the number of degrees of
- * freedom on each level.
- */
- vector<unsigned int> mg_used_dofs;
-
-#if (__GNUC__==2) && (__GNUC_MINOR__ < 95)
- // this seems to be disallowed
- // by the standard, so gcc2.95
- // does not accept it
- friend class MGDoFObjectAccessor<1, dim>;
- friend class MGDoFObjectAccessor<2, dim>;
- friend class MGDoFObjectAccessor<3, dim>;
-#else
- // this, however, may grant
- // access to too many classes,
- // but ...
- template <int dim1, int dim2> friend class MGDoFObjectAccessor;
-#endif
-};
-
-
-/* ----------------------- Inline functions of MGVertexDoFs -------------------*/
-
-
-template <int dim>
-inline
-void MGDoFHandler<dim>::MGVertexDoFs::set_index (const unsigned int level,
- const unsigned int dof_number,
- const unsigned int dofs_per_vertex,
- const unsigned int index) {
- Assert ((level >= coarsest_level) && (level <= finest_level),
- ExcInvalidLevel(level));
- Assert (dof_number < dofs_per_vertex,
- ExcIndexRange(dof_number, 0, dofs_per_vertex));
-
- indices[(level-coarsest_level)*dofs_per_vertex + dof_number] = index;
-};
-
-
-template <int dim>
-inline
-unsigned int
-MGDoFHandler<dim>::MGVertexDoFs::get_index (const unsigned int level,
- const unsigned int dof_number,
- const unsigned int dofs_per_vertex) const {
- Assert ((level >= coarsest_level) && (level <= finest_level),
- ExcInvalidLevel(level));
- Assert (dof_number < dofs_per_vertex,
- ExcIndexRange (dof_number, 0, dofs_per_vertex));
-
- return indices[(level-coarsest_level)*dofs_per_vertex + dof_number];
-};
-
-
-/*---------------------------- mg_dof.h ---------------------------*/
-
-#endif
-/*---------------------------- mg_dof.h ---------------------------*/
+++ /dev/null
-//---------------------------- mg_dof_tools.h ---------------------------
-// $Id$
-// Version: $Name$
-//
-// Copyright (C) 1998, 1999, 2000 by the deal.II authors
-//
-// This file is subject to QPL and may not be distributed
-// without copyright and license information. Please refer
-// to the file deal.II/doc/license.html for the text and
-// further information on this license.
-//
-//---------------------------- mg_dof_tools.h ---------------------------
-#ifndef __deal2__mg_dof_tools_h
-#define __deal2__mg_dof_tools_h
-
-
-#include <lac/forward_declarations.h>
-#include <grid/forward_declarations.h>
-
-//TODO: Consider incorporating these functions in DoFTools
-
-/**
- * This is a collection of functions operating on, and manipulating
- * the numbers of degrees of freedom in a multilevel triangulation. It
- * is similar in purpose and function to the #DoFTools# class, but
- * operates on #MGDoFHandler# objects instead of #DoFHandler#
- * objects. See there and the documentation of the member functions
- * for more information.
- *
- * All member functions are static, so there is no need to create an
- * object of class #MGDoFTools#.
- *
- * @author Wolfgang Bangerth and others, 1999
- */
-class MGDoFTools
-{
- public:
- /**
- * Write the sparsity structure
- * of the matrix belonging to the
- * specified #level#. The sparsity pattern
- * is not compressed, so before
- * creating the actual matrix
- * you have to compress the
- * matrix yourself, using
- * #SparseMatrixStruct::compress()#.
- *
- * There is no need to consider
- * hanging nodes here, since only
- * one level is considered.
- */
- template <int dim>
- static void
- make_sparsity_pattern (const MGDoFHandler<dim> &dof_handler,
- SparsityPattern &sparsity,
- const unsigned int level);
-
- /**
- * Make a sparsity pattern including fluxes
- * of discontinuous Galerkin methods.
- * @see make_sparsity_pattern
- * $see DoFTools
- */
- template <int dim>
- static void
- make_flux_sparsity_pattern (const MGDoFHandler<dim> &dof_handler,
- SparsityPattern &sparsity,
- const unsigned int level);
-
-
-};
-
-
-#endif
+++ /dev/null
-//---------------------------- mg_smoother.h ---------------------------
-// $Id$
-// Version: $Name$
-//
-// Copyright (C) 1998, 1999, 2000 by the deal.II authors
-//
-// This file is subject to QPL and may not be distributed
-// without copyright and license information. Please refer
-// to the file deal.II/doc/license.html for the text and
-// further information on this license.
-//
-//---------------------------- mg_smoother.h ---------------------------
-#ifndef __deal2__mg_smoother_h
-#define __deal2__mg_smoother_h
-
-
-#include <lac/forward_declarations.h>
-#include <multigrid/mg_base.h>
-#include <grid/forward_declarations.h>
-#include <base/smartpointer.h>
-#include <vector>
-
-
-/**
- * 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. Does
- * nothing in particular, but
- * needs to be declared anyhow.
- */
- virtual ~MGSmootherBase();
-
- /**
- * Smooth the residual of #u# on
- * the given level. If $S$ is the
- * smoothing operator, then this
- * function should do the
- * following operation:
- * $u += S (rhs - Au)$, where #u# and
- * #rhs# are the input parameters.
- *
- * 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<double> &u,
- const Vector<double> &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
- */
-class MGSmootherIdentity : public MGSmootherBase
-{
- public:
- /**
- * Implementation of the
- * interface in #MGSmootherBase#.
- * This function does nothing,
- * which by comparison with the
- * definition of this function
- * means that the the smoothing
- * operator equals the null
- * operator.
- */
- virtual void smooth (const unsigned int level,
- Vector<double> &u,
- const Vector<double> &rhs) const;
-};
-
-
-/**
- * Base class for multigrid smoothers. It declares the interface
- * to smoothers by inheriting #MGSmootherBase# and implements some functionality for setting
- * the values
- * of vectors at interior boundaries (i.e. boundaries between differing
- * levels of the triangulation) to zero, by building a list of these degrees
- * of freedom's indices at construction time.
- *
- * @author Wolfgang Bangerth, Guido Kanschat, 1999
- */
-class MGSmoother : public MGSmootherBase
-{
- private:
- /**
- * Default constructor. Made private
- * to prevent it being called, which
- * is necessary since this could
- * circumpass the set-up of the list
- * if interior boundary dofs.
- */
- MGSmoother ();
-
- public:
-
- /**
- * Constructor. This one collects
- * the indices of the degrees of freedom
- * on the interior boundaries between
- * the different levels, which are
- * needed by the function
- * #set_zero_interior_boundaries#.
- *
- * Since this function is implemented
- * a bit different in 1d (there are no
- * faces of cells, just vertices),
- * there are actually two sets of
- * constructors, namely this one for 1d
- * and the following one for all other
- * dimensions.
- */
- MGSmoother (const MGDoFHandler<1> &mg_dof,
- unsigned int steps);
-
- /**
- * Constructor. This one collects
- * the indices of the degrees of freedom
- * on the interior boundaries between
- * the different levels, which are
- * needed by the function
- * #set_zero_interior_boundaries#.
- *
- * The parameter steps indicates the number of smoothing
- * steps to be executed by #smooth#.
- */
- template <int dim>
- MGSmoother (const MGDoFHandler<dim> &mg_dof,
- unsigned int steps);
-
- /**
- * Reset the values of the degrees of
- * freedom on interior boundaries between
- * different levels to zero in the given
- * data vector #u#.
- *
- * Since the coarsest level (#level==0#)
- * has no interior boundaries, this
- * function does nothing in this case.
- */
- void set_zero_interior_boundary (const unsigned int level,
- Vector<double> &u) const;
- /**
- * Modify the number of smoothing steps.
- */
- void set_steps(unsigned int steps);
- /**
- * How many steps should be used?
- */
- unsigned int get_steps() const;
-
- private:
- /**
- * Number of smoothing steps.
- */
- unsigned int steps;
- /**
- * For each level, we store a list of
- * degree of freedom indices which are
- * located on interior boundaries between
- * differing levels of the triangulation.
- * Since the coarsest level has no
- * interior boundary dofs, the first
- * entry refers to the second level.
- *
- * These arrays are set by the constructor.
- * The entries for each level are sorted ascendingly.
- */
- vector<vector<unsigned int> > interior_boundary_dofs;
-};
-
-
-/**
- * Implementation of a smoother using matrix builtin relaxation methods.
- *
- * @author Wolfgang Bangerth, Guido Kanschat, 1999
- */
-template<typename number>
-class MGSmootherRelaxation : public MGSmoother
-{
- public:
- /**
- * Type of the smoothing
- * function of the matrix.
- */
- typedef void
- (SparseMatrix<number>::* function_ptr)(Vector<double>&, const Vector<double>&,
- typename SparseMatrix<number>::value_type) const;
-
- /**
- * Constructor.
- * The constructor uses an #MGDoFHandler# to initialize
- * data structures for interior level boundary handling
- * in #MGSmoother#.
- *
- * Furthermore, it takes a pointer to the
- * level matrices and their smoothing function.
- * This function must perform one relaxation step
- * like #SparseMatrix<number>::SOR_step# does. Do not
- * use the preconditioning methods because they apply
- * a preconditioning matrix to the residual.
- *
- * The final two parameters are the number of relaxation
- * steps and the relaxation parameter.
- */
-
- template<int dim>
- MGSmootherRelaxation(const MGDoFHandler<dim> &mg_dof,
- const MGLevelObject<SparseMatrix<number> > &matrix,
- const function_ptr function,
- const unsigned int steps,
- const double omega = 1.);
-
- /**
- * Implementation of the interface in #MGSmootherBase#.
- * We use the SOR method in #SparseMatrix# for the real work
- * and find a way to keep the boundary values.
- */
- virtual void smooth (const unsigned int level,
- Vector<double> &u,
- const Vector<double> &rhs) const;
- private:
- /**
- * Pointer to the matrices.
- */
- SmartPointer<const MGLevelObject<SparseMatrix<number> > > matrix;
-
- /**
- * Pointer to the relaxation function.
- */
- function_ptr relaxation;
-
- /**
- * Relaxation parameter.
- */
- double omega;
-
- /**
- * Auxiliary vector.
- */
- mutable Vector<double> h1;
-
- /**
- * Auxiliary vector.
- */
- mutable Vector<double> h2;
-};
-
-
-/* ------------------------------- Inline functions -------------------------- */
-
-inline
-void
-MGSmoother::set_steps(unsigned int i)
-{
- steps = i;
-}
-
-inline
-unsigned int
-MGSmoother::get_steps() const
-{
- return steps;
-}
-
-
-#endif
-
+++ /dev/null
-// $Id$
-
-
-template<typename number>
-template<int dim>
-MGSmootherRelaxation<number>::MGSmootherRelaxation (const MGDoFHandler<dim> &mg_dof,
- const MGLevelObject<SparseMatrix<number> > &matrix,
- const function_ptr relaxation,
- const unsigned int steps,
- const double omega)
- :
- MGSmoother(mg_dof, steps),
- matrix(&matrix),
- relaxation(relaxation),
- omega(omega)
-{};
-
-
-template<typename number>
-void
-MGSmootherRelaxation<number>::smooth (const unsigned int level,
- Vector<double> &u,
- const Vector<double> &rhs) const
-{
- h1.reinit(u);
- h2.reinit(u);
- for(unsigned i=0;i<get_steps();++i)
- {
- (*matrix)[level].residual(h1, u, rhs);
- ((*matrix)[level].*relaxation)(h2,h1,omega);
- set_zero_interior_boundary(level,h2);
- u.add(h2);
- }
-}
-
+++ /dev/null
-//---------------------------- multigrid.h ---------------------------
-// $Id$
-// Version: $Name$
-//
-// Copyright (C) 1998, 1999, 2000 by the deal.II authors
-//
-// This file is subject to QPL and may not be distributed
-// without copyright and license information. Please refer
-// to the file deal.II/doc/license.html for the text and
-// further information on this license.
-//
-//---------------------------- multigrid.h ---------------------------
-#ifndef __deal2__multigrid_h
-#define __deal2__multigrid_h
-
-
-#include <base/subscriptor.h>
-#include <base/smartpointer.h>
-#include <lac/forward_declarations.h>
-#include <grid/forward_declarations.h>
-#include <lac/sparse_matrix.h>
-#include <lac/vector.h>
-#include <multigrid/mg_base.h>
-
-#include <vector>
-
-
-/**
- * Implementation of multigrid with matrices.
- * While #MGBase# 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 Multigrid : public MGBase
-{
- public:
- /**
- * Constructor. Take the object
- * describing the multilevel
- * hierarchy of degrees of
- * freedom, an object describing
- * constraints to degrees of
- * freedom on the global grid
- * (for example for hanging
- * nodes), a set of matrices
- * which describe the equation
- * discretized on each level,
- * their respective sparsity
- * patterns, and an object
- * mediating the transfer of
- * solutions between different
- * levels.
- *
- * This function already
- * initializes the vectors which
- * will be used later on in the
- * course of the
- * computations. You should
- * therefore create objects of
- * this type as late as possible.
- */
-//TODO: minlevel, maxlevel?
- Multigrid(const MGDoFHandler<dim> &mg_dof_handler,
- const ConstraintMatrix &constraints,
- const MGLevelObject<SparsityPattern> &level_sparsities,
- const MGLevelObject<SparseMatrix<double> > &level_matrices,
- const MGTransferBase &transfer,
- const unsigned int minlevel = 0,
- const unsigned int maxlevel = 1000000);
-
- /**
- * Transfer from a vector on the
- * global grid to vectors defined
- * on each of the levels
- * separately, i.a. an #MGVector#.
- */
- template<typename number>
- void copy_to_mg (const Vector<number> &src);
-
- /**
- * Transfer from multi-level vector to
- * normal vector.
- *
- * Copies data from active
- * portions of an MGVector into
- * the respective positions of a
- * #Vector<double>#. In order to
- * keep the result consistent,
- * constrained degrees of freedom
- * are set to zero.
- */
- template<typename number>
- void copy_from_mg (Vector<number> &dst) const;
-
- /**
- * Negative #vmult# on a level.
-//TODO: what does this function do?
- * @see MGBase.
- */
- virtual void level_vmult (const unsigned int level,
- Vector<double> &dest,
- const Vector<double> &src,
- const Vector<double> &rhs);
-
- /**
- * Read-only access to level matrices.
- */
- const SparseMatrix<double>& get_matrix (const unsigned int level) const;
-
- private:
- /**
- * Associated #MGDoFHandler#.
- */
- SmartPointer<const MGDoFHandler<dim> > mg_dof_handler;
-
-
-/**
- * Sparsity patterns for each level.
- */
- SmartPointer<const MGLevelObject<SparsityPattern> > level_sparsities;
-
- /**
- * Matrices for each level. The
- * matrices are prepared by the
- * constructor of #Multigrid# and
- * can be accessed for
- * assembling.
- */
- SmartPointer<const MGLevelObject<SparseMatrix<double> > > level_matrices;
-
- /**
- * Pointer to the object holding
- * constraints.
- */
- SmartPointer<const ConstraintMatrix> constraints;
-};
-
-
-/**
- * Implementation of the #MGTransferBase# interface for which the transfer
- * operations are prebuilt upon construction of the object of this class as
- * matrices. This is the fast way, since it only needs to build the operation
- * once by looping over all cells and storing the result in a matrix for
- * each level, but requires additional memory.
- *
- * @author Wolfgang Bangerth, Guido Kanschat, 1999
- */
-class MGTransferPrebuilt : public MGTransferBase
-{
- public:
- /**
- * Destructor.
- */
- virtual ~MGTransferPrebuilt ();
-
- /**
- * Actually build the prolongation
- * matrices for each level.
- */
- template <int dim>
- void build_matrices (const MGDoFHandler<dim> &mg_dof);
-
- /**
- * Prolongate a vector from level
- * #to_level-1# to level
- * #to_level#. The previous
- * content of #dst# is
- * overwritten.
- *
- * #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<double> &dst,
- const Vector<double> &src) const;
-
- /**
- * Restrict a vector from level
- * #from_level# to level
- * #from_level-1# and add this
- * restriction to
- * #dst#. Obviously, if the
- * refined region on level
- * #from_level# is smaller than
- * that on level #from_level-1#,
- * some degrees of freedom in
- * #dst# are not covered and will
- * not be altered. For the other
- * degress of freedom, the result
- * of the restriction is added.
- *
- * #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_and_add (const unsigned int from_level,
- Vector<double> &dst,
- const Vector<double> &src) const;
-
- private:
-
- vector<SparsityPattern> prolongation_sparsities;
-
- /**
- * The actual prolongation matrix.
- * column indices belong to the
- * dof indices of the mother cell,
- * i.e. the coarse level.
- * while row indices belong to the
- * child cell, i.e. the fine level.
- */
- vector<SparseMatrix<float> > prolongation_matrices;
-};
-
-
-/* --------------------------- inline functions --------------------- */
-
-
-template<int dim>
-inline
-const SparseMatrix<double>&
-Multigrid<dim>::get_matrix (const unsigned int level) const
-{
- Assert((level>=minlevel) && (level<maxlevel),
- ExcIndexRange(level, minlevel, maxlevel));
-
- return (*level_matrices)[level];
-}
-
-
-#endif
+++ /dev/null
-//---------------------------- multigrid.templates.h ---------------------------
-// $Id$
-// Version: $Name$
-//
-// Copyright (C) 1998, 1999, 2000 by the deal.II authors
-//
-// This file is subject to QPL and may not be distributed
-// without copyright and license information. Please refer
-// to the file deal.II/doc/license.html for the text and
-// further information on this license.
-//
-//---------------------------- multigrid.templates.h ---------------------------
-#ifndef __deal2__multigrid_templates_h
-#define __deal2__multigrid_templates_h
-
-
-#include <dofs/dof_constraints.h>
-#include <multigrid/multigrid.h>
-#include <algorithm>
-#include <fstream>
-
-#include <lac/sparse_matrix.h>
-
-
-//TODO: This function needs to be specially implemented, since in 2d mode we use faces
-#if deal_II_dimension == 1
-
-template <int dim>
-template <typename number>
-void
-Multigrid<dim>::copy_to_mg (const Vector<number>&)
-{
- Assert(false, ExcNotImplemented());
-}
-
-#else
-
-template <int dim>
-template <typename number>
-void
-Multigrid<dim>::copy_to_mg (const Vector<number>& src)
-{
- const unsigned int dofs_per_cell = mg_dof_handler->get_fe().dofs_per_cell;
- const unsigned int dofs_per_face = mg_dof_handler->get_fe().dofs_per_face;
-
- // set the elements of the vectors
- // on all levels to zero
- defect.clear();
-// constraints->condense(src);
-
- vector<unsigned int> global_dof_indices (dofs_per_cell);
- vector<unsigned int> level_dof_indices (dofs_per_cell);
- vector<unsigned int> level_face_indices (dofs_per_face);
-
- // initialize the objects with
- // their correct size
- for (unsigned int level=minlevel; level<=maxlevel ; ++level)
- {
- const unsigned int system_size = (*level_matrices)[level].m();
-
- solution[level].reinit(system_size);
- defect[level].reinit(system_size);
- };
-
- // traverse the grid top-down
- // (i.e. starting with the most
- // refined grid). this way, we can
- // always get that part of one
- // level of the output vector which
- // corresponds to a region which is
- // more refined, by restriction of
- // the respective vector on the
- // next finer level, which we then
- // already have built.
- for (int level=maxlevel; level>=static_cast<int>(minlevel); --level)
- {
- DoFHandler<dim>::active_cell_iterator
- global_cell = mg_dof_handler->DoFHandler<dim>::begin_active(level);
- MGDoFHandler<dim>::active_cell_iterator
- level_cell = mg_dof_handler->begin_active(level);
- const MGDoFHandler<dim>::active_cell_iterator
- level_end = mg_dof_handler->end_active(level);
-
- for (; level_cell!=level_end; ++level_cell, ++global_cell)
- {
- // get the dof numbers of
- // this cell for the global
- // and the level-wise
- // numbering
- global_cell->get_dof_indices(global_dof_indices);
- level_cell->get_mg_dof_indices (level_dof_indices);
-
- // transfer the global
- // defect in the vector
- // into the level-wise one
- for (unsigned int i=0; i<dofs_per_cell; ++i)
- defect[level](level_dof_indices[i]) = src(global_dof_indices[i]);
-
-//TODO: what happens here?
- // Delete values on refinement edge
- for (unsigned int face_n=0; face_n<GeometryInfo<dim>::faces_per_cell; ++face_n)
- {
- const MGDoFHandler<dim>::face_iterator face = level_cell->face(face_n);
- if (face->has_children())
- {
- face->get_mg_dof_indices(level_face_indices);
- for (unsigned int i=0; i<dofs_per_face; ++i)
- defect[level](level_face_indices[i]) = 0.;
- };
- };
- };
-
- // for that part of the level
- // which is further refined:
- // get the defect by
- // restriction of the defect on
- // one level higher
- if (static_cast<unsigned int>(level) < maxlevel)
- transfer->restrict_and_add (level+1, defect[level], defect[level+1]);
- };
-};
-
-#endif
-
-
-template <int dim>
-template <typename number>
-void
-Multigrid<dim>::copy_from_mg(Vector<number> &dst) const
-{
- const unsigned int dofs_per_cell = mg_dof_handler->get_fe().dofs_per_cell;
-
- vector<unsigned int> global_dof_indices (dofs_per_cell);
- vector<unsigned int> level_dof_indices (dofs_per_cell);
-
- DoFHandler<dim>::active_cell_iterator
- global_cell = mg_dof_handler->DoFHandler<dim>::begin_active();
- MGDoFHandler<dim>::active_cell_iterator
- level_cell = mg_dof_handler->begin_active();
- const MGDoFHandler<dim>::active_cell_iterator
- endc = mg_dof_handler->end();
-
- // traverse all cells and copy the
- // data appropriately to the output
- // vector
- for (; level_cell != endc; ++level_cell, ++global_cell)
- {
- const unsigned int level = level_cell->level();
-
- // get the dof numbers of
- // this cell for the global
- // and the level-wise
- // numbering
- global_cell->get_dof_indices (global_dof_indices);
- level_cell->get_mg_dof_indices(level_dof_indices);
-
- // copy level-wise data to
- // global vector
- for (unsigned int i=0; i<dofs_per_cell; ++i)
- dst(global_dof_indices[i]) = solution[level](level_dof_indices[i]);
- };
-
- // clear constrained nodes
- constraints->set_zero(dst);
-}
-
-
-#endif
+++ /dev/null
-//---------------------------- mg_base.cc ---------------------------
-// $Id$
-// Version: $Name$
-//
-// Copyright (C) 1998, 1999, 2000 by the deal.II authors
-//
-// This file is subject to QPL and may not be distributed
-// without copyright and license information. Please refer
-// to the file deal.II/doc/license.html for the text and
-// further information on this license.
-//
-//---------------------------- mg_base.cc ---------------------------
-
-
-#include <multigrid/mg_base.h>
-#include <multigrid/mg_smoother.h>
-#include <iostream>
-#include <cmath>
-
-
-//TODO: this function is only for debugging purposes and should be removed sometimes
-template<class VECTOR>
-static
-void print_vector(ostream& s, const VECTOR& v, const char* text)
-{
- const unsigned int n = (unsigned int)(sqrt(v.size())+.3);
- unsigned int k=0;
-
- s << endl << "splot '-' title '" << text << "'" << endl;
-
- // write the vector entries in a
- // kind of square
- for (unsigned int i=0;i<n;++i)
- {
- for (unsigned int j=0;j<n;++j)
- s << '\n' << v(k++);
- s << endl;
- }
- s << "e\npause -1" << endl;
-}
-
-
-//////////////////////////////////////////////////////////////////////
-
-
-MGBase::~MGBase ()
-{};
-
-
-MGBase::MGBase(const MGTransferBase &transfer,
- const unsigned minlevel,
- const unsigned maxlevel)
- :
- maxlevel(maxlevel),
- minlevel(minlevel),
- defect(minlevel,maxlevel),
- solution(minlevel,maxlevel),
- transfer(&transfer)
-{
- Assert(minlevel <= maxlevel,
- ExcSwitchedLevels(minlevel, maxlevel));
-};
-
-
-void
-MGBase::vcycle(const MGSmootherBase &pre_smooth,
- const MGSmootherBase &post_smooth,
- const MGCoarseGridSolver &coarse_grid_solver)
-{
-// static int k=0;
-// cout << "set title 'cycle " << ++k << "'\n";
-
- level_mgstep(maxlevel, pre_smooth, post_smooth, coarse_grid_solver);
-};
-
-
-void
-MGBase::level_mgstep(const unsigned int level,
- const MGSmootherBase &pre_smooth,
- const MGSmootherBase &post_smooth,
- const MGCoarseGridSolver &coarse_grid_solver)
-{
- solution[level] = 0.;
-
- if (level == minlevel)
- {
- coarse_grid_solver(level, solution[level], defect[level]);
- return;
- }
-
- // smoothing of the residual by modifying s
- pre_smooth.smooth(level, solution[level], defect[level]);
- // t = d-As
- t.reinit(solution[level].size());
- level_vmult(level, t, solution[level], defect[level]);
-// print_vector(cout,t,"T");
-
- // make t rhs of lower level
-//TODO: this function adds the restricted t to defect[level-1].
-//TODO: why don't we have to clear it before?
- transfer->restrict_and_add (level, defect[level-1], t);
-// print_vector(cout,defect[level-1],"Dl-1");
-
- // do recursion
- level_mgstep(level-1, pre_smooth, post_smooth, coarse_grid_solver);
-
- // reset size of the auxiliary
- // vector, since it has been
- // resized in the recursive call to
- // level_mgstep directly above
- t.reinit(solution[level].size());
-
- // do coarse grid correction
-// print_vector(cout,solution[level-1],"Sl-1");
- transfer->prolongate(level, t, solution[level-1]);
-// print_vector(cout,t,"T");
- solution[level] += t;
-
- // smoothing (modify solution again)
-//TODO: what happens here? smooth overwrites the solution[level],
-//TODO: so the previous two statements should have no effect. No?
- post_smooth.smooth(level, solution[level], defect[level]);
-};
-
-
-//////////////////////////////////////////////////////////////////////
-
-MGCoarseGridSolver::~MGCoarseGridSolver()
-{};
-
-
-//////////////////////////////////////////////////////////////////////
-
-MGTransferBase::~MGTransferBase()
-{};
-
-
+++ /dev/null
-//---------------------------- mg_dof_accessor.cc ---------------------------
-// $Id$
-// Version: $Name$
-//
-// Copyright (C) 1998, 1999, 2000 by the deal.II authors
-//
-// This file is subject to QPL and may not be distributed
-// without copyright and license information. Please refer
-// to the file deal.II/doc/license.html for the text and
-// further information on this license.
-//
-//---------------------------- mg_dof_accessor.cc ---------------------------
-
-
-#include <dofs/dof_levels.h>
-#include <multigrid/mg_dof_accessor.h>
-#include <multigrid/mg_dof_handler.h>
-#include <grid/tria_iterator.h>
-#include <grid/tria_iterator.templates.h>
-#include <fe/fe.h>
-
-#include <lac/vector.h>
-#include <lac/full_matrix.h>
-#include <lac/sparse_matrix.h>
-
-
-/* ------------------------ MGDoFLineAccessor --------------------------- */
-
-template <int dim>
-MGDoFObjectAccessor<1, dim>::MGDoFObjectAccessor (Triangulation<dim> *tria,
- const int level,
- const int index,
- const AccessorData *local_data) :
- MGDoFAccessor<dim> (local_data),
- MGDoFObjectAccessor_Inheritance<1,dim>::BaseClass(tria,level,index,local_data)
-{};
-
-
-template <int dim>
-unsigned int MGDoFObjectAccessor<1, dim>::mg_dof_index (const unsigned int i) const
-{
- Assert (dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (mg_dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- // make sure a FE has been selected
- // and enough room was reserved
- Assert (&dof_handler->get_fe() != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (i<dof_handler->get_fe().dofs_per_line,
- ExcIndexRange (i, 0, dof_handler->get_fe().dofs_per_line));
-
- return mg_dof_handler->mg_levels[present_level]
- ->line_dofs[present_index*dof_handler->get_fe().dofs_per_line+i];
-};
-
-
-template <int dim>
-void MGDoFObjectAccessor<1, dim>::set_mg_dof_index (const unsigned int i,
- const unsigned int index) const
-{
- Assert (dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (mg_dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- // make sure a FE has been selected
- // and enough room was reserved
- Assert (&dof_handler->get_fe() != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (i<dof_handler->get_fe().dofs_per_line,
- ExcIndexRange (i, 0, dof_handler->get_fe().dofs_per_line));
-
- mg_dof_handler->mg_levels[present_level]
- ->line_dofs[present_index*dof_handler->get_fe().dofs_per_line+i] = index;
-};
-
-
-template <int dim>
-unsigned int MGDoFObjectAccessor<1, dim>::mg_vertex_dof_index (const unsigned int vertex,
- const unsigned int i) const
-{
- Assert (dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (mg_dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (&dof_handler->get_fe() != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (vertex<2, ExcIndexRange (i,0,2));
- Assert (i<dof_handler->get_fe().dofs_per_vertex,
- ExcIndexRange (i, 0, dof_handler->get_fe().dofs_per_vertex));
-
- return (mg_dof_handler->mg_vertex_dofs[vertex_index(vertex)]
- .get_index (present_level, i, dof_handler->get_fe().dofs_per_vertex));
-};
-
-
-template <int dim>
-void MGDoFObjectAccessor<1, dim>::set_mg_vertex_dof_index (const unsigned int vertex,
- const unsigned int i,
- const unsigned int index) const
-{
- Assert (dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (mg_dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (&dof_handler->get_fe() != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (vertex<2, ExcIndexRange (i,0,2));
- Assert (i<dof_handler->get_fe().dofs_per_vertex,
- ExcIndexRange (i, 0, dof_handler->get_fe().dofs_per_vertex));
-
- mg_dof_handler->mg_vertex_dofs[vertex_index(vertex)]
- .set_index (present_level, i, dof_handler->get_fe().dofs_per_vertex, index);
-};
-
-
-template <int dim>
-void
-MGDoFObjectAccessor<1, dim>::get_mg_dof_indices (vector<unsigned int> &dof_indices) const
-{
- Assert (dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (mg_dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (&dof_handler->get_fe() != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (dof_indices.size() == (2*dof_handler->get_fe().dofs_per_vertex +
- dof_handler->get_fe().dofs_per_line),
- DoFAccessor<dim>::ExcVectorDoesNotMatch());
-
- const unsigned int dofs_per_vertex = dof_handler->get_fe().dofs_per_vertex,
- dofs_per_line = dof_handler->get_fe().dofs_per_line;
- vector<unsigned int>::iterator next = dof_indices.begin();
- for (unsigned int vertex=0; vertex<2; ++vertex)
- for (unsigned int d=0; d<dofs_per_vertex; ++d)
- *next++ = mg_vertex_dof_index(vertex,d);
- for (unsigned int d=0; d<dofs_per_line; ++d)
- *next++ = mg_dof_index(d);
-
- Assert (next == dof_indices.end(),
- ExcInternalError());
-};
-
-
-template <int dim>
-template <typename number>
-void
-MGDoFObjectAccessor<1,dim>::get_mg_dof_values (const Vector<number> &values,
- Vector<number> &dof_values) const
-{
- Assert (dof_handler != 0, DoFAccessor<1>::DoFAccessor<1>::ExcInvalidObject());
- Assert (mg_dof_handler != 0, DoFAccessor<1>::DoFAccessor<1>::ExcInvalidObject());
- Assert (&dof_handler->get_fe() != 0, DoFAccessor<1>::DoFAccessor<1>::ExcInvalidObject());
- Assert (dof_values.size() == dof_handler->get_fe().dofs_per_cell,
- DoFAccessor<1>::ExcVectorDoesNotMatch());
- Assert (values.size() == dof_handler->n_dofs(),
- DoFAccessor<1>::ExcVectorDoesNotMatch());
-
- const unsigned int dofs_per_vertex = dof_handler->get_fe().dofs_per_vertex,
- dofs_per_line = dof_handler->get_fe().dofs_per_line;
- vector<number>::iterator next_dof_value=dof_values.begin();
- for (unsigned int vertex=0; vertex<2; ++vertex)
- for (unsigned int d=0; d<dofs_per_vertex; ++d)
- *next_dof_value++ = values(mg_vertex_dof_index(vertex,d));
- for (unsigned int d=0; d<dofs_per_line; ++d)
- *next_dof_value++ = values(mg_dof_index(d));
-
- Assert (next_dof_value == dof_values.end(),
- ExcInternalError());
-};
-
-
-template <int dim>
-TriaIterator<dim,MGDoFObjectAccessor<1, dim> >
-MGDoFObjectAccessor<1, dim>::child (const unsigned int i) const
-{
- TriaIterator<dim,MGDoFObjectAccessor<1, dim> > q (tria,
- present_level+1,
- child_index (i),
- mg_dof_handler);
-
-#ifdef DEBUG
- if (q.state() != past_the_end)
- Assert (q->used(), typename TriaAccessor<dim>::ExcUnusedCellAsChild());
-#endif
- return q;
-};
-
-
-template <int dim>
-void
-MGDoFObjectAccessor<1, dim>::copy_from (const MGDoFObjectAccessor<1, dim> &a)
-{
- DoFObjectAccessor<1, dim>::copy_from (a);
- set_mg_dof_handler (a.mg_dof_handler);
-};
-
-
-/* ------------------------ MGDoFQuadAccessor --------------------------- */
-
-template <int dim>
-MGDoFObjectAccessor<2, dim>::MGDoFObjectAccessor (Triangulation<dim> *tria,
- const int level,
- const int index,
- const AccessorData *local_data) :
- MGDoFAccessor<dim> (local_data),
- MGDoFObjectAccessor_Inheritance<2,dim>::BaseClass(tria,level,index,local_data)
-{};
-
-
-template <int dim>
-inline
-unsigned int MGDoFObjectAccessor<2, dim>::mg_dof_index (const unsigned int i) const
-{
- Assert (dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (mg_dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- // make sure a FE has been selected
- // and enough room was reserved
- Assert (&dof_handler->get_fe() != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (i<dof_handler->get_fe().dofs_per_quad,
- ExcIndexRange (i, 0, dof_handler->get_fe().dofs_per_quad));
-
- return mg_dof_handler->mg_levels[present_level]
- ->quad_dofs[present_index*dof_handler->get_fe().dofs_per_quad+i];
-};
-
-
-template <int dim>
-void MGDoFObjectAccessor<2, dim>::set_mg_dof_index (const unsigned int i,
- const unsigned int index) const
-{
- Assert (dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (mg_dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- // make sure a FE has been selected
- // and enough room was reserved
- Assert (&dof_handler->get_fe() != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (i<dof_handler->get_fe().dofs_per_quad,
- ExcIndexRange (i, 0, dof_handler->get_fe().dofs_per_quad));
-
- mg_dof_handler->mg_levels[present_level]
- ->quad_dofs[present_index*dof_handler->get_fe().dofs_per_quad+i] = index;
-};
-
-
-template <int dim>
-inline
-unsigned int MGDoFObjectAccessor<2, dim>::mg_vertex_dof_index (const unsigned int vertex,
- const unsigned int i) const
-{
- Assert (dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (mg_dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (&dof_handler->get_fe() != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (vertex<4, ExcIndexRange (i,0,4));
- Assert (i<dof_handler->get_fe().dofs_per_vertex,
- ExcIndexRange (i, 0, dof_handler->get_fe().dofs_per_vertex));
-
- return (mg_dof_handler->mg_vertex_dofs[vertex_index(vertex)]
- .get_index (present_level, i, dof_handler->get_fe().dofs_per_vertex));
-};
-
-
-template <int dim>
-void MGDoFObjectAccessor<2, dim>::set_mg_vertex_dof_index (const unsigned int vertex,
- const unsigned int i,
- const unsigned int index) const
-{
- Assert (dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (mg_dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (&dof_handler->get_fe() != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (vertex<4, ExcIndexRange (i,0,4));
- Assert (i<dof_handler->get_fe().dofs_per_vertex,
- ExcIndexRange (i, 0, dof_handler->get_fe().dofs_per_vertex));
-
- mg_dof_handler->mg_vertex_dofs[vertex_index(vertex)]
- .set_index (present_level, i, dof_handler->get_fe().dofs_per_vertex, index);
-};
-
-
-template <int dim>
-void
-MGDoFObjectAccessor<2, dim>::get_mg_dof_indices (vector<unsigned int> &dof_indices) const
-{
- Assert (dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (mg_dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (&dof_handler->get_fe() != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (dof_indices.size() == (4*dof_handler->get_fe().dofs_per_vertex +
- 4*dof_handler->get_fe().dofs_per_line +
- dof_handler->get_fe().dofs_per_quad),
- DoFAccessor<2>::ExcVectorDoesNotMatch());
-
- const unsigned int dofs_per_vertex = dof_handler->get_fe().dofs_per_vertex,
- dofs_per_line = dof_handler->get_fe().dofs_per_line,
- dofs_per_quad = dof_handler->get_fe().dofs_per_quad;
- vector<unsigned int>::iterator next = dof_indices.begin();
- for (unsigned int vertex=0; vertex<4; ++vertex)
- for (unsigned int d=0; d<dofs_per_vertex; ++d)
- *next++ = mg_vertex_dof_index(vertex,d);
- for (unsigned int line=0; line<4; ++line)
- for (unsigned int d=0; d<dofs_per_line; ++d)
- *next++ = this->line(line)->mg_dof_index(d);
- for (unsigned int d=0; d<dofs_per_quad; ++d)
- *next++ = mg_dof_index(d);
-
- Assert (next == dof_indices.end(),
- ExcInternalError());
-};
-
-
-template <int dim>
-template <typename number>
-void
-MGDoFObjectAccessor<2,dim>::get_mg_dof_values (const Vector<number> &values,
- Vector<number> &dof_values) const
-{
- Assert (dof_handler != 0, DoFAccessor<2>::ExcInvalidObject());
- Assert (mg_dof_handler != 0, DoFAccessor<2>::ExcInvalidObject());
- Assert (&dof_handler->get_fe() != 0, DoFAccessor<2>::ExcInvalidObject());
- Assert (dof_values.size() == dof_handler->get_fe().dofs_per_cell,
- DoFAccessor<2>::ExcVectorDoesNotMatch());
- Assert (values.size() == mg_dof_handler->n_dofs(present_level),
- DoFAccessor<2>::ExcVectorDoesNotMatch());
-
- const unsigned int dofs_per_vertex = dof_handler->get_fe().dofs_per_vertex,
- dofs_per_line = dof_handler->get_fe().dofs_per_line,
- dofs_per_quad = dof_handler->get_fe().dofs_per_quad;
- vector<number>::iterator next_dof_value=dof_values.begin();
- for (unsigned int vertex=0; vertex<4; ++vertex)
- for (unsigned int d=0; d<dofs_per_vertex; ++d)
- *next_dof_value++ = values(mg_vertex_dof_index(vertex,d));
- for (unsigned int line=0; line<4; ++line)
- for (unsigned int d=0; d<dofs_per_line; ++d)
- *next_dof_value++ = values(this->line(line)->mg_dof_index(d));
- for (unsigned int d=0; d<dofs_per_quad; ++d)
- *next_dof_value++ = values(mg_dof_index(d));
-
- Assert (next_dof_value == dof_values.end(),
- ExcInternalError());
-};
-
-
-template <int dim>
-TriaIterator<dim,MGDoFObjectAccessor<1, dim> >
-MGDoFObjectAccessor<2, dim>::line (const unsigned int i) const
-{
- Assert (i<4, ExcIndexRange (i, 0, 4));
-
- return TriaIterator<dim,MGDoFObjectAccessor<1, dim> >
- (
- tria,
- present_level,
- line_index (i),
- mg_dof_handler
- );
-};
-
-
-template <int dim>
-TriaIterator<dim,MGDoFObjectAccessor<2, dim> >
-MGDoFObjectAccessor<2, dim>::child (const unsigned int i) const
-{
- TriaIterator<dim,MGDoFObjectAccessor<2, dim> > q (tria,
- present_level+1,
- child_index (i),
- mg_dof_handler);
-
-#ifdef DEBUG
- if (q.state() != past_the_end)
- Assert (q->used(), typename TriaAccessor<dim>::ExcUnusedCellAsChild());
-#endif
- return q;
-};
-
-
-template <int dim>
-void
-MGDoFObjectAccessor<2, dim>::copy_from (const MGDoFObjectAccessor<2, dim> &a)
-{
- DoFObjectAccessor<2, dim>::copy_from (a);
- set_mg_dof_handler (a.mg_dof_handler);
-};
-
-
-/* ------------------------ MGDoFHexAccessor --------------------------- */
-
-template <int dim>
-MGDoFObjectAccessor<3, dim>::MGDoFObjectAccessor (Triangulation<dim> *tria,
- const int level,
- const int index,
- const AccessorData *local_data) :
- MGDoFAccessor<dim> (local_data),
- MGDoFObjectAccessor_Inheritance<3,dim>::BaseClass(tria,level,index,local_data)
-{};
-
-
-template <int dim>
-inline
-unsigned int MGDoFObjectAccessor<3, dim>::mg_dof_index (const unsigned int i) const
-{
- Assert (dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (mg_dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- // make sure a FE has been selected
- // and enough room was reserved
- Assert (&dof_handler->get_fe() != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (i<dof_handler->get_fe().dofs_per_hex,
- ExcIndexRange (i, 0, dof_handler->get_fe().dofs_per_hex));
-
- return mg_dof_handler->mg_levels[present_level]
- ->hex_dofs[present_index*dof_handler->get_fe().dofs_per_hex+i];
-};
-
-
-template <int dim>
-void MGDoFObjectAccessor<3, dim>::set_mg_dof_index (const unsigned int i,
- const unsigned int index) const
-{
- Assert (dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (mg_dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- // make sure a FE has been selected
- // and enough room was reserved
- Assert (&dof_handler->get_fe() != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (i<dof_handler->get_fe().dofs_per_hex,
- ExcIndexRange (i, 0, dof_handler->get_fe().dofs_per_hex));
-
- mg_dof_handler->mg_levels[present_level]
- ->hex_dofs[present_index*dof_handler->get_fe().dofs_per_hex+i] = index;
-};
-
-
-template <int dim>
-inline
-unsigned int MGDoFObjectAccessor<3, dim>::mg_vertex_dof_index (const unsigned int vertex,
- const unsigned int i) const
-{
- Assert (dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (mg_dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (&dof_handler->get_fe() != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (vertex<8, ExcIndexRange (i,0,8));
- Assert (i<dof_handler->get_fe().dofs_per_vertex,
- ExcIndexRange (i, 0, dof_handler->get_fe().dofs_per_vertex));
-
- return (mg_dof_handler->mg_vertex_dofs[vertex_index(vertex)]
- .get_index (present_level, i, dof_handler->get_fe().dofs_per_vertex));
-};
-
-
-template <int dim>
-void MGDoFObjectAccessor<3, dim>::set_mg_vertex_dof_index (const unsigned int vertex,
- const unsigned int i,
- const unsigned int index) const
-{
- Assert (dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (mg_dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (&dof_handler->get_fe() != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (vertex<4, ExcIndexRange (i,0,4));
- Assert (i<dof_handler->get_fe().dofs_per_vertex,
- ExcIndexRange (i, 0, dof_handler->get_fe().dofs_per_vertex));
-
- mg_dof_handler->mg_vertex_dofs[vertex_index(vertex)]
- .set_index (present_level, i, dof_handler->get_fe().dofs_per_vertex, index);
-};
-
-
-template <int dim>
-void
-MGDoFObjectAccessor<3, dim>::get_mg_dof_indices (vector<unsigned int> &dof_indices) const
-{
- Assert (dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (mg_dof_handler != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (&dof_handler->get_fe() != 0, DoFAccessor<dim>::ExcInvalidObject());
- Assert (dof_indices.size() == (8*dof_handler->get_fe().dofs_per_vertex +
- 12*dof_handler->get_fe().dofs_per_line +
- 6*dof_handler->get_fe().dofs_per_quad +
- dof_handler->get_fe().dofs_per_hex),
- DoFAccessor<3>::ExcVectorDoesNotMatch());
-
- const unsigned int dofs_per_vertex = dof_handler->get_fe().dofs_per_vertex,
- dofs_per_line = dof_handler->get_fe().dofs_per_line,
- dofs_per_quad = dof_handler->get_fe().dofs_per_quad,
- dofs_per_hex = dof_handler->get_fe().dofs_per_hex;
- vector<unsigned int>::iterator next = dof_indices.begin();
- for (unsigned int vertex=0; vertex<8; ++vertex)
- for (unsigned int d=0; d<dofs_per_vertex; ++d)
- *next++ = mg_vertex_dof_index(vertex,d);
- for (unsigned int line=0; line<12; ++line)
- for (unsigned int d=0; d<dofs_per_line; ++d)
- *next++ = this->line(line)->mg_dof_index(d);
- for (unsigned int quad=0; quad<12; ++quad)
- for (unsigned int d=0; d<dofs_per_quad; ++d)
- *next++ = this->quad(quad)->mg_dof_index(d);
- for (unsigned int d=0; d<dofs_per_hex; ++d)
- *next++ = mg_dof_index(d);
-
- Assert (next == dof_indices.end(),
- ExcInternalError());
-};
-
-
-template <int dim>
-template <typename number>
-void
-MGDoFObjectAccessor<3,dim>::get_mg_dof_values (const Vector<number> &values,
- Vector<number> &dof_values) const
-{
- Assert (dof_handler != 0, DoFAccessor<2>::ExcInvalidObject());
- Assert (mg_dof_handler != 0, DoFAccessor<2>::ExcInvalidObject());
- Assert (&dof_handler->get_fe() != 0, DoFAccessor<2>::ExcInvalidObject());
- Assert (dof_values.size() == dof_handler->get_fe().dofs_per_cell,
- DoFAccessor<3>::ExcVectorDoesNotMatch());
- Assert (values.size() == mg_dof_handler->n_dofs(present_level),
- DoFAccessor<3>::ExcVectorDoesNotMatch());
-
- const unsigned int dofs_per_vertex = dof_handler->get_fe().dofs_per_vertex,
- dofs_per_line = dof_handler->get_fe().dofs_per_line,
- dofs_per_quad = dof_handler->get_fe().dofs_per_quad,
- dofs_per_hex = dof_handler->get_fe().dofs_per_hex;
- vector<number>::iterator next_dof_value=dof_values.begin();
- for (unsigned int vertex=0; vertex<8; ++vertex)
- for (unsigned int d=0; d<dofs_per_vertex; ++d)
- *next_dof_value++ = values(mg_vertex_dof_index(vertex,d));
- for (unsigned int line=0; line<12; ++line)
- for (unsigned int d=0; d<dofs_per_line; ++d)
- *next_dof_value++ = values(this->line(line)->mg_dof_index(d));
- for (unsigned int quad=0; quad<6; ++quad)
- for (unsigned int d=0; d<dofs_per_quad; ++d)
- *next_dof_value++ = values(this->quad(quad)->mg_dof_index(d));
- for (unsigned int d=0; d<dofs_per_hex; ++d)
- *next_dof_value++ = values(dof_index(d));
-
- Assert (next_dof_value == dof_values.end(),
- ExcInternalError());
-};
-
-
-template <int dim>
-TriaIterator<dim,MGDoFObjectAccessor<1, dim> >
-MGDoFObjectAccessor<3, dim>::line (const unsigned int i) const {
- Assert (i<12, ExcIndexRange (i, 0, 12));
-
- return TriaIterator<dim,MGDoFObjectAccessor<1, dim> >
- (
- tria,
- present_level,
- line_index (i),
- mg_dof_handler
- );
-};
-
-
-template <int dim>
-TriaIterator<dim,MGDoFObjectAccessor<2, dim> >
-MGDoFObjectAccessor<3, dim>::quad (const unsigned int i) const {
- Assert (i<12, ExcIndexRange (i, 0, 6));
-
- return TriaIterator<dim,MGDoFObjectAccessor<2, dim> >
- (
- tria,
- present_level,
- quad_index (i),
- mg_dof_handler
- );
-};
-
-
-template <int dim>
-TriaIterator<dim,MGDoFObjectAccessor<3, dim> >
-MGDoFObjectAccessor<3, dim>::child (const unsigned int i) const {
- TriaIterator<dim,MGDoFObjectAccessor<3, dim> > q (tria,
- present_level+1,
- child_index (i),
- mg_dof_handler);
-
-#ifdef DEBUG
- if (q.state() != past_the_end)
- Assert (q->used(), typename TriaAccessor<dim>::ExcUnusedCellAsChild());
-#endif
- return q;
-};
-
-
-template <int dim>
-void
-MGDoFObjectAccessor<3, dim>::copy_from (const MGDoFObjectAccessor<3, dim> &a) {
- DoFObjectAccessor<3, dim>::copy_from (a);
- set_mg_dof_handler (a.mg_dof_handler);
-};
-
-
-/*------------------------- Functions: MGDoFCellAccessor -----------------------*/
-
-template <int dim>
-TriaIterator<dim,MGDoFCellAccessor<dim> >
-MGDoFCellAccessor<dim>::neighbor (const unsigned int i) const {
- TriaIterator<dim,MGDoFCellAccessor<dim> > q (tria,
- neighbor_level (i),
- neighbor_index (i),
- mg_dof_handler);
-
-#ifdef DEBUG
- if (q.state() != past_the_end)
- Assert (q->used(), typename TriaAccessor<dim>::ExcUnusedCellAsNeighbor());
-#endif
- return q;
-};
-
-
-template <int dim>
-TriaIterator<dim,MGDoFCellAccessor<dim> >
-MGDoFCellAccessor<dim>::child (const unsigned int i) const {
- TriaIterator<dim,MGDoFCellAccessor<dim> > q (tria,
- present_level+1,
- child_index (i),
- mg_dof_handler);
-
-#ifdef DEBUG
- if (q.state() != past_the_end)
- Assert (q->used(), typename TriaAccessor<dim>::ExcUnusedCellAsChild());
-#endif
- return q;
-};
-
-
-#if deal_II_dimension == 1
-
-template <>
-MGDoFCellAccessor<1>::face_iterator
-MGDoFCellAccessor<1>::face (const unsigned int) const
-{
- Assert (false, ExcNotUsefulForThisDimension());
- return face_iterator();
-};
-
-#endif
-
-
-#if deal_II_dimension == 2
-
-template <>
-MGDoFCellAccessor<2>::face_iterator
-MGDoFCellAccessor<2>::face (const unsigned int i) const
-{
- return line(i);
-};
-
-#endif
-
-
-#if deal_II_dimension == 3
-
-template <>
-MGDoFCellAccessor<3>::face_iterator
-MGDoFCellAccessor<3>::face (const unsigned int i) const
-{
- return quad(i);
-};
-
-#endif
-
-
-// explicit instantiations
-
-template
-void
-MGDoFObjectAccessor<1,deal_II_dimension>::
-get_mg_dof_values (const Vector<double> &values,
- Vector<double> &dof_values) const;
-
-template
-void
-MGDoFObjectAccessor<1,deal_II_dimension>::
-get_mg_dof_values (const Vector<float> &values,
- Vector<float> &dof_values) const;
-
-
-#if deal_II_dimension >= 2
-
-template
-void
-MGDoFObjectAccessor<2,deal_II_dimension>::
-get_mg_dof_values (const Vector<double> &values,
- Vector<double> &dof_values) const;
-
-template
-void
-MGDoFObjectAccessor<2,deal_II_dimension>::
-get_mg_dof_values (const Vector<float> &values,
- Vector<float> &dof_values) const;
-
-#endif
-
-
-#if deal_II_dimension >= 3
-
-template
-void
-MGDoFObjectAccessor<3,deal_II_dimension>::
-get_mg_dof_values (const Vector<double> &values,
- Vector<double> &dof_values) const;
-
-template
-void
-MGDoFObjectAccessor<3,deal_II_dimension>::
-get_mg_dof_values (const Vector<float> &values,
- Vector<float> &dof_values) const;
-
-#endif
-
-
-#if deal_II_dimension == 1
-template class MGDoFObjectAccessor<1, 1>;
-template class MGDoFCellAccessor<1>;
-
-template class TriaRawIterator<1,MGDoFCellAccessor<1> >;
-template class TriaIterator<1,MGDoFCellAccessor<1> >;
-template class TriaActiveIterator<1,MGDoFCellAccessor<1> >;
-#endif
-
-#if deal_II_dimension == 2
-template class MGDoFObjectAccessor<1, 2>;
-template class MGDoFObjectAccessor<2, 2>;
-template class MGDoFCellAccessor<2>;
-
-template class TriaRawIterator <2,MGDoFObjectAccessor<1, 2> >;
-template class TriaIterator <2,MGDoFObjectAccessor<1, 2> >;
-template class TriaActiveIterator<2,MGDoFObjectAccessor<1, 2> >;
-
-template class TriaRawIterator <2,MGDoFCellAccessor<2> >;
-template class TriaIterator <2,MGDoFCellAccessor<2> >;
-template class TriaActiveIterator<2,MGDoFCellAccessor<2> >;
-#endif
-
-
-#if deal_II_dimension == 3
-template class MGDoFObjectAccessor<1, 3>;
-template class MGDoFObjectAccessor<2, 3>;
-template class MGDoFObjectAccessor<3, 3>;
-template class MGDoFCellAccessor<3>;
-
-template class TriaRawIterator <3,MGDoFObjectAccessor<1, 3> >;
-template class TriaIterator <3,MGDoFObjectAccessor<1, 3> >;
-template class TriaActiveIterator<3,MGDoFObjectAccessor<1, 3> >;
-
-template class TriaRawIterator <3,MGDoFObjectAccessor<2, 3> >;
-template class TriaIterator <3,MGDoFObjectAccessor<2, 3> >;
-template class TriaActiveIterator<3,MGDoFObjectAccessor<2, 3> >;
-
-template class TriaRawIterator <3,MGDoFCellAccessor<3> >;
-template class TriaIterator <3,MGDoFCellAccessor<3> >;
-template class TriaActiveIterator<3,MGDoFCellAccessor<3> >;
-#endif
-
-
+++ /dev/null
-//---------------------------- mg_dof_handler.cc ---------------------------
-// $Id$
-// Version: $Name$
-//
-// Copyright (C) 1998, 1999, 2000 by the deal.II authors
-//
-// This file is subject to QPL and may not be distributed
-// without copyright and license information. Please refer
-// to the file deal.II/doc/license.html for the text and
-// further information on this license.
-//
-//---------------------------- mg_dof_handler.cc ---------------------------
-
-
-#include <dofs/dof_levels.h>
-#include <dofs/dof_constraints.h>
-#include <multigrid/mg_dof_handler.h>
-#include <multigrid/mg_dof_accessor.h>
-#include <grid/tria_levels.h>
-#include <grid/tria_accessor.h>
-#include <grid/tria_iterator.h>
-#include <grid/tria.h>
-#include <grid/geometry_info.h>
-#include <fe/fe.h>
-#include <lac/sparse_matrix.h>
-
-#include <algorithm>
-
-
-/* ------------------------ MGVertexDoFs ----------------------------------- */
-
-template <int dim>
-MGDoFHandler<dim>::MGVertexDoFs::MGVertexDoFs () :
- coarsest_level (1<<30),
- finest_level (0),
- indices (0)
-{};
-
-
-template <int dim>
-void MGDoFHandler<dim>::MGVertexDoFs::init (const unsigned int cl,
- const unsigned int fl,
- const unsigned int dofs_per_vertex) {
- Assert (indices == 0, ExcInternalError());
-
- coarsest_level = cl;
- finest_level = fl;
-
- const unsigned int n_levels = finest_level-coarsest_level + 1;
-
- indices = new unsigned int[n_levels * dofs_per_vertex];
- Assert (indices != 0, ExcNoMemory ());
-
- for (unsigned int i=0; i<n_levels*dofs_per_vertex; ++i)
- indices[i] = DoFHandler<dim>::invalid_dof_index;
-};
-
-
-template <int dim>
-MGDoFHandler<dim>::MGVertexDoFs::~MGVertexDoFs () {
- delete[] indices;
-};
-
-
-template <int dim>
-inline
-unsigned int MGDoFHandler<dim>::MGVertexDoFs::get_coarsest_level () const {
- return coarsest_level;
-};
-
-
-template <int dim>
-inline
-unsigned int MGDoFHandler<dim>::MGVertexDoFs::get_finest_level () const {
- return finest_level;
-};
-
-
-/* ------------------------ MGDoFHandler ------------------------------------- */
-
-template <int dim>
-MGDoFHandler<dim>::MGDoFHandler (Triangulation<dim> &tria) :
- DoFHandler<dim> (tria)
-{};
-
-
-template <int dim>
-MGDoFHandler<dim>::~MGDoFHandler () {};
-
-
-#if deal_II_dimension == 1
-
-template <>
-MGDoFHandler<1>::raw_cell_iterator
-MGDoFHandler<1>::begin_raw (const unsigned int level) const {
- return begin_raw_line (level);
-};
-
-
-template <>
-MGDoFHandler<1>::cell_iterator
-MGDoFHandler<1>::begin (const unsigned int level) const {
- return begin_line (level);
-};
-
-
-template <>
-MGDoFHandler<1>::active_cell_iterator
-MGDoFHandler<1>::begin_active (const unsigned int level) const {
- return begin_active_line (level);
-};
-
-
-template <>
-MGDoFHandler<1>::raw_cell_iterator
-MGDoFHandler<1>::end () const {
- return end_line ();
-};
-
-
-template <>
-MGDoFHandler<1>::raw_cell_iterator
-MGDoFHandler<1>::last_raw () const {
- return last_raw_line ();
-};
-
-
-template <>
-MGDoFHandler<1>::raw_cell_iterator
-MGDoFHandler<1>::last_raw (const unsigned int level) const {
- return last_raw_line (level);
-};
-
-
-template <>
-MGDoFHandler<1>::cell_iterator
-MGDoFHandler<1>::last () const {
- return last_line ();
-};
-
-
-template <>
-MGDoFHandler<1>::cell_iterator
-MGDoFHandler<1>::last (const unsigned int level) const {
- return last_line (level);
-};
-
-
-template <>
-MGDoFHandler<1>::active_cell_iterator
-MGDoFHandler<1>::last_active () const {
- return last_active_line ();
-};
-
-
-template <>
-MGDoFHandler<1>::active_cell_iterator
-MGDoFHandler<1>::last_active (const unsigned int level) const {
- return last_active_line (level);
-};
-
-
-template <>
-MGDoFDimensionInfo<1>::raw_face_iterator
-MGDoFHandler<1>::begin_raw_face (const unsigned int) const {
- Assert (false, ExcFunctionNotUseful());
- return 0;
-};
-
-
-template <>
-MGDoFDimensionInfo<1>::face_iterator
-MGDoFHandler<1>::begin_face (const unsigned int) const {
- Assert (false, ExcFunctionNotUseful());
- return 0;
-};
-
-
-template <>
-MGDoFDimensionInfo<1>::active_face_iterator
-MGDoFHandler<1>::begin_active_face (const unsigned int) const {
- Assert (false, ExcFunctionNotUseful());
- return 0;
-};
-
-
-template <>
-MGDoFDimensionInfo<1>::raw_face_iterator
-MGDoFHandler<1>::end_face () const {
- Assert (false, ExcFunctionNotUseful());
- return 0;
-};
-
-
-template <>
-MGDoFDimensionInfo<1>::raw_face_iterator
-MGDoFHandler<1>::last_raw_face () const {
- Assert (false, ExcFunctionNotUseful());
- return 0;
-};
-
-
-template <>
-MGDoFDimensionInfo<1>::raw_face_iterator
-MGDoFHandler<1>::last_raw_face (const unsigned int) const {
- Assert (false, ExcFunctionNotUseful());
- return 0;
-};
-
-
-template <>
-MGDoFDimensionInfo<1>::face_iterator
-MGDoFHandler<1>::last_face () const {
- Assert (false, ExcFunctionNotUseful());
- return 0;
-};
-
-
-template <>
-MGDoFDimensionInfo<1>::face_iterator
-MGDoFHandler<1>::last_face (const unsigned int) const {
- Assert (false, ExcFunctionNotUseful());
- return 0;
-};
-
-
-template <>
-MGDoFDimensionInfo<1>::active_face_iterator
-MGDoFHandler<1>::last_active_face () const {
- Assert (false, ExcFunctionNotUseful());
- return 0;
-};
-
-
-template <>
-MGDoFDimensionInfo<1>::active_face_iterator
-MGDoFHandler<1>::last_active_face (const unsigned int) const {
- Assert (false, ExcFunctionNotUseful());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::raw_quad_iterator
-MGDoFHandler<1>::begin_raw_quad (const unsigned int) const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::quad_iterator
-MGDoFHandler<1>::begin_quad (const unsigned int) const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::active_quad_iterator
-MGDoFHandler<1>::begin_active_quad (const unsigned int) const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::raw_quad_iterator
-MGDoFHandler<1>::end_quad () const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::raw_quad_iterator
-MGDoFHandler<1>::last_raw_quad (const unsigned int) const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::quad_iterator
-MGDoFHandler<1>::last_quad (const unsigned int) const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::active_quad_iterator
-MGDoFHandler<1>::last_active_quad (const unsigned int) const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::raw_quad_iterator
-MGDoFHandler<1>::last_raw_quad () const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::quad_iterator
-MGDoFHandler<1>::last_quad () const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::active_quad_iterator
-MGDoFHandler<1>::last_active_quad () const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::raw_hex_iterator
-MGDoFHandler<1>::begin_raw_hex (const unsigned int) const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::hex_iterator
-MGDoFHandler<1>::begin_hex (const unsigned int) const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::active_hex_iterator
-MGDoFHandler<1>::begin_active_hex (const unsigned int) const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::raw_hex_iterator
-MGDoFHandler<1>::end_hex () const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::raw_hex_iterator
-MGDoFHandler<1>::last_raw_hex (const unsigned int) const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::hex_iterator
-MGDoFHandler<1>::last_hex (const unsigned int) const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::active_hex_iterator
-MGDoFHandler<1>::last_active_hex (const unsigned int) const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::raw_hex_iterator
-MGDoFHandler<1>::last_raw_hex () const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::hex_iterator
-MGDoFHandler<1>::last_hex () const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<1>::active_hex_iterator
-MGDoFHandler<1>::last_active_hex () const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-#endif
-
-
-#if deal_II_dimension == 2
-
-template <>
-MGDoFHandler<2>::raw_cell_iterator
-MGDoFHandler<2>::begin_raw (const unsigned int level) const {
- return begin_raw_quad (level);
-};
-
-
-template <>
-MGDoFHandler<2>::cell_iterator
-MGDoFHandler<2>::begin (const unsigned int level) const {
- return begin_quad (level);
-};
-
-
-template <>
-MGDoFHandler<2>::active_cell_iterator
-MGDoFHandler<2>::begin_active (const unsigned int level) const {
- return begin_active_quad (level);
-};
-
-
-template <>
-MGDoFHandler<2>::raw_cell_iterator
-MGDoFHandler<2>::end () const {
- return end_quad ();
-};
-
-
-template <>
-MGDoFHandler<2>::raw_cell_iterator
-MGDoFHandler<2>::last_raw () const {
- return last_raw_quad ();
-};
-
-
-template <>
-MGDoFHandler<2>::raw_cell_iterator
-MGDoFHandler<2>::last_raw (const unsigned int level) const {
- return last_raw_quad (level);
-};
-
-
-template <>
-MGDoFHandler<2>::cell_iterator
-MGDoFHandler<2>::last () const {
- return last_quad ();
-};
-
-
-template <>
-MGDoFHandler<2>::cell_iterator
-MGDoFHandler<2>::last (const unsigned int level) const {
- return last_quad (level);
-};
-
-
-template <>
-MGDoFHandler<2>::active_cell_iterator
-MGDoFHandler<2>::last_active () const {
- return last_active_quad ();
-};
-
-
-template <>
-MGDoFHandler<2>::active_cell_iterator
-MGDoFHandler<2>::last_active (const unsigned int level) const {
- return last_active_quad (level);
-};
-
-
-template <>
-MGDoFDimensionInfo<2>::raw_face_iterator
-MGDoFHandler<2>::begin_raw_face (const unsigned int level) const {
- return begin_raw_line (level);
-};
-
-
-template <>
-MGDoFDimensionInfo<2>::face_iterator
-MGDoFHandler<2>::begin_face (const unsigned int level) const {
- return begin_line (level);
-};
-
-
-template <>
-MGDoFDimensionInfo<2>::active_face_iterator
-MGDoFHandler<2>::begin_active_face (const unsigned int level) const {
- return begin_active_line (level);
-};
-
-
-template <>
-MGDoFDimensionInfo<2>::raw_face_iterator
-MGDoFHandler<2>::end_face () const {
- return end_line ();
-};
-
-
-template <>
-MGDoFDimensionInfo<2>::raw_face_iterator
-MGDoFHandler<2>::last_raw_face () const {
- return last_raw_line ();
-};
-
-
-template <>
-MGDoFDimensionInfo<2>::raw_face_iterator
-MGDoFHandler<2>::last_raw_face (const unsigned int level) const {
- return last_raw_line (level);
-};
-
-
-template <>
-MGDoFDimensionInfo<2>::face_iterator
-MGDoFHandler<2>::last_face () const {
- return last_line ();
-};
-
-
-template <>
-MGDoFDimensionInfo<2>::face_iterator
-MGDoFHandler<2>::last_face (const unsigned int level) const {
- return last_line (level);
-};
-
-
-template <>
-MGDoFDimensionInfo<2>::active_face_iterator
-MGDoFHandler<2>::last_active_face () const {
- return last_active_line ();
-};
-
-
-template <>
-MGDoFDimensionInfo<2>::active_face_iterator
-MGDoFHandler<2>::last_active_face (const unsigned int level) const {
- return last_active_line (level);
-};
-
-
-template <>
-MGDoFHandler<2>::raw_hex_iterator
-MGDoFHandler<2>::begin_raw_hex (const unsigned int) const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<2>::hex_iterator
-MGDoFHandler<2>::begin_hex (const unsigned int) const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<2>::active_hex_iterator
-MGDoFHandler<2>::begin_active_hex (const unsigned int) const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<2>::raw_hex_iterator
-MGDoFHandler<2>::end_hex () const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<2>::raw_hex_iterator
-MGDoFHandler<2>::last_raw_hex (const unsigned int) const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<2>::hex_iterator
-MGDoFHandler<2>::last_hex (const unsigned int) const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<2>::active_hex_iterator
-MGDoFHandler<2>::last_active_hex (const unsigned int) const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<2>::raw_hex_iterator
-MGDoFHandler<2>::last_raw_hex () const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<2>::hex_iterator
-MGDoFHandler<2>::last_hex () const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-template <>
-MGDoFHandler<2>::active_hex_iterator
-MGDoFHandler<2>::last_active_hex () const {
- Assert (false, ExcNotImplemented());
- return 0;
-};
-
-
-#endif
-
-
-#if deal_II_dimension == 3
-
-template <>
-MGDoFHandler<3>::raw_cell_iterator
-MGDoFHandler<3>::begin_raw (const unsigned int level) const {
- return begin_raw_hex (level);
-};
-
-
-template <>
-MGDoFHandler<3>::cell_iterator
-MGDoFHandler<3>::begin (const unsigned int level) const {
- return begin_hex (level);
-};
-
-
-template <>
-MGDoFHandler<3>::active_cell_iterator
-MGDoFHandler<3>::begin_active (const unsigned int level) const {
- return begin_active_hex (level);
-};
-
-
-template <>
-MGDoFHandler<3>::raw_cell_iterator
-MGDoFHandler<3>::end () const {
- return end_hex ();
-};
-
-
-template <>
-MGDoFHandler<3>::raw_cell_iterator
-MGDoFHandler<3>::last_raw () const {
- return last_raw_hex ();
-};
-
-
-template <>
-MGDoFHandler<3>::raw_cell_iterator
-MGDoFHandler<3>::last_raw (const unsigned int level) const {
- return last_raw_hex (level);
-};
-
-
-template <>
-MGDoFHandler<3>::cell_iterator
-MGDoFHandler<3>::last () const {
- return last_hex ();
-};
-
-
-template <>
-MGDoFHandler<3>::cell_iterator
-MGDoFHandler<3>::last (const unsigned int level) const {
- return last_hex (level);
-};
-
-
-template <>
-MGDoFHandler<3>::active_cell_iterator
-MGDoFHandler<3>::last_active () const {
- return last_active_hex ();
-};
-
-
-template <>
-MGDoFHandler<3>::active_cell_iterator
-MGDoFHandler<3>::last_active (const unsigned int level) const {
- return last_active_hex (level);
-};
-
-
-template <>
-MGDoFDimensionInfo<3>::raw_face_iterator
-MGDoFHandler<3>::begin_raw_face (const unsigned int level) const {
- return begin_raw_quad (level);
-};
-
-
-template <>
-MGDoFDimensionInfo<3>::face_iterator
-MGDoFHandler<3>::begin_face (const unsigned int level) const {
- return begin_quad (level);
-};
-
-
-template <>
-MGDoFDimensionInfo<3>::active_face_iterator
-MGDoFHandler<3>::begin_active_face (const unsigned int level) const {
- return begin_active_quad (level);
-};
-
-
-template <>
-MGDoFDimensionInfo<3>::raw_face_iterator
-MGDoFHandler<3>::end_face () const {
- return end_quad ();
-};
-
-
-template <>
-MGDoFDimensionInfo<3>::raw_face_iterator
-MGDoFHandler<3>::last_raw_face () const {
- return last_raw_quad ();
-};
-
-
-template <>
-MGDoFDimensionInfo<3>::raw_face_iterator
-MGDoFHandler<3>::last_raw_face (const unsigned int level) const {
- return last_raw_quad (level);
-};
-
-
-template <>
-MGDoFDimensionInfo<3>::face_iterator
-MGDoFHandler<3>::last_face () const {
- return last_quad ();
-};
-
-
-template <>
-MGDoFDimensionInfo<3>::face_iterator
-MGDoFHandler<3>::last_face (const unsigned int level) const {
- return last_quad (level);
-};
-
-
-template <>
-MGDoFDimensionInfo<3>::active_face_iterator
-MGDoFHandler<3>::last_active_face () const {
- return last_active_quad ();
-};
-
-
-template <>
-MGDoFDimensionInfo<3>::active_face_iterator
-MGDoFHandler<3>::last_active_face (const unsigned int level) const {
- return last_active_quad (level);
-};
-
-#endif
-
-
-template <int dim>
-typename MGDoFHandler<dim>::raw_line_iterator
-MGDoFHandler<dim>::begin_raw_line (const unsigned int level) const {
- return raw_line_iterator (tria,
- tria->begin_raw_line(level)->level(),
- tria->begin_raw_line(level)->index(),
- this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::line_iterator
-MGDoFHandler<dim>::begin_line (const unsigned int level) const {
- return line_iterator (tria,
- tria->begin_line(level)->level(),
- tria->begin_line(level)->index(),
- this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::active_line_iterator
-MGDoFHandler<dim>::begin_active_line (const unsigned int level) const {
- return active_line_iterator (tria,
- tria->begin_active_line(level)->level(),
- tria->begin_active_line(level)->index(),
- this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::raw_quad_iterator
-MGDoFHandler<dim>::begin_raw_quad (const unsigned int level) const {
- return raw_quad_iterator (tria,
- tria->begin_raw_quad(level)->level(),
- tria->begin_raw_quad(level)->index(),
- this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::quad_iterator
-MGDoFHandler<dim>::begin_quad (const unsigned int level) const {
- return quad_iterator (tria,
- tria->begin_quad(level)->level(),
- tria->begin_quad(level)->index(),
- this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::active_quad_iterator
-MGDoFHandler<dim>::begin_active_quad (const unsigned int level) const {
- return active_quad_iterator (tria,
- tria->begin_active_quad(level)->level(),
- tria->begin_active_quad(level)->index(),
- this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::raw_hex_iterator
-MGDoFHandler<dim>::begin_raw_hex (const unsigned int level) const {
- return raw_hex_iterator (tria,
- tria->begin_raw_hex(level)->level(),
- tria->begin_raw_hex(level)->index(),
- this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::hex_iterator
-MGDoFHandler<dim>::begin_hex (const unsigned int level) const {
- return hex_iterator (tria,
- tria->begin_hex(level)->level(),
- tria->begin_hex(level)->index(),
- this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::active_hex_iterator
-MGDoFHandler<dim>::begin_active_hex (const unsigned int level) const {
- return active_hex_iterator (tria,
- tria->begin_active_hex(level)->level(),
- tria->begin_active_hex(level)->index(),
- this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::raw_line_iterator
-MGDoFHandler<dim>::end_line () const {
- return raw_line_iterator (tria, -1, -1, this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::raw_quad_iterator
-MGDoFHandler<dim>::end_quad () const {
- return raw_quad_iterator (tria, -1, -1, this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::raw_hex_iterator
-MGDoFHandler<dim>::end_hex () const {
- return raw_hex_iterator (tria, -1, -1, this);
-};
-
-
-template <int dim>
-typename MGDoFDimensionInfo<dim>::raw_cell_iterator
-MGDoFHandler<dim>::end_raw (const unsigned int level) const {
- return (level == mg_levels.size()-1 ?
- end() :
- begin_raw (level+1));
-};
-
-
-template <int dim>
-typename MGDoFDimensionInfo<dim>::cell_iterator
-MGDoFHandler<dim>::end (const unsigned int level) const {
- return (level == mg_levels.size()-1 ?
- cell_iterator(end()) :
- begin (level+1));
-};
-
-
-template <int dim>
-typename MGDoFDimensionInfo<dim>::active_cell_iterator
-MGDoFHandler<dim>::end_active (const unsigned int level) const {
- return (level == mg_levels.size()-1 ?
- active_cell_iterator(end()) :
- begin_active (level+1));
-};
-
-
-template <int dim>
-typename MGDoFDimensionInfo<dim>::raw_face_iterator
-MGDoFHandler<dim>::end_raw_face (const unsigned int level) const {
- return (level == mg_levels.size()-1 ?
- end_face() :
- begin_raw_face (level+1));
-};
-
-
-template <int dim>
-typename MGDoFDimensionInfo<dim>::face_iterator
-MGDoFHandler<dim>::end_face (const unsigned int level) const {
- return (level == mg_levels.size()-1 ?
- face_iterator(end_face()) :
- begin_face (level+1));
-};
-
-
-template <int dim>
-typename MGDoFDimensionInfo<dim>::active_face_iterator
-MGDoFHandler<dim>::end_active_face (const unsigned int level) const {
- return (level == mg_levels.size()-1 ?
- active_face_iterator(end_face()) :
- begin_active_face (level+1));
-};
-
-
-template <int dim>
-typename MGDoFDimensionInfo<dim>::raw_line_iterator
-MGDoFHandler<dim>::end_raw_line (const unsigned int level) const {
- return (level == mg_levels.size()-1 ?
- end_line() :
- begin_raw_line (level+1));
-};
-
-
-template <int dim>
-typename MGDoFDimensionInfo<dim>::line_iterator
-MGDoFHandler<dim>::end_line (const unsigned int level) const {
- return (level == mg_levels.size()-1 ?
- line_iterator(end_line()) :
- begin_line (level+1));
-};
-
-
-template <int dim>
-typename MGDoFDimensionInfo<dim>::active_line_iterator
-MGDoFHandler<dim>::end_active_line (const unsigned int level) const {
- return (level == mg_levels.size()-1 ?
- active_line_iterator(end_line()) :
- begin_active_line (level+1));
-};
-
-
-template <int dim>
-typename MGDoFDimensionInfo<dim>::raw_quad_iterator
-MGDoFHandler<dim>::end_raw_quad (const unsigned int level) const {
- return (level == mg_levels.size()-1 ?
- end_quad() :
- begin_raw_quad (level+1));
-};
-
-
-template <int dim>
-typename MGDoFDimensionInfo<dim>::quad_iterator
-MGDoFHandler<dim>::end_quad (const unsigned int level) const {
- return (level == mg_levels.size()-1 ?
- quad_iterator(end_quad()) :
- begin_quad (level+1));
-};
-
-
-template <int dim>
-typename MGDoFDimensionInfo<dim>::active_quad_iterator
-MGDoFHandler<dim>::end_active_quad (const unsigned int level) const {
- return (level == mg_levels.size()-1 ?
- active_quad_iterator(end_quad()) :
- begin_active_quad (level+1));
-};
-
-
-template <int dim>
-typename MGDoFDimensionInfo<dim>::raw_hex_iterator
-MGDoFHandler<dim>::end_raw_hex (const unsigned int level) const {
- return (level == mg_levels.size()-1 ?
- end_hex() :
- begin_raw_hex (level+1));
-};
-
-
-template <int dim>
-typename MGDoFDimensionInfo<dim>::hex_iterator
-MGDoFHandler<dim>::end_hex (const unsigned int level) const {
- return (level == mg_levels.size()-1 ?
- hex_iterator(end_hex()) :
- begin_hex (level+1));
-};
-
-
-template <int dim>
-typename MGDoFDimensionInfo<dim>::active_hex_iterator
-MGDoFHandler<dim>::end_active_hex (const unsigned int level) const {
- return (level == mg_levels.size()-1 ?
- active_hex_iterator(end_hex()) :
- begin_active_hex (level+1));
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::raw_line_iterator
-MGDoFHandler<dim>::last_raw_line (const unsigned int level) const {
- return raw_line_iterator (tria,
- tria->last_raw_line(level)->level(),
- tria->last_raw_line(level)->index(),
- this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::line_iterator
-MGDoFHandler<dim>::last_line (const unsigned int level) const {
- return line_iterator (tria,
- tria->last_line(level)->level(),
- tria->last_line(level)->index(),
- this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::active_line_iterator
-MGDoFHandler<dim>::last_active_line (const unsigned int level) const {
- return active_line_iterator (tria,
- tria->last_active_line(level)->level(),
- tria->last_active_line(level)->index(),
- this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::raw_quad_iterator
-MGDoFHandler<dim>::last_raw_quad (const unsigned int level) const {
- return raw_quad_iterator (tria,
- tria->last_raw_quad(level)->level(),
- tria->last_raw_quad(level)->index(),
- this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::quad_iterator
-MGDoFHandler<dim>::last_quad (const unsigned int level) const {
- return quad_iterator (tria,
- tria->last_quad(level)->level(),
- tria->last_quad(level)->index(),
- this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::active_quad_iterator
-MGDoFHandler<dim>::last_active_quad (const unsigned int level) const {
- return active_quad_iterator (tria,
- tria->last_active_quad(level)->level(),
- tria->last_active_quad(level)->index(),
- this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::raw_hex_iterator
-MGDoFHandler<dim>::last_raw_hex (const unsigned int level) const {
- return raw_hex_iterator (tria,
- tria->last_raw_hex(level)->level(),
- tria->last_raw_hex(level)->index(),
- this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::hex_iterator
-MGDoFHandler<dim>::last_hex (const unsigned int level) const {
- return hex_iterator (tria,
- tria->last_hex(level)->level(),
- tria->last_hex(level)->index(),
- this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::active_hex_iterator
-MGDoFHandler<dim>::last_active_hex (const unsigned int level) const {
- return active_hex_iterator (tria,
- tria->last_active_hex(level)->level(),
- tria->last_active_hex(level)->index(),
- this);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::raw_line_iterator
-MGDoFHandler<dim>::last_raw_line () const {
- return last_raw_line (mg_levels.size()-1);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::raw_quad_iterator
-MGDoFHandler<dim>::last_raw_quad () const {
- return last_raw_quad (mg_levels.size()-1);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::raw_hex_iterator
-MGDoFHandler<dim>::last_raw_hex () const {
- return last_raw_hex (mg_levels.size()-1);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::line_iterator
-MGDoFHandler<dim>::last_line () const {
- return last_line (mg_levels.size()-1);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::quad_iterator
-MGDoFHandler<dim>::last_quad () const {
- return last_quad (mg_levels.size()-1);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::hex_iterator
-MGDoFHandler<dim>::last_hex () const {
- return last_hex (mg_levels.size()-1);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::active_line_iterator
-MGDoFHandler<dim>::last_active_line () const {
- return last_active_line (mg_levels.size()-1);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::active_quad_iterator
-MGDoFHandler<dim>::last_active_quad () const {
- return last_active_quad (mg_levels.size()-1);
-};
-
-
-template <int dim>
-typename MGDoFHandler<dim>::active_hex_iterator
-MGDoFHandler<dim>::last_active_hex () const {
- return last_active_hex (mg_levels.size()-1);
-};
-
-
-//------------------------------------------------------------------
-
-
-template <int dim>
-void MGDoFHandler<dim>::distribute_dofs (const FiniteElement<dim> &fe,
- unsigned int offset)
-{
- // first distribute global dofs
- DoFHandler<dim>::distribute_dofs (fe, offset);
-
-
-// reserve space for the MG dof numbers
- reserve_space ();
- mg_used_dofs.resize (tria->n_levels(), 0);
-
- // clear user flags because we will
- // need them
- tria->clear_user_flags ();
-
- // now distribute indices on each level
- // separately
- for (unsigned int level=0; level<tria->n_levels(); ++level)
- {
- unsigned int next_free_dof = offset;
- cell_iterator cell = begin(level),
- endc = end(level);
-
- for (; cell != endc; ++cell)
- next_free_dof = distribute_dofs_on_cell (cell, next_free_dof);
-
- mg_used_dofs[level] = next_free_dof;
- };
-};
-
-
-#if deal_II_dimension == 1
-
-template <>
-unsigned int
-MGDoFHandler<1>::distribute_dofs_on_cell (cell_iterator &cell,
- unsigned int next_free_dof) {
-
- // distribute dofs of vertices
- if (selected_fe->dofs_per_vertex > 0)
- for (unsigned int v=0; v<GeometryInfo<1>::vertices_per_cell; ++v)
- {
- cell_iterator neighbor = cell->neighbor(v);
-
- if (neighbor.state() == valid)
- {
- // has neighbor already been processed?
- if (neighbor->user_flag_set() &&
- (neighbor->level() == cell->level()))
- // copy dofs if the neighbor is on
- // the same level (only then are
- // mg dofs the same)
- {
- if (v==0)
- for (unsigned int d=0; d<selected_fe->dofs_per_vertex; ++d)
- cell->set_mg_vertex_dof_index (0, d,
- neighbor->mg_vertex_dof_index (1, d));
- else
- for (unsigned int d=0; d<selected_fe->dofs_per_vertex; ++d)
- cell->set_mg_vertex_dof_index (1, d,
- neighbor->mg_vertex_dof_index (0, d));
-
- // next neighbor
- continue;
- };
- };
-
- // otherwise: create dofs newly
- for (unsigned int d=0; d<selected_fe->dofs_per_vertex; ++d)
- cell->set_mg_vertex_dof_index (v, d, next_free_dof++);
- };
-
- // dofs of line
- if (selected_fe->dofs_per_line > 0)
- for (unsigned int d=0; d<selected_fe->dofs_per_line; ++d)
- cell->set_mg_dof_index (d, next_free_dof++);
-
- // note that this cell has been processed
- cell->set_user_flag ();
-
- return next_free_dof;
-};
-
-#endif
-
-
-#if deal_II_dimension == 2
-
-template <>
-unsigned int
-MGDoFHandler<2>::distribute_dofs_on_cell (cell_iterator &cell,
- unsigned int next_free_dof) {
- if (selected_fe->dofs_per_vertex > 0)
- // number dofs on vertices
- for (unsigned int vertex=0; vertex<GeometryInfo<2>::vertices_per_cell; ++vertex)
- // check whether dofs for this
- // vertex have been distributed
- // (only check the first dof)
- if (cell->mg_vertex_dof_index(vertex, 0) == DoFHandler<2>::invalid_dof_index)
- for (unsigned int d=0; d<selected_fe->dofs_per_vertex; ++d)
- cell->set_mg_vertex_dof_index (vertex, d, next_free_dof++);
-
- // for the four sides
- if (selected_fe->dofs_per_line > 0)
- for (unsigned int side=0; side<GeometryInfo<2>::faces_per_cell; ++side)
- {
- line_iterator line = cell->line(side);
-
- // distribute dofs if necessary:
- // check whether line dof is already
- // numbered (check only first dof)
- if (line->mg_dof_index(0) == DoFHandler<2>::invalid_dof_index)
- // if not: distribute dofs
- for (unsigned int d=0; d<selected_fe->dofs_per_line; ++d)
- line->set_mg_dof_index (d, next_free_dof++);
- };
-
-
-// dofs of quad
- if (selected_fe->dofs_per_quad > 0)
- for (unsigned int d=0; d<selected_fe->dofs_per_quad; ++d)
- cell->set_mg_dof_index (d, next_free_dof++);
-
-
-// note that this cell has been processed
- cell->set_user_flag ();
-
- return next_free_dof;
-};
-
-#endif
-
-
-#if deal_II_dimension == 3
-
-template <>
-unsigned int
-MGDoFHandler<3>::distribute_dofs_on_cell (cell_iterator &cell,
- unsigned int next_free_dof) {
- if (selected_fe->dofs_per_vertex > 0)
- // number dofs on vertices
- for (unsigned int vertex=0; vertex<GeometryInfo<3>::vertices_per_cell; ++vertex)
- // check whether dofs for this
- // vertex have been distributed
- // (only check the first dof)
- if (cell->mg_vertex_dof_index(vertex, 0) == DoFHandler<3>::invalid_dof_index)
- for (unsigned int d=0; d<selected_fe->dofs_per_vertex; ++d)
- cell->set_mg_vertex_dof_index (vertex, d, next_free_dof++);
-
- // for the lines
- if (selected_fe->dofs_per_line > 0)
- for (unsigned int l=0; l<GeometryInfo<3>::lines_per_cell; ++l)
- {
- line_iterator line = cell->line(l);
-
- // distribute dofs if necessary:
- // check whether line dof is already
- // numbered (check only first dof)
- if (line->mg_dof_index(0) == DoFHandler<3>::invalid_dof_index)
- // if not: distribute dofs
- for (unsigned int d=0; d<selected_fe->dofs_per_line; ++d)
- line->set_mg_dof_index (d, next_free_dof++);
- };
-
- // for the quads
- if (selected_fe->dofs_per_quad > 0)
- for (unsigned int q=0; q<GeometryInfo<3>::lines_per_cell; ++q)
- {
- quad_iterator quad = cell->quad(q);
-
- // distribute dofs if necessary:
- // check whether line dof is already
- // numbered (check only first dof)
- if (quad->mg_dof_index(0) == DoFHandler<3>::invalid_dof_index)
- // if not: distribute dofs
- for (unsigned int d=0; d<selected_fe->dofs_per_line; ++d)
- quad->set_mg_dof_index (d, next_free_dof++);
- };
-
-
-// dofs of cell
- if (selected_fe->dofs_per_hex > 0)
- for (unsigned int d=0; d<selected_fe->dofs_per_hex; ++d)
- cell->set_mg_dof_index (d, next_free_dof++);
-
-
-// note that this cell has been processed
- cell->set_user_flag ();
-
- return next_free_dof;
-};
-
-#endif
-
-
-template <int dim>
-unsigned int MGDoFHandler<dim>::n_dofs (const unsigned int level) const {
- Assert (level < mg_used_dofs.size(), ExcInvalidLevel(level));
-
- return mg_used_dofs[level];
-};
-
-
-#if deal_II_dimension == 1
-
-template <>
-void MGDoFHandler<1>::renumber_dofs (const unsigned int level,
- const vector<unsigned int> &new_numbers) {
- Assert (new_numbers.size() == n_dofs(level), ExcRenumberingIncomplete());
-
- // note that we can not use cell iterators
- // in this function since then we would
- // renumber the dofs on the interface of
- // two cells more than once. Anyway, this
- // ways it's not only more correct but also
- // faster
- for (vector<MGVertexDoFs>::iterator i=mg_vertex_dofs.begin();
- i!=mg_vertex_dofs.end(); ++i)
- // if the present vertex lives on
- // the present level
- if ((i->get_coarsest_level() <= level) &&
- (i->get_finest_level() >= level))
- for (unsigned int d=0; d<selected_fe->dofs_per_vertex; ++d)
- i->set_index (level, d, selected_fe->dofs_per_vertex,
- new_numbers[i->get_index (level, d,
- selected_fe->dofs_per_vertex)]);
-
- for (vector<unsigned int>::iterator i=mg_levels[level]->line_dofs.begin();
- i!=mg_levels[level]->line_dofs.end(); ++i)
- {
- Assert (*i != DoFHandler<1>::invalid_dof_index, ExcInternalError());
- *i = new_numbers[*i];
- };
-};
-
-#endif
-
-
-#if deal_II_dimension == 2
-
-template <>
-void MGDoFHandler<2>::renumber_dofs (const unsigned int level,
- const vector<unsigned int> &new_numbers) {
- Assert (new_numbers.size() == n_dofs(level), ExcRenumberingIncomplete());
-
- for (vector<MGVertexDoFs>::iterator i=mg_vertex_dofs.begin();
- i!=mg_vertex_dofs.end(); ++i)
- // if the present vertex lives on
- // the present level
- if ((i->get_coarsest_level() <= level) &&
- (i->get_finest_level() >= level))
- for (unsigned int d=0; d<selected_fe->dofs_per_vertex; ++d)
- i->set_index (level, d, selected_fe->dofs_per_vertex,
- new_numbers[i->get_index (level, d,
- selected_fe->dofs_per_vertex)]);
-
- for (vector<unsigned int>::iterator i=mg_levels[level]->line_dofs.begin();
- i!=mg_levels[level]->line_dofs.end(); ++i)
- {
- Assert (*i != DoFHandler<2>::invalid_dof_index, ExcInternalError());
- *i = new_numbers[*i];
- };
-
- for (vector<unsigned int>::iterator i=mg_levels[level]->quad_dofs.begin();
- i!=mg_levels[level]->quad_dofs.end(); ++i)
- {
- Assert (*i != DoFHandler<2>::invalid_dof_index, ExcInternalError());
- *i = new_numbers[*i];
- };
-};
-
-#endif
-
-
-#if deal_II_dimension == 3
-
-template <>
-void MGDoFHandler<3>::renumber_dofs (const unsigned int level,
- const vector<unsigned int> &new_numbers) {
- Assert (new_numbers.size() == n_dofs(level), ExcRenumberingIncomplete());
-
- for (vector<MGVertexDoFs>::iterator i=mg_vertex_dofs.begin();
- i!=mg_vertex_dofs.end(); ++i)
- // if the present vertex lives on
- // the present level
- if ((i->get_coarsest_level() <= level) &&
- (i->get_finest_level() >= level))
- for (unsigned int d=0; d<selected_fe->dofs_per_vertex; ++d)
- i->set_index (level, d, selected_fe->dofs_per_vertex,
- new_numbers[i->get_index (level, d,
- selected_fe->dofs_per_vertex)]);
-
- for (vector<unsigned int>::iterator i=mg_levels[level]->line_dofs.begin();
- i!=mg_levels[level]->line_dofs.end(); ++i)
- {
- Assert (*i != DoFHandler<3>::invalid_dof_index, ExcInternalError());
- *i = new_numbers[*i];
- };
-
- for (vector<unsigned int>::iterator i=mg_levels[level]->quad_dofs.begin();
- i!=mg_levels[level]->quad_dofs.end(); ++i)
- {
- Assert (*i != DoFHandler<3>::invalid_dof_index, ExcInternalError());
- *i = new_numbers[*i];
- };
-
- for (vector<unsigned int>::iterator i=mg_levels[level]->hex_dofs.begin();
- i!=mg_levels[level]->hex_dofs.end(); ++i)
- {
- Assert (*i != DoFHandler<3>::invalid_dof_index, ExcInternalError());
- *i = new_numbers[*i];
- };
-};
-
-#endif
-
-
-#if deal_II_dimension == 1
-
-template <>
-void MGDoFHandler<1>::reserve_space () {
- const unsigned int dim = 1;
-
- Assert (selected_fe != 0, ExcNoFESelected());
- Assert (tria->n_levels() > 0, ExcInvalidTriangulation());
-
- //////////////////////////
- // DESTRUCTION
-
- // delete all levels and set them up
- // newly, since vectors are
- // troublesome if you want to change
- // their size
- for (unsigned int i=0; i<mg_levels.size(); ++i)
- delete mg_levels[i];
- mg_levels.resize (0);
-
- // also delete vector of vertex indices
- // this calls the destructor which
- // must free the space
- mg_vertex_dofs.resize (0);
-
- ////////////////////////////
- // CONSTRUCTION
-
- // first allocate space for the
- // lines on each level
- for (unsigned int i=0; i<tria->n_levels(); ++i)
- {
- mg_levels.push_back (new DoFLevel<1>);
-
- mg_levels.back()->line_dofs = vector<unsigned int>(tria->levels[i]->lines.lines.size() *
- selected_fe->dofs_per_line,
- DoFHandler<1>::invalid_dof_index);
- };
-
- // now allocate space for the
- // vertices. To this end, we need
- // to construct as many objects as
- // there are vertices and let them
- // allocate enough space for their
- // vertex indices on the levels they
- // live on. We need therefore to
- // count to how many levels a cell
- // belongs to, which we do by looping
- // over all cells and storing the
- // maximum and minimum level each
- // vertex we pass by belongs to
- mg_vertex_dofs.resize (tria->vertices.size());
-
- vector<unsigned int> min_level (tria->vertices.size(), tria->n_levels());
- vector<unsigned int> max_level (tria->vertices.size(), 0);
-
- Triangulation<dim>::cell_iterator cell = tria->begin(),
- endc = tria->end();
- for (; cell!=endc; ++cell)
- for (unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell;
- ++vertex)
- {
- const unsigned int vertex_index = cell->vertex_index(vertex);
- if (min_level[vertex_index] > static_cast<unsigned int>(cell->level()))
- min_level[vertex_index] = cell->level();
- if (max_level[vertex_index] < static_cast<unsigned int>(cell->level()))
- max_level[vertex_index] = cell->level();
- };
-
- // now allocate the needed space
- for (unsigned int vertex=0; vertex<tria->vertices.size(); ++vertex)
- {
- Assert (min_level[vertex] < tria->n_levels(), ExcInternalError());
- Assert (max_level[vertex] >= min_level[vertex], ExcInternalError());
-
- mg_vertex_dofs[vertex].init (min_level[vertex],
- max_level[vertex],
- selected_fe->dofs_per_vertex);
- };
-};
-
-#endif
-
-
-#if deal_II_dimension == 2
-
-template <>
-void MGDoFHandler<2>::reserve_space () {
- const unsigned int dim = 2;
-
- Assert (selected_fe != 0, ExcNoFESelected());
- Assert (tria->n_levels() > 0, ExcInvalidTriangulation());
-
- ////////////////////////////
- // DESTRUCTION
-
- // delete all levels and set them up
- // newly, since vectors are
- // troublesome if you want to change
- // their size
- for (unsigned int i=0; i<mg_levels.size(); ++i)
- delete mg_levels[i];
- mg_levels.clear ();
-
- // also delete vector of vertex indices
- // this calls the destructor which
- // must free the space
- mg_vertex_dofs.clear ();
-
-
-////////////////////////////
- // CONSTRUCTION
-
- // first allocate space for the
- // lines and quads on each level
- for (unsigned int i=0; i<tria->n_levels(); ++i)
- {
- mg_levels.push_back (new DoFLevel<2>);
-
- mg_levels.back()->line_dofs = vector<unsigned int> (tria->levels[i]->lines.lines.size() *
- selected_fe->dofs_per_line,
- DoFHandler<2>::invalid_dof_index);
- mg_levels.back()->quad_dofs = vector<unsigned int> (tria->levels[i]->quads.quads.size() *
- selected_fe->dofs_per_quad,
- DoFHandler<2>::invalid_dof_index);
- };
-
-
-// now allocate space for the
- // vertices. To this end, we need
- // to construct as many objects as
- // there are vertices and let them
- // allocate enough space for their
- // vertex indices on the levels they
- // live on. We need therefore to
- // count to how many levels a cell
- // belongs to, which we do by looping
- // over all cells and storing the
- // maximum and minimum level each
- // vertex we pass by belongs to
- mg_vertex_dofs.resize (tria->vertices.size());
-
- vector<unsigned int> min_level (tria->vertices.size(), tria->n_levels());
- vector<unsigned int> max_level (tria->vertices.size(), 0);
-
- Triangulation<dim>::cell_iterator cell = tria->begin(),
- endc = tria->end();
- for (; cell!=endc; ++cell)
- for (unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell;
- ++vertex)
- {
- const unsigned int vertex_index = cell->vertex_index(vertex);
- if (min_level[vertex_index] > static_cast<unsigned int>(cell->level()))
- min_level[vertex_index] = cell->level();
- if (max_level[vertex_index] < static_cast<unsigned int>(cell->level()))
- max_level[vertex_index] = cell->level();
- };
-
-
-// now allocate the needed space
- for (unsigned int vertex=0; vertex<tria->vertices.size(); ++vertex)
- {
- Assert (min_level[vertex] < tria->n_levels(), ExcInternalError());
- Assert (max_level[vertex] >= min_level[vertex], ExcInternalError());
-
- mg_vertex_dofs[vertex].init (min_level[vertex],
- max_level[vertex],
- selected_fe->dofs_per_vertex);
- };
-};
-
-#endif
-
-
-#if deal_II_dimension == 3
-
-template <>
-void MGDoFHandler<3>::reserve_space () {
- const unsigned int dim = 3;
-
- Assert (selected_fe != 0, ExcNoFESelected());
- Assert (tria->n_levels() > 0, ExcInvalidTriangulation());
-
- ////////////////////////////
- // DESTRUCTION
-
- // delete all levels and set them up
- // newly, since vectors are
- // troublesome if you want to change
- // their size
- for (unsigned int i=0; i<mg_levels.size(); ++i)
- delete mg_levels[i];
- mg_levels.clear ();
-
- // also delete vector of vertex indices
- // this calls the destructor which
- // must free the space
- mg_vertex_dofs.clear ();
-
-
-////////////////////////////
- // CONSTRUCTION
-
- // first allocate space for the
- // lines and quads on each level
- for (unsigned int i=0; i<tria->n_levels(); ++i)
- {
- mg_levels.push_back (new DoFLevel<3>);
-
- mg_levels.back()->line_dofs = vector<unsigned int> (tria->levels[i]->lines.lines.size() *
- selected_fe->dofs_per_line,
- DoFHandler<3>::invalid_dof_index);
- mg_levels.back()->quad_dofs = vector<unsigned int> (tria->levels[i]->quads.quads.size() *
- selected_fe->dofs_per_quad,
- DoFHandler<3>::invalid_dof_index);
- mg_levels.back()->hex_dofs = vector<unsigned int> (tria->levels[i]->hexes.hexes.size() *
- selected_fe->dofs_per_hex,
- DoFHandler<3>::invalid_dof_index);
- };
-
-
-// now allocate space for the
- // vertices. To this end, we need
- // to construct as many objects as
- // there are vertices and let them
- // allocate enough space for their
- // vertex indices on the levels they
- // live on. We need therefore to
- // count to how many levels a cell
- // belongs to, which we do by looping
- // over all cells and storing the
- // maximum and minimum level each
- // vertex we pass by belongs to
- mg_vertex_dofs.resize (tria->vertices.size());
-
- vector<unsigned int> min_level (tria->vertices.size(), tria->n_levels());
- vector<unsigned int> max_level (tria->vertices.size(), 0);
-
- Triangulation<dim>::cell_iterator cell = tria->begin(),
- endc = tria->end();
- for (; cell!=endc; ++cell)
- for (unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell;
- ++vertex)
- {
- const unsigned int vertex_index = cell->vertex_index(vertex);
- if (min_level[vertex_index] > static_cast<unsigned int>(cell->level()))
- min_level[vertex_index] = cell->level();
- if (max_level[vertex_index] < static_cast<unsigned int>(cell->level()))
- max_level[vertex_index] = cell->level();
- };
-
-
-// now allocate the needed space
- for (unsigned int vertex=0; vertex<tria->vertices.size(); ++vertex)
- {
- Assert (min_level[vertex] < tria->n_levels(), ExcInternalError());
- Assert (max_level[vertex] >= min_level[vertex], ExcInternalError());
-
- mg_vertex_dofs[vertex].init (min_level[vertex],
- max_level[vertex],
- selected_fe->dofs_per_vertex);
- };
-};
-
-#endif
-
-
-// explicite instantiations
-template class MGDoFHandler<deal_II_dimension>;
+++ /dev/null
-//---------------------------- mg_dof_tools.cc ---------------------------
-// $Id$
-// Version: $Name$
-//
-// Copyright (C) 1998, 1999, 2000 by the deal.II authors
-//
-// This file is subject to QPL and may not be distributed
-// without copyright and license information. Please refer
-// to the file deal.II/doc/license.html for the text and
-// further information on this license.
-//
-//---------------------------- mg_dof_tools.cc ---------------------------
-
-
-#include <lac/sparsity_pattern.h>
-#include <multigrid/mg_dof_handler.h>
-#include <multigrid/mg_dof_accessor.h>
-#include <grid/tria_iterator.h>
-#include <multigrid/mg_dof_tools.h>
-#include <fe/fe.h>
-
-
-template <int dim>
-void MGDoFTools::make_sparsity_pattern (const MGDoFHandler<dim> &dof,
- SparsityPattern &sparsity,
- const unsigned int level)
-{
- const unsigned int n_dofs = dof.n_dofs(level);
-
- Assert (sparsity.n_rows() == n_dofs,
- ExcDimensionMismatch (sparsity.n_rows(), n_dofs));
- Assert (sparsity.n_cols() == n_dofs,
- ExcDimensionMismatch (sparsity.n_cols(), n_dofs));
-
- const unsigned int dofs_per_cell = dof.get_fe().dofs_per_cell;
- vector<unsigned int> dofs_on_this_cell(dofs_per_cell);
- MGDoFHandler<dim>::cell_iterator cell = dof.begin(level),
- endc = dof.end(level);
- for (; cell!=endc; ++cell)
- {
- cell->get_mg_dof_indices (dofs_on_this_cell);
- // make sparsity pattern for this cell
- for (unsigned int i=0; i<dofs_per_cell; ++i)
- for (unsigned int j=0; j<dofs_per_cell; ++j)
- sparsity.add (dofs_on_this_cell[i],
- dofs_on_this_cell[j]);
- }
-}
-
-template<int dim>
-void
-MGDoFTools::make_flux_sparsity_pattern (const MGDoFHandler<dim> &dof,
- SparsityPattern &sparsity,
- const unsigned int level)
-{
- const unsigned int n_dofs = dof.n_dofs(level);
-
- Assert (sparsity.n_rows() == n_dofs,
- ExcDimensionMismatch (sparsity.n_rows(), n_dofs));
- Assert (sparsity.n_cols() == n_dofs,
- ExcDimensionMismatch (sparsity.n_cols(), n_dofs));
-
- const unsigned int dofs_per_cell = dof.get_fe().dofs_per_cell;
- vector<unsigned int> dofs_on_this_cell(dofs_per_cell);
- vector<unsigned int> dofs_on_other_cell(dofs_per_cell);
- MGDoFHandler<dim>::cell_iterator cell = dof.begin(level),
- endc = dof.end(level);
- for (; cell!=endc; ++cell)
- {
- cell->get_mg_dof_indices (dofs_on_this_cell);
- // make sparsity pattern for this cell
- for (unsigned int i=0; i<dofs_per_cell; ++i)
- for (unsigned int j=0; j<dofs_per_cell; ++j)
- sparsity.add (dofs_on_this_cell[i],
- dofs_on_this_cell[j]);
-
- // Loop over all interior neighbors
- for (unsigned int face = 0;
- face < GeometryInfo<dim>::faces_per_cell;
- ++face)
- {
- if ( (! cell->at_boundary(face)) &&
- (static_cast<unsigned int>(cell->neighbor_level(face)) == level) )
- {
- MGDoFHandler<dim>::cell_iterator neighbor = cell->neighbor(face);
- neighbor->get_mg_dof_indices (dofs_on_other_cell);
- // only add one direction
- // The other is taken care of
- // by neighbor.
- for (unsigned int i=0; i<dofs_per_cell; ++i)
- {
- for (unsigned int j=0; j<dofs_per_cell; ++j)
- {
- sparsity.add (dofs_on_this_cell[i],
- dofs_on_other_cell[j]);
- }
- }
- }
- }
- }
-}
-
-
-// explicit instantiations
-template void MGDoFTools::make_sparsity_pattern (const MGDoFHandler<deal_II_dimension> &,
- SparsityPattern &,
- const unsigned int);
-template void MGDoFTools::make_flux_sparsity_pattern (const MGDoFHandler<deal_II_dimension> &,
- SparsityPattern &,
- const unsigned int);
-
+++ /dev/null
-//---------------------------- mg_smoother.cc ---------------------------
-// $Id$
-// Version: $Name$
-//
-// Copyright (C) 1998, 1999, 2000 by the deal.II authors
-//
-// This file is subject to QPL and may not be distributed
-// without copyright and license information. Please refer
-// to the file deal.II/doc/license.html for the text and
-// further information on this license.
-//
-//---------------------------- mg_smoother.cc ---------------------------
-
-
-#include <grid/tria.h>
-#include <multigrid/mg_dof_handler.h>
-#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>
-
-
-//////////////////////////////////////////////////////////////////////
-
-
-MGSmootherBase::~MGSmootherBase()
-{};
-
-
-//////////////////////////////////////////////////////////////////////
-
-
-void
-MGSmootherIdentity::smooth (const unsigned int,
- Vector<double> &,
- const Vector<double> &) const
-{};
-
-
-#if deal_II_dimension == 1
-
-MGSmoother::MGSmoother (const MGDoFHandler<1> &/*mg_dof*/, unsigned int steps)
- :
- steps(steps)
-{
- Assert (false, ExcNotImplemented());
-};
-
-#endif
-
-
-#if deal_II_dimension > 1
-
-template <int dim>
-MGSmoother::MGSmoother (const MGDoFHandler<dim> &mg_dof, unsigned int steps)
- :
- steps(steps)
-{
- const unsigned int n_levels = mg_dof.get_tria().n_levels();
-
- // allocate the right number of
- // elements
- interior_boundary_dofs.resize (n_levels-1);
-
- // use a temporary to store the
- // indices. this allows to not allocate
- // to much memory by later copying the
- // content of this vector to its final
- // destination
- vector<unsigned int> boundary_dofs;
-
- // temporary to hold the dof indices
- // on a face between to levels
- vector<unsigned int> dofs_on_face (mg_dof.get_fe().dofs_per_face);
-
- for (unsigned int level=1; level<n_levels; ++level)
- {
- boundary_dofs.clear ();
-
- // for each cell on this level:
- // find out whether a face is
- // at the boundary of this level's
- // cells and if so add the dofs
- // to the interior boundary dofs
- for (MGDoFHandler<dim>::cell_iterator cell=mg_dof.begin(level);
- cell != mg_dof.end(level); ++cell)
- for (unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
- if ((cell->neighbor(face).state() == valid) &&
- (static_cast<unsigned int>(cell->neighbor(face)->level())
- == level-1))
- {
- // get indices of this face
- cell->face(face)->get_mg_dof_indices (dofs_on_face);
- // append them to the levelwise
- // list
- boundary_dofs.insert (boundary_dofs.end(),
- dofs_on_face.begin(),
- dofs_on_face.end());
- };
-
- // now sort the list of interior boundary
- // dofs and eliminate duplicates
- sort (boundary_dofs.begin(), boundary_dofs.end());
- boundary_dofs.erase (unique (boundary_dofs.begin(),
- boundary_dofs.end()),
- boundary_dofs.end());
-
- // now finally copy the result
- // for this level its destination
- interior_boundary_dofs[level-1] = boundary_dofs;
- };
-};
-
-#endif
-
-//////////////////////////////////////////////////////////////////////
-
-void
-MGSmoother::set_zero_interior_boundary (const unsigned int level,
- Vector<double> &u) const
-{
- if (level==0)
- return;
- else
- for (vector<unsigned int>::const_iterator p=interior_boundary_dofs[level-1].begin();
- p!=interior_boundary_dofs[level-1].end(); ++p)
- u(*p) = 0;
-};
-
-//////////////////////////////////////////////////////////////////////
-
-
-// explicit instantiations
-// don't do the following instantiation in 1d, since there is a specialized
-// function there
-#if deal_II_dimension > 1
-template MGSmoother::MGSmoother (const MGDoFHandler<deal_II_dimension>&, unsigned int);
-#endif
-
-template
-MGSmootherRelaxation<float>
-::MGSmootherRelaxation(const MGDoFHandler<deal_II_dimension> &mg_dof,
- const MGLevelObject<SparseMatrix<float> > &matrix,
- const function_ptr relaxation,
- const unsigned int steps,
- const double omega);
-
-template
-MGSmootherRelaxation<double>
-::MGSmootherRelaxation(const MGDoFHandler<deal_II_dimension> &mg_dof,
- const MGLevelObject<SparseMatrix<double> > &matrix,
- const function_ptr relaxation,
- const unsigned int steps,
- const double omega);
-
+++ /dev/null
-//---------------------------- multigrid.cc ---------------------------
-// $Id$
-// Version: $Name$
-//
-// Copyright (C) 1998, 1999, 2000 by the deal.II authors
-//
-// This file is subject to QPL and may not be distributed
-// without copyright and license information. Please refer
-// to the file deal.II/doc/license.html for the text and
-// further information on this license.
-//
-//---------------------------- multigrid.cc ---------------------------
-
-
-#include <grid/tria.h>
-#include <multigrid/mg_dof_handler.h>
-#include <multigrid/mg_dof_accessor.h>
-#include <grid/tria_iterator.h>
-#include <fe/fe.h>
-#include <multigrid/multigrid.h>
-#include <multigrid/multigrid.templates.h>
-#include <multigrid/mg_smoother.h>
-#include <lac/vector.h>
-
-
-/* --------------------------------- MG ------------------------------------ */
-
-template <int dim>
-Multigrid<dim>::Multigrid (const MGDoFHandler<dim> &mg_dof_handler,
- const ConstraintMatrix &constraints,
- const MGLevelObject<SparsityPattern> &level_sparsities,
- const MGLevelObject<SparseMatrix<double> > &level_matrices,
- const MGTransferBase &transfer,
- const unsigned int minlevel,
- const unsigned int maxlevel)
- :
- MGBase(transfer, minlevel,
- min(mg_dof_handler.get_tria().n_levels()-1,
- maxlevel)),
- mg_dof_handler(&mg_dof_handler),
- level_sparsities(&level_sparsities),
- level_matrices(&level_matrices),
- constraints(&constraints)
-{};
-
-
-template <int dim>
-void
-Multigrid<dim>::level_vmult (const unsigned int level,
- Vector<double> &result,
- const Vector<double> &u,
- const Vector<double> &/* rhs */)
-{
- (*level_matrices)[level].vmult(result,u);
- result.scale(-1.);
-}
-
-
-/* ----------------------------- MGTransferPrebuilt ------------------------ */
-
-
-MGTransferPrebuilt::~MGTransferPrebuilt ()
-{};
-
-
-template <int dim>
-void MGTransferPrebuilt::build_matrices (const MGDoFHandler<dim> &mg_dof)
-{
- const unsigned int n_levels = mg_dof.get_tria().n_levels();
- const unsigned int dofs_per_cell = mg_dof.get_fe().dofs_per_cell;
-
- // reset the size of the array of
- // matrices
- prolongation_sparsities.clear ();
- prolongation_matrices.clear ();
- prolongation_sparsities.reserve (n_levels-1);
- prolongation_matrices.reserve (n_levels-1);
-
-
-// two fields which will store the
- // indices of the multigrid dofs
- // for a cell and one of its children
- vector<unsigned int> dof_indices_mother (dofs_per_cell);
- vector<unsigned int> dof_indices_child (dofs_per_cell);
-
- // for each level: first build the sparsity
- // pattern of the matrices and then build the
- // matrices themselves. note that we only
- // need to take care of cells on the coarser
- // level which have children
- for (unsigned int level=0; level<n_levels-1; ++level)
- {
- prolongation_sparsities.push_back (SparsityPattern());
- // reset the dimension of the structure.
- // note that for the number of entries
- // per row, the number of mother dofs
- // coupling to a child dof is
- // necessary. this, of course, is the
- // number of degrees of freedom per
- // cell
- prolongation_sparsities.back().reinit (mg_dof.n_dofs(level+1),
- mg_dof.n_dofs(level),
-//TODO: evil hack, must be corrected!
- dofs_per_cell+1);
-
- for (typename MGDoFHandler<dim>::cell_iterator cell=mg_dof.begin(level);
- cell != mg_dof.end(level); ++cell)
- if (cell->has_children())
- {
- cell->get_mg_dof_indices (dof_indices_mother);
-
- for (unsigned int child=0;
- child<GeometryInfo<dim>::children_per_cell; ++child)
- {
- // set an alias to the
- // prolongation matrix for
- // this child
- const FullMatrix<double> &prolongation
- = mg_dof.get_fe().prolongate(child);
-
- cell->child(child)->get_mg_dof_indices (dof_indices_child);
-
- // now tag the entries in the
- // matrix which will be used
- // for this pair of mother/child
- for (unsigned int i=0; i<dofs_per_cell; ++i)
- for (unsigned int j=0; j<dofs_per_cell; ++j)
- if (prolongation(i,j) != 0)
- {
- prolongation_sparsities[level].add
- (dof_indices_child[i],
- dof_indices_mother[j]);
- };
- };
- };
- prolongation_sparsities[level].compress ();
-
- prolongation_matrices.push_back (SparseMatrix<float>());
- prolongation_matrices[level].reinit (prolongation_sparsities[level]);
-
-
-// now actually build the matrices
- for (MGDoFHandler<dim>::cell_iterator cell=mg_dof.begin(level);
- cell != mg_dof.end(level); ++cell)
- if (cell->has_children())
- {
- cell->get_mg_dof_indices (dof_indices_mother);
-
- for (unsigned int child=0;
- child<GeometryInfo<dim>::children_per_cell; ++child)
- {
- // set an alias to the
- // prolongation matrix for
- // this child
- const FullMatrix<double> &prolongation
- = mg_dof.get_fe().prolongate(child);
-
- cell->child(child)->get_mg_dof_indices (dof_indices_child);
-
- // now set the entries in the
- // matrix
- for (unsigned int i=0; i<dofs_per_cell; ++i)
- for (unsigned int j=0; j<dofs_per_cell; ++j)
- if (prolongation(i,j) != 0)
- {
- prolongation_matrices[level].set
- (dof_indices_child[i],
- dof_indices_mother[j],
- prolongation(i,j));
- };
- };
- };
- };
-};
-
-
-void MGTransferPrebuilt::prolongate (const unsigned int to_level,
- Vector<double> &dst,
- const Vector<double> &src) const
-{
- Assert ((to_level >= 1) && (to_level<=prolongation_matrices.size()),
- ExcIndexRange (to_level, 1, prolongation_matrices.size()+1));
-
- prolongation_matrices[to_level-1].vmult (dst, src);
-};
-
-
-void MGTransferPrebuilt::restrict_and_add (const unsigned int from_level,
- Vector<double> &dst,
- const Vector<double> &src) const
-{
- Assert ((from_level >= 1) && (from_level<=prolongation_matrices.size()),
- ExcIndexRange (from_level, 1, prolongation_matrices.size()+1));
-
- prolongation_matrices[from_level-1].Tvmult_add (dst, src);
-};
-
-
-// explicit instatations
-template class Multigrid<deal_II_dimension>;
-template
-void MGTransferPrebuilt::build_matrices (const MGDoFHandler<deal_II_dimension> &mg_dof);
-
-template
-void Multigrid<deal_II_dimension>::copy_to_mg (const Vector<float>& src);
-
-template
-void Multigrid<deal_II_dimension>::copy_to_mg (const Vector<double>& src);
-
-template
-void Multigrid<deal_II_dimension>::copy_from_mg (Vector<float>& src) const;
-
-template
-void Multigrid<deal_II_dimension>::copy_from_mg (Vector<double>& src) const;
-
-