]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Redo some documentation to use manifolds.
authorDavid Wells <wellsd2@rpi.edu>
Thu, 5 Apr 2018 21:39:50 +0000 (17:39 -0400)
committerDavid Wells <wellsd2@rpi.edu>
Fri, 20 Apr 2018 21:09:06 +0000 (17:09 -0400)
We should describe things in terms of the Manifold classes since the
Boundary classes are deprecated.

doc/doxygen/headers/glossary.h
doc/doxygen/headers/manifold.h
examples/step-34/doc/intro.dox
include/deal.II/grid/tria.h
include/deal.II/grid/tria_accessor.h

index bac9d73eb6d0a89bb2edcb18375a1553e678eff4..1c9834c868c0d2639e91b87125284103e3e6199d 100644 (file)
  * particular boundary indicator. This method is still supported, and
  * it allows the Triangulation object to use a different method of
  * finding new points on faces and edges to be refined; the default is
- * to use a StraightBoundary object for all faces and edges. The
+ * to use a FlatManifold object for all faces and edges. The
  * results section of step-49 has a worked example that shows all of
  * this in action.
  *
index ca4831509839c709887d571081d1b880d2c1c675..23c3c2645a2eba64ba9f34f017268ad6c0da50d1 100644 (file)
  * reasonable implementations. More complicated examples can be described
  * using the techniques shown in step-53 and step-54.
  *
- * The boundary of a Triangulation is a special case of Manifold, for
- * which additional information can be useful in user codes, such as
- * normal vectors to surfaces or to curves. If your coarse mesh is reasonably
- * shaped, you might be interested in only attaching a manifold
- * description to boundary portion of your domain. This can be done
- * using the Triangulation::set_boundary() function, which take as arguments a Boundary
- * object (derived from Manifold). Notice that Triangulation uses only
- * the Manifold interface, not the Boundary interface. Other tools,
- * however, might need to compute exact normals at quadrature points,
- * and therefore a wrapper to query Boundary objects is provided.
- *
- *
  * <h3>An example</h3>
  *
  * A simple example why dealing with curved geometries is already provided
  *  GridGenerator::hyper_shell (triangulation,
  *                              center, inner_radius, outer_radius,
  *                              3);    // three circumferential cells
- *  const HyperShellBoundary<2> boundary_description(center);
+ *  const SphericalManifold<2> boundary_description(center);
  *  triangulation.set_all_manifold_ids_on_boundary(0);
  *  triangulation.set_manifold (0, boundary_description);
  *
index 7fc4aa661ae35e1e45dcd437ce5255018de3045d..c3c135c1d2842438bb09200f54eb7d7842b0a553 100644 (file)
@@ -670,7 +670,7 @@ derivatives.
 
 The testcase we will be solving is for a circular (in 2d) or spherical
 (in 3d) obstacle. Meshes for these geometries will be read in from
-files in the current directory and an object of type HyperBallBoundary
+files in the current directory and an object of type SphericalManifold
 will then be attached to the triangulation to allow mesh refinement
 that respects the continuous geometry behind the discrete initial
 mesh.
index c6a0eea3dafdc2167f4df6ede47fd92ddb945605..0759fd3169e333436ce1111e8bebb407c88c4498 100644 (file)
@@ -25,7 +25,6 @@
 #include <deal.II/base/iterator_range.h>
 #include <deal.II/grid/tria_iterator_selector.h>
 
-// Ignore deprecation warnings for auto_ptr.
 #include <boost/signals2.hpp>
 #include <boost/serialization/vector.hpp>
 #include <boost/serialization/map.hpp>
@@ -946,54 +945,81 @@ namespace internal
  * been called.
  *
  *
- * <h3>%Boundary approximation</h3>
+ * <h3>Describing curved geometries</h3>
  *
- * You can specify a class describing the boundary for each boundary component.
- * If a new vertex is created on a side or face at the boundary, this class is
- * used to compute where it will be placed. The manifold indicator of the face
- * will be used to determine the proper component. See Manifold for the details.
- * Usage with the Triangulation object is then like this (let @p Ball be a
- * class derived from Manifold<tt><2></tt>):
+ * deal.II implements all geometries (curved and otherwise) with classes
+ * inheriting from Manifold; see the documentation of Manifold, step-49, or
+ * the
+ * @ref manifold
+ * module for examples and a complete description of the algorithms. By
+ * default, all cells in a Triangulation have a flat geometry, meaning that
+ * all lines in the Triangulation are assumed to be straight. If a cell has a
+ * manifold_id that is not equal to numbers::flat_manifold_id then the
+ * Triangulation uses the associated Manifold object for computations on that
+ * cell (e.g., cell refinement). Here is a quick example, taken from the
+ * implementation of GridGenerator::hyper_ball(), that sets up a polar grid:
  *
- *   @code
- *     int main ()
- *     {
- *       Triangulation<2> tria;
- *                                        // set the boundary function
- *                                        // for all boundaries with
- *                                        // boundary indicator 0
- *       Ball ball;
- *       tria.set_manifold (0, ball);
- *
- *   // read some coarse grid
- *
- *
- *   Triangulation<2>::active_cell_iterator cell, endc;
- *   for (int i=0; i<8; ++i)
+ * @code
+ * int main ()
+ * {
+ *   Triangulation<2> triangulation;
+ *   const Point<2> vertices[8] = {Point<2>(-1,-1),
+ *                                 Point<2>(+1,-1),
+ *                                 Point<2>(-1,-1)*0.5,
+ *                                 Point<2>(+1,-1)*0.5,
+ *                                 Point<2>(-1,+1)*0.5,
+ *                                 Point<2>(+1,+1)*0.5,
+ *                                 Point<2>(-1,+1),
+ *                                 Point<2>(+1,+1)};
+ *   const int cell_vertices[5][4] = {{0, 1, 2, 3},
+ *                                    {0, 2, 6, 4},
+ *                                    {2, 3, 4, 5},
+ *                                    {1, 7, 3, 5},
+ *                                    {6, 4, 7, 5}};
+ *
+ *   std::vector<CellData<2>> cells(5, CellData<2>());
+ *   for (unsigned int i=0; i<5; ++i)
+ *     for (unsigned int j=0; j<4; ++j)
+ *       cells[i].vertices[j] = cell_vertices[i][j];
+ *
+ *   triangulation.create_triangulation ({std::begin(vertices), std::end(vertices)},
+ *                                       cells,
+ *                                       SubCellData());
+ *   triangulation.set_all_manifold_ids_on_boundary(42);
+ *   // set_manifold stores a copy of its second argument, so a temporary is okay
+ *   triangulation.set_manifold(42, PolarManifold<2>());
+ *   for (unsigned int i = 0; i < 4; ++i)
  *     {
- *       cell = tria.begin_active();
- *       endc = tria.end();
- *
  *       // refine all boundary cells
- *       for (; cell!=endc; ++cell)
+ *       for (const auto &cell : triangulation.active_cell_iterators())
  *         if (cell->at_boundary())
  *           cell->set_refine_flag();
  *
- *       tria.execute_coarsening_and_refinement();
+ *       triangulation.execute_coarsening_and_refinement();
  *     }
  * }
  * @endcode
  *
+ * This will set up a grid where the boundary lines will be refined by
+ * performing calculations in polar coordinates. When the mesh is refined the
+ * cells adjacent to the boundary will use this new line midpoint (as well as
+ * the other three midpoints and original cell vertices) to calculate the cell
+ * midpoint with a transfinite interpolation: this propagates the curved
+ * boundary into the interior in a smooth way. It is possible to generate a
+ * better grid (which interpolates across all cells between two different
+ * Manifold descriptions, instead of just going one cell at a time) by using
+ * TransfiniteInterpolationManifold; see the documentation of that class for
+ * more information.
+ *
  * You should take note of one caveat: if you have concave boundaries, you
  * must make sure that a new boundary vertex does not lie too much inside the
  * cell which is to be refined. The reason is that the center vertex is placed
- * at the point which is the arithmetic mean of the vertices of the original
- * cell. Therefore if your new boundary vertex is too near the center of the
- * old quadrilateral or hexahedron, the distance to the midpoint vertex will
- * become too small, thus generating distorted cells. This issue is discussed
- * extensively in
- * @ref GlossDistorted "distorted cells".
- *
+ * at the point which is a weighted average of the vertices of the original
+ * cell, new face midpoints, and (in 3D) new line midpoints. Therefore if your
+ * new boundary vertex is too near the center of the old quadrilateral or
+ * hexahedron, the distance to the midpoint vertex will become too small, thus
+ * generating distorted cells. This issue is discussed extensively in @ref
+ * GlossDistorted "distorted cells".
  *
  * <h3>Getting notice when a triangulation changes</h3>
  *
index c98b798e3548ae90bd9569f68a21fb45e0f4aed0..7846428484815898d756ca63b6be60370e834ae2 100644 (file)
@@ -1031,11 +1031,11 @@ public:
   /**
    * Return a constant reference to the manifold object used for this object.
    *
-   * As explained in
-   * @ref boundary "Boundary and manifold description for triangulations",
-   * the process involved in finding the appropriate manifold description
-   * involves querying both the manifold or boundary indicators. See there for
-   * more information.
+   * As explained in the
+   * @ref manifold
+   * module, the process involved in finding the appropriate manifold
+   * description involves querying both the manifold or boundary
+   * indicators. See there for more information.
    */
   const Manifold<dim,spacedim> &get_manifold () const;
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.