Adjust test range and fix singleaccuracy tests

This changes the limit of exp() tests by one ulp
to account for rounding that leads to tests failing
for a single value when the CPU uses different
denormals-to-zero modes for SIMD and non-SIMD (ARMv7).

It also fixes the singleaccuracy tests to
check the correct functions for log2/exp2/pow2/cbrt
instead of the double-precision flavors, and the
implementation of powSingleAccuracy has been
copied from the (faster) single precision type
instead of just calling the double precision version.

Change-Id: I67d027cb916c4f78f36719f9ff554e9957879763

diff --git a/src/gromacs/simd/simd_math.h b/src/gromacs/simd/simd_math.h
index 727a654f484a045895f5beaa3bb7bf462062ed92..dab7a5adb54edf940429f29794e957bde8071deb 100644 (file)
--- a/src/gromacs/simd/simd_math.h
+++ b/src/gromacs/simd/simd_math.h
@@ -3972,7 +3972,27 @@ template <MathOptimization opt = MathOptimization::Safe>
 static inline SimdDouble gmx_simdcall
 powSingleAccuracy(SimdDouble x, SimdDouble y)
 {
-    return pow<opt>(x, y);
+    SimdDouble xcorr;
+
+    if (opt == MathOptimization::Safe)
+    {
+        xcorr = max(x, SimdDouble(std::numeric_limits<double>::min()));
+    }
+    else
+    {
+        xcorr = x;
+    }
+
+    SimdDouble result = exp2SingleAccuracy<opt>(y * log2SingleAccuracy(xcorr));
+
+    if (opt == MathOptimization::Safe)
+    {
+        // if x==0 and y>0 we explicitly set the result to 0.0
+        // For any x with y==0, the result will already be 1.0 since we multiply by y (0.0) and call exp().
+        result = blend(result, setZero(), x == setZero() && setZero() < y );
+    }
+
+    return result;
 }
 
 /*! \brief SIMD erf(x). Double precision SIMD data, single accuracy.

diff --git a/src/gromacs/simd/tests/simd_math.cpp b/src/gromacs/simd/tests/simd_math.cpp
index b8ed46f28b8a727fb98c727542dfb553c21e3319..ef2b2c1e9da13ad4e509cdb6fdcc39cf4907dd1e 100644 (file)
--- a/src/gromacs/simd/tests/simd_math.cpp
+++ b/src/gromacs/simd/tests/simd_math.cpp
@@ -747,7 +747,10 @@ TEST_F(SimdMathTest, exp)
     // Relevant threshold points. See the exp2 test for more details about the values; these are simply
     // scaled by log(2) due to the difference between exp2 and exp.
     const real      lowestReal                       = -std::numeric_limits<real>::max();
-    const real      lowestRealThatProducesNormal     = (std::numeric_limits<real>::min_exponent - 1)*std::log(2.0);
+    // In theory the smallest value should be (min_exponent-1)*log(2), but rounding after the multiplication will cause this
+    // value to be a single ulp too low. This might cause failed tests on CPUs that use different DTZ modes for SIMD vs.
+    // non-SIMD arithmetics (ARM v7), so multiply by (1.0-eps) to increase it by a single ulp.
+    const real      lowestRealThatProducesNormal     = (std::numeric_limits<real>::min_exponent - 1)*std::log(2.0)*(1-std::numeric_limits<real>::epsilon());
     const real      lowestRealThatProducesDenormal   = lowestRealThatProducesNormal - std::numeric_limits<real>::digits*std::log(2.0);
     const real      highestRealThatProducesNormal    = (std::numeric_limits<real>::max_exponent - 1)*std::log(2.0);
     CompareSettings settings;
@@ -1131,7 +1134,7 @@ TEST_F(SimdMathTest, cbrtSingleAccuracy)
     CompareSettings settings {
         Range(low, high), ulpTol_, absTol_, MatchRule::Normal
     };
-    GMX_EXPECT_SIMD_FUNC_NEAR(std::cbrt, cbrt, settings);
+    GMX_EXPECT_SIMD_FUNC_NEAR(std::cbrt, cbrtSingleAccuracy, settings);
 }
 
 TEST_F(SimdMathTest, invcbrtSingleAccuracy)
@@ -1146,13 +1149,13 @@ TEST_F(SimdMathTest, invcbrtSingleAccuracy)
     CompareSettings settings {
         Range(low, high), ulpTol_, absTol_, MatchRule::Normal
     };
-    GMX_EXPECT_SIMD_FUNC_NEAR(refInvCbrt, invcbrt, settings);
+    GMX_EXPECT_SIMD_FUNC_NEAR(refInvCbrt, invcbrtSingleAccuracy, settings);
 
     // Positive values
     low      = std::numeric_limits<real>::min();
     high     = std::numeric_limits<real>::max();
     settings = { Range(low, high), ulpTol_, absTol_, MatchRule::Normal };
-    GMX_EXPECT_SIMD_FUNC_NEAR(refInvCbrt, invcbrt, settings);
+    GMX_EXPECT_SIMD_FUNC_NEAR(refInvCbrt, invcbrtSingleAccuracy, settings);
 }
 
 TEST_F(SimdMathTest, log2SingleAccuracy)
@@ -1197,15 +1200,15 @@ TEST_F(SimdMathTest, exp2SingleAccuracy)
 
     // Below subnormal range all results should be zero
     settings = { Range(lowestReal, lowestRealThatProducesDenormal), ulpTol_, absTol_, MatchRule::Normal };
-    GMX_EXPECT_SIMD_FUNC_NEAR(std::exp2, exp2, settings);
+    GMX_EXPECT_SIMD_FUNC_NEAR(std::exp2, exp2SingleAccuracy, settings);
 
     // Subnormal range, require matching, but DTZ is fine
     settings = { Range(lowestRealThatProducesDenormal, lowestRealThatProducesNormal), ulpTol_, absTol_, MatchRule::Dtz };
-    GMX_EXPECT_SIMD_FUNC_NEAR(std::exp2, exp2, settings);
+    GMX_EXPECT_SIMD_FUNC_NEAR(std::exp2, exp2SingleAccuracy, settings);
 
     // Normal range, standard result expected
     settings = { Range(lowestRealThatProducesNormal, highestRealThatProducesNormal), ulpTol_, absTol_, MatchRule::Normal };
-    GMX_EXPECT_SIMD_FUNC_NEAR(std::exp2, exp2, settings);
+    GMX_EXPECT_SIMD_FUNC_NEAR(std::exp2, exp2SingleAccuracy, settings);
 }
 
 TEST_F(SimdMathTest, exp2SingleAccuracyUnsafe)
@@ -1227,22 +1230,28 @@ TEST_F(SimdMathTest, expSingleAccuracy)
 {
     // See threshold point comments in normal exp() test
     const real      lowestReal                       = -std::numeric_limits<real>::max();
-    const real      lowestRealThatProducesNormal     = (std::numeric_limits<real>::min_exponent - 1)*std::log(2.0);
+    // In theory the smallest value should be (min_exponent-1)*log(2), but rounding after the multiplication will cause this
+    // value to be a single ulp too low. This might cause failed tests on CPUs that use different DTZ modes for SIMD vs.
+    // non-SIMD arithmetics (ARM v7), so multiply by (1.0-eps) to increase it by a single ulp.
+    const real      lowestRealThatProducesNormal     =  (std::numeric_limits<real>::min_exponent - 1)*std::log(2.0)*(1.0-std::numeric_limits<real>::epsilon());
     const real      lowestRealThatProducesDenormal   = lowestRealThatProducesNormal - std::numeric_limits<real>::digits*std::log(2.0);
     const real      highestRealThatProducesNormal    = (std::numeric_limits<real>::max_exponent - 1)*std::log(2.0);
     CompareSettings settings;
 
+    // Increase the allowed error by the difference between the actual precision and single
+    setUlpTolSingleAccuracy(ulpTol_);
+
     // Below subnormal range all results should be zero (so, match the reference)
     settings = { Range(lowestReal, lowestRealThatProducesDenormal), ulpTol_, absTol_, MatchRule::Normal };
-    GMX_EXPECT_SIMD_FUNC_NEAR(std::exp, exp, settings);
+    GMX_EXPECT_SIMD_FUNC_NEAR(std::exp, expSingleAccuracy, settings);
 
     // Subnormal range, require matching, but DTZ is fine
     settings = { Range(lowestRealThatProducesDenormal, lowestRealThatProducesNormal), ulpTol_, absTol_, MatchRule::Dtz };
-    GMX_EXPECT_SIMD_FUNC_NEAR(std::exp, exp, settings);
+    GMX_EXPECT_SIMD_FUNC_NEAR(std::exp, expSingleAccuracy, settings);
 
     // Normal range, standard result expected
     settings = { Range(lowestRealThatProducesNormal, highestRealThatProducesNormal), ulpTol_, absTol_, MatchRule::Normal };
-    GMX_EXPECT_SIMD_FUNC_NEAR(std::exp, exp, settings);
+    GMX_EXPECT_SIMD_FUNC_NEAR(std::exp, expSingleAccuracy, settings);
 }
 
 TEST_F(SimdMathTest, expSingleAccuracyUnsafe)
@@ -1267,6 +1276,9 @@ TEST_F(SimdMathTest, powSingleAccuracy)
     // simple single-line function, so here we just test a handful of values to catch typos
     // and then some special values.
 
+    // Increase the allowed error by the difference between the actual precision and single
+    setUlpTolSingleAccuracy(ulpTol_);
+
     GMX_EXPECT_SIMD_REAL_NEAR(setSimdRealFrom3R(std::pow(c0, c3), std::pow(c1, c4), std::pow(c2, c5)),
                               powSingleAccuracy(rSimd_c0c1c2, rSimd_c3c4c5));
 
@@ -1274,7 +1286,7 @@ TEST_F(SimdMathTest, powSingleAccuracy)
                               powSingleAccuracy(rSimd_c0c1c2, rSimd_m3m0m4));
 
     // 0^0 = 1 , 0^c1=0, -c1^0=1
-    GMX_EXPECT_SIMD_REAL_NEAR(setSimdRealFrom3R(1.0, 0.0, 1.0), pow(setSimdRealFrom3R(0, 0.0, -c1), setSimdRealFrom3R(0.0, c1, 0.0)));
+    GMX_EXPECT_SIMD_REAL_NEAR(setSimdRealFrom3R(1.0, 0.0, 1.0), powSingleAccuracy(setSimdRealFrom3R(0, 0.0, -c1), setSimdRealFrom3R(0.0, c1, 0.0)));
 }
 
 TEST_F(SimdMathTest, powSingleAccuracyUnsafe)
@@ -1283,6 +1295,9 @@ TEST_F(SimdMathTest, powSingleAccuracyUnsafe)
     // simple single-line function, so here we just test a handful of values to catch typos
     // and then some special values.
 
+    // Increase the allowed error by the difference between the actual precision and single
+    setUlpTolSingleAccuracy(ulpTol_);
+
     GMX_EXPECT_SIMD_REAL_NEAR(setSimdRealFrom3R(std::pow(c0, c3), std::pow(c1, c4), std::pow(c2, c5)),
                               powSingleAccuracy<MathOptimization::Unsafe>(rSimd_c0c1c2, rSimd_c3c4c5));