200 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#if __CLC_FPSIZE == 3210 11_CLC_OVERLOAD _CLC_DEF __CLC_GENTYPE __clc_cosh(__CLC_GENTYPE x) {12 // After dealing with special cases the computation is split into regions as13 // follows. abs(x) >= max_cosh_arg: cosh(x) = sign(x)*Inf abs(x) >=14 // small_threshold: cosh(x) = sign(x)*exp(abs(x))/2 computed using the15 // splitexp and scaleDouble functions as for exp_amd().16 // abs(x) < small_threshold:17 // compute p = exp(y) - 1 and then z = 0.5*(p+(p/(p+1.0)))18 // cosh(x) is then z.19 20 const __CLC_GENTYPE max_cosh_arg = 0x1.65a9fap+6f;21 const __CLC_GENTYPE small_threshold = 0x1.0a2b24p+3f;22 23 __CLC_UINTN ux = __CLC_AS_UINTN(x);24 __CLC_GENTYPE y = __clc_fabs(x);25 __CLC_UINTN aux = __CLC_AS_UINTN(y);26 27 // Find the integer part y0 of y and the increment dy = y - y0. We then28 // compute z = sinh(y) = sinh(y0)cosh(dy) + cosh(y0)sinh(dy) z = cosh(y) =29 // cosh(y0)cosh(dy) + sinh(y0)sinh(dy) where sinh(y0) and cosh(y0) are30 // tabulated above.31 32 __CLC_INTN ind = __CLC_CONVERT_INTN(y);33 ind = __CLC_CONVERT_UINTN(ind) > 36U ? 0 : ind;34 35 __CLC_GENTYPE dy = y - __CLC_CONVERT_GENTYPE(ind);36 __CLC_GENTYPE dy2 = dy * dy;37 38 __CLC_GENTYPE sdy = __clc_mad(39 dy2,40 __clc_mad(41 dy2,42 __clc_mad(43 dy2,44 __clc_mad(45 dy2,46 __clc_mad(dy2,47 __clc_mad(dy2, 0.7746188980094184251527126e-12f,48 0.160576793121939886190847e-9f),49 0.250521176994133472333666e-7f),50 0.275573191913636406057211e-5f),51 0.198412698413242405162014e-3f),52 0.833333333333329931873097e-2f),53 0.166666666666666667013899e0f);54 sdy = __clc_mad(sdy, dy * dy2, dy);55 56 __CLC_GENTYPE cdy = __clc_mad(57 dy2,58 __clc_mad(59 dy2,60 __clc_mad(61 dy2,62 __clc_mad(63 dy2,64 __clc_mad(dy2,65 __clc_mad(dy2, 0.1163921388172173692062032e-10f,66 0.208744349831471353536305e-8f),67 0.275573350756016588011357e-6f),68 0.248015872460622433115785e-4f),69 0.138888888889814854814536e-2f),70 0.416666666666660876512776e-1f),71 0.500000000000000005911074e0f);72 cdy = __clc_mad(cdy, dy2, 1.0f);73 74 __CLC_GENTYPE sinhcoshh = __CLC_USE_TABLE(sinhcosh_tbl_head, ind);75 __CLC_GENTYPE sinhcosht = __CLC_USE_TABLE(sinhcosh_tbl_tail, ind);76 __CLC_GENTYPE z = __clc_mad(sinhcoshh, sdy, sinhcosht * cdy);77 78 // When exp(-x) is insignificant compared to exp(x), return exp(x)/279 __CLC_GENTYPE t = __clc_exp(y - 0x1.62e500p-1f);80 __CLC_GENTYPE zsmall = __clc_mad(0x1.a0210ep-18f, t, t);81 z = y >= small_threshold ? zsmall : z;82 83 // Corner cases84 z = y >= max_cosh_arg ? __CLC_AS_GENTYPE((__CLC_UINTN)PINFBITPATT_SP32) : z;85 z = aux > PINFBITPATT_SP32 ? __CLC_GENTYPE_NAN : z;86 z = aux < 0x38800000 ? 1.0f : z;87 88 return z;89}90 91#elif __CLC_FPSIZE == 6492 93_CLC_OVERLOAD _CLC_DEF __CLC_GENTYPE __clc_cosh(__CLC_GENTYPE x) {94 // After dealing with special cases the computation is split into95 // regions as follows:96 //97 // abs(x) >= max_cosh_arg:98 // cosh(x) = sign(x)*Inf99 //100 // abs(x) >= small_threshold:101 // cosh(x) = sign(x)*exp(abs(x))/2 computed using the102 // splitexp and scaleDouble functions as for exp_amd().103 //104 // abs(x) < small_threshold:105 // compute p = exp(y) - 1 and then z = 0.5*(p+(p/(p+1.0)))106 // cosh(x) is then sign(x)*z.107 108 // This is ln(2^1025) = 0x408633ce8fb9f87e109 const __CLC_GENTYPE max_cosh_arg = 7.10475860073943977113e+02;110 111 // This is where exp(-x) is insignificant compared to exp(x) = ln(2^27)112 const __CLC_GENTYPE small_threshold = 0x1.2b708872320e2p+4;113 114 __CLC_GENTYPE y = __clc_fabs(x);115 116 // In this range we find the integer part y0 of y117 // and the increment dy = y - y0. We then compute118 // z = cosh(y) = cosh(y0)cosh(dy) + sinh(y0)sinh(dy)119 // where sinh(y0) and cosh(y0) are tabulated above.120 121 __CLC_INTN ind = __clc_min(__CLC_CONVERT_INTN(y), 36);122 __CLC_GENTYPE dy = y - __CLC_CONVERT_GENTYPE(ind);123 __CLC_GENTYPE dy2 = dy * dy;124 125 __CLC_GENTYPE sdy =126 dy * dy2 *127 __clc_fma(128 dy2,129 __clc_fma(130 dy2,131 __clc_fma(132 dy2,133 __clc_fma(134 dy2,135 __clc_fma(dy2,136 __clc_fma(dy2, 0.7746188980094184251527126e-12,137 0.160576793121939886190847e-9),138 0.250521176994133472333666e-7),139 0.275573191913636406057211e-5),140 0.198412698413242405162014e-3),141 0.833333333333329931873097e-2),142 0.166666666666666667013899e0);143 144 __CLC_GENTYPE cdy =145 dy2 *146 __clc_fma(147 dy2,148 __clc_fma(149 dy2,150 __clc_fma(151 dy2,152 __clc_fma(153 dy2,154 __clc_fma(dy2,155 __clc_fma(dy2, 0.1163921388172173692062032e-10,156 0.208744349831471353536305e-8),157 0.275573350756016588011357e-6),158 0.248015872460622433115785e-4),159 0.138888888889814854814536e-2),160 0.416666666666660876512776e-1),161 0.500000000000000005911074e0);162 163 // At this point sinh(dy) is approximated by dy + sdy,164 // and cosh(dy) is approximated by 1 + cdy.165 __CLC_GENTYPE cl = __CLC_USE_TABLE(cosh_tbl_head, ind);166 __CLC_GENTYPE ct = __CLC_USE_TABLE(cosh_tbl_tail, ind);167 __CLC_GENTYPE sl = __CLC_USE_TABLE(sinh_tbl_head, ind);168 __CLC_GENTYPE st = __CLC_USE_TABLE(sinh_tbl_tail, ind);169 170 __CLC_GENTYPE z =171 __clc_fma(172 sl, dy,173 __clc_fma(sl, sdy,174 __clc_fma(cl, cdy,175 __clc_fma(st, dy, __clc_fma(st, sdy, ct * cdy)) +176 ct))) +177 cl;178 179 // Other cases180 z = y < 0x1.0p-28 ? 1.0 : z;181 182 __CLC_GENTYPE t = __clc_exp(y - 0x1.62e42fefa3800p-1);183 t = __clc_fma(t, -0x1.ef35793c76641p-45, t);184 z = y >= small_threshold ? t : z;185 186 z = y >= max_cosh_arg ? __CLC_AS_GENTYPE((__CLC_ULONGN)PINFBITPATT_DP64) : z;187 188 z = __clc_isinf(x) || __clc_isnan(x) ? y : z;189 190 return z;191}192 193#elif __CLC_FPSIZE == 16194 195_CLC_OVERLOAD _CLC_DEF __CLC_GENTYPE __clc_cosh(__CLC_GENTYPE x) {196 return __CLC_CONVERT_GENTYPE(__clc_cosh(__CLC_CONVERT_FLOATN(x)));197}198 199#endif200