]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix step-33. 15437/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Thu, 22 Jun 2023 16:47:33 +0000 (10:47 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Thu, 22 Jun 2023 16:47:33 +0000 (10:47 -0600)
examples/step-33/doc/intro.dox
examples/step-33/step-33.cc

index 58ee54687b03fc5aeb61dc34200ecb25952147e7..4f5b09cceceb26de668529596824d915589702e0 100644 (file)
@@ -304,7 +304,7 @@ the continuous and DG methodologies. It also simplifies the generation
 of the Jacobian because we do not have to track constrained degrees of
 freedom through the automatic differentiation used to compute it.
 
-@note Whereas this program was written in 2008, we were unaware of any
+@note When this program was written in 2008, we were unaware of any
 publication that would actually have used this approach. However, a
 more recent paper by A. Dedner, R. Kl&ouml;fkorn, and M. Kr&auml;nkel
 ("Continuous Finite-Elements on Non-Conforming Grids Using
@@ -312,6 +312,26 @@ Discontinuous Galerkin Stabilization", Proceedings of Finite Volumes
 for Complex Applications VII - Methods and Theoretical Aspects,
 Springer, 2014) comes close.
 
+@note In hindsight, dealing with hanging nodes in this way is perhaps
+  not the best choice. deal.II contains many places that assume that
+  an element that has degrees of freedom on vertices, edges, and faces
+  (like the FE_Q element used here) represents finite element functions
+  on a triangulation that have certain continuity properties. Calling
+  DoFTools::make_hanging_node_constraints() and putting the corresponding
+  output into an AffineConstraints object allows one to enforce
+  these assumptions when one calls AffineConstraints::distribute()
+  at the end of the function that solves linear systems -- this call
+  can be found in all other tutorial programs that use meshes with
+  hanging nodes, and ensures that the resulting function's value at
+  hanging nodes is correct with regard to the values of the nodes on
+  the adjacent parent cell. On the other hand, this program does not do
+  this, and this runs afoul of an assertion that was later added to the
+  library and that tested that the function is indeed continuous.
+  To address the resulting failure due to the assertion not being
+  satisfied, we *do* after all enforce continuity at hanging nodes
+  just before we transfer the solution from one mesh to another,
+  in the `ConservationLaw::refine_grid()` function.
+
 Further, we enforce a maximum number of refinement levels to keep refinement under check.  It is the
 author's experience that for adaptivity for a time dependent problem, refinement can easily lead the simulation to
 a screeching halt, because of time step restrictions if the mesh
index 0e5b75d23f5f0e89a971f2b35dd4b373f9659084..6401e674eb861dc29d48a61f69bb69d90317ba2c 100644 (file)
@@ -2258,6 +2258,29 @@ namespace Step33
           cell->set_coarsen_flag();
       }
 
+    // The next step addresses a problem mentioned in a remark in the
+    // introduction: The SolutionTransfer class we want to use later on tests
+    // the assumption that the solution function is continuous at hanging
+    // nodes. This is not actually the case in this program because we chose
+    // (perhaps unwisely) to enforce hanging node constraints in a weak way,
+    // as one would for example do with discontinuous elements. But the elements
+    // we use here are continuous (namely, multiple copies of FE_Q), and so
+    // the assertion would fail and the program abort. To avoid the issue
+    // (without having to rewrite the whole program), we simply ensure that the
+    // solution *does* satisfy the hanging node constraints, but creating
+    // an AffineConstraint object that contains the hanging node constraints and
+    // applying the constraints to the two solution vectors we want the
+    // SolutionTransfer class to transfer to the next mesh:
+    {
+      AffineConstraints<double> hanging_node_constraints;
+      DoFTools::make_hanging_node_constraints(dof_handler,
+                                              hanging_node_constraints);
+      hanging_node_constraints.close();
+
+      hanging_node_constraints.distribute(old_solution);
+      hanging_node_constraints.distribute(predictor);
+    }
+
     // Then we need to transfer the various solution vectors from the old to
     // the new grid while we do the refinement. The SolutionTransfer class is
     // our friend here; it has a fairly extensive documentation, including

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.