1 //==- lib/Support/ScaledNumber.cpp - Support for scaled numbers -*- C++ -*-===// 2 // 3 // The LLVM Compiler Infrastructure 4 // 5 // This file is distributed under the University of Illinois Open Source 6 // License. See LICENSE.TXT for details. 7 // 8 //===----------------------------------------------------------------------===// 9 // 10 // Implementation of some scaled number algorithms. 11 // 12 //===----------------------------------------------------------------------===// 13 14 #include "llvm/Support/ScaledNumber.h" 15 #include "llvm/ADT/APFloat.h" 16 #include "llvm/ADT/ArrayRef.h" 17 #include "llvm/Support/Debug.h" 18 #include "llvm/Support/raw_ostream.h" 19 20 using namespace llvm; 21 using namespace llvm::ScaledNumbers; 22 23 std::pair<uint64_t, int16_t> ScaledNumbers::multiply64(uint64_t LHS, 24 uint64_t RHS) { 25 // Separate into two 32-bit digits (U.L). 26 auto getU = [](uint64_t N) { return N >> 32; }; 27 auto getL = [](uint64_t N) { return N & UINT32_MAX; }; 28 uint64_t UL = getU(LHS), LL = getL(LHS), UR = getU(RHS), LR = getL(RHS); 29 30 // Compute cross products. 31 uint64_t P1 = UL * UR, P2 = UL * LR, P3 = LL * UR, P4 = LL * LR; 32 33 // Sum into two 64-bit digits. 34 uint64_t Upper = P1, Lower = P4; 35 auto addWithCarry = [&](uint64_t N) { 36 uint64_t NewLower = Lower + (getL(N) << 32); 37 Upper += getU(N) + (NewLower < Lower); 38 Lower = NewLower; 39 }; 40 addWithCarry(P2); 41 addWithCarry(P3); 42 43 // Check whether the upper digit is empty. 44 if (!Upper) 45 return std::make_pair(Lower, 0); 46 47 // Shift as little as possible to maximize precision. 48 unsigned LeadingZeros = countLeadingZeros(Upper); 49 int Shift = 64 - LeadingZeros; 50 if (LeadingZeros) 51 Upper = Upper << LeadingZeros | Lower >> Shift; 52 return getRounded(Upper, Shift, 53 Shift && (Lower & UINT64_C(1) << (Shift - 1))); 54 } 55 56 static uint64_t getHalf(uint64_t N) { return (N >> 1) + (N & 1); } 57 58 std::pair<uint32_t, int16_t> ScaledNumbers::divide32(uint32_t Dividend, 59 uint32_t Divisor) { 60 assert(Dividend && "expected non-zero dividend"); 61 assert(Divisor && "expected non-zero divisor"); 62 63 // Use 64-bit math and canonicalize the dividend to gain precision. 64 uint64_t Dividend64 = Dividend; 65 int Shift = 0; 66 if (int Zeros = countLeadingZeros(Dividend64)) { 67 Shift -= Zeros; 68 Dividend64 <<= Zeros; 69 } 70 uint64_t Quotient = Dividend64 / Divisor; 71 uint64_t Remainder = Dividend64 % Divisor; 72 73 // If Quotient needs to be shifted, leave the rounding to getAdjusted(). 74 if (Quotient > UINT32_MAX) 75 return getAdjusted<uint32_t>(Quotient, Shift); 76 77 // Round based on the value of the next bit. 78 return getRounded<uint32_t>(Quotient, Shift, Remainder >= getHalf(Divisor)); 79 } 80 81 std::pair<uint64_t, int16_t> ScaledNumbers::divide64(uint64_t Dividend, 82 uint64_t Divisor) { 83 assert(Dividend && "expected non-zero dividend"); 84 assert(Divisor && "expected non-zero divisor"); 85 86 // Minimize size of divisor. 87 int Shift = 0; 88 if (int Zeros = countTrailingZeros(Divisor)) { 89 Shift -= Zeros; 90 Divisor >>= Zeros; 91 } 92 93 // Check for powers of two. 94 if (Divisor == 1) 95 return std::make_pair(Dividend, Shift); 96 97 // Maximize size of dividend. 98 if (int Zeros = countLeadingZeros(Dividend)) { 99 Shift -= Zeros; 100 Dividend <<= Zeros; 101 } 102 103 // Start with the result of a divide. 104 uint64_t Quotient = Dividend / Divisor; 105 Dividend %= Divisor; 106 107 // Continue building the quotient with long division. 108 while (!(Quotient >> 63) && Dividend) { 109 // Shift Dividend and check for overflow. 110 bool IsOverflow = Dividend >> 63; 111 Dividend <<= 1; 112 --Shift; 113 114 // Get the next bit of Quotient. 115 Quotient <<= 1; 116 if (IsOverflow || Divisor <= Dividend) { 117 Quotient |= 1; 118 Dividend -= Divisor; 119 } 120 } 121 122 return getRounded(Quotient, Shift, Dividend >= getHalf(Divisor)); 123 } 124 125 int ScaledNumbers::compareImpl(uint64_t L, uint64_t R, int ScaleDiff) { 126 assert(ScaleDiff >= 0 && "wrong argument order"); 127 assert(ScaleDiff < 64 && "numbers too far apart"); 128 129 uint64_t L_adjusted = L >> ScaleDiff; 130 if (L_adjusted < R) 131 return -1; 132 if (L_adjusted > R) 133 return 1; 134 135 return L > L_adjusted << ScaleDiff ? 1 : 0; 136 } 137 138 static void appendDigit(std::string &Str, unsigned D) { 139 assert(D < 10); 140 Str += '0' + D % 10; 141 } 142 143 static void appendNumber(std::string &Str, uint64_t N) { 144 while (N) { 145 appendDigit(Str, N % 10); 146 N /= 10; 147 } 148 } 149 150 static bool doesRoundUp(char Digit) { 151 switch (Digit) { 152 case '5': 153 case '6': 154 case '7': 155 case '8': 156 case '9': 157 return true; 158 default: 159 return false; 160 } 161 } 162 163 static std::string toStringAPFloat(uint64_t D, int E, unsigned Precision) { 164 assert(E >= ScaledNumbers::MinScale); 165 assert(E <= ScaledNumbers::MaxScale); 166 167 // Find a new E, but don't let it increase past MaxScale. 168 int LeadingZeros = ScaledNumberBase::countLeadingZeros64(D); 169 int NewE = std::min(ScaledNumbers::MaxScale, E + 63 - LeadingZeros); 170 int Shift = 63 - (NewE - E); 171 assert(Shift <= LeadingZeros); 172 assert(Shift == LeadingZeros || NewE == ScaledNumbers::MaxScale); 173 assert(Shift >= 0 && Shift < 64 && "undefined behavior"); 174 D <<= Shift; 175 E = NewE; 176 177 // Check for a denormal. 178 unsigned AdjustedE = E + 16383; 179 if (!(D >> 63)) { 180 assert(E == ScaledNumbers::MaxScale); 181 AdjustedE = 0; 182 } 183 184 // Build the float and print it. 185 uint64_t RawBits[2] = {D, AdjustedE}; 186 APFloat Float(APFloat::x87DoubleExtended, APInt(80, RawBits)); 187 SmallVector<char, 24> Chars; 188 Float.toString(Chars, Precision, 0); 189 return std::string(Chars.begin(), Chars.end()); 190 } 191 192 static std::string stripTrailingZeros(const std::string &Float) { 193 size_t NonZero = Float.find_last_not_of('0'); 194 assert(NonZero != std::string::npos && "no . in floating point string"); 195 196 if (Float[NonZero] == '.') 197 ++NonZero; 198 199 return Float.substr(0, NonZero + 1); 200 } 201 202 std::string ScaledNumberBase::toString(uint64_t D, int16_t E, int Width, 203 unsigned Precision) { 204 if (!D) 205 return "0.0"; 206 207 // Canonicalize exponent and digits. 208 uint64_t Above0 = 0; 209 uint64_t Below0 = 0; 210 uint64_t Extra = 0; 211 int ExtraShift = 0; 212 if (E == 0) { 213 Above0 = D; 214 } else if (E > 0) { 215 if (int Shift = std::min(int16_t(countLeadingZeros64(D)), E)) { 216 D <<= Shift; 217 E -= Shift; 218 219 if (!E) 220 Above0 = D; 221 } 222 } else if (E > -64) { 223 Above0 = D >> -E; 224 Below0 = D << (64 + E); 225 } else if (E == -64) { 226 // Special case: shift by 64 bits is undefined behavior. 227 Below0 = D; 228 } else if (E > -120) { 229 Below0 = D >> (-E - 64); 230 Extra = D << (128 + E); 231 ExtraShift = -64 - E; 232 } 233 234 // Fall back on APFloat for very small and very large numbers. 235 if (!Above0 && !Below0) 236 return toStringAPFloat(D, E, Precision); 237 238 // Append the digits before the decimal. 239 std::string Str; 240 size_t DigitsOut = 0; 241 if (Above0) { 242 appendNumber(Str, Above0); 243 DigitsOut = Str.size(); 244 } else 245 appendDigit(Str, 0); 246 std::reverse(Str.begin(), Str.end()); 247 248 // Return early if there's nothing after the decimal. 249 if (!Below0) 250 return Str + ".0"; 251 252 // Append the decimal and beyond. 253 Str += '.'; 254 uint64_t Error = UINT64_C(1) << (64 - Width); 255 256 // We need to shift Below0 to the right to make space for calculating 257 // digits. Save the precision we're losing in Extra. 258 Extra = (Below0 & 0xf) << 56 | (Extra >> 8); 259 Below0 >>= 4; 260 size_t SinceDot = 0; 261 size_t AfterDot = Str.size(); 262 do { 263 if (ExtraShift) { 264 --ExtraShift; 265 Error *= 5; 266 } else 267 Error *= 10; 268 269 Below0 *= 10; 270 Extra *= 10; 271 Below0 += (Extra >> 60); 272 Extra = Extra & (UINT64_MAX >> 4); 273 appendDigit(Str, Below0 >> 60); 274 Below0 = Below0 & (UINT64_MAX >> 4); 275 if (DigitsOut || Str.back() != '0') 276 ++DigitsOut; 277 ++SinceDot; 278 } while (Error && (Below0 << 4 | Extra >> 60) >= Error / 2 && 279 (!Precision || DigitsOut <= Precision || SinceDot < 2)); 280 281 // Return early for maximum precision. 282 if (!Precision || DigitsOut <= Precision) 283 return stripTrailingZeros(Str); 284 285 // Find where to truncate. 286 size_t Truncate = 287 std::max(Str.size() - (DigitsOut - Precision), AfterDot + 1); 288 289 // Check if there's anything to truncate. 290 if (Truncate >= Str.size()) 291 return stripTrailingZeros(Str); 292 293 bool Carry = doesRoundUp(Str[Truncate]); 294 if (!Carry) 295 return stripTrailingZeros(Str.substr(0, Truncate)); 296 297 // Round with the first truncated digit. 298 for (std::string::reverse_iterator I(Str.begin() + Truncate), E = Str.rend(); 299 I != E; ++I) { 300 if (*I == '.') 301 continue; 302 if (*I == '9') { 303 *I = '0'; 304 continue; 305 } 306 307 ++*I; 308 Carry = false; 309 break; 310 } 311 312 // Add "1" in front if we still need to carry. 313 return stripTrailingZeros(std::string(Carry, '1') + Str.substr(0, Truncate)); 314 } 315 316 raw_ostream &ScaledNumberBase::print(raw_ostream &OS, uint64_t D, int16_t E, 317 int Width, unsigned Precision) { 318 return OS << toString(D, E, Width, Precision); 319 } 320 321 void ScaledNumberBase::dump(uint64_t D, int16_t E, int Width) { 322 print(dbgs(), D, E, Width, 0) << "[" << Width << ":" << D << "*2^" << E 323 << "]"; 324 } 325