]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Make the tolerance in step-1 relative. 10923/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Wed, 16 Sep 2020 17:41:19 +0000 (11:41 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Wed, 16 Sep 2020 18:46:50 +0000 (12:46 -0600)
examples/step-1/step-1.cc
examples/step-2/step-2.cc

index 38eadd26f1edacce0faf0c45d6e124127ff55fc5..76ea18a660053cb9df429b7f2a219ac5af68f472 100644 (file)
@@ -210,14 +210,39 @@ void second_grid()
               // If this cell is at the inner boundary, then at least one of its
               // vertices must sit on the inner ring and therefore have a radial
               // distance from the center of exactly 0.5, up to floating point
-              // accuracy. Compute this distance, and if we have found a vertex
-              // with this property flag this cell for later refinement. We can
-              // then also break the loop over all vertices and move on to the
-              // next cell.
+              // accuracy. So we compute this distance, and if we find a vertex
+              // with this property, we flag this cell for later refinement. We
+              // can then also break the loop over all vertices and move on to
+              // the next cell.
+              //
+              // Because the distance from the center is computed as a floating
+              // point number, we have to expect that whatever we compute is
+              // only accurate to within
+              // [round-off](https://en.wikipedia.org/wiki/Round-off_error). As
+              // a consequence, we can never expect to compare the distance
+              // with the inner radius by equality: A statement such as
+              // `if (distance_from_center == inner_radius)` will fail
+              // unless we get exceptionally lucky. Rather, we need to do this
+              // comparison with a certain tolerance, and the usual way to do
+              // this is to write it as `if (std::abs(distance_from_center -
+              // inner_radius) <= tolerance)`
+              // where `tolerance` is some small number larger
+              // than round-off. The question is how to choose it: We could just
+              // pick, say, `1e-10`, but this is only appropriate if the objects
+              // we compare are of size one. If we had created a mesh with cells
+              // of size `1e+10`, then `1e-10` would be far lower than round-off
+              // and, as before, the comparison will only succeed if we get
+              // exceptionally lucky. Rather, it is almost always useful to make
+              // the tolerance *relative* to a typical "scale" of the objects
+              // being compared. Here, the "scale" would be the inner radius, or
+              // maybe the diameter of cells. We choose the former and set the
+              // tolerance equal to $10^{-6}$ times the inner radius of the
+              // annulus.
               const double distance_from_center =
                 center.distance(cell->vertex(v));
 
-              if (std::fabs(distance_from_center - inner_radius) < 1e-10)
+              if (std::fabs(distance_from_center - inner_radius) <=
+                  1e-6 * inner_radius)
                 {
                   cell->set_refine_flag();
                   break;
index e8a25a9d04c14b57d1ef0f62bb3cd113ae3fbc9a..74ea05a6f5cb1cbdb020c9de9ad7194818d597d7 100644 (file)
@@ -80,7 +80,8 @@ void make_grid(Triangulation<2> &triangulation)
             const double distance_from_center =
               center.distance(cell->vertex(v));
 
-            if (std::fabs(distance_from_center - inner_radius) < 1e-10)
+            if (std::fabs(distance_from_center - inner_radius) <=
+                1e-6 * inner_radius)
               {
                 cell->set_refine_flag();
                 break;

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.