From e333a1713efe3ad93aba69d39becab13e066dc8f Mon Sep 17 00:00:00 2001 From: "Toby D. Young" Date: Sun, 13 Jan 2013 00:08:43 +0000 Subject: [PATCH] Generate a parallelepiped grid in 1d, 2d, and 3d. git-svn-id: https://svn.dealii.org/trunk@28033 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/grid/grid_generator.h | 30 +++++++- deal.II/source/grid/grid_generator.cc | 77 ++++++++++++++++++- deal.II/source/grid/grid_generator.inst.in | 7 +- 3 files changed, 107 insertions(+), 7 deletions(-) diff --git a/deal.II/include/deal.II/grid/grid_generator.h b/deal.II/include/deal.II/grid/grid_generator.h index d9a15f9b5f..9f7b344ed2 100644 --- a/deal.II/include/deal.II/grid/grid_generator.h +++ b/deal.II/include/deal.II/grid/grid_generator.h @@ -303,11 +303,33 @@ public: * function. */ template - static void - parallelogram(Triangulation &tria, - const Tensor<2,dim> &corners, - const bool colorize=false); + static + void + parallelogram(Triangulation &tria, + const Tensor<2,dim> &corners, + const bool colorize=false); + /** + * A parallelepiped. The first corner point is the origin. The + * dim adjacent points are the rank one subtensors of the + * tensor provided and additional points will be sums of those + * dim vectors. Thus in 1d this is a line, in 2d this is a + * parallelogram, and in 3d a parallelepiped. Colorizing is done + * according to hyper_rectangle(). + * + * @note This function silently reorders the vertices on the cells + * to lexiographic ordering if they are not already ordered that way + * (see GridReordering::reorder_grid()). + * + * @note The triangulation needs to be void upon calling this + * function. + */ + template + static + void + parallelepiped (Triangulation &tria, + const Tensor<2, dim> &corners, + const bool colorize = false); /** * Hypercube with a layer of diff --git a/deal.II/source/grid/grid_generator.cc b/deal.II/source/grid/grid_generator.cc index b976ca453a..2d3e21fab2 100644 --- a/deal.II/source/grid/grid_generator.cc +++ b/deal.II/source/grid/grid_generator.cc @@ -373,9 +373,9 @@ GridGenerator::torus (Triangulation<2,3> &tria, template<> void GridGenerator::parallelogram ( - Triangulation<2> &tria, + Triangulation<2> &tria, const Tensor<2,2> &corners, - const bool colorize) + const bool colorize) { std::vector > vertices (GeometryInfo<2>::vertices_per_cell); @@ -397,6 +397,79 @@ GridGenerator::parallelogram ( +// Parallelepiped implementation in 1d, 2d, and 3d. @note The +// implementation in 1d is equivalent to hyper_rectangle, and in 2d is +// equivalent to parallelogram. Reordering of vertices to lexiographic +// occurs if needed. For this reason the triangulation is explicitly +// constructed for each dim here. +template +void +GridGenerator::parallelepiped (Triangulation &tria, + const Tensor<2, dim> &corners, + const bool colorize) +{ + // Check that none of the user defined vertices overlap + for (unsigned int i=0; i > vertices (GeometryInfo::vertices_per_cell); + + switch (dim) + { + // A line (1d parallelepiped) + case 1: + vertices[1] = corners[0]; + break; + + // A parallelogram (2d parallelepiped) + case 2: + // assign corners to vertices: + vertices[1] = corners[0]; + vertices[2] = corners[1]; + + // compose the remaining vertex: + vertices[3] = vertices[1] + vertices[2]; + break; + + // A parallelepiped (3d parallelepiped) + case 3: + // assign corners to vertices: + vertices[1] = corners[0]; + vertices[2] = corners[1]; + vertices[4] = corners[2]; + + // compose the remaining vertices: + vertices[3] = vertices[1] + vertices[2]; + vertices[5] = vertices[1] + vertices[4]; + vertices[6] = vertices[2] + vertices[4]; + vertices[7] = vertices[1] + vertices[2] + vertices[4]; + break; + + default: + Assert (false, ExcNotImplemented()); + } + + // Prepare cell data and wipe material identity + std::vector > cells (1); + for (unsigned int i=0; i::vertices_per_cell; ++i) + cells[0].vertices[i] = i; + cells[0].material_id = 0; + + // Check ordering of vertices and create triangulation + GridReordering::reorder_cells (cells); + tria.create_triangulation (vertices, cells, SubCellData()); + + // Finally assign boundary indicatorsaccording to hyper_rectangle + if (colorize) + colorize_hyper_rectangle (tria); +} + + + + template void GridGenerator::subdivided_hyper_cube (Triangulation &tria, diff --git a/deal.II/source/grid/grid_generator.inst.in b/deal.II/source/grid/grid_generator.inst.in index f16f7005fd..0561f0a4a4 100644 --- a/deal.II/source/grid/grid_generator.inst.in +++ b/deal.II/source/grid/grid_generator.inst.in @@ -64,7 +64,12 @@ for (deal_II_dimension : DIMENSIONS) const Point &, const bool ); - + + template void + GridGenerator::parallelepiped ( + Triangulation &, + const Tensor<2, deal_II_dimension> &, + const bool); #if deal_II_dimension > 1 -- 2.39.5