]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Avoid the use of infinities.
authorWolfgang Bangerth <bangerth@colostate.edu>
Wed, 18 Dec 2024 04:18:32 +0000 (21:18 -0700)
committerWolfgang Bangerth <bangerth@colostate.edu>
Wed, 18 Dec 2024 04:37:17 +0000 (21:37 -0700)
include/deal.II/fe/mapping.h
include/deal.II/fe/mapping_q_internal.h
source/fe/mapping.cc
source/fe/mapping_q.cc

index 3978d43e61d878eb08706e48334c9ec24d5e0e7b..e578c853547e264f1c2c95e7015c51e48fc29466 100644 (file)
@@ -498,7 +498,8 @@ public:
    * MappingQ. The only difference in behavior is that this function
    * will never throw an ExcTransformationFailed() exception. If the
    * transformation fails for `real_points[i]`, the returned `unit_points[i]`
-   * contains std::numeric_limits<double>::infinity() as the first entry.
+   * contains std::numeric_limits<double>::lowest() as the first component
+   * of the point, marking this one point as invalid.
    */
   virtual void
   transform_points_real_to_unit_cell(
index e4299356bcfe8d17f9b7448312aa5b6bb4b582be..210f90565f77de5abeda875441f2764f61be1508 100644 (file)
@@ -552,7 +552,7 @@ namespace internal
       const unsigned int newton_iteration_limit = 20;
 
       Point<dim, Number> invalid_point;
-      invalid_point[0]                = std::numeric_limits<double>::infinity();
+      invalid_point[0]                = std::numeric_limits<double>::lowest();
       bool tried_project_to_unit_cell = false;
 
       unsigned int newton_iteration            = 0;
index 314a2d18a887774841d7b89ca4a82e8d32364a24..c7fb0891e249cc456e42baab7cda7dd23491098f 100644 (file)
@@ -142,8 +142,10 @@ Mapping<dim, spacedim>::transform_points_real_to_unit_cell(
         }
       catch (typename Mapping<dim>::ExcTransformationFailed &)
         {
+          // If the transformation for this one point failed, mark it
+          // as invalid as described in the documentation.
           unit_points[i]    = Point<dim>();
-          unit_points[i][0] = std::numeric_limits<double>::infinity();
+          unit_points[i][0] = std::numeric_limits<double>::lowest();
         }
     }
 }
index 60bdd76ccea7bdddd8609abab3978bc420980404..734fb98646513b288fccabc479efb0c5e4f54d6d 100644 (file)
@@ -615,7 +615,7 @@ MappingQ<dim, spacedim>::transform_real_to_unit_cell(
   // statement may throw an exception, which we simply pass up to the caller
   const Point<dim> p_unit =
     this->transform_real_to_unit_cell_internal(cell, p, initial_p_unit);
-  AssertThrow(numbers::is_finite(p_unit[0]),
+  AssertThrow(p_unit[0] != std::numeric_limits<double>::lowest(),
               (typename Mapping<dim, spacedim>::ExcTransformationFailed()));
   return p_unit;
 }
@@ -690,7 +690,7 @@ MappingQ<dim, spacedim>::transform_points_real_to_unit_cell(
         // determinants) from other SIMD lanes. Repeat the computation in this
         // unlikely case with scalar arguments.
         for (unsigned int j = 0; j < n_lanes && i + j < n_points; ++j)
-          if (numbers::is_finite(unit_point[0][j]))
+          if (unit_point[0][j] != std::numeric_limits<double>::lowest())
             for (unsigned int d = 0; d < dim; ++d)
               unit_points[i + j][d] = unit_point[d][j];
           else

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.