contrib/mpfr/set_d64.c

   1 /* mpfr_set_decimal64 -- convert a IEEE 754r decimal64 float to
   2                          a multiple precision floating-point number
   3
   4 See http://gcc.gnu.org/ml/gcc/2006-06/msg00691.html,
   5 http://gcc.gnu.org/onlinedocs/gcc/Decimal-Float.html,
   6 and TR 24732 <http://www.open-std.org/jtc1/sc22/wg14/www/projects#24732>.
   7
   8 Copyright 2006, 2007, 2008, 2009 Free Software Foundation, Inc.
   9 Contributed by the Arenaire and Cacao projects, INRIA.
  10
  11 This file is part of the GNU MPFR Library.
  12
  13 The GNU MPFR Library is free software; you can redistribute it and/or modify
  14 it under the terms of the GNU Lesser General Public License as published by
  15 the Free Software Foundation; either version 2.1 of the License, or (at your
  16 option) any later version.
  17
  18 The GNU MPFR Library is distributed in the hope that it will be useful, but
  19 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
  20 or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
  21 License for more details.
  22
  23 You should have received a copy of the GNU Lesser General Public License
  24 along with the GNU MPFR Library; see the file COPYING.LIB.  If not, write to
  25 the Free Software Foundation, Inc., 51 Franklin St, Fifth Floor, Boston,
  26 MA 02110-1301, USA. */
  27
  28 #define MPFR_NEED_LONGLONG_H
  29 #include "mpfr-impl.h"
  30
  31 #ifdef MPFR_WANT_DECIMAL_FLOATS
  32
  33 #ifdef DPD_FORMAT
  34   /* conversion 10-bits to 3 digits */
  35 static unsigned int T[1024] = {
  36   0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 80, 81, 800, 801, 880, 881, 10, 11, 12, 13,
  37   14, 15, 16, 17, 18, 19, 90, 91, 810, 811, 890, 891, 20, 21, 22, 23, 24, 25,
  38   26, 27, 28, 29, 82, 83, 820, 821, 808, 809, 30, 31, 32, 33, 34, 35, 36, 37,
  39   38, 39, 92, 93, 830, 831, 818, 819, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49,
  40   84, 85, 840, 841, 88, 89, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 94, 95,
  41   850, 851, 98, 99, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 86, 87, 860, 861,
  42   888, 889, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 96, 97, 870, 871, 898,
  43   899, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 180, 181, 900, 901,
  44   980, 981, 110, 111, 112, 113, 114, 115, 116, 117, 118, 119, 190, 191, 910,
  45   911, 990, 991, 120, 121, 122, 123, 124, 125, 126, 127, 128, 129, 182, 183,
  46   920, 921, 908, 909, 130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 192,
  47   193, 930, 931, 918, 919, 140, 141, 142, 143, 144, 145, 146, 147, 148, 149,
  48   184, 185, 940, 941, 188, 189, 150, 151, 152, 153, 154, 155, 156, 157, 158,
  49   159, 194, 195, 950, 951, 198, 199, 160, 161, 162, 163, 164, 165, 166, 167,
  50   168, 169, 186, 187, 960, 961, 988, 989, 170, 171, 172, 173, 174, 175, 176,
  51   177, 178, 179, 196, 197, 970, 971, 998, 999, 200, 201, 202, 203, 204, 205,
  52   206, 207, 208, 209, 280, 281, 802, 803, 882, 883, 210, 211, 212, 213, 214,
  53   215, 216, 217, 218, 219, 290, 291, 812, 813, 892, 893, 220, 221, 222, 223,
  54   224, 225, 226, 227, 228, 229, 282, 283, 822, 823, 828, 829, 230, 231, 232,
  55   233, 234, 235, 236, 237, 238, 239, 292, 293, 832, 833, 838, 839, 240, 241,
  56   242, 243, 244, 245, 246, 247, 248, 249, 284, 285, 842, 843, 288, 289, 250,
  57   251, 252, 253, 254, 255, 256, 257, 258, 259, 294, 295, 852, 853, 298, 299,
  58   260, 261, 262, 263, 264, 265, 266, 267, 268, 269, 286, 287, 862, 863, 888,
  59   889, 270, 271, 272, 273, 274, 275, 276, 277, 278, 279, 296, 297, 872, 873,
  60   898, 899, 300, 301, 302, 303, 304, 305, 306, 307, 308, 309, 380, 381, 902,
  61   903, 982, 983, 310, 311, 312, 313, 314, 315, 316, 317, 318, 319, 390, 391,
  62   912, 913, 992, 993, 320, 321, 322, 323, 324, 325, 326, 327, 328, 329, 382,
  63   383, 922, 923, 928, 929, 330, 331, 332, 333, 334, 335, 336, 337, 338, 339,
  64   392, 393, 932, 933, 938, 939, 340, 341, 342, 343, 344, 345, 346, 347, 348,
  65   349, 384, 385, 942, 943, 388, 389, 350, 351, 352, 353, 354, 355, 356, 357,
  66   358, 359, 394, 395, 952, 953, 398, 399, 360, 361, 362, 363, 364, 365, 366,
  67   367, 368, 369, 386, 387, 962, 963, 988, 989, 370, 371, 372, 373, 374, 375,
  68   376, 377, 378, 379, 396, 397, 972, 973, 998, 999, 400, 401, 402, 403, 404,
  69   405, 406, 407, 408, 409, 480, 481, 804, 805, 884, 885, 410, 411, 412, 413,
  70   414, 415, 416, 417, 418, 419, 490, 491, 814, 815, 894, 895, 420, 421, 422,
  71   423, 424, 425, 426, 427, 428, 429, 482, 483, 824, 825, 848, 849, 430, 431,
  72   432, 433, 434, 435, 436, 437, 438, 439, 492, 493, 834, 835, 858, 859, 440,
  73   441, 442, 443, 444, 445, 446, 447, 448, 449, 484, 485, 844, 845, 488, 489,
  74   450, 451, 452, 453, 454, 455, 456, 457, 458, 459, 494, 495, 854, 855, 498,
  75   499, 460, 461, 462, 463, 464, 465, 466, 467, 468, 469, 486, 487, 864, 865,
  76   888, 889, 470, 471, 472, 473, 474, 475, 476, 477, 478, 479, 496, 497, 874,
  77   875, 898, 899, 500, 501, 502, 503, 504, 505, 506, 507, 508, 509, 580, 581,
  78   904, 905, 984, 985, 510, 511, 512, 513, 514, 515, 516, 517, 518, 519, 590,
  79   591, 914, 915, 994, 995, 520, 521, 522, 523, 524, 525, 526, 527, 528, 529,
  80   582, 583, 924, 925, 948, 949, 530, 531, 532, 533, 534, 535, 536, 537, 538,
  81   539, 592, 593, 934, 935, 958, 959, 540, 541, 542, 543, 544, 545, 546, 547,
  82   548, 549, 584, 585, 944, 945, 588, 589, 550, 551, 552, 553, 554, 555, 556,
  83   557, 558, 559, 594, 595, 954, 955, 598, 599, 560, 561, 562, 563, 564, 565,
  84   566, 567, 568, 569, 586, 587, 964, 965, 988, 989, 570, 571, 572, 573, 574,
  85   575, 576, 577, 578, 579, 596, 597, 974, 975, 998, 999, 600, 601, 602, 603,
  86   604, 605, 606, 607, 608, 609, 680, 681, 806, 807, 886, 887, 610, 611, 612,
  87   613, 614, 615, 616, 617, 618, 619, 690, 691, 816, 817, 896, 897, 620, 621,
  88   622, 623, 624, 625, 626, 627, 628, 629, 682, 683, 826, 827, 868, 869, 630,
  89   631, 632, 633, 634, 635, 636, 637, 638, 639, 692, 693, 836, 837, 878, 879,
  90   640, 641, 642, 643, 644, 645, 646, 647, 648, 649, 684, 685, 846, 847, 688,
  91   689, 650, 651, 652, 653, 654, 655, 656, 657, 658, 659, 694, 695, 856, 857,
  92   698, 699, 660, 661, 662, 663, 664, 665, 666, 667, 668, 669, 686, 687, 866,
  93   867, 888, 889, 670, 671, 672, 673, 674, 675, 676, 677, 678, 679, 696, 697,
  94   876, 877, 898, 899, 700, 701, 702, 703, 704, 705, 706, 707, 708, 709, 780,
  95   781, 906, 907, 986, 987, 710, 711, 712, 713, 714, 715, 716, 717, 718, 719,
  96   790, 791, 916, 917, 996, 997, 720, 721, 722, 723, 724, 725, 726, 727, 728,
  97   729, 782, 783, 926, 927, 968, 969, 730, 731, 732, 733, 734, 735, 736, 737,
  98   738, 739, 792, 793, 936, 937, 978, 979, 740, 741, 742, 743, 744, 745, 746,
  99   747, 748, 749, 784, 785, 946, 947, 788, 789, 750, 751, 752, 753, 754, 755,
 100   756, 757, 758, 759, 794, 795, 956, 957, 798, 799, 760, 761, 762, 763, 764,
 101   765, 766, 767, 768, 769, 786, 787, 966, 967, 988, 989, 770, 771, 772, 773,
 102   774, 775, 776, 777, 778, 779, 796, 797, 976, 977, 998, 999 };
 103 #endif
 104
 105 /* Convert d to a decimal string (one-to-one correspondence, no rounding).
 106    The string s needs to have at least 23 characters.
 107  */
 108 static void
 109 decimal64_to_string (char *s, _Decimal64 d)
 110 {
 111   union ieee_double_extract x;
 112   union ieee_double_decimal64 y;
 113   char *t;
 114   unsigned int Gh; /* most 5 significant bits from combination field */
 115   int exp; /* exponent */
 116   mp_limb_t rp[2];
 117   mp_size_t rn = 2;
 118   unsigned int i;
 119 #ifdef DPD_FORMAT
 120   unsigned int d0, d1, d2, d3, d4, d5;
 121 #endif
 122
 123   /* now convert BID or DPD to string */
 124   y.d64 = d;
 125   x.d = y.d;
 126   Gh = x.s.exp >> 6;
 127   if (Gh == 31)
 128     {
 129       sprintf (s, "NaN");
 130       return;
 131     }
 132   else if (Gh == 30)
 133     {
 134       if (x.s.sig == 0)
 135         sprintf (s, "Inf");
 136       else
 137         sprintf (s, "-Inf");
 138       return;
 139     }
 140   t = s;
 141   if (x.s.sig)
 142     *t++ = '-';
 143
 144 #ifdef DPD_FORMAT
 145   if (Gh < 24)
 146     {
 147       exp = (x.s.exp >> 1) & 768;
 148       d0 = Gh & 7;
 149     }
 150   else
 151     {
 152       exp = (x.s.exp & 384) << 1;
 153       d0 = 8 | (Gh & 1);
 154     }
 155   exp |= (x.s.exp & 63) << 2;
 156   exp |= x.s.manh >> 18;
 157   d1 = (x.s.manh >> 8) & 1023;
 158   d2 = ((x.s.manh << 2) | (x.s.manl >> 30)) & 1023;
 159   d3 = (x.s.manl >> 20) & 1023;
 160   d4 = (x.s.manl >> 10) & 1023;
 161   d5 = x.s.manl & 1023;
 162   sprintf (t, "%1u%3u%3u%3u%3u%3u", d0, T[d1], T[d2], T[d3], T[d4], T[d5]);
 163   /* Warning: some characters may be blank */
 164   for (i = 0; i < 16; i++)
 165     if (t[i] == ' ')
 166       t[i] = '0';
 167   t += 16;
 168 #else /* BID */
 169   if (Gh < 24)
 170     {
 171       /* the biased exponent E is formed from G[0] to G[9] and the
 172          significand from bits G[10] through the end of the decoding */
 173       exp = x.s.exp >> 1;
 174       /* manh has 20 bits, manl has 32 bits */
 175       rp[1] = ((x.s.exp & 1) << 20) | x.s.manh;
 176       rp[0] = x.s.manl;
 177     }
 178   else
 179     {
 180       /* the biased exponent is formed from G[2] to G[11] */
 181       exp = (x.s.exp & 511) << 1;
 182       rp[1] = x.s.manh;
 183       rp[0] = x.s.manl;
 184       exp |= rp[1] >> 19;
 185       rp[1] &= 524287; /* 2^19-1: cancel G[11] */
 186       rp[1] |= 2097152; /* add 2^21 */
 187     }
 188 #if BITS_PER_MP_LIMB >= 54
 189   rp[0] |= rp[1] << 32;
 190   rn = 1;
 191 #endif
 192   while (rn > 0 && rp[rn - 1] == 0)
 193     rn --;
 194   if (rn == 0)
 195     {
 196       *t = 0;
 197       i = 1;
 198     }
 199   else
 200     {
 201       i = mpn_get_str ((unsigned char*)t, 10, rp, rn);
 202     }
 203   while (i-- > 0)
 204     *t++ += '0';
 205 #endif /* DPD or BID */
 206
 207   exp -= 398; /* unbiased exponent */
 208   t += sprintf (t, "E%d", exp);
 209 }
 210
 211 int
 212 mpfr_set_decimal64 (mpfr_ptr r, _Decimal64 d, mp_rnd_t rnd_mode)
 213 {
 214   char s[23]; /* need 1 character for sign,
 215                      16 characters for mantissa,
 216                       1 character for exponent,
 217                       4 characters for exponent (including sign),
 218                       1 character for terminating \0. */
 219
 220   decimal64_to_string (s, d);
 221   return mpfr_set_str (r, s, 10, rnd_mode);
 222 }
 223
 224 #endif /* MPFR_WANT_DECIMAL_FLOATS */