From de5b04cb9b40ec472133741d07ab38241e0405e7 Mon Sep 17 00:00:00 2001 From: Guido Kanschat Date: Wed, 28 Dec 2005 10:58:55 +0000 Subject: [PATCH] reverted to templates in order to keep documentation readable git-svn-id: https://svn.dealii.org/trunk@11930 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/grid/grid_generator.h | 222 ++------ deal.II/deal.II/source/grid/grid_generator.cc | 536 ++++++++++++------ 2 files changed, 395 insertions(+), 363 deletions(-) diff --git a/deal.II/deal.II/include/grid/grid_generator.h b/deal.II/deal.II/include/grid/grid_generator.h index 84a47c5439..6b95d3d88c 100644 --- a/deal.II/deal.II/include/grid/grid_generator.h +++ b/deal.II/deal.II/include/grid/grid_generator.h @@ -301,26 +301,13 @@ class GridGenerator * void upon calling this * function. */ - static void enclosed_hyper_cube (Triangulation<2> &tria, - const double left = 0., - const double right= 1., - const double thickness = 1., - const bool colorize = false); - - /** - * Declaration of same function - * for different space dimension. - * - * @note The triangulation needs to be - * void upon calling this - * function. - */ - static void enclosed_hyper_cube (Triangulation<3> &tria, + template + static void enclosed_hyper_cube (Triangulation &tria, const double left = 0., const double right= 1., const double thickness = 1., const bool colorize = false); - + /** * Initialize the given * triangulation with a @@ -346,32 +333,9 @@ class GridGenerator * void upon calling this * function. */ - static void hyper_ball (Triangulation<1> &tria, - const Point<1> ¢er = Point<1>(), - const double radius = 1.); - - /** - * Declaration of same function - * for different space dimension. - * - * @note The triangulation needs to be - * void upon calling this - * function. - */ - static void hyper_ball (Triangulation<2> &tria, - const Point<2> ¢er = Point<2>(), - const double radius = 1.); - - /** - * Declaration of same function - * for different space dimension. - * - * @note The triangulation needs to be - * void upon calling this - * function. - */ - static void hyper_ball (Triangulation<3> &tria, - const Point<3> ¢er = Point<3>(), + template + static void hyper_ball (Triangulation &tria, + const Point ¢er = Point(), const double radius = 1.); /** @@ -383,6 +347,13 @@ class GridGenerator * @p yz-plane is a circle of * radius @p radius. * + * In two dimensions, the + * cylinder is a rectangle from + * x=-half_length to + * x=+half_length and + * from y=-radius to + * y=radius. + * * The boundaries are colored * according to the following * scheme: 0 for the hull of the @@ -394,40 +365,11 @@ class GridGenerator * void upon calling this * function. */ - static void cylinder (Triangulation<3> &tria, - const double radius = 1., - const double half_length = 1.); - - /** - * Projection of the - * three-dimensional cylinder - * into the @p xy-plane. - * Therefore, this is simply a square. - * - * Coloring is like in 3D. - * - * @note The triangulation needs to be - * void upon calling this - * function. - */ - static void cylinder (Triangulation<2> &tria, + template + static void cylinder (Triangulation &tria, const double radius = 1., const double half_length = 1.); - - /** - * Non-implemented dummy for compilation - * purposes. - * - * @note The triangulation needs to be - * void upon calling this - * function. - */ - static void cylinder (Triangulation<1> &tria, - const double radius, - const double half_length); - - - + /** * Initialize the given * triangulation with a hyper-L @@ -452,31 +394,8 @@ class GridGenerator * void upon calling this * function. */ - static void hyper_L (Triangulation<1> &tria, - const double left = -1., - const double right= 1.); - - /** - * Declaration of same function - * for different space dimension. - * - * @note The triangulation needs to be - * void upon calling this - * function. - */ - static void hyper_L (Triangulation<2> &tria, - const double left = -1., - const double right= 1.); - - /** - * Declaration of same function - * for different space dimension. - * - * @note The triangulation needs to be - * void upon calling this - * function. - */ - static void hyper_L (Triangulation<3> &tria, + template + static void hyper_L (Triangulation &tria, const double left = -1., const double right= 1.); @@ -510,46 +429,12 @@ class GridGenerator * void upon calling this * function. */ - static void hyper_cube_slit (Triangulation<1> &tria, - const double left = 0., - const double right= 1., - const bool colorize = false); - - /** - * Declaration of same function - * for different space dimension. - * - * If colorize is selected, then - * the two edges forming the slit - * carry numbers 1 and 2, while - * the outer boundary has number - * 0. - * - * @note The triangulation needs to be - * void upon calling this - * function. - */ - static void hyper_cube_slit (Triangulation<2> &tria, - const double left = 0., - const double right= 1., - const bool colorize = false); - - /** - * Declaration of same function - * for different space dimension. - * - * Colorization is not - * implemented in 3D. - * - * @note The triangulation needs to be - * void upon calling this - * function. - */ - static void hyper_cube_slit (Triangulation<3> &tria, + template + static void hyper_cube_slit (Triangulation &tria, const double left = 0., const double right= 1., const bool colorize = false); - + /** * Produce a hyper-shell, * i.e. the space between two @@ -581,40 +466,13 @@ class GridGenerator * void upon calling this * function. */ - static void hyper_shell (Triangulation<1> &tria, - const Point<1> ¢er, - const double inner_radius, - const double outer_radius, - const unsigned int n_cells = 0); - - /** - * Declaration of same function - * for different space dimension. - * - * @note The triangulation needs to be - * void upon calling this - * function. - */ - static void hyper_shell (Triangulation<2> &tria, - const Point<2> ¢er, - const double inner_radius, - const double outer_radius, - const unsigned int n_cells = 0); - - /** - * Declaration of same function - * for different space dimension. - * - * @note The triangulation needs to be - * void upon calling this - * function. - */ - static void hyper_shell (Triangulation<3> &tria, - const Point<3> ¢er, + template + static void hyper_shell (Triangulation &tria, + const Point ¢er, const double inner_radius, const double outer_radius, const unsigned int n_cells = 0); - + /** * Produce a domain that is the space * between two cylinders in 3d, with @@ -628,14 +486,17 @@ class GridGenerator * ratio. The same holds for @p * n_axial_cells. * - * The triangulation needs to be void - * upon calling this function. + * @note Although this function + * is declared as a template, it + * does not make sense in 1D and + * 2D. * - * @note The triangulation needs to be - * void upon calling this + * @note The triangulation needs + * to be void upon calling this * function. */ - static void cylinder_shell (Triangulation<3> &tria, + template + static void cylinder_shell (Triangulation &tria, const double length, const double inner_radius, const double outer_radius, @@ -659,8 +520,9 @@ class GridGenerator * void upon calling this * function. */ - static void half_hyper_ball (Triangulation<2> &tria, - const Point<2> ¢er = Point<2>(), + template + static void half_hyper_ball (Triangulation &tria, + const Point ¢er = Point(), const double radius = 1.); /** @@ -705,8 +567,9 @@ class GridGenerator * void upon calling this * function. */ - static void half_hyper_shell (Triangulation<2> &tria, - const Point<2> ¢er, + template + static void half_hyper_shell (Triangulation &tria, + const Point ¢er, const double inner_radius, const double outer_radius, const unsigned int n_cells = 0); @@ -738,13 +601,6 @@ class GridGenerator static void laplace_transformation (Triangulation &tria, const std::map > &new_points); - /** - * Declaration of same function - * for different space dimension. - */ - static void laplace_transformation (Triangulation<1> &tria, - const std::map > &new_points); - /** * Exception */ diff --git a/deal.II/deal.II/source/grid/grid_generator.cc b/deal.II/deal.II/source/grid/grid_generator.cc index 2d8a9a6697..13d080b579 100644 --- a/deal.II/deal.II/source/grid/grid_generator.cc +++ b/deal.II/deal.II/source/grid/grid_generator.cc @@ -101,6 +101,7 @@ GridGenerator::hyper_rectangle (Triangulation &tria, #if deal_II_dimension == 1 +// Implementation for 1D only void GridGenerator::colorize_hyper_rectangle (Triangulation<1> &) { @@ -110,6 +111,7 @@ GridGenerator::colorize_hyper_rectangle (Triangulation<1> &) #else +// Implementation for dimensions except 1 template void GridGenerator::colorize_hyper_rectangle (Triangulation &tria) @@ -140,6 +142,7 @@ void GridGenerator::hyper_cube (Triangulation &tria, } +// Implementation for 2D only template void GridGenerator::parallelogram ( @@ -407,6 +410,7 @@ subdivided_hyper_rectangle (Triangulation &tria, #if deal_II_dimension == 1 +// Implementation for 1D only void GridGenerator::colorize_subdivided_hyper_rectangle (Triangulation<1> &, const Point<1> &, @@ -420,6 +424,7 @@ GridGenerator::colorize_subdivided_hyper_rectangle (Triangulation<1> &, #else +// Implementation for dimensions except 1 template void GridGenerator::colorize_subdivided_hyper_rectangle (Triangulation &tria, @@ -471,7 +476,10 @@ GridGenerator::colorize_subdivided_hyper_rectangle (Triangulation &tria, #if deal_II_dimension == 1 -void GridGenerator::hyper_cube_slit (Triangulation<1> &, +// The following implementations are for 1D only. See below for other +// dimensions. +template +void GridGenerator::hyper_cube_slit (Triangulation &, const double, const double, const bool) @@ -481,7 +489,22 @@ void GridGenerator::hyper_cube_slit (Triangulation<1> &, -void GridGenerator::hyper_L (Triangulation<1> &, +// Implementation for 1D only +template +void GridGenerator::enclosed_hyper_cube (Triangulation&, + const double, + const double, + const double, + const bool) +{ + Assert (false, ExcNotImplemented()); +} + + + +// Implementation for 1D only +template +void GridGenerator::hyper_L (Triangulation &, const double, const double) { @@ -490,8 +513,10 @@ void GridGenerator::hyper_L (Triangulation<1> &, -void GridGenerator::hyper_ball (Triangulation<1> &, - const Point<1> &, +// Implementation for 1D only +template +void GridGenerator::hyper_ball (Triangulation &, + const Point &, const double) { Assert (false, ExcNotImplemented()); @@ -499,7 +524,9 @@ void GridGenerator::hyper_ball (Triangulation<1> &, -void GridGenerator::cylinder (Triangulation<1> &, +// Implementation for 1D only +template +void GridGenerator::cylinder (Triangulation &, const double, const double) { @@ -508,8 +535,46 @@ void GridGenerator::cylinder (Triangulation<1> &, -void GridGenerator::hyper_shell (Triangulation<1> &, - const Point<1> &, +// Implementation for 1D only +template +void GridGenerator::hyper_shell (Triangulation &, + const Point &, + const double, + const double, + const unsigned int) +{ + Assert (false, ExcNotImplemented()); +} + +// Implementation for 1D only +template +void GridGenerator::cylinder_shell (Triangulation&, + const double, + const double, + const double, + const unsigned int, + const unsigned int) +{ + Assert (false, ExcNotImplemented()); +} + + +// Implementation for 1D only +template +void +GridGenerator::half_hyper_ball (Triangulation&, + const Point&, + const double) +{ + Assert (false, ExcNotImplemented()); +} + + +// Implementation for 1D only +template +void +GridGenerator::half_hyper_shell (Triangulation&, + const Point&, const double, const double, const unsigned int) @@ -523,13 +588,15 @@ void GridGenerator::hyper_shell (Triangulation<1> &, #if deal_II_dimension == 2 -void GridGenerator::enclosed_hyper_cube (Triangulation<2> &tria, +// Implementation for 2D only +template +void GridGenerator::enclosed_hyper_cube (Triangulation &tria, const double l, const double r, const double d, const bool colorize) { - std::vector > vertices(16); + std::vector > vertices(16); double coords[4]; coords[0] = l-d; coords[1] = l; @@ -539,14 +606,14 @@ void GridGenerator::enclosed_hyper_cube (Triangulation<2> &tria, unsigned int k=0; for (unsigned int i0=0;i0<4;++i0) for (unsigned int i1=0;i1<4;++i1) - vertices[k++] = Point<2>(coords[i1], coords[i0]); + vertices[k++] = Point(coords[i1], coords[i0]); const unsigned char materials[9] = { 5, 4, 6, 1, 0, 2, 9, 8,10 }; - std::vector > cells(9); + std::vector > cells(9); k = 0; for (unsigned int i0=0;i0<3;++i0) for (unsigned int i1=0;i1<3;++i1) @@ -566,28 +633,30 @@ void GridGenerator::enclosed_hyper_cube (Triangulation<2> &tria, +// Implementation for 2D only +template void -GridGenerator::hyper_cube_slit (Triangulation<2> &tria, +GridGenerator::hyper_cube_slit (Triangulation &tria, const double left, const double right, const bool colorize) { const double rl2=(right+left)/2; - const Point<2> vertices[10] = { Point<2>(left, left ), - Point<2>(rl2, left ), - Point<2>(rl2, rl2 ), - Point<2>(left, rl2 ), - Point<2>(right,left ), - Point<2>(right,rl2 ), - Point<2>(rl2, right), - Point<2>(left, right), - Point<2>(right,right), - Point<2>(rl2, left ) }; + const Point vertices[10] = { Point(left, left ), + Point(rl2, left ), + Point(rl2, rl2 ), + Point(left, rl2 ), + Point(right,left ), + Point(right,rl2 ), + Point(rl2, right), + Point(left, right), + Point(right,right), + Point(rl2, left ) }; const int cell_vertices[4][4] = { { 0,1,3,2 }, { 9,4,2,5 }, { 3,2,7,6 }, { 2,5,6,8 } }; - std::vector > cells (4, CellData<2>()); + std::vector > cells (4, CellData()); for (unsigned int i=0; i<4; ++i) { for (unsigned int j=0; j<4; ++j) @@ -595,13 +664,13 @@ GridGenerator::hyper_cube_slit (Triangulation<2> &tria, cells[i].material_id = 0; }; tria.create_triangulation ( - std::vector >(&vertices[0], &vertices[10]), + std::vector >(&vertices[0], &vertices[10]), cells, SubCellData()); // no boundary information if (colorize) { - Triangulation<2>::cell_iterator cell = tria.begin(); + typename Triangulation::cell_iterator cell = tria.begin(); cell->face(1)->set_boundary_indicator(1); ++cell; cell->face(3)->set_boundary_indicator(2); @@ -610,12 +679,13 @@ GridGenerator::hyper_cube_slit (Triangulation<2> &tria, //TODO: Colorize edges as circumference, left and right radius +// Implementation for 2D only +template void -GridGenerator::hyper_L (Triangulation<2> &tria, +GridGenerator::hyper_L (Triangulation &tria, const double a, const double b) { - const unsigned int dim=2; const Point vertices[8] = { Point (a,a), Point ((a+b)/2,a), Point (b,a), @@ -628,7 +698,7 @@ GridGenerator::hyper_L (Triangulation<2> &tria, {1, 2, 4, 5}, {3, 4, 6, 7}}; - std::vector > cells (3, CellData<2>()); + std::vector > cells (3, CellData()); for (unsigned int i=0; i<3; ++i) { @@ -645,23 +715,25 @@ GridGenerator::hyper_L (Triangulation<2> &tria, +// Implementation for 2D only +template void -GridGenerator::hyper_ball (Triangulation<2> &tria, - const Point<2> &p, +GridGenerator::hyper_ball (Triangulation &tria, + const Point &p, const double radius) { // equilibrate cell sizes at // transition from the inner part // to the radial cells const double a = 1./(1+std::sqrt(2.0)); - const Point<2> vertices[8] = { p+Point<2>(-1,-1)*(radius/std::sqrt(2.0)), - p+Point<2>(+1,-1)*(radius/std::sqrt(2.0)), - p+Point<2>(-1,-1)*(radius/std::sqrt(2.0)*a), - p+Point<2>(+1,-1)*(radius/std::sqrt(2.0)*a), - p+Point<2>(-1,+1)*(radius/std::sqrt(2.0)*a), - p+Point<2>(+1,+1)*(radius/std::sqrt(2.0)*a), - p+Point<2>(-1,+1)*(radius/std::sqrt(2.0)), - p+Point<2>(+1,+1)*(radius/std::sqrt(2.0)) }; + const Point vertices[8] = { p+Point(-1,-1)*(radius/std::sqrt(2.0)), + p+Point(+1,-1)*(radius/std::sqrt(2.0)), + p+Point(-1,-1)*(radius/std::sqrt(2.0)*a), + p+Point(+1,-1)*(radius/std::sqrt(2.0)*a), + p+Point(-1,+1)*(radius/std::sqrt(2.0)*a), + p+Point(+1,+1)*(radius/std::sqrt(2.0)*a), + p+Point(-1,+1)*(radius/std::sqrt(2.0)), + p+Point(+1,+1)*(radius/std::sqrt(2.0)) }; const int cell_vertices[5][4] = {{0, 1, 2, 3}, {0, 2, 6, 4}, @@ -669,7 +741,7 @@ GridGenerator::hyper_ball (Triangulation<2> &tria, {1, 7, 3, 5}, {6, 4, 7, 5}}; - std::vector > cells (5, CellData<2>()); + std::vector > cells (5, CellData()); for (unsigned int i=0; i<5; ++i) { @@ -679,15 +751,17 @@ GridGenerator::hyper_ball (Triangulation<2> &tria, }; tria.create_triangulation ( - std::vector >(&vertices[0], &vertices[8]), + std::vector >(&vertices[0], &vertices[8]), cells, SubCellData()); // no boundary information } -void GridGenerator::hyper_shell (Triangulation<2> &tria, - const Point<2> ¢er, +// Implementation for 2D only +template +void GridGenerator::hyper_shell (Triangulation &tria, + const Point ¢er, const double inner_radius, const double outer_radius, const unsigned int n_cells) @@ -718,10 +792,10 @@ void GridGenerator::hyper_shell (Triangulation<2> &tria, // first N ones are on the // outer one, and all are // numbered counter-clockwise - std::vector > vertices(2*N); + std::vector > vertices(2*N); for (unsigned int i=0; i( std::cos(2*pi * i/N), + vertices[i] = Point( std::cos(2*pi * i/N), std::sin(2*pi * i/N)) * outer_radius; vertices[i+N] = vertices[i] * (inner_radius/outer_radius); @@ -729,7 +803,7 @@ void GridGenerator::hyper_shell (Triangulation<2> &tria, vertices[i+N] += center; }; - std::vector > cells (N, CellData<2>()); + std::vector > cells (N, CellData()); for (unsigned int i=0; i &tria, +// Implementation for 2D only +template void -GridGenerator::cylinder (Triangulation<2> &tria, +GridGenerator::cylinder (Triangulation &tria, const double radius, const double half_length) { - Point<2> p1 (-half_length, -radius); - Point<2> p2 (half_length, radius); + Point p1 (-half_length, -radius); + Point p2 (half_length, radius); hyper_rectangle(tria, p1, p2, true); - Triangulation<2>::face_iterator f = tria.begin_face(); - Triangulation<2>::face_iterator end = tria.end_face(); + typename Triangulation::face_iterator f = tria.begin_face(); + typename Triangulation::face_iterator end = tria.end_face(); while (f != end) { switch (f->boundary_indicator()) @@ -779,31 +855,46 @@ GridGenerator::cylinder (Triangulation<2> &tria, +// Implementation for 2D only +template +void GridGenerator::cylinder_shell (Triangulation&, + const double, + const double, + const double, + const unsigned int, + const unsigned int) +{ + Assert (false, ExcNotImplemented()); +} + + //TODO: Colorize edges as circumference and cut plane +// Implementation for 2D only +template void -GridGenerator::half_hyper_ball (Triangulation<2> &tria, - const Point<2> &p, +GridGenerator::half_hyper_ball (Triangulation &tria, + const Point &p, const double radius) { // equilibrate cell sizes at // transition from the inner part // to the radial cells const double a = 1./(1+std::sqrt(2.0)); - const Point<2> vertices[8] = { p+Point<2>(0,-1)*radius, - p+Point<2>(+1,-1)*(radius/std::sqrt(2.0)), - p+Point<2>(0,-1)*(radius/std::sqrt(2.0)*a), - p+Point<2>(+1,-1)*(radius/std::sqrt(2.0)*a), - p+Point<2>(0,+1)*(radius/std::sqrt(2.0)*a), - p+Point<2>(+1,+1)*(radius/std::sqrt(2.0)*a), - p+Point<2>(0,+1)*radius, - p+Point<2>(+1,+1)*(radius/std::sqrt(2.0)) }; + const Point vertices[8] = { p+Point(0,-1)*radius, + p+Point(+1,-1)*(radius/std::sqrt(2.0)), + p+Point(0,-1)*(radius/std::sqrt(2.0)*a), + p+Point(+1,-1)*(radius/std::sqrt(2.0)*a), + p+Point(0,+1)*(radius/std::sqrt(2.0)*a), + p+Point(+1,+1)*(radius/std::sqrt(2.0)*a), + p+Point(0,+1)*radius, + p+Point(+1,+1)*(radius/std::sqrt(2.0)) }; const int cell_vertices[5][4] = {{0, 1, 2, 3}, {2, 3, 4, 5}, {1, 7, 3, 5}, {6, 4, 7, 5}}; - std::vector > cells (4, CellData<2>()); + std::vector > cells (4, CellData()); for (unsigned int i=0; i<4; ++i) { @@ -813,16 +904,18 @@ GridGenerator::half_hyper_ball (Triangulation<2> &tria, }; tria.create_triangulation ( - std::vector >(&vertices[0], &vertices[8]), + std::vector >(&vertices[0], &vertices[8]), cells, SubCellData()); // no boundary information } +// Implementation for 2D only +template void -GridGenerator::half_hyper_shell (Triangulation<2> &tria, - const Point<2> ¢er, +GridGenerator::half_hyper_shell (Triangulation &tria, + const Point ¢er, const double inner_radius, const double outer_radius, const unsigned int n_cells) @@ -852,7 +945,7 @@ GridGenerator::half_hyper_shell (Triangulation<2> &tria, // first N+1 ones are on the // outer one, and all are // numbered counter-clockwise - std::vector > vertices(2*(N+1)); + std::vector > vertices(2*(N+1)); for (unsigned int i=0; i<=N; ++i) { // enforce that the x-coordinates @@ -861,7 +954,7 @@ GridGenerator::half_hyper_shell (Triangulation<2> &tria, // zero (contrary to what we may // compute using the imprecise // value of pi) - vertices[i] = Point<2>( ( (i==0) || (i==N) ? + vertices[i] = Point( ( (i==0) || (i==N) ? 0 : std::cos(pi * i/N - pi/2) ), std::sin(pi * i/N - pi/2)) * outer_radius; @@ -872,7 +965,7 @@ GridGenerator::half_hyper_shell (Triangulation<2> &tria, }; - std::vector > cells (N, CellData<2>()); + std::vector > cells (N, CellData()); for (unsigned int i=0; i &tria, #if deal_II_dimension == 3 -void GridGenerator::hyper_cube_slit (Triangulation<3> &, +// Implementation for 3D only +template +void GridGenerator::hyper_cube_slit (Triangulation &, const double, const double, const bool) @@ -905,13 +1000,15 @@ void GridGenerator::hyper_cube_slit (Triangulation<3> &, -void GridGenerator::enclosed_hyper_cube (Triangulation<3> &tria, +// Implementation for 3D only +template +void GridGenerator::enclosed_hyper_cube (Triangulation &tria, const double l, const double r, const double d, const bool colorize) { - std::vector > vertices(64); + std::vector > vertices(64); double coords[4]; coords[0] = l-d; coords[1] = l; @@ -922,7 +1019,7 @@ void GridGenerator::enclosed_hyper_cube (Triangulation<3> &tria, for (unsigned int i0=0;i0<4;++i0) for (unsigned int i1=0;i1<4;++i1) for (unsigned int i2=0;i2<4;++i2) - vertices[k++] = Point<3>(coords[i2], coords[i1], coords[i0]); + vertices[k++] = Point(coords[i2], coords[i1], coords[i0]); const unsigned char materials[27] = { 21,20,22, @@ -936,7 +1033,7 @@ void GridGenerator::enclosed_hyper_cube (Triangulation<3> &tria, 41,40,42 }; - std::vector > cells(27); + std::vector > cells(27); k = 0; for (unsigned int i0=0;i0<3;++i0) for (unsigned int i1=0;i1<3;++i1) @@ -962,12 +1059,13 @@ void GridGenerator::enclosed_hyper_cube (Triangulation<3> &tria, +// Implementation for 3D only +template void -GridGenerator::hyper_L (Triangulation<3> &tria, +GridGenerator::hyper_L (Triangulation &tria, const double a, const double b) { - const unsigned int dim=3; // we slice out the top back right // part of the cube const Point vertices[26] @@ -1011,7 +1109,7 @@ GridGenerator::hyper_L (Triangulation<3> &tria, {10, 11, 19, 20, 13, 14, 22, 23}, {12, 13, 21, 22, 15, 16, 24, 25}}; - std::vector > cells (7, CellData<3>()); + std::vector > cells (7, CellData()); for (unsigned int i=0; i<7; ++i) { @@ -1028,37 +1126,39 @@ GridGenerator::hyper_L (Triangulation<3> &tria, +// Implementation for 3D only +template void -GridGenerator::hyper_ball (Triangulation<3> &tria, - const Point<3> &p, +GridGenerator::hyper_ball (Triangulation &tria, + const Point &p, const double radius) { const double a = 1./(1+std::sqrt(3.0)); // equilibrate cell sizes at transition // from the inner part to the radial // cells const unsigned int n_vertices = 16; - const Point<3> vertices[n_vertices] + const Point vertices[n_vertices] = { // first the vertices of the inner // cell - p+Point<3>(-1,-1,-1)*(radius/std::sqrt(3.0)*a), - p+Point<3>(+1,-1,-1)*(radius/std::sqrt(3.0)*a), - p+Point<3>(+1,-1,+1)*(radius/std::sqrt(3.0)*a), - p+Point<3>(-1,-1,+1)*(radius/std::sqrt(3.0)*a), - p+Point<3>(-1,+1,-1)*(radius/std::sqrt(3.0)*a), - p+Point<3>(+1,+1,-1)*(radius/std::sqrt(3.0)*a), - p+Point<3>(+1,+1,+1)*(radius/std::sqrt(3.0)*a), - p+Point<3>(-1,+1,+1)*(radius/std::sqrt(3.0)*a), + p+Point(-1,-1,-1)*(radius/std::sqrt(3.0)*a), + p+Point(+1,-1,-1)*(radius/std::sqrt(3.0)*a), + p+Point(+1,-1,+1)*(radius/std::sqrt(3.0)*a), + p+Point(-1,-1,+1)*(radius/std::sqrt(3.0)*a), + p+Point(-1,+1,-1)*(radius/std::sqrt(3.0)*a), + p+Point(+1,+1,-1)*(radius/std::sqrt(3.0)*a), + p+Point(+1,+1,+1)*(radius/std::sqrt(3.0)*a), + p+Point(-1,+1,+1)*(radius/std::sqrt(3.0)*a), // now the eight vertices at // the outer sphere - p+Point<3>(-1,-1,-1)*(radius/std::sqrt(3.0)), - p+Point<3>(+1,-1,-1)*(radius/std::sqrt(3.0)), - p+Point<3>(+1,-1,+1)*(radius/std::sqrt(3.0)), - p+Point<3>(-1,-1,+1)*(radius/std::sqrt(3.0)), - p+Point<3>(-1,+1,-1)*(radius/std::sqrt(3.0)), - p+Point<3>(+1,+1,-1)*(radius/std::sqrt(3.0)), - p+Point<3>(+1,+1,+1)*(radius/std::sqrt(3.0)), - p+Point<3>(-1,+1,+1)*(radius/std::sqrt(3.0)), + p+Point(-1,-1,-1)*(radius/std::sqrt(3.0)), + p+Point(+1,-1,-1)*(radius/std::sqrt(3.0)), + p+Point(+1,-1,+1)*(radius/std::sqrt(3.0)), + p+Point(-1,-1,+1)*(radius/std::sqrt(3.0)), + p+Point(-1,+1,-1)*(radius/std::sqrt(3.0)), + p+Point(+1,+1,-1)*(radius/std::sqrt(3.0)), + p+Point(+1,+1,+1)*(radius/std::sqrt(3.0)), + p+Point(-1,+1,+1)*(radius/std::sqrt(3.0)), }; // one needs to draw the seven cubes to @@ -1072,57 +1172,59 @@ GridGenerator::hyper_ball (Triangulation<3> &tria, {8, 9, 0, 1, 11, 10, 3, 2}, // front {12, 4, 13, 5, 15, 7, 14, 6}}; // back - std::vector > cells (n_cells, CellData<3>()); + std::vector > cells (n_cells, CellData()); for (unsigned int i=0; i::vertices_per_cell; ++j) + for (unsigned int j=0; j::vertices_per_cell; ++j) cells[i].vertices[j] = cell_vertices[i][j]; cells[i].material_id = 0; }; tria.create_triangulation ( - std::vector >(&vertices[0], &vertices[n_vertices]), + std::vector >(&vertices[0], &vertices[n_vertices]), cells, SubCellData()); // no boundary information } +// Implementation for 3D only +template void -GridGenerator::cylinder (Triangulation<3> &tria, +GridGenerator::cylinder (Triangulation &tria, const double radius, const double half_length) { - // Copy the base from hyper_ball<2> + // Copy the base from hyper_ball // and transform it to yz const double d = radius/std::sqrt(2.0); const double a = d/(1+std::sqrt(2.0)); - Point<3> vertices[24] = { - Point<3>(-d, -half_length,-d), - Point<3>( d, -half_length,-d), - Point<3>(-a, -half_length,-a), - Point<3>( a, -half_length,-a), - Point<3>(-a, -half_length, a), - Point<3>( a, -half_length, a), - Point<3>(-d, -half_length, d), - Point<3>( d, -half_length, d), - Point<3>(-d, 0,-d), - Point<3>( d, 0,-d), - Point<3>(-a, 0,-a), - Point<3>( a, 0,-a), - Point<3>(-a, 0, a), - Point<3>( a, 0, a), - Point<3>(-d, 0, d), - Point<3>( d, 0, d), - Point<3>(-d, half_length,-d), - Point<3>( d, half_length,-d), - Point<3>(-a, half_length,-a), - Point<3>( a, half_length,-a), - Point<3>(-a, half_length, a), - Point<3>( a, half_length, a), - Point<3>(-d, half_length, d), - Point<3>( d, half_length, d), + Point vertices[24] = { + Point(-d, -half_length,-d), + Point( d, -half_length,-d), + Point(-a, -half_length,-a), + Point( a, -half_length,-a), + Point(-a, -half_length, a), + Point( a, -half_length, a), + Point(-d, -half_length, d), + Point( d, -half_length, d), + Point(-d, 0,-d), + Point( d, 0,-d), + Point(-a, 0,-a), + Point( a, 0,-a), + Point(-a, 0, a), + Point( a, 0, a), + Point(-d, 0, d), + Point( d, 0, d), + Point(-d, half_length,-d), + Point( d, half_length,-d), + Point(-a, half_length,-a), + Point( a, half_length,-a), + Point(-a, half_length, a), + Point( a, half_length, a), + Point(-d, half_length, d), + Point( d, half_length, d), }; // Turn cylinder such that y->x for (unsigned int i=0;i<24;++i) @@ -1143,7 +1245,7 @@ GridGenerator::cylinder (Triangulation<3> &tria, for (unsigned int j=0;j<8;++j) cell_vertices[i+5][j] = cell_vertices[i][j]+8; - std::vector > cells (10, CellData<3>()); + std::vector > cells (10, CellData()); for (unsigned int i=0; i<10; ++i) { @@ -1153,16 +1255,16 @@ GridGenerator::cylinder (Triangulation<3> &tria, }; tria.create_triangulation ( - std::vector >(&vertices[0], &vertices[24]), + std::vector >(&vertices[0], &vertices[24]), cells, SubCellData()); // no boundary information - Triangulation<3>::cell_iterator cell = tria.begin(); - Triangulation<3>::cell_iterator end = tria.end(); + typename Triangulation::cell_iterator cell = tria.begin(); + typename Triangulation::cell_iterator end = tria.end(); while (cell != end) { - for (unsigned int i=0;i::faces_per_cell;++i) + for (unsigned int i=0;i::faces_per_cell;++i) { if (cell->face(i)->boundary_indicator() == 255) continue; @@ -1178,18 +1280,47 @@ GridGenerator::cylinder (Triangulation<3> &tria, -void GridGenerator::hyper_shell (Triangulation<3> &, - const Point<3> &, - const double , - const double , - const unsigned int ) +// Implementation for 3D only +template +void +GridGenerator::half_hyper_ball (Triangulation&, + const Point&, + const double) { Assert (false, ExcNotImplemented()); } -void GridGenerator::cylinder_shell (Triangulation<3> &tria, +// Implementation for 3D only +template +void GridGenerator::hyper_shell (Triangulation&, + const Point&, + const double, + const double, + const unsigned int) +{ + Assert (false, ExcNotImplemented()); +} + + +// Implementation for 3D only +template +void +GridGenerator::half_hyper_shell (Triangulation&, + const Point&, + const double, + const double, + const unsigned int) +{ + Assert (false, ExcNotImplemented()); +} + + + +// Implementation for 3D only +template +void GridGenerator::cylinder_shell (Triangulation &tria, const double length, const double inner_radius, const double outer_radius, @@ -1227,26 +1358,26 @@ void GridGenerator::cylinder_shell (Triangulation<3> &tria, // first N ones are on the // outer one, and all are // numbered counter-clockwise - std::vector > vertices_2d(2*N_r); + std::vector > vertices_2d(2*N_r); for (unsigned int i=0; i( std::cos(2*pi * i/N_r), + vertices_2d[i] = Point( std::cos(2*pi * i/N_r), std::sin(2*pi * i/N_r)) * outer_radius; vertices_2d[i+N_r] = vertices_2d[i] * (inner_radius/outer_radius); }; - std::vector > vertices_3d; + std::vector > vertices_3d; vertices_3d.reserve (2*N_r*(N_z+1)); for (unsigned int j=0; j<=N_z; ++j) for (unsigned int i=0; i<2*N_r; ++i) { - const Point<3> v (vertices_2d[i][0], + const Point v (vertices_2d[i][0], vertices_2d[i][1], j*length/N_z); vertices_3d.push_back (v); } - std::vector > cells (N_r*N_z, CellData<3>()); + std::vector > cells (N_r*N_z, CellData()); for (unsigned int j=0; j &tria, } + #endif +// make the following function inline. this is so that it becomes marked +// internal/weak for the linker and we don't get multiply defined errors +// when linking with more than one dimension at a time. Usually we used +// the trick of putting these functions in a .all_dimensions.cc file, but +// this is not necessary here as this is an internal only function. +inline +void GridGenerator::laplace_solve (const SparseMatrix &S, + const std::map &m, + Vector &u) +{ + const unsigned int n_dofs=S.n(); + FilteredMatrix > SF (S); + SolverControl control (1000, 1.e-10, false, false); + PrimitiveVectorMemory > mem; + SolverCG > solver (control, mem); + PreconditionJacobi > > prec; + Vector f(n_dofs); + + SF.add_constraints(m); + prec.initialize (SF); + SF.apply_constraints (f, true); + solver.solve(SF, u, f, prec); +} + + #if deal_II_dimension == 1 -void GridGenerator::laplace_transformation (Triangulation<1> &, - const std::map > &) +// Implementation for 1D only +template +void GridGenerator::laplace_transformation (Triangulation &, + const std::map > &) { Assert(false, ExcNotImplemented()); } #else +// Implementation for dimensions except 1 template void GridGenerator::laplace_transformation (Triangulation &tria, const std::map > &new_points) @@ -1367,40 +1527,29 @@ void GridGenerator::laplace_transformation (Triangulation &tria, #endif -// make the following function inline. this is so that it becomes marked -// internal/weak for the linker and we don't get multiply defined errors -// when linking with more than one dimension at a time. Usually we used -// the trick of putting these functions in a .all_dimensions.cc file, but -// this is not necessary here as this is an internal only function. -inline -void GridGenerator::laplace_solve (const SparseMatrix &S, - const std::map &m, - Vector &u) -{ - const unsigned int n_dofs=S.n(); - FilteredMatrix > SF (S); - SolverControl control (1000, 1.e-10, false, false); - PrimitiveVectorMemory > mem; - SolverCG > solver (control, mem); - PreconditionJacobi > > prec; - Vector f(n_dofs); - - SF.add_constraints(m); - prec.initialize (SF); - SF.apply_constraints (f, true); - solver.solve(SF, u, f, prec); -} +// explicit instantiations +template void +GridGenerator::hyper_cube ( + Triangulation &, const double, const double); +template void +GridGenerator::subdivided_hyper_cube ( + Triangulation &, + const unsigned int, const double, const double); -// explicit instantiations template void GridGenerator::hyper_rectangle ( Triangulation &, const Point&, const Point&, const bool); + template void -GridGenerator::hyper_cube ( - Triangulation &, const double, const double); +GridGenerator::subdivided_hyper_rectangle +(Triangulation &, + const std::vector&, + const Point&, + const Point&, bool); + template void GridGenerator::parallelogram ( Triangulation &, @@ -1408,23 +1557,50 @@ GridGenerator::parallelogram ( const bool); template void -GridGenerator::subdivided_hyper_cube ( - Triangulation &, - const unsigned int, const double, const double); +GridGenerator::enclosed_hyper_cube ( + Triangulation&, double, double, double, bool); template void -GridGenerator::subdivided_hyper_rectangle -(Triangulation &, - const std::vector&, - const Point&, - const Point&, - const bool); +GridGenerator::hyper_ball ( + Triangulation&, + const Point&, double); + +template void +GridGenerator::cylinder ( + Triangulation&, double, double); + + +template void +GridGenerator::hyper_L ( + Triangulation&, double, double); + +template void +GridGenerator::hyper_cube_slit ( + Triangulation&, double, double, bool); + +template void +GridGenerator::hyper_shell ( + Triangulation&, + const Point&, double, double, unsigned int); + + +template void +GridGenerator::cylinder_shell ( + Triangulation&, + double, double, double, unsigned int, unsigned int); + +template void +GridGenerator::half_hyper_ball ( + Triangulation&, const Point&, double); + +template void +GridGenerator::half_hyper_shell ( + Triangulation&, + const Point&, double, double, unsigned int); + -#if deal_II_dimension != 1 template void GridGenerator:: laplace_transformation (Triangulation &, const std::map > &); - -#endif -- 2.39.5