1 // Tencent is pleased to support the open source community by making RapidJSON available. 2 // 3 // Copyright (C) 2015 THL A29 Limited, a Tencent company, and Milo Yip. All rights reserved. 4 // 5 // Licensed under the MIT License (the "License"); you may not use this file except 6 // in compliance with the License. You may obtain a copy of the License at 7 // 8 // http://opensource.org/licenses/MIT 9 // 10 // Unless required by applicable law or agreed to in writing, software distributed 11 // under the License is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR 12 // CONDITIONS OF ANY KIND, either express or implied. See the License for the 13 // specific language governing permissions and limitations under the License. 14 15 #ifndef RAPIDJSON_STRTOD_ 16 #define RAPIDJSON_STRTOD_ 17 18 #include "../rapidjson.h" 19 #include "ieee754.h" 20 #include "biginteger.h" 21 #include "diyfp.h" 22 #include "pow10.h" 23 24 RAPIDJSON_NAMESPACE_BEGIN 25 namespace internal { 26 27 inline double FastPath(double significand, int exp) { 28 if (exp < -308) 29 return 0.0; 30 else if (exp >= 0) 31 return significand * internal::Pow10(exp); 32 else 33 return significand / internal::Pow10(-exp); 34 } 35 36 inline double StrtodNormalPrecision(double d, int p) { 37 if (p < -308) { 38 // Prevent expSum < -308, making Pow10(p) = 0 39 d = FastPath(d, -308); 40 d = FastPath(d, p + 308); 41 } 42 else 43 d = FastPath(d, p); 44 return d; 45 } 46 47 template <typename T> 48 inline T Min3(T a, T b, T c) { 49 T m = a; 50 if (m > b) m = b; 51 if (m > c) m = c; 52 return m; 53 } 54 55 inline int CheckWithinHalfULP(double b, const BigInteger& d, int dExp) { 56 const Double db(b); 57 const uint64_t bInt = db.IntegerSignificand(); 58 const int bExp = db.IntegerExponent(); 59 const int hExp = bExp - 1; 60 61 int dS_Exp2 = 0, dS_Exp5 = 0, bS_Exp2 = 0, bS_Exp5 = 0, hS_Exp2 = 0, hS_Exp5 = 0; 62 63 // Adjust for decimal exponent 64 if (dExp >= 0) { 65 dS_Exp2 += dExp; 66 dS_Exp5 += dExp; 67 } 68 else { 69 bS_Exp2 -= dExp; 70 bS_Exp5 -= dExp; 71 hS_Exp2 -= dExp; 72 hS_Exp5 -= dExp; 73 } 74 75 // Adjust for binary exponent 76 if (bExp >= 0) 77 bS_Exp2 += bExp; 78 else { 79 dS_Exp2 -= bExp; 80 hS_Exp2 -= bExp; 81 } 82 83 // Adjust for half ulp exponent 84 if (hExp >= 0) 85 hS_Exp2 += hExp; 86 else { 87 dS_Exp2 -= hExp; 88 bS_Exp2 -= hExp; 89 } 90 91 // Remove common power of two factor from all three scaled values 92 int common_Exp2 = Min3(dS_Exp2, bS_Exp2, hS_Exp2); 93 dS_Exp2 -= common_Exp2; 94 bS_Exp2 -= common_Exp2; 95 hS_Exp2 -= common_Exp2; 96 97 BigInteger dS = d; 98 dS.MultiplyPow5(static_cast<unsigned>(dS_Exp5)) <<= static_cast<unsigned>(dS_Exp2); 99 100 BigInteger bS(bInt); 101 bS.MultiplyPow5(static_cast<unsigned>(bS_Exp5)) <<= static_cast<unsigned>(bS_Exp2); 102 103 BigInteger hS(1); 104 hS.MultiplyPow5(static_cast<unsigned>(hS_Exp5)) <<= static_cast<unsigned>(hS_Exp2); 105 106 BigInteger delta(0); 107 dS.Difference(bS, &delta); 108 109 return delta.Compare(hS); 110 } 111 112 inline bool StrtodFast(double d, int p, double* result) { 113 // Use fast path for string-to-double conversion if possible 114 // see http://www.exploringbinary.com/fast-path-decimal-to-floating-point-conversion/ 115 if (p > 22 && p < 22 + 16) { 116 // Fast Path Cases In Disguise 117 d *= internal::Pow10(p - 22); 118 p = 22; 119 } 120 121 if (p >= -22 && p <= 22 && d <= 9007199254740991.0) { // 2^53 - 1 122 *result = FastPath(d, p); 123 return true; 124 } 125 else 126 return false; 127 } 128 129 // Compute an approximation and see if it is within 1/2 ULP 130 inline bool StrtodDiyFp(const char* decimals, size_t length, size_t decimalPosition, int exp, double* result) { 131 uint64_t significand = 0; 132 size_t i = 0; // 2^64 - 1 = 18446744073709551615, 1844674407370955161 = 0x1999999999999999 133 for (; i < length; i++) { 134 if (significand > RAPIDJSON_UINT64_C2(0x19999999, 0x99999999) || 135 (significand == RAPIDJSON_UINT64_C2(0x19999999, 0x99999999) && decimals[i] > '5')) 136 break; 137 significand = significand * 10u + static_cast<unsigned>(decimals[i] - '0'); 138 } 139 140 if (i < length && decimals[i] >= '5') // Rounding 141 significand++; 142 143 size_t remaining = length - i; 144 const unsigned kUlpShift = 3; 145 const unsigned kUlp = 1 << kUlpShift; 146 int error = (remaining == 0) ? 0 : kUlp / 2; 147 148 DiyFp v(significand, 0); 149 v = v.Normalize(); 150 error <<= -v.e; 151 152 const int dExp = (int)decimalPosition - (int)i + exp; 153 154 int actualExp; 155 DiyFp cachedPower = GetCachedPower10(dExp, &actualExp); 156 if (actualExp != dExp) { 157 static const DiyFp kPow10[] = { 158 DiyFp(RAPIDJSON_UINT64_C2(0xa0000000, 00000000), -60), // 10^1 159 DiyFp(RAPIDJSON_UINT64_C2(0xc8000000, 00000000), -57), // 10^2 160 DiyFp(RAPIDJSON_UINT64_C2(0xfa000000, 00000000), -54), // 10^3 161 DiyFp(RAPIDJSON_UINT64_C2(0x9c400000, 00000000), -50), // 10^4 162 DiyFp(RAPIDJSON_UINT64_C2(0xc3500000, 00000000), -47), // 10^5 163 DiyFp(RAPIDJSON_UINT64_C2(0xf4240000, 00000000), -44), // 10^6 164 DiyFp(RAPIDJSON_UINT64_C2(0x98968000, 00000000), -40) // 10^7 165 }; 166 int adjustment = dExp - actualExp - 1; 167 RAPIDJSON_ASSERT(adjustment >= 0 && adjustment < 7); 168 v = v * kPow10[adjustment]; 169 if (length + static_cast<unsigned>(adjustment)> 19u) // has more digits than decimal digits in 64-bit 170 error += kUlp / 2; 171 } 172 173 v = v * cachedPower; 174 175 error += kUlp + (error == 0 ? 0 : 1); 176 177 const int oldExp = v.e; 178 v = v.Normalize(); 179 error <<= oldExp - v.e; 180 181 const unsigned effectiveSignificandSize = Double::EffectiveSignificandSize(64 + v.e); 182 unsigned precisionSize = 64 - effectiveSignificandSize; 183 if (precisionSize + kUlpShift >= 64) { 184 unsigned scaleExp = (precisionSize + kUlpShift) - 63; 185 v.f >>= scaleExp; 186 v.e += scaleExp; 187 error = (error >> scaleExp) + 1 + static_cast<int>(kUlp); 188 precisionSize -= scaleExp; 189 } 190 191 DiyFp rounded(v.f >> precisionSize, v.e + static_cast<int>(precisionSize)); 192 const uint64_t precisionBits = (v.f & ((uint64_t(1) << precisionSize) - 1)) * kUlp; 193 const uint64_t halfWay = (uint64_t(1) << (precisionSize - 1)) * kUlp; 194 if (precisionBits >= halfWay + static_cast<unsigned>(error)) { 195 rounded.f++; 196 if (rounded.f & (DiyFp::kDpHiddenBit << 1)) { // rounding overflows mantissa (issue #340) 197 rounded.f >>= 1; 198 rounded.e++; 199 } 200 } 201 202 *result = rounded.ToDouble(); 203 204 return halfWay - static_cast<unsigned>(error) >= precisionBits || precisionBits >= halfWay + static_cast<unsigned>(error); 205 } 206 207 inline double StrtodBigInteger(double approx, const char* decimals, size_t length, size_t decimalPosition, int exp) { 208 const BigInteger dInt(decimals, length); 209 const int dExp = (int)decimalPosition - (int)length + exp; 210 Double a(approx); 211 int cmp = CheckWithinHalfULP(a.Value(), dInt, dExp); 212 if (cmp < 0) 213 return a.Value(); // within half ULP 214 else if (cmp == 0) { 215 // Round towards even 216 if (a.Significand() & 1) 217 return a.NextPositiveDouble(); 218 else 219 return a.Value(); 220 } 221 else // adjustment 222 return a.NextPositiveDouble(); 223 } 224 225 inline double StrtodFullPrecision(double d, int p, const char* decimals, size_t length, size_t decimalPosition, int exp) { 226 RAPIDJSON_ASSERT(d >= 0.0); 227 RAPIDJSON_ASSERT(length >= 1); 228 229 double result; 230 if (StrtodFast(d, p, &result)) 231 return result; 232 233 // Trim leading zeros 234 while (*decimals == '0' && length > 1) { 235 length--; 236 decimals++; 237 decimalPosition--; 238 } 239 240 // Trim trailing zeros 241 while (decimals[length - 1] == '0' && length > 1) { 242 length--; 243 decimalPosition--; 244 exp++; 245 } 246 247 // Trim right-most digits 248 const int kMaxDecimalDigit = 780; 249 if ((int)length > kMaxDecimalDigit) { 250 int delta = (int(length) - kMaxDecimalDigit); 251 exp += delta; 252 decimalPosition -= static_cast<unsigned>(delta); 253 length = kMaxDecimalDigit; 254 } 255 256 // If too small, underflow to zero 257 if (int(length) + exp < -324) 258 return 0.0; 259 260 if (StrtodDiyFp(decimals, length, decimalPosition, exp, &result)) 261 return result; 262 263 // Use approximation from StrtodDiyFp and make adjustment with BigInteger comparison 264 return StrtodBigInteger(result, decimals, length, decimalPosition, exp); 265 } 266 267 } // namespace internal 268 RAPIDJSON_NAMESPACE_END 269 270 #endif // RAPIDJSON_STRTOD_ 271
