## root / src / common / bigdecimal.cc @ e1750c09

History | View | Annotate | Download (13.4 KB)

1 | 01873262 | Georg Kunz | ```
//==========================================================================
``` |
---|---|---|---|

2 | ```
// BIGDECIMAL.CC - part of
``` |
||

3 | ```
// OMNeT++/OMNEST
``` |
||

4 | ```
// Discrete System Simulation in C++
``` |
||

5 | ```
//
``` |
||

6 | ```
// Author: Tamas Borbely
``` |
||

7 | ```
//
``` |
||

8 | ```
//==========================================================================
``` |
||

9 | |||

10 | ```
/*--------------------------------------------------------------*
``` |
||

11 | ```
Copyright (C) 2006-2008 OpenSim Ltd.
``` |
||

12 | |||

13 | ```
This file is distributed WITHOUT ANY WARRANTY. See the file
``` |
||

14 | ```
`license' for details on this and other legal matters.
``` |
||

15 | ```
*--------------------------------------------------------------*/
``` |
||

16 | |||

17 | #include <sstream> |
||

18 | #include <assert.h> |
||

19 | #include <string.h> |
||

20 | #include "opp_ctype.h" |
||

21 | #include "platmisc.h" |
||

22 | #include "bigdecimal.h" |
||

23 | #include "commonutil.h" //NaN & friends |
||

24 | |||

25 | USING_NAMESPACE |
||

26 | |||

27 | ```
// helpers
``` |
||

28 | static inline int64 max(int64 x, int64 y) { return x > y ? x : y; } |
||

29 | static inline int64 min(int64 x, int64 y) { return x < y ? x : y; } |
||

30 | static inline int64 abs(int64 x) { return x >= 0 ? x : -x; } |
||

31 | static inline int sgn(int64 x) { return (x > 0 ? 1 : (x < 0 ? -1 : 0)); } |
||

32 | |||

33 | BigDecimal BigDecimal::Zero(0, 0); |
||

34 | BigDecimal BigDecimal::One(1, 0); |
||

35 | BigDecimal BigDecimal::MinusOne(-1, 0); |
||

36 | ```
BigDecimal BigDecimal::NaN(0, INT_MAX);
``` |
||

37 | ```
BigDecimal BigDecimal::PositiveInfinity(1, INT_MAX);
``` |
||

38 | ```
BigDecimal BigDecimal::NegativeInfinity(-1, INT_MAX);
``` |
||

39 | BigDecimal BigDecimal::Nil; |
||

40 | |||

41 | static int64 powersOfTen[21]; |
||

42 | static double negativePowersOfTen[21]; |
||

43 | |||

44 | class PowersOfTenInitializer |
||

45 | { |
||

46 | ```
public:
``` |
||

47 | PowersOfTenInitializer(); |
||

48 | }; |
||

49 | |||

50 | PowersOfTenInitializer initializer; |
||

51 | |||

52 | PowersOfTenInitializer::PowersOfTenInitializer() |
||

53 | { |
||

54 | ```
int64 power = 1;
``` |
||

55 | for (unsigned int i = 0; i < sizeof(powersOfTen) / sizeof(*powersOfTen); i++) { |
||

56 | powersOfTen[i] = power; |
||

57 | ```
power *= 10;
``` |
||

58 | } |
||

59 | |||

60 | double negativePower = 1; |
||

61 | for (unsigned int i = 0; i < sizeof(negativePowersOfTen) / sizeof(*negativePowersOfTen); i++) { |
||

62 | negativePowersOfTen[i] = negativePower; |
||

63 | negativePower /= 10.0; |
||

64 | } |
||

65 | } |
||

66 | |||

67 | ```
void BigDecimal::normalize()
``` |
||

68 | { |
||

69 | ```
// special values
``` |
||

70 | if (scale == INT_MAX && (intVal == -1 || intVal == 0 || intVal == 1)) |
||

71 | ```
return;
``` |
||

72 | |||

73 | ```
// zero
``` |
||

74 | if (intVal == 0) |
||

75 | { |
||

76 | ```
intVal = 0;
``` |
||

77 | ```
scale = 0;
``` |
||

78 | ```
return;
``` |
||

79 | } |
||

80 | |||

81 | ```
// underflow, XXX should throw an exception?
``` |
||

82 | ```
if (scale < minScale - INT64_MAX_DIGITS)
``` |
||

83 | { |
||

84 | ```
intVal = 0;
``` |
||

85 | ```
scale = 0;
``` |
||

86 | ```
return;
``` |
||

87 | } |
||

88 | |||

89 | ```
// overflow
``` |
||

90 | ```
if (scale > maxScale + INT64_MAX_DIGITS)
``` |
||

91 | throw opp_runtime_error("BigDecimal::normalize(): scale too big: %d.", scale); // XXX should be +-Infinity? |
||

92 | |||

93 | ```
// transform scale between minScale and maxScale
``` |
||

94 | ```
if (scale < minScale)
``` |
||

95 | { |
||

96 | ```
while (scale < minScale)
``` |
||

97 | { |
||

98 | ```
intVal /= 10;
``` |
||

99 | scale++; |
||

100 | |||

101 | if (intVal == 0) |
||

102 | { |
||

103 | ```
scale = 0;
``` |
||

104 | ```
return;
``` |
||

105 | } |
||

106 | } |
||

107 | } |
||

108 | else if (scale > maxScale) |
||

109 | { |
||

110 | ```
while (scale > maxScale)
``` |
||

111 | { |
||

112 | if (intVal > INT64_MAX/10) |
||

113 | throw opp_runtime_error("BigDecimal::normalize(): arithmetic overflow"); |
||

114 | |||

115 | ```
intVal *= 10;
``` |
||

116 | scale--; |
||

117 | } |
||

118 | } |
||

119 | |||

120 | ```
// strip trailing zeros
``` |
||

121 | while ((intVal % 10 == 0) && scale < maxScale) |
||

122 | { |
||

123 | ```
intVal /= 10;
``` |
||

124 | scale++; |
||

125 | } |
||

126 | } |
||

127 | |||

128 | int64 BigDecimal::getDigits(int scale, int numDigits) const |
||

129 | { |
||

130 | ```
assert(!this->isSpecial());
``` |
||

131 | |||

132 | int start = max(scale, this->scale); // inclusive |
||

133 | int end = scale+numDigits; // exclusive |
||

134 | |||

135 | ```
if (start >= end)
``` |
||

136 | return 0; |
||

137 | |||

138 | ```
int64 val = abs(this->intVal);
``` |
||

139 | for (int i = this->scale; i < start; ++i) |
||

140 | ```
val /= 10;
``` |
||

141 | |||

142 | if (val == 0) |
||

143 | return 0; |
||

144 | |||

145 | ```
int64 result = 0;
``` |
||

146 | ```
int digit;
``` |
||

147 | ```
int64 multiplier = 1;
``` |
||

148 | for (int i = start; i < end; ++i) |
||

149 | { |
||

150 | ```
digit = val % 10;
``` |
||

151 | ```
val /= 10;
``` |
||

152 | result += multiplier*digit; |
||

153 | ```
multiplier *= 10;
``` |
||

154 | } |
||

155 | |||

156 | for (int i = 0; i < (start-scale); ++i) |
||

157 | ```
result *= 10;
``` |
||

158 | |||

159 | ```
return result;
``` |
||

160 | } |
||

161 | |||

162 | const BigDecimal& BigDecimal::operator=(double d) |
||

163 | { |
||

164 | ```
// check NaN and infinity
``` |
||

165 | ```
if (::isNaN(d))
``` |
||

166 | return *this = NaN; |
||

167 | else if (::isPositiveInfinity(d)) |
||

168 | return *this = PositiveInfinity; |
||

169 | else if (::isNegativeInfinity(d)) |
||

170 | return *this = NegativeInfinity; |
||

171 | |||

172 | int sign = 1; |
||

173 | if (d < 0.0) |
||

174 | { |
||

175 | ```
sign = -1;
``` |
||

176 | d = -d; |
||

177 | } |
||

178 | |||

179 | ```
int exponent;
``` |
||

180 | double mantissa = frexp(d, &exponent); // d = mantissa * 2 ^ exponent |
||

181 | |||

182 | for (int i=0; i < 52; ++i) |
||

183 | { |
||

184 | mantissa *= 2.0; |
||

185 | exponent--; |
||

186 | } |
||

187 | |||

188 | ```
int64 intVal = (int64)mantissa; // d = intVal * 2 ^ exponent
``` |
||

189 | |||

190 | |||

191 | ```
// special case for 0.0, next loop would not terminate
``` |
||

192 | if (intVal == 0) |
||

193 | { |
||

194 | this->intVal = 0; |
||

195 | this->scale = 0; |
||

196 | return *this; |
||

197 | } |
||

198 | |||

199 | ```
// normalize
``` |
||

200 | while ((intVal & 1) == 0) |
||

201 | { |
||

202 | ```
intVal >>= 1;
``` |
||

203 | exponent++; |
||

204 | } |
||

205 | |||

206 | ```
//printf("intVal=%I64d, exponent=%d\n", intVal, exponent);
``` |
||

207 | |||

208 | ```
int scale;
``` |
||

209 | if (exponent < 0) |
||

210 | { |
||

211 | scale = exponent; |
||

212 | for (int i = exponent; i < 0; ++i) |
||

213 | { |
||

214 | if (intVal <= INT64_MAX/5) |
||

215 | { |
||

216 | ```
intVal *= 5;
``` |
||

217 | } |
||

218 | ```
else
``` |
||

219 | { |
||

220 | ```
intVal /= 2;
``` |
||

221 | scale++; |
||

222 | } |
||

223 | } |
||

224 | } |
||

225 | ```
else
``` |
||

226 | { |
||

227 | ```
scale = 0;
``` |
||

228 | for (int i = 0; i < exponent; ++i) |
||

229 | { |
||

230 | if (intVal <= INT64_MAX/2) |
||

231 | { |
||

232 | ```
intVal *= 2;
``` |
||

233 | } |
||

234 | ```
else
``` |
||

235 | { |
||

236 | ```
intVal /= 5;
``` |
||

237 | scale++; |
||

238 | } |
||

239 | } |
||

240 | } |
||

241 | |||

242 | ```
this->intVal = sign * intVal;
``` |
||

243 | ```
this->scale = scale;
``` |
||

244 | ```
this->normalize();
``` |
||

245 | return *this; |
||

246 | } |
||

247 | |||

248 | bool BigDecimal::operator<(const BigDecimal &x) const |
||

249 | { |
||

250 | ```
if (isSpecial() || x.isSpecial())
``` |
||

251 | { |
||

252 | ```
if (isNil() || x.isNil())
``` |
||

253 | throw opp_runtime_error("BigDecimal::operator<() received Nil."); |
||

254 | else if (isNaN() || x.isNaN()) |
||

255 | return false; |
||

256 | else if (x == PositiveInfinity) |
||

257 | return *this != PositiveInfinity; |
||

258 | else if (*this == NegativeInfinity) |
||

259 | ```
return x != NegativeInfinity;
``` |
||

260 | ```
else
``` |
||

261 | return false; |
||

262 | } |
||

263 | |||

264 | ```
if (scale == x.scale)
``` |
||

265 | ```
return intVal < x.intVal;
``` |
||

266 | ```
if (sgn(intVal) < sgn(x.intVal))
``` |
||

267 | return true; |
||

268 | ```
if (sgn(intVal) > sgn(x.intVal))
``` |
||

269 | return false; |
||

270 | |||

271 | assert((intVal < 0 && x.intVal < 0) || (intVal > 0 && x.intVal > 0)); |
||

272 | bool negatives = intVal < 0; |
||

273 | |||

274 | ```
// compare absolute values by comparing most significant digits first
``` |
||

275 | bool result = false; |
||

276 | for (int s = max(scale,x.scale); s > min(scale,x.scale)-18; s-=18) |
||

277 | { |
||

278 | int64 digits = this->getDigits(s, 18); |
||

279 | ```
int64 digitsX = x.getDigits(s, 18);
``` |
||

280 | ```
if (digits < digitsX)
``` |
||

281 | { |
||

282 | ```
result = true;
``` |
||

283 | ```
break;
``` |
||

284 | } |
||

285 | else if (digits > digitsX) |
||

286 | { |
||

287 | ```
result = false;
``` |
||

288 | ```
break;
``` |
||

289 | } |
||

290 | } |
||

291 | |||

292 | ```
return negatives ? !result : result;
``` |
||

293 | } |
||

294 | |||

295 | int64 BigDecimal::getMantissaForScale(int reqScale) const |
||

296 | { |
||

297 | ```
if (isSpecial())
``` |
||

298 | throw opp_runtime_error("BigDecimal: cannot return mantissa for Nil, NaN or +/-Inf value"); |
||

299 | checkScale(reqScale); |
||

300 | ```
int scaleDiff = scale - reqScale;
``` |
||

301 | if (scaleDiff == 0) |
||

302 | ```
return intVal;
``` |
||

303 | else if (scaleDiff < 0) |
||

304 | ```
return intVal / powersOfTen[-scaleDiff];
``` |
||

305 | ```
else
``` |
||

306 | ```
return intVal * powersOfTen[scaleDiff];
``` |
||

307 | } |
||

308 | |||

309 | double BigDecimal::dbl() const |
||

310 | { |
||

311 | ```
if (isSpecial())
``` |
||

312 | { |
||

313 | ```
if (isNaN())
``` |
||

314 | ```
return ::NaN;
``` |
||

315 | else if (*this == PositiveInfinity) |
||

316 | ```
return POSITIVE_INFINITY;
``` |
||

317 | else if (*this == NegativeInfinity) |
||

318 | ```
return NEGATIVE_INFINITY;
``` |
||

319 | else // Nil |
||

320 | throw opp_runtime_error("BigDecimal::dbl(): received Nil."); // XXX should return NaN? |
||

321 | } |
||

322 | |||

323 | return (double)intVal * negativePowersOfTen[-scale]; |
||

324 | } |
||

325 | |||

326 | std::string BigDecimal::str() const |
||

327 | { |
||

328 | ```
// delegate to operator<<
``` |
||

329 | std::stringstream out; |
||

330 | ```
out << *this;
``` |
||

331 | ```
return out.str();
``` |
||

332 | } |
||

333 | |||

334 | char *BigDecimal::ttoa(char *buf, const BigDecimal &x, char *&endp) |
||

335 | { |
||

336 | ```
// special values
``` |
||

337 | ```
if (x.isSpecial())
``` |
||

338 | { |
||

339 | ```
if (x.isNaN())
``` |
||

340 | { |
||

341 | ```
strcpy(buf, "NaN");
``` |
||

342 | ```
endp = buf+3;
``` |
||

343 | } |
||

344 | else if (x == PositiveInfinity) |
||

345 | { |
||

346 | ```
strcpy(buf, "+Inf");
``` |
||

347 | ```
endp = buf+4;
``` |
||

348 | } |
||

349 | else if (x == NegativeInfinity) |
||

350 | { |
||

351 | ```
strcpy(buf, "-Inf");
``` |
||

352 | ```
endp = buf+4;
``` |
||

353 | } |
||

354 | else // Nil |
||

355 | throw opp_runtime_error("BigDecimal::ttoa(): received Nil."); |
||

356 | ```
return buf;
``` |
||

357 | } |
||

358 | |||

359 | int64 intVal = x.getIntValue(); |
||

360 | ```
int scale = x.getScale();
``` |
||

361 | |||

362 | ```
// prepare for conversion
``` |
||

363 | endp = buf+63; //19+19+5 should be enough, but play it safe |
||

364 | ```
*endp = '\0';
``` |
||

365 | ```
char *s = endp;
``` |
||

366 | if (intVal==0) |
||

367 | {*--s = '0'; return s;} |
||

368 | |||

369 | ```
// convert digits
``` |
||

370 | bool negative = intVal<0; |
||

371 | ```
if (negative) intVal = -intVal;
``` |
||

372 | |||

373 | bool skipzeros = true; |
||

374 | ```
int decimalplace = scale;
``` |
||

375 | ```
do {
``` |
||

376 | ```
int64 res = intVal / 10;
``` |
||

377 | int digit = intVal - (10*res); |
||

378 | |||

379 | if (skipzeros && (digit!=0 || decimalplace>=0)) |
||

380 | ```
skipzeros = false;
``` |
||

381 | if (decimalplace++==0 && s!=endp) |
||

382 | ```
*--s = '.';
``` |
||

383 | ```
if (!skipzeros)
``` |
||

384 | ```
*--s = '0'+digit;
``` |
||

385 | intVal = res; |
||

386 | ```
} while (intVal);
``` |
||

387 | |||

388 | ```
// add leading zeros, decimal point, etc if needed
``` |
||

389 | if (decimalplace<=0) |
||

390 | { |
||

391 | while (decimalplace++ < 0) |
||

392 | ```
*--s = '0';
``` |
||

393 | ```
*--s = '.';
``` |
||

394 | ```
*--s = '0';
``` |
||

395 | } |
||

396 | |||

397 | if (negative) *--s = '-'; |
||

398 | ```
return s;
``` |
||

399 | } |
||

400 | |||

401 | const BigDecimal BigDecimal::parse(const char *s) |
||

402 | { |
||

403 | const char *endp; |
||

404 | ```
return parse(s, endp);
``` |
||

405 | } |
||

406 | |||

407 | #define OVERFLOW_CHECK(c,s) if (!(c)) throw opp_runtime_error("BigDecimal::parse(\"%s\"): arithmetic overflow", (s)); |
||

408 | |||

409 | |||

410 | const BigDecimal BigDecimal::parse(const char *s, const char *&endp) |
||

411 | { |
||

412 | ```
int64 intVal = 0;
``` |
||

413 | ```
int digit;
``` |
||

414 | int digits = 0; |
||

415 | int scale = 0; |
||

416 | int sign = 1; |
||

417 | const char *p = s; |
||

418 | |||

419 | ```
// skip leading spaces
``` |
||

420 | ```
while (opp_isspace(*p))
``` |
||

421 | ++p; |
||

422 | |||

423 | ```
// optional signs
``` |
||

424 | if (*p == '-') |
||

425 | { |
||

426 | ```
sign = -1;
``` |
||

427 | ++p; |
||

428 | } |
||

429 | else if (*p == '+') |
||

430 | ++p; |
||

431 | |||

432 | ```
// parse special numbers
``` |
||

433 | ```
if (opp_isalpha(*p))
``` |
||

434 | { |
||

435 | if (strncasecmp(p, "nan", 3) == 0) |
||

436 | { |
||

437 | ```
endp = p+3;
``` |
||

438 | ```
return NaN;
``` |
||

439 | } |
||

440 | else if (strncasecmp(p, "inf", 3) == 0) // inf and infinity |
||

441 | { |
||

442 | ```
endp = p+3;
``` |
||

443 | if (strncasecmp(endp, "inity", 5) == 0) |
||

444 | ```
endp += 5;
``` |
||

445 | return sign > 0 ? PositiveInfinity : NegativeInfinity; |
||

446 | } |
||

447 | ```
else
``` |
||

448 | { |
||

449 | endp = p; |
||

450 | return Zero; // XXX should return Nil? |
||

451 | } |
||

452 | } |
||

453 | else if (*p=='1' && *(p+1)=='.' && *(p+2)=='#') |
||

454 | { |
||

455 | if (strncasecmp(p+3, "ind", 3) == 0) |
||

456 | { |
||

457 | ```
endp = p+6;
``` |
||

458 | ```
return NaN;
``` |
||

459 | } |
||

460 | else if (strncasecmp(p+3, "inf", 6) == 0) |
||

461 | { |
||

462 | ```
endp = p+6;
``` |
||

463 | return sign > 0 ? PositiveInfinity : NegativeInfinity; |
||

464 | } |
||

465 | } |
||

466 | |||

467 | ```
// digits before decimal
``` |
||

468 | ```
while (opp_isdigit(*p))
``` |
||

469 | { |
||

470 | ```
OVERFLOW_CHECK(intVal <= INT64_MAX / 10, s);
``` |
||

471 | ```
intVal *= 10;
``` |
||

472 | ```
digit = ((*p++)-'0');
``` |
||

473 | OVERFLOW_CHECK(intVal <= INT64_MAX - digit, s); |
||

474 | intVal += digit; |
||

475 | digits++; |
||

476 | } |
||

477 | |||

478 | if (digits == 0 && *p != '.') |
||

479 | { |
||

480 | throw opp_runtime_error("BigDecimal::parse(\"%s\"): missing digits", s); |
||

481 | } |
||

482 | |||

483 | ```
// digits after decimal
``` |
||

484 | if (*p == '.') |
||

485 | { |
||

486 | p++; |
||

487 | ```
while (opp_isdigit(*p))
``` |
||

488 | { |
||

489 | ```
OVERFLOW_CHECK(intVal <= INT64_MAX / 10, s);
``` |
||

490 | ```
intVal *= 10;
``` |
||

491 | ```
digit = ((*p++)-'0');
``` |
||

492 | OVERFLOW_CHECK(intVal <= INT64_MAX - digit, s); |
||

493 | intVal += digit; |
||

494 | digits++; |
||

495 | scale--; |
||

496 | } |
||

497 | } |
||

498 | |||

499 | endp = p; |
||

500 | ```
return BigDecimal(sign*intVal, scale);
``` |
||

501 | } |
||

502 | |||

503 | const BigDecimal operator+(const BigDecimal& x, const BigDecimal& y) |
||

504 | { |
||

505 | ```
// 1. try to add exactly
``` |
||

506 | ```
int scale = std::min(x.scale, y.scale);
``` |
||

507 | ```
int xm = x.scale - scale;
``` |
||

508 | ```
int ym = y.scale - scale;
``` |
||

509 | |||

510 | const int NUMPOWERS = sizeof(powersOfTen) / sizeof(*powersOfTen); |
||

511 | |||

512 | if (!x.isSpecial() && !y.isSpecial() && 0 <= xm && xm < NUMPOWERS && 0 <= ym && ym < NUMPOWERS) |
||

513 | { |
||

514 | int64 xmp = powersOfTen[xm]; |
||

515 | int64 xv = x.intVal * xmp; |
||

516 | |||

517 | ```
if (xv / xmp == x.intVal) {
``` |
||

518 | int64 ymp = powersOfTen[ym]; |
||

519 | int64 yv = y.intVal * ymp; |
||

520 | |||

521 | ```
if (yv / ymp == y.intVal) {
``` |
||

522 | ```
bool sameSign = haveSameSign(xv, yv);
``` |
||

523 | int64 intVal = xv + yv; |
||

524 | |||

525 | ```
if (!sameSign || haveSameSign(intVal, yv))
``` |
||

526 | ```
return BigDecimal(intVal, scale);
``` |
||

527 | } |
||

528 | } |
||

529 | } |
||

530 | |||

531 | ```
// 2. add with precision loss
``` |
||

532 | ```
return BigDecimal(x.dbl()+y.dbl());
``` |
||

533 | } |
||

534 | |||

535 | const BigDecimal operator-(const BigDecimal& x, const BigDecimal& y) |
||

536 | { |
||

537 | ```
// 1. try to subtract exactly
``` |
||

538 | ```
int scale = std::min(x.scale, y.scale);
``` |
||

539 | ```
int xm = x.scale - scale;
``` |
||

540 | ```
int ym = y.scale - scale;
``` |
||

541 | |||

542 | const int NUMPOWERS = sizeof(powersOfTen) / sizeof(*powersOfTen); |
||

543 | |||

544 | if (!x.isSpecial() && !y.isSpecial() && 0 <= xm && xm < NUMPOWERS && 0 <= ym && ym < NUMPOWERS) |
||

545 | { |
||

546 | int64 xmp = powersOfTen[xm]; |
||

547 | int64 xv = x.intVal * xmp; |
||

548 | |||

549 | ```
if (xv / xmp == x.intVal) {
``` |
||

550 | int64 ymp = powersOfTen[ym]; |
||

551 | int64 yv = y.intVal * ymp; |
||

552 | |||

553 | ```
if (yv / ymp == y.intVal) {
``` |
||

554 | ```
bool differentSign = !haveSameSign(xv, yv);
``` |
||

555 | int64 intVal = xv - yv; |
||

556 | |||

557 | ```
if (!differentSign || !haveSameSign(intVal, yv))
``` |
||

558 | ```
return BigDecimal(intVal, scale);
``` |
||

559 | } |
||

560 | } |
||

561 | } |
||

562 | |||

563 | ```
// 2. subtract with precision loss
``` |
||

564 | ```
return BigDecimal(x.dbl()-y.dbl());
``` |
||

565 | } |