]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add more documentation.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 25 Feb 2013 04:21:05 +0000 (04:21 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 25 Feb 2013 04:21:05 +0000 (04:21 +0000)
git-svn-id: https://svn.dealii.org/trunk@28543 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/examples/step-49/step-49.cc

index e8b9937687d7341a2eb0aefc83bce55b11977bfc..499aff5b378fd30f7dca81d003daabc3f16833cc 100644 (file)
@@ -1,4 +1,4 @@
-/* $Id:  $
+/* $Id$
  *
  * Copyright (C) 2013 by the deal.II authors
  *
@@ -8,6 +8,15 @@
  * further information on this license.
  */
 
+
+// This tutorial program is odd in the sense that, unlike for most other
+// steps, the introduction already provides most of the information on how to
+// use the various strategies to generate meshes. Consequently, there is
+// little that remains to be commented on here, and we intersperse the code
+// with relatively little text. In essence, the code here simply provides a
+// reference implementation of what has already been described in the
+// introduction.
+
 // @sect3{Include files}
 
 #include <deal.II/grid/tria.h>
 
 using namespace dealii;
 
-// @sect3{Generate some output for a given mesh}
-
+// @sect3{Generating output for a given mesh}
+
+// The following function generates some output for any of the meshes we will
+// be generating in the remainder of this program. In particular, it generates
+// the following information:
+//
+// - Some general information about the number of space dimensions in which
+//   this mesh lives and its number of cells.
+// - The number of boundary faces that use each boundary indicator, so that
+//   it can be compared with what we expect.
+//
+// Finally, the function outputs the mesh in encapsulated postscript (EPS)
+// format that can easily be visualized in the same way as was done in step-1.
 template<int dim>
-void mesh_info(Triangulation<dim> &tria, const char *filename)
+void mesh_info(const Triangulation<dim> &tria,
+              const std::string        &filename)
 {
-  // some general info:
   std::cout << "Mesh info:" << std::endl
             << " dimension: " << dim << std::endl
             << " no. of cells: " << tria.n_active_cells() << std::endl;
 
-  // loop over all the cells and find how often each boundary indicator is used:
+  // Next loop over all faces of all cells and find how often each boundary
+  // indicator is used:
   {
     std::map<unsigned int, unsigned int> boundary_count;
     typename Triangulation<dim>::active_cell_iterator
@@ -60,9 +81,9 @@ void mesh_info(Triangulation<dim> &tria, const char *filename)
     std::cout  << std::endl;
   }
 
-  // Now we want to write a graphical representation of the mesh to an output
-  // file.
-  std::ofstream out (filename);
+  // Finally, produce a graphical representation of the mesh to an output
+  // file:
+  std::ofstream out (filename.c_str());
   GridOut grid_out;
   grid_out.write_eps (tria, out);
   std::cout << " written to " << filename << std::endl;
@@ -70,8 +91,12 @@ void mesh_info(Triangulation<dim> &tria, const char *filename)
 
 // @sect3{Main routines}
 
+// @sect4{grid_1: Loading a mesh generated by gmsh}
 
-// load from a msh generated by gmsh
+// In this first example, we show how to load the mesh for which we have
+// discussed in the introduction how to generate it. This follows the same
+// pattern as used in step-5 to load a mesh, although there it was written in
+// a different file format (UCD instead of MSH).
 void grid_1 ()
 {
   Triangulation<2> triangulation;
@@ -86,28 +111,52 @@ void grid_1 ()
 }
 
 
-// merge triangulations
-// note: the vertices have to be exactly on top of each other!
+// @sect4{grid_2: Merging triangulations}
+
+// Here, we first create two triangulations and then merge them into one.  As
+// discussed in the introduction, it is important to ensure that the vertices
+// at the common interface are located at the same coordinates.
 void grid_2 ()
 {
-  Triangulation<2> triangulation;
-
-  Triangulation<2> tria1, tria2;
-
-
+  Triangulation<2> tria1;
   GridGenerator::hyper_cube_with_cylindrical_hole (tria1, 0.25, 1.0);
+
+  Triangulation<2> tria2;
   std::vector< unsigned int > repetitions(2);
   repetitions[0]=3;
   repetitions[1]=2;
   GridGenerator::subdivided_hyper_rectangle (tria2, repetitions,
-                                             Point<2>(1.0,-1.0), Point<2>(4.0,1.0));
+                                             Point<2>(1.0,-1.0),
+                                            Point<2>(4.0,1.0));
 
+  Triangulation<2> triangulation;
   GridGenerator::merge_triangulations (tria1, tria2, triangulation);
 
   mesh_info(triangulation, "grid-2.eps");
 }
 
-// move vertices
+
+// @sect4{grid_3: Moving vertices}
+
+// In this function, we move vertices of a mesh. This is simpler than one
+// usually expects: if you ask a cell using <code>cell-@>vertex(i)</code> for
+// the coordinates of its <code>i</code>th vertex, it doesn't just provide the
+// location of this vertex but in fact a reference to the location where these
+// coordinates are stored. We can then modify the value stored there.
+//
+// So this is what we do in the first part of this function: We create a
+// square of geometry $[-1,1]^2$ with a circular hole with radius 0.25 located
+// at the origin. We then loop over all cells and all vertices and if a vertex
+// has a $y$ coordinate equal to one, we move it upward by 0.5.
+//
+// Note that this sort of procedure does not usually work this way because one
+// will typically encounter the same vertices multiple times and may move them
+// more than once. It works here because we select the vertices we want to use
+// based on their geometric location, and a vertex moved once will fail this
+// test in the future. A more general approach to this problem would have been
+// to keep a std::set of of those vertex indices that we have already moved
+// (which we can obtain using <code>cell-@>vertex_index(i)</code> and only
+// move those vertices whose index isn't in the set yet.
 void grid_3 ()
 {
   Triangulation<2> triangulation;
@@ -122,74 +171,110 @@ void grid_3 ()
         {
           Point<2> &v = cell->vertex(i);
           if (std::abs(v(1)-1.0)<1e-5)
-            v(1)+=0.5;
+            v(1) += 0.5;
         }
-
     }
 
-  // here we are going to set a boundary descriptor for the round hole
-  // and refine the mesh twice.
-  const HyperBallBoundary<2> boundary_description(Point<2>(0,0),0.25);
+  // In the second step we will refine the mesh twice. To do this correctly,
+  // we have to associate a geometry object with the boundary of the hole;
+  // since the boundary of the hole has boundary indicator 1 (see the
+  // documentation of the function that generates the mesh), we need to create
+  // an object that describes a circle (i.e., a hyper ball) with appropriate
+  // center and radius and assign it to the triangulation. We can then refine
+  // twice:
+  const HyperBallBoundary<2> boundary_description(Point<2>(0,0), 0.25);
   triangulation.set_boundary (1, boundary_description);
   triangulation.refine_global(2);
 
-  mesh_info(triangulation, "grid-3.eps");
-
-  // remove boundary object from Triangulation again:
+  // The mesh so generated is then passed to the function that generates
+  // output. In a final step we remove the boundary object again so that it is
+  // no longer in use by the triangulation when it is destroyed (the boundary
+  // object is destroyed first in this function since it was declared after
+  // the triangulation).
+  mesh_info (triangulation, "grid-3.eps");
   triangulation.set_boundary (1);
 }
 
-// demonstrate extrude_triangulation
+
+// @sect4{grid_4: Demonstrating extrude_triangulation}
+
+// This example takes the initial grid from the previous function and simply extrudes it into the third space dimension:
 void grid_4()
 {
   Triangulation<2> triangulation;
   Triangulation<3> out;
   GridGenerator::hyper_cube_with_cylindrical_hole (triangulation, 0.25, 1.0);
 
-  GridGenerator::extrude_triangulation(triangulation, 3, 2.0, out);
+  GridGenerator::extrude_triangulation (triangulation, 3, 2.0, out);
   mesh_info(out, "grid-4.eps");
-  
+
 }
 
+
+// @sect4{grid_5: Demonstrating GridTools::transform, part 1}
+
+// This and the next example first create a mesh and then transform it by
+// moving every node of the mesh according to a function that takes a point
+// and returns a mapped point. In this case, we transform $(x,y) \mapsto
+// (x,y+\sin(\pi x/5))$.
+//
+// GridTools::transform takes a triangulation and any kind of object that can
+// be called like a function as arguments. This function-like argument can be
+// simply the address of a function as in the current case, or an object that
+// has an <code>operator()</code> as in the next example, or for example a
+// <code>std::function@<Point@<2@>(const Point@<2@>)</code> object one can get
+// via <code>std::bind</code> in more complex cases.
 struct Grid5Func
 {
     Point<2> operator() (const Point<2> & in) const
       {
-       return Point<2>(in(0), in(1)+sin(in(0)/5.0*3.14159));
-      }      
+       return Point<2>(in(0),
+                       in(1) + std::sin(in(0)/5.0*3.14159));
+      }
 };
 
 
-// demonstrate GridTools::transform
 void grid_5()
 {
   Triangulation<2> tria;
-  std::vector< unsigned int > repetitions(2);
-  repetitions[0]=14;
-  repetitions[1]=2;
+  std::vector<unsigned int> repetitions(2);
+  repetitions[0] = 14;
+  repetitions[1] = 2;
   GridGenerator::subdivided_hyper_rectangle (tria, repetitions,
-                                             Point<2>(0.0,0.0), Point<2>(10.0,1.0));
+                                             Point<2>(0.0,0.0),
+                                            Point<2>(10.0,1.0));
 
-  
   GridTools::transform(Grid5Func(), tria);
-  mesh_info(tria, "grid-5.eps");  
+  mesh_info(tria, "grid-5.eps");
 }
 
+
+
+// @sect4{grid_6: Demonstrating GridTools::transform, part 2}
+
+// In this second example of transforming points from an original to a new
+// mesh, we will use the mapping $(x,y) \mapsto (x,\tanh(2y)/\tanh(2))$. To
+// make things more interesting, rather than doing so in a single function as
+// in the previous example, we here create an object with an
+// <code>operator()</code> that will be called by GridTools::transform. Of
+// course, this object may in reality be much more complex: the object may
+// have member variables that play a role in computing the new locations of
+// vertices.
 struct Grid6Func
 {
-    double trans(double x) const
+    double trans(const double x) const
       {
-       return tanh(2*x)/tanh(2);
+       return std::tanh(2*x)/tanh(2);
       }
-    
+
     Point<2> operator() (const Point<2> & in) const
       {
-       return Point<2>((in(0)), trans(in(1)));
-      }      
+       return Point<2> (in(0),
+                        trans(in(1)));
+      }
 };
-  
 
-// demonstrate GridTools::transform
+
 void grid_6()
 {
   Triangulation<2> tria;
@@ -200,25 +285,32 @@ void grid_6()
                                              Point<2>(0.0,0.0), Point<2>(1.0,1.0));
 
   GridTools::transform(Grid6Func(), tria);
-  mesh_info(tria, "grid-6.eps");  
+  mesh_info(tria, "grid-6.eps");
 }
 
 
+// @sect4{grid_7: Demonstrating distort_random}
 
-//demonstrate distort_random
+// In this last example, we create a mesh and then distort its (interior)
+// vertices by a random perturbation. This is not something you want to do for
+// production computations, but it is a useful tool for testing
+// discretizations and codes to make sure they don't work just by accident
+// because the mesh happens to be uniformly structured and supporting
+// super-convergence properties.
 void grid_7()
 {
   Triangulation<2> tria;
-  std::vector< unsigned int > repetitions(2);
-  repetitions[0]=16;
-  repetitions[1]=16;
+  std::vector<unsigned int> repetitions(2);
+  repetitions[0] = repetitions[1] = 16;
   GridGenerator::subdivided_hyper_rectangle (tria, repetitions,
-                                             Point<2>(0.0,0.0), Point<2>(1.0,1.0));
+                                             Point<2>(0.0,0.0),
+                                            Point<2>(1.0,1.0));
 
   tria.distort_random(0.3, true);
-  mesh_info(tria, "grid-7.eps");  
+  mesh_info(tria, "grid-7.eps");
 }
 
+
 // @sect3{The main function}
 
 // Finally, the main function. There isn't much to do here, only to call the

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.