Annotation of OpenXM_contrib/gmp/mpq/set_d.c, Revision 1.1.1.2
1.1 maekawa 1: /* mpq_set_d(mpq_t q, double d) -- Set q to d without rounding.
2:
1.1.1.2 ! ohara 3: Copyright 2000, 2002 Free Software Foundation, Inc.
1.1 maekawa 4:
5: This file is part of the GNU MP Library.
6:
7: The GNU MP Library is free software; you can redistribute it and/or modify
8: it under the terms of the GNU Lesser General Public License as published by
9: the Free Software Foundation; either version 2.1 of the License, or (at your
10: option) any later version.
11:
12: The GNU MP Library is distributed in the hope that it will be useful, but
13: WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
14: or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
15: License for more details.
16:
17: You should have received a copy of the GNU Lesser General Public License
18: along with the GNU MP Library; see the file COPYING.LIB. If not, write to
19: the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
20: MA 02111-1307, USA. */
21:
22: #include "gmp.h"
23: #include "gmp-impl.h"
24: #include "longlong.h"
25:
26: #if BITS_PER_MP_LIMB != 32 && BITS_PER_MP_LIMB != 64
1.1.1.2 ! ohara 27: choke me
1.1 maekawa 28: #endif
29:
30: void
31: mpq_set_d (mpq_ptr dest, double d)
32: {
33: int negative;
34: mp_exp_t exp;
1.1.1.2 ! ohara 35: mp_limb_t tp[LIMBS_PER_DOUBLE];
1.1 maekawa 36: mp_ptr np, dp;
37: mp_size_t nn, dn;
38: int c;
39:
40: negative = d < 0;
41: d = ABS (d);
42:
43: exp = __gmp_extract_double (tp, d);
44:
45: /* There are two main version of the conversion. The `then' arm handles
46: things that have a fractional part, while the `else' part handles
47: only integers. */
48: #if BITS_PER_MP_LIMB == 32
49: if (exp <= 1 || (exp == 2 && tp[0] != 0))
50: #else
51: if (exp <= 1)
52: #endif
53: {
54: if (d == 0.0)
55: {
56: SIZ(&(dest->_mp_num)) = 0;
57: SIZ(&(dest->_mp_den)) = 1;
58: PTR(&(dest->_mp_den))[0] = 1;
59: return;
60: }
61:
62: dn = -exp;
63: if (dest->_mp_num._mp_alloc < 3)
64: _mpz_realloc (&(dest->_mp_num), 3);
65: np = PTR(&(dest->_mp_num));
66: #if BITS_PER_MP_LIMB == 32
67: if ((tp[0] | tp[1]) == 0)
68: np[0] = tp[2], nn = 1;
69: else if (tp[0] == 0)
70: np[1] = tp[2], np[0] = tp[1], nn = 2;
71: else
72: np[2] = tp[2], np[1] = tp[1], np[0] = tp[0], nn = 3;
73: #else
74: if (tp[0] == 0)
75: np[0] = tp[1], nn = 1;
76: else
77: np[1] = tp[1], np[0] = tp[0], nn = 2;
78: #endif
79: dn += nn + 1;
80: if (dest->_mp_den._mp_alloc < dn)
81: _mpz_realloc (&(dest->_mp_den), dn);
82: dp = PTR(&(dest->_mp_den));
83: MPN_ZERO (dp, dn - 1);
84: dp[dn - 1] = 1;
85: count_trailing_zeros (c, np[0] | dp[0]);
86: if (c != 0)
87: {
88: mpn_rshift (np, np, nn, c);
89: nn -= np[nn - 1] == 0;
90: mpn_rshift (dp, dp, dn, c);
91: dn -= dp[dn - 1] == 0;
92: }
93: SIZ(&(dest->_mp_den)) = dn;
94: SIZ(&(dest->_mp_num)) = negative ? -nn : nn;
95: }
96: else
97: {
98: nn = exp;
99: if (dest->_mp_num._mp_alloc < nn)
100: _mpz_realloc (&(dest->_mp_num), nn);
101: np = PTR(&(dest->_mp_num));
102: #if BITS_PER_MP_LIMB == 32
103: switch (nn)
104: {
105: default:
106: MPN_ZERO (np, nn - 3);
107: np += nn - 3;
108: /* fall through */
109: case 3:
110: np[2] = tp[2], np[1] = tp[1], np[0] = tp[0];
111: break;
112: case 2:
113: np[1] = tp[2], np[0] = tp[1];
114: break;
115: }
116: #else
117: switch (nn)
118: {
119: default:
120: MPN_ZERO (np, nn - 2);
121: np += nn - 2;
122: /* fall through */
123: case 2:
124: np[1] = tp[1], np[0] = tp[0];
125: break;
126: }
127: #endif
128: dp = PTR(&(dest->_mp_den));
129: dp[0] = 1;
130: SIZ(&(dest->_mp_den)) = 1;
131: SIZ(&(dest->_mp_num)) = negative ? -nn : nn;
132: }
133: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>