source/libs/gmp/gmp-src/mpf/get_str.c

   1 /* mpf_get_str (digit_ptr, exp, base, n_digits, a) -- Convert the floating
   2    point number A to a base BASE number and store N_DIGITS raw digits at
   3    DIGIT_PTR, and the base BASE exponent in the word pointed to by EXP.  For
   4    example, the number 3.1416 would be returned as "31416" in DIGIT_PTR and
   5    1 in EXP.
   6
   7 Copyright 1993-1997, 2000-2003, 2005, 2006, 2011, 2015 Free Software
   8 Foundation, Inc.
   9
  10 This file is part of the GNU MP Library.
  11
  12 The GNU MP Library is free software; you can redistribute it and/or modify
  13 it under the terms of either:
  14
  15   * the GNU Lesser General Public License as published by the Free
  16     Software Foundation; either version 3 of the License, or (at your
  17     option) any later version.
  18
  19 or
  20
  21   * the GNU General Public License as published by the Free Software
  22     Foundation; either version 2 of the License, or (at your option) any
  23     later version.
  24
  25 or both in parallel, as here.
  26
  27 The GNU MP Library is distributed in the hope that it will be useful, but
  28 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
  29 or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
  30 for more details.
  31
  32 You should have received copies of the GNU General Public License and the
  33 GNU Lesser General Public License along with the GNU MP Library.  If not,
  34 see https://www.gnu.org/licenses/.  */
  35
  36 #include <stdlib.h>             /* for NULL */
  37 #include "gmp.h"
  38 #include "gmp-impl.h"
  39 #include "longlong.h"           /* for count_leading_zeros */
  40
  41 /* Could use some more work.
  42
  43    1. Allocation is excessive.  Try to combine areas.  Perhaps use result
  44       string area for temp limb space?
  45    2. We generate up to two limbs of extra digits.  This is because we don't
  46       check the exact number of bits in the input operand, and from that
  47       compute an accurate exponent (variable e in the code).  It would be
  48       cleaner and probably somewhat faster to change this.
  49 */
  50
  51 /* Compute base^exp and return the most significant prec limbs in rp[].
  52    Put the count of omitted low limbs in *ign.
  53    Return the actual size (which might be less than prec).
  54    Allocation of rp[] and the temporary tp[] should be 2*prec+2 limbs.  */
  55 static mp_size_t
  56 mpn_pow_1_highpart (mp_ptr rp, mp_size_t *ignp,
  57                     mp_limb_t base, unsigned long exp,
  58                     mp_size_t prec, mp_ptr tp)
  59 {
  60   mp_size_t ign;                /* counts number of ignored low limbs in r */
  61   mp_size_t off;                /* keeps track of offset where value starts */
  62   mp_ptr passed_rp = rp;
  63   mp_size_t rn;
  64   int cnt;
  65   int i;
  66
  67   if (exp == 0)
  68     {
  69       rp[0] = 1;
  70       *ignp = 0;
  71       return 1;
  72     }
  73
  74   rp[0] = base;
  75   rn = 1;
  76   off = 0;
  77   ign = 0;
  78   count_leading_zeros (cnt, exp);
  79   for (i = GMP_LIMB_BITS - cnt - 2; i >= 0; i--)
  80     {
  81       mpn_sqr (tp, rp + off, rn);
  82       rn = 2 * rn;
  83       rn -= tp[rn - 1] == 0;
  84       ign <<= 1;
  85
  86       off = 0;
  87       if (rn > prec)
  88         {
  89           ign += rn - prec;
  90           off = rn - prec;
  91           rn = prec;
  92         }
  93       MP_PTR_SWAP (rp, tp);
  94
  95       if (((exp >> i) & 1) != 0)
  96         {
  97           mp_limb_t cy;
  98           cy = mpn_mul_1 (rp, rp + off, rn, base);
  99           rp[rn] = cy;
 100           rn += cy != 0;
 101           off = 0;
 102         }
 103     }
 104
 105   if (rn > prec)
 106     {
 107       ASSERT (rn == prec + 1);
 108
 109       ign += rn - prec;
 110       rp += rn - prec;
 111       rn = prec;
 112     }
 113
 114   /* With somewhat less than 50% probability, we can skip this copy.  */
 115   if (passed_rp != rp + off)
 116     MPN_COPY_INCR (passed_rp, rp + off, rn);
 117   *ignp = ign;
 118   return rn;
 119 }
 120
 121 char *
 122 mpf_get_str (char *dbuf, mp_exp_t *exp, int base, size_t n_digits, mpf_srcptr u)
 123 {
 124   mp_exp_t ue;
 125   mp_size_t n_limbs_needed;
 126   size_t max_digits;
 127   mp_ptr up, pp, tp;
 128   mp_size_t un, pn, tn;
 129   unsigned char *tstr;
 130   mp_exp_t exp_in_base;
 131   size_t n_digits_computed;
 132   mp_size_t i;
 133   const char *num_to_text;
 134   size_t alloc_size = 0;
 135   char *dp;
 136   TMP_DECL;
 137
 138   up = PTR(u);
 139   un = ABSIZ(u);
 140   ue = EXP(u);
 141
 142   if (base >= 0)
 143     {
 144       num_to_text = "0123456789abcdefghijklmnopqrstuvwxyz";
 145       if (base <= 1)
 146         base = 10;
 147       else if (base > 36)
 148         {
 149           num_to_text = "0123456789ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz";
 150           if (base > 62)
 151             return NULL;
 152         }
 153     }
 154   else
 155     {
 156       base = -base;
 157       if (base <= 1)
 158         base = 10;
 159       else if (base > 36)
 160         return NULL;
 161       num_to_text = "0123456789ABCDEFGHIJKLMNOPQRSTUVWXYZ";
 162     }
 163
 164   MPF_SIGNIFICANT_DIGITS (max_digits, base, PREC(u));
 165   if (n_digits == 0 || n_digits > max_digits)
 166     n_digits = max_digits;
 167
 168   if (dbuf == 0)
 169     {
 170       /* We didn't get a string from the user.  Allocate one (and return
 171          a pointer to it) with space for `-' and terminating null.  */
 172       alloc_size = n_digits + 2;
 173       dbuf = (char *) (*__gmp_allocate_func) (n_digits + 2);
 174     }
 175
 176   if (un == 0)
 177     {
 178       *exp = 0;
 179       *dbuf = 0;
 180       n_digits = 0;
 181       goto done;
 182     }
 183
 184   TMP_MARK;
 185
 186   /* Allocate temporary digit space.  We can't put digits directly in the user
 187      area, since we generate more digits than requested.  (We allocate
 188      2 * GMP_LIMB_BITS extra bytes because of the digit block nature of the
 189      conversion.)  */
 190   tstr = (unsigned char *) TMP_ALLOC (n_digits + 2 * GMP_LIMB_BITS + 3);
 191
 192   LIMBS_PER_DIGIT_IN_BASE (n_limbs_needed, n_digits, base);
 193
 194   if (un > n_limbs_needed)
 195     {
 196       up += un - n_limbs_needed;
 197       un = n_limbs_needed;
 198     }
 199
 200   TMP_ALLOC_LIMBS_2 (pp, 2 * n_limbs_needed + 4,
 201                      tp, 2 * n_limbs_needed + 4);
 202
 203   if (ue <= n_limbs_needed)
 204     {
 205       /* We need to multiply number by base^n to get an n_digits integer part.  */
 206       mp_size_t n_more_limbs_needed, ign, off;
 207       unsigned long e;
 208
 209       n_more_limbs_needed = n_limbs_needed - ue;
 210       DIGITS_IN_BASE_PER_LIMB (e, n_more_limbs_needed, base);
 211
 212       pn = mpn_pow_1_highpart (pp, &ign, (mp_limb_t) base, e, n_limbs_needed + 1, tp);
 213       if (un > pn)
 214         mpn_mul (tp, up, un, pp, pn);   /* FIXME: mpn_mul_highpart */
 215       else
 216         mpn_mul (tp, pp, pn, up, un);   /* FIXME: mpn_mul_highpart */
 217       tn = un + pn;
 218       tn -= tp[tn - 1] == 0;
 219       off = un - ue - ign;
 220       if (off < 0)
 221         {
 222           MPN_COPY_DECR (tp - off, tp, tn);
 223           MPN_ZERO (tp, -off);
 224           tn -= off;
 225           off = 0;
 226         }
 227       n_digits_computed = mpn_get_str (tstr, base, tp + off, tn - off);
 228
 229       exp_in_base = n_digits_computed - e;
 230     }
 231   else
 232     {
 233       /* We need to divide number by base^n to get an n_digits integer part.  */
 234       mp_size_t n_less_limbs_needed, ign, off, xn;
 235       unsigned long e;
 236       mp_ptr dummyp, xp;
 237
 238       n_less_limbs_needed = ue - n_limbs_needed;
 239       DIGITS_IN_BASE_PER_LIMB (e, n_less_limbs_needed, base);
 240
 241       pn = mpn_pow_1_highpart (pp, &ign, (mp_limb_t) base, e, n_limbs_needed + 1, tp);
 242
 243       xn = n_limbs_needed + (n_less_limbs_needed-ign);
 244       xp = TMP_ALLOC_LIMBS (xn);
 245       off = xn - un;
 246       MPN_ZERO (xp, off);
 247       MPN_COPY (xp + off, up, un);
 248
 249       dummyp = TMP_ALLOC_LIMBS (pn);
 250       mpn_tdiv_qr (tp, dummyp, (mp_size_t) 0, xp, xn, pp, pn);
 251       tn = xn - pn + 1;
 252       tn -= tp[tn - 1] == 0;
 253       n_digits_computed = mpn_get_str (tstr, base, tp, tn);
 254
 255       exp_in_base = n_digits_computed + e;
 256     }
 257
 258   /* We should normally have computed too many digits.  Round the result
 259      at the point indicated by n_digits.  */
 260   if (n_digits_computed > n_digits)
 261     {
 262       size_t i;
 263       /* Round the result.  */
 264       if (tstr[n_digits] * 2 >= base)
 265         {
 266           n_digits_computed = n_digits;
 267           for (i = n_digits - 1;; i--)
 268             {
 269               unsigned int x;
 270               x = ++(tstr[i]);
 271               if (x != base)
 272                 break;
 273               n_digits_computed--;
 274               if (i == 0)
 275                 {
 276                   /* We had something like `bbbbbbb...bd', where 2*d >= base
 277                      and `b' denotes digit with significance base - 1.
 278                      This rounds up to `1', increasing the exponent.  */
 279                   tstr[0] = 1;
 280                   n_digits_computed = 1;
 281                   exp_in_base++;
 282                   break;
 283                 }
 284             }
 285         }
 286     }
 287
 288   /* We might have fewer digits than requested as a result of rounding above,
 289      (i.e. 0.999999 => 1.0) or because we have a number that simply doesn't
 290      need many digits in this base (e.g., 0.125 in base 10).  */
 291   if (n_digits > n_digits_computed)
 292     n_digits = n_digits_computed;
 293
 294   /* Remove trailing 0.  There can be many zeros.  */
 295   while (n_digits != 0 && tstr[n_digits - 1] == 0)
 296     n_digits--;
 297
 298   dp = dbuf + (SIZ(u) < 0);
 299
 300   /* Translate to ASCII and copy to result string.  */
 301   for (i = 0; i < n_digits; i++)
 302     dp[i] = num_to_text[tstr[i]];
 303   dp[n_digits] = 0;
 304
 305   *exp = exp_in_base;
 306
 307   if (SIZ(u) < 0)
 308     {
 309       dbuf[0] = '-';
 310       n_digits++;
 311     }
 312
 313   TMP_FREE;
 314
 315  done:
 316   /* If the string was alloced then resize it down to the actual space
 317      required.  */
 318   if (alloc_size != 0)
 319     {
 320       __GMP_REALLOCATE_FUNC_MAYBE_TYPE (dbuf, alloc_size, n_digits + 1, char);
 321     }
 322
 323   return dbuf;
 324 }