--- /dev/null
+/*---------------------------- intergrid_map.h ---------------------------*/
+/* $Id$ */
+#ifndef __intergrid_map_H
+#define __intergrid_map_H
+/*---------------------------- intergrid_map.h ---------------------------*/
+
+
+
+/**
+ * This class provides a map between two grids which are derived from
+ * the same coarse grid. For each cell iterator of the source map, it provides
+ * the respective cell iterator on the destination map, through its
+ * #operator []#.
+ *
+ * Usually, the two grids will be refined differently. Then, the value
+ * returned for an iterator on the source grid will be either:
+ * \begin{itemize}
+ * \item The same cell on the destination grid, if it exists there;
+ * \item The most refined cell of the destination grid from which the
+ * pendant of the source cell could be obtained by refinement. This
+ * cell is always active and has a refinement level less than that
+ * of the source cell.
+ * \end{itemize}
+ * Keys for this map are all cells on the source grid, whether active or
+ * not.
+ *
+ * The implementation of this class is such that not only cell iterators
+ * into triangulations can be mapped, but also iterators into objects of
+ * type #DoFHandler# and #MGDoFHandler#. The extension to other classes
+ * offering iterator functions and some minor additional requirements is
+ * simple.
+ *
+ * In practice, use of this class is as follows:
+ * \begin{verbatim}
+ * // have two grids, which are derived from the
+ * // same coarse grid
+ * Triangulation<dim> tria1, tria2;
+ * DoFHandler<dim> dof_handler_1(tria1), dof_handler_2(tria2);
+ * ...
+ * // do something with these objects, e.g.
+ * // refine the triangulations differently,
+ * // distribute degrees of freedom, etc
+ * ...
+ * // create the mapping
+ * InterGridMap<DoFHandler,dim> grid_1_to_2_map;
+ * grid_1_to_2_map.make_mapping (dof_handler_1,
+ * dof_handler_2);
+ * ...
+ * typename DoFHandler<dim>::cell_iterator cell = dof_handler_1.begin(),
+ * endc = dof_handler_1.end();
+ * for (; cell!=endc; ++cell)
+ * // now do something with the cell of dof_handler_2
+ * // corresponding to #cell# (which is one of
+ * // dof_handler_1
+ * f( grid_1_to_2_map[cell]);
+ * \end{verbatim}
+ *
+ * Note that the template parameters to this class have to be given as
+ * #InterGridMap<DoFHandler,2>#, i.e. the dimension is given explicitely and
+ * no dimension is attributed to the first parameter, which here is
+ * #DoFHandler# (and could equally well be #Triangulation# or #MGDoFHandler#).
+ *
+ * @author Wolfgang Bangerth, 1999
+ */
+template <template <int> class GridClass, int dim>
+class InterGridMap
+{
+ public:
+ /**
+ * Typedef to the iterator type of
+ * the grid class under consideration.
+ */
+ typedef typename GridClass<dim>::cell_iterator cell_iterator;
+
+ /**
+ * Create the mapping between the two
+ * grids.
+ */
+ void make_mapping (const GridClass<dim> &source_grid,
+ const GridClass<dim> &destination_grid);
+
+ /**
+ * Access operator: give a cell
+ * on the source grid and receive
+ * the respective cell on the
+ * other grid, or if that does not
+ * exist, the most refined cell
+ * of which the source cell would
+ * be created if it were further
+ * refined.
+ */
+ cell_iterator operator [] (const cell_iterator &source_cell) const;
+
+ /**
+ * Delete all data of this class.
+ */
+ void clear ();
+
+ /**
+ * Exception
+ */
+ DeclException1 (ExcInvalidKey,
+ cell_iterator,
+ << "The iterator " << arg1 << " is not valid as key for "
+ << "this map.");
+ /**
+ * Exception
+ */
+ DeclException0 (ExcIncompatibleGrids);
+
+ private:
+ /**
+ * The actual data. Hold one iterator
+ * for each cell on each level.
+ */
+ vector<vector<cell_iterator> > mapping;
+
+ /**
+ * Set the mapping for the pair of
+ * cells given. These shall match
+ * in level of refinement and all
+ * other properties.
+ */
+ void set_mapping (const cell_iterator &src_cell,
+ const cell_iterator &dst_cell);
+
+ /**
+ * Set the value of the key #src_cell#
+ * to #dst_cell#. Do so as well for
+ * all the children and their children
+ * of #src_cell#. This function is
+ * used for cells which are more
+ * refined on #src_grid# than on
+ * #dst_grid#; then all values of
+ * the hierarchy of cells and their
+ * children point to one cell on the
+ * #dst_grid#.
+ */
+ void set_entries_to_cell (const cell_iterator &src_cell,
+ const cell_iterator &dst_cell);
+};
+
+
+
+/*---------------------------- intergrid_map.h ---------------------------*/
+/* end of #ifndef __intergrid_map_H */
+#endif
+/*---------------------------- intergrid_map.h ---------------------------*/
--- /dev/null
+/* $Id$ */
+
+#include <grid/tria.h>
+#include <grid/dof.h>
+#include <grid/mg_dof.h>
+#include <grid/tria_accessor.h>
+#include <grid/dof_accessor.h>
+#include <grid/mg_dof_accessor.h>
+#include <grid/tria_iterator.h>
+#include <grid//intergrid_map.h>
+
+
+
+
+// helper function to aquire the number of levels within a grid
+template <class GridClass>
+unsigned int
+get_n_levels (const GridClass &grid)
+{
+ // all objects we deal with are able
+ // to deliver a pointer to the
+ // underlying triangulation.
+ //
+ // for the triangulation as GridClass
+ // of this object, there is a
+ // specialization of this function
+ return grid.get_tria().n_levels();
+};
+
+
+// specialization for grid==tria
+template <int dim>
+unsigned int
+get_n_levels (const Triangulation<dim> &grid)
+{
+ // if GridClass==Triangulation, then
+ // we can ask directly.
+ return grid.n_levels();
+};
+
+
+
+
+
+template <template <int> class GridClass, int dim>
+void InterGridMap<GridClass,dim>::make_mapping (const GridClass<dim> &source_grid,
+ const GridClass<dim> &destination_grid)
+{
+ // first delete all contents
+ clear ();
+
+
+ // then set up the containers from
+ // scratch and fill them with end-iterators
+ const unsigned int n_levels = get_n_levels(source_grid);
+ mapping.resize (n_levels);
+ for (unsigned int level=0; level<n_levels; ++level)
+ {
+ // first find out about the highest
+ // index used on this level. We could
+ // in principle ask the triangulation
+ // about this, but we would have to
+ // know the underlying data structure
+ // for this and we would like to
+ // avoid such knowledge here
+ unsigned int n_cells = 0;
+ cell_iterator cell = source_grid.begin(level),
+ endc = source_grid.end(level);
+ for (; cell!=endc; ++cell)
+ if (static_cast<unsigned int>(cell->index()) > n_cells)
+ n_cells = cell->index();
+
+ // note: n_cells is now the largest
+ // zero-based index, but we need the
+ // number of cells, which is one larger
+ mapping[level].resize (n_cells+1, destination_grid.end());
+ };
+
+ // now make up the mapping
+ // loop over all cells and set the user
+ // pointers as well as the contents of
+ // the two arrays. note that the function
+ // takes a *reference* to the int and
+ // this may change it
+ cell_iterator src_cell = source_grid.begin(0),
+ dst_cell = destination_grid.begin(0),
+ endc = source_grid.end(0);
+ for (; src_cell!=endc; ++src_cell, ++dst_cell)
+ set_mapping (src_cell, dst_cell);
+
+ // little assertion that the two grids
+ // are indeed related:
+ Assert (dst_cell == destination_grid.end(0),
+ ExcIncompatibleGrids ());
+};
+
+
+
+template <template <int> class GridClass, int dim>
+void
+InterGridMap<GridClass,dim>::set_mapping (const cell_iterator &src_cell,
+ const cell_iterator &dst_cell)
+{
+ // first set the map for this cell
+ mapping[src_cell->level()][src_cell->index()] = dst_cell;
+
+ // if both cells have children, we may
+ // recurse further into the hierarchy
+ if (src_cell->has_children() && dst_cell->has_children())
+ for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c)
+ set_mapping (src_cell->child(c),
+ dst_cell->child(c));
+ else
+ if (src_cell->has_children() &&
+ !dst_cell->has_children())
+ // src grid is more refined here.
+ // set entries for all children
+ // of this cell to the one
+ // dst_cell
+ for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c)
+ set_entries_to_cell (src_cell->child(c),
+ dst_cell);
+ // else (no cell is refined or
+ // dst_cell is refined): no pointers
+ // to be set
+};
+
+
+
+template <template <int> class GridClass, int dim>
+void
+InterGridMap<GridClass,dim>::set_entries_to_cell (const cell_iterator &src_cell,
+ const cell_iterator &dst_cell)
+{
+ // first set the map for this cell
+ mapping[src_cell->level()][src_cell->index()] = dst_cell;
+
+ // then do so for the children as well
+ for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c)
+ set_entries_to_cell (src_cell->child(c),
+ dst_cell);
+};
+
+
+
+template <template <int> class GridClass, int dim>
+typename InterGridMap<GridClass,dim>::cell_iterator
+InterGridMap<GridClass,dim>::operator [] (const cell_iterator &source_cell) const
+{
+ Assert (source_cell.state() == valid,
+ ExcInvalidKey (source_cell));
+ Assert (source_cell->level() <= static_cast<int>(mapping.size()),
+ ExcInvalidKey (source_cell));
+ Assert (source_cell->index() <= static_cast<int>(mapping[source_cell->level()].size()),
+ ExcInvalidKey (source_cell));
+
+ return mapping[source_cell->level()][source_cell->index()];
+};
+
+
+
+template <template <int> class GridClass, int dim>
+void InterGridMap<GridClass,dim>::clear ()
+{
+ mapping.clear ();
+};
+
+
+
+
+// explicit instantiations
+template class InterGridMap<Triangulation, deal_II_dimension>;
+template class InterGridMap<DoFHandler, deal_II_dimension>;
+template class InterGridMap<MGDoFHandler, deal_II_dimension>;
+
+