From 0c92994d44aa2bc8bdab9f1ca6a7ec13f2b46d99 Mon Sep 17 00:00:00 2001 From: David Wells Date: Sun, 15 Jul 2018 19:32:43 -0400 Subject: [PATCH] Add GridGenerator::concentric_hyper_shells. This function generates a grid comprised of (as the name states) concentric hyper shells, where the shells are (controlled by the skewness parameter) clustered near the inner boundary. --- .../images/concentric_hyper_shells_2d.svg | 737 ++++++++++++++++++ doc/news/changes/minor/20180715DavidWells | 5 + include/deal.II/grid/grid_generator.h | 82 ++ source/grid/grid_generator.cc | 121 +++ source/grid/grid_generator.inst.in | 9 + tests/grid/concentric_hyper_shell_01.cc | 90 +++ tests/grid/concentric_hyper_shell_01.output | 69 ++ tests/grid/concentric_hyper_shell_02.cc | 88 +++ tests/grid/concentric_hyper_shell_02.output | 69 ++ 9 files changed, 1270 insertions(+) create mode 100644 doc/doxygen/images/concentric_hyper_shells_2d.svg create mode 100644 doc/news/changes/minor/20180715DavidWells create mode 100644 tests/grid/concentric_hyper_shell_01.cc create mode 100644 tests/grid/concentric_hyper_shell_01.output create mode 100644 tests/grid/concentric_hyper_shell_02.cc create mode 100644 tests/grid/concentric_hyper_shell_02.output diff --git a/doc/doxygen/images/concentric_hyper_shells_2d.svg b/doc/doxygen/images/concentric_hyper_shells_2d.svg new file mode 100644 index 0000000000..ac9383d959 --- /dev/null +++ b/doc/doxygen/images/concentric_hyper_shells_2d.svg @@ -0,0 +1,737 @@ + + + + + + image/svg+xml + + Qt SVG Document + + + + + Qt SVG Document + Generated with Qt + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/doc/news/changes/minor/20180715DavidWells b/doc/news/changes/minor/20180715DavidWells new file mode 100644 index 0000000000..b589f5e14b --- /dev/null +++ b/doc/news/changes/minor/20180715DavidWells @@ -0,0 +1,5 @@ +New: A new function GridGenerator::concentric_hyper_shells has been added. This +function generates a grid comprised of shells of varying widths, where the +shells near the center are much thinner than those near the outer boundary. +
+(David Wells, 2018/07/15) diff --git a/include/deal.II/grid/grid_generator.h b/include/deal.II/grid/grid_generator.h index e470279f47..56da7be560 100644 --- a/include/deal.II/grid/grid_generator.h +++ b/include/deal.II/grid/grid_generator.h @@ -991,6 +991,88 @@ namespace GridGenerator const unsigned int repetitions = 1, const bool colorize = false); + /** + * Produce a grid consisting of concentric shells. The primary difference + * between this function and GridGenerator::hyper_shell is that this + * function permits unevenly spaced (in the radial direction) coarse level + * cells. + * + * The parameters @p center, @p inner_radius, and @p outer_radius behave in + * the same way as the first three arguments to + * GridGenerator::hyper_shell. @p n_shells gives the total number of shells + * to use (i.e., the number of cells in the radial direction). The outer + * radius of the $k$th shell is given by + * + * @f[ + * r = r_{\text{inner}} + (r_\text{outer} - r_\text{inner}) + * \frac{1 - \tanh(\text{skewness}(1 - k/\text{n_shells}))} + * {\tanh(\text{skewness})} + * @f] + * + * where @p skewness is a parameter controlling the shell spacing in the + * radial direction: values of @p skewness close to zero correspond to even + * spacing, while larger values of @p skewness (such as $2$ or $3$) + * correspond to shells biased to the inner radius. + * + * @p n_cells_per_shell is the same as in GridGenerator::hyper_shell: in 2d + * the default choice of zero will result in 8 cells per shell (and 12 in + * 3d). The only valid values in 3d are 6 (the default), 12, and 96 cells: + * see the documentation of GridGenerator::hyper_shell for more information. + * + * If @p colorize is true then the outer boundary of the merged + * shells has a boundary id of $1$ and the inner boundary has a boundary id + * of $0$. + * + * Example: The following code (see, e.g., step-10 for instructions on how + * to visualize GNUPLOT output) + * + * @code + * #include + * + * #include + * #include + * #include + * + * #include + * + * int main() + * { + * using namespace dealii; + * + * Triangulation<2> triangulation; + * GridGenerator::concentric_hyper_shells(triangulation, + * Point<2>(), + * 1.0, + * 2.0, + * 5u, + * 2.0); + * + * GridOut grid_out; + * GridOutFlags::Gnuplot gnuplot_flags(false, 10, true); + * grid_out.set_flags(gnuplot_flags); + * + * const MappingQGeneric<2> mapping(3); + * std::ofstream out("out.gpl"); + * grid_out.write_gnuplot(triangulation, out, &mapping); + * } + * @endcode + * + * generates the following output: + * + * @image html concentric_hyper_shells_2d.svg + * + */ + template + void + concentric_hyper_shells(Triangulation &triangulation, + const Point & center, + const double inner_radius = 0.125, + const double outer_radius = 0.25, + const unsigned int n_shells = 1, + const double skewness = 0.1, + const unsigned int n_cells_per_shell = 0, + const bool colorize = false); + /** * Produce a ring of cells in 3d that is cut open, twisted and glued * together again. This results in a kind of moebius-loop. diff --git a/source/grid/grid_generator.cc b/source/grid/grid_generator.cc index 313a3302f0..a05207881c 100644 --- a/source/grid/grid_generator.cc +++ b/source/grid/grid_generator.cc @@ -4701,6 +4701,127 @@ namespace GridGenerator + template + void + concentric_hyper_shells(Triangulation &triangulation, + const Point & center, + const double inner_radius, + const double outer_radius, + const unsigned int n_shells, + const double skewness, + const unsigned int n_cells, + const bool colorize) + { + Assert(dim == 2 || dim == 3, ExcNotImplemented()); + (void)colorize; + (void)n_cells; + Assert(inner_radius < outer_radius, + ExcMessage("outer_radius has to be bigger than inner_radius.")); + if (n_shells == 0) + return; // empty Triangulation + + std::vector radii; + radii.push_back(inner_radius); + for (unsigned int shell_n = 1; shell_n < n_shells; ++shell_n) + if (skewness == 0.0) + // same as below, but works in the limiting case of zero skewness + radii.push_back(inner_radius + + (outer_radius - inner_radius) * + (1.0 - (1.0 - double(shell_n) / n_shells))); + else + radii.push_back( + inner_radius + + (outer_radius - inner_radius) * + (1.0 - std::tanh(skewness * (1.0 - double(shell_n) / n_shells)) / + std::tanh(skewness))); + radii.push_back(outer_radius); + + auto min_line_length = [](const Triangulation &tria) -> double { + double length = std::numeric_limits::max(); + for (const auto cell : tria.active_cell_iterators()) + for (unsigned int n = 0; n < GeometryInfo::lines_per_cell; ++n) + length = std::min(length, cell->line(n)->diameter()); + return length; + }; + + double grid_vertex_tolerance = 0.0; + for (unsigned int shell_n = 0; shell_n < radii.size() - 1; ++shell_n) + { + Triangulation current_shell; + GridGenerator::hyper_shell(current_shell, + center, + radii[shell_n], + radii[shell_n + 1], + n_cells == 0 ? (dim == 2 ? 8 : 12) : + n_cells); + + // The innermost shell has the smallest cells: use that to set the + // vertex merging tolerance + if (grid_vertex_tolerance == 0.0) + grid_vertex_tolerance = 0.5 * min_line_length(current_shell); + + Triangulation temp(std::move(triangulation)); + triangulation.clear(); + GridGenerator::merge_triangulations(current_shell, + temp, + triangulation, + grid_vertex_tolerance); + } + + const types::manifold_id manifold_id = 0; + triangulation.set_all_manifold_ids(manifold_id); + if (dim == 2) + triangulation.set_manifold(manifold_id, PolarManifold(center)); + else if (dim == 3) + triangulation.set_manifold(manifold_id, SphericalManifold(center)); + + // We use boundary vertex positions to see if things are on the inner or + // outer boundary. + constexpr double radial_vertex_tolerance = + 100.0 * std::numeric_limits::epsilon(); + auto assert_vertex_distance_within_tolerance = + [center, radial_vertex_tolerance]( + const TriaIterator> face, + const double radius) { + (void)face; + (void)radius; + for (unsigned int vertex_n = 0; + vertex_n < GeometryInfo::vertices_per_face; + ++vertex_n) + { + Assert(std::abs((face->vertex(vertex_n) - center).norm() - radius) < + (center.norm() + radius) * radial_vertex_tolerance, + ExcInternalError()); + } + }; + if (colorize) + for (const auto &cell : triangulation.active_cell_iterators()) + for (unsigned int face_n = 0; + face_n < GeometryInfo::faces_per_cell; + ++face_n) + { + auto face = cell->face(face_n); + if (face->at_boundary()) + { + if (((face->vertex(0) - center).norm() - inner_radius) < + (center.norm() + inner_radius) * radial_vertex_tolerance) + { + // we must be at an inner face, but check + assert_vertex_distance_within_tolerance(face, inner_radius); + face->set_all_boundary_ids(0); + } + else + { + // we must be at an outer face, but check + assert_vertex_distance_within_tolerance(face, outer_radius); + face->set_all_boundary_ids(1); + } + } + } + } + + + template <> void hyper_cube_with_cylindrical_hole(Triangulation<3> & triangulation, const double inner_radius, diff --git a/source/grid/grid_generator.inst.in b/source/grid/grid_generator.inst.in index b64c8f60a4..6ac4262276 100644 --- a/source/grid/grid_generator.inst.in +++ b/source/grid/grid_generator.inst.in @@ -183,6 +183,15 @@ for (deal_II_dimension : DIMENSIONS) const bool); #endif + template void + concentric_hyper_shells(Triangulation &, + const Point &, + const double, + const double, + const unsigned int, + const double, + const unsigned int, + const bool); \} } diff --git a/tests/grid/concentric_hyper_shell_01.cc b/tests/grid/concentric_hyper_shell_01.cc new file mode 100644 index 0000000000..916298b7fa --- /dev/null +++ b/tests/grid/concentric_hyper_shell_01.cc @@ -0,0 +1,90 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2018 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + +#include +#include + +#include "../tests.h" + +// Test GridGenerator::concentric_hyper_shell with skewness = 2.0 + +int +main() +{ + initlog(); + + // 2D test + { + Triangulation<2> triangulation; + const Point<2> center(0.0, 1.0); + GridGenerator::concentric_hyper_shells( + triangulation, center, 1.0, 4.0, 3u, 2.0, 4, true); + + for (const auto &cell : triangulation.active_cell_iterators()) + { + deallog << "vertices: " << cell->vertex(0) << ", " << cell->vertex(1) + << ", " << cell->vertex(2) << ", " << cell->vertex(3) + << std::endl; + + bool manifold_ids_are_zero = cell->manifold_id() == 0; + for (unsigned int face_n = 0; face_n < GeometryInfo<2>::faces_per_cell; + ++face_n) + manifold_ids_are_zero &= cell->face(face_n)->manifold_id() == 0; + AssertThrow(manifold_ids_are_zero, ExcInternalError()); + + for (unsigned int face_n = 0; face_n < GeometryInfo<2>::faces_per_cell; + ++face_n) + if (cell->face(face_n)->at_boundary()) + deallog << "boundary face center distance to origin: " + << (cell->face(face_n)->center(/*respect_manifold*/ true) - + center) + .norm() + << " boundary id: " << cell->face(face_n)->boundary_id() + << std::endl; + } + } + + // 3D test + { + Triangulation<3> triangulation; + const Point<3> center(0.0, 1.0, 2.0); + GridGenerator::concentric_hyper_shells( + triangulation, center, 1.0, 5.0, 2u, 2.0, 0, true); + + for (const auto &cell : triangulation.active_cell_iterators()) + { + deallog << "vertices: " << cell->vertex(0) << ", " << cell->vertex(1) + << ", " << cell->vertex(2) << ", " << cell->vertex(3) << ", " + << cell->vertex(4) << ", " << cell->vertex(5) << ", " + << cell->vertex(6) << ", " << cell->vertex(7) << std::endl; + + bool manifold_ids_are_zero = cell->manifold_id() == 0; + for (unsigned int face_n = 0; face_n < GeometryInfo<3>::faces_per_cell; + ++face_n) + manifold_ids_are_zero &= cell->face(face_n)->manifold_id() == 0; + AssertThrow(manifold_ids_are_zero, ExcInternalError()); + + for (unsigned int face_n = 0; face_n < GeometryInfo<3>::faces_per_cell; + ++face_n) + if (cell->face(face_n)->at_boundary()) + deallog << "boundary face center distance to origin: " + << (cell->face(face_n)->center(/*respect_manifold*/ true) - + center) + .norm() + << " boundary id: " << cell->face(face_n)->boundary_id() + << std::endl; + } + } +} diff --git a/tests/grid/concentric_hyper_shell_01.output b/tests/grid/concentric_hyper_shell_01.output new file mode 100644 index 0000000000..adc606bda6 --- /dev/null +++ b/tests/grid/concentric_hyper_shell_01.output @@ -0,0 +1,69 @@ + +DEAL::vertices: 4.00000 1.00000, 2.44929e-16 5.00000, 2.18641 1.00000, 1.33879e-16 3.18641 +DEAL::boundary face center distance to origin: 4.00000 boundary id: 1 +DEAL::vertices: 2.44929e-16 5.00000, -4.00000 1.00000, 1.33879e-16 3.18641, -2.18641 1.00000 +DEAL::boundary face center distance to origin: 4.00000 boundary id: 1 +DEAL::vertices: -4.00000 1.00000, -7.34788e-16 -3.00000, -2.18641 1.00000, -4.01637e-16 -1.18641 +DEAL::boundary face center distance to origin: 4.00000 boundary id: 1 +DEAL::vertices: -7.34788e-16 -3.00000, 4.00000 1.00000, -4.01637e-16 -1.18641, 2.18641 1.00000 +DEAL::boundary face center distance to origin: 4.00000 boundary id: 1 +DEAL::vertices: 2.18641 1.00000, 1.33879e-16 3.18641, 1.29242 1.00000, 7.91377e-17 2.29242 +DEAL::vertices: 1.33879e-16 3.18641, -2.18641 1.00000, 7.91377e-17 2.29242, -1.29242 1.00000 +DEAL::vertices: -2.18641 1.00000, -4.01637e-16 -1.18641, -1.29242 1.00000, -2.37413e-16 -0.292417 +DEAL::vertices: -4.01637e-16 -1.18641, 2.18641 1.00000, -2.37413e-16 -0.292417, 1.29242 1.00000 +DEAL::vertices: 1.29242 1.00000, 7.91377e-17 2.29242, 1.00000 1.00000, 6.12323e-17 2.00000 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 7.91377e-17 2.29242, -1.29242 1.00000, 6.12323e-17 2.00000, -1.00000 1.00000 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: -1.29242 1.00000, -2.37413e-16 -0.292417, -1.00000 1.00000, -1.83697e-16 0.00000 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: -2.37413e-16 -0.292417, 1.29242 1.00000, -1.83697e-16 0.00000, 1.00000 1.00000 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 0.00000 -0.839949 2.00000, -1.06229 -0.0622949 3.06229, -1.06229 -0.0622949 0.937705, -1.83995 1.00000 2.00000, 0.00000 -4.00000 2.00000, -2.88675 -1.88675 4.88675, -2.88675 -1.88675 -0.886751, -5.00000 1.00000 2.00000 +DEAL::boundary face center distance to origin: 5.00000 boundary id: 1 +DEAL::vertices: -1.06229 -0.0622949 3.06229, 0.00000 1.00000 3.83995, -1.83995 1.00000 2.00000, -1.06229 2.06229 3.06229, -2.88675 -1.88675 4.88675, 0.00000 1.00000 7.00000, -5.00000 1.00000 2.00000, -2.88675 3.88675 4.88675 +DEAL::boundary face center distance to origin: 5.00000 boundary id: 1 +DEAL::vertices: 0.00000 -0.839949 2.00000, 1.06229 -0.0622949 3.06229, -1.06229 -0.0622949 3.06229, 0.00000 1.00000 3.83995, 0.00000 -4.00000 2.00000, 2.88675 -1.88675 4.88675, -2.88675 -1.88675 4.88675, 0.00000 1.00000 7.00000 +DEAL::boundary face center distance to origin: 5.00000 boundary id: 1 +DEAL::vertices: 1.06229 -0.0622949 0.937705, 1.83995 1.00000 2.00000, 0.00000 -0.839949 2.00000, 1.06229 -0.0622949 3.06229, 2.88675 -1.88675 -0.886751, 5.00000 1.00000 2.00000, 0.00000 -4.00000 2.00000, 2.88675 -1.88675 4.88675 +DEAL::boundary face center distance to origin: 5.00000 boundary id: 1 +DEAL::vertices: 1.83995 1.00000 2.00000, 1.06229 2.06229 3.06229, 1.06229 -0.0622949 3.06229, 0.00000 1.00000 3.83995, 5.00000 1.00000 2.00000, 2.88675 3.88675 4.88675, 2.88675 -1.88675 4.88675, 0.00000 1.00000 7.00000 +DEAL::boundary face center distance to origin: 5.00000 boundary id: 1 +DEAL::vertices: 1.06229 2.06229 3.06229, 0.00000 2.83995 2.00000, 0.00000 1.00000 3.83995, -1.06229 2.06229 3.06229, 2.88675 3.88675 4.88675, 0.00000 6.00000 2.00000, 0.00000 1.00000 7.00000, -2.88675 3.88675 4.88675 +DEAL::boundary face center distance to origin: 5.00000 boundary id: 1 +DEAL::vertices: 1.83995 1.00000 2.00000, 1.06229 2.06229 0.937705, 1.06229 2.06229 3.06229, 0.00000 2.83995 2.00000, 5.00000 1.00000 2.00000, 2.88675 3.88675 -0.886751, 2.88675 3.88675 4.88675, 0.00000 6.00000 2.00000 +DEAL::boundary face center distance to origin: 5.00000 boundary id: 1 +DEAL::vertices: 1.06229 -0.0622949 0.937705, 0.00000 1.00000 0.160051, 1.83995 1.00000 2.00000, 1.06229 2.06229 0.937705, 2.88675 -1.88675 -0.886751, 0.00000 1.00000 -3.00000, 5.00000 1.00000 2.00000, 2.88675 3.88675 -0.886751 +DEAL::boundary face center distance to origin: 5.00000 boundary id: 1 +DEAL::vertices: 0.00000 1.00000 0.160051, -1.06229 2.06229 0.937705, 1.06229 2.06229 0.937705, 0.00000 2.83995 2.00000, 0.00000 1.00000 -3.00000, -2.88675 3.88675 -0.886751, 2.88675 3.88675 -0.886751, 0.00000 6.00000 2.00000 +DEAL::boundary face center distance to origin: 5.00000 boundary id: 1 +DEAL::vertices: -1.06229 2.06229 0.937705, -1.83995 1.00000 2.00000, 0.00000 2.83995 2.00000, -1.06229 2.06229 3.06229, -2.88675 3.88675 -0.886751, -5.00000 1.00000 2.00000, 0.00000 6.00000 2.00000, -2.88675 3.88675 4.88675 +DEAL::boundary face center distance to origin: 5.00000 boundary id: 1 +DEAL::vertices: 0.00000 1.00000 0.160051, -1.06229 -0.0622949 0.937705, -1.06229 2.06229 0.937705, -1.83995 1.00000 2.00000, 0.00000 1.00000 -3.00000, -2.88675 -1.88675 -0.886751, -2.88675 3.88675 -0.886751, -5.00000 1.00000 2.00000 +DEAL::boundary face center distance to origin: 5.00000 boundary id: 1 +DEAL::vertices: 1.06229 -0.0622949 0.937705, 0.00000 -0.839949 2.00000, 0.00000 1.00000 0.160051, -1.06229 -0.0622949 0.937705, 2.88675 -1.88675 -0.886751, 0.00000 -4.00000 2.00000, 0.00000 1.00000 -3.00000, -2.88675 -1.88675 -0.886751 +DEAL::boundary face center distance to origin: 5.00000 boundary id: 1 +DEAL::vertices: 0.00000 0.00000 2.00000, -0.577350 0.422650 2.57735, -0.577350 0.422650 1.42265, -1.00000 1.00000 2.00000, 0.00000 -0.839949 2.00000, -1.06229 -0.0622949 3.06229, -1.06229 -0.0622949 0.937705, -1.83995 1.00000 2.00000 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: -0.577350 0.422650 2.57735, 0.00000 1.00000 3.00000, -1.00000 1.00000 2.00000, -0.577350 1.57735 2.57735, -1.06229 -0.0622949 3.06229, 0.00000 1.00000 3.83995, -1.83995 1.00000 2.00000, -1.06229 2.06229 3.06229 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 0.00000 0.00000 2.00000, 0.577350 0.422650 2.57735, -0.577350 0.422650 2.57735, 0.00000 1.00000 3.00000, 0.00000 -0.839949 2.00000, 1.06229 -0.0622949 3.06229, -1.06229 -0.0622949 3.06229, 0.00000 1.00000 3.83995 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 0.577350 0.422650 1.42265, 1.00000 1.00000 2.00000, 0.00000 0.00000 2.00000, 0.577350 0.422650 2.57735, 1.06229 -0.0622949 0.937705, 1.83995 1.00000 2.00000, 0.00000 -0.839949 2.00000, 1.06229 -0.0622949 3.06229 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 1.00000 1.00000 2.00000, 0.577350 1.57735 2.57735, 0.577350 0.422650 2.57735, 0.00000 1.00000 3.00000, 1.83995 1.00000 2.00000, 1.06229 2.06229 3.06229, 1.06229 -0.0622949 3.06229, 0.00000 1.00000 3.83995 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 0.577350 1.57735 2.57735, 0.00000 2.00000 2.00000, 0.00000 1.00000 3.00000, -0.577350 1.57735 2.57735, 1.06229 2.06229 3.06229, 0.00000 2.83995 2.00000, 0.00000 1.00000 3.83995, -1.06229 2.06229 3.06229 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 1.00000 1.00000 2.00000, 0.577350 1.57735 1.42265, 0.577350 1.57735 2.57735, 0.00000 2.00000 2.00000, 1.83995 1.00000 2.00000, 1.06229 2.06229 0.937705, 1.06229 2.06229 3.06229, 0.00000 2.83995 2.00000 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 0.577350 0.422650 1.42265, 0.00000 1.00000 1.00000, 1.00000 1.00000 2.00000, 0.577350 1.57735 1.42265, 1.06229 -0.0622949 0.937705, 0.00000 1.00000 0.160051, 1.83995 1.00000 2.00000, 1.06229 2.06229 0.937705 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 0.00000 1.00000 1.00000, -0.577350 1.57735 1.42265, 0.577350 1.57735 1.42265, 0.00000 2.00000 2.00000, 0.00000 1.00000 0.160051, -1.06229 2.06229 0.937705, 1.06229 2.06229 0.937705, 0.00000 2.83995 2.00000 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: -0.577350 1.57735 1.42265, -1.00000 1.00000 2.00000, 0.00000 2.00000 2.00000, -0.577350 1.57735 2.57735, -1.06229 2.06229 0.937705, -1.83995 1.00000 2.00000, 0.00000 2.83995 2.00000, -1.06229 2.06229 3.06229 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 0.00000 1.00000 1.00000, -0.577350 0.422650 1.42265, -0.577350 1.57735 1.42265, -1.00000 1.00000 2.00000, 0.00000 1.00000 0.160051, -1.06229 -0.0622949 0.937705, -1.06229 2.06229 0.937705, -1.83995 1.00000 2.00000 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 0.577350 0.422650 1.42265, 0.00000 0.00000 2.00000, 0.00000 1.00000 1.00000, -0.577350 0.422650 1.42265, 1.06229 -0.0622949 0.937705, 0.00000 -0.839949 2.00000, 0.00000 1.00000 0.160051, -1.06229 -0.0622949 0.937705 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 diff --git a/tests/grid/concentric_hyper_shell_02.cc b/tests/grid/concentric_hyper_shell_02.cc new file mode 100644 index 0000000000..e8d4a70852 --- /dev/null +++ b/tests/grid/concentric_hyper_shell_02.cc @@ -0,0 +1,88 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2018 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + +#include +#include + +#include "../tests.h" + +// Test GridGenerator::concentric_hyper_shell with skewness = 0.0 + +int +main() +{ + initlog(); + + // 2D test + { + Triangulation<2> triangulation; + const Point<2> center(0.0, 1.0); + GridGenerator::concentric_hyper_shells( + triangulation, center, 1.0, 4.0, 3u, 0.0, 4, true); + + for (const auto &cell : triangulation.active_cell_iterators()) + { + deallog << "vertices: " << cell->vertex(0) << ", " << cell->vertex(1) + << ", " << cell->vertex(2) << ", " << cell->vertex(3) + << std::endl; + + bool manifold_ids_are_zero = cell->manifold_id() == 0; + for (unsigned int face_n = 0; face_n < GeometryInfo<2>::faces_per_cell; + ++face_n) + manifold_ids_are_zero &= cell->face(face_n)->manifold_id() == 0; + AssertThrow(manifold_ids_are_zero, ExcInternalError()); + + for (unsigned int face_n = 0; face_n < GeometryInfo<2>::faces_per_cell; + ++face_n) + if (cell->face(face_n)->at_boundary()) + deallog << "boundary face center distance to origin: " + << (cell->face(face_n)->center(true) - center).norm() + << " boundary id: " << cell->face(face_n)->boundary_id() + << std::endl; + } + } + + // 3D test + { + Triangulation<3> triangulation; + const Point<3> center(0.0, 1.0, 2.0); + GridGenerator::concentric_hyper_shells( + triangulation, center, 1.0, 20.0, 2u, 0.0, 0, true); + + for (const auto &cell : triangulation.active_cell_iterators()) + { + deallog << "vertices: " << cell->vertex(0) << ", " << cell->vertex(1) + << ", " << cell->vertex(2) << ", " << cell->vertex(3) << ", " + << cell->vertex(4) << ", " << cell->vertex(5) << ", " + << cell->vertex(6) << ", " << cell->vertex(7) << std::endl; + + bool manifold_ids_are_zero = cell->manifold_id() == 0; + for (unsigned int face_n = 0; face_n < GeometryInfo<3>::faces_per_cell; + ++face_n) + manifold_ids_are_zero &= cell->face(face_n)->manifold_id() == 0; + AssertThrow(manifold_ids_are_zero, ExcInternalError()); + + for (unsigned int face_n = 0; face_n < GeometryInfo<3>::faces_per_cell; + ++face_n) + if (cell->face(face_n)->at_boundary()) + { + deallog << "boundary face center distance to origin: " + << (cell->face(face_n)->center(true) - center).norm() + << " boundary id: " << cell->face(face_n)->boundary_id() + << std::endl; + } + } + } +} diff --git a/tests/grid/concentric_hyper_shell_02.output b/tests/grid/concentric_hyper_shell_02.output new file mode 100644 index 0000000000..ac78edd140 --- /dev/null +++ b/tests/grid/concentric_hyper_shell_02.output @@ -0,0 +1,69 @@ + +DEAL::vertices: 4.00000 1.00000, 2.44929e-16 5.00000, 3.00000 1.00000, 1.83697e-16 4.00000 +DEAL::boundary face center distance to origin: 4.00000 boundary id: 1 +DEAL::vertices: 2.44929e-16 5.00000, -4.00000 1.00000, 1.83697e-16 4.00000, -3.00000 1.00000 +DEAL::boundary face center distance to origin: 4.00000 boundary id: 1 +DEAL::vertices: -4.00000 1.00000, -7.34788e-16 -3.00000, -3.00000 1.00000, -5.51091e-16 -2.00000 +DEAL::boundary face center distance to origin: 4.00000 boundary id: 1 +DEAL::vertices: -7.34788e-16 -3.00000, 4.00000 1.00000, -5.51091e-16 -2.00000, 3.00000 1.00000 +DEAL::boundary face center distance to origin: 4.00000 boundary id: 1 +DEAL::vertices: 3.00000 1.00000, 1.83697e-16 4.00000, 2.00000 1.00000, 1.22465e-16 3.00000 +DEAL::vertices: 1.83697e-16 4.00000, -3.00000 1.00000, 1.22465e-16 3.00000, -2.00000 1.00000 +DEAL::vertices: -3.00000 1.00000, -5.51091e-16 -2.00000, -2.00000 1.00000, -3.67394e-16 -1.00000 +DEAL::vertices: -5.51091e-16 -2.00000, 3.00000 1.00000, -3.67394e-16 -1.00000, 2.00000 1.00000 +DEAL::vertices: 2.00000 1.00000, 1.22465e-16 3.00000, 1.00000 1.00000, 6.12323e-17 2.00000 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 1.22465e-16 3.00000, -2.00000 1.00000, 6.12323e-17 2.00000, -1.00000 1.00000 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: -2.00000 1.00000, -3.67394e-16 -1.00000, -1.00000 1.00000, -1.83697e-16 0.00000 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: -3.67394e-16 -1.00000, 2.00000 1.00000, -1.83697e-16 0.00000, 1.00000 1.00000 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 0.00000 -9.50000 2.00000, -6.06218 -5.06218 8.06218, -6.06218 -5.06218 -4.06218, -10.5000 1.00000 2.00000, 0.00000 -19.0000 2.00000, -11.5470 -10.5470 13.5470, -11.5470 -10.5470 -9.54701, -20.0000 1.00000 2.00000 +DEAL::boundary face center distance to origin: 20.0000 boundary id: 1 +DEAL::vertices: -6.06218 -5.06218 8.06218, 0.00000 1.00000 12.5000, -10.5000 1.00000 2.00000, -6.06218 7.06218 8.06218, -11.5470 -10.5470 13.5470, 0.00000 1.00000 22.0000, -20.0000 1.00000 2.00000, -11.5470 12.5470 13.5470 +DEAL::boundary face center distance to origin: 20.0000 boundary id: 1 +DEAL::vertices: 0.00000 -9.50000 2.00000, 6.06218 -5.06218 8.06218, -6.06218 -5.06218 8.06218, 0.00000 1.00000 12.5000, 0.00000 -19.0000 2.00000, 11.5470 -10.5470 13.5470, -11.5470 -10.5470 13.5470, 0.00000 1.00000 22.0000 +DEAL::boundary face center distance to origin: 20.0000 boundary id: 1 +DEAL::vertices: 6.06218 -5.06218 -4.06218, 10.5000 1.00000 2.00000, 0.00000 -9.50000 2.00000, 6.06218 -5.06218 8.06218, 11.5470 -10.5470 -9.54701, 20.0000 1.00000 2.00000, 0.00000 -19.0000 2.00000, 11.5470 -10.5470 13.5470 +DEAL::boundary face center distance to origin: 20.0000 boundary id: 1 +DEAL::vertices: 10.5000 1.00000 2.00000, 6.06218 7.06218 8.06218, 6.06218 -5.06218 8.06218, 0.00000 1.00000 12.5000, 20.0000 1.00000 2.00000, 11.5470 12.5470 13.5470, 11.5470 -10.5470 13.5470, 0.00000 1.00000 22.0000 +DEAL::boundary face center distance to origin: 20.0000 boundary id: 1 +DEAL::vertices: 6.06218 7.06218 8.06218, 0.00000 11.5000 2.00000, 0.00000 1.00000 12.5000, -6.06218 7.06218 8.06218, 11.5470 12.5470 13.5470, 0.00000 21.0000 2.00000, 0.00000 1.00000 22.0000, -11.5470 12.5470 13.5470 +DEAL::boundary face center distance to origin: 20.0000 boundary id: 1 +DEAL::vertices: 10.5000 1.00000 2.00000, 6.06218 7.06218 -4.06218, 6.06218 7.06218 8.06218, 0.00000 11.5000 2.00000, 20.0000 1.00000 2.00000, 11.5470 12.5470 -9.54701, 11.5470 12.5470 13.5470, 0.00000 21.0000 2.00000 +DEAL::boundary face center distance to origin: 20.0000 boundary id: 1 +DEAL::vertices: 6.06218 -5.06218 -4.06218, 0.00000 1.00000 -8.50000, 10.5000 1.00000 2.00000, 6.06218 7.06218 -4.06218, 11.5470 -10.5470 -9.54701, 0.00000 1.00000 -18.0000, 20.0000 1.00000 2.00000, 11.5470 12.5470 -9.54701 +DEAL::boundary face center distance to origin: 20.0000 boundary id: 1 +DEAL::vertices: 0.00000 1.00000 -8.50000, -6.06218 7.06218 -4.06218, 6.06218 7.06218 -4.06218, 0.00000 11.5000 2.00000, 0.00000 1.00000 -18.0000, -11.5470 12.5470 -9.54701, 11.5470 12.5470 -9.54701, 0.00000 21.0000 2.00000 +DEAL::boundary face center distance to origin: 20.0000 boundary id: 1 +DEAL::vertices: -6.06218 7.06218 -4.06218, -10.5000 1.00000 2.00000, 0.00000 11.5000 2.00000, -6.06218 7.06218 8.06218, -11.5470 12.5470 -9.54701, -20.0000 1.00000 2.00000, 0.00000 21.0000 2.00000, -11.5470 12.5470 13.5470 +DEAL::boundary face center distance to origin: 20.0000 boundary id: 1 +DEAL::vertices: 0.00000 1.00000 -8.50000, -6.06218 -5.06218 -4.06218, -6.06218 7.06218 -4.06218, -10.5000 1.00000 2.00000, 0.00000 1.00000 -18.0000, -11.5470 -10.5470 -9.54701, -11.5470 12.5470 -9.54701, -20.0000 1.00000 2.00000 +DEAL::boundary face center distance to origin: 20.0000 boundary id: 1 +DEAL::vertices: 6.06218 -5.06218 -4.06218, 0.00000 -9.50000 2.00000, 0.00000 1.00000 -8.50000, -6.06218 -5.06218 -4.06218, 11.5470 -10.5470 -9.54701, 0.00000 -19.0000 2.00000, 0.00000 1.00000 -18.0000, -11.5470 -10.5470 -9.54701 +DEAL::boundary face center distance to origin: 20.0000 boundary id: 1 +DEAL::vertices: 0.00000 0.00000 2.00000, -0.577350 0.422650 2.57735, -0.577350 0.422650 1.42265, -1.00000 1.00000 2.00000, 0.00000 -9.50000 2.00000, -6.06218 -5.06218 8.06218, -6.06218 -5.06218 -4.06218, -10.5000 1.00000 2.00000 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: -0.577350 0.422650 2.57735, 0.00000 1.00000 3.00000, -1.00000 1.00000 2.00000, -0.577350 1.57735 2.57735, -6.06218 -5.06218 8.06218, 0.00000 1.00000 12.5000, -10.5000 1.00000 2.00000, -6.06218 7.06218 8.06218 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 0.00000 0.00000 2.00000, 0.577350 0.422650 2.57735, -0.577350 0.422650 2.57735, 0.00000 1.00000 3.00000, 0.00000 -9.50000 2.00000, 6.06218 -5.06218 8.06218, -6.06218 -5.06218 8.06218, 0.00000 1.00000 12.5000 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 0.577350 0.422650 1.42265, 1.00000 1.00000 2.00000, 0.00000 0.00000 2.00000, 0.577350 0.422650 2.57735, 6.06218 -5.06218 -4.06218, 10.5000 1.00000 2.00000, 0.00000 -9.50000 2.00000, 6.06218 -5.06218 8.06218 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 1.00000 1.00000 2.00000, 0.577350 1.57735 2.57735, 0.577350 0.422650 2.57735, 0.00000 1.00000 3.00000, 10.5000 1.00000 2.00000, 6.06218 7.06218 8.06218, 6.06218 -5.06218 8.06218, 0.00000 1.00000 12.5000 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 0.577350 1.57735 2.57735, 0.00000 2.00000 2.00000, 0.00000 1.00000 3.00000, -0.577350 1.57735 2.57735, 6.06218 7.06218 8.06218, 0.00000 11.5000 2.00000, 0.00000 1.00000 12.5000, -6.06218 7.06218 8.06218 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 1.00000 1.00000 2.00000, 0.577350 1.57735 1.42265, 0.577350 1.57735 2.57735, 0.00000 2.00000 2.00000, 10.5000 1.00000 2.00000, 6.06218 7.06218 -4.06218, 6.06218 7.06218 8.06218, 0.00000 11.5000 2.00000 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 0.577350 0.422650 1.42265, 0.00000 1.00000 1.00000, 1.00000 1.00000 2.00000, 0.577350 1.57735 1.42265, 6.06218 -5.06218 -4.06218, 0.00000 1.00000 -8.50000, 10.5000 1.00000 2.00000, 6.06218 7.06218 -4.06218 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 0.00000 1.00000 1.00000, -0.577350 1.57735 1.42265, 0.577350 1.57735 1.42265, 0.00000 2.00000 2.00000, 0.00000 1.00000 -8.50000, -6.06218 7.06218 -4.06218, 6.06218 7.06218 -4.06218, 0.00000 11.5000 2.00000 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: -0.577350 1.57735 1.42265, -1.00000 1.00000 2.00000, 0.00000 2.00000 2.00000, -0.577350 1.57735 2.57735, -6.06218 7.06218 -4.06218, -10.5000 1.00000 2.00000, 0.00000 11.5000 2.00000, -6.06218 7.06218 8.06218 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 0.00000 1.00000 1.00000, -0.577350 0.422650 1.42265, -0.577350 1.57735 1.42265, -1.00000 1.00000 2.00000, 0.00000 1.00000 -8.50000, -6.06218 -5.06218 -4.06218, -6.06218 7.06218 -4.06218, -10.5000 1.00000 2.00000 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 +DEAL::vertices: 0.577350 0.422650 1.42265, 0.00000 0.00000 2.00000, 0.00000 1.00000 1.00000, -0.577350 0.422650 1.42265, 6.06218 -5.06218 -4.06218, 0.00000 -9.50000 2.00000, 0.00000 1.00000 -8.50000, -6.06218 -5.06218 -4.06218 +DEAL::boundary face center distance to origin: 1.00000 boundary id: 0 -- 2.39.5