1 /*
2 * Single-precision cosh(x) function.
3 *
4 * Copyright (c) 2022-2023, Arm Limited.
5 * SPDX-License-Identifier: MIT OR Apache-2.0 WITH LLVM-exception
6 */
7
8 #include "math_config.h"
9 #include "pl_sig.h"
10 #include "pl_test.h"
11
12 #define AbsMask 0x7fffffff
13 #define TinyBound 0x20000000 /* 0x1p-63: Round to 1 below this. */
14 #define SpecialBound \
15 0x42ad496c /* 0x1.5a92d8p+6: expf overflows above this, so have to use \
16 special case. */
17
18 float
19 optr_aor_exp_f32 (float);
20
21 static NOINLINE float
specialcase(float x,uint32_t iax)22 specialcase (float x, uint32_t iax)
23 {
24 if (iax == 0x7f800000)
25 return INFINITY;
26 if (iax > 0x7f800000)
27 return __math_invalidf (x);
28 if (iax <= TinyBound)
29 /* For tiny x, avoid underflow by just returning 1. */
30 return 1;
31 /* Otherwise SpecialBound <= |x| < Inf. x is too large to calculate exp(x)
32 without overflow, so use exp(|x|/2) instead. For large x cosh(x) is
33 dominated by exp(x), so return:
34 cosh(x) ~= (exp(|x|/2))^2 / 2. */
35 float t = optr_aor_exp_f32 (asfloat (iax) / 2);
36 return (0.5 * t) * t;
37 }
38
39 /* Approximation for single-precision cosh(x) using exp.
40 cosh(x) = (exp(x) + exp(-x)) / 2.
41 The maximum error is 1.89 ULP, observed for |x| > SpecialBound:
42 coshf(0x1.65898cp+6) got 0x1.f00aep+127 want 0x1.f00adcp+127.
43 The maximum error observed for TinyBound < |x| < SpecialBound is 1.02 ULP:
44 coshf(0x1.50a3cp+0) got 0x1.ff21dcp+0 want 0x1.ff21dap+0. */
45 float
coshf(float x)46 coshf (float x)
47 {
48 uint32_t ix = asuint (x);
49 uint32_t iax = ix & AbsMask;
50 float ax = asfloat (iax);
51
52 if (unlikely (iax <= TinyBound || iax >= SpecialBound))
53 {
54 /* x is tiny, large or special. */
55 return specialcase (x, iax);
56 }
57
58 /* Compute cosh using the definition:
59 coshf(x) = exp(x) / 2 + exp(-x) / 2. */
60 float t = optr_aor_exp_f32 (ax);
61 return 0.5f * t + 0.5f / t;
62 }
63
64 PL_SIG (S, F, 1, cosh, -10.0, 10.0)
65 PL_TEST_ULP (coshf, 1.89)
66 PL_TEST_INTERVAL (coshf, 0, 0x1p-63, 100)
67 PL_TEST_INTERVAL (coshf, 0, 0x1.5a92d8p+6, 80000)
68 PL_TEST_INTERVAL (coshf, 0x1.5a92d8p+6, inf, 2000)
69 PL_TEST_INTERVAL (coshf, -0, -0x1p-63, 100)
70 PL_TEST_INTERVAL (coshf, -0, -0x1.5a92d8p+6, 80000)
71 PL_TEST_INTERVAL (coshf, -0x1.5a92d8p+6, -inf, 2000)
72