src/System.Private.CoreLib/shared/System/Number.Grisu3.cs

   1 // Licensed to the .NET Foundation under one or more agreements.
   2 // The .NET Foundation licenses this file to you under the MIT license.
   3 // See the LICENSE file in the project root for more information.
   4
   5 using System.Diagnostics;
   6
   7 namespace System
   8 {
   9     internal static partial class Number
  10     {
  11         internal static unsafe class Grisu3
  12         {
  13             private const int Alpha = -59;
  14             private const double D1Log210 = 0.301029995663981195;
  15             private const int Gamma = -32;
  16             private const int PowerDecimalExponentDistance = 8;
  17             private const int PowerMinDecimalExponent = -348;
  18             private const int PowerMaxDecimalExponent = 340;
  19             private const int PowerOffset = -PowerMinDecimalExponent;
  20             private const uint Ten4 = 10000;
  21             private const uint Ten5 = 100000;
  22             private const uint Ten6 = 1000000;
  23             private const uint Ten7 = 10000000;
  24             private const uint Ten8 = 100000000;
  25             private const uint Ten9 = 1000000000;
  26
  27             private static readonly short[] s_CachedPowerBinaryExponents = new short[]
  28             {
  29                 -1220,
  30                 -1193,
  31                 -1166,
  32                 -1140,
  33                 -1113,
  34                 -1087,
  35                 -1060,
  36                 -1034,
  37                 -1007,
  38                 -980,
  39                 -954,
  40                 -927,
  41                 -901,
  42                 -874,
  43                 -847,
  44                 -821,
  45                 -794,
  46                 -768,
  47                 -741,
  48                 -715,
  49                 -688,
  50                 -661,
  51                 -635,
  52                 -608,
  53                 -582,
  54                 -555,
  55                 -529,
  56                 -502,
  57                 -475,
  58                 -449,
  59                 -422,
  60                 -396,
  61                 -369,
  62                 -343,
  63                 -316,
  64                 -289,
  65                 -263,
  66                 -236,
  67                 -210,
  68                 -183,
  69                 -157,
  70                 -130,
  71                 -103,
  72                 -77,
  73                 -50,
  74                 -24,
  75                 3,
  76                 30,
  77                 56,
  78                 83,
  79                 109,
  80                 136,
  81                 162,
  82                 189,
  83                 216,
  84                 242,
  85                 269,
  86                 295,
  87                 322,
  88                 348,
  89                 375,
  90                 402,
  91                 428,
  92                 455,
  93                 481,
  94                 508,
  95                 534,
  96                 561,
  97                 588,
  98                 614,
  99                 641,
 100                 667,
 101                 694,
 102                 720,
 103                 747,
 104                 774,
 105                 800,
 106                 827,
 107                 853,
 108                 880,
 109                 907,
 110                 933,
 111                 960,
 112                 986,
 113                 1013,
 114                 1039,
 115                 1066,
 116             };
 117
 118             private static readonly short[] s_CachedPowerDecimalExponents = new short[]
 119             {
 120                 PowerMinDecimalExponent,
 121                 -340,
 122                 -332,
 123                 -324,
 124                 -316,
 125                 -308,
 126                 -300,
 127                 -292,
 128                 -284,
 129                 -276,
 130                 -268,
 131                 -260,
 132                 -252,
 133                 -244,
 134                 -236,
 135                 -228,
 136                 -220,
 137                 -212,
 138                 -204,
 139                 -196,
 140                 -188,
 141                 -180,
 142                 -172,
 143                 -164,
 144                 -156,
 145                 -148,
 146                 -140,
 147                 -132,
 148                 -124,
 149                 -116,
 150                 -108,
 151                 -100,
 152                 -92,
 153                 -84,
 154                 -76,
 155                 -68,
 156                 -60,
 157                 -52,
 158                 -44,
 159                 -36,
 160                 -28,
 161                 -20,
 162                 -12,
 163                 -4,
 164                 4,
 165                 12,
 166                 20,
 167                 28,
 168                 36,
 169                 44,
 170                 52,
 171                 60,
 172                 68,
 173                 76,
 174                 84,
 175                 92,
 176                 100,
 177                 108,
 178                 116,
 179                 124,
 180                 132,
 181                 140,
 182                 148,
 183                 156,
 184                 164,
 185                 172,
 186                 180,
 187                 188,
 188                 196,
 189                 204,
 190                 212,
 191                 220,
 192                 228,
 193                 236,
 194                 244,
 195                 252,
 196                 260,
 197                 268,
 198                 276,
 199                 284,
 200                 292,
 201                 300,
 202                 308,
 203                 316,
 204                 324,
 205                 332,
 206                 PowerMaxDecimalExponent,
 207             };
 208
 209             private static readonly uint[] s_CachedPowerOfTen = new uint[]
 210             {
 211                 1,          // 10^0
 212                 10,         // 10^1
 213                 100,        // 10^2
 214                 1000,       // 10^3
 215                 10000,      // 10^4
 216                 100000,     // 10^5
 217                 1000000,    // 10^6
 218                 10000000,   // 10^7
 219                 100000000,  // 10^8
 220                 1000000000, // 10^9
 221             };
 222
 223             private static readonly ulong[] s_CachedPowerSignificands = new ulong[]
 224             {
 225                 0xFA8FD5A0081C0288,
 226                 0xBAAEE17FA23EBF76,
 227                 0x8B16FB203055AC76,
 228                 0xCF42894A5DCE35EA,
 229                 0x9A6BB0AA55653B2D,
 230                 0xE61ACF033D1A45DF,
 231                 0xAB70FE17C79AC6CA,
 232                 0xFF77B1FCBEBCDC4F,
 233                 0xBE5691EF416BD60C,
 234                 0x8DD01FAD907FFC3C,
 235                 0xD3515C2831559A83,
 236                 0x9D71AC8FADA6C9B5,
 237                 0xEA9C227723EE8BCB,
 238                 0xAECC49914078536D,
 239                 0x823C12795DB6CE57,
 240                 0xC21094364DFB5637,
 241                 0x9096EA6F3848984F,
 242                 0xD77485CB25823AC7,
 243                 0xA086CFCD97BF97F4,
 244                 0xEF340A98172AACE5,
 245                 0xB23867FB2A35B28E,
 246                 0x84C8D4DFD2C63F3B,
 247                 0xC5DD44271AD3CDBA,
 248                 0x936B9FCEBB25C996,
 249                 0xDBAC6C247D62A584,
 250                 0xA3AB66580D5FDAF6,
 251                 0xF3E2F893DEC3F126,
 252                 0xB5B5ADA8AAFF80B8,
 253                 0x87625F056C7C4A8B,
 254                 0xC9BCFF6034C13053,
 255                 0x964E858C91BA2655,
 256                 0xDFF9772470297EBD,
 257                 0xA6DFBD9FB8E5B88F,
 258                 0xF8A95FCF88747D94,
 259                 0xB94470938FA89BCF,
 260                 0x8A08F0F8BF0F156B,
 261                 0xCDB02555653131B6,
 262                 0x993FE2C6D07B7FAC,
 263                 0xE45C10C42A2B3B06,
 264                 0xAA242499697392D3,
 265                 0xFD87B5F28300CA0E,
 266                 0xBCE5086492111AEB,
 267                 0x8CBCCC096F5088CC,
 268                 0xD1B71758E219652C,
 269                 0x9C40000000000000,
 270                 0xE8D4A51000000000,
 271                 0xAD78EBC5AC620000,
 272                 0x813F3978F8940984,
 273                 0xC097CE7BC90715B3,
 274                 0x8F7E32CE7BEA5C70,
 275                 0xD5D238A4ABE98068,
 276                 0x9F4F2726179A2245,
 277                 0xED63A231D4C4FB27,
 278                 0xB0DE65388CC8ADA8,
 279                 0x83C7088E1AAB65DB,
 280                 0xC45D1DF942711D9A,
 281                 0x924D692CA61BE758,
 282                 0xDA01EE641A708DEA,
 283                 0xA26DA3999AEF774A,
 284                 0xF209787BB47D6B85,
 285                 0xB454E4A179DD1877,
 286                 0x865B86925B9BC5C2,
 287                 0xC83553C5C8965D3D,
 288                 0x952AB45CFA97A0B3,
 289                 0xDE469FBD99A05FE3,
 290                 0xA59BC234DB398C25,
 291                 0xF6C69A72A3989F5C,
 292                 0xB7DCBF5354E9BECE,
 293                 0x88FCF317F22241E2,
 294                 0xCC20CE9BD35C78A5,
 295                 0x98165AF37B2153DF,
 296                 0xE2A0B5DC971F303A,
 297                 0xA8D9D1535CE3B396,
 298                 0xFB9B7CD9A4A7443C,
 299                 0xBB764C4CA7A44410,
 300                 0x8BAB8EEFB6409C1A,
 301                 0xD01FEF10A657842C,
 302                 0x9B10A4E5E9913129,
 303                 0xE7109BFBA19C0C9D,
 304                 0xAC2820D9623BF429,
 305                 0x80444B5E7AA7CF85,
 306                 0xBF21E44003ACDD2D,
 307                 0x8E679C2F5E44FF8F,
 308                 0xD433179D9C8CB841,
 309                 0x9E19DB92B4E31BA9,
 310                 0xEB96BF6EBADF77D9,
 311                 0xAF87023B9BF0EE6B,
 312             };
 313
 314             public static bool Run(double value, int precision, ref NumberBuffer number)
 315             {
 316                 // ========================================================================================================================================
 317                 // This implementation is based on the paper: http://www.cs.tufts.edu/~nr/cs257/archive/florian-loitsch/printf.pdf
 318                 // You must read this paper to fully understand the code.
 319                 //
 320                 // Deviation: Instead of generating shortest digits, we generate the digits according to the input count.
 321                 // Therefore, we do not need m+ and m- which are used to determine the exact range of values.
 322                 // ========================================================================================================================================
 323                 //
 324                 // Overview:
 325                 //
 326                 // The idea of Grisu3 is to leverage additional bits and cached power of ten to produce the digits.
 327                 // We need to create a handmade floating point data structure DiyFp to extend the bits of double.
 328                 // We also need to cache the powers of ten for digits generation. By choosing the correct index of powers
 329                 // we need to start with, we can eliminate the expensive big num divide operation.
 330                 //
 331                 // Grisu3 is imprecision for some numbers. Fortunately, the algorithm itself can determine that and give us
 332                 // a success/fail flag. We may fall back to other algorithms (For instance, Dragon4) if it fails.
 333                 //
 334                 // w: the normalized DiyFp from the input value.
 335                 // mk: The index of the cached powers.
 336                 // cmk: The cached power.
 337                 // D: Product: w * cmk.
 338                 // kappa: A factor used for generating digits. See step 5 of the Grisu3 procedure in the paper.
 339
 340                 // Handle sign bit.
 341                 if (double.IsNegative(value))
 342                 {
 343                     value = -value;
 344                     number.sign = true;
 345                 }
 346                 else
 347                 {
 348                     number.sign = false;
 349                 }
 350
 351                 // Step 1: Determine the normalized DiyFp w.
 352
 353                 DiyFp.GenerateNormalizedDiyFp(value, out DiyFp w);
 354
 355                 // Step 2: Find the cached power of ten.
 356
 357                 // Compute the proper index mk.
 358                 int mk = KComp(w.e + DiyFp.SignificandLength);
 359
 360                 // Retrieve the cached power of ten.
 361                 CachedPower(mk, out DiyFp cmk, out int decimalExponent);
 362
 363                 // Step 3: Scale the w with the cached power of ten.
 364
 365                 DiyFp.Multiply(ref w, ref cmk, out DiyFp D);
 366
 367                 // Step 4: Generate digits.
 368
 369                 bool isSuccess = DigitGen(ref D, precision, number.GetDigitsPointer(), out int length, out int kappa);
 370
 371                 if (isSuccess)
 372                 {
 373                     number.digits[precision] = '\0';
 374                     number.scale = (length - decimalExponent + kappa);
 375                 }
 376
 377                 return isSuccess;
 378             }
 379
 380             // Returns the biggest power of ten that is less than or equal to the given number.
 381             static void BiggestPowerTenLessThanOrEqualTo(uint number, int bits, out uint power, out int exponent)
 382             {
 383                 switch (bits)
 384                 {
 385                     case 32:
 386                     case 31:
 387                     case 30:
 388                     {
 389                         if (Ten9 <= number)
 390                         {
 391                             power = Ten9;
 392                             exponent = 9;
 393                             break;
 394                         }
 395
 396                         goto case 29;
 397                     }
 398
 399                     case 29:
 400                     case 28:
 401                     case 27:
 402                     {
 403                         if (Ten8 <= number)
 404                         {
 405                             power = Ten8;
 406                             exponent = 8;
 407                             break;
 408                         }
 409
 410                         goto case 26;
 411                     }
 412
 413                     case 26:
 414                     case 25:
 415                     case 24:
 416                     {
 417                         if (Ten7 <= number)
 418                         {
 419                             power = Ten7;
 420                             exponent = 7;
 421                             break;
 422                         }
 423
 424                         goto case 23;
 425                     }
 426
 427                     case 23:
 428                     case 22:
 429                     case 21:
 430                     case 20:
 431                     {
 432                         if (Ten6 <= number)
 433                         {
 434                             power = Ten6;
 435                             exponent = 6;
 436                             break;
 437                         }
 438
 439                         goto case 19;
 440                     }
 441
 442                     case 19:
 443                     case 18:
 444                     case 17:
 445                     {
 446                         if (Ten5 <= number)
 447                         {
 448                             power = Ten5;
 449                             exponent = 5;
 450                             break;
 451                         }
 452
 453                         goto case 16;
 454                     }
 455
 456                     case 16:
 457                     case 15:
 458                     case 14:
 459                     {
 460                         if (Ten4 <= number)
 461                         {
 462                             power = Ten4;
 463                             exponent = 4;
 464                             break;
 465                         }
 466
 467                         goto case 13;
 468                     }
 469
 470                     case 13:
 471                     case 12:
 472                     case 11:
 473                     case 10:
 474                     {
 475                         if (1000 <= number)
 476                         {
 477                             power = 1000;
 478                             exponent = 3;
 479                             break;
 480                         }
 481
 482                         goto case 9;
 483                     }
 484
 485                     case 9:
 486                     case 8:
 487                     case 7:
 488                     {
 489                         if (100 <= number)
 490                         {
 491                             power = 100;
 492                             exponent = 2;
 493                             break;
 494                         }
 495
 496                         goto case 6;
 497                     }
 498
 499                     case 6:
 500                     case 5:
 501                     case 4:
 502                     {
 503                         if (10 <= number)
 504                         {
 505                             power = 10;
 506                             exponent = 1;
 507                             break;
 508                         }
 509
 510                         goto case 3;
 511                     }
 512
 513                     case 3:
 514                     case 2:
 515                     case 1:
 516                     {
 517                         if (1 <= number)
 518                         {
 519                             power = 1;
 520                             exponent = 0;
 521                             break;
 522                         }
 523
 524                         goto case 0;
 525                     }
 526
 527                     case 0:
 528                     {
 529                         power = 0;
 530                         exponent = -1;
 531                         break;
 532                     }
 533
 534                     default:
 535                     {
 536                         power = 0;
 537                         exponent = 0;
 538
 539                         Debug.Fail("unreachable");
 540                         break;
 541                     }
 542                 }
 543             }
 544
 545             private static void CachedPower(int k, out DiyFp cmk, out int decimalExponent)
 546             {
 547                 int index = ((PowerOffset + k - 1) / PowerDecimalExponentDistance) + 1;
 548                 cmk = new DiyFp(s_CachedPowerSignificands[index], s_CachedPowerBinaryExponents[index]);
 549                 decimalExponent = s_CachedPowerDecimalExponents[index];
 550             }
 551
 552             private static bool DigitGen(ref DiyFp mp, int precision, char* digits, out int length, out int k)
 553             {
 554                 // Split the input mp to two parts. Part 1 is integral. Part 2 can be used to calculate
 555                 // fractional.
 556                 //
 557                 // mp: the input DiyFp scaled by cached power.
 558                 // K: final kappa.
 559                 // p1: part 1.
 560                 // p2: part 2.
 561
 562                 Debug.Assert(precision > 0);
 563                 Debug.Assert(digits != null);
 564                 Debug.Assert(mp.e >= Alpha);
 565                 Debug.Assert(mp.e <= Gamma);
 566
 567                 ulong mpF = mp.f;
 568                 int mpE = mp.e;
 569
 570                 var one = new DiyFp(1UL << -mpE, mpE);
 571
 572                 ulong oneF = one.f;
 573                 int oneNegE = -one.e;
 574
 575                 ulong ulp = 1;
 576
 577                 uint p1 = (uint)(mpF >> oneNegE);
 578                 ulong p2 = mpF & (oneF - 1);
 579
 580                 // When p2 (fractional part) is zero, we can predicate if p1 is good to produce the numbers in requested digit count:
 581                 //
 582                 // - When requested digit count >= 11, p1 is not be able to exhaust the count as 10^(11 - 1) > uint.MaxValue >= p1.
 583                 // - When p1 < 10^(count - 1), p1 is not be able to exhaust the count.
 584                 // - Otherwise, p1 may have chance to exhaust the count.
 585                 if ((p2 == 0) && ((precision >= 11) || (p1 < s_CachedPowerOfTen[precision - 1])))
 586                 {
 587                     length = 0;
 588                     k = 0;
 589                     return false;
 590                 }
 591
 592                 // Note: The code in the paper simply assigns div to Ten9 and kappa to 10 directly.
 593                 // That means we need to check if any leading zero of the generated
 594                 // digits during the while loop, which hurts the performance.
 595                 //
 596                 // Now if we can estimate what the div and kappa, we do not need to check the leading zeros.
 597                 // The idea is to find the biggest power of 10 that is less than or equal to the given number.
 598                 // Then we don't need to worry about the leading zeros and we can get 10% performance gain.
 599                 int index = 0;
 600                 BiggestPowerTenLessThanOrEqualTo(p1, (DiyFp.SignificandLength - oneNegE), out uint div, out int kappa);
 601                 kappa++;
 602
 603                 // Produce integral.
 604                 while (kappa > 0)
 605                 {
 606                     int d = (int)(Math.DivRem(p1, div, out p1));
 607                     digits[index] = (char)('0' + d);
 608
 609                     index++;
 610                     precision--;
 611                     kappa--;
 612
 613                     if (precision == 0)
 614                     {
 615                         break;
 616                     }
 617
 618                     div /= 10;
 619                 }
 620
 621                 // End up here if we already exhausted the digit count.
 622                 if (precision == 0)
 623                 {
 624                     ulong rest = ((ulong)(p1) << oneNegE) + p2;
 625
 626                     length = index;
 627                     k = kappa;
 628
 629                     return RoundWeed(
 630                         digits,
 631                         index,
 632                         rest,
 633                         ((ulong)(div)) << oneNegE,
 634                         ulp,
 635                         ref k
 636                     );
 637                 }
 638
 639                 // We have to generate digits from part2 if we have requested digit count left
 640                 // and part2 is greater than ulp.
 641                 while ((precision > 0) && (p2 > ulp))
 642                 {
 643                     p2 *= 10;
 644
 645                     int d = (int)(p2 >> oneNegE);
 646                     digits[index] = (char)('0' + d);
 647
 648                     index++;
 649                     precision--;
 650                     kappa--;
 651
 652                     p2 &= (oneF - 1);
 653
 654                     ulp *= 10;
 655                 }
 656
 657                 // If we haven't exhausted the requested digit counts, the Grisu3 algorithm fails.
 658                 if (precision != 0)
 659                 {
 660                     length = 0;
 661                     k = 0;
 662                     return false;
 663                 }
 664
 665                 length = index;
 666                 k = kappa;
 667
 668                 return RoundWeed(digits, index, p2, oneF, ulp, ref k);
 669             }
 670
 671             private static int KComp(int e)
 672             {
 673                 return (int)(Math.Ceiling((Alpha - e + DiyFp.SignificandLength - 1) * D1Log210));
 674             }
 675
 676             private static bool RoundWeed(char* buffer, int len, ulong rest, ulong tenKappa, ulong ulp, ref int kappa)
 677             {
 678                 Debug.Assert(rest < tenKappa);
 679
 680                 // 1. tenKappa <= ulp: we don't have an idea which way to round.
 681                 // 2. Even if tenKappa > ulp, but if tenKappa <= 2 * ulp we cannot find the way to round.
 682                 // Note: to prevent overflow, we need to use tenKappa - ulp <= ulp.
 683                 if ((tenKappa <= ulp) || ((tenKappa - ulp) <= ulp))
 684                 {
 685                     return false;
 686                 }
 687
 688                 // tenKappa >= 2 * (rest + ulp). We should round down.
 689                 // Note: to prevent overflow, we need to check if tenKappa > 2 * rest as a prerequisite.
 690                 if (((tenKappa - rest) > rest) && ((tenKappa - (2 * rest)) >= (2 * ulp)))
 691                 {
 692                     return true;
 693                 }
 694
 695                 // tenKappa <= 2 * (rest - ulp). We should round up.
 696                 // Note: to prevent overflow, we need to check if rest > ulp as a prerequisite.
 697                 if ((rest > ulp) && (tenKappa <= (rest - ulp) || ((tenKappa - (rest - ulp)) <= (rest - ulp))))
 698                 {
 699                     // Find all 9s from end to start.
 700                     buffer[len - 1]++;
 701                     for (int i = len - 1; i > 0; i--)
 702                     {
 703                         if (buffer[i] != (char)('0' + 10))
 704                         {
 705                             // We end up a number less than 9.
 706                             break;
 707                         }
 708
 709                         // Current number becomes 0 and add the promotion to the next number.
 710                         buffer[i] = '0';
 711                         buffer[i - 1]++;
 712                     }
 713
 714                     if (buffer[0] == (char)('0' + 10))
 715                     {
 716                         // First number is '0' + 10 means all numbers are 9.
 717                         // We simply make the first number to 1 and increase the kappa.
 718                         buffer[0] = '1';
 719                         kappa++;
 720                     }
 721
 722                     return true;
 723                 }
 724
 725                 return false;
 726             }
 727         }
 728     }
 729 }