Stephen Hines | 2d1fdb2 | 2014-05-28 23:58:16 -0700 | [diff] [blame] | 1 | //===-- lib/fp_lib.h - Floating-point utilities -------------------*- C -*-===// |
| 2 | // |
| 3 | // The LLVM Compiler Infrastructure |
| 4 | // |
| 5 | // This file is dual licensed under the MIT and the University of Illinois Open |
| 6 | // Source Licenses. See LICENSE.TXT for details. |
| 7 | // |
| 8 | //===----------------------------------------------------------------------===// |
| 9 | // |
| 10 | // This file is a configuration header for soft-float routines in compiler-rt. |
| 11 | // This file does not provide any part of the compiler-rt interface, but defines |
| 12 | // many useful constants and utility routines that are used in the |
| 13 | // implementation of the soft-float routines in compiler-rt. |
| 14 | // |
| 15 | // Assumes that float, double and long double correspond to the IEEE-754 |
| 16 | // binary32, binary64 and binary 128 types, respectively, and that integer |
| 17 | // endianness matches floating point endianness on the target platform. |
| 18 | // |
| 19 | //===----------------------------------------------------------------------===// |
| 20 | |
| 21 | #ifndef FP_LIB_HEADER |
| 22 | #define FP_LIB_HEADER |
| 23 | |
| 24 | #include <stdint.h> |
| 25 | #include <stdbool.h> |
| 26 | #include <limits.h> |
| 27 | #include "int_lib.h" |
| 28 | |
Stephen Hines | 6a211c5 | 2014-07-21 00:49:56 -0700 | [diff] [blame] | 29 | // x86_64 FreeBSD prior v9.3 define fixed-width types incorrectly in |
| 30 | // 32-bit mode. |
| 31 | #if defined(__FreeBSD__) && defined(__i386__) |
| 32 | # include <sys/param.h> |
| 33 | # if __FreeBSD_version < 903000 // v9.3 |
| 34 | # define uint64_t unsigned long long |
| 35 | # define int64_t long long |
| 36 | # undef UINT64_C |
| 37 | # define UINT64_C(c) (c ## ULL) |
| 38 | # endif |
| 39 | #endif |
| 40 | |
Stephen Hines | 2d1fdb2 | 2014-05-28 23:58:16 -0700 | [diff] [blame] | 41 | #if defined SINGLE_PRECISION |
| 42 | |
| 43 | typedef uint32_t rep_t; |
| 44 | typedef int32_t srep_t; |
| 45 | typedef float fp_t; |
| 46 | #define REP_C UINT32_C |
| 47 | #define significandBits 23 |
| 48 | |
Pirama Arumuga Nainar | 799172d | 2016-03-03 15:50:30 -0800 | [diff] [blame^] | 49 | static __inline int rep_clz(rep_t a) { |
Stephen Hines | 2d1fdb2 | 2014-05-28 23:58:16 -0700 | [diff] [blame] | 50 | return __builtin_clz(a); |
| 51 | } |
| 52 | |
| 53 | // 32x32 --> 64 bit multiply |
Pirama Arumuga Nainar | 799172d | 2016-03-03 15:50:30 -0800 | [diff] [blame^] | 54 | static __inline void wideMultiply(rep_t a, rep_t b, rep_t *hi, rep_t *lo) { |
Stephen Hines | 2d1fdb2 | 2014-05-28 23:58:16 -0700 | [diff] [blame] | 55 | const uint64_t product = (uint64_t)a*b; |
| 56 | *hi = product >> 32; |
| 57 | *lo = product; |
| 58 | } |
| 59 | COMPILER_RT_ABI fp_t __addsf3(fp_t a, fp_t b); |
| 60 | |
| 61 | #elif defined DOUBLE_PRECISION |
| 62 | |
| 63 | typedef uint64_t rep_t; |
| 64 | typedef int64_t srep_t; |
| 65 | typedef double fp_t; |
| 66 | #define REP_C UINT64_C |
| 67 | #define significandBits 52 |
| 68 | |
Pirama Arumuga Nainar | 799172d | 2016-03-03 15:50:30 -0800 | [diff] [blame^] | 69 | static __inline int rep_clz(rep_t a) { |
Stephen Hines | 2d1fdb2 | 2014-05-28 23:58:16 -0700 | [diff] [blame] | 70 | #if defined __LP64__ |
| 71 | return __builtin_clzl(a); |
| 72 | #else |
| 73 | if (a & REP_C(0xffffffff00000000)) |
| 74 | return __builtin_clz(a >> 32); |
| 75 | else |
| 76 | return 32 + __builtin_clz(a & REP_C(0xffffffff)); |
| 77 | #endif |
| 78 | } |
| 79 | |
| 80 | #define loWord(a) (a & 0xffffffffU) |
| 81 | #define hiWord(a) (a >> 32) |
| 82 | |
| 83 | // 64x64 -> 128 wide multiply for platforms that don't have such an operation; |
| 84 | // many 64-bit platforms have this operation, but they tend to have hardware |
| 85 | // floating-point, so we don't bother with a special case for them here. |
Pirama Arumuga Nainar | 799172d | 2016-03-03 15:50:30 -0800 | [diff] [blame^] | 86 | static __inline void wideMultiply(rep_t a, rep_t b, rep_t *hi, rep_t *lo) { |
Stephen Hines | 2d1fdb2 | 2014-05-28 23:58:16 -0700 | [diff] [blame] | 87 | // Each of the component 32x32 -> 64 products |
| 88 | const uint64_t plolo = loWord(a) * loWord(b); |
| 89 | const uint64_t plohi = loWord(a) * hiWord(b); |
| 90 | const uint64_t philo = hiWord(a) * loWord(b); |
| 91 | const uint64_t phihi = hiWord(a) * hiWord(b); |
| 92 | // Sum terms that contribute to lo in a way that allows us to get the carry |
| 93 | const uint64_t r0 = loWord(plolo); |
| 94 | const uint64_t r1 = hiWord(plolo) + loWord(plohi) + loWord(philo); |
| 95 | *lo = r0 + (r1 << 32); |
| 96 | // Sum terms contributing to hi with the carry from lo |
| 97 | *hi = hiWord(plohi) + hiWord(philo) + hiWord(r1) + phihi; |
| 98 | } |
| 99 | #undef loWord |
| 100 | #undef hiWord |
| 101 | |
| 102 | COMPILER_RT_ABI fp_t __adddf3(fp_t a, fp_t b); |
| 103 | |
| 104 | #elif defined QUAD_PRECISION |
| 105 | #if __LDBL_MANT_DIG__ == 113 |
| 106 | #define CRT_LDBL_128BIT |
| 107 | typedef __uint128_t rep_t; |
| 108 | typedef __int128_t srep_t; |
| 109 | typedef long double fp_t; |
| 110 | #define REP_C (__uint128_t) |
| 111 | // Note: Since there is no explicit way to tell compiler the constant is a |
| 112 | // 128-bit integer, we let the constant be casted to 128-bit integer |
| 113 | #define significandBits 112 |
| 114 | |
Pirama Arumuga Nainar | 799172d | 2016-03-03 15:50:30 -0800 | [diff] [blame^] | 115 | static __inline int rep_clz(rep_t a) { |
Stephen Hines | 2d1fdb2 | 2014-05-28 23:58:16 -0700 | [diff] [blame] | 116 | const union |
| 117 | { |
| 118 | __uint128_t ll; |
| 119 | #if _YUGA_BIG_ENDIAN |
| 120 | struct { uint64_t high, low; } s; |
| 121 | #else |
| 122 | struct { uint64_t low, high; } s; |
| 123 | #endif |
| 124 | } uu = { .ll = a }; |
| 125 | |
| 126 | uint64_t word; |
| 127 | uint64_t add; |
| 128 | |
| 129 | if (uu.s.high){ |
| 130 | word = uu.s.high; |
| 131 | add = 0; |
| 132 | } |
| 133 | else{ |
| 134 | word = uu.s.low; |
| 135 | add = 64; |
| 136 | } |
| 137 | return __builtin_clzll(word) + add; |
| 138 | } |
| 139 | |
| 140 | #define Word_LoMask UINT64_C(0x00000000ffffffff) |
| 141 | #define Word_HiMask UINT64_C(0xffffffff00000000) |
| 142 | #define Word_FullMask UINT64_C(0xffffffffffffffff) |
| 143 | #define Word_1(a) (uint64_t)((a >> 96) & Word_LoMask) |
| 144 | #define Word_2(a) (uint64_t)((a >> 64) & Word_LoMask) |
| 145 | #define Word_3(a) (uint64_t)((a >> 32) & Word_LoMask) |
| 146 | #define Word_4(a) (uint64_t)(a & Word_LoMask) |
| 147 | |
| 148 | // 128x128 -> 256 wide multiply for platforms that don't have such an operation; |
| 149 | // many 64-bit platforms have this operation, but they tend to have hardware |
| 150 | // floating-point, so we don't bother with a special case for them here. |
Pirama Arumuga Nainar | 799172d | 2016-03-03 15:50:30 -0800 | [diff] [blame^] | 151 | static __inline void wideMultiply(rep_t a, rep_t b, rep_t *hi, rep_t *lo) { |
Stephen Hines | 2d1fdb2 | 2014-05-28 23:58:16 -0700 | [diff] [blame] | 152 | |
| 153 | const uint64_t product11 = Word_1(a) * Word_1(b); |
| 154 | const uint64_t product12 = Word_1(a) * Word_2(b); |
| 155 | const uint64_t product13 = Word_1(a) * Word_3(b); |
| 156 | const uint64_t product14 = Word_1(a) * Word_4(b); |
| 157 | const uint64_t product21 = Word_2(a) * Word_1(b); |
| 158 | const uint64_t product22 = Word_2(a) * Word_2(b); |
| 159 | const uint64_t product23 = Word_2(a) * Word_3(b); |
| 160 | const uint64_t product24 = Word_2(a) * Word_4(b); |
| 161 | const uint64_t product31 = Word_3(a) * Word_1(b); |
| 162 | const uint64_t product32 = Word_3(a) * Word_2(b); |
| 163 | const uint64_t product33 = Word_3(a) * Word_3(b); |
| 164 | const uint64_t product34 = Word_3(a) * Word_4(b); |
| 165 | const uint64_t product41 = Word_4(a) * Word_1(b); |
| 166 | const uint64_t product42 = Word_4(a) * Word_2(b); |
| 167 | const uint64_t product43 = Word_4(a) * Word_3(b); |
| 168 | const uint64_t product44 = Word_4(a) * Word_4(b); |
| 169 | |
| 170 | const __uint128_t sum0 = (__uint128_t)product44; |
| 171 | const __uint128_t sum1 = (__uint128_t)product34 + |
| 172 | (__uint128_t)product43; |
| 173 | const __uint128_t sum2 = (__uint128_t)product24 + |
| 174 | (__uint128_t)product33 + |
| 175 | (__uint128_t)product42; |
| 176 | const __uint128_t sum3 = (__uint128_t)product14 + |
| 177 | (__uint128_t)product23 + |
| 178 | (__uint128_t)product32 + |
| 179 | (__uint128_t)product41; |
| 180 | const __uint128_t sum4 = (__uint128_t)product13 + |
| 181 | (__uint128_t)product22 + |
| 182 | (__uint128_t)product31; |
| 183 | const __uint128_t sum5 = (__uint128_t)product12 + |
| 184 | (__uint128_t)product21; |
| 185 | const __uint128_t sum6 = (__uint128_t)product11; |
| 186 | |
| 187 | const __uint128_t r0 = (sum0 & Word_FullMask) + |
| 188 | ((sum1 & Word_LoMask) << 32); |
| 189 | const __uint128_t r1 = (sum0 >> 64) + |
| 190 | ((sum1 >> 32) & Word_FullMask) + |
| 191 | (sum2 & Word_FullMask) + |
| 192 | ((sum3 << 32) & Word_HiMask); |
| 193 | |
| 194 | *lo = r0 + (r1 << 64); |
| 195 | *hi = (r1 >> 64) + |
| 196 | (sum1 >> 96) + |
| 197 | (sum2 >> 64) + |
| 198 | (sum3 >> 32) + |
| 199 | sum4 + |
| 200 | (sum5 << 32) + |
| 201 | (sum6 << 64); |
| 202 | } |
| 203 | #undef Word_1 |
| 204 | #undef Word_2 |
| 205 | #undef Word_3 |
| 206 | #undef Word_4 |
| 207 | #undef Word_HiMask |
| 208 | #undef Word_LoMask |
| 209 | #undef Word_FullMask |
| 210 | #endif // __LDBL_MANT_DIG__ == 113 |
| 211 | #else |
| 212 | #error SINGLE_PRECISION, DOUBLE_PRECISION or QUAD_PRECISION must be defined. |
| 213 | #endif |
| 214 | |
| 215 | #if defined(SINGLE_PRECISION) || defined(DOUBLE_PRECISION) || defined(CRT_LDBL_128BIT) |
| 216 | #define typeWidth (sizeof(rep_t)*CHAR_BIT) |
| 217 | #define exponentBits (typeWidth - significandBits - 1) |
| 218 | #define maxExponent ((1 << exponentBits) - 1) |
| 219 | #define exponentBias (maxExponent >> 1) |
| 220 | |
| 221 | #define implicitBit (REP_C(1) << significandBits) |
| 222 | #define significandMask (implicitBit - 1U) |
| 223 | #define signBit (REP_C(1) << (significandBits + exponentBits)) |
| 224 | #define absMask (signBit - 1U) |
| 225 | #define exponentMask (absMask ^ significandMask) |
| 226 | #define oneRep ((rep_t)exponentBias << significandBits) |
| 227 | #define infRep exponentMask |
| 228 | #define quietBit (implicitBit >> 1) |
| 229 | #define qnanRep (exponentMask | quietBit) |
| 230 | |
Pirama Arumuga Nainar | 799172d | 2016-03-03 15:50:30 -0800 | [diff] [blame^] | 231 | static __inline rep_t toRep(fp_t x) { |
Stephen Hines | 2d1fdb2 | 2014-05-28 23:58:16 -0700 | [diff] [blame] | 232 | const union { fp_t f; rep_t i; } rep = {.f = x}; |
| 233 | return rep.i; |
| 234 | } |
| 235 | |
Pirama Arumuga Nainar | 799172d | 2016-03-03 15:50:30 -0800 | [diff] [blame^] | 236 | static __inline fp_t fromRep(rep_t x) { |
Stephen Hines | 2d1fdb2 | 2014-05-28 23:58:16 -0700 | [diff] [blame] | 237 | const union { fp_t f; rep_t i; } rep = {.i = x}; |
| 238 | return rep.f; |
| 239 | } |
| 240 | |
Pirama Arumuga Nainar | 799172d | 2016-03-03 15:50:30 -0800 | [diff] [blame^] | 241 | static __inline int normalize(rep_t *significand) { |
Stephen Hines | 2d1fdb2 | 2014-05-28 23:58:16 -0700 | [diff] [blame] | 242 | const int shift = rep_clz(*significand) - rep_clz(implicitBit); |
| 243 | *significand <<= shift; |
| 244 | return 1 - shift; |
| 245 | } |
| 246 | |
Pirama Arumuga Nainar | 799172d | 2016-03-03 15:50:30 -0800 | [diff] [blame^] | 247 | static __inline void wideLeftShift(rep_t *hi, rep_t *lo, int count) { |
Stephen Hines | 2d1fdb2 | 2014-05-28 23:58:16 -0700 | [diff] [blame] | 248 | *hi = *hi << count | *lo >> (typeWidth - count); |
| 249 | *lo = *lo << count; |
| 250 | } |
| 251 | |
Pirama Arumuga Nainar | 799172d | 2016-03-03 15:50:30 -0800 | [diff] [blame^] | 252 | static __inline void wideRightShiftWithSticky(rep_t *hi, rep_t *lo, unsigned int count) { |
Stephen Hines | 2d1fdb2 | 2014-05-28 23:58:16 -0700 | [diff] [blame] | 253 | if (count < typeWidth) { |
| 254 | const bool sticky = *lo << (typeWidth - count); |
| 255 | *lo = *hi << (typeWidth - count) | *lo >> count | sticky; |
| 256 | *hi = *hi >> count; |
| 257 | } |
| 258 | else if (count < 2*typeWidth) { |
| 259 | const bool sticky = *hi << (2*typeWidth - count) | *lo; |
| 260 | *lo = *hi >> (count - typeWidth) | sticky; |
| 261 | *hi = 0; |
| 262 | } else { |
| 263 | const bool sticky = *hi | *lo; |
| 264 | *lo = sticky; |
| 265 | *hi = 0; |
| 266 | } |
| 267 | } |
| 268 | #endif |
| 269 | |
| 270 | #endif // FP_LIB_HEADER |