fixed bug concerning escape grammage in FlatExponential

Environment/FlatExponential.h

-
 /*
- * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ * (c) Copyright 2019 CORSIKA Project, corsika-project@lists.kit.edu
  *
  * See file AUTHORS for a list of contributors.
  *
@@ -21,6 +20,7 @@
 #include <corsika/units/PhysicalUnits.h>
 
 #include <cassert>
+#include <limits>
 
 /**
  *
@@ -76,7 +76,14 @@ namespace corsika::environment {
       if (vDotA == 0) {
         return pGrammage / GetMassDensity(line.GetR0());
       } else {
-        return fLambda / vDotA * log(pGrammage * vDotA / (fRho0 * fLambda) + 1);
+        auto const logArg = pGrammage * vDotA / (fRho0 * fLambda) + 1;
+        if (logArg > 0) {
+          return fLambda / vDotA * log(pGrammage * vDotA / (fRho0 * fLambda) + 1);
+        } else {
+          return std::numeric_limits<typename decltype(
+                     pGrammage)::value_type>::infinity() *
+                 corsika::units::si::meter;
+        }
       }
     }
   };

Environment/testEnvironment.cc

@@ -75,7 +75,19 @@ TEST_CASE("FlatExponential") {
     REQUIRE((medium.ArclengthFromGrammage(trajectory, exact) / length) == Approx(1));
   }
 
-  SECTION("vertical") {
+  SECTION("escape grammage") {
+    Line const line(gOrigin, Vector<SpeedType::dimension_type>(
+                                 gCS, {0_m / second, 0_m / second, -5_m / second}));
+    Trajectory<Line> const trajectory(line, tEnd);
+
+    GrammageType const escapeGrammage = rho0 * lambda;
+
+    REQUIRE(trajectory.NormalizedDirection().dot(axis).magnitude() < 0);
+    REQUIRE(medium.ArclengthFromGrammage(trajectory, 1.2 * escapeGrammage) ==
+            std::numeric_limits<typename GrammageType::value_type>::infinity() * 1_m);
+  }
+
+  SECTION("inclined") {
     Line const line(gOrigin, Vector<SpeedType::dimension_type>(
                                  gCS, {0_m / second, 5_m / second, 5_m / second}));
     Trajectory<Line> const trajectory(line, tEnd);
@@ -83,7 +95,6 @@ TEST_CASE("FlatExponential") {
     LengthType const length = 2 * lambda;
     GrammageType const exact =
         rho0 * lambda * (exp(cosTheta * length / lambda) - 1) / cosTheta;
-
     REQUIRE((medium.IntegratedGrammage(trajectory, length) / exact) == Approx(1));
     REQUIRE((medium.ArclengthFromGrammage(trajectory, exact) / length) == Approx(1));
   }