diff options
author | Rasmus Munk Larsen <rmlarsen@google.com> | 2021-02-18 20:49:18 +0000 |
---|---|---|
committer | Rasmus Munk Larsen <rmlarsen@google.com> | 2021-02-18 20:49:18 +0000 |
commit | 7f09d3487de97882585a72819ec7e0f9100c7121 (patch) | |
tree | a6d9877365028e9870ff5f3d02705d7145060869 /Eigen/src | |
parent | 12fd3dd655e37ba26e7ab236d32163e0aa35da39 (diff) |
Use the Cephes double subtraction trick in pexp<float> even when FMA is available. Otherwise the accuracy drops from 1 ulp to 3 ulp.
Diffstat (limited to 'Eigen/src')
-rw-r--r-- | Eigen/src/Core/arch/Default/GenericPacketMathFunctions.h | 15 |
1 files changed, 4 insertions, 11 deletions
diff --git a/Eigen/src/Core/arch/Default/GenericPacketMathFunctions.h b/Eigen/src/Core/arch/Default/GenericPacketMathFunctions.h index abdbcdbb9..e4c7e6223 100644 --- a/Eigen/src/Core/arch/Default/GenericPacketMathFunctions.h +++ b/Eigen/src/Core/arch/Default/GenericPacketMathFunctions.h @@ -1,4 +1,3 @@ - // This file is part of Eigen, a lightweight C++ template library // for linear algebra. // @@ -468,16 +467,10 @@ Packet pexp_float(const Packet _x) // Get r = x - m*ln(2). If no FMA instructions are available, m*ln(2) is // subtracted out in two parts, m*C1+m*C2 = m*ln(2), to avoid accumulating // truncation errors. - Packet r; -#ifdef EIGEN_HAS_SINGLE_INSTRUCTION_MADD - const Packet cst_nln2 = pset1<Packet>(-0.6931471805599453f); - r = pmadd(m, cst_nln2, x); -#else - const Packet cst_cephes_exp_C1 = pset1<Packet>(0.693359375f); - const Packet cst_cephes_exp_C2 = pset1<Packet>(-2.12194440e-4f); - r = psub(x, pmul(m, cst_cephes_exp_C1)); - r = psub(r, pmul(m, cst_cephes_exp_C2)); -#endif + const Packet cst_cephes_exp_C1 = pset1<Packet>(-0.693359375f); + const Packet cst_cephes_exp_C2 = pset1<Packet>(2.12194440e-4f); + Packet r = pmadd(m, cst_cephes_exp_C1, x); + r = pmadd(m, cst_cephes_exp_C2, r); Packet r2 = pmul(r, r); Packet r3 = pmul(r2, r); |