1 //===-- Implementation of expm1f function ---------------------------------===//
2 //
3 // Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
4 // See https://llvm.org/LICENSE.txt for license information.
5 // SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
6 //
7 //===----------------------------------------------------------------------===//
8 
9 #include "src/math/expm1f.h"
10 #include "src/__support/FPUtil/BasicOperations.h"
11 #include "src/__support/FPUtil/PolyEval.h"
12 #include "src/__support/common.h"
13 #include "src/math/expf.h"
14 
15 namespace __llvm_libc {
16 
17 // When |x| > Ln2, catastrophic cancellation does not occur with the
18 // subtraction expf(x) - 1.0f, so we use it to compute expm1f(x).
19 //
20 // We divide [-Ln2; Ln2] into 3 subintervals [-Ln2; -1/8], [-1/8; 1/8],
21 // [1/8; Ln2]. And we use a degree-6 polynomial to approximate exp(x) - 1 in
22 // each interval. The coefficients were generated by Sollya's fpminmax.
23 //
24 // See libc/utils/mathtools/expm1f.sollya for more detail.
25 INLINE_FMA
26 LLVM_LIBC_FUNCTION(float, expm1f, (float x)) {
27   const float ln2 =
28       0.69314718055994530941723212145817656807550013436025f; // For C++17:
29                                                              // 0x1.62e'42ffp-1
30   float abs_x = __llvm_libc::fputil::abs(x);
31 
32   if (abs_x <= ln2) {
33     if (abs_x <= 0.125f) {
34       return x * __llvm_libc::fputil::polyeval(
35                      x, 1.0f, 0.5f, 0.16666664183139801025390625f,
36                      4.1666664183139801025390625e-2f,
37                      8.3379410207271575927734375e-3f,
38                      1.3894210569560527801513671875e-3f);
39     }
40     if (x > 0.125f) {
41       return __llvm_libc::fputil::polyeval(
42           x, 1.23142086749794543720781803131103515625e-7f,
43           0.9999969005584716796875f, 0.500031292438507080078125f,
44           0.16650259494781494140625f, 4.21491153538227081298828125e-2f,
45           7.53940828144550323486328125e-3f,
46           2.05591344274580478668212890625e-3f);
47     }
48     return __llvm_libc::fputil::polyeval(
49         x, -6.899231408397099585272371768951416015625e-8f,
50         0.999998271465301513671875f, 0.4999825656414031982421875f,
51         0.16657467186450958251953125f, 4.1390590369701385498046875e-2f,
52         7.856394164264202117919921875e-3f,
53         9.380675037391483783721923828125e-4f);
54   }
55   return expf(x) - 1.0f;
56 }
57 
58 } // namespace __llvm_libc
59