version 1.1.1.1, 2000/09/09 14:12:16 |
version 1.1.1.2, 2003/08/25 16:06:00 |
|
|
/* __gmp_extract_double -- convert from double to array of mp_limb_t. |
/* __gmp_extract_double -- convert from double to array of mp_limb_t. |
|
|
Copyright (C) 1996, 1999, 2000 Free Software Foundation, Inc. |
Copyright 1996, 1999, 2000, 2001, 2002 Free Software Foundation, Inc. |
|
|
This file is part of the GNU MP Library. |
This file is part of the GNU MP Library. |
|
|
Line 30 MA 02111-1307, USA. */ |
|
Line 30 MA 02111-1307, USA. */ |
|
#define _GMP_IEEE_FLOATS 0 |
#define _GMP_IEEE_FLOATS 0 |
#endif |
#endif |
|
|
|
#define BITS_IN_MANTISSA 53 |
|
|
/* Extract a non-negative double in d. */ |
/* Extract a non-negative double in d. */ |
|
|
int |
int |
#if __STDC__ |
|
__gmp_extract_double (mp_ptr rp, double d) |
__gmp_extract_double (mp_ptr rp, double d) |
#else |
|
__gmp_extract_double (rp, d) |
|
mp_ptr rp; |
|
double d; |
|
#endif |
|
{ |
{ |
long exp; |
long exp; |
unsigned sc; |
unsigned sc; |
mp_limb_t manh, manl; |
#ifdef _LONG_LONG_LIMB |
|
#define BITS_PER_PART 64 /* somewhat bogus */ |
|
unsigned long long int manl; |
|
#else |
|
#define BITS_PER_PART GMP_LIMB_BITS |
|
unsigned long int manh, manl; |
|
#endif |
|
|
/* BUGS |
/* BUGS |
|
|
1. Should handle Inf and NaN in IEEE specific code. |
1. Should handle Inf and NaN in IEEE specific code. |
2. Handle Inf and NaN also in default code, to avoid hangs. |
2. Handle Inf and NaN also in default code, to avoid hangs. |
3. Generalize to handle all BITS_PER_MP_LIMB >= 32. |
3. Generalize to handle all GMP_LIMB_BITS >= 32. |
4. This lits is incomplete and misspelled. |
4. This lits is incomplete and misspelled. |
*/ |
*/ |
|
|
|
ASSERT (d >= 0.0); |
|
|
if (d == 0.0) |
if (d == 0.0) |
{ |
{ |
rp[0] = 0; |
MPN_ZERO (rp, LIMBS_PER_DOUBLE); |
rp[1] = 0; |
|
#if BITS_PER_MP_LIMB == 32 |
|
rp[2] = 0; |
|
#endif |
|
return 0; |
return 0; |
} |
} |
|
|
Line 72 __gmp_extract_double (rp, d) |
|
Line 72 __gmp_extract_double (rp, d) |
|
union ieee_double_extract x; |
union ieee_double_extract x; |
x.d = d; |
x.d = d; |
exp = x.s.exp; |
exp = x.s.exp; |
#if BITS_PER_MP_LIMB == 64 |
#if BITS_PER_PART == 64 /* generalize this to BITS_PER_PART > BITS_IN_MANTISSA */ |
manl = (((mp_limb_t) 1 << 63) |
manl = (((mp_limb_t) 1 << 63) |
| ((mp_limb_t) x.s.manh << 43) | ((mp_limb_t) x.s.manl << 11)); |
| ((mp_limb_t) x.s.manh << 43) | ((mp_limb_t) x.s.manl << 11)); |
if (exp == 0) |
if (exp == 0) |
Line 87 __gmp_extract_double (rp, d) |
|
Line 87 __gmp_extract_double (rp, d) |
|
} |
} |
while ((mp_limb_signed_t) manl >= 0); |
while ((mp_limb_signed_t) manl >= 0); |
} |
} |
#else |
#endif |
|
#if BITS_PER_PART == 32 |
manh = ((mp_limb_t) 1 << 31) | (x.s.manh << 11) | (x.s.manl >> 21); |
manh = ((mp_limb_t) 1 << 31) | (x.s.manh << 11) | (x.s.manl >> 21); |
manl = x.s.manl << 11; |
manl = x.s.manl << 11; |
if (exp == 0) |
if (exp == 0) |
Line 104 __gmp_extract_double (rp, d) |
|
Line 105 __gmp_extract_double (rp, d) |
|
while ((mp_limb_signed_t) manh >= 0); |
while ((mp_limb_signed_t) manh >= 0); |
} |
} |
#endif |
#endif |
|
#if BITS_PER_PART != 32 && BITS_PER_PART != 64 |
|
You need to generalize the code above to handle this. |
|
#endif |
exp -= 1022; /* Remove IEEE bias. */ |
exp -= 1022; /* Remove IEEE bias. */ |
} |
} |
#else |
#else |
Line 140 __gmp_extract_double (rp, d) |
|
Line 144 __gmp_extract_double (rp, d) |
|
} |
} |
} |
} |
|
|
d *= MP_BASE_AS_DOUBLE; |
d *= (4.0 * ((unsigned long int) 1 << (BITS_PER_PART - 2))); |
#if BITS_PER_MP_LIMB == 64 |
#if BITS_PER_PART == 64 |
manl = d; |
manl = d; |
#else |
#else |
manh = d; |
manh = d; |
manl = (d - manh) * MP_BASE_AS_DOUBLE; |
manl = (d - manh) * (4.0 * ((unsigned long int) 1 << (BITS_PER_PART - 2))); |
#endif |
#endif |
} |
} |
#endif |
#endif /* IEEE */ |
|
|
sc = (unsigned) exp % BITS_PER_MP_LIMB; |
/* Up until here, we have ignored the actual limb size. Remains |
|
to split manh,,manl into an array of LIMBS_PER_DOUBLE limbs. |
|
*/ |
|
|
|
sc = (unsigned) (exp + 64 * GMP_NUMB_BITS) % GMP_NUMB_BITS; |
|
|
/* We add something here to get rounding right. */ |
/* We add something here to get rounding right. */ |
exp = (exp + 2048) / BITS_PER_MP_LIMB - 2048 / BITS_PER_MP_LIMB + 1; |
exp = (exp + 64 * GMP_NUMB_BITS) / GMP_NUMB_BITS - 64 * GMP_NUMB_BITS / GMP_NUMB_BITS + 1; |
|
|
#if BITS_PER_MP_LIMB == 64 |
#if LIMBS_PER_DOUBLE == 2 |
|
#if GMP_NAIL_BITS == 0 |
if (sc != 0) |
if (sc != 0) |
{ |
{ |
rp[1] = manl >> (BITS_PER_MP_LIMB - sc); |
rp[1] = manl >> (GMP_LIMB_BITS - sc); |
rp[0] = manl << sc; |
rp[0] = manl << sc; |
} |
} |
else |
else |
Line 168 __gmp_extract_double (rp, d) |
|
Line 177 __gmp_extract_double (rp, d) |
|
exp--; |
exp--; |
} |
} |
#else |
#else |
|
if (sc > GMP_NAIL_BITS) |
|
{ |
|
rp[1] = manl >> (GMP_LIMB_BITS - sc); |
|
rp[0] = (manl << (sc - GMP_NAIL_BITS)) & GMP_NUMB_MASK; |
|
} |
|
else |
|
{ |
|
if (sc == 0) |
|
{ |
|
rp[1] = manl >> GMP_NAIL_BITS; |
|
rp[0] = (manl << GMP_NUMB_BITS - GMP_NAIL_BITS) & GMP_NUMB_MASK; |
|
exp--; |
|
} |
|
else |
|
{ |
|
rp[1] = manl >> (GMP_LIMB_BITS - sc); |
|
rp[0] = (manl >> (GMP_NAIL_BITS - sc)) & GMP_NUMB_MASK; |
|
} |
|
} |
|
#endif |
|
#endif |
|
|
|
#if LIMBS_PER_DOUBLE == 3 |
|
#if GMP_NAIL_BITS == 0 |
if (sc != 0) |
if (sc != 0) |
{ |
{ |
rp[2] = manh >> (BITS_PER_MP_LIMB - sc); |
rp[2] = manh >> (GMP_LIMB_BITS - sc); |
rp[1] = (manl >> (BITS_PER_MP_LIMB - sc)) | (manh << sc); |
rp[1] = (manh << sc) | (manl >> (GMP_LIMB_BITS - sc)); |
rp[0] = manl << sc; |
rp[0] = manl << sc; |
} |
} |
else |
else |
Line 181 __gmp_extract_double (rp, d) |
|
Line 214 __gmp_extract_double (rp, d) |
|
rp[0] = 0; |
rp[0] = 0; |
exp--; |
exp--; |
} |
} |
|
#else |
|
if (sc > GMP_NAIL_BITS) |
|
{ |
|
rp[2] = (manh >> (GMP_LIMB_BITS - sc)); |
|
rp[1] = ((manh << (sc - GMP_NAIL_BITS)) | |
|
(manl >> (GMP_LIMB_BITS - sc + GMP_NAIL_BITS))) & GMP_NUMB_MASK; |
|
if (sc >= 2 * GMP_NAIL_BITS) |
|
rp[0] = (manl << sc - 2 * GMP_NAIL_BITS) & GMP_NUMB_MASK; |
|
else |
|
rp[0] = manl >> (2 * GMP_NAIL_BITS - sc) & GMP_NUMB_MASK; |
|
} |
|
else |
|
{ |
|
if (sc == 0) |
|
{ |
|
rp[2] = manh >> GMP_NAIL_BITS; |
|
rp[1] = ((manh << GMP_NUMB_BITS - GMP_NAIL_BITS) | (manl >> 2 * GMP_NAIL_BITS)) & GMP_NUMB_MASK; |
|
rp[0] = 0; |
|
exp--; |
|
} |
|
else |
|
{ |
|
rp[2] = (manh >> (GMP_LIMB_BITS - sc)); |
|
rp[1] = (manh >> (GMP_NAIL_BITS - sc)) & GMP_NUMB_MASK; |
|
rp[0] = ((manh << (GMP_NUMB_BITS - GMP_NAIL_BITS + sc)) |
|
| (manl >> (GMP_LIMB_BITS - (GMP_NUMB_BITS - GMP_NAIL_BITS + sc)))) & GMP_NUMB_MASK; |
|
} |
|
} |
|
#endif |
|
#endif |
|
|
|
#if LIMBS_PER_DOUBLE > 3 |
|
/* Insert code for splitting manh,,manl into LIMBS_PER_DOUBLE |
|
mp_limb_t's at rp. */ |
|
if (sc != 0) |
|
{ |
|
/* This is not perfect, and would fail for GMP_LIMB_BITS == 16. |
|
The ASSERT_ALWAYS should catch the problematic cases. */ |
|
ASSERT_ALWAYS ((manl << sc) == 0); |
|
manl = (manh << sc) | (manl >> (GMP_LIMB_BITS - sc)); |
|
manh = manh >> (GMP_LIMB_BITS - sc); |
|
} |
|
{ |
|
int i; |
|
for (i = LIMBS_PER_DOUBLE - 1; i >= 0; i--) |
|
{ |
|
rp[i] = manh >> (BITS_PER_LONGINT - GMP_LIMB_BITS); |
|
manh = ((manh << GMP_LIMB_BITS) |
|
| (manl >> (BITS_PER_LONGINT - GMP_LIMB_BITS))); |
|
manl = manl << GMP_LIMB_BITS; |
|
} |
|
} |
#endif |
#endif |
|
|
return exp; |
return exp; |