Annotation of OpenXM_contrib/gmp/mpf/get_str.c, Revision 1.1.1.4
1.1 maekawa 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:
1.1.1.4 ! ohara 7: Copyright 1993, 1994, 1995, 1996, 1997, 2000, 2001, 2002 Free Software
! 8: Foundation, Inc.
1.1 maekawa 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
1.1.1.2 maekawa 13: it under the terms of the GNU Lesser General Public License as published by
14: the Free Software Foundation; either version 2.1 of the License, or (at your
1.1 maekawa 15: option) any later version.
16:
17: The GNU MP Library is distributed in the hope that it will be useful, but
18: WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
1.1.1.2 maekawa 19: or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
1.1 maekawa 20: License for more details.
21:
1.1.1.2 maekawa 22: You should have received a copy of the GNU Lesser General Public License
1.1 maekawa 23: along with the GNU MP Library; see the file COPYING.LIB. If not, write to
24: the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
25: MA 02111-1307, USA. */
26:
1.1.1.4 ! ohara 27: #include <string.h> /* for strlen */
1.1 maekawa 28: #include "gmp.h"
29: #include "gmp-impl.h"
30: #include "longlong.h"
31:
32: /*
1.1.1.2 maekawa 33: The conversion routine works like this:
34:
35: 1. If U >= 1, compute U' = U / base**n, where n is chosen such that U' is
36: the largest number smaller than 1.
37: 2. Else, if U < 1, compute U' = U * base**n, where n is chosen such that U'
38: is the largest number smaller than 1.
39: 3. Convert U' (by repeatedly multiplying it by base). This process can
40: easily be interrupted when the needed number of digits are generated.
1.1 maekawa 41: */
42:
43: char *
44: mpf_get_str (char *digit_ptr, mp_exp_t *exp, int base, size_t n_digits, mpf_srcptr u)
45: {
1.1.1.2 maekawa 46: mp_ptr up;
1.1 maekawa 47: mp_size_t usize;
48: mp_exp_t uexp;
1.1.1.2 maekawa 49: mp_size_t prec;
1.1 maekawa 50: unsigned char *str;
51: char *num_to_text;
52: mp_ptr rp;
1.1.1.2 maekawa 53: mp_size_t rsize;
1.1 maekawa 54: mp_limb_t big_base;
55: size_t digits_computed_so_far;
56: int dig_per_u;
57: unsigned char *tstr;
58: mp_exp_t exp_in_base;
1.1.1.2 maekawa 59: int cnt;
1.1.1.4 ! ohara 60: size_t alloc_size = 0;
1.1 maekawa 61: TMP_DECL (marker);
62:
63: TMP_MARK (marker);
64: usize = u->_mp_size;
65: uexp = u->_mp_exp;
1.1.1.2 maekawa 66: prec = u->_mp_prec + 1;
1.1 maekawa 67:
68: if (base >= 0)
69: {
70: if (base == 0)
71: base = 10;
72: num_to_text = "0123456789abcdefghijklmnopqrstuvwxyz";
73: }
74: else
75: {
76: base = -base;
77: num_to_text = "0123456789ABCDEFGHIJKLMNOPQRSTUVWXYZ";
78: }
79:
80: /* Don't compute more digits than U can accurately represent.
81: Also, if 0 digits were requested, give *exactly* as many digits
82: as can be accurately represented. */
83: {
1.1.1.4 ! ohara 84: size_t max_digits;
! 85: MPF_SIGNIFICANT_DIGITS (max_digits, base, prec-1);
1.1 maekawa 86: if (n_digits == 0 || n_digits > max_digits)
87: n_digits = max_digits;
1.1.1.2 maekawa 88: #if 0
89: /* This seems to work, but check it better before enabling it. */
90: else
91: /* Limit the computation precision if only a limited digits are
92: desired. We could probably decrease both this, and avoid the +1
93: for setting prec above. */
94: prec = 2 + (mp_size_t)
1.1.1.4 ! ohara 95: (n_digits / (GMP_NUMB_BITS * __mp_bases[base].chars_per_bit_exactly));
1.1.1.2 maekawa 96: #endif
1.1 maekawa 97: }
98:
99: if (digit_ptr == 0)
100: {
101: /* We didn't get a string from the user. Allocate one (and return
102: a pointer to it) with space for `-' and terminating null. */
1.1.1.4 ! ohara 103: alloc_size = n_digits + 2;
! 104: digit_ptr = (char *) (*__gmp_allocate_func) (n_digits + 2);
1.1 maekawa 105: }
106:
107: if (usize == 0)
108: {
109: *exp = 0;
110: *digit_ptr = 0;
1.1.1.4 ! ohara 111: goto done;
1.1 maekawa 112: }
113:
114: str = (unsigned char *) digit_ptr;
115:
116: if (usize < 0)
117: {
118: *digit_ptr = '-';
119: str++;
120: usize = -usize;
121: }
122:
1.1.1.2 maekawa 123: up = PTR (u);
1.1 maekawa 124:
1.1.1.2 maekawa 125: if (uexp > 0)
1.1 maekawa 126: {
1.1.1.2 maekawa 127: /* U >= 1. Compute U' = U / base**n, where n is chosen such that U' < 1. */
128: mp_size_t ralloc;
129: mp_ptr tp;
130: int i;
1.1 maekawa 131:
1.1.1.2 maekawa 132: /* Limit the number of digits to develop for small integers. */
133: #if 0
134: if (exp_in_base < n_digits)
135: n_digits = exp_in_base;
136: #endif
1.1 maekawa 137:
1.1.1.2 maekawa 138: count_leading_zeros (cnt, up[usize - 1]);
1.1.1.4 ! ohara 139: exp_in_base = (((double) uexp * GMP_NUMB_BITS - cnt + GMP_NAIL_BITS)
1.1.1.2 maekawa 140: * __mp_bases[base].chars_per_bit_exactly);
141: exp_in_base += 1;
142:
143: ralloc = (prec + 1) * 2;
144: rp = (mp_ptr) TMP_ALLOC (ralloc * BYTES_PER_MP_LIMB);
145: tp = (mp_ptr) TMP_ALLOC (ralloc * BYTES_PER_MP_LIMB);
146:
147: rp[0] = base;
148: rsize = 1;
149: count_leading_zeros (cnt, exp_in_base);
1.1.1.4 ! ohara 150: for (i = GMP_LIMB_BITS - cnt - 2; i >= 0; i--)
1.1.1.2 maekawa 151: {
1.1.1.4 ! ohara 152: mpn_sqr_n (tp, rp, rsize);
1.1.1.2 maekawa 153: rsize = 2 * rsize;
154: rsize -= tp[rsize - 1] == 0;
1.1 maekawa 155:
1.1.1.2 maekawa 156: if (rsize > prec)
157: {
158: MPN_COPY (rp, tp + rsize - prec, prec + 1);
159: rsize = prec;
160: }
161: else
162: MP_PTR_SWAP (rp, tp);
1.1 maekawa 163:
1.1.1.2 maekawa 164: if (((exp_in_base >> i) & 1) != 0)
165: {
166: mp_limb_t cy;
167: cy = mpn_mul_1 (rp, rp, rsize, (mp_limb_t) base);
168: rp[rsize] = cy;
169: rsize += cy != 0;
170: }
171: }
172:
173: count_leading_zeros (cnt, rp[rsize - 1]);
1.1.1.4 ! ohara 174: cnt -= GMP_NAIL_BITS;
1.1.1.2 maekawa 175: if (cnt != 0)
1.1 maekawa 176: {
1.1.1.2 maekawa 177: mpn_lshift (rp, rp, rsize, cnt);
178:
179: if (usize < rsize)
180: {
181: /* Pad out U to the size of R while shifting it.
182: (Reuse temporary space at tp.) */
183: mp_limb_t cy;
184:
185: MPN_ZERO (tp, rsize - usize);
186: cy = mpn_lshift (tp + rsize - usize, up, usize, cnt);
187: up = tp;
188: usize = rsize;
189: if (cy)
190: up[usize++] = cy;
191: ASSERT_ALWAYS (usize <= ralloc); /* sufficient space? */
192: }
193: else
194: {
195: /* Copy U to temporary space. */
196: /* FIXME: Allocate more space for tp above, and reuse it here. */
197: mp_limb_t cy;
198: mp_ptr tup = (mp_ptr) TMP_ALLOC ((usize + 1) * BYTES_PER_MP_LIMB);
199:
200: cy = mpn_lshift (tup, up, usize, cnt);
201: up = tup;
202: if (cy)
203: up[usize++] = cy;
204: }
1.1 maekawa 205: }
206: else
207: {
1.1.1.2 maekawa 208: if (usize < rsize)
209: {
210: /* Pad out U to the size of R. (Reuse temporary space at tp.) */
211: MPN_ZERO (tp, rsize - usize);
212: MPN_COPY (tp + rsize - usize, up, usize);
213: up = tp;
214: usize = rsize;
215: }
216: else
217: {
218: /* Copy U to temporary space. */
219: mp_ptr tmp = (mp_ptr) TMP_ALLOC (usize * BYTES_PER_MP_LIMB);
220: MPN_COPY (tmp, up, usize);
221: up = tmp;
222: }
223: }
224:
225: {
226: mp_ptr qp;
227: qp = (mp_ptr) TMP_ALLOC (prec * BYTES_PER_MP_LIMB);
228: mpn_divrem (qp, prec - (usize - rsize), up, usize, rp, rsize);
229: rsize = prec;
230: rp = qp;
231: }
232: }
233: else
234: {
235: /* U < 1. Compute U' = U * base**n, where n is chosen such that U' is
236: the greatest number that still satisfies U' < 1. */
237: mp_size_t ralloc;
238: mp_ptr tp;
239: int i;
240:
241: uexp = -uexp;
242: count_leading_zeros (cnt, up[usize - 1]);
1.1.1.4 ! ohara 243: exp_in_base = (((double) uexp * GMP_NUMB_BITS + cnt - GMP_NAIL_BITS - 1)
1.1.1.2 maekawa 244: * __mp_bases[base].chars_per_bit_exactly);
245: if (exp_in_base < 0)
246: exp_in_base = 0;
247:
248: if (exp_in_base != 0)
249: {
250: ralloc = (prec + 1) * 2;
251: rp = (mp_ptr) TMP_ALLOC (ralloc * BYTES_PER_MP_LIMB);
252: tp = (mp_ptr) TMP_ALLOC (ralloc * BYTES_PER_MP_LIMB);
1.1 maekawa 253:
254: rp[0] = base;
255: rsize = 1;
256: count_leading_zeros (cnt, exp_in_base);
1.1.1.4 ! ohara 257: for (i = GMP_LIMB_BITS - cnt - 2; i >= 0; i--)
1.1 maekawa 258: {
1.1.1.4 ! ohara 259: mpn_sqr_n (tp, rp, rsize);
1.1 maekawa 260: rsize = 2 * rsize;
261: rsize -= tp[rsize - 1] == 0;
1.1.1.2 maekawa 262: if (rsize > prec)
263: {
264: MPN_COPY (rp, tp + rsize - prec, prec + 1);
265: rsize = prec;
266: }
267: else
268: MP_PTR_SWAP (rp, tp);
1.1 maekawa 269:
270: if (((exp_in_base >> i) & 1) != 0)
271: {
1.1.1.2 maekawa 272: mp_limb_t cy;
1.1 maekawa 273: cy = mpn_mul_1 (rp, rp, rsize, (mp_limb_t) base);
274: rp[rsize] = cy;
275: rsize += cy != 0;
276: }
277: }
278:
279: {
1.1.1.2 maekawa 280: mp_limb_t cy;
281: tp = (mp_ptr) TMP_ALLOC ((rsize + usize) * BYTES_PER_MP_LIMB);
282: if (rsize > usize)
283: cy = mpn_mul (tp, rp, rsize, up, usize);
1.1 maekawa 284: else
1.1.1.2 maekawa 285: cy = mpn_mul (tp, up, usize, rp, rsize);
286: rsize += usize;
287: rsize -= cy == 0;
288: rp = tp;
1.1 maekawa 289: }
1.1.1.2 maekawa 290: exp_in_base = -exp_in_base;
1.1 maekawa 291: }
1.1.1.2 maekawa 292: else
1.1 maekawa 293: {
1.1.1.2 maekawa 294: rp = (mp_ptr) TMP_ALLOC (usize * BYTES_PER_MP_LIMB);
295: MPN_COPY (rp, up, usize);
296: rsize = usize;
1.1 maekawa 297: }
298: }
299:
1.1.1.2 maekawa 300: big_base = __mp_bases[base].big_base;
301: dig_per_u = __mp_bases[base].chars_per_limb;
1.1 maekawa 302:
1.1.1.2 maekawa 303: /* Hack for correctly (although not optimally) converting to bases that are
304: powers of 2. If we deem it important, we could handle powers of 2 by
305: shifting and masking (just like mpn_get_str). */
306: if (big_base < 10) /* logarithm of base when power of two */
1.1 maekawa 307: {
1.1.1.2 maekawa 308: int logbase = big_base;
1.1.1.4 ! ohara 309: if (dig_per_u * logbase == GMP_NUMB_BITS)
1.1.1.2 maekawa 310: dig_per_u--;
311: big_base = (mp_limb_t) 1 << (dig_per_u * logbase);
312: /* fall out to general code... */
313: }
1.1 maekawa 314:
1.1.1.2 maekawa 315: /* Now that we have normalized the number, develop the digits, essentially by
316: multiplying it by BASE. We initially develop at least 3 extra digits,
317: since the two leading digits might become zero, and we need one extra for
318: rounding the output properly. */
1.1 maekawa 319:
1.1.1.2 maekawa 320: /* Allocate temporary digit space. We can't put digits directly in the user
321: area, since we generate more digits than requested. (We allocate
1.1.1.4 ! ohara 322: GMP_LIMB_BITS extra bytes because of the digit block nature of the
1.1.1.2 maekawa 323: conversion.) */
1.1.1.4 ! ohara 324: tstr = (unsigned char *) TMP_ALLOC (n_digits + GMP_LIMB_BITS + 3);
1.1 maekawa 325:
1.1.1.2 maekawa 326: for (digits_computed_so_far = 0; digits_computed_so_far < n_digits + 3;
327: digits_computed_so_far += dig_per_u)
328: {
329: mp_limb_t cy;
1.1.1.3 maekawa 330: /* For speed: skip trailing zero limbs. */
1.1.1.2 maekawa 331: if (rp[0] == 0)
1.1 maekawa 332: {
1.1.1.2 maekawa 333: rp++;
334: rsize--;
335: if (rsize == 0)
1.1.1.3 maekawa 336: break;
1.1 maekawa 337: }
338:
1.1.1.2 maekawa 339: cy = mpn_mul_1 (rp, rp, rsize, big_base);
1.1 maekawa 340:
1.1.1.4 ! ohara 341: if (! POW2_P (GMP_NUMB_BITS))
! 342: {
! 343: if (digits_computed_so_far == 0 && cy == 0)
! 344: cy = mpn_mul_1 (rp, rp, rsize, big_base);
! 345: }
! 346:
1.1.1.2 maekawa 347: ASSERT_ALWAYS (! (digits_computed_so_far == 0 && cy == 0));
1.1 maekawa 348:
1.1.1.2 maekawa 349: /* Convert N1 from BIG_BASE to a string of digits in BASE
350: using single precision operations. */
1.1 maekawa 351: {
1.1.1.2 maekawa 352: int i;
353: unsigned char *s = tstr + digits_computed_so_far + dig_per_u;
354: for (i = dig_per_u - 1; i >= 0; i--)
1.1 maekawa 355: {
1.1.1.2 maekawa 356: *--s = cy % base;
357: cy /= base;
1.1 maekawa 358: }
359: }
1.1.1.2 maekawa 360: }
1.1 maekawa 361:
1.1.1.2 maekawa 362: /* We can have at most two leading 0. Remove them. */
1.1 maekawa 363: if (tstr[0] == 0)
364: {
365: tstr++;
366: digits_computed_so_far--;
367: exp_in_base--;
368:
1.1.1.2 maekawa 369: if (tstr[0] == 0)
1.1 maekawa 370: {
1.1.1.2 maekawa 371: tstr++;
372: digits_computed_so_far--;
373: exp_in_base--;
374:
1.1.1.3 maekawa 375: ASSERT_ALWAYS (tstr[0] != 0);
1.1 maekawa 376: }
377: }
378:
1.1.1.3 maekawa 379: /* We should normally have computed too many digits. Round the result
380: at the point indicated by n_digits. */
381: if (digits_computed_so_far > n_digits)
382: {
383: size_t i;
384: /* Round the result. */
385: if (tstr[n_digits] * 2 >= base)
386: {
387: digits_computed_so_far = n_digits;
388: for (i = n_digits - 1;; i--)
389: {
390: unsigned int x;
391: x = ++(tstr[i]);
392: if (x != base)
393: break;
394: digits_computed_so_far--;
395: if (i == 0)
396: {
397: /* We had something like `bbbbbbb...bd', where 2*d >= base
398: and `b' denotes digit with significance base - 1.
399: This rounds up to `1', increasing the exponent. */
400: tstr[0] = 1;
401: digits_computed_so_far = 1;
402: exp_in_base++;
1.1.1.2 maekawa 403: break;
1.1.1.3 maekawa 404: }
405: }
406: }
407: }
1.1.1.2 maekawa 408:
1.1.1.3 maekawa 409: /* We might have fewer digits than requested as a result of rounding above,
410: (i.e. 0.999999 => 1.0) or because we have a number that simply doesn't
411: need many digits in this base (e.g., 0.125 in base 10). */
412: if (n_digits > digits_computed_so_far)
413: n_digits = digits_computed_so_far;
414:
415: /* Remove trailing 0. There can be many zeros. */
416: while (n_digits != 0 && tstr[n_digits - 1] == 0)
417: n_digits--;
1.1.1.2 maekawa 418:
1.1.1.3 maekawa 419: /* Translate to ASCII and copy to result string. */
420: while (n_digits != 0)
421: {
422: *str++ = num_to_text[*tstr++];
1.1.1.2 maekawa 423: n_digits--;
1.1.1.3 maekawa 424: }
1.1.1.2 maekawa 425:
1.1 maekawa 426: *str = 0;
427: *exp = exp_in_base;
1.1.1.4 ! ohara 428:
! 429: done:
1.1 maekawa 430: TMP_FREE (marker);
1.1.1.4 ! ohara 431:
! 432: /* If the string was alloced then resize it down to the actual space
! 433: required. */
! 434: if (alloc_size != 0)
! 435: {
! 436: size_t actual_size = strlen (digit_ptr) + 1;
! 437: __GMP_REALLOCATE_FUNC_MAYBE_TYPE (digit_ptr, alloc_size, actual_size,
! 438: char);
! 439: }
! 440:
1.1 maekawa 441: return digit_ptr;
442: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to get_str.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / gmp / mpf