Convert step-61 to use face loops with cell->face_iterators().

This functionality was previously implemented in e2d7a2ff171a111480307d97dce2aa4db063a953
by @masterleinad but had crossed paths with 3f96459cacd673c231357e87c15e95315857e749
(in #8999).

diff --git a/examples/step-61/step-61.cc b/examples/step-61/step-61.cc
index ff77d4ecd134e1da0c156617faa5efac9fa5dfe4..97c0b84a6605dc40d020a64104a0f8b1f0f53785 100644 (file)
--- a/examples/step-61/step-61.cc
+++ b/examples/step-61/step-61.cc
@@ -442,8 +442,8 @@ namespace Step61
     // We need <code>FEValuesExtractors</code> to access the @p interior and
     // @p face component of the shape functions.
     const FEValuesExtractors::Vector velocities(0);
-    const FEValuesExtractors::Scalar interior(0);
-    const FEValuesExtractors::Scalar face(1);
+    const FEValuesExtractors::Scalar pressure_interior(0);
+    const FEValuesExtractors::Scalar pressure_face(1);
 
     // This finally gets us in position to loop over all cells. On
     // each cell, we will first calculate the various cell matrices
@@ -505,7 +505,8 @@ namespace Step61
                 fe_values_dgrt[velocities].divergence(i, q);
               for (unsigned int j = 0; j < dofs_per_cell; ++j)
                 {
-                  const double phi_j_interior = fe_values[interior].value(j, q);
+                  const double phi_j_interior =
+                    fe_values[pressure_interior].value(j, q);
 
                   cell_matrix_G(i, j) -=
                     (div_v_i * phi_j_interior * fe_values.JxW(q));
@@ -518,12 +519,10 @@ namespace Step61
         // of the polynomial space and the dot product of a basis function of
         // the Raviart-Thomas space and the normal vector. So we loop over all
         // the faces of the element and obtain the normal vector.
-        for (unsigned int face_n = 0;
-             face_n < GeometryInfo<dim>::faces_per_cell;
-             ++face_n)
+        for (const auto &face : cell->face_iterators())
           {
-            fe_face_values.reinit(cell, face_n);
-            fe_face_values_dgrt.reinit(cell_dgrt, face_n);
+            fe_face_values.reinit(cell, face);
+            fe_face_values_dgrt.reinit(cell_dgrt, face);
 
             for (unsigned int q = 0; q < n_face_q_points; ++q)
               {
@@ -536,7 +535,7 @@ namespace Step61
                     for (unsigned int j = 0; j < dofs_per_cell; ++j)
                       {
                         const double phi_j_face =
-                          fe_face_values[face].value(j, q);
+                          fe_face_values[pressure_face].value(j, q);
 
                         cell_matrix_G(i, j) +=
                           ((v_i * normal) * phi_j_face * fe_face_values.JxW(q));
@@ -582,7 +581,7 @@ namespace Step61
         for (unsigned int q = 0; q < n_q_points; ++q)
           for (unsigned int i = 0; i < dofs_per_cell; ++i)
             {
-              cell_rhs(i) += (fe_values[interior].value(i, q) *
+              cell_rhs(i) += (fe_values[pressure_interior].value(i, q) *
                               right_hand_side_values[q] * fe_values.JxW(q));
             }
 
@@ -694,9 +693,8 @@ namespace Step61
     std::vector<Tensor<2, dim>> coefficient_values(n_q_points_dgrt);
 
     const FEValuesExtractors::Vector velocities(0);
-    const FEValuesExtractors::Scalar pressure(dim);
-    const FEValuesExtractors::Scalar interior(0);
-    const FEValuesExtractors::Scalar face(1);
+    const FEValuesExtractors::Scalar pressure_interior(0);
+    const FEValuesExtractors::Scalar pressure_face(1);
 
     // In the introduction, we explained how to calculate the numerical velocity
     // on the cell. We need the pressure solution values on each cell,
@@ -760,19 +758,18 @@ namespace Step61
                 fe_values_dgrt[velocities].divergence(i, q);
               for (unsigned int j = 0; j < dofs_per_cell; ++j)
                 {
-                  const double phi_j_interior = fe_values[interior].value(j, q);
+                  const double phi_j_interior =
+                    fe_values[pressure_interior].value(j, q);
 
                   cell_matrix_G(i, j) -=
                     (div_v_i * phi_j_interior * fe_values.JxW(q));
                 }
             }
 
-        for (unsigned int face_n = 0;
-             face_n < GeometryInfo<dim>::faces_per_cell;
-             ++face_n)
+        for (const auto &face : cell->face_iterators())
           {
-            fe_face_values.reinit(cell, face_n);
-            fe_face_values_dgrt.reinit(cell_dgrt, face_n);
+            fe_face_values.reinit(cell, face);
+            fe_face_values_dgrt.reinit(cell_dgrt, face);
 
             for (unsigned int q = 0; q < n_face_q_points; ++q)
               {
@@ -785,7 +782,7 @@ namespace Step61
                     for (unsigned int j = 0; j < dofs_per_cell; ++j)
                       {
                         const double phi_j_face =
-                          fe_face_values[face].value(j, q);
+                          fe_face_values[pressure_face].value(j, q);
 
                         cell_matrix_G(i, j) +=
                           ((v_i * normal) * phi_j_face * fe_face_values.JxW(q));