From: David Wells Date: Sat, 7 Apr 2018 23:21:14 +0000 (-0400) Subject: Fix some examples in the manifold module. X-Git-Tag: v9.0.0-rc1~142^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=15784ae30dfb5951bd93e70384bff63092fb9c33;p=dealii.git Fix some examples in the manifold module. The second to last example was wrong: it set manifold ids on the interior, but should not have. Otherwise, this commit: 1. uses a consistent initialization order 2. refers to SphericalManifold as 'manifold' 3. uses the relatively new Triangulation::set_all_manifold_ids function 4. explicitly remove the manifolds attached to the Triangulations so that the effect of each manifold is clear. to tidy things up. --- diff --git a/doc/doxygen/headers/manifold.h b/doc/doxygen/headers/manifold.h index 23c3c2645a..c668389315 100644 --- a/doc/doxygen/headers/manifold.h +++ b/doc/doxygen/headers/manifold.h @@ -81,20 +81,27 @@ * *

An example

* - * A simple example why dealing with curved geometries is already provided - * by step-1, though it is not elaborated there. For example, consider this + * A simple example why dealing with curved geometries is already provided by + * step-1, though it is not elaborated there. By default, the functions in + * GridGenerator will attach manifolds to meshes when needed. In each code + * snippet below we call Triangulation::reset_all_manifolds() to remove these + * manifolds and handle all Manifold attachment in the example itself to make + * the impact of the choice of Manifold clear. + * + * Consider this * small variation of the second_grid() function shown there, - * where we simply refine every cell several times and do not deal - * with boundaries at all: + * where we simply refine every cell several times: * @code - * Triangulation<2> triangulation; - * * const Point<2> center (1,0); * const double inner_radius = 0.5, * outer_radius = 1.0; + * Triangulation<2> triangulation; * GridGenerator::hyper_shell (triangulation, * center, inner_radius, outer_radius, * 10); + * // as noted above: disable all non-Cartesian manifolds + * // for demonstration purposes: + * triangulation.reset_all_manifolds(); * * triangulation.refine_global (3); * @endcode @@ -110,20 +117,20 @@ * in the middle of the existing ones, regardless of what we may have intended * (but omitted to describe in code). * - * This is easily remedied. step-1 already shows how to do this. Consider this - * code: + * This is easily remedied. Consider this code: * @code - * Triangulation<2> triangulation; - * * const Point<2> center (1,0); + * const SphericalManifold<2> manifold(center); * const double inner_radius = 0.5, * outer_radius = 1.0; + * Triangulation<2> triangulation; * GridGenerator::hyper_shell (triangulation, * center, inner_radius, outer_radius, * 10); - * const SphericalManifold<2> boundary_description(center); + * // again disable all manifolds for demonstration purposes + * triangulation.reset_all_manifolds(); * triangulation.set_all_manifold_ids_on_boundary(0); - * triangulation.set_boundary (0, boundary_description); + * triangulation.set_manifold (0, manifold); * * triangulation.refine_global (3); * @endcode @@ -146,26 +153,23 @@ * This can be remedied by assigning a manifold description not only to * the lines along the boundary, but also to the radial lines and cells (which, * in turn, will inherit it to the new lines that are created upon mesh - * refinement). This code achieves this: + * refinement). This is exactly what GridGenerator::hyper_shell() does by default. + * For demonstration purposes, we disable the default Manifold behavior and then + * duplicate it manually: * @code - * Triangulation<2> triangulation; - * * const Point<2> center (1,0); + * const SphericalManifold<2> manifold(center); * const double inner_radius = 0.5, * outer_radius = 1.0; + * Triangulation<2> triangulation; * GridGenerator::hyper_shell (triangulation, * center, inner_radius, outer_radius, * 10); - * const SphericalManifold<2> boundary_description(center); - * triangulation.set_all_manifold_ids_on_boundary(0); - * triangulation.set_manifold (0, boundary_description); - * - * Triangulation<2>::active_cell_iterator - * cell = triangulation.begin_active(), - * endc = triangulation.end(); - * for (; cell!=endc; ++cell) - * cell->set_all_manifold_ids (0); - * + * // again disable all manifolds for demonstration purposes + * triangulation.reset_all_manifolds(); + * // reenable the + * triangulation.set_all_manifold_ids(0); + * triangulation.set_manifold (0, manifold); * triangulation.refine_global (3); * @endcode * This leads to the following mesh: @@ -188,17 +192,17 @@ * not so much of an issue for the meshes shown above, but is sometimes an * issue. Consider, for example, the following code and mesh: * @code - * Triangulation<2> triangulation; - * * const Point<2> center (1,0); + * const SphericalManifold<2> manifold(center); * const double inner_radius = 0.5, * outer_radius = 1.0; + * Triangulation<2> triangulation; * GridGenerator::hyper_shell (triangulation, * center, inner_radius, outer_radius, * 3); // three circumferential cells - * const SphericalManifold<2> boundary_description(center); + * triangulation.reset_all_manifolds(); * triangulation.set_all_manifold_ids_on_boundary(0); - * triangulation.set_manifold (0, boundary_description); + * triangulation.set_manifold (0, manifold); * * triangulation.refine_global (3); * @endcode @@ -217,23 +221,17 @@ * generators), we observe the following: * * @code - * Triangulation<2> triangulation; - * * const Point<2> center (1,0); + * const SphericalManifold<2> manifold(center); * const double inner_radius = 0.8, * outer_radius = 1.0; + * Triangulation<2> triangulation; * GridGenerator::hyper_shell (triangulation, * center, inner_radius, outer_radius, * 3); // three circumferential cells - * const SphericalManifold<2> boundary_description(center); + * triangulation.reset_all_manifolds(); * triangulation.set_all_manifold_ids_on_boundary(0); - * triangulation.set_manifold (0, boundary_description); - * - * Triangulation<2>::active_cell_iterator - * cell = triangulation.begin_active(), - * endc = triangulation.end(); - * for (; cell!=endc; ++cell) - * cell->set_all_manifold_ids (0); + * triangulation.set_manifold (0, manifold); * * triangulation.refine_global (3); * @endcode @@ -247,29 +245,23 @@ * not only to the boundary but also to interior cells and edges, using * the same code as above: * @code - * Triangulation<2> triangulation; - * * const Point<2> center (1,0); * const double inner_radius = 0.8, * outer_radius = 1.0; + * Triangulation<2> triangulation; * GridGenerator::hyper_shell (triangulation, * center, inner_radius, outer_radius, * 3); // three circumferential cells - * const SphericalManifold<2> boundary_description(center); - * triangulation.set_all_manifold_ids(0); - * triangulation.set_manifold (0, boundary_description); - * - * Triangulation<2>::active_cell_iterator - * cell = triangulation.begin_active(), - * endc = triangulation.end(); - * for (; cell!=endc; ++cell) - * cell->set_all_manifold_ids (0); * * triangulation.refine_global (3); * @endcode * * @image html hypershell-all-3.png "" * + * In this last example we finally let GridGenerator do its job and we keep + * the default manifold configuration, which is a SphericalManifold on every + * cell and face. + * * Here, even starting with an initial, inappropriately chosen mesh retains * our ability to adequately refine the mesh into one that will serve us * well. This example may be manufactured here, but it is relevant, for example