]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add GridGenerator for a general cell. 3155/head
authorBruno Turcksin <bruno.turcksin@gmail.com>
Thu, 22 Sep 2016 21:39:05 +0000 (17:39 -0400)
committerBruno Turcksin <bruno.turcksin@gmail.com>
Mon, 3 Oct 2016 12:43:09 +0000 (08:43 -0400)
include/deal.II/grid/grid_generator.h
source/grid/grid_generator.cc
source/grid/grid_generator.inst.in
tests/grid/grid_generator_general_cell.cc [new file with mode: 0644]
tests/grid/grid_generator_general_cell.output [new file with mode: 0644]

index 267d57e8ba31bb081319c3fca97844e8322c4e6c..77db7ec137d9e16f3b13670e62357390946e7cd9 100644 (file)
@@ -52,11 +52,12 @@ namespace GridGenerator
    * in 2D, etc) consisting of exactly one cell. The hypercube volume is the
    * tensor product interval $[left,right]^{\text{dim}}$ in the present number
    * of dimensions, where the limits are given as arguments. They default to
-   * zero and unity, then producing the unit hypercube. If the argument @p
-   * colorize is false, all boundary indicators are set to zero ("not
-   * colorized") for 2d and 3d. If it is true, the boundary is colorized as in
-   * hyper_rectangle(). In 1d the indicators are always colorized, see
-   * hyper_rectangle().
+   * zero and unity, then producing the unit hypercube.
+   *
+   * If the argument @p colorize is false, all boundary indicators are set to
+   * zero ("not colorized") for 2d and 3d. If it is true, the boundary is
+   * colorized as in hyper_rectangle(). In 1d the indicators are always
+   * colorized, see hyper_rectangle().
    *
    * @image html hyper_cubes.png
    *
@@ -291,6 +292,26 @@ namespace GridGenerator
   cheese (Triangulation<dim, spacedim> &tria,
           const std::vector<unsigned int> &holes);
 
+  /**
+   * A general quadrilateral in 2d or a general hexahedron in 3d. It is the
+   * responsibility of the user to provide the vertices in the right order (see
+   * the documentation of the GeometryInfo class) because the vertices are stored
+   * in the same order as they are given. It is also important to make sure that
+   * the volume of the cell is positive.
+   *
+   * If the argument @p colorize is false, all boundary indicators are set to
+   * zero ("not colorized") for 2d and 3d. If it is true, the boundary is
+   * colorized as in hyper_rectangle(). In 1d the indicators are always
+   * colorized, see hyper_rectangle().
+   *
+   * @author Bruno Turcksin
+   */
+  template <int dim>
+  void
+  general_cell(Triangulation<dim> &tria,
+               const std::vector<Point<dim> > &vertices,
+               const bool colorize = false);
+
   /**
    * A parallelogram. The first corner point is the origin. The @p dim
    * adjacent points are the ones given in the second argument and the fourth
index b16b56b59bae9eaa07904788f4ddea5ad6fcaea3..a4a7aa7ca802622437f3363309ecb3f97ca675ce 100644 (file)
@@ -749,6 +749,33 @@ namespace GridGenerator
     tria.set_all_manifold_ids_on_boundary(0);
   }
 
+
+
+  template <int dim>
+  void
+  general_cell (Triangulation<dim> &tria,
+                const std::vector<Point<dim> > &vertices,
+                const bool colorize)
+  {
+    Assert (vertices.size() == dealii::GeometryInfo<dim>::vertices_per_cell,
+            ExcMessage("Wrong number of vertices."));
+
+    // First create a hyper_rectangle and then deform it.
+    hyper_cube(tria, 0, 1, colorize);
+
+    typename Triangulation<dim>::active_cell_iterator cell = tria.begin_active();
+    for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
+      cell->vertex(i) = vertices[i];
+
+    // Check that the order of the vertices makes sense, i.e., the volume of the
+    // cell is positive.
+    Assert(GridTools::volume(tria) > 0.,
+           ExcMessage("The volume of the cell is not greater than zero. "
+                      "This could be due to the wrong ordering of the vertices."));
+  }
+
+
+
   template<>
   void
   parallelogram (Triangulation<3> &,
index 3d60241d21e185ba53e20c53bc601af2036c5238..f2af7dd898b71a44d645a6308c4313ad0148b6e0 100644 (file)
@@ -122,6 +122,12 @@ for (deal_II_dimension : DIMENSIONS)
         const Point<deal_II_dimension>                &,
         const bool                       );
 
+    template
+    void
+    general_cell (Triangulation<deal_II_dimension> &,
+                  const std::vector<Point<deal_II_dimension> > &,
+                  const bool);
+
     template void
     simplex<deal_II_dimension> (
         Triangulation<deal_II_dimension, deal_II_dimension> &,
diff --git a/tests/grid/grid_generator_general_cell.cc b/tests/grid/grid_generator_general_cell.cc
new file mode 100644 (file)
index 0000000..edbd294
--- /dev/null
@@ -0,0 +1,84 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 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 output for GridGenerator::general_cell()
+
+#include "../tests.h"
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_out.h>
+
+void dim2(std::ostream &os)
+{
+  std::vector<Point<2> > vertices(4);
+  vertices[0](0) = -1.;
+  vertices[0](1) = -1.;
+  vertices[1](0) = 1.;
+  vertices[1](1) = -1.5;
+  vertices[2](0) = 1.5;
+  vertices[2](1) = 1.5;
+  vertices[3](0) = 2.;
+  vertices[3](1) = 0.5;
+
+  Triangulation<2> tria;
+  GridGenerator::general_cell<2> (tria, vertices);
+
+  GridOut gout;
+  gout.write_vtk(tria, os);
+}
+
+void dim3(std::ostream &os)
+{
+  std::vector<Point<3> > vertices(8);
+  vertices[0](0) = -1.;
+  vertices[0](1) = -1.;
+  vertices[0](2) = -1.;
+  vertices[1](0) = 1.;
+  vertices[1](1) = -1.5;
+  vertices[1](2) = -1.5;
+  vertices[2](0) = 2.;
+  vertices[2](1) = 1.5;
+  vertices[2](2) = -2.;
+  vertices[3](0) = 2.5;
+  vertices[3](1) = 0.5;
+  vertices[3](2) = -3.;
+  vertices[4](0) = -1.;
+  vertices[4](1) = -1.;
+  vertices[4](2) = 1.;
+  vertices[5](0) = 1.;
+  vertices[5](1) = -1.5;
+  vertices[5](2) = 1.5;
+  vertices[6](0) = 2.;
+  vertices[6](1) = 1.5;
+  vertices[6](2) = 2.;
+  vertices[7](0) = 2.;
+  vertices[7](1) = 0.5;
+  vertices[7](2) = 3.;
+
+  Triangulation<3> tria;
+  GridGenerator::general_cell<3> (tria, vertices);
+
+  GridOut gout;
+  gout.write_vtk(tria, os);
+}
+
+
+int main()
+{
+  initlog(true);
+  std::ostream &logfile = deallog.get_file_stream();
+  dim2(logfile);
+  dim3(logfile);
+}
diff --git a/tests/grid/grid_generator_general_cell.output b/tests/grid/grid_generator_general_cell.output
new file mode 100644 (file)
index 0000000..edf1cc7
--- /dev/null
@@ -0,0 +1,69 @@
+
+# vtk DataFile Version 3.0
+#This file was generated 
+ASCII
+DATASET UNSTRUCTURED_GRID
+
+POINTS 4 double
+-1.00000 -1.00000 0
+1.00000 -1.50000 0
+1.50000 1.50000 0
+2.00000 0.500000 0
+
+CELLS 1 5
+4      0       1       3       2
+
+CELL_TYPES 1
+ 9
+POINT_DATA 4
+SCALARS level double 1
+LOOKUP_TABLE default
+0.00000 0.00000 0.00000 0.00000 
+SCALARS manifold double 1
+LOOKUP_TABLE default
+-1.00000 -1.00000 -1.00000 -1.00000 
+SCALARS material double 1
+LOOKUP_TABLE default
+0.00000 0.00000 0.00000 0.00000 
+SCALARS subdomain double 1
+LOOKUP_TABLE default
+0.00000 0.00000 0.00000 0.00000 
+SCALARS level_subdomain double 1
+LOOKUP_TABLE default
+0.00000 0.00000 0.00000 0.00000 
+# vtk DataFile Version 3.0
+#This file was generated 
+ASCII
+DATASET UNSTRUCTURED_GRID
+
+POINTS 8 double
+-1.00000 -1.00000 -1.00000
+1.00000 -1.50000 -1.50000
+2.00000 1.50000 -2.00000
+2.50000 0.500000 -3.00000
+-1.00000 -1.00000 1.00000
+1.00000 -1.50000 1.50000
+2.00000 1.50000 2.00000
+2.00000 0.500000 3.00000
+
+CELLS 1 9
+8      0       1       3       2       4       5       7       6
+
+CELL_TYPES 1
+ 12
+POINT_DATA 8
+SCALARS level double 1
+LOOKUP_TABLE default
+0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
+SCALARS manifold double 1
+LOOKUP_TABLE default
+-1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 
+SCALARS material double 1
+LOOKUP_TABLE default
+0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
+SCALARS subdomain double 1
+LOOKUP_TABLE default
+0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
+SCALARS level_subdomain double 1
+LOOKUP_TABLE default
+0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 

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.