19 #ifndef RAPIDJSON_DIYFP_H_
20 #define RAPIDJSON_DIYFP_H_
22 #include "../rapidjson.h"
26 #if defined(_MSC_VER) && defined(_M_AMD64) && !defined(__INTEL_COMPILER)
28 #if !defined(_ARM64EC_)
29 #pragma intrinsic(_umul128)
31 #pragma comment(lib,"softintrin")
35 RAPIDJSON_NAMESPACE_BEGIN
40 RAPIDJSON_DIAG_OFF(effc++)
45 RAPIDJSON_DIAG_OFF(padded)
51 DiyFp(uint64_t fp,
int exp) : f(fp), e(exp) {}
53 explicit DiyFp(
double d) {
59 int biased_e =
static_cast<int>((u.u64 & kDpExponentMask) >> kDpSignificandSize);
60 uint64_t significand = (u.u64 & kDpSignificandMask);
62 f = significand + kDpHiddenBit;
63 e = biased_e - kDpExponentBias;
67 e = kDpMinExponent + 1;
71 DiyFp operator-(
const DiyFp& rhs)
const {
72 return DiyFp(f - rhs.f, e);
75 DiyFp operator*(
const DiyFp& rhs)
const {
76 #if defined(_MSC_VER) && defined(_M_AMD64)
78 uint64_t l = _umul128(f, rhs.f, &h);
79 if (l & (uint64_t(1) << 63))
81 return DiyFp(h, e + rhs.e + 64);
82 #elif defined(__GNUC__) && (__GNUC__ > 4 || (__GNUC__ == 4 && __GNUC_MINOR__ >= 6)) && defined(__x86_64__)
83 __extension__
typedef unsigned __int128 uint128;
84 uint128 p =
static_cast<uint128
>(f) * static_cast<uint128>(rhs.f);
85 uint64_t h =
static_cast<uint64_t
>(p >> 64);
86 uint64_t l =
static_cast<uint64_t
>(p);
87 if (l & (uint64_t(1) << 63))
89 return DiyFp(h, e + rhs.e + 64);
91 const uint64_t M32 = 0xFFFFFFFF;
92 const uint64_t a = f >> 32;
93 const uint64_t b = f & M32;
94 const uint64_t c = rhs.f >> 32;
95 const uint64_t d = rhs.f & M32;
96 const uint64_t ac = a * c;
97 const uint64_t bc = b * c;
98 const uint64_t ad = a * d;
99 const uint64_t bd = b * d;
100 uint64_t tmp = (bd >> 32) + (ad & M32) + (bc & M32);
102 return DiyFp(ac + (ad >> 32) + (bc >> 32) + (tmp >> 32), e + rhs.e + 64);
106 DiyFp Normalize()
const {
107 int s =
static_cast<int>(clzll(f));
108 return DiyFp(f << s, e - s);
111 DiyFp NormalizeBoundary()
const {
113 while (!(res.f & (kDpHiddenBit << 1))) {
117 res.f <<= (kDiySignificandSize - kDpSignificandSize - 2);
118 res.e = res.e - (kDiySignificandSize - kDpSignificandSize - 2);
122 void NormalizedBoundaries(DiyFp* minus, DiyFp* plus)
const {
123 DiyFp pl = DiyFp((f << 1) + 1, e - 1).NormalizeBoundary();
124 DiyFp mi = (f == kDpHiddenBit) ? DiyFp((f << 2) - 1, e - 2) : DiyFp((f << 1) - 1, e - 1);
125 mi.f <<= mi.e - pl.e;
131 double ToDouble()
const {
137 if (e < kDpDenormalExponent) {
141 if (e >= kDpMaxExponent) {
143 return std::numeric_limits<double>::infinity();
145 const uint64_t be = (e == kDpDenormalExponent && (f & kDpHiddenBit) == 0) ? 0 :
146 static_cast<uint64_t
>(e + kDpExponentBias);
147 u.u64 = (f & kDpSignificandMask) | (be << kDpSignificandSize);
151 static const int kDiySignificandSize = 64;
152 static const int kDpSignificandSize = 52;
153 static const int kDpExponentBias = 0x3FF + kDpSignificandSize;
154 static const int kDpMaxExponent = 0x7FF - kDpExponentBias;
155 static const int kDpMinExponent = -kDpExponentBias;
156 static const int kDpDenormalExponent = -kDpExponentBias + 1;
165 inline DiyFp GetCachedPowerByIndex(
size_t index) {
167 static const uint64_t kCachedPowers_F[] = {
213 static const int16_t kCachedPowers_E[] = {
214 -1220, -1193, -1166, -1140, -1113, -1087, -1060, -1034, -1007, -980,
215 -954, -927, -901, -874, -847, -821, -794, -768, -741, -715,
216 -688, -661, -635, -608, -582, -555, -529, -502, -475, -449,
217 -422, -396, -369, -343, -316, -289, -263, -236, -210, -183,
218 -157, -130, -103, -77, -50, -24, 3, 30, 56, 83,
219 109, 136, 162, 189, 216, 242, 269, 295, 322, 348,
220 375, 402, 428, 455, 481, 508, 534, 561, 588, 614,
221 641, 667, 694, 720, 747, 774, 800, 827, 853, 880,
222 907, 933, 960, 986, 1013, 1039, 1066
225 return DiyFp(kCachedPowers_F[index], kCachedPowers_E[index]);
228 inline DiyFp GetCachedPower(
int e,
int* K) {
231 double dk = (-61 - e) * 0.30102999566398114 + 347;
232 int k =
static_cast<int>(dk);
236 unsigned index =
static_cast<unsigned>((k >> 3) + 1);
237 *K = -(-348 +
static_cast<int>(index << 3));
239 return GetCachedPowerByIndex(index);
242 inline DiyFp GetCachedPower10(
int exp,
int *outExp) {
244 unsigned index =
static_cast<unsigned>(exp + 348) / 8u;
245 *outExp = -348 +
static_cast<int>(index) * 8;
246 return GetCachedPowerByIndex(index);
255 RAPIDJSON_DIAG_OFF(padded)
259 RAPIDJSON_NAMESPACE_END
261 #endif // RAPIDJSON_DIYFP_H_
#define RAPIDJSON_UINT64_C2(high32, low32)
Construct a 64-bit literal by a pair of 32-bit integer.
Definition: rapidjson.h:320
#define RAPIDJSON_ASSERT(x)
Assertion.
Definition: rapidjson.h:437