From 99721eee6eceb94a5c3c51b3f24e21745b84f6af Mon Sep 17 00:00:00 2001 From: wolf Date: Mon, 7 Apr 2003 15:41:14 +0000 Subject: [PATCH] subdivided_hyper_cube. git-svn-id: https://svn.dealii.org/trunk@7365 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/grid/grid_generator.h | 34 ++++++- deal.II/deal.II/source/grid/grid_generator.cc | 91 +++++++++++++++++++ 2 files changed, 123 insertions(+), 2 deletions(-) diff --git a/deal.II/deal.II/include/grid/grid_generator.h b/deal.II/deal.II/include/grid/grid_generator.h index e7414a49d3..69af332390 100644 --- a/deal.II/deal.II/include/grid/grid_generator.h +++ b/deal.II/deal.II/include/grid/grid_generator.h @@ -38,7 +38,13 @@ template class SparseMatrix; * hypercube domain triangulated with exactly one element. You can * get tensor product meshes by successive refinement of this cell. * - * @item Rectangular coordinate parallel domains as a generalization + * If you want the hypercube subdivided a certain number of + * times (and if this is not achievable by hierarchic + * refinement, for example 3 times), then the + * @p{subdivided_hyper_cube} function may be what you are + * looking for. + * + * @item Rectangular coordinate-parallel domains as a generalization * of hypercubes are generated by * @ref{GridGenerator}@p{::hyper_rectangle (tria, p1, p2)}, with two * opposite corner points @p{p1} and @p{p2}. @@ -125,7 +131,7 @@ template class SparseMatrix; * transform (simple-shaped) grids to a more complicated ones, like a * shell onto a grid of an airfoil, for example. * - * @author Wolfgang Bangerth, Ralf Hartmann, Guido Kanschat, Stefan Nauber, 1998, 1999, 2000, 2001, 2002. + * @author Wolfgang Bangerth, Ralf Hartmann, Guido Kanschat, Stefan Nauber, 1998, 1999, 2000, 2001, 2002, 2003. */ class GridGenerator { @@ -154,6 +160,23 @@ class GridGenerator const double left = 0., const double right= 1.); + /** + * Same as @ref{hyper_cube}, but + * with the difference that not + * only one cell is created but + * each coordinate direction is + * subdivided into + * @p{repetitions} cells. Thus, + * the number of cells filling + * the given volume is + * @p{repetitions^d}. + */ + template + static void subdivided_hyper_cube (Triangulation &tria, + const unsigned int repetitions, + const double left = 0., + const double right= 1.); + /** * Create a coordinate-parallel * parallelepiped from the two @@ -561,6 +584,13 @@ class GridGenerator * Exception */ DeclException0 (ExcInvalidRadii); + /** + * Exception + */ + DeclException1 (ExcInvalidRepetitions, + int, + << "The number of repetitions " << arg1 + << " must be >=1."); private: /** diff --git a/deal.II/deal.II/source/grid/grid_generator.cc b/deal.II/deal.II/source/grid/grid_generator.cc index 9ce2b3521f..87144cebde 100644 --- a/deal.II/deal.II/source/grid/grid_generator.cc +++ b/deal.II/deal.II/source/grid/grid_generator.cc @@ -158,6 +158,90 @@ void GridGenerator::hyper_cube (Triangulation &tria, } +template +void +GridGenerator::subdivided_hyper_cube (Triangulation &tria, + const unsigned int repetitions, + const double left, + const double right) +{ + Assert (repetitions >= 1, ExcInvalidRepetitions(repetitions)); + + // first generate the necessary + // points + const double delta = (right-left)/repetitions; + std::vector > points; + switch (dim) + { + case 1: + for (unsigned int x=0; x<=repetitions; ++x) + points.push_back (Point (left+x*delta)); + break; + + case 2: + for (unsigned int y=0; y<=repetitions; ++y) + for (unsigned int x=0; x<=repetitions; ++x) + points.push_back (Point (left+x*delta, + left+y*delta)); + break; + + case 3: + for (unsigned int z=0; z<=repetitions; ++z) + for (unsigned int y=0; y<=repetitions; ++y) + for (unsigned int x=0; x<=repetitions; ++x) + points.push_back (Point (left+x*delta, + left+y*delta, + left+z*delta)); + break; + + default: + Assert (false, ExcNotImplemented()); + } + + // next create the cells + // Prepare cell data + std::vector > cells; + switch (dim) + { + case 1: + cells.resize (repetitions); + for (unsigned int x=0; x &, @@ -1014,6 +1098,13 @@ GridGenerator::hyper_cube (Triangulation & const double, const double); +template void +GridGenerator::subdivided_hyper_cube (Triangulation &, + const unsigned int, + const double, + const double); + + #if deal_II_dimension != 1 template void GridGenerator:: -- 2.39.5