From: Bruno Turcksin Date: Thu, 22 Sep 2016 21:39:05 +0000 (-0400) Subject: Add GridGenerator for a general cell. X-Git-Tag: v8.5.0-rc1~618^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=d095fcf82b288c1b23fb0325ce169b026df14568;p=dealii.git Add GridGenerator for a general cell. --- diff --git a/include/deal.II/grid/grid_generator.h b/include/deal.II/grid/grid_generator.h index 267d57e8ba..77db7ec137 100644 --- a/include/deal.II/grid/grid_generator.h +++ b/include/deal.II/grid/grid_generator.h @@ -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 &tria, const std::vector &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 + void + general_cell(Triangulation &tria, + const std::vector > &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 diff --git a/source/grid/grid_generator.cc b/source/grid/grid_generator.cc index b16b56b59b..a4a7aa7ca8 100644 --- a/source/grid/grid_generator.cc +++ b/source/grid/grid_generator.cc @@ -749,6 +749,33 @@ namespace GridGenerator tria.set_all_manifold_ids_on_boundary(0); } + + + template + void + general_cell (Triangulation &tria, + const std::vector > &vertices, + const bool colorize) + { + Assert (vertices.size() == dealii::GeometryInfo::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::active_cell_iterator cell = tria.begin_active(); + for (unsigned int i=0; i::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> &, diff --git a/source/grid/grid_generator.inst.in b/source/grid/grid_generator.inst.in index 3d60241d21..f2af7dd898 100644 --- a/source/grid/grid_generator.inst.in +++ b/source/grid/grid_generator.inst.in @@ -122,6 +122,12 @@ for (deal_II_dimension : DIMENSIONS) const Point &, const bool ); + template + void + general_cell (Triangulation &, + const std::vector > &, + const bool); + template void simplex ( Triangulation &, diff --git a/tests/grid/grid_generator_general_cell.cc b/tests/grid/grid_generator_general_cell.cc new file mode 100644 index 0000000000..edbd294928 --- /dev/null +++ b/tests/grid/grid_generator_general_cell.cc @@ -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 +#include +#include + +void dim2(std::ostream &os) +{ + std::vector > 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 > 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 index 0000000000..edf1cc7562 --- /dev/null +++ b/tests/grid/grid_generator_general_cell.output @@ -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