]> https://gitweb.dealii.org/ - dealii.git/commitdiff
step-30: Remove references to an old version of step-12.
authorDavid Wells <drwells@email.unc.edu>
Wed, 8 May 2019 22:10:42 +0000 (18:10 -0400)
committerDavid Wells <drwells@email.unc.edu>
Thu, 9 May 2019 14:47:18 +0000 (10:47 -0400)
examples/step-30/doc/intro.dox
examples/step-30/doc/results.dox
examples/step-30/step-30.cc

index 81690a84407d011125783554bb13860949d53062..e38bcb05b1bd575c4b8448de1f7862d0ad4aa4cd 100644 (file)
@@ -483,8 +483,4 @@ extended to cover $[-1,1]\times[0,1]$ in 2D, where the flow field $\beta$ descri
 counterclockwise quarter circle around the origin in the right half of the
 domain and is parallel to the x-axis in the left part of the domain. The inflow
 boundary is again located at $x=1$ and along the positive part of the x-axis,
-and the boundary conditions are chosen as in step-12. Compared to step-12 we
-only use the more effective second assembling technique. In order to make
-comparisons more effective, we decided to keep function names like @p
-assemble_system2 even if there is only one of these routines in this tutorial
-program.
+and the boundary conditions are chosen as in step-12.
index 63d73aea48063ffcdf6c36bf33e5b9e836fea7f9..a71e3208e0eae33fa640642b7c7df7fad56f3d77 100644 (file)
@@ -9,42 +9,42 @@ Performing a 2D run with isotropic refinement...
 Cycle 0:
    Number of active cells:       128
    Number of degrees of freedom: 512
-Time of assemble_system2: 0.040003
+Time of assemble_system: 0.040003
 Writing grid to <grid-0.iso.eps>...
 Writing grid to <grid-0.iso.gnuplot>...
 Writing solution to <sol-0.iso.gnuplot>...
 Cycle 1:
    Number of active cells:       239
    Number of degrees of freedom: 956
-Time of assemble_system2: 0.072005
+Time of assemble_system: 0.072005
 Writing grid to <grid-1.iso.eps>...
 Writing grid to <grid-1.iso.gnuplot>...
 Writing solution to <sol-1.iso.gnuplot>...
 Cycle 2:
    Number of active cells:       491
    Number of degrees of freedom: 1964
-Time of assemble_system2: 0.144009
+Time of assemble_system: 0.144009
 Writing grid to <grid-2.iso.eps>...
 Writing grid to <grid-2.iso.gnuplot>...
 Writing solution to <sol-2.iso.gnuplot>...
 Cycle 3:
    Number of active cells:       1031
    Number of degrees of freedom: 4124
-Time of assemble_system2: 0.296019
+Time of assemble_system: 0.296019
 Writing grid to <grid-3.iso.eps>...
 Writing grid to <grid-3.iso.gnuplot>...
 Writing solution to <sol-3.iso.gnuplot>...
 Cycle 4:
    Number of active cells:       2027
    Number of degrees of freedom: 8108
-Time of assemble_system2: 0.576036
+Time of assemble_system: 0.576036
 Writing grid to <grid-4.iso.eps>...
 Writing grid to <grid-4.iso.gnuplot>...
 Writing solution to <sol-4.iso.gnuplot>...
 Cycle 5:
    Number of active cells:       4019
    Number of degrees of freedom: 16076
-Time of assemble_system2: 1.13607
+Time of assemble_system: 1.13607
 Writing grid to <grid-5.iso.eps>...
 Writing grid to <grid-5.iso.gnuplot>...
 Writing solution to <sol-5.iso.gnuplot>...
@@ -54,42 +54,42 @@ Performing a 2D run with anisotropic refinement...
 Cycle 0:
    Number of active cells:       128
    Number of degrees of freedom: 512
-Time of assemble_system2: 0.040003
+Time of assemble_system: 0.040003
 Writing grid to <grid-0.aniso.eps>...
 Writing grid to <grid-0.aniso.gnuplot>...
 Writing solution to <sol-0.aniso.gnuplot>...
 Cycle 1:
    Number of active cells:       171
    Number of degrees of freedom: 684
-Time of assemble_system2: 0.048003
+Time of assemble_system: 0.048003
 Writing grid to <grid-1.aniso.eps>...
 Writing grid to <grid-1.aniso.gnuplot>...
 Writing solution to <sol-1.aniso.gnuplot>...
 Cycle 2:
    Number of active cells:       255
    Number of degrees of freedom: 1020
-Time of assemble_system2: 0.072005
+Time of assemble_system: 0.072005
 Writing grid to <grid-2.aniso.eps>...
 Writing grid to <grid-2.aniso.gnuplot>...
 Writing solution to <sol-2.aniso.gnuplot>...
 Cycle 3:
    Number of active cells:       397
    Number of degrees of freedom: 1588
-Time of assemble_system2: 0.16401
+Time of assemble_system: 0.16401
 Writing grid to <grid-3.aniso.eps>...
 Writing grid to <grid-3.aniso.gnuplot>...
 Writing solution to <sol-3.aniso.gnuplot>...
 Cycle 4:
    Number of active cells:       658
    Number of degrees of freedom: 2632
-Time of assemble_system2: 0.192012
+Time of assemble_system: 0.192012
 Writing grid to <grid-4.aniso.eps>...
 Writing grid to <grid-4.aniso.gnuplot>...
 Writing solution to <sol-4.aniso.gnuplot>...
 Cycle 5:
    Number of active cells:       1056
    Number of degrees of freedom: 4224
-Time of assemble_system2: 0.304019
+Time of assemble_system: 0.304019
 Writing grid to <grid-5.aniso.eps>...
 Writing grid to <grid-5.aniso.gnuplot>...
 Writing solution to <sol-5.aniso.gnuplot>...
@@ -148,4 +148,3 @@ problems. However, that is not always the case. Considering boundary layers in
 compressible viscous flows, for example, the mesh is always aligned with the
 anisotropic features, thus anisotropic refinement will almost always increase the
 efficiency of computations on adapted grids for these cases.
-
index 0bdf2a3028566d45512706b937f8dba4b37b4181..fc4659da2f286f2b6c8b454116c74890385de550 100644 (file)
@@ -152,8 +152,7 @@ namespace Step30
   // @sect3{Class: DGTransportEquation}
   //
   // This declaration of this class is utterly unaffected by our current
-  // changes.  The only substantial change is that we use only the second
-  // assembly scheme described in step-12.
+  // changes.
   template <int dim>
   class DGTransportEquation
   {
@@ -168,12 +167,12 @@ namespace Step30
                                 FullMatrix<double> &     ui_vi_matrix,
                                 Vector<double> &         cell_vector) const;
 
-    void assemble_face_term2(const FEFaceValuesBase<dim> &fe_v,
-                             const FEFaceValuesBase<dim> &fe_v_neighbor,
-                             FullMatrix<double> &         ui_vi_matrix,
-                             FullMatrix<double> &         ue_vi_matrix,
-                             FullMatrix<double> &         ui_ve_matrix,
-                             FullMatrix<double> &         ue_ve_matrix) const;
+    void assemble_face_term(const FEFaceValuesBase<dim> &fe_v,
+                            const FEFaceValuesBase<dim> &fe_v_neighbor,
+                            FullMatrix<double> &         ui_vi_matrix,
+                            FullMatrix<double> &         ue_vi_matrix,
+                            FullMatrix<double> &         ui_ve_matrix,
+                            FullMatrix<double> &         ue_ve_matrix) const;
 
   private:
     const Beta<dim>           beta_function;
@@ -258,7 +257,7 @@ namespace Step30
 
 
   template <int dim>
-  void DGTransportEquation<dim>::assemble_face_term2(
+  void DGTransportEquation<dim>::assemble_face_term(
     const FEFaceValuesBase<dim> &fe_v,
     const FEFaceValuesBase<dim> &fe_v_neighbor,
     FullMatrix<double> &         ui_vi_matrix,
@@ -309,10 +308,8 @@ namespace Step30
 
   // @sect3{Class: DGMethod}
   //
-  // Even the main class of this program stays more or less the same. We omit
-  // one of the assembly routines and use only the second, more effective one
-  // of the two presented in step-12. However, we introduce a new routine
-  // (set_anisotropic_flags) and modify another one (refine_grid).
+  // This declaration is much like that of step-12. However, we introduce a
+  // new routine (set_anisotropic_flags) and modify another one (refine_grid).
   template <int dim>
   class DGMethod
   {
@@ -324,8 +321,7 @@ namespace Step30
 
   private:
     void setup_system();
-    void assemble_system1();
-    void assemble_system2();
+    void assemble_system();
     void solve(Vector<double> &solution);
     void refine_grid();
     void set_anisotropic_flags();
@@ -413,18 +409,19 @@ namespace Step30
   }
 
 
-  // @sect4{Function: assemble_system2}
+  // @sect4{Function: assemble_system}
   //
-  // We proceed with the <code>assemble_system2</code> function that
-  // implements the DG discretization in its second version. This function is
-  // very similar to the <code>assemble_system2</code> function from step-12,
-  // even the four cases considered for the neighbor-relations of a cell are
-  // the same, namely a) cell is at the boundary, b) there are finer
-  // neighboring cells, c) the neighbor is neither coarser nor finer and d)
-  // the neighbor is coarser.  However, the way in which we decide upon which
-  // case we have are modified in the way described in the introduction.
+  // We proceed with the <code>assemble_system</code> function that implements
+  // the DG discretization. This function does the same thing as the
+  // <code>assemble_system</code> function from step-12 (but without
+  // MeshWorker).  The four cases considered for the neighbor-relations of a
+  // cell are the same as the isotropic case, namely a) cell is at the
+  // boundary, b) there are finer neighboring cells, c) the neighbor is
+  // neither coarser nor finer and d) the neighbor is coarser.  However, the
+  // way in which we decide upon which case we have are modified in the way
+  // described in the introduction.
   template <int dim>
-  void DGMethod<dim>::assemble_system2()
+  void DGMethod<dim>::assemble_system()
   {
     const unsigned int dofs_per_cell = dof_handler.get_fe().dofs_per_cell;
     std::vector<types::global_dof_index> dofs(dofs_per_cell);
@@ -463,10 +460,7 @@ namespace Step30
 
     Vector<double> cell_vector(dofs_per_cell);
 
-    typename 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())
       {
         ui_vi_matrix = 0;
         cell_vector  = 0;
@@ -481,7 +475,7 @@ namespace Step30
              face_no < GeometryInfo<dim>::faces_per_cell;
              ++face_no)
           {
-            typename DoFHandler<dim>::face_iterator face = cell->face(face_no);
+            const auto face = cell->face(face_no);
 
             // Case (a): The face is at the boundary.
             if (face->at_boundary())
@@ -494,8 +488,7 @@ namespace Step30
               {
                 Assert(cell->neighbor(face_no).state() == IteratorState::valid,
                        ExcInternalError());
-                typename DoFHandler<dim>::cell_iterator neighbor =
-                  cell->neighbor(face_no);
+                const auto neighbor = cell->neighbor(face_no);
 
                 // Case (b): This is an internal face and the neighbor
                 // is refined (which we can test by asking whether the
@@ -536,7 +529,7 @@ namespace Step30
                         // use the @p neighbor_child_on_subface function. it
                         // takes care of all the complicated situations of
                         // anisotropic refinement and non-standard faces.
-                        typename DoFHandler<dim>::cell_iterator neighbor_child =
+                        const auto neighbor_child =
                           cell->neighbor_child_on_subface(face_no, subface_no);
                         Assert(!neighbor_child->has_children(),
                                ExcInternalError());
@@ -549,12 +542,12 @@ namespace Step30
                         fe_v_subface.reinit(cell, face_no, subface_no);
                         fe_v_face_neighbor.reinit(neighbor_child, neighbor2);
 
-                        dg.assemble_face_term2(fe_v_subface,
-                                               fe_v_face_neighbor,
-                                               ui_vi_matrix,
-                                               ue_vi_matrix,
-                                               ui_ve_matrix,
-                                               ue_ve_matrix);
+                        dg.assemble_face_term(fe_v_subface,
+                                              fe_v_face_neighbor,
+                                              ui_vi_matrix,
+                                              ue_vi_matrix,
+                                              ui_ve_matrix,
+                                              ue_ve_matrix);
 
                         neighbor_child->get_dof_indices(dofs_neighbor);
 
@@ -622,12 +615,12 @@ namespace Step30
                         fe_v_face.reinit(cell, face_no);
                         fe_v_face_neighbor.reinit(neighbor, neighbor2);
 
-                        dg.assemble_face_term2(fe_v_face,
-                                               fe_v_face_neighbor,
-                                               ui_vi_matrix,
-                                               ue_vi_matrix,
-                                               ui_ve_matrix,
-                                               ue_ve_matrix);
+                        dg.assemble_face_term(fe_v_face,
+                                              fe_v_face_neighbor,
+                                              ui_vi_matrix,
+                                              ue_vi_matrix,
+                                              ui_ve_matrix,
+                                              ue_ve_matrix);
 
                         neighbor->get_dof_indices(dofs_neighbor);
 
@@ -697,11 +690,8 @@ namespace Step30
                                                   gradient_indicator);
 
     // and scale it to obtain an error indicator.
-    typename DoFHandler<dim>::active_cell_iterator cell =
-                                                     dof_handler.begin_active(),
-                                                   endc = dof_handler.end();
-    for (unsigned int cell_no = 0; cell != endc; ++cell, ++cell_no)
-      gradient_indicator(cell_no) *=
+    for (const auto &cell : triangulation.active_cell_iterators())
+      gradient_indicator[cell->active_cell_index()] *=
         std::pow(cell->diameter(), 1 + 1.0 * dim / 2);
     // Then we use this indicator to flag the 30 percent of the cells with
     // highest error indicator to be refined.
@@ -749,11 +739,7 @@ namespace Step30
                                          update_values);
 
     // Now we need to loop over all active cells.
-    typename 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())
       // We only need to consider cells which are flagged for refinement.
       if (cell->refine_flag_set())
         {
@@ -764,16 +750,14 @@ namespace Step30
                face_no < GeometryInfo<dim>::faces_per_cell;
                ++face_no)
             {
-              typename DoFHandler<dim>::face_iterator face =
-                cell->face(face_no);
+              const auto face = cell->face(face_no);
 
               if (!face->at_boundary())
                 {
                   Assert(cell->neighbor(face_no).state() ==
                            IteratorState::valid,
                          ExcInternalError());
-                  typename DoFHandler<dim>::cell_iterator neighbor =
-                    cell->neighbor(face_no);
+                  const auto neighbor = cell->neighbor(face_no);
 
                   std::vector<double> u(fe_v_face.n_quadrature_points);
                   std::vector<double> u_neighbor(fe_v_face.n_quadrature_points);
@@ -795,10 +779,9 @@ namespace Step30
                         {
                           // get an iterator pointing to the cell behind the
                           // present subface...
-                          typename DoFHandler<dim>::cell_iterator
-                            neighbor_child =
-                              cell->neighbor_child_on_subface(face_no,
-                                                              subface_no);
+                          const auto neighbor_child =
+                            cell->neighbor_child_on_subface(face_no,
+                                                            subface_no);
                           Assert(!neighbor_child->has_children(),
                                  ExcInternalError());
                           // ... and reinit the respective FEFaceValues and
@@ -1032,8 +1015,8 @@ namespace Step30
                   << std::endl;
 
         Timer assemble_timer;
-        assemble_system2();
-        std::cout << "Time of assemble_system2: " << assemble_timer.cpu_time()
+        assemble_system();
+        std::cout << "Time of assemble_system: " << assemble_timer.cpu_time()
                   << std::endl;
         solve(solution2);
 

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.