]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
First version.
authorcvs <cvs@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 8 Oct 1999 19:06:49 +0000 (19:06 +0000)
committercvs <cvs@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 8 Oct 1999 19:06:49 +0000 (19:06 +0000)
git-svn-id: https://svn.dealii.org/trunk@1751 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/grid/intergrid_map.h [new file with mode: 0644]
deal.II/deal.II/source/grid/intergrid_map.cc [new file with mode: 0644]

diff --git a/deal.II/deal.II/include/grid/intergrid_map.h b/deal.II/deal.II/include/grid/intergrid_map.h
new file mode 100644 (file)
index 0000000..ed15c1b
--- /dev/null
@@ -0,0 +1,148 @@
+/*----------------------------   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     ---------------------------*/
diff --git a/deal.II/deal.II/source/grid/intergrid_map.cc b/deal.II/deal.II/source/grid/intergrid_map.cc
new file mode 100644 (file)
index 0000000..a0af3fd
--- /dev/null
@@ -0,0 +1,176 @@
+/* $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>;
+
+  

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.