]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix range-based for loops over the cells of triangulations. 15614/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Mon, 3 Jul 2023 14:52:56 +0000 (08:52 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Mon, 3 Jul 2023 16:30:47 +0000 (10:30 -0600)
source/grid/grid_generator.cc
tests/grid/grid_generator_airfoil_01.cc
tests/matrix_free/pbc_orientation_01.cc
tests/matrix_free/pbc_orientation_02.cc
tests/simplex/matrix_free_range_iteration_01.cc

index e0a2ed86994372e17341e886efae037c7aa78902..30a3fe72b1cdb667d79adf066d0da6b4dae8ddcf 100644 (file)
@@ -978,21 +978,21 @@ namespace GridGenerator
                                                 std::sin(gamma) * edge_length);
 
           // loop over vertices of all cells
-          for (auto &cell : tria)
+          for (auto &cell : tria.cell_iterators())
             for (const unsigned int v : GeometryInfo<2>::vertex_indices())
               {
                 // vertex has been already processed: nothing to do
-                if (vertex_processed[cell.vertex_index(v)])
+                if (vertex_processed[cell->vertex_index(v)])
                   continue;
 
                 // mark vertex as processed
-                vertex_processed[cell.vertex_index(v)] = true;
+                vertex_processed[cell->vertex_index(v)] = true;
 
-                auto &node = cell.vertex(v);
+                auto &node = cell->vertex(v);
 
                 // distinguish blocks
-                if (cell.material_id() == id_block_1 ||
-                    cell.material_id() == id_block_4) // block 1 and 4
+                if (cell->material_id() == id_block_1 ||
+                    cell->material_id() == id_block_4) // block 1 and 4
                   {
                     // step 1: rotate block 1 clockwise by gamma and move block
                     // 1 so that A(0) is on y-axis so that faces AD and BC are
@@ -1001,7 +1001,7 @@ namespace GridGenerator
                     // positive) Move trapeze to be in first quadrant by adding
                     // trapeze_offset
                     Point<2, double> node_;
-                    if (cell.material_id() == id_block_1)
+                    if (cell->material_id() == id_block_1)
                       {
                         node_ = Point<2, double>(rotation_matrix_1 *
                                                    (node - horizontal_offset) +
@@ -1010,7 +1010,7 @@ namespace GridGenerator
                     // step 1: rotate block 4 counterclockwise and move down so
                     // that trapeze is located in fourth quadrant (subtracting
                     // trapeze_offset)
-                    else if (cell.material_id() == id_block_4)
+                    else if (cell->material_id() == id_block_4)
                       {
                         node_ = Point<2, double>(rotation_matrix_2 *
                                                    (node - horizontal_offset) -
@@ -1054,19 +1054,19 @@ namespace GridGenerator
                         bias_alpha(1 - (1.0 * iy) / n_cells_y);
                       const double   theta = node_(0);
                       const Point<2> p(-height * std::cos(theta) + center_mesh,
-                                       ((cell.material_id() == id_block_1) ?
+                                       ((cell->material_id() == id_block_1) ?
                                           (height) :
                                           (-height)) *
                                          std::sin(theta));
-                      node =
-                        airfoil_1D[(
-                          (cell.material_id() == id_block_1) ? (0) : (1))][ix] *
-                          alpha +
-                        p * (1 - alpha);
+                      node = airfoil_1D[(
+                               (cell->material_id() == id_block_1) ? (0) : (1))]
+                                       [ix] *
+                               alpha +
+                             p * (1 - alpha);
                     }
                   }
-                else if (cell.material_id() == id_block_2 ||
-                         cell.material_id() == id_block_5) // block 2 and 5
+                else if (cell->material_id() == id_block_2 ||
+                         cell->material_id() == id_block_5) // block 2 and 5
                   {
                     // geometric parameters and indices for interpolation
                     Assert(
@@ -1094,19 +1094,19 @@ namespace GridGenerator
                     const Point<2> p(ix * dx + center_mesh +
                                        incline_factor * length_b2 * ix /
                                          n_cells_x_1,
-                                     ((cell.material_id() == id_block_2) ?
+                                     ((cell->material_id() == id_block_2) ?
                                         (height) :
                                         (-height)));
                     // interpolate between y = height and upper airfoil points
                     // (block2) or y = -height and lower airfoil points (block5)
                     node = airfoil_1D[(
-                             (cell.material_id() == id_block_2) ? (0) : (1))]
+                             (cell->material_id() == id_block_2) ? (0) : (1))]
                                      [n_cells_x_0 + ix] *
                              alpha +
                            p * (1 - alpha);
                   }
-                else if (cell.material_id() == id_block_3 ||
-                         cell.material_id() == id_block_6) // block 3 and 6
+                else if (cell->material_id() == id_block_3 ||
+                         cell->material_id() == id_block_6) // block 3 and 6
                   {
                     // compute indices ix and iy
                     const double dx = length_b2 / n_cells_x_2;
@@ -1124,7 +1124,7 @@ namespace GridGenerator
                     // points G and H to the right
                     const Point<2> p1(J(0) - (1 - incline_factor) * length_b2 *
                                                (alpha_x),
-                                      ((cell.material_id() == id_block_3) ?
+                                      ((cell->material_id() == id_block_3) ?
                                          (height) :
                                          (-height)));
                     // define points on HJ but use tail_y as y-coordinate, in
@@ -1135,7 +1135,7 @@ namespace GridGenerator
                 else
                   {
                     Assert(false,
-                           ExcIndexRange(cell.material_id(),
+                           ExcIndexRange(cell->material_id(),
                                          id_block_1,
                                          id_block_6));
                   }
@@ -8000,22 +8000,22 @@ namespace GridGenerator
     // (ii) create new midpoint vertex locations for each face (and record their
     // new indices in the 'face_to_new_vertex_indices' vector),
     // (iii) create new midpoint vertex locations for each cell (dim = 2 only)
-    for (const auto &cell : ref_tria)
+    for (const auto &cell : ref_tria.cell_iterators())
       {
         // temporary array storing the global indices of each cell entity in the
         // sequence: vertices, edges/faces, cell
         std::array<unsigned int, dim == 2 ? 9 : 14> local_vertex_indices;
 
         // (i) copy the existing vertex locations
-        for (const auto v : cell.vertex_indices())
+        for (const auto v : cell->vertex_indices())
           {
-            const auto v_global = cell.vertex_index(v);
+            const auto v_global = cell->vertex_index(v);
 
             if (old_to_new_vertex_indices[v_global] ==
                 numbers::invalid_unsigned_int)
               {
                 old_to_new_vertex_indices[v_global] = vertices.size();
-                vertices.push_back(cell.vertex(v));
+                vertices.push_back(cell->vertex(v));
               }
 
             AssertIndexRange(v, local_vertex_indices.size());
@@ -8023,32 +8023,32 @@ namespace GridGenerator
           }
 
         // (ii) create new midpoint vertex locations for each face
-        for (const auto f : cell.face_indices())
+        for (const auto f : cell->face_indices())
           {
-            const auto f_global = cell.face_index(f);
+            const auto f_global = cell->face_index(f);
 
             if (face_to_new_vertex_indices[f_global] ==
                 numbers::invalid_unsigned_int)
               {
                 face_to_new_vertex_indices[f_global] = vertices.size();
                 vertices.push_back(
-                  cell.face(f)->center(/*respect_manifold*/ true));
+                  cell->face(f)->center(/*respect_manifold*/ true));
               }
 
-            AssertIndexRange(cell.n_vertices() + f,
+            AssertIndexRange(cell->n_vertices() + f,
                              local_vertex_indices.size());
-            local_vertex_indices[cell.n_vertices() + f] =
+            local_vertex_indices[cell->n_vertices() + f] =
               face_to_new_vertex_indices[f_global];
           }
 
         // (iii) create new midpoint vertex locations for each cell
         if (dim == 2)
           {
-            AssertIndexRange(cell.n_vertices() + cell.n_faces(),
+            AssertIndexRange(cell->n_vertices() + cell->n_faces(),
                              local_vertex_indices.size());
-            local_vertex_indices[cell.n_vertices() + cell.n_faces()] =
+            local_vertex_indices[cell->n_vertices() + cell->n_faces()] =
               vertices.size();
-            vertices.push_back(cell.center(/*respect_manifold*/ true));
+            vertices.push_back(cell->center(/*respect_manifold*/ true));
           }
 
         // helper function for creating cells and subcells
@@ -8139,13 +8139,13 @@ namespace GridGenerator
             }
         };
 
-        const auto material_id_cell = cell.material_id();
+        const auto material_id_cell = cell->material_id();
 
         // create cells one by one
         if (dim == 2)
           {
             // get cell-manifold id from current quad cell
-            const auto manifold_id_cell = cell.manifold_id();
+            const auto manifold_id_cell = cell->manifold_id();
             // inherit cell manifold
             for (const auto &cell_vertices : table_2D_cell)
               add_cell(dim, cell_vertices, material_id_cell, manifold_id_cell);
@@ -8163,7 +8163,7 @@ namespace GridGenerator
         else if (dim == 3)
           {
             // get cell-manifold id from current quad cell
-            const auto manifold_id_cell = cell.manifold_id();
+            const auto manifold_id_cell = cell->manifold_id();
             // inherit cell manifold
             for (const auto &cell_vertices : vertex_ids_for_cells_3d)
               add_cell(dim, cell_vertices, material_id_cell, manifold_id_cell);
@@ -8188,10 +8188,10 @@ namespace GridGenerator
           Assert(false, ExcNotImplemented());
 
         // Set up sub-cell data.
-        for (const auto f : cell.face_indices())
+        for (const auto f : cell->face_indices())
           {
-            const auto bid = cell.face(f)->boundary_id();
-            const auto mid = cell.face(f)->manifold_id();
+            const auto bid = cell->face(f)->boundary_id();
+            const auto mid = cell->face(f)->manifold_id();
 
             // process boundary-faces: set boundary and manifold ids
             if (dim == 2) // 2d boundary-faces
@@ -8221,9 +8221,9 @@ namespace GridGenerator
         // triangulation.
         if (dim == 3)
           {
-            for (const auto e : cell.line_indices())
+            for (const auto e : cell->line_indices())
               {
-                auto edge = cell.line(e);
+                auto edge = cell->line(e);
                 // Rather than use add_cell(), which does additional index
                 // translation, just add edges directly into subcell_data since
                 // we already know the correct global vertex indices.
index 9aed786a01e320b727183ba9669674328b3894b3..592210cd3bdde852d89f8bf9256c597ed280983d 100644 (file)
@@ -32,17 +32,18 @@ template <int dim>
 void
 print_triangulation(const Triangulation<dim, dim> &tria)
 {
-  for (const auto &cell : tria)
+  for (const auto &cell : tria.cell_iterators())
     {
-      deallog << cell.material_id() << ' ';
+      deallog << cell->material_id() << ' ';
 
       for (const unsigned int f : GeometryInfo<dim>::face_indices())
-        deallog << (cell.face(f)->at_boundary() ? cell.face(f)->boundary_id() :
-                                                  -1)
+        deallog << (cell->face(f)->at_boundary() ?
+                      cell->face(f)->boundary_id() :
+                      -1)
                 << ' ';
 
       for (const unsigned int v : GeometryInfo<dim>::vertex_indices())
-        deallog << cell.vertex(v) << ' ';
+        deallog << cell->vertex(v) << ' ';
 
       deallog << std::endl;
     }
index 22c895bda28092516d46b28a67b51c28af1fcc73..c8c981ce27cbd1e819fcf97275d009662f2dbebe 100644 (file)
@@ -120,8 +120,8 @@ test()
       DoFHandler<3> dof_handler(tria);
       dof_handler.distribute_dofs(fe);
 
-      for (auto &cell : tria)
-        for (auto &face : cell.face_iterators())
+      for (auto &cell : tria.cell_iterators())
+        for (auto &face : cell->face_iterators())
           {
             if (std::abs(face->center()[2] - 0.0) < 10e-6)
               face->set_boundary_id(4);
index 1b5f150c19b4fce2ec3203ec0705af2d64471235..28759c519e15bcd6f2cf3ca6f2e59d30385efad8 100644 (file)
@@ -120,8 +120,8 @@ test()
       DoFHandler<3> dof_handler(tria);
       dof_handler.distribute_dofs(fe);
 
-      for (auto &cell : tria)
-        for (auto &face : cell.face_iterators())
+      for (auto &cell : tria.cell_iterators())
+        for (auto &face : cell->face_iterators())
           {
             if (std::abs(face->center()[2] - 0.0) < 10e-6)
               face->set_boundary_id(4);
index 7be3bc738cea91dd200d2f68a3f0377486adb378..be3a78d8be1b635857de8aa5ac3ec54b67c781f6 100644 (file)
@@ -54,9 +54,9 @@ test()
 
   unsigned int counter = 0;
 
-  for (auto &cell : dof_handler)
+  for (auto &cell : dof_handler.cell_iterators())
     if (counter++ < tria.n_cells() / 2)
-      cell.set_active_fe_index(1);
+      cell->set_active_fe_index(1);
 
   dof_handler.distribute_dofs(fe);
 

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.