diff options
author | Joe Ramsay <Joe.Ramsay@arm.com> | 2023-10-05 17:10:49 +0100 |
---|---|---|
committer | Szabolcs Nagy <szabolcs.nagy@arm.com> | 2023-10-23 15:00:45 +0100 |
commit | b39e9db5e305365db0c347fd308b7c0d86e3507b (patch) | |
tree | e7fccd5c92b1769601497ce2b83f97405dee54f4 /sysdeps/aarch64/fpu | |
parent | f554334c05a95c6b4df532ddc88cd3e72dc7d04c (diff) | |
download | glibc-b39e9db5e305365db0c347fd308b7c0d86e3507b.zip glibc-b39e9db5e305365db0c347fd308b7c0d86e3507b.tar.gz glibc-b39e9db5e305365db0c347fd308b7c0d86e3507b.tar.bz2 |
aarch64: Add vector implementations of exp2 routines
Some routines reuse table from v_exp_data.c
Diffstat (limited to 'sysdeps/aarch64/fpu')
-rw-r--r-- | sysdeps/aarch64/fpu/Makefile | 1 | ||||
-rw-r--r-- | sysdeps/aarch64/fpu/Versions | 4 | ||||
-rw-r--r-- | sysdeps/aarch64/fpu/bits/math-vector.h | 4 | ||||
-rw-r--r-- | sysdeps/aarch64/fpu/exp2_advsimd.c | 128 | ||||
-rw-r--r-- | sysdeps/aarch64/fpu/exp2_sve.c | 111 | ||||
-rw-r--r-- | sysdeps/aarch64/fpu/exp2f_advsimd.c | 124 | ||||
-rw-r--r-- | sysdeps/aarch64/fpu/exp2f_sve.c | 75 | ||||
-rw-r--r-- | sysdeps/aarch64/fpu/test-double-advsimd-wrappers.c | 1 | ||||
-rw-r--r-- | sysdeps/aarch64/fpu/test-double-sve-wrappers.c | 1 | ||||
-rw-r--r-- | sysdeps/aarch64/fpu/test-float-advsimd-wrappers.c | 1 | ||||
-rw-r--r-- | sysdeps/aarch64/fpu/test-float-sve-wrappers.c | 1 |
11 files changed, 451 insertions, 0 deletions
diff --git a/sysdeps/aarch64/fpu/Makefile b/sysdeps/aarch64/fpu/Makefile index a1bbc9b..9c7c768 100644 --- a/sysdeps/aarch64/fpu/Makefile +++ b/sysdeps/aarch64/fpu/Makefile @@ -1,5 +1,6 @@ libmvec-supported-funcs = cos \ exp \ + exp2 \ log \ sin \ tan diff --git a/sysdeps/aarch64/fpu/Versions b/sysdeps/aarch64/fpu/Versions index f0ca094..05de432 100644 --- a/sysdeps/aarch64/fpu/Versions +++ b/sysdeps/aarch64/fpu/Versions @@ -18,6 +18,10 @@ libmvec { _ZGVsMxv_sinf; } GLIBC_2.39 { + _ZGVnN4v_exp2f; + _ZGVnN2v_exp2; + _ZGVsMxv_exp2f; + _ZGVsMxv_exp2; _ZGVnN4v_tanf; _ZGVnN2v_tan; _ZGVsMxv_tanf; diff --git a/sysdeps/aarch64/fpu/bits/math-vector.h b/sysdeps/aarch64/fpu/bits/math-vector.h index 6193213..50921b2 100644 --- a/sysdeps/aarch64/fpu/bits/math-vector.h +++ b/sysdeps/aarch64/fpu/bits/math-vector.h @@ -51,12 +51,14 @@ typedef __SVBool_t __sv_bool_t; __vpcs __f32x4_t _ZGVnN4v_cosf (__f32x4_t); __vpcs __f32x4_t _ZGVnN4v_expf (__f32x4_t); +__vpcs __f32x4_t _ZGVnN4v_exp2f (__f32x4_t); __vpcs __f32x4_t _ZGVnN4v_logf (__f32x4_t); __vpcs __f32x4_t _ZGVnN4v_sinf (__f32x4_t); __vpcs __f32x4_t _ZGVnN4v_tanf (__f32x4_t); __vpcs __f64x2_t _ZGVnN2v_cos (__f64x2_t); __vpcs __f64x2_t _ZGVnN2v_exp (__f64x2_t); +__vpcs __f64x2_t _ZGVnN2v_exp2 (__f64x2_t); __vpcs __f64x2_t _ZGVnN2v_log (__f64x2_t); __vpcs __f64x2_t _ZGVnN2v_sin (__f64x2_t); __vpcs __f64x2_t _ZGVnN2v_tan (__f64x2_t); @@ -68,12 +70,14 @@ __vpcs __f64x2_t _ZGVnN2v_tan (__f64x2_t); __sv_f32_t _ZGVsMxv_cosf (__sv_f32_t, __sv_bool_t); __sv_f32_t _ZGVsMxv_expf (__sv_f32_t, __sv_bool_t); +__sv_f32_t _ZGVsMxv_exp2f (__sv_f32_t, __sv_bool_t); __sv_f32_t _ZGVsMxv_logf (__sv_f32_t, __sv_bool_t); __sv_f32_t _ZGVsMxv_sinf (__sv_f32_t, __sv_bool_t); __sv_f32_t _ZGVsMxv_tanf (__sv_f32_t, __sv_bool_t); __sv_f64_t _ZGVsMxv_cos (__sv_f64_t, __sv_bool_t); __sv_f64_t _ZGVsMxv_exp (__sv_f64_t, __sv_bool_t); +__sv_f64_t _ZGVsMxv_exp2 (__sv_f64_t, __sv_bool_t); __sv_f64_t _ZGVsMxv_log (__sv_f64_t, __sv_bool_t); __sv_f64_t _ZGVsMxv_sin (__sv_f64_t, __sv_bool_t); __sv_f64_t _ZGVsMxv_tan (__sv_f64_t, __sv_bool_t); diff --git a/sysdeps/aarch64/fpu/exp2_advsimd.c b/sysdeps/aarch64/fpu/exp2_advsimd.c new file mode 100644 index 0000000..8e26033 --- /dev/null +++ b/sysdeps/aarch64/fpu/exp2_advsimd.c @@ -0,0 +1,128 @@ +/* Double-precision vector (AdvSIMD) exp2 function + + Copyright (C) 2023 Free Software Foundation, Inc. + This file is part of the GNU C Library. + + The GNU C Library is free software; you can redistribute it and/or + modify it under the terms of the GNU Lesser General Public + License as published by the Free Software Foundation; either + version 2.1 of the License, or (at your option) any later version. + + The GNU C Library is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + Lesser General Public License for more details. + + You should have received a copy of the GNU Lesser General Public + License along with the GNU C Library; if not, see + <https://www.gnu.org/licenses/>. */ + +#include "v_math.h" +#include "poly_advsimd_f64.h" + +#define N (1 << V_EXP_TABLE_BITS) +#define IndexMask (N - 1) +#define BigBound 1022.0 +#define UOFlowBound 1280.0 + +static const struct data +{ + float64x2_t poly[4]; + float64x2_t shift, scale_big_bound, scale_uoflow_bound; +} data = { + /* Coefficients are computed using Remez algorithm with + minimisation of the absolute error. */ + .poly = { V2 (0x1.62e42fefa3686p-1), V2 (0x1.ebfbdff82c241p-3), + V2 (0x1.c6b09b16de99ap-5), V2 (0x1.3b2abf5571ad8p-7) }, + .shift = V2 (0x1.8p52 / N), + .scale_big_bound = V2 (BigBound), + .scale_uoflow_bound = V2 (UOFlowBound), +}; + +static inline uint64x2_t +lookup_sbits (uint64x2_t i) +{ + return (uint64x2_t){ __v_exp_data[i[0] & IndexMask], + __v_exp_data[i[1] & IndexMask] }; +} + +#if WANT_SIMD_EXCEPT + +# define TinyBound 0x2000000000000000 /* asuint64(0x1p-511). */ +# define Thres 0x2080000000000000 /* asuint64(512.0) - TinyBound. */ + +/* Call scalar exp2 as a fallback. */ +static float64x2_t VPCS_ATTR NOINLINE +special_case (float64x2_t x) +{ + return v_call_f64 (exp2, x, x, v_u64 (0xffffffffffffffff)); +} + +#else + +# define SpecialOffset 0x6000000000000000 /* 0x1p513. */ +/* SpecialBias1 + SpecialBias1 = asuint(1.0). */ +# define SpecialBias1 0x7000000000000000 /* 0x1p769. */ +# define SpecialBias2 0x3010000000000000 /* 0x1p-254. */ + +static float64x2_t VPCS_ATTR +special_case (float64x2_t s, float64x2_t y, float64x2_t n, + const struct data *d) +{ + /* 2^(n/N) may overflow, break it up into s1*s2. */ + uint64x2_t b = vandq_u64 (vclezq_f64 (n), v_u64 (SpecialOffset)); + float64x2_t s1 = vreinterpretq_f64_u64 (vsubq_u64 (v_u64 (SpecialBias1), b)); + float64x2_t s2 = vreinterpretq_f64_u64 ( + vaddq_u64 (vsubq_u64 (vreinterpretq_u64_f64 (s), v_u64 (SpecialBias2)), b)); + uint64x2_t cmp = vcagtq_f64 (n, d->scale_uoflow_bound); + float64x2_t r1 = vmulq_f64 (s1, s1); + float64x2_t r0 = vmulq_f64 (vfmaq_f64 (s2, s2, y), s1); + return vbslq_f64 (cmp, r1, r0); +} + +#endif + +/* Fast vector implementation of exp2. + Maximum measured error is 1.65 ulp. + _ZGVnN2v_exp2(-0x1.4c264ab5b559bp-6) got 0x1.f8db0d4df721fp-1 + want 0x1.f8db0d4df721dp-1. */ +VPCS_ATTR +float64x2_t V_NAME_D1 (exp2) (float64x2_t x) +{ + const struct data *d = ptr_barrier (&data); + uint64x2_t cmp; +#if WANT_SIMD_EXCEPT + uint64x2_t ia = vreinterpretq_u64_f64 (vabsq_f64 (x)); + cmp = vcgeq_u64 (vsubq_u64 (ia, v_u64 (TinyBound)), v_u64 (Thres)); + /* If any special case (inf, nan, small and large x) is detected, + fall back to scalar for all lanes. */ + if (__glibc_unlikely (v_any_u64 (cmp))) + return special_case (x); +#else + cmp = vcagtq_f64 (x, d->scale_big_bound); +#endif + + /* n = round(x/N). */ + float64x2_t z = vaddq_f64 (d->shift, x); + uint64x2_t u = vreinterpretq_u64_f64 (z); + float64x2_t n = vsubq_f64 (z, d->shift); + + /* r = x - n/N. */ + float64x2_t r = vsubq_f64 (x, n); + + /* s = 2^(n/N). */ + uint64x2_t e = vshlq_n_u64 (u, 52 - V_EXP_TABLE_BITS); + u = lookup_sbits (u); + float64x2_t s = vreinterpretq_f64_u64 (vaddq_u64 (u, e)); + + /* y ~ exp2(r) - 1. */ + float64x2_t r2 = vmulq_f64 (r, r); + float64x2_t y = v_pairwise_poly_3_f64 (r, r2, d->poly); + y = vmulq_f64 (r, y); + +#if !WANT_SIMD_EXCEPT + if (__glibc_unlikely (v_any_u64 (cmp))) + return special_case (s, y, n, d); +#endif + return vfmaq_f64 (s, s, y); +} diff --git a/sysdeps/aarch64/fpu/exp2_sve.c b/sysdeps/aarch64/fpu/exp2_sve.c new file mode 100644 index 0000000..c9314a9 --- /dev/null +++ b/sysdeps/aarch64/fpu/exp2_sve.c @@ -0,0 +1,111 @@ +/* Double-precision vector (SVE) exp2 function + + Copyright (C) 2023 Free Software Foundation, Inc. + This file is part of the GNU C Library. + + The GNU C Library is free software; you can redistribute it and/or + modify it under the terms of the GNU Lesser General Public + License as published by the Free Software Foundation; either + version 2.1 of the License, or (at your option) any later version. + + The GNU C Library is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + Lesser General Public License for more details. + + You should have received a copy of the GNU Lesser General Public + License along with the GNU C Library; if not, see + <https://www.gnu.org/licenses/>. */ + +#include "sv_math.h" +#include "poly_sve_f64.h" + +#define N (1 << V_EXP_TABLE_BITS) + +#define BigBound 1022 +#define UOFlowBound 1280 + +static const struct data +{ + double poly[4]; + double shift, big_bound, uoflow_bound; +} data = { + /* Coefficients are computed using Remez algorithm with + minimisation of the absolute error. */ + .poly = { 0x1.62e42fefa3686p-1, 0x1.ebfbdff82c241p-3, 0x1.c6b09b16de99ap-5, + 0x1.3b2abf5571ad8p-7 }, + .shift = 0x1.8p52 / N, + .uoflow_bound = UOFlowBound, + .big_bound = BigBound, +}; + +#define SpecialOffset 0x6000000000000000 /* 0x1p513. */ +/* SpecialBias1 + SpecialBias1 = asuint(1.0). */ +#define SpecialBias1 0x7000000000000000 /* 0x1p769. */ +#define SpecialBias2 0x3010000000000000 /* 0x1p-254. */ + +/* Update of both special and non-special cases, if any special case is + detected. */ +static inline svfloat64_t +special_case (svbool_t pg, svfloat64_t s, svfloat64_t y, svfloat64_t n, + const struct data *d) +{ + /* s=2^n may overflow, break it up into s=s1*s2, + such that exp = s + s*y can be computed as s1*(s2+s2*y) + and s1*s1 overflows only if n>0. */ + + /* If n<=0 then set b to 0x6, 0 otherwise. */ + svbool_t p_sign = svcmple (pg, n, 0.0); /* n <= 0. */ + svuint64_t b = svdup_u64_z (p_sign, SpecialOffset); + + /* Set s1 to generate overflow depending on sign of exponent n. */ + svfloat64_t s1 = svreinterpret_f64 (svsubr_x (pg, b, SpecialBias1)); + /* Offset s to avoid overflow in final result if n is below threshold. */ + svfloat64_t s2 = svreinterpret_f64 ( + svadd_x (pg, svsub_x (pg, svreinterpret_u64 (s), SpecialBias2), b)); + + /* |n| > 1280 => 2^(n) overflows. */ + svbool_t p_cmp = svacgt (pg, n, d->uoflow_bound); + + svfloat64_t r1 = svmul_x (pg, s1, s1); + svfloat64_t r2 = svmla_x (pg, s2, s2, y); + svfloat64_t r0 = svmul_x (pg, r2, s1); + + return svsel (p_cmp, r1, r0); +} + +/* Fast vector implementation of exp2. + Maximum measured error is 1.65 ulp. + _ZGVsMxv_exp2(-0x1.4c264ab5b559bp-6) got 0x1.f8db0d4df721fp-1 + want 0x1.f8db0d4df721dp-1. */ +svfloat64_t SV_NAME_D1 (exp2) (svfloat64_t x, svbool_t pg) +{ + const struct data *d = ptr_barrier (&data); + svbool_t no_big_scale = svacle (pg, x, d->big_bound); + svbool_t special = svnot_z (pg, no_big_scale); + + /* Reduce x to k/N + r, where k is integer and r in [-1/2N, 1/2N]. */ + svfloat64_t shift = sv_f64 (d->shift); + svfloat64_t kd = svadd_x (pg, x, shift); + svuint64_t ki = svreinterpret_u64 (kd); + /* kd = k/N. */ + kd = svsub_x (pg, kd, shift); + svfloat64_t r = svsub_x (pg, x, kd); + + /* scale ~= 2^(k/N). */ + svuint64_t idx = svand_x (pg, ki, N - 1); + svuint64_t sbits = svld1_gather_index (pg, __v_exp_data, idx); + /* This is only a valid scale when -1023*N < k < 1024*N. */ + svuint64_t top = svlsl_x (pg, ki, 52 - V_EXP_TABLE_BITS); + svfloat64_t scale = svreinterpret_f64 (svadd_x (pg, sbits, top)); + + /* Approximate exp2(r) using polynomial. */ + svfloat64_t r2 = svmul_x (pg, r, r); + svfloat64_t p = sv_pairwise_poly_3_f64_x (pg, r, r2, d->poly); + svfloat64_t y = svmul_x (pg, r, p); + + /* Assemble exp2(x) = exp2(r) * scale. */ + if (__glibc_unlikely (svptest_any (pg, special))) + return special_case (pg, scale, y, kd, d); + return svmla_x (pg, scale, scale, y); +} diff --git a/sysdeps/aarch64/fpu/exp2f_advsimd.c b/sysdeps/aarch64/fpu/exp2f_advsimd.c new file mode 100644 index 0000000..70b3ab6 --- /dev/null +++ b/sysdeps/aarch64/fpu/exp2f_advsimd.c @@ -0,0 +1,124 @@ +/* Single-precision vector (AdvSIMD) exp2 function + + Copyright (C) 2023 Free Software Foundation, Inc. + This file is part of the GNU C Library. + + The GNU C Library is free software; you can redistribute it and/or + modify it under the terms of the GNU Lesser General Public + License as published by the Free Software Foundation; either + version 2.1 of the License, or (at your option) any later version. + + The GNU C Library is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + Lesser General Public License for more details. + + You should have received a copy of the GNU Lesser General Public + License along with the GNU C Library; if not, see + <https://www.gnu.org/licenses/>. */ + +#include "v_math.h" + +static const struct data +{ + float32x4_t poly[5]; + uint32x4_t exponent_bias; +#if !WANT_SIMD_EXCEPT + float32x4_t special_bound, scale_thresh; +#endif +} data = { + /* maxerr: 1.962 ulp. */ + .poly = { V4 (0x1.59977ap-10f), V4 (0x1.3ce9e4p-7f), V4 (0x1.c6bd32p-5f), + V4 (0x1.ebf9bcp-3f), V4 (0x1.62e422p-1f) }, + .exponent_bias = V4 (0x3f800000), +#if !WANT_SIMD_EXCEPT + .special_bound = V4 (126.0f), + .scale_thresh = V4 (192.0f), +#endif +}; + +#define C(i) d->poly[i] + +#if WANT_SIMD_EXCEPT + +# define TinyBound v_u32 (0x20000000) /* asuint (0x1p-63). */ +# define BigBound v_u32 (0x42800000) /* asuint (0x1p6). */ +# define SpecialBound v_u32 (0x22800000) /* BigBound - TinyBound. */ + +static float32x4_t VPCS_ATTR NOINLINE +special_case (float32x4_t x, float32x4_t y, uint32x4_t cmp) +{ + /* If fenv exceptions are to be triggered correctly, fall back to the scalar + routine for special lanes. */ + return v_call_f32 (exp2f, x, y, cmp); +} + +#else + +# define SpecialOffset v_u32 (0x82000000) +# define SpecialBias v_u32 (0x7f000000) + +static float32x4_t VPCS_ATTR NOINLINE +special_case (float32x4_t poly, float32x4_t n, uint32x4_t e, uint32x4_t cmp1, + float32x4_t scale, const struct data *d) +{ + /* 2^n may overflow, break it up into s1*s2. */ + uint32x4_t b = vandq_u32 (vclezq_f32 (n), SpecialOffset); + float32x4_t s1 = vreinterpretq_f32_u32 (vaddq_u32 (b, SpecialBias)); + float32x4_t s2 = vreinterpretq_f32_u32 (vsubq_u32 (e, b)); + uint32x4_t cmp2 = vcagtq_f32 (n, d->scale_thresh); + float32x4_t r2 = vmulq_f32 (s1, s1); + float32x4_t r1 = vmulq_f32 (vfmaq_f32 (s2, poly, s2), s1); + /* Similar to r1 but avoids double rounding in the subnormal range. */ + float32x4_t r0 = vfmaq_f32 (scale, poly, scale); + float32x4_t r = vbslq_f32 (cmp1, r1, r0); + return vbslq_f32 (cmp2, r2, r); +} + +#endif + +float32x4_t VPCS_ATTR V_NAME_F1 (exp2) (float32x4_t x) +{ + const struct data *d = ptr_barrier (&data); + float32x4_t n, r, r2, scale, p, q, poly; + uint32x4_t cmp, e; + +#if WANT_SIMD_EXCEPT + /* asuint(|x|) - TinyBound >= BigBound - TinyBound. */ + uint32x4_t ia = vreinterpretq_u32_f32 (vabsq_f32 (x)); + cmp = vcgeq_u32 (vsubq_u32 (ia, TinyBound), SpecialBound); + float32x4_t xm = x; + /* If any lanes are special, mask them with 1 and retain a copy of x to allow + special_case to fix special lanes later. This is only necessary if fenv + exceptions are to be triggered correctly. */ + if (__glibc_unlikely (v_any_u32 (cmp))) + x = vbslq_f32 (cmp, v_f32 (1), x); +#endif + + /* exp2(x) = 2^n (1 + poly(r)), with 1 + poly(r) in [1/sqrt(2),sqrt(2)] + x = n + r, with r in [-1/2, 1/2]. */ + n = vrndaq_f32 (x); + r = vsubq_f32 (x, n); + e = vshlq_n_u32 (vreinterpretq_u32_s32 (vcvtaq_s32_f32 (x)), 23); + scale = vreinterpretq_f32_u32 (vaddq_u32 (e, d->exponent_bias)); + +#if !WANT_SIMD_EXCEPT + cmp = vcagtq_f32 (n, d->special_bound); +#endif + + r2 = vmulq_f32 (r, r); + p = vfmaq_f32 (C (1), C (0), r); + q = vfmaq_f32 (C (3), C (2), r); + q = vfmaq_f32 (q, p, r2); + p = vmulq_f32 (C (4), r); + poly = vfmaq_f32 (p, q, r2); + + if (__glibc_unlikely (v_any_u32 (cmp))) +#if WANT_SIMD_EXCEPT + return special_case (xm, vfmaq_f32 (scale, poly, scale), cmp); +#else + return special_case (poly, n, e, cmp, scale, d); +#endif + + return vfmaq_f32 (scale, poly, scale); +} diff --git a/sysdeps/aarch64/fpu/exp2f_sve.c b/sysdeps/aarch64/fpu/exp2f_sve.c new file mode 100644 index 0000000..88d10da --- /dev/null +++ b/sysdeps/aarch64/fpu/exp2f_sve.c @@ -0,0 +1,75 @@ +/* Single-precision vector (SVE) exp2 function + + Copyright (C) 2023 Free Software Foundation, Inc. + This file is part of the GNU C Library. + + The GNU C Library is free software; you can redistribute it and/or + modify it under the terms of the GNU Lesser General Public + License as published by the Free Software Foundation; either + version 2.1 of the License, or (at your option) any later version. + + The GNU C Library is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + Lesser General Public License for more details. + + You should have received a copy of the GNU Lesser General Public + License along with the GNU C Library; if not, see + <https://www.gnu.org/licenses/>. */ + +#include "sv_math.h" +#include "poly_sve_f32.h" + +static const struct data +{ + float poly[5]; + float shift, thres; +} data = { + /* Coefficients copied from the polynomial in AdvSIMD variant, reversed for + compatibility with polynomial helpers. */ + .poly = { 0x1.62e422p-1f, 0x1.ebf9bcp-3f, 0x1.c6bd32p-5f, 0x1.3ce9e4p-7f, + 0x1.59977ap-10f }, + /* 1.5*2^17 + 127. */ + .shift = 0x1.903f8p17f, + /* Roughly 87.3. For x < -Thres, the result is subnormal and not handled + correctly by FEXPA. */ + .thres = 0x1.5d5e2ap+6f, +}; + +static svfloat32_t NOINLINE +special_case (svfloat32_t x, svfloat32_t y, svbool_t special) +{ + return sv_call_f32 (exp2f, x, y, special); +} + +/* Single-precision SVE exp2f routine. Implements the same algorithm + as AdvSIMD exp2f. + Worst case error is 1.04 ULPs. + SV_NAME_F1 (exp2)(0x1.943b9p-1) got 0x1.ba7eb2p+0 + want 0x1.ba7ebp+0. */ +svfloat32_t SV_NAME_F1 (exp2) (svfloat32_t x, const svbool_t pg) +{ + const struct data *d = ptr_barrier (&data); + /* exp2(x) = 2^n (1 + poly(r)), with 1 + poly(r) in [1/sqrt(2),sqrt(2)] + x = n + r, with r in [-1/2, 1/2]. */ + svfloat32_t shift = sv_f32 (d->shift); + svfloat32_t z = svadd_x (pg, x, shift); + svfloat32_t n = svsub_x (pg, z, shift); + svfloat32_t r = svsub_x (pg, x, n); + + svbool_t special = svacgt (pg, x, d->thres); + svfloat32_t scale = svexpa (svreinterpret_u32 (z)); + + /* Polynomial evaluation: poly(r) ~ exp2(r)-1. + Evaluate polynomial use hybrid scheme - offset ESTRIN by 1 for + coefficients 1 to 4, and apply most significant coefficient directly. */ + svfloat32_t r2 = svmul_x (pg, r, r); + svfloat32_t p14 = sv_pairwise_poly_3_f32_x (pg, r, r2, d->poly + 1); + svfloat32_t p0 = svmul_x (pg, r, d->poly[0]); + svfloat32_t poly = svmla_x (pg, p0, r2, p14); + + if (__glibc_unlikely (svptest_any (pg, special))) + return special_case (x, svmla_x (pg, scale, scale, poly), special); + + return svmla_x (pg, scale, scale, poly); +} diff --git a/sysdeps/aarch64/fpu/test-double-advsimd-wrappers.c b/sysdeps/aarch64/fpu/test-double-advsimd-wrappers.c index f76726e..b2b36fd 100644 --- a/sysdeps/aarch64/fpu/test-double-advsimd-wrappers.c +++ b/sysdeps/aarch64/fpu/test-double-advsimd-wrappers.c @@ -25,6 +25,7 @@ VPCS_VECTOR_WRAPPER (cos_advsimd, _ZGVnN2v_cos) VPCS_VECTOR_WRAPPER (exp_advsimd, _ZGVnN2v_exp) +VPCS_VECTOR_WRAPPER (exp2_advsimd, _ZGVnN2v_exp2) VPCS_VECTOR_WRAPPER (log_advsimd, _ZGVnN2v_log) VPCS_VECTOR_WRAPPER (sin_advsimd, _ZGVnN2v_sin) VPCS_VECTOR_WRAPPER (tan_advsimd, _ZGVnN2v_tan) diff --git a/sysdeps/aarch64/fpu/test-double-sve-wrappers.c b/sysdeps/aarch64/fpu/test-double-sve-wrappers.c index 5a9e5b5..88b76ed 100644 --- a/sysdeps/aarch64/fpu/test-double-sve-wrappers.c +++ b/sysdeps/aarch64/fpu/test-double-sve-wrappers.c @@ -34,6 +34,7 @@ SVE_VECTOR_WRAPPER (cos_sve, _ZGVsMxv_cos) SVE_VECTOR_WRAPPER (exp_sve, _ZGVsMxv_exp) +SVE_VECTOR_WRAPPER (exp2_sve, _ZGVsMxv_exp2) SVE_VECTOR_WRAPPER (log_sve, _ZGVsMxv_log) SVE_VECTOR_WRAPPER (sin_sve, _ZGVsMxv_sin) SVE_VECTOR_WRAPPER (tan_sve, _ZGVsMxv_tan) diff --git a/sysdeps/aarch64/fpu/test-float-advsimd-wrappers.c b/sysdeps/aarch64/fpu/test-float-advsimd-wrappers.c index 66c6227..02ab609 100644 --- a/sysdeps/aarch64/fpu/test-float-advsimd-wrappers.c +++ b/sysdeps/aarch64/fpu/test-float-advsimd-wrappers.c @@ -25,6 +25,7 @@ VPCS_VECTOR_WRAPPER (cosf_advsimd, _ZGVnN4v_cosf) VPCS_VECTOR_WRAPPER (expf_advsimd, _ZGVnN4v_expf) +VPCS_VECTOR_WRAPPER (exp2f_advsimd, _ZGVnN4v_exp2f) VPCS_VECTOR_WRAPPER (logf_advsimd, _ZGVnN4v_logf) VPCS_VECTOR_WRAPPER (sinf_advsimd, _ZGVnN4v_sinf) VPCS_VECTOR_WRAPPER (tanf_advsimd, _ZGVnN4v_tanf) diff --git a/sysdeps/aarch64/fpu/test-float-sve-wrappers.c b/sysdeps/aarch64/fpu/test-float-sve-wrappers.c index 3cbf8b4..fa41ce0 100644 --- a/sysdeps/aarch64/fpu/test-float-sve-wrappers.c +++ b/sysdeps/aarch64/fpu/test-float-sve-wrappers.c @@ -34,6 +34,7 @@ SVE_VECTOR_WRAPPER (cosf_sve, _ZGVsMxv_cosf) SVE_VECTOR_WRAPPER (expf_sve, _ZGVsMxv_expf) +SVE_VECTOR_WRAPPER (exp2f_sve, _ZGVsMxv_exp2f) SVE_VECTOR_WRAPPER (logf_sve, _ZGVsMxv_logf) SVE_VECTOR_WRAPPER (sinf_sve, _ZGVsMxv_sinf) SVE_VECTOR_WRAPPER (tanf_sve, _ZGVsMxv_tanf) |