]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Generalizing GridGenerator::general_cell to tria<dim,spacedim> 7654/head
authorGiovanni Alzetta <giovannialzetta@hotmail.it>
Tue, 29 Jan 2019 15:29:15 +0000 (16:29 +0100)
committerGiovanni Alzetta <giovannialzetta@hotmail.it>
Wed, 30 Jan 2019 08:50:55 +0000 (09:50 +0100)
doc/news/changes/minor/20190130GiovanniAlzetta [new file with mode: 0644]
include/deal.II/grid/grid_generator.h
source/grid/grid_generator.cc
source/grid/grid_generator.inst.in
tests/grid/grid_generator_general_cell_01.cc [moved from tests/grid/grid_generator_general_cell.cc with 100% similarity]
tests/grid/grid_generator_general_cell_01.output [moved from tests/grid/grid_generator_general_cell.output with 100% similarity]
tests/grid/grid_generator_general_cell_02.cc [new file with mode: 0644]
tests/grid/grid_generator_general_cell_02.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20190130GiovanniAlzetta b/doc/news/changes/minor/20190130GiovanniAlzetta
new file mode 100644 (file)
index 0000000..a9b6441
--- /dev/null
@@ -0,0 +1,5 @@
+Changed: Function GridGenerator::general_cell() can now generate
+a cell of dimension `dim` inside a space of dimension `spacedim`
+with `dim <= spacedim`
+<br>
+(Giovanni Alzetta, 2019/01/30)
index 354037d6f20ac9f909e472193388754e988852cd..0a4d52818f7d2d375e3a7b347fdb789efc3b24fd 100644 (file)
@@ -462,11 +462,13 @@ namespace GridGenerator
                         const bool          colorize           = false);
 
   /**
-   * 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.
+   * A general @p dim -dimensional cell (a segment if dim is 1, a quadrilateral
+   * if @p dim is 2, or a hexahedron if @p dim is 3) immersed in a
+   * @p spacedim -dimensional space. 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, then all boundary indicators are
    * set to zero for 2d and 3d. If it is true, the boundary is colorized as in
@@ -474,13 +476,17 @@ namespace GridGenerator
    * @ref GlossColorization "the glossary entry on colorization"). In 1d the
    * indicators are always colorized, see hyper_rectangle().
    *
+   * @param tria The triangulation that will be created
+   * @param vertices The 2^dim vertices of the cell
+   * @param colorize If true, set different boundary ids.
+   *
    * @author Bruno Turcksin
    */
-  template <int dim>
+  template <int dim, int spacedim>
   void
-  general_cell(Triangulation<dim> &           tria,
-               const std::vector<Point<dim>> &vertices,
-               const bool                     colorize = false);
+  general_cell(Triangulation<dim, spacedim> &      tria,
+               const std::vector<Point<spacedim>> &vertices,
+               const bool                          colorize = false);
 
   /**
    * A parallelogram. The first corner point is the origin. The @p dim
index a51cbe29a31239f7b02c27141db7dcc7514fa70a..9738d79372872203dfd3a6eccdfb84e73f1c7437 100644 (file)
@@ -749,11 +749,11 @@ namespace GridGenerator
 
 
 
-  template <int dim>
+  template <int dim, int spacedim>
   void
-  general_cell(Triangulation<dim> &           tria,
-               const std::vector<Point<dim>> &vertices,
-               const bool                     colorize)
+  general_cell(Triangulation<dim, spacedim> &      tria,
+               const std::vector<Point<spacedim>> &vertices,
+               const bool                          colorize)
   {
     Assert(vertices.size() == dealii::GeometryInfo<dim>::vertices_per_cell,
            ExcMessage("Wrong number of vertices."));
@@ -761,7 +761,7 @@ namespace GridGenerator
     // First create a hyper_rectangle and then deform it.
     hyper_cube(tria, 0, 1, colorize);
 
-    typename Triangulation<dim>::active_cell_iterator cell =
+    typename Triangulation<dim, spacedim>::active_cell_iterator cell =
       tria.begin_active();
     for (unsigned int i = 0; i < GeometryInfo<dim>::vertices_per_cell; ++i)
       cell->vertex(i) = vertices[i];
index 0084a285121079791286b0ae5da760438277f3fc..2253edee987a25c2af44745954618d2e69fbee9b 100644 (file)
@@ -106,6 +106,10 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
           &cells_to_remove,
         Triangulation<deal_II_dimension, deal_II_space_dimension> &result);
 
+      template void
+      general_cell(Triangulation<deal_II_dimension, deal_II_space_dimension> &,
+                   const std::vector<Point<deal_II_space_dimension>> &,
+                   const bool);
 #endif
     \}
   }
@@ -153,11 +157,6 @@ 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_02.cc b/tests/grid/grid_generator_general_cell_02.cc
new file mode 100644 (file)
index 0000000..cd60e7a
--- /dev/null
@@ -0,0 +1,95 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 - 2019 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+// Test output for GridGenerator::general_cell() for dim != spacedim
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/tria.h>
+
+#include "../tests.h"
+
+void
+dim_2_3(std::ostream &os)
+{
+  std::vector<Point<3>> vertices(4);
+  vertices[0](0) = -1.;
+  vertices[0](1) = -1.;
+  vertices[0](0) = 0.;
+
+  vertices[1](0) = 1.;
+  vertices[1](1) = -1.5;
+  vertices[1](2) = 0.;
+
+  vertices[2](0) = 1.5;
+  vertices[2](1) = 1.5;
+  vertices[2](2) = 0.;
+
+  vertices[3](0) = 2.;
+  vertices[3](1) = 0.5;
+  vertices[3](2) = 0.;
+
+  Triangulation<2, 3> tria;
+  GridGenerator::general_cell<2, 3>(tria, vertices);
+
+  GridOut gout;
+  gout.write_vtk(tria, os);
+}
+
+void
+dim_1_3(std::ostream &os)
+{
+  std::vector<Point<3>> vertices(2);
+  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;
+
+  Triangulation<1, 3> tria;
+  GridGenerator::general_cell<1, 3>(tria, vertices);
+
+  GridOut gout;
+  gout.write_vtk(tria, os);
+}
+
+void
+dim_1_2(std::ostream &os)
+{
+  std::vector<Point<2>> vertices(2);
+  vertices[0](0) = -1.;
+  vertices[0](1) = -1.;
+
+  vertices[1](0) = 1.;
+  vertices[1](1) = -1.5;
+
+  Triangulation<1, 2> tria;
+  GridGenerator::general_cell<1, 2>(tria, vertices);
+
+  GridOut gout;
+  gout.write_vtk(tria, os);
+}
+
+
+int
+main()
+{
+  initlog(true);
+  std::ostream &logfile = deallog.get_file_stream();
+  dim_2_3(logfile);
+  dim_1_3(logfile);
+  dim_1_2(logfile);
+}
diff --git a/tests/grid/grid_generator_general_cell_02.output b/tests/grid/grid_generator_general_cell_02.output
new file mode 100644 (file)
index 0000000..8a44d13
--- /dev/null
@@ -0,0 +1,87 @@
+
+# vtk DataFile Version 3.0
+Triangulation generated with deal.II
+ASCII
+DATASET UNSTRUCTURED_GRID
+POINTS 4 double
+0.00000 -1.00000 0.00000
+1.00000 -1.50000 0.00000
+1.50000 1.50000 0.00000
+2.00000 0.500000 0.00000
+
+CELLS 1 5
+4 0 1 3 2
+
+CELL_TYPES 1
+9 
+
+
+
+CELL_DATA 1
+SCALARS MaterialID int 1
+LOOKUP_TABLE default
+0 
+
+
+
+SCALARS ManifoldID int 1
+LOOKUP_TABLE default
+-1 
+
+
+# vtk DataFile Version 3.0
+Triangulation generated with deal.II
+ASCII
+DATASET UNSTRUCTURED_GRID
+POINTS 2 double
+-1.00000 -1.00000 -1.00000
+1.00000 -1.50000 -1.50000
+
+CELLS 1 3
+2 0 1
+
+CELL_TYPES 1
+3 
+
+
+
+CELL_DATA 1
+SCALARS MaterialID int 1
+LOOKUP_TABLE default
+0 
+
+
+
+SCALARS ManifoldID int 1
+LOOKUP_TABLE default
+-1 
+
+
+# vtk DataFile Version 3.0
+Triangulation generated with deal.II
+ASCII
+DATASET UNSTRUCTURED_GRID
+POINTS 2 double
+-1.00000 -1.00000 0
+1.00000 -1.50000 0
+
+CELLS 1 3
+2 0 1
+
+CELL_TYPES 1
+3 
+
+
+
+CELL_DATA 1
+SCALARS MaterialID int 1
+LOOKUP_TABLE default
+0 
+
+
+
+SCALARS ManifoldID int 1
+LOOKUP_TABLE default
+-1 
+
+

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.