]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
New function GridGenerator::merge_triangulations.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 17 Jun 2011 19:55:07 +0000 (19:55 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 17 Jun 2011 19:55:07 +0000 (19:55 +0000)
git-svn-id: https://svn.dealii.org/trunk@23834 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/grid/grid_generator.h
deal.II/source/grid/grid_generator.cc
deal.II/source/grid/grid_generator.inst.in

index c0602a73fd78a852d7826d4cd68fbecfa392743d..6ed08bb76d3a7b4e72771816cb3954d58a4dbef7 100644 (file)
@@ -198,6 +198,12 @@ should be fixed now.
 <h3>Specific improvements</h3>
 
 <ol>
+<li> New: The new function GridGenerator::merge_triangulations can be used to compose
+coarse meshes from simpler ones by merging their cells into a single
+triangulation object.
+<br>
+(Wolfgang Bangerth 2011/06/17)
+
 <li> Fixed: If an FEValues object was kept around until after the triangulation
 on which it works has been refined or coarsened, and is then reinitialized
 with a cell from the refined triangulation, it could compute wrong results or
index cd8a19e46ada8c6321c79ddd3701e87ba3d449ed..4c236005e3bccde11c4a5bc34f27b9dd3a09fdf3 100644 (file)
@@ -50,7 +50,7 @@ template <typename number> class SparseMatrix;
  * @ingroup grid
  * @author Wolfgang Bangerth, Ralf Hartmann, Guido Kanschat, Stefan
  * Nauber, Joerg Weimar, Yaqi Wang, Luca Heltai, 1998, 1999, 2000, 2001, 2002,
- * 2003, 2006, 2007, 2008.
+ * 2003, 2006, 2007, 2008, 2009, 2010, 2011.
  */
 class GridGenerator
 {
@@ -772,6 +772,55 @@ class GridGenerator
                         const double         r);
 
                                     /**
+                                     * Given the two triangulations
+                                     * specified as the first two
+                                     * arguments, create the
+                                     * triangulation that contains
+                                     * the cells of both
+                                     * triangulation and store it in
+                                     * the third parameter. Previous
+                                     * content of @p result will be
+                                     * deleted.
+                                     * 
+                                     * This function is most often used 
+                                     * to compose meshes for more
+                                     * complicated geometries if the
+                                     * geometry can be composed of
+                                     * simpler parts for which functions
+                                     * exist to generate coarse meshes.
+                                     * For example, the channel mesh used
+                                     * in step-35 could in principle be
+                                     * created using a mesh created by the 
+                                     * GridGenerator::hyper_cube_with_cylindrical_hole
+                                     * function and several rectangles,
+                                     * and merging them using the current
+                                     * function. The rectangles will
+                                     * have to be translated to the
+                                     * right for this, a task that can
+                                     * be done using the GridTools::shift
+                                     * function (other tools to transform
+                                     * individual mesh building blocks are
+                                     * GridTools::transform, GridTools::rotate,
+                                     * and GridTools::scale).
+                                     * 
+                                     * @note The two input triangulations
+                                     * must be coarse meshes that have
+                                     * no refined cells. 
+                                     * 
+                                     * @note for a related operation
+                                     * on refined meshes when both
+                                     * meshes are derived from the
+                                     * same coarse mesh, see
+                                     * GridTools::create_union_triangulation .
+                                     */
+    template <int dim, int spacedim>
+    static
+    void
+    merge_triangulations (const Triangulation<dim, spacedim> &triangulation_1,
+                         const Triangulation<dim, spacedim> &triangulation_2,
+                         Triangulation<dim, spacedim>       &result);
+
+                                     /**
                                      * This function transformes the
                                      * @p Triangulation @p tria
                                      * smoothly to a domain that is
index d9199b09f2e05fc637708a67aacab4055d1c3d18..08df26b1e60f7633eaa36bf863d86c7eaeb8e1f3 100644 (file)
@@ -2685,6 +2685,59 @@ void GridGenerator::cylinder_shell (Triangulation<3>   &tria,
 
 
 
+template <int dim, int spacedim>
+void
+GridGenerator::
+merge_triangulations (const Triangulation<dim, spacedim> &triangulation_1,
+                     const Triangulation<dim, spacedim> &triangulation_2,
+                     Triangulation<dim, spacedim>       &result)
+{
+  Assert (triangulation_1.n_levels() == 1,
+         ExcMessage ("The input triangulations must be coarse meshes."));
+  Assert (triangulation_2.n_levels() == 1,
+         ExcMessage ("The input triangulations must be coarse meshes."));
+  
+  // get the union of the set of vertices
+  std::vector<Point<spacedim> > vertices = triangulation_1.get_vertices();
+  vertices.insert (vertices.end(),
+                  triangulation_2.get_vertices().begin(),
+                  triangulation_2.get_vertices().end());
+  
+  // now form the union of the set of cells
+  std::vector<CellData<dim> > cells;
+  cells.reserve (triangulation_1.n_cells() + triangulation_2.n_cells());
+  for (typename Triangulation<dim,spacedim>::cell_iterator
+    cell = triangulation_1.begin(); cell != triangulation_1.end(); ++cell)
+  {
+    CellData<dim> this_cell;
+    for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
+      this_cell.vertices[v] = cell->vertex_index(v);
+    this_cell.material_id = cell->material_id();
+    cells.push_back (this_cell);
+  }
+  
+  // now do the same for the other other mesh. note that we have to 
+  // translate the vertex indices
+  for (typename Triangulation<dim,spacedim>::cell_iterator
+    cell = triangulation_2.begin(); cell != triangulation_2.end(); ++cell)
+  {
+    CellData<dim> this_cell;
+    for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
+      this_cell.vertices[v] = cell->vertex_index(v) + triangulation_1.n_vertices();
+    this_cell.material_id = cell->material_id();
+    cells.push_back (this_cell);
+  }
+  
+  // throw out duplicated vertices from the two meshes
+  // and create the triangulation
+  SubCellData subcell_data;
+  std::vector<unsigned int> considered_vertices;
+  GridTools::delete_duplicated_vertices (vertices, cells, subcell_data, considered_vertices);
+  result.clear ();
+  result.create_triangulation (vertices, cells, subcell_data);
+}
+
+
 
 // make the following function inline. this is so that it becomes marked
 // internal/weak for the linker and we don't get multiply defined errors
index 076f9b4a6a81687ea98bf1cc7b5877cd7c25ca68..0ec9986dafabbe5f6a4ccc285793859e4bb4058b 100644 (file)
@@ -51,6 +51,20 @@ for (deal_II_dimension : DIMENSIONS)
       const bool                       );
 
 
+    template
+    void
+    GridGenerator::merge_triangulations (const Triangulation<deal_II_dimension,deal_II_dimension> &triangulation_1,
+                     const Triangulation<deal_II_dimension,deal_II_dimension> &triangulation_2,
+                     Triangulation<deal_II_dimension,deal_II_dimension>       &result);
+
+#if deal_II_dimension < 3
+    template
+    void
+    GridGenerator::merge_triangulations (const Triangulation<deal_II_dimension,deal_II_dimension+1> &triangulation_1,
+                     const Triangulation<deal_II_dimension,deal_II_dimension+1> &triangulation_2,
+                     Triangulation<deal_II_dimension,deal_II_dimension+1>       &result);
+#endif
+    
 #if deal_II_dimension > 1
   template void
     GridGenerator::

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.