]> https://gitweb.dealii.org/ - dealii.git/commitdiff
step-70: Remove quadrature_formula class attributes 12818/head
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Tue, 12 Oct 2021 09:02:01 +0000 (11:02 +0200)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Tue, 12 Oct 2021 09:02:01 +0000 (11:02 +0200)
examples/step-70/step-70.cc

index 79212bec0bed7a93246002a2b48a75adfe6b7b97..3e2d0360e9cec1f3e58ad7b4ac70afcc03c35172 100644 (file)
@@ -622,12 +622,12 @@ namespace Step70
     parallel::distributed::Triangulation<spacedim>      fluid_tria;
     parallel::distributed::Triangulation<dim, spacedim> solid_tria;
 
-    // Next come descriptions of the finite elements in use, along with
-    // appropriate quadrature formulas and the corresponding DoFHandler objects.
-    // For the current implementation, only `fluid_fe` is really necessary. For
-    // completeness, and to allow easy extension, we also keep the `solid_fe`
-    // around, which is however initialized to a FE_Nothing finite element
-    // space, i.e., one that has no degrees of freedom.
+    // Next come descriptions of the finite elements in use, along with the
+    // corresponding DoFHandler objects. For the current implementation, only
+    // `fluid_fe` is really necessary. For completeness, and to allow easy
+    // extension, we also keep the `solid_fe` around, which is however
+    // initialized to a FE_Nothing finite element space, i.e., one that has no
+    // degrees of freedom.
     //
     // We declare both finite element spaces as `std::unique_ptr` objects rather
     // than regular member variables, to allow their generation after
@@ -636,9 +636,6 @@ namespace Step70
     std::unique_ptr<FiniteElement<spacedim>>      fluid_fe;
     std::unique_ptr<FiniteElement<dim, spacedim>> solid_fe;
 
-    std::unique_ptr<Quadrature<spacedim>> fluid_quadrature_formula;
-    std::unique_ptr<Quadrature<dim>>      solid_quadrature_formula;
-
     DoFHandler<spacedim>      fluid_dh;
     DoFHandler<dim, spacedim> solid_dh;
 
@@ -1202,11 +1199,6 @@ namespace Step70
 
     solid_fe = std::make_unique<FE_Nothing<dim, spacedim>>();
     solid_dh.distribute_dofs(*solid_fe);
-
-    fluid_quadrature_formula =
-      std::make_unique<QGauss<spacedim>>(par.velocity_degree + 1);
-    solid_quadrature_formula =
-      std::make_unique<QGauss<dim>>(par.velocity_degree + 1);
   }
 
 
@@ -1335,14 +1327,15 @@ namespace Step70
 
     TimerOutput::Scope t(computing_timer, "Assemble Stokes terms");
 
+    QGauss<spacedim>   quadrature_formula(fluid_fe->degree + 1);
     FEValues<spacedim> fe_values(*fluid_fe,
-                                 *fluid_quadrature_formula,
+                                 quadrature_formula,
                                  update_values | update_gradients |
                                    update_quadrature_points |
                                    update_JxW_values);
 
     const unsigned int dofs_per_cell = fluid_fe->n_dofs_per_cell();
-    const unsigned int n_q_points    = fluid_quadrature_formula->size();
+    const unsigned int n_q_points    = quadrature_formula.size();
 
     FullMatrix<double> cell_matrix(dofs_per_cell, dofs_per_cell);
     FullMatrix<double> cell_matrix2(dofs_per_cell, dofs_per_cell);

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.