From 6ec0977f6fb87e292cdeffe6f61f744996c85f96 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Mon, 19 Apr 2004 21:09:40 +0000 Subject: [PATCH] Implement a cylinder shell in 3d. git-svn-id: https://svn.dealii.org/trunk@9066 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/grid/grid_generator.h | 22 +++++ deal.II/deal.II/source/grid/grid_generator.cc | 80 +++++++++++++++++++ tests/bits/Makefile | 3 +- tests/bits/cylinder_shell_01.cc | 44 ++++++++++ 4 files changed, 148 insertions(+), 1 deletion(-) create mode 100644 tests/bits/cylinder_shell_01.cc diff --git a/deal.II/deal.II/include/grid/grid_generator.h b/deal.II/deal.II/include/grid/grid_generator.h index 206d93254a..bd07bcf03e 100644 --- a/deal.II/deal.II/include/grid/grid_generator.h +++ b/deal.II/deal.II/include/grid/grid_generator.h @@ -536,6 +536,28 @@ class GridGenerator const double outer_radius, const unsigned int n_cells = 0); + /** + * Produce a domain that is the space + * between two cylinders in 3d, with + * given length, inner and outer radius + * and a given number of elements for + * this initial triangulation. If @p + * n_radial_cells is zero (as is the + * default), then it is computed + * adaptively such that the resulting + * elements have the least aspect + * ratio. The same holds for @p + * n_axial_cells. + * + * The triangulation needs to be void + * upon calling this function. + */ + static void cylinder_shell (Triangulation<3> &tria, + const double length, + const double inner_radius, + const double outer_radius, + const unsigned int n_radial_cells = 0, + const unsigned int n_axial_cells = 0); /** * This class produces a half hyper-ball, diff --git a/deal.II/deal.II/source/grid/grid_generator.cc b/deal.II/deal.II/source/grid/grid_generator.cc index 55e392be21..8b1b7a1006 100644 --- a/deal.II/deal.II/source/grid/grid_generator.cc +++ b/deal.II/deal.II/source/grid/grid_generator.cc @@ -1172,6 +1172,86 @@ void GridGenerator::hyper_shell (Triangulation<3> &, } + +void GridGenerator::cylinder_shell (Triangulation<3> &tria, + const double length, + const double inner_radius, + const double outer_radius, + const unsigned int n_radial_cells, + const unsigned int n_axial_cells) +{ + Assert ((inner_radius > 0) && (inner_radius < outer_radius), + ExcInvalidRadii ()); + + const double pi = deal_II_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_r = (n_radial_cells == 0 ? + static_cast + (std::ceil((2*pi* (outer_radius + inner_radius)/2) / + (outer_radius - inner_radius))) : + n_radial_cells); + const unsigned int N_z = (n_axial_cells == 0 ? + static_cast + (std::ceil (length / + (2*pi*(outer_radius + inner_radius)/2/N_r))) : + n_axial_cells); + + // set up N vertices on the + // outer and N vertices on + // the inner circle. the + // first N ones are on the + // outer one, and all are + // numbered counter-clockwise + std::vector > vertices_2d(2*N_r); + for (unsigned int i=0; i( 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; + 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], + vertices_2d[i][1], + j*length/N_z); + vertices_3d.push_back (v); + } + + std::vector > cells (N_r*N_z, CellData<3>()); + + for (unsigned int j=0; j works as +// expected + + +#include "../tests.h" +#include +#include +#include +#include +#include +#include + + + +int main () +{ + std::ofstream logfile("cylinder_shell_01.output"); + deallog.attach(logfile); + deallog.depth_console(0); + logfile.precision (2); + + // generate a hyperball in 3d + Triangulation<3> tria; + GridGenerator::cylinder_shell (tria, 1, .8, 1); + + // make sure that all cells have positive + // volume + for (Triangulation<3>::active_cell_iterator cell=tria.begin_active(); + cell!=tria.end(); ++cell) + deallog << cell << ' ' << cell->measure () << std::endl; +} -- 2.39.5