From 568e55196178c1cd820a8d9225d32f7ab18356f2 Mon Sep 17 00:00:00 2001 From: heister Date: Thu, 7 Jul 2011 19:40:51 +0000 Subject: [PATCH] Add GridGenerator::quarter_hyper_shell. Only implemented in 2d for now. git-svn-id: https://svn.dealii.org/trunk@23931 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/grid/grid_generator.h | 33 ++++++ deal.II/source/grid/grid_generator.cc | 106 +++++++++++++++++- 2 files changed, 138 insertions(+), 1 deletion(-) diff --git a/deal.II/include/deal.II/grid/grid_generator.h b/deal.II/include/deal.II/grid/grid_generator.h index 52dccacb41..188f02aad6 100644 --- a/deal.II/include/deal.II/grid/grid_generator.h +++ b/deal.II/include/deal.II/grid/grid_generator.h @@ -685,6 +685,39 @@ class GridGenerator const double outer_radius, const unsigned int n_cells = 0); + + /** + * Produce a quarter hyper-shell, + * i.e. the space between two + * circles in two space + * dimensions and the region + * between two spheres in 3d, + * with given inner and outer + * radius and a given number of + * elements for this initial + * triangulation. All components are + * restricted to be positive, so the + * opening angle is 90 degrees. + * + * If the number of + * initial cells is zero (as is + * the default), then it is + * computed adaptively such that + * the resulting elements have + * the least aspect ratio. + * + * @note The triangulation needs to be + * void upon calling this + * function. Only implemented in 2d so far. + */ + template + static void quarter_hyper_shell (Triangulation &tria, + const Point ¢er, + const double inner_radius, + const double outer_radius, + const unsigned int n_cells = 0, + const bool colorize = false); + /** * Produce a domain that is the space * between two cylinders in 3d, with diff --git a/deal.II/source/grid/grid_generator.cc b/deal.II/source/grid/grid_generator.cc index 0bf9d67358..623006c4e3 100644 --- a/deal.II/source/grid/grid_generator.cc +++ b/deal.II/source/grid/grid_generator.cc @@ -1213,7 +1213,16 @@ GridGenerator::half_hyper_shell (Triangulation<1>&, Assert (false, ExcNotImplemented()); } - +template <> +void GridGenerator::quarter_hyper_shell (Triangulation<1> &tria, + const Point<1> ¢er, + const double inner_radius, + const double outer_radius, + const unsigned int n_cells, + const bool colorize) +{ + Assert (false, ExcNotImplemented()); +} template <> void GridGenerator::enclosed_hyper_cube (Triangulation<2> &tria, @@ -1689,6 +1698,90 @@ GridGenerator::half_hyper_shell (Triangulation<2> &tria, } +template <> +void GridGenerator::quarter_hyper_shell (Triangulation<2> &tria, + const Point<2> ¢er, + const double inner_radius, + const double outer_radius, + const unsigned int n_cells, + const bool colorize) + { + Assert ((inner_radius > 0) && (inner_radius < outer_radius), + ExcInvalidRadii ()); + + const double pi = numbers::PI; + // determine the number of cells + // for the grid. if not provided by + // the user determine it such that + // the length of each cell on the + // median (in the middle between + // the two circles) is equal to its + // radial extent (which is the + // difference between the two + // radii) + const unsigned int N = (n_cells == 0 ? + static_cast + (std::ceil((pi* (outer_radius + inner_radius)/4) / + (outer_radius - inner_radius))) : + n_cells); + + // set up N+1 vertices on the + // outer and N+1 vertices on + // the inner circle. the + // first N+1 ones are on the + // outer one, and all are + // numbered counter-clockwise + std::vector > vertices(2*(N+1)); + for (unsigned int i=0; i<=N; ++i) + { + // enforce that the x-coordinates + // of the last point is exactly + // zero (contrary to what we may + // compute using the imprecise + // value of pi) + vertices[i] = Point<2>( ( (i==N) ? + 0 : + std::cos(pi * i/N/2) ), + std::sin(pi * i/N/2)) * outer_radius; + vertices[i+N+1] = vertices[i] * (inner_radius/outer_radius); + + vertices[i] += center; + vertices[i+N+1] += center; + }; + + + std::vector > cells (N, CellData<2>()); + + for (unsigned int i=0; i::cell_iterator cell = tria.begin(); + int i=0; + for (;cell!=tria.end();++cell) + { + std::cout << i++ << ": "; + for (int j=0;j<4;++j) + std::cout << cell->face(j)->at_boundary(); + std::cout << std::endl; + cell->face(2)->set_boundary_indicator(1); + } + tria.begin()->face(0)->set_boundary_indicator(3); + + tria.last()->face(1)->set_boundary_indicator(2); + } +} + // Implementation for 3D only @@ -2603,6 +2696,17 @@ GridGenerator::half_hyper_shell (Triangulation<3>& tria, } +template <> +void GridGenerator::quarter_hyper_shell (Triangulation<3> &tria, + const Point<3> ¢er, + const double inner_radius, + const double outer_radius, + const unsigned int n_cells, + const bool colorize) +{ + Assert (false, ExcNotImplemented()); +} + // Implementation for 3D only template <> -- 2.39.5