]> https://gitweb.dealii.org/ - dealii.git/commitdiff
examples/step-27: Update indenting and modernize
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Sun, 1 Jul 2018 21:49:09 +0000 (23:49 +0200)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Sun, 1 Jul 2018 21:49:09 +0000 (23:49 +0200)
examples/step-27/doc/intro.dox
examples/step-27/step-27.cc

index cc1fa99ef8f64f97c0d5c4e02387a3356fc08f0b..8e9de94d7c1573d2bd85fea0577171b232f25681 100644 (file)
@@ -149,9 +149,7 @@ inactive on it. The general outline of this reads like this:
 
 @code
   hp::DoFHandler<dim> dof_handler (triangulation);
-  for (typename hp::DoFHandler<dim>::active_cell_iterator
-         cell = dof_handler.begin_active();
-       cell != dof_handler.end(); ++cell)
+  for (auto &cell: dof_handler.active_cell_iterators())
     cell->set_active_fe_index (...);
   dof_handler.distribute_dofs (fe_collection);
 @endcode
@@ -171,9 +169,8 @@ finite element field to ensure that it is continuous. This is
 conceptually very similar to how we compute hanging node constraints,
 and in fact the code looks exactly the same:
 @code
-  ConstraintMatrix constraints;
-  DoFTools::make_hanging_node_constraints (dof_handler,
-                                          constraints);
+  AffineConstraints<double> constraints;
+  DoFTools::make_hanging_node_constraints (dof_handler, constraints);
 @endcode
 In other words, the DoFTools::make_hanging_node_constraints deals not
 only with hanging node constraints, but also with $hp$ constraints at
@@ -218,10 +215,7 @@ i.e. the interesting part of the loop over all cells would look like this:
                                  update_values    |  update_gradients |
                                  update_q_points  |  update_JxW_values);
 
-  typename hp::DoFHandler<dim>::active_cell_iterator
-    cell = dof_handler.begin_active(),
-    endc = dof_handler.end();
-  for (; cell!=endc; ++cell)
+  for (const auto &cell: dof_handler.active_cell_iterators())
     {
       hp_fe_values.reinit (cell,
                            cell->active_fe_index(),
@@ -624,7 +618,7 @@ challenge.
 Most programs built on deal.II use the DoFTools::make_sparsity_pattern
 function to allocate the sparsity pattern of a matrix, and later add a few
 more entries necessary to handle constrained degrees of freedom using
-ConstraintMatrix::condense. The sparsity pattern is then compressed using
+AffineConstraints::condense. The sparsity pattern is then compressed using
 SparsityPattern::compress. This method is explained in step-6 and used in
 most tutorial programs. In order to work, it needs an initial upper estimate
 for the maximal number of nonzero entries per row, something that can be had
@@ -671,7 +665,7 @@ of freedom, but is constrained against many more degrees of freedom.
 
 It turns out that the strategy presented first in step-6 to eliminate the
 constraints while computing the element matrices and vectors with
-ConstraintMatrix::distribute_local_to_global is the most efficient approach
+AffineConstraints::distribute_local_to_global is the most efficient approach
 also for this case. The alternative strategy to first build the matrix without
 constraints and then "condensing" away constrained degrees of freedom is
 considerably more expensive. It turns out that building the sparsity pattern
@@ -685,7 +679,7 @@ In our program, we will also treat the boundary conditions as (possibly
 inhomogeneous) constraints and eliminate the matrix rows and columns to
 those as well. All we have to do for this is to call the function that
 interpolates the Dirichlet boundary conditions already in the setup phase in
-order to tell the ConstraintMatrix object about them, and then do the
+order to tell the AffineConstraints object about them, and then do the
 transfer from local to global data on matrix and vector simultaneously. This
 is exactly what we've shown in step-6.
 
index c0ea4491b7eb33c60b8087734ce0d2bf186e0dd7..81120240d1c1a903212325d59488d9bb66b5bf68 100644 (file)
@@ -115,7 +115,7 @@ namespace Step27
     std::vector<double>                     ln_k;
     Table<dim, std::complex<double>>        fourier_coefficients;
 
-    ConstraintMatrix constraints;
+    AffineConstraints<double> constraints;
 
     SparsityPattern      sparsity_pattern;
     SparseMatrix<double> system_matrix;
@@ -328,10 +328,7 @@ namespace Step27
 
     std::vector<types::global_dof_index> local_dof_indices;
 
-    typename hp::DoFHandler<dim>::active_cell_iterator cell = dof_handler
-                                                                .begin_active(),
-                                                       endc = dof_handler.end();
-    for (; cell != endc; ++cell)
+    for (const auto &cell : dof_handler.active_cell_iterators())
       {
         const unsigned int dofs_per_cell = cell->get_fe().dofs_per_cell;
 
@@ -354,11 +351,13 @@ namespace Step27
             {
               for (unsigned int j = 0; j < dofs_per_cell; ++j)
                 cell_matrix(i, j) +=
-                  (fe_values.shape_grad(i, q_point) *
-                   fe_values.shape_grad(j, q_point) * fe_values.JxW(q_point));
+                  (fe_values.shape_grad(i, q_point) * // grad phi_i(x_q)
+                   fe_values.shape_grad(j, q_point) * // grad phi_j(x_q)
+                   fe_values.JxW(q_point));           // dx
 
-              cell_rhs(i) += (fe_values.shape_value(i, q_point) *
-                              rhs_values[q_point] * fe_values.JxW(q_point));
+              cell_rhs(i) += (fe_values.shape_value(i, q_point) * // phi_i(x_q)
+                              rhs_values[q_point] *               // f(x_q)
+                              fe_values.JxW(q_point));            // dx
             }
 
         local_dof_indices.resize(dofs_per_cell);
@@ -438,14 +437,9 @@ namespace Step27
     // all integers, so that it what we use:
     {
       Vector<float> fe_degrees(triangulation.n_active_cells());
-      {
-        typename hp::DoFHandler<dim>::active_cell_iterator
-          cell = dof_handler.begin_active(),
-          endc = dof_handler.end();
-        for (; cell != endc; ++cell)
-          fe_degrees(cell->active_cell_index()) =
-            fe_collection[cell->active_fe_index()].degree;
-      }
+      for (const auto &cell : dof_handler.active_cell_iterators())
+        fe_degrees(cell->active_cell_index()) =
+          fe_collection[cell->active_fe_index()].degree;
 
       // With now all data vectors available -- solution, estimated errors and
       // smoothness indicators, and finite element degrees --, we create a
@@ -499,21 +493,16 @@ namespace Step27
                                                smoothness_indicators.end()),
             min_smoothness = *std::max_element(smoothness_indicators.begin(),
                                                smoothness_indicators.end());
-      {
-        typename hp::DoFHandler<dim>::active_cell_iterator
-          cell = dof_handler.begin_active(),
-          endc = dof_handler.end();
-        for (; cell != endc; ++cell)
-          if (cell->refine_flag_set())
-            {
-              max_smoothness =
-                std::max(max_smoothness,
-                         smoothness_indicators(cell->active_cell_index()));
-              min_smoothness =
-                std::min(min_smoothness,
-                         smoothness_indicators(cell->active_cell_index()));
-            }
-      }
+      for (const auto &cell : dof_handler.active_cell_iterators())
+        if (cell->refine_flag_set())
+          {
+            max_smoothness =
+              std::max(max_smoothness,
+                       smoothness_indicators(cell->active_cell_index()));
+            min_smoothness =
+              std::min(min_smoothness,
+                       smoothness_indicators(cell->active_cell_index()));
+          }
       const float threshold_smoothness = (max_smoothness + min_smoothness) / 2;
 
       // With this, we can go back, loop over all cells again, and for those
@@ -524,20 +513,15 @@ namespace Step27
       // degree and in return remove the flag indicating that the cell should
       // undergo bisection. For all other cells, the refinement flags remain
       // untouched:
-      {
-        typename hp::DoFHandler<dim>::active_cell_iterator
-          cell = dof_handler.begin_active(),
-          endc = dof_handler.end();
-        for (; cell != endc; ++cell)
-          if (cell->refine_flag_set() &&
-              (smoothness_indicators(cell->active_cell_index()) >
-               threshold_smoothness) &&
-              (cell->active_fe_index() + 1 < fe_collection.size()))
-            {
-              cell->clear_refine_flag();
-              cell->set_active_fe_index(cell->active_fe_index() + 1);
-            }
-      }
+      for (const auto &cell : dof_handler.active_cell_iterators())
+        if (cell->refine_flag_set() &&
+            (smoothness_indicators(cell->active_cell_index()) >
+             threshold_smoothness) &&
+            (cell->active_fe_index() + 1 < fe_collection.size()))
+          {
+            cell->clear_refine_flag();
+            cell->set_active_fe_index(cell->active_fe_index() + 1);
+          }
 
       // At the end of this procedure, we then refine the mesh. During this
       // process, children of cells undergoing bisection inherit their mother
@@ -561,43 +545,28 @@ namespace Step27
   {
     const unsigned int dim = 2;
 
-    static const Point<2> vertices_1[] = {
-      Point<2>(-1., -1.),      Point<2>(-1. / 2, -1.),
-      Point<2>(0., -1.),       Point<2>(+1. / 2, -1.),
-      Point<2>(+1, -1.),
-
-      Point<2>(-1., -1. / 2.), Point<2>(-1. / 2, -1. / 2.),
-      Point<2>(0., -1. / 2.),  Point<2>(+1. / 2, -1. / 2.),
-      Point<2>(+1, -1. / 2.),
-
-      Point<2>(-1., 0.),       Point<2>(-1. / 2, 0.),
-      Point<2>(+1. / 2, 0.),   Point<2>(+1, 0.),
-
-      Point<2>(-1., 1. / 2.),  Point<2>(-1. / 2, 1. / 2.),
-      Point<2>(0., 1. / 2.),   Point<2>(+1. / 2, 1. / 2.),
-      Point<2>(+1, 1. / 2.),
-
-      Point<2>(-1., 1.),       Point<2>(-1. / 2, 1.),
-      Point<2>(0., 1.),        Point<2>(+1. / 2, 1.),
-      Point<2>(+1, 1.)};
-    const unsigned int n_vertices = sizeof(vertices_1) / sizeof(vertices_1[0]);
-    const std::vector<Point<dim>> vertices(&vertices_1[0],
-                                           &vertices_1[n_vertices]);
-    static const int cell_vertices[][GeometryInfo<dim>::vertices_per_cell] = {
-      {0, 1, 5, 6},
-      {1, 2, 6, 7},
-      {2, 3, 7, 8},
-      {3, 4, 8, 9},
-      {5, 6, 10, 11},
-      {8, 9, 12, 13},
-      {10, 11, 14, 15},
-      {12, 13, 17, 18},
-      {14, 15, 19, 20},
-      {15, 16, 20, 21},
-      {16, 17, 21, 22},
-      {17, 18, 22, 23}};
-    const unsigned int n_cells =
-      sizeof(cell_vertices) / sizeof(cell_vertices[0]);
+    const std::vector<Point<2>> vertices = {
+      {-1.0, -1.0}, {-0.5, -1.0}, {+0.0, -1.0}, {+0.5, -1.0}, {+1.0, -1.0}, //
+      {-1.0, -0.5}, {-0.5, -0.5}, {+0.0, -0.5}, {+0.5, -0.5}, {+1.0, -0.5}, //
+      {-1.0, +0.0}, {-0.5, +0.0}, {+0.5, +0.0}, {+1.0, +0.0},               //
+      {-1.0, +0.5}, {-0.5, +0.5}, {+0.0, +0.5}, {+0.5, +0.5}, {+1.0, +0.5}, //
+      {-1.0, +1.0}, {-0.5, +1.0}, {+0.0, +1.0}, {+0.5, +1.0}, {+1.0, +1.0}};
+
+    const std::vector<std::array<int, GeometryInfo<dim>::vertices_per_cell>>
+      cell_vertices = {{0, 1, 5, 6},
+                       {1, 2, 6, 7},
+                       {2, 3, 7, 8},
+                       {3, 4, 8, 9},
+                       {5, 6, 10, 11},
+                       {8, 9, 12, 13},
+                       {10, 11, 14, 15},
+                       {12, 13, 17, 18},
+                       {14, 15, 19, 20},
+                       {15, 16, 20, 21},
+                       {16, 17, 21, 22},
+                       {17, 18, 22, 23}};
+
+    const unsigned int n_cells = cell_vertices.size();
 
     std::vector<CellData<dim>> cells(n_cells, CellData<dim>());
     for (unsigned int i = 0; i < n_cells; ++i)
@@ -698,10 +667,7 @@ namespace Step27
     Vector<double> local_dof_values;
 
     // Then here is the loop:
-    typename hp::DoFHandler<dim>::active_cell_iterator cell = dof_handler
-                                                                .begin_active(),
-                                                       endc = dof_handler.end();
-    for (; cell != endc; ++cell)
+    for (const auto &cell : dof_handler.active_cell_iterators())
       {
         // Inside the loop, we first need to get the values of the local
         // degrees of freedom (which we put into the
@@ -752,8 +718,8 @@ namespace Step27
         // We have to calculate the logarithms of absolute
         // values of coefficients and use it in linear regression fit to
         // obtain $\mu$.
-        for (unsigned int f = 0; f < res.second.size(); f++)
-          res.second[f] = std::log(res.second[f]);
+        for (double &f : res.second)
+          f = std::log(f);
 
         std::pair<double, double> fit =
           FESeries::linear_regression(ln_k, res.second);

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.