Home | History | Annotate | Line # | Download | only in libbid
      1 /* Copyright (C) 2007-2024 Free Software Foundation, Inc.
      2 
      3 This file is part of GCC.
      4 
      5 GCC is free software; you can redistribute it and/or modify it under
      6 the terms of the GNU General Public License as published by the Free
      7 Software Foundation; either version 3, or (at your option) any later
      8 version.
      9 
     10 GCC is distributed in the hope that it will be useful, but WITHOUT ANY
     11 WARRANTY; without even the implied warranty of MERCHANTABILITY or
     12 FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
     13 for more details.
     14 
     15 Under Section 7 of GPL version 3, you are granted additional
     16 permissions described in the GCC Runtime Library Exception, version
     17 3.1, as published by the Free Software Foundation.
     18 
     19 You should have received a copy of the GNU General Public License and
     20 a copy of the GCC Runtime Library Exception along with this program;
     21 see the files COPYING3 and COPYING.RUNTIME respectively.  If not, see
     22 <http://www.gnu.org/licenses/>.  */
     23 
     24 #include <ctype.h>
     25 #include "bid_internal.h"
     26 #include "bid128_2_str.h"
     27 #include "bid128_2_str_macros.h"
     28 
     29 #define MAX_FORMAT_DIGITS     16
     30 #define DECIMAL_EXPONENT_BIAS 398
     31 #define MAX_DECIMAL_EXPONENT  767
     32 
     33 #if DECIMAL_CALL_BY_REFERENCE
     34 
     35 void
     36 bid64_to_string (char *ps, UINT64 * px
     37 		 _EXC_FLAGS_PARAM _EXC_MASKS_PARAM _EXC_INFO_PARAM) {
     38   UINT64 x;
     39 #else
     40 
     41 void
     42 bid64_to_string (char *ps, UINT64 x
     43 		 _EXC_FLAGS_PARAM _EXC_MASKS_PARAM _EXC_INFO_PARAM) {
     44 #endif
     45 // the destination string (pointed to by ps) must be pre-allocated
     46   UINT64 sign_x, coefficient_x, D, ER10;
     47   int istart, exponent_x, j, digits_x, bin_expon_cx;
     48   int_float tempx;
     49   UINT32 MiDi[12], *ptr;
     50   UINT64 HI_18Dig, LO_18Dig, Tmp;
     51   char *c_ptr_start, *c_ptr;
     52   int midi_ind, k_lcv, len;
     53   unsigned int save_fpsf;
     54 
     55 #if DECIMAL_CALL_BY_REFERENCE
     56   x = *px;
     57 #endif
     58 
     59   save_fpsf = *pfpsf; // place holder only
     60   // unpack arguments, check for NaN or Infinity
     61   if (!unpack_BID64 (&sign_x, &exponent_x, &coefficient_x, x)) {
     62     // x is Inf. or NaN or 0
     63 
     64     // Inf or NaN?
     65     if ((x & 0x7800000000000000ull) == 0x7800000000000000ull) {
     66       if ((x & 0x7c00000000000000ull) == 0x7c00000000000000ull) {
     67     ps[0] = (sign_x) ? '-' : '+';
     68     ps[1] = ((x & MASK_SNAN) == MASK_SNAN)? 'S':'Q';
     69 	ps[2] = 'N';
     70 	ps[3] = 'a';
     71 	ps[4] = 'N';
     72 	ps[5] = 0;
     73 	return;
     74       }
     75       // x is Inf
     76       ps[0] = (sign_x) ? '-' : '+';
     77       ps[1] = 'I';
     78       ps[2] = 'n';
     79       ps[3] = 'f';
     80       ps[4] = 0;
     81       return;
     82     }
     83     // 0
     84     istart = 0;
     85     if (sign_x) {
     86       ps[istart++] = '-';
     87     }
     88 
     89     ps[istart++] = '0';
     90     ps[istart++] = 'E';
     91 
     92     exponent_x -= 398;
     93     if (exponent_x < 0) {
     94       ps[istart++] = '-';
     95       exponent_x = -exponent_x;
     96     } else
     97       ps[istart++] = '+';
     98 
     99     if (exponent_x) {
    100       // get decimal digits in coefficient_x
    101       tempx.d = (float) exponent_x;
    102       bin_expon_cx = ((tempx.i >> 23) & 0xff) - 0x7f;
    103       digits_x = estimate_decimal_digits[bin_expon_cx];
    104       if ((UINT64)exponent_x >= power10_table_128[digits_x].w[0])
    105 	digits_x++;
    106 
    107       j = istart + digits_x - 1;
    108       istart = j + 1;
    109 
    110       // 2^32/10
    111       ER10 = 0x1999999a;
    112 
    113       while (exponent_x > 9) {
    114 	D = (UINT64) exponent_x *ER10;
    115 	D >>= 32;
    116 	exponent_x = exponent_x - (D << 1) - (D << 3);
    117 
    118 	ps[j--] = '0' + (char) exponent_x;
    119 	exponent_x = D;
    120       }
    121       ps[j] = '0' + (char) exponent_x;
    122     } else {
    123       ps[istart++] = '0';
    124     }
    125 
    126     ps[istart] = 0;
    127 
    128     return;
    129   }
    130   // convert expon, coeff to ASCII
    131   exponent_x -= DECIMAL_EXPONENT_BIAS;
    132 
    133   ER10 = 0x1999999a;
    134 
    135   istart = 0;
    136   if (sign_x) {
    137     ps[0] = '-';
    138     istart = 1;
    139   }
    140   // if zero or non-canonical, set coefficient to '0'
    141   if ((coefficient_x > 9999999999999999ull) ||	// non-canonical
    142       ((coefficient_x == 0))	// significand is zero
    143     ) {
    144     ps[istart++] = '0';
    145   } else {
    146     /* ****************************************************
    147        This takes a bid coefficient in C1.w[1],C1.w[0]
    148        and put the converted character sequence at location
    149        starting at &(str[k]). The function returns the number
    150        of MiDi returned. Note that the character sequence
    151        does not have leading zeros EXCEPT when the input is of
    152        zero value. It will then output 1 character '0'
    153        The algorithm essentailly tries first to get a sequence of
    154        Millenial Digits "MiDi" and then uses table lookup to get the
    155        character strings of these MiDis.
    156        **************************************************** */
    157     /* Algorithm first decompose possibly 34 digits in hi and lo
    158        18 digits. (The high can have at most 16 digits). It then
    159        uses macro that handle 18 digit portions.
    160        The first step is to get hi and lo such that
    161        2^(64) C1.w[1] + C1.w[0] = hi * 10^18  + lo,   0 <= lo < 10^18.
    162        We use a table lookup method to obtain the hi and lo 18 digits.
    163        [C1.w[1],C1.w[0]] = c_8 2^(107) + c_7 2^(101) + ... + c_0 2^(59) + d
    164        where 0 <= d < 2^59 and each c_j has 6 bits. Because d fits in
    165        18 digits,  we set hi = 0, and lo = d to begin with.
    166        We then retrieve from a table, for j = 0, 1, ..., 8
    167        that gives us A and B where c_j 2^(59+6j) = A * 10^18 + B.
    168        hi += A ; lo += B; After each accumulation into lo, we normalize
    169        immediately. So at the end, we have the decomposition as we need. */
    170 
    171     Tmp = coefficient_x >> 59;
    172     LO_18Dig = (coefficient_x << 5) >> 5;
    173     HI_18Dig = 0;
    174     k_lcv = 0;
    175 
    176     while (Tmp) {
    177       midi_ind = (int) (Tmp & 0x000000000000003FLL);
    178       midi_ind <<= 1;
    179       Tmp >>= 6;
    180       HI_18Dig += mod10_18_tbl[k_lcv][midi_ind++];
    181       LO_18Dig += mod10_18_tbl[k_lcv++][midi_ind];
    182       __L0_Normalize_10to18 (HI_18Dig, LO_18Dig);
    183     }
    184 
    185     ptr = MiDi;
    186     __L1_Split_MiDi_6_Lead (LO_18Dig, ptr);
    187     len = ptr - MiDi;
    188     c_ptr_start = &(ps[istart]);
    189     c_ptr = c_ptr_start;
    190 
    191     /* now convert the MiDi into character strings */
    192     __L0_MiDi2Str_Lead (MiDi[0], c_ptr);
    193     for (k_lcv = 1; k_lcv < len; k_lcv++) {
    194       __L0_MiDi2Str (MiDi[k_lcv], c_ptr);
    195     }
    196     istart = istart + (c_ptr - c_ptr_start);
    197   }
    198 
    199   ps[istart++] = 'E';
    200 
    201   if (exponent_x < 0) {
    202     ps[istart++] = '-';
    203     exponent_x = -exponent_x;
    204   } else
    205     ps[istart++] = '+';
    206 
    207   if (exponent_x) {
    208     // get decimal digits in coefficient_x
    209     tempx.d = (float) exponent_x;
    210     bin_expon_cx = ((tempx.i >> 23) & 0xff) - 0x7f;
    211     digits_x = estimate_decimal_digits[bin_expon_cx];
    212     if ((UINT64)exponent_x >= power10_table_128[digits_x].w[0])
    213       digits_x++;
    214 
    215     j = istart + digits_x - 1;
    216     istart = j + 1;
    217 
    218     // 2^32/10
    219     ER10 = 0x1999999a;
    220 
    221     while (exponent_x > 9) {
    222       D = (UINT64) exponent_x *ER10;
    223       D >>= 32;
    224       exponent_x = exponent_x - (D << 1) - (D << 3);
    225 
    226       ps[j--] = '0' + (char) exponent_x;
    227       exponent_x = D;
    228     }
    229     ps[j] = '0' + (char) exponent_x;
    230   } else {
    231     ps[istart++] = '0';
    232   }
    233 
    234   ps[istart] = 0;
    235 
    236   return;
    237 
    238 }
    239 
    240 
    241 #if DECIMAL_CALL_BY_REFERENCE
    242 void
    243 bid64_from_string (UINT64 * pres, char *ps
    244 		   _RND_MODE_PARAM _EXC_FLAGS_PARAM
    245                    _EXC_MASKS_PARAM _EXC_INFO_PARAM) {
    246 #else
    247 UINT64
    248 bid64_from_string (char *ps
    249 		   _RND_MODE_PARAM _EXC_FLAGS_PARAM
    250                    _EXC_MASKS_PARAM _EXC_INFO_PARAM) {
    251 #endif
    252   UINT64 sign_x, coefficient_x = 0, rounded = 0, res;
    253   int expon_x = 0, sgn_expon, ndigits, add_expon = 0, midpoint =
    254     0, rounded_up = 0;
    255   int dec_expon_scale = 0, right_radix_leading_zeros = 0, rdx_pt_enc =
    256     0;
    257   unsigned fpsc;
    258   char c;
    259   unsigned int save_fpsf;
    260 
    261 #if DECIMAL_CALL_BY_REFERENCE
    262 #if !DECIMAL_GLOBAL_ROUNDING
    263   _IDEC_round rnd_mode = *prnd_mode;
    264 #endif
    265 #endif
    266 
    267   save_fpsf = *pfpsf; // place holder only
    268   // eliminate leading whitespace
    269   while (((*ps == ' ') || (*ps == '\t')) && (*ps))
    270     ps++;
    271 
    272   // get first non-whitespace character
    273   c = *ps;
    274 
    275   // detect special cases (INF or NaN)
    276   if (!c || (c != '.' && c != '-' && c != '+' && (c < '0' || c > '9'))) {
    277     // Infinity?
    278     if ((tolower_macro (ps[0]) == 'i' && tolower_macro (ps[1]) == 'n' &&
    279         tolower_macro (ps[2]) == 'f') && (!ps[3] ||
    280         (tolower_macro (ps[3]) == 'i' &&
    281         tolower_macro (ps[4]) == 'n' && tolower_macro (ps[5]) == 'i' &&
    282         tolower_macro (ps[6]) == 't' && tolower_macro (ps[7]) == 'y' &&
    283         !ps[8]))) {
    284       res = 0x7800000000000000ull;
    285       BID_RETURN (res);
    286     }
    287     // return sNaN
    288     if (tolower_macro (ps[0]) == 's' && tolower_macro (ps[1]) == 'n' &&
    289         tolower_macro (ps[2]) == 'a' && tolower_macro (ps[3]) == 'n') {
    290         // case insensitive check for snan
    291       res = 0x7e00000000000000ull;
    292       BID_RETURN (res);
    293     } else {
    294       // return qNaN
    295       res = 0x7c00000000000000ull;
    296       BID_RETURN (res);
    297     }
    298   }
    299   // detect +INF or -INF
    300   if ((tolower_macro (ps[1]) == 'i' && tolower_macro (ps[2]) == 'n' &&
    301       tolower_macro (ps[3]) == 'f') && (!ps[4] ||
    302       (tolower_macro (ps[4]) == 'i' && tolower_macro (ps[5]) == 'n' &&
    303       tolower_macro (ps[6]) == 'i' && tolower_macro (ps[7]) == 't' &&
    304       tolower_macro (ps[8]) == 'y' && !ps[9]))) {
    305     if (c == '+')
    306       res = 0x7800000000000000ull;
    307     else if (c == '-')
    308       res = 0xf800000000000000ull;
    309     else
    310       res = 0x7c00000000000000ull;
    311     BID_RETURN (res);
    312   }
    313   // if +sNaN, +SNaN, -sNaN, or -SNaN
    314   if (tolower_macro (ps[1]) == 's' && tolower_macro (ps[2]) == 'n'
    315       && tolower_macro (ps[3]) == 'a' && tolower_macro (ps[4]) == 'n') {
    316     if (c == '-')
    317       res = 0xfe00000000000000ull;
    318     else
    319       res = 0x7e00000000000000ull;
    320     BID_RETURN (res);
    321   }
    322   // determine sign
    323   if (c == '-')
    324     sign_x = 0x8000000000000000ull;
    325   else
    326     sign_x = 0;
    327 
    328   // get next character if leading +/- sign
    329   if (c == '-' || c == '+') {
    330     ps++;
    331     c = *ps;
    332   }
    333   // if c isn't a decimal point or a decimal digit, return NaN
    334   if (c != '.' && (c < '0' || c > '9')) {
    335     // return NaN
    336     res = 0x7c00000000000000ull | sign_x;
    337     BID_RETURN (res);
    338   }
    339 
    340   rdx_pt_enc = 0;
    341 
    342   // detect zero (and eliminate/ignore leading zeros)
    343   if (*(ps) == '0' || *(ps) == '.') {
    344 
    345     if (*(ps) == '.') {
    346       rdx_pt_enc = 1;
    347       ps++;
    348     }
    349     // if all numbers are zeros (with possibly 1 radix point, the number is zero
    350     // should catch cases such as: 000.0
    351     while (*ps == '0') {
    352       ps++;
    353       // for numbers such as 0.0000000000000000000000000000000000001001,
    354       // we want to count the leading zeros
    355       if (rdx_pt_enc) {
    356 	right_radix_leading_zeros++;
    357       }
    358       // if this character is a radix point, make sure we haven't already
    359       // encountered one
    360       if (*(ps) == '.') {
    361 	if (rdx_pt_enc == 0) {
    362 	  rdx_pt_enc = 1;
    363 	  // if this is the first radix point, and the next character is NULL,
    364           // we have a zero
    365 	  if (!*(ps + 1)) {
    366 	    res =
    367 	      ((UINT64) (398 - right_radix_leading_zeros) << 53) |
    368 	      sign_x;
    369 	    BID_RETURN (res);
    370 	  }
    371 	  ps = ps + 1;
    372 	} else {
    373 	  // if 2 radix points, return NaN
    374 	  res = 0x7c00000000000000ull | sign_x;
    375 	  BID_RETURN (res);
    376 	}
    377       } else if (!*(ps)) {
    378 	//pres->w[1] = 0x3040000000000000ull | sign_x;
    379 	res =
    380 	  ((UINT64) (398 - right_radix_leading_zeros) << 53) | sign_x;
    381 	BID_RETURN (res);
    382       }
    383     }
    384   }
    385 
    386   c = *ps;
    387 
    388   ndigits = 0;
    389   while ((c >= '0' && c <= '9') || c == '.') {
    390     if (c == '.') {
    391       if (rdx_pt_enc) {
    392 	// return NaN
    393 	res = 0x7c00000000000000ull | sign_x;
    394 	BID_RETURN (res);
    395       }
    396       rdx_pt_enc = 1;
    397       ps++;
    398       c = *ps;
    399       continue;
    400     }
    401     dec_expon_scale += rdx_pt_enc;
    402 
    403     ndigits++;
    404     if (ndigits <= 16) {
    405       coefficient_x = (coefficient_x << 1) + (coefficient_x << 3);
    406       coefficient_x += (UINT64) (c - '0');
    407     } else if (ndigits == 17) {
    408       // coefficient rounding
    409 		switch(rnd_mode){
    410 	case ROUNDING_TO_NEAREST:
    411       midpoint = (c == '5' && !(coefficient_x & 1)) ? 1 : 0;
    412           // if coefficient is even and c is 5, prepare to round up if
    413           // subsequent digit is nonzero
    414       // if str[MAXDIG+1] > 5, we MUST round up
    415       // if str[MAXDIG+1] == 5 and coefficient is ODD, ROUND UP!
    416       if (c > '5' || (c == '5' && (coefficient_x & 1))) {
    417 	coefficient_x++;
    418 	rounded_up = 1;
    419 	break;
    420 
    421 	case ROUNDING_DOWN:
    422 		if(sign_x) { coefficient_x++; rounded_up=1; }
    423 		break;
    424 	case ROUNDING_UP:
    425 		if(!sign_x) { coefficient_x++; rounded_up=1; }
    426 		break;
    427 	case ROUNDING_TIES_AWAY:
    428 		if(c>='5') { coefficient_x++; rounded_up=1; }
    429 		break;
    430 	  }
    431 	if (coefficient_x == 10000000000000000ull) {
    432 	  coefficient_x = 1000000000000000ull;
    433 	  add_expon = 1;
    434 	}
    435       }
    436       if (c > '0')
    437 	rounded = 1;
    438       add_expon += 1;
    439     } else { // ndigits > 17
    440       add_expon++;
    441       if (midpoint && c > '0') {
    442 	coefficient_x++;
    443 	midpoint = 0;
    444 	rounded_up = 1;
    445       }
    446       if (c > '0')
    447 	rounded = 1;
    448     }
    449     ps++;
    450     c = *ps;
    451   }
    452 
    453   add_expon -= (dec_expon_scale + right_radix_leading_zeros);
    454 
    455   if (!c) {
    456     res =
    457       fast_get_BID64_check_OF (sign_x,
    458 			       add_expon + DECIMAL_EXPONENT_BIAS,
    459 			       coefficient_x, 0, &fpsc);
    460     BID_RETURN (res);
    461   }
    462 
    463   if (c != 'E' && c != 'e') {
    464     // return NaN
    465     res = 0x7c00000000000000ull | sign_x;
    466     BID_RETURN (res);
    467   }
    468   ps++;
    469   c = *ps;
    470   sgn_expon = (c == '-') ? 1 : 0;
    471   if (c == '-' || c == '+') {
    472     ps++;
    473     c = *ps;
    474   }
    475   if (!c || c < '0' || c > '9') {
    476     // return NaN
    477     res = 0x7c00000000000000ull | sign_x;
    478     BID_RETURN (res);
    479   }
    480 
    481   while (c >= '0' && c <= '9') {
    482     expon_x = (expon_x << 1) + (expon_x << 3);
    483     expon_x += (int) (c - '0');
    484 
    485     ps++;
    486     c = *ps;
    487   }
    488 
    489   if (c) {
    490     // return NaN
    491     res = 0x7c00000000000000ull | sign_x;
    492     BID_RETURN (res);
    493   }
    494 
    495   if (sgn_expon)
    496     expon_x = -expon_x;
    497 
    498   expon_x += add_expon + DECIMAL_EXPONENT_BIAS;
    499 
    500   if (expon_x < 0) {
    501     if (rounded_up)
    502       coefficient_x--;
    503     rnd_mode = 0;
    504     res =
    505       get_BID64_UF (sign_x, expon_x, coefficient_x, rounded, rnd_mode,
    506 		    &fpsc);
    507     BID_RETURN (res);
    508   }
    509   res = get_BID64 (sign_x, expon_x, coefficient_x, rnd_mode, &fpsc);
    510   BID_RETURN (res);
    511 }
    512