1 /* 2 * Single-precision vector 2^x function. 3 * 4 * Copyright (c) 2019, Arm Limited. 5 * SPDX-License-Identifier: MIT 6 */ 7 8 #include "mathlib.h" 9 #include "v_math.h" 10 #if V_SUPPORTED 11 12 static const float Poly[] = { 13 /* maxerr: 0.878 ulp. */ 14 0x1.416b5ep-13f, 0x1.5f082ep-10f, 0x1.3b2dep-7f, 0x1.c6af7cp-5f, 0x1.ebfbdcp-3f, 0x1.62e43p-1f 15 }; 16 #define C0 v_f32 (Poly[0]) 17 #define C1 v_f32 (Poly[1]) 18 #define C2 v_f32 (Poly[2]) 19 #define C3 v_f32 (Poly[3]) 20 #define C4 v_f32 (Poly[4]) 21 #define C5 v_f32 (Poly[5]) 22 23 #define Shift v_f32 (0x1.8p23f) 24 #define InvLn2 v_f32 (0x1.715476p+0f) 25 #define Ln2hi v_f32 (0x1.62e4p-1f) 26 #define Ln2lo v_f32 (0x1.7f7d1cp-20f) 27 28 VPCS_ATTR 29 static v_f32_t 30 specialcase (v_f32_t poly, v_f32_t n, v_u32_t e, v_f32_t absn) 31 { 32 /* 2^n may overflow, break it up into s1*s2. */ 33 v_u32_t b = v_cond_u32 (n <= v_f32 (0.0f)) & v_u32 (0x83000000); 34 v_f32_t s1 = v_as_f32_u32 (v_u32 (0x7f000000) + b); 35 v_f32_t s2 = v_as_f32_u32 (e - b); 36 v_u32_t cmp = v_cond_u32 (absn > v_f32 (192.0f)); 37 v_f32_t r1 = s1 * s1; 38 v_f32_t r0 = poly * s1 * s2; 39 return v_as_f32_u32 ((cmp & v_as_u32_f32 (r1)) | (~cmp & v_as_u32_f32 (r0))); 40 } 41 42 VPCS_ATTR 43 v_f32_t 44 V_NAME(exp2f_1u) (v_f32_t x) 45 { 46 v_f32_t n, r, scale, poly, absn; 47 v_u32_t cmp, e; 48 49 /* exp2(x) = 2^n * poly(r), with poly(r) in [1/sqrt(2),sqrt(2)] 50 x = n + r, with r in [-1/2, 1/2]. */ 51 #if 0 52 v_f32_t z; 53 z = x + Shift; 54 n = z - Shift; 55 r = x - n; 56 e = v_as_u32_f32 (z) << 23; 57 #else 58 n = v_round_f32 (x); 59 r = x - n; 60 e = v_as_u32_s32 (v_round_s32 (x)) << 23; 61 #endif 62 scale = v_as_f32_u32 (e + v_u32 (0x3f800000)); 63 absn = v_abs_f32 (n); 64 cmp = v_cond_u32 (absn > v_f32 (126.0f)); 65 poly = v_fma_f32 (C0, r, C1); 66 poly = v_fma_f32 (poly, r, C2); 67 poly = v_fma_f32 (poly, r, C3); 68 poly = v_fma_f32 (poly, r, C4); 69 poly = v_fma_f32 (poly, r, C5); 70 poly = v_fma_f32 (poly, r, v_f32 (1.0f)); 71 if (unlikely (v_any_u32 (cmp))) 72 return specialcase (poly, n, e, absn); 73 return scale * poly; 74 } 75 #endif 76