]> https://gitweb.dealii.org/ - dealii.git/commitdiff
step-34: C++ modernization.
authorDavid Wells <drwells@email.unc.edu>
Sat, 11 May 2019 17:19:20 +0000 (13:19 -0400)
committerDavid Wells <drwells@email.unc.edu>
Sun, 12 May 2019 02:58:22 +0000 (22:58 -0400)
1. Use dof_handler instead of dh to match the other tutorial programs.
2. Use std::to_string
3. Avoid using deprecated features (e.g., update_cell_normal_vectors)
4. Ranged-for loops
5. Remove a circular builds-on dependency with step-38

examples/step-34/doc/builds-on
examples/step-34/step-34.cc

index 239775109a5c9181270c2b38653b05e034e0fd51..48a0f738761717ac5d7b461e0b04951a515849bc 100644 (file)
@@ -1 +1 @@
-step-4 step-38
+step-4
index 4922584df712a07762432a819f2e8b0acf163eb3..2855d065b501a8266c16b0d93078612a6f5b2485 100644 (file)
@@ -238,7 +238,7 @@ namespace Step34
 
     Triangulation<dim - 1, dim> tria;
     FE_Q<dim - 1, dim>          fe;
-    DoFHandler<dim - 1, dim>    dh;
+    DoFHandler<dim - 1, dim>    dof_handler;
     MappingQ<dim - 1, dim>      mapping;
 
     // In BEM methods, the matrix that is generated is dense. Depending on the
@@ -317,7 +317,7 @@ namespace Step34
   BEMProblem<dim>::BEMProblem(const unsigned int fe_degree,
                               const unsigned int mapping_degree)
     : fe(fe_degree)
-    , dh(tria)
+    , dof_handler(tria)
     , mapping(mapping_degree, true)
     , wind(dim)
     , singular_quadrature_order(5)
@@ -436,15 +436,13 @@ namespace Step34
     }
     prm.leave_subsection();
 
-    prm.enter_subsection(std::string("Wind function ") +
-                         Utilities::int_to_string(dim) + std::string("d"));
+    prm.enter_subsection("Wind function " + std::to_string(dim) + "d");
     {
       wind.parse_parameters(prm);
     }
     prm.leave_subsection();
 
-    prm.enter_subsection(std::string("Exact solution ") +
-                         Utilities::int_to_string(dim) + std::string("d"));
+    prm.enter_subsection("Exact solution " + std::to_string(dim) + "d");
     {
       exact_solution.parse_parameters(prm);
     }
@@ -460,7 +458,7 @@ namespace Step34
     // the two simulations, we could do this by setting the corresponding "Run
     // 2d simulation" or "Run 3d simulation" flag to false:
     run_in_this_dimension =
-      prm.get_bool("Run " + Utilities::int_to_string(dim) + "d simulation");
+      prm.get_bool("Run " + std::to_string(dim) + "d simulation");
   }
 
 
@@ -535,9 +533,9 @@ namespace Step34
   {
     tria.refine_global(1);
 
-    dh.distribute_dofs(fe);
+    dof_handler.distribute_dofs(fe);
 
-    const unsigned int n_dofs = dh.n_dofs();
+    const unsigned int n_dofs = dof_handler.n_dofs();
 
     system_matrix.reinit(n_dofs, n_dofs);
 
@@ -561,7 +559,7 @@ namespace Step34
     FEValues<dim - 1, dim> fe_v(mapping,
                                 fe,
                                 *quadrature,
-                                update_values | update_cell_normal_vectors |
+                                update_values | update_normal_vectors |
                                   update_quadrature_points | update_JxW_values);
 
     const unsigned int n_q_points = fe_v.n_quadrature_points;
@@ -586,9 +584,9 @@ namespace Step34
 
     // We construct a vector of support points which will be used in the local
     // integrations:
-    std::vector<Point<dim>> support_points(dh.n_dofs());
+    std::vector<Point<dim>> support_points(dof_handler.n_dofs());
     DoFTools::map_dofs_to_support_points<dim - 1, dim>(mapping,
-                                                       dh,
+                                                       dof_handler,
                                                        support_points);
 
 
@@ -596,11 +594,7 @@ namespace Step34
     // we first initialize the FEValues object and get the values of
     // $\mathbf{\tilde v}$ at the quadrature points (this vector field should
     // be constant, but it doesn't hurt to be more general):
-    typename DoFHandler<dim - 1, dim>::active_cell_iterator cell =
-                                                              dh.begin_active(),
-                                                            endc = dh.end();
-
-    for (cell = dh.begin_active(); cell != endc; ++cell)
+    for (const auto &cell : dof_handler.active_cell_iterators())
       {
         fe_v.reinit(cell);
         cell->get_dof_indices(local_dof_indices);
@@ -616,7 +610,7 @@ namespace Step34
         // of the local degrees of freedom is the same as the support point
         // $i$. A the beginning of the loop we therefore check whether this is
         // the case, and we store which one is the singular index:
-        for (unsigned int i = 0; i < dh.n_dofs(); ++i)
+        for (unsigned int i = 0; i < dof_handler.n_dofs(); ++i)
           {
             local_matrix_row_i = 0;
 
@@ -681,8 +675,8 @@ namespace Step34
                   mapping,
                   fe,
                   singular_quadrature,
-                  update_jacobians | update_values |
-                    update_cell_normal_vectors | update_quadrature_points);
+                  update_jacobians | update_values | update_normal_vectors |
+                    update_quadrature_points);
 
                 fe_v_singular.reinit(cell);
 
@@ -737,12 +731,12 @@ namespace Step34
     // of the alpha angles, or solid angles (see the formula in the
     // introduction for this). The result is then added back onto the system
     // matrix object to yield the final form of the matrix:
-    Vector<double> ones(dh.n_dofs());
+    Vector<double> ones(dof_handler.n_dofs());
     ones.add(-1.);
 
     system_matrix.vmult(alpha, ones);
     alpha.add(1);
-    for (unsigned int i = 0; i < dh.n_dofs(); ++i)
+    for (unsigned int i = 0; i < dof_handler.n_dofs(); ++i)
       system_matrix(i, i) += alpha(i);
   }
 
@@ -768,7 +762,7 @@ namespace Step34
   {
     Vector<float> difference_per_cell(tria.n_active_cells());
     VectorTools::integrate_difference(mapping,
-                                      dh,
+                                      dof_handler,
                                       phi,
                                       exact_solution,
                                       difference_per_cell,
@@ -788,7 +782,7 @@ namespace Step34
 
     const double       alpha_error    = difference_per_node.linfty_norm();
     const unsigned int n_active_cells = tria.n_active_cells();
-    const unsigned int n_dofs         = dh.n_dofs();
+    const unsigned int n_dofs         = dof_handler.n_dofs();
 
     deallog << "Cycle " << cycle << ':' << std::endl
             << "   Number of active cells:       " << n_active_cells
@@ -934,15 +928,10 @@ namespace Step34
     external_dh.distribute_dofs(external_fe);
     external_phi.reinit(external_dh.n_dofs());
 
-    typename DoFHandler<dim - 1, dim>::active_cell_iterator cell =
-                                                              dh.begin_active(),
-                                                            endc = dh.end();
-
-
     FEValues<dim - 1, dim> fe_v(mapping,
                                 fe,
                                 *quadrature,
-                                update_values | update_cell_normal_vectors |
+                                update_values | update_normal_vectors |
                                   update_quadrature_points | update_JxW_values);
 
     const unsigned int n_q_points = fe_v.n_quadrature_points;
@@ -958,7 +947,7 @@ namespace Step34
                                               external_dh,
                                               external_support_points);
 
-    for (cell = dh.begin_active(); cell != endc; ++cell)
+    for (const auto &cell : dof_handler.active_cell_iterators())
       {
         fe_v.reinit(cell);
 
@@ -996,9 +985,8 @@ namespace Step34
     data_out.add_data_vector(external_phi, "external_phi");
     data_out.build_patches();
 
-    const std::string filename =
-      Utilities::int_to_string(dim) + "d_external.vtk";
-    std::ofstream file(filename);
+    const std::string filename = std::to_string(dim) + "d_external.vtk";
+    std::ofstream     file(filename);
 
     data_out.write_vtk(file);
   }
@@ -1013,7 +1001,7 @@ namespace Step34
   {
     DataOut<dim - 1, DoFHandler<dim - 1, dim>> dataout;
 
-    dataout.attach_dof_handler(dh);
+    dataout.attach_dof_handler(dof_handler);
     dataout.add_data_vector(
       phi, "phi", DataOut<dim - 1, DoFHandler<dim - 1, dim>>::type_dof_data);
     dataout.add_data_vector(
@@ -1025,9 +1013,8 @@ namespace Step34
       mapping.get_degree(),
       DataOut<dim - 1, DoFHandler<dim - 1, dim>>::curved_inner_cells);
 
-    std::string filename =
-      (Utilities::int_to_string(dim) + "d_boundary_solution_" +
-       Utilities::int_to_string(cycle) + ".vtk");
+    const std::string filename = std::to_string(dim) + "d_boundary_solution_" +
+                                 std::to_string(cycle) + ".vtk";
     std::ofstream file(filename);
 
     dataout.write_vtk(file);

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.