* 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
*
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
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> &,
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> &,
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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);
+}
--- /dev/null
+
+# 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