From 44b4a358b529b60e50647b64d4489666468c9df1 Mon Sep 17 00:00:00 2001 From: kanschat Date: Wed, 6 Jun 2007 21:09:30 +0000 Subject: [PATCH] start hyper_shell in 3D git-svn-id: https://svn.dealii.org/trunk@14762 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/source/grid/grid_generator.cc | 86 +++++++++++++++++-- 1 file changed, 80 insertions(+), 6 deletions(-) diff --git a/deal.II/deal.II/source/grid/grid_generator.cc b/deal.II/deal.II/source/grid/grid_generator.cc index 0426ea62f0..0741e91c57 100644 --- a/deal.II/deal.II/source/grid/grid_generator.cc +++ b/deal.II/deal.II/source/grid/grid_generator.cc @@ -42,6 +42,39 @@ DEAL_II_NAMESPACE_OPEN +namespace +{ +#if deal_II_dimension == 3 + + // Corner points of the cube [-1,1]^3 + const Point<3> hexagon[8] = + { + Point<3>(-1,-1,-1), + Point<3>(+1,-1,-1), + Point<3>(-1,+1,-1), + Point<3>(+1,+1,-1), + Point<3>(-1,-1,+1), + Point<3>(+1,-1,+1), + Point<3>(-1,+1,+1), + Point<3>(+1,+1,+1) + }; + + // Octahedron inscribed in the cube + // [-1,1]^3 + const Point<3> octahedron[6] = + { + Point<3>(-1, 0, 0), + Point<3>( 1, 0, 0), + Point<3>( 0,-1, 0), + Point<3>( 0, 1, 0), + Point<3>( 0, 0,-1), + Point<3>( 0, 0, 1) + }; + +#endif +} + + template void GridGenerator::hyper_rectangle (Triangulation &tria, @@ -2039,13 +2072,54 @@ GridGenerator::half_hyper_ball (Triangulation& tria, // Implementation for 3D only template -void GridGenerator::hyper_shell (Triangulation&, - const Point&, - const double, - const double, - const unsigned int) +void GridGenerator::hyper_shell (Triangulation& tria, + const Point& p, + const double inner_radius, + const double outer_radius, + const unsigned int n) { - Assert (false, ExcNotImplemented()); + Assert ((inner_radius > 0) && (inner_radius < outer_radius), + ExcInvalidRadii ()); + + const double irad = inner_radius/std::sqrt(3.0); + const double orad = outer_radius/std::sqrt(3.0); + std::vector > vertices; + std::vector > cells; + + // Start with the shell bounded by + // two nested cubes + if (n <= 6) + { + for (unsigned int i=0;i<8;++i) + vertices.push_back(p+hexagon[i]*irad); + for (unsigned int i=0;i<8;++i) + vertices.push_back(p+hexagon[i]*orad); + // one needs to draw the seven cubes to + // understand what's going on here + const unsigned int n_cells = 6; + const int cell_vertices[n_cells][8] = {{8, 9, 10, 11, 0, 1, 2, 3}, // bottom + {9, 11, 1, 3, 13, 15, 5, 7}, // right + {12, 13, 4, 5, 14, 15, 6, 7}, // top + {8, 0, 10, 2, 12, 4, 14, 6}, // left + {8, 9, 0, 1, 12, 13, 4, 5}, // front + {10, 2, 11, 3, 14, 6, 15, 7}}; // back + + cells.resize(n_cells, CellData()); + + for (unsigned int i=0; i::vertices_per_cell; ++j) + cells[i].vertices[j] = cell_vertices[i][j]; + cells[i].material_id = 0; + } + } + else + { + Assert(false, ExcIndexRange(n, 1, 7)); + } + + tria.create_triangulation (vertices, cells, + SubCellData()); // no boundary information } -- 2.39.5