]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add GridGenerator::create_triangulation_with_removed_cells().
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 13 Mar 2015 13:17:41 +0000 (08:17 -0500)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 13 Mar 2015 13:17:41 +0000 (08:17 -0500)
doc/news/changes.h
include/deal.II/grid/grid_generator.h
source/grid/grid_generator.cc
source/grid/grid_generator.inst.in
tests/grid/grid_generator_07.cc [new file with mode: 0644]
tests/grid/grid_generator_07.output [new file with mode: 0644]

index 909ba4aa0891d7012c2bd368c6b584a0c13a8733..28ad3ff801393c85dbb51b791b5d073f6df5967f 100644 (file)
@@ -376,6 +376,12 @@ inconvenience this causes.
 <h3>Specific improvements</h3>
 
 <ol>
+  <li> New: GridGenerator::create_triangulation_with_removed_cells() creates
+  a new mesh out of an existing one by dropping individual cells.
+  <br>
+  (Wolfgang Bangerth, 2015/03/13)
+  </li>
+
   <li> New: Add MueLu preconditioner from Trilinos. This is a new algebraic 
   multigrid package. The input parameters are almost the same as the ones 
   from ML so that the two preconditioners can be easily swapped.
index 11d691e452e78bc8e31b900e2d138824f15955bb..6b29f7f0182bef67e2ae86949d129d557642a683 100644 (file)
@@ -705,6 +705,44 @@ namespace GridGenerator
                               const Triangulation<dim, spacedim> &triangulation_2,
                               Triangulation<dim, spacedim>       &result);
 
+  /**
+   * This function creates a triangulation that consists of the same cells
+   * as are present in the first argument, except those cells that are listed
+   * in the second argument. The purpose of the function is to generate
+   * geometries <i>subtractively</i> from the geometry described by an
+   * existing triangulation. A prototypical case is a 2d domain with
+   * rectangular holes. This can be achieved by first meshing the entire
+   * domain and then using this function to get rid of the cells that
+   * are located at the holes. Likewise, you could create the mesh
+   * that GridGenerator::hyper_L() produces by starting with a
+   * GridGenerator::hyper_cube(), refining it once, and then calling
+   * the current function with a single cell in the second argument.
+   *
+   * @param[in] input_triangulation The original triangulation that serves as
+   *   the template from which the new one is to be created.
+   * @param[in] cells_to_remove A list of cells of the triangulation provided
+   *   as first argument that should be removed (i.e., that should not
+   *   show up in the result.
+   * @param[out] result The resulting triangulation that consists of the same
+   *   cells as are in @p input_triangulation, with the exception of the cells
+   *   listed in @p cells_to_remove.
+   *
+   * @pre Because we cannot create triangulations de novo that contain
+   *   adaptively refined cells, the input triangulation needs to have all
+   *   of its cells on the same level. Oftentimes, this will in fact be
+   *   the coarsest level, but it is allowed to pass in a triangulation
+   *   that has been refined <i>globally</i> a number of times. The output
+   *   triangulation will in that case simply be a mesh with only one
+   *   level that consists of the active cells of the input minus the
+   *   ones listed in the second argument. However, the input triangulation
+   *   must not have been <i>adaptively</i> refined.
+   */
+  template <int dim, int spacedim>
+  void
+  create_triangulation_with_removed_cells (const Triangulation<dim, spacedim> &input_triangulation,
+                                           const std::set<typename Triangulation<dim, spacedim>::active_cell_iterator> &cells_to_remove,
+                                           Triangulation<dim, spacedim>       &result);
+
 
   /**
    * Take a 2d Triangulation that is being extruded in z direction by the
index 6ff40765ffccddd4ba2daf781b88fff6e27d96d8..f6048ba87067f5357e8c021dc4c4b1db4daceb0c 100644 (file)
@@ -3438,6 +3438,51 @@ namespace GridGenerator
 
 
 
+  template <int dim, int spacedim>
+  void
+  create_triangulation_with_removed_cells (const Triangulation<dim, spacedim> &input_triangulation,
+                                           const std::set<typename Triangulation<dim, spacedim>::active_cell_iterator> &cells_to_remove,
+                                           Triangulation<dim, spacedim>       &result)
+  {
+    // simply copy the vertices; we will later strip those
+    // that turn out to be unused
+    std::vector<Point<spacedim> > vertices = input_triangulation.get_vertices();
+
+    // the loop through the cells and copy stuff, excluding
+    // the ones we are to remove
+    std::vector<CellData<dim> > cells;
+    for (typename Triangulation<dim,spacedim>::active_cell_iterator
+         cell = input_triangulation.begin_active(); cell != input_triangulation.end(); ++cell)
+      if (cells_to_remove.find(cell) == cells_to_remove.end())
+        {
+          Assert (static_cast<unsigned int>(cell->level()) == input_triangulation.n_levels()-1,
+                  ExcMessage ("Your input triangulation appears to have "
+                              "adaptively refined cells. This is not allowed. You can "
+                              "only call this function on a triangulation in which "
+                              "all cells are on the same refinement level."));
+
+          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);
+        }
+
+    // throw out duplicated vertices from the two meshes, reorder vertices as
+    // necessary and create the triangulation
+    SubCellData subcell_data;
+    std::vector<unsigned int> considered_vertices;
+    GridTools::delete_duplicated_vertices (vertices, cells,
+                                           subcell_data,
+                                           considered_vertices);
+
+    // then clear the old triangulation and create the new one
+    result.clear ();
+    result.create_triangulation (vertices, cells, subcell_data);
+  }
+
+
+
   void
   extrude_triangulation(const Triangulation<2, 2> &input,
                         const unsigned int n_slices,
index eb1e25278c1baa981926443c18e685c2673bfd7f..f17220b8d1a5383b005070daff6c3696af0a845f 100644 (file)
@@ -42,6 +42,12 @@ namespace GridGenerator
       (const Triangulation<deal_II_dimension, deal_II_space_dimension> &triangulation_1,
        const Triangulation<deal_II_dimension, deal_II_space_dimension> &triangulation_2,
        Triangulation<deal_II_dimension, deal_II_space_dimension>       &result);
+
+    template
+      void
+      create_triangulation_with_removed_cells (const Triangulation<deal_II_dimension, deal_II_space_dimension> &input_triangulation,
+                                               const std::set<typename Triangulation<deal_II_dimension, deal_II_space_dimension>::active_cell_iterator> &cells_to_remove,
+                                               Triangulation<deal_II_dimension, deal_II_space_dimension>       &result);
       
 #endif
 \}  
diff --git a/tests/grid/grid_generator_07.cc b/tests/grid/grid_generator_07.cc
new file mode 100644 (file)
index 0000000..e6c0976
--- /dev/null
@@ -0,0 +1,65 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2005 - 2015 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+
+// Test GridGenerator::create_triangulation_with_removed_cells
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/base/tensor.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/fe/fe_q.h>
+
+#include <fstream>
+#include <iomanip>
+
+
+template<int dim>
+void test(std::ostream &out)
+{
+  Triangulation<dim> triangulation;
+  Triangulation<dim> tr;
+  GridGenerator::hyper_cube (triangulation);
+  triangulation.refine_global (3);
+
+  // remove all cells but the first. this is the hardest case to handle as it
+  // makes a bunch of vertices unused
+  std::set<typename Triangulation<dim>::active_cell_iterator>
+    cells_to_remove;
+  for (typename Triangulation<dim>::active_cell_iterator
+        cell = ++triangulation.begin_active();
+       cell != triangulation.end(); ++cell)
+    cells_to_remove.insert (cell);
+  
+  GridGenerator::create_triangulation_with_removed_cells(triangulation,
+                                                        cells_to_remove, tr);
+  GridOut go;
+  go.write_gnuplot(tr, out);
+}
+
+
+int main()
+{
+  std::ofstream logfile("output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  test<2>(logfile);
+  test<3>(logfile);
+}
diff --git a/tests/grid/grid_generator_07.output b/tests/grid/grid_generator_07.output
new file mode 100644 (file)
index 0000000..e4ce372
--- /dev/null
@@ -0,0 +1,32 @@
+
+0.00000 0.00000 0 0
+0.125000 0.00000 0 0
+0.125000 0.125000 0 0
+0.00000 0.125000 0 0
+0.00000 0.00000 0 0
+
+
+0.00000 0.00000 0.00000 0 0
+0.125000 0.00000 0.00000 0 0
+0.125000 0.00000 0.125000 0 0
+0.00000 0.00000 0.125000 0 0
+0.00000 0.00000 0.00000 0 0
+
+0.00000 0.125000 0.00000 0 0
+0.125000 0.125000 0.00000 0 0
+0.125000 0.125000 0.125000 0 0
+0.00000 0.125000 0.125000 0 0
+0.00000 0.125000 0.00000 0 0
+
+0.00000 0.00000 0.00000 0 0
+0.00000 0.125000 0.00000 0 0
+
+0.125000 0.00000 0.00000 0 0
+0.125000 0.125000 0.00000 0 0
+
+0.125000 0.00000 0.125000 0 0
+0.125000 0.125000 0.125000 0 0
+
+0.00000 0.00000 0.125000 0 0
+0.00000 0.125000 0.125000 0 0
+

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.