]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Mention range-based for loops in the tutorials. 79/head
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Thu, 14 Aug 2014 13:28:08 +0000 (08:28 -0500)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Thu, 14 Aug 2014 13:45:29 +0000 (08:45 -0500)
examples/step-1/step-1.cc
examples/step-3/step-3.cc

index 3653accf2415d809566fcd2f68c75614bf2879a3..04e8fb1ca904fdb52a4527cba80052693c5dc6ef 100644 (file)
@@ -130,51 +130,80 @@ void second_grid ()
   // refine the grid in five steps towards the inner circle of the domain:
   for (unsigned int step=0; step<5; ++step)
     {
-      // Next, we need an iterator which points to a cell and which we will
-      // move over all active cells one by one (active cells are those that
-      // are not further refined, and the only ones that can be marked for
-      // further refinement, obviously). By convention, we almost always use
-      // the names <code>cell</code> and <code>endc</code> for the iterator
-      // pointing to the present cell and to the <code>one-past-the-end</code>
-      // iterator:
+      // Next, we need an iterator that points to a cell and which we will
+      // move over all active cells one by one. In a sense, you can think of a
+      // triangulation as a collection of cells. If it was an array, you would
+      // just get a pointer that you move from one to the next. In
+      // triangulations, cells aren't stored as an array, so simple pointers
+      // do not work, but one can generalize pointers to iterators (see <a
+      // href="http://en.wikipedia.org/wiki/Iterator#C.2B.2B">this wikipedia
+      // link</a> for more information). We will then get an iterator to the
+      // first cell and iterate over all of the cells until we hit the last
+      // one.
+      //
+      // The second important piece is that we only need the active cells.
+      // Active cells are those that are not further refined, and the only
+      // ones that can be marked for further refinement, obviously. deal.II
+      // provides iterator categories that allow us to iterate over <i>all</i>
+      // cells (including the parent cells of active ones) or only over the
+      // active cells. Because we want the latter, we need to choose
+      // Triangulation::active_cell_iterator as data type.
+      //
+      // Finally, by convention, we almost always use the names
+      // <code>cell</code> and <code>endc</code> for the iterator pointing to
+      // the present cell and to the "one-past-the-end" iterator. This is, in
+      // a sense a misnomer, because the object is not really a "cell": it is
+      // an iterator/pointer to a cell. We should really have started to call
+      // these objects <code>cell_iterator</code> when deal.II started in
+      // 1998, but it is what it is.
+      //
+      // After declaring the iterator variable, the loop over all cells is
+      // then rather trivial, and looks like any loop involving pointers
+      // instead of iterators:
       Triangulation<2>::active_cell_iterator
       cell = triangulation.begin_active(),
       endc = triangulation.end();
-
-      // The loop over all cells is then rather trivial, and looks like any
-      // loop involving pointers instead of iterators:
       for (; cell!=endc; ++cell)
-        // Next, we want to loop over all vertices of the cells. Since we are
-        // in 2d, we know that each cell has exactly four vertices. However,
-        // instead of penning down a 4 in the loop bound, we make a first
-        // attempt at writing it in a dimension-independent way by which we
-        // find out about the number of vertices of a cell. Using the
-        // GeometryInfo class, we will later have an easier time getting the
-        // program to also run in 3d: we only have to change all occurrences
-        // of <code>&lt;2&gt;</code> to <code>&lt;3&gt;</code>, and do not
-        // have to audit our code for the hidden appearance of magic numbers
-        // like a 4 that needs to be replaced by an 8:
-        for (unsigned int v=0;
-             v < GeometryInfo<2>::vertices_per_cell;
-             ++v)
-          {
-            // If this cell is at the inner boundary, then at least one of its
-            // vertices must sit on the inner ring and therefore have a radial
-            // distance from the center of exactly 0.5, up to floating point
-            // accuracy. Compute this distance, and if we have found a vertex
-            // with this property flag this cell for later refinement. We can
-            // then also break the loop over all vertices and move on to the
-            // next cell.
-            const double distance_from_center
-              = center.distance (cell->vertex(v));
-
-            if (std::fabs(distance_from_center - inner_radius) < 1e-10)
-              {
-                cell->set_refine_flag ();
-                break;
-              }
-          }
-
+       {
+         // @note Writing a loop like this requires a lot of typing, but it
+         // is the only way of doing it in C++98 and C++03. However, if you
+         // have a C++11-compliant compiler, you can also use the C++11
+         // range-based for loop style that requires significantly less
+         // typing. Take a look at @ref CPP11 "the deal.II C++11 page" to see
+         // how this works.
+         //
+         // Next, we want to loop over all vertices of the cells. Since we are
+         // in 2d, we know that each cell has exactly four vertices. However,
+         // instead of penning down a 4 in the loop bound, we make a first
+         // attempt at writing it in a dimension-independent way by which we
+         // find out about the number of vertices of a cell. Using the
+         // GeometryInfo class, we will later have an easier time getting the
+         // program to also run in 3d: we only have to change all occurrences
+         // of <code>&lt;2&gt;</code> to <code>&lt;3&gt;</code>, and do not
+         // have to audit our code for the hidden appearance of magic numbers
+         // like a 4 that needs to be replaced by an 8:
+         for (unsigned int v=0;
+              v < GeometryInfo<2>::vertices_per_cell;
+              ++v)
+           {
+             // If this cell is at the inner boundary, then at least one of its
+             // vertices must sit on the inner ring and therefore have a radial
+             // distance from the center of exactly 0.5, up to floating point
+             // accuracy. Compute this distance, and if we have found a vertex
+             // with this property flag this cell for later refinement. We can
+             // then also break the loop over all vertices and move on to the
+             // next cell.
+             const double distance_from_center
+               = center.distance (cell->vertex(v));
+
+             if (std::fabs(distance_from_center - inner_radius) < 1e-10)
+               {
+                 cell->set_refine_flag ();
+                 break;
+               }
+           }
+       }
+      
       // Now that we have marked all the cells that we want refined, we let
       // the triangulation actually do this refinement. The function that does
       // so owes its long name to the fact that one can also mark cells for
index c692a1f58f3faee3715832e4410741ae982315af..0ca76a4654d55fa738a87087ed31012bf4bec186 100644 (file)
@@ -356,13 +356,18 @@ void Step3::assemble_system ()
   std::vector<types::global_dof_index> local_dof_indices (dofs_per_cell);
 
   // Now for the loop over all cells. We have seen before how this works, so
-  // this should be familiar including the conventional names for these
-  // variables:
+  // the following code should be familiar including the conventional names
+  // for these variables:
   DoFHandler<2>::active_cell_iterator
   cell = dof_handler.begin_active(),
   endc = dof_handler.end();
   for (; cell!=endc; ++cell)
     {
+      // @note As already mentioned in step-1, there is a more convenient way
+      // of writing such loops if your compiler supports the C++11
+      // standard. See @ref CPP11 "the deal.II C++11 page" to see
+      // how this works.
+      //
       // We are now sitting on one cell, and we would like the values and
       // gradients of the shape functions be computed, as well as the
       // determinants of the Jacobian matrices of the mapping between

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.