156 lines · plain
1//===----------------------------------------------------------------------===//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-exception6//7//===----------------------------------------------------------------------===//8//9// Algorithm:10//11// e^x = 2^(x/ln(2)) = 2^(x*(64/ln(2))/64)12//13// x*(64/ln(2)) = n + f, |f| <= 0.5, n is integer14// n = 64*m + j, 0 <= j < 6415//16// e^x = 2^((64*m + j + f)/64)17// = (2^m) * (2^(j/64)) * 2^(f/64)18// = (2^m) * (2^(j/64)) * e^(f*(ln(2)/64))19//20// f = x*(64/ln(2)) - n21// r = f*(ln(2)/64) = x - n*(ln(2)/64)22//23// e^x = (2^m) * (2^(j/64)) * e^r24//25// (2^(j/64)) is precomputed26//27// e^r = 1 + r + (r^2)/2! + (r^3)/3! + (r^4)/4! + (r^5)/5!28// e^r = 1 + q29//30// q = r + (r^2)/2! + (r^3)/3! + (r^4)/4! + (r^5)/5!31//32// e^x = (2^m) * ( (2^(j/64)) + q*(2^(j/64)) )33//34//===----------------------------------------------------------------------===//35 36#if __CLC_FPSIZE == 3237 38_CLC_DEF _CLC_OVERLOAD __CLC_GENTYPE __clc_exp10(__CLC_GENTYPE x) {39 // 128*log2/log10 : 38.5318394449895940 const __CLC_GENTYPE X_MAX = 0x1.344134p+5f;41 // -149*log2/log10 : -44.853469353933242 const __CLC_GENTYPE X_MIN = -0x1.66d3e8p+5f;43 // 64*log10/log2 : 212.603398072791244 const __CLC_GENTYPE R_64_BY_LOG10_2 = 0x1.a934f0p+7f;45 // log2/(64 * log10) lead : 0.00469970746 const __CLC_GENTYPE R_LOG10_2_BY_64_LD = 0x1.340000p-8f;47 // log2/(64 * log10) tail : 0.0000038866505748 const __CLC_GENTYPE R_LOG10_2_BY_64_TL = 0x1.04d426p-18f;49 const __CLC_GENTYPE R_LN10 = 0x1.26bb1cp+1f;50 51 __CLC_INTN return_nan = __clc_isnan(x);52 __CLC_INTN return_inf = x > X_MAX;53 __CLC_INTN return_zero = x < X_MIN;54 55 __CLC_INTN n = __CLC_CONVERT_INTN(x * R_64_BY_LOG10_2);56 57 __CLC_GENTYPE fn = __CLC_CONVERT_GENTYPE(n);58 __CLC_INTN j = n & 0x3f;59 __CLC_INTN m = n >> 6;60 __CLC_INTN m2 = m << EXPSHIFTBITS_SP32;61 __CLC_GENTYPE r;62 63 r = R_LN10 *64 __clc_mad(fn, -R_LOG10_2_BY_64_TL, __clc_mad(fn, -R_LOG10_2_BY_64_LD, x));65 66 // Truncated Taylor series for e^r67 __CLC_GENTYPE z2 =68 __clc_mad(__clc_mad(__clc_mad(r, 0x1.555556p-5f, 0x1.555556p-3f), r,69 0x1.000000p-1f),70 r * r, r);71 72 __CLC_GENTYPE two_to_jby64 = __CLC_USE_TABLE(exp_tbl, j);73 z2 = __clc_mad(two_to_jby64, z2, two_to_jby64);74 75 __CLC_GENTYPE z2s = z2 * __CLC_AS_GENTYPE((__CLC_UINTN)0x1 << (m + 149));76 __CLC_GENTYPE z2n = __CLC_AS_GENTYPE(__CLC_AS_INTN(z2) + m2);77 z2 = m <= -126 ? z2s : z2n;78 79 z2 = return_inf ? __CLC_AS_GENTYPE((__CLC_UINTN)PINFBITPATT_SP32) : z2;80 z2 = return_zero ? 0.0f : z2;81 z2 = return_nan ? x : z2;82 return z2;83}84 85#elif __CLC_FPSIZE == 6486 87_CLC_DEF _CLC_OVERLOAD __CLC_GENTYPE __clc_exp10(__CLC_GENTYPE x) {88 // 1024*ln(2)/ln(10)89 const __CLC_GENTYPE X_MAX = 0x1.34413509f79ffp+8;90 // -1074*ln(2)/ln(10)91 const __CLC_GENTYPE X_MIN = -0x1.434e6420f4374p+8;92 // 64*ln(10)/ln(2)93 const __CLC_GENTYPE R_64_BY_LOG10_2 = 0x1.a934f0979a371p+7;94 // head ln(2)/(64*ln(10))95 const __CLC_GENTYPE R_LOG10_2_BY_64_LD = 0x1.3441350000000p-8;96 // tail ln(2)/(64*ln(10))97 const __CLC_GENTYPE R_LOG10_2_BY_64_TL = 0x1.3ef3fde623e25p-37;98 // ln(10)99 const __CLC_GENTYPE R_LN10 = 0x1.26bb1bbb55516p+1;100 101 __CLC_INTN n = __CLC_CONVERT_INTN(x * R_64_BY_LOG10_2);102 103 __CLC_GENTYPE dn = __CLC_CONVERT_GENTYPE(n);104 105 __CLC_INTN j = n & 0x3f;106 __CLC_INTN m = n >> 6;107 108 __CLC_GENTYPE r = R_LN10 * __clc_fma(-R_LOG10_2_BY_64_TL, dn,109 __clc_fma(-R_LOG10_2_BY_64_LD, dn, x));110 111 // 6 term tail of Taylor expansion of e^r112 __CLC_GENTYPE z2 =113 r * __clc_fma(114 r,115 __clc_fma(r,116 __clc_fma(r,117 __clc_fma(r,118 __clc_fma(r, 0x1.6c16c16c16c17p-10,119 0x1.1111111111111p-7),120 0x1.5555555555555p-5),121 0x1.5555555555555p-3),122 0x1.0000000000000p-1),123 1.0);124 125 __CLC_GENTYPE tv0 = __CLC_USE_TABLE(two_to_jby64_ep_tbl_head, j);126 __CLC_GENTYPE tv1 = __CLC_USE_TABLE(two_to_jby64_ep_tbl_tail, j);127 z2 = __clc_fma(tv0 + tv1, z2, tv1) + tv0;128 129 __CLC_INTN small_value =130 (m < -1022) || ((m == -1022) && __CLC_CONVERT_INTN(z2 < 1.0));131 132 __CLC_INTN n1 = m >> 2;133 __CLC_INTN n2 = m - n1;134 __CLC_GENTYPE z3 =135 z2 * __CLC_AS_GENTYPE((__CLC_CONVERT_LONGN(n1) + 1023) << 52);136 z3 *= __CLC_AS_GENTYPE((__CLC_CONVERT_LONGN(n2) + 1023) << 52);137 138 z2 = __clc_ldexp(z2, m);139 z2 = __CLC_CONVERT_LONGN(small_value) ? z3 : z2;140 141 z2 = __clc_isnan(x) ? x : z2;142 143 z2 = x > X_MAX ? __CLC_AS_GENTYPE((__CLC_ULONGN)PINFBITPATT_DP64) : z2;144 z2 = x < X_MIN ? 0.0 : z2;145 146 return z2;147}148 149#elif __CLC_FPSIZE == 16150 151_CLC_OVERLOAD _CLC_DEF __CLC_GENTYPE __clc_exp10(__CLC_GENTYPE x) {152 return __CLC_CONVERT_GENTYPE(__clc_exp10(__CLC_CONVERT_FLOATN(x)));153}154 155#endif156