added new functions mpfr_set_binary32 and mpfr_get_binary32

fixed bug in mpfr_get_d and mpfr_get_decimal64 for RNDA


git-svn-id: svn://scm.gforge.inria.fr/svn/mpfr/trunk@6424 280ebfd0-de03-0410-8827-d642c229c3f4

11 files changed, 331 insertions, 132 deletions

diff --git a/Makefile.am b/Makefile.am
index 31382b1e6..744cb6dbe 100644
--- a/Makefile.am
+++ b/Makefile.am
@@ -17,13 +17,13 @@ SUBDIRS = tests
 nobase_dist_doc_DATA = AUTHORS BUGS COPYING COPYING.LESSER FAQ.html NEWS TODO \
   examples/ReadMe examples/divworst.c examples/rndo-add.c examples/sample.c
 
-EXTRA_DIST = PATCHES VERSION get_patches.sh round_raw_generic.c gen_inverse.h jyn_asympt.c
+EXTRA_DIST = PATCHES VERSION get_patches.sh round_raw_generic.c gen_inverse.h jyn_asympt.c scale2.c ieee_floats.h
 
 include_HEADERS = mpfr.h mpf2mpfr.h
 
 lib_LTLIBRARIES = libmpfr.la
 
-libmpfr_la_SOURCES = mpfr.h mpf2mpfr.h mpfr-gmp.h mpfr-impl.h mpfr-longlong.h mpfr-thread.h exceptions.c extract.c uceil_exp2.c uceil_log2.c ufloor_log2.c add.c add1.c add_ui.c agm.c clear.c cmp.c cmp_abs.c cmp_si.c cmp_ui.c comparisons.c div_2exp.c div_2si.c div_2ui.c div.c div_ui.c dump.c eq.c exp10.c exp2.c exp3.c exp.c frac.c get_d.c get_exp.c get_str.c init.c inp_str.c isinteger.c isinf.c isnan.c isnum.c const_log2.c log.c modf.c mul_2exp.c mul_2si.c mul_2ui.c mul.c mul_ui.c neg.c next.c out_str.c printf.c vasprintf.c const_pi.c pow.c pow_si.c pow_ui.c print_raw.c print_rnd_mode.c reldiff.c round_prec.c set.c setmax.c setmin.c set_d.c set_dfl_prec.c set_exp.c set_rnd.c set_f.c set_prc_raw.c set_prec.c set_q.c set_si.c set_str.c set_str_raw.c set_ui.c set_z.c sqrt.c sqrt_ui.c sub.c sub1.c sub_ui.c rint.c ui_div.c ui_sub.c urandomb.c get_z_exp.c swap.c  factorial.c cosh.c sinh.c tanh.c sinh_cosh.c acosh.c asinh.c atanh.c atan.c cmp2.c exp_2.c asin.c const_euler.c cos.c sin.c tan.c fma.c fms.c hypot.c log1p.c expm1.c log2.c log10.c ui_pow.c ui_pow_ui.c minmax.c dim.c signbit.c copysign.c setsign.c gmp_op.c init2.c acos.c sin_cos.c set_nan.c set_inf.c set_zero.c powerof2.c gamma.c set_ld.c get_ld.c cbrt.c volatile.c fits_s.h fits_sshort.c fits_sint.c fits_slong.c fits_u.h fits_ushort.c fits_uint.c fits_ulong.c fits_uintmax.c fits_intmax.c get_si.c get_ui.c zeta.c cmp_d.c erf.c inits.c inits2.c clears.c sgn.c check.c sub1sp.c version.c mpn_exp.c mpfr-gmp.c mp_clz_tab.c sum.c add1sp.c free_cache.c si_op.c cmp_ld.c set_ui_2exp.c set_si_2exp.c set_uj.c set_sj.c get_sj.c get_uj.c get_z.c iszero.c cache.c sqr.c int_ceil_log2.c isqrt.c strtofr.c pow_z.c logging.c mulders.c get_f.c round_p.c erfc.c atan2.c subnormal.c const_catalan.c root.c gen_inverse.h sec.c csc.c cot.c eint.c sech.c csch.c coth.c round_near_x.c constant.c abort_prec_max.c stack_interface.c lngamma.c zeta_ui.c set_d64.c get_d64.c jn.c yn.c rem1.c get_patches.c add_d.c sub_d.c d_sub.c mul_d.c div_d.c d_div.c li2.c rec_sqrt.c min_prec.c buildopt.c digamma.c bernoulli.c isregular.c
+libmpfr_la_SOURCES = mpfr.h mpf2mpfr.h mpfr-gmp.h mpfr-impl.h mpfr-longlong.h mpfr-thread.h exceptions.c extract.c uceil_exp2.c uceil_log2.c ufloor_log2.c add.c add1.c add_ui.c agm.c clear.c cmp.c cmp_abs.c cmp_si.c cmp_ui.c comparisons.c div_2exp.c div_2si.c div_2ui.c div.c div_ui.c dump.c eq.c exp10.c exp2.c exp3.c exp.c frac.c get_d.c get_exp.c get_str.c init.c inp_str.c isinteger.c isinf.c isnan.c isnum.c const_log2.c log.c modf.c mul_2exp.c mul_2si.c mul_2ui.c mul.c mul_ui.c neg.c next.c out_str.c printf.c vasprintf.c const_pi.c pow.c pow_si.c pow_ui.c print_raw.c print_rnd_mode.c reldiff.c round_prec.c set.c setmax.c setmin.c set_d.c set_dfl_prec.c set_exp.c set_rnd.c set_f.c set_prc_raw.c set_prec.c set_q.c set_si.c set_str.c set_str_raw.c set_ui.c set_z.c sqrt.c sqrt_ui.c sub.c sub1.c sub_ui.c rint.c ui_div.c ui_sub.c urandomb.c get_z_exp.c swap.c  factorial.c cosh.c sinh.c tanh.c sinh_cosh.c acosh.c asinh.c atanh.c atan.c cmp2.c exp_2.c asin.c const_euler.c cos.c sin.c tan.c fma.c fms.c hypot.c log1p.c expm1.c log2.c log10.c ui_pow.c ui_pow_ui.c minmax.c dim.c signbit.c copysign.c setsign.c gmp_op.c init2.c acos.c sin_cos.c set_nan.c set_inf.c set_zero.c powerof2.c gamma.c set_ld.c get_ld.c cbrt.c volatile.c fits_s.h fits_sshort.c fits_sint.c fits_slong.c fits_u.h fits_ushort.c fits_uint.c fits_ulong.c fits_uintmax.c fits_intmax.c get_si.c get_ui.c zeta.c cmp_d.c erf.c inits.c inits2.c clears.c sgn.c check.c sub1sp.c version.c mpn_exp.c mpfr-gmp.c mp_clz_tab.c sum.c add1sp.c free_cache.c si_op.c cmp_ld.c set_ui_2exp.c set_si_2exp.c set_uj.c set_sj.c get_sj.c get_uj.c get_z.c iszero.c cache.c sqr.c int_ceil_log2.c isqrt.c strtofr.c pow_z.c logging.c mulders.c get_f.c round_p.c erfc.c atan2.c subnormal.c const_catalan.c root.c gen_inverse.h sec.c csc.c cot.c eint.c sech.c csch.c coth.c round_near_x.c constant.c abort_prec_max.c stack_interface.c lngamma.c zeta_ui.c set_d64.c get_d64.c jn.c yn.c rem1.c get_patches.c add_d.c sub_d.c d_sub.c mul_d.c div_d.c d_div.c li2.c rec_sqrt.c min_prec.c buildopt.c digamma.c bernoulli.c isregular.c set_binary32.c get_binary32.c
 
 libmpfr_la_LIBADD = @LIBOBJS@
 
diff --git a/get_d.c b/get_d.c
index 99b7f92af..943f96fdc 100644
--- a/get_d.c
+++ b/get_d.c
@@ -26,125 +26,8 @@ http://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc.,
 #define MPFR_NEED_LONGLONG_H
 #include "mpfr-impl.h"
 
-/* "double" NaN and infinities are written as explicit bytes to be sure of
-   getting what we want, and to be sure of not depending on libm.
-
-   Could use 4-byte "float" values and let the code convert them, but it
-   seems more direct to give exactly what we want.  Certainly for gcc 3.0.2
-   on alphaev56-unknown-freebsd4.3 the NaN must be 8-bytes, since that
-   compiler+system was seen incorrectly converting from a "float" NaN.  */
-
-#if _GMP_IEEE_FLOATS
-
-/* The "d" field guarantees alignment to a suitable boundary for a double.
-   Could use a union instead, if we checked the compiler supports union
-   initializers.  */
-struct dbl_bytes {
-  unsigned char b[8];
-  double d;
-};
-
-#define MPFR_DBL_INFP  (* (const double *) dbl_infp.b)
-#define MPFR_DBL_INFM  (* (const double *) dbl_infm.b)
-#define MPFR_DBL_NAN   (* (const double *) dbl_nan.b)
-
-#if HAVE_DOUBLE_IEEE_LITTLE_ENDIAN
-static const struct dbl_bytes dbl_infp =
-  { { 0, 0, 0, 0, 0, 0, 0xF0, 0x7F }, 0.0 };
-static const struct dbl_bytes dbl_infm =
-  { { 0, 0, 0, 0, 0, 0, 0xF0, 0xFF }, 0.0 };
-static const struct dbl_bytes dbl_nan  =
-  { { 0, 0, 0, 0, 0, 0, 0xF8, 0x7F }, 0.0 };
-#endif
-#if HAVE_DOUBLE_IEEE_LITTLE_SWAPPED
-static const struct dbl_bytes dbl_infp =
-  { { 0, 0, 0xF0, 0x7F, 0, 0, 0, 0 }, 0.0 };
-static const struct dbl_bytes dbl_infm =
-  { { 0, 0, 0xF0, 0xFF, 0, 0, 0, 0 }, 0.0 };
-static const struct dbl_bytes dbl_nan  =
-  { { 0, 0, 0xF8, 0x7F, 0, 0, 0, 0 }, 0.0 };
-#endif
-#if HAVE_DOUBLE_IEEE_BIG_ENDIAN
-static const struct dbl_bytes dbl_infp =
-  { { 0x7F, 0xF0, 0, 0, 0, 0, 0, 0 }, 0.0 };
-static const struct dbl_bytes dbl_infm =
-  { { 0xFF, 0xF0, 0, 0, 0, 0, 0, 0 }, 0.0 };
-static const struct dbl_bytes dbl_nan  =
-  { { 0x7F, 0xF8, 0, 0, 0, 0, 0, 0 }, 0.0 };
-#endif
-
-#else /* _GMP_IEEE_FLOATS */
-
-#define MPFR_DBL_INFP DBL_POS_INF
-#define MPFR_DBL_INFM DBL_NEG_INF
-#define MPFR_DBL_NAN DBL_NAN
-
-#endif /* _GMP_IEEE_FLOATS */
-
-
-/* multiplies 1/2 <= d <= 1 by 2^exp */
-static double
-mpfr_scale2 (double d, int exp)
-{
-#if _GMP_IEEE_FLOATS
-  {
-    union ieee_double_extract x;
-
-    if (MPFR_UNLIKELY (d == 1.0))
-      {
-        d = 0.5;
-        exp ++;
-      }
-
-    /* now 1/2 <= d < 1 */
-
-    /* infinities and zeroes have already been checked */
-    MPFR_ASSERTD (-1073 <= exp && exp <= 1025);
-
-    x.d = d;
-    if (MPFR_UNLIKELY (exp < -1021)) /* subnormal case */
-      {
-        x.s.exp += exp + 52;
-        x.d *= DBL_EPSILON;
-      }
-    else /* normalized case */
-      {
-        x.s.exp += exp;
-      }
-    return x.d;
-  }
-#else /* _GMP_IEEE_FLOATS */
-  {
-    double factor;
-
-    /* An overflow may occurs (example: 0.5*2^1024) */
-    if (d < 1.0)
-      {
-        d += d;
-        exp--;
-      }
-    /* Now 1.0 <= d < 2.0 */
-
-    if (exp < 0)
-      {
-        factor = 0.5;
-        exp = -exp;
-      }
-    else
-      {
-        factor = 2.0;
-      }
-    while (exp != 0)
-      {
-        if ((exp & 1) != 0)
-          d *= factor;
-        exp >>= 1;
-        factor *= factor;
-      }
-    return d;
-  }
-#endif
-}
+#include "ieee_floats.h"
+#include "scale2.c"
 
 /* Assumes IEEE-754 double precision; otherwise, only an approximated
    result will be returned, without any guaranty (and special cases
@@ -174,6 +57,9 @@ mpfr_get_d (mpfr_srcptr src, mpfr_rnd_t rnd_mode)
   e = MPFR_GET_EXP (src);
   negative = MPFR_IS_NEG (src);
 
+  if (MPFR_UNLIKELY(rnd_mode == MPFR_RNDA))
+    rnd_mode = negative ? MPFR_RNDD : MPFR_RNDU;
+
   /* the smallest normalized number is 2^(-1022)=0.1e-1021, and the smallest
      subnormal is 2^(-1074)=0.1e-1073 */
   if (MPFR_UNLIKELY (e < -1073))
@@ -188,7 +74,8 @@ mpfr_get_d (mpfr_srcptr src, mpfr_rnd_t rnd_mode)
         (rnd_mode == MPFR_RNDU ||
          (rnd_mode == MPFR_RNDN && mpfr_cmp_si_2exp(src, 1, -1075) > 0)
          ? DBL_MIN : 0.0);
-      if (d != 0.0)
+      if (d != 0.0) /* we multiply DBL_MIN = 2^(-1022) by DBL_EPSILON = 2^(-52)
+                       to get +-2^(-1074) */
         d *= DBL_EPSILON;
     }
   /* the largest normalized number is 2^1024*(1-2^(-53))=0.111...111e1024 */
diff --git a/get_d64.c b/get_d64.c
index 331ca1d8f..11f5ab32d 100644
--- a/get_d64.c
+++ b/get_d64.c
@@ -317,6 +317,9 @@ mpfr_get_decimal64 (mpfr_srcptr src, mpfr_rnd_t rnd_mode)
   e = MPFR_GET_EXP (src);
   negative = MPFR_IS_NEG (src);
 
+  if (MPFR_UNLIKELY(rnd_mode == MPFR_RNDA))
+    rnd_mode = negative ? MPFR_RNDD : MPFR_RNDU;
+
   /* the smallest decimal64 number is 10^(-398),
      with 2^(-1323) < 10^(-398) < 2^(-1322) */
   if (MPFR_UNLIKELY (e < -1323)) /* src <= 2^(-1324) < 1/2*10^(-398) */
diff --git a/ieee_floats.h b/ieee_floats.h
new file mode 100644
index 000000000..08400c522
--- /dev/null
+++ b/ieee_floats.h
@@ -0,0 +1,76 @@
+/* auxiliary data to generate special IEEE floats (NaN, +Inf, -Inf)
+
+Copyright 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 Free Software Foundation, Inc.
+Contributed by the Arenaire and Cacao projects, INRIA.
+
+This file is part of the GNU MPFR Library.
+
+The GNU MPFR Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Lesser General Public License as published by
+the Free Software Foundation; either version 3 of the License, or (at your
+option) any later version.
+
+The GNU MPFR Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
+License for more details.
+
+You should have received a copy of the GNU Lesser General Public License
+along with the GNU MPFR Library; see the file COPYING.LESSER.  If not, see
+http://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc.,
+51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */
+
+/* "double" NaN and infinities are written as explicit bytes to be sure of
+   getting what we want, and to be sure of not depending on libm.
+
+   Could use 4-byte "float" values and let the code convert them, but it
+   seems more direct to give exactly what we want.  Certainly for gcc 3.0.2
+   on alphaev56-unknown-freebsd4.3 the NaN must be 8-bytes, since that
+   compiler+system was seen incorrectly converting from a "float" NaN.  */
+
+#if _GMP_IEEE_FLOATS
+
+/* The "d" field guarantees alignment to a suitable boundary for a double.
+   Could use a union instead, if we checked the compiler supports union
+   initializers.  */
+struct dbl_bytes {
+  unsigned char b[8];
+  double d;
+};
+
+#define MPFR_DBL_INFP  (* (const double *) dbl_infp.b)
+#define MPFR_DBL_INFM  (* (const double *) dbl_infm.b)
+#define MPFR_DBL_NAN   (* (const double *) dbl_nan.b)
+
+#if HAVE_DOUBLE_IEEE_LITTLE_ENDIAN
+static const struct dbl_bytes dbl_infp =
+  { { 0, 0, 0, 0, 0, 0, 0xF0, 0x7F }, 0.0 };
+static const struct dbl_bytes dbl_infm =
+  { { 0, 0, 0, 0, 0, 0, 0xF0, 0xFF }, 0.0 };
+static const struct dbl_bytes dbl_nan  =
+  { { 0, 0, 0, 0, 0, 0, 0xF8, 0x7F }, 0.0 };
+#endif
+#if HAVE_DOUBLE_IEEE_LITTLE_SWAPPED
+static const struct dbl_bytes dbl_infp =
+  { { 0, 0, 0xF0, 0x7F, 0, 0, 0, 0 }, 0.0 };
+static const struct dbl_bytes dbl_infm =
+  { { 0, 0, 0xF0, 0xFF, 0, 0, 0, 0 }, 0.0 };
+static const struct dbl_bytes dbl_nan  =
+  { { 0, 0, 0xF8, 0x7F, 0, 0, 0, 0 }, 0.0 };
+#endif
+#if HAVE_DOUBLE_IEEE_BIG_ENDIAN
+static const struct dbl_bytes dbl_infp =
+  { { 0x7F, 0xF0, 0, 0, 0, 0, 0, 0 }, 0.0 };
+static const struct dbl_bytes dbl_infm =
+  { { 0xFF, 0xF0, 0, 0, 0, 0, 0, 0 }, 0.0 };
+static const struct dbl_bytes dbl_nan  =
+  { { 0x7F, 0xF8, 0, 0, 0, 0, 0, 0 }, 0.0 };
+#endif
+
+#else /* _GMP_IEEE_FLOATS */
+
+#define MPFR_DBL_INFP DBL_POS_INF
+#define MPFR_DBL_INFM DBL_NEG_INF
+#define MPFR_DBL_NAN DBL_NAN
+
+#endif /* _GMP_IEEE_FLOATS */
diff --git a/mpfr-impl.h b/mpfr-impl.h
index 46573231b..e1636b7c7 100644
--- a/mpfr-impl.h
+++ b/mpfr-impl.h
@@ -382,6 +382,11 @@ __MPFR_DECLSPEC extern const mpfr_t __gmpfr_four;
 #endif
 #define MPFR_LIMBS_PER_DOUBLE ((IEEE_DBL_MANT_DIG-1)/BITS_PER_MP_LIMB+1)
 
+#ifndef IEEE_FLT_MANT_DIG
+#define IEEE_FLT_MANT_DIG 24
+#endif
+#define MPFR_LIMBS_PER_BINARY32 ((IEEE_FLT_MANT_DIG-1)/BITS_PER_MP_LIMB+1)
+
 /* Visual C++ doesn't support +1.0/.00, -1.0/0.0 and 0.0/0.0
    at compile time. */
 #if defined(_MSC_VER) && defined(_WIN32) && (_MSC_VER >= 1200)
diff --git a/mpfr.h b/mpfr.h
index 5cbb47d86..0de810bb4 100644
--- a/mpfr.h
+++ b/mpfr.h
@@ -281,6 +281,7 @@ __MPFR_DECLSPEC void mpfr_set_default_prec _MPFR_PROTO((mpfr_prec_t));
 __MPFR_DECLSPEC mp_prec_t mpfr_get_default_prec _MPFR_PROTO((void));
 
 __MPFR_DECLSPEC int mpfr_set_d _MPFR_PROTO ((mpfr_ptr, double, mpfr_rnd_t));
+__MPFR_DECLSPEC int mpfr_set_binary32 _MPFR_PROTO ((mpfr_ptr, float, mpfr_rnd_t));
 #ifdef MPFR_WANT_DECIMAL_FLOATS
 __MPFR_DECLSPEC int mpfr_set_decimal64 _MPFR_PROTO ((mpfr_ptr, _Decimal64,
                                                      mpfr_rnd_t));
@@ -341,6 +342,7 @@ __MPFR_DECLSPEC uintmax_t mpfr_get_uj _MPFR_PROTO ((mpfr_srcptr, mpfr_rnd_t));
 #endif
 
 __MPFR_DECLSPEC mp_exp_t mpfr_get_z_exp _MPFR_PROTO ((mpz_ptr, mpfr_srcptr));
+__MPFR_DECLSPEC float mpfr_get_binary32 _MPFR_PROTO ((mpfr_srcptr, mpfr_rnd_t));
 __MPFR_DECLSPEC double mpfr_get_d _MPFR_PROTO ((mpfr_srcptr, mpfr_rnd_t));
 #ifdef MPFR_WANT_DECIMAL_FLOATS
 __MPFR_DECLSPEC _Decimal64 mpfr_get_decimal64 _MPFR_PROTO ((mpfr_srcptr,
diff --git a/mpfr.texi b/mpfr.texi
index 42accea66..734c0d59f 100644
--- a/mpfr.texi
+++ b/mpfr.texi
@@ -1068,6 +1068,7 @@ These functions assign new values to already initialized floats
 @deftypefunx int mpfr_set_si (mpfr_t @var{rop}, long int @var{op}, mpfr_rnd_t @var{rnd})
 @deftypefunx int mpfr_set_uj (mpfr_t @var{rop}, uintmax_t @var{op}, mpfr_rnd_t @var{rnd})
 @deftypefunx int mpfr_set_sj (mpfr_t @var{rop}, intmax_t @var{op}, mpfr_rnd_t @var{rnd})
+@deftypefunx int mpfr_set_binary32 (mpfr_t @var{rop}, float @var{op}, mpfr_rnd_t @var{rnd})
 @deftypefunx int mpfr_set_d (mpfr_t @var{rop}, double @var{op}, mpfr_rnd_t @var{rnd})
 @deftypefunx int mpfr_set_ld (mpfr_t @var{rop}, long double @var{op}, mpfr_rnd_t @var{rnd})
 @deftypefunx int mpfr_set_decimal64 (mpfr_t @var{rop}, _Decimal64 @var{op}, mpfr_rnd_t @var{rnd})
@@ -1080,8 +1081,8 @@ Note that the input 0 is converted to +0 by @code{mpfr_set_ui},
 @code{mpfr_set_si}, @code{mpfr_set_sj}, @code{mpfr_set_uj},
 @code{mpfr_set_z}, @code{mpfr_set_q} and
 @code{mpfr_set_f}, regardless of the rounding mode.
-If the system does not support the IEEE 754 standard, @code{mpfr_set_d},
-@code{mpfr_set_ld} and
+If the system does not support the IEEE 754 standard,
+@code{mpfr_set_binary32}, @code{mpfr_set_d}, @code{mpfr_set_ld} and
 @code{mpfr_set_decimal64} might not preserve the signed zeros.
 The @code{mpfr_set_decimal64} function is built only with the configure
 option @samp{--enable-decimal-float}, which also requires
@@ -1094,7 +1095,8 @@ denominator) can not be representable as a @code{mpfr_t}.
 
 Note: If you want to store a floating-point constant to a @code{mpfr_t},
 you should use @code{mpfr_set_str} (or one of the MPFR constant functions,
-such as @code{mpfr_const_pi} for @m{\pi,Pi}) instead of @code{mpfr_set_d},
+such as @code{mpfr_const_pi} for @m{\pi,Pi}) instead of 
+@code{mpfr_set_binary32}, @code{mpfr_set_d},
 @code{mpfr_set_ld} or @code{mpfr_set_decimal64}.
 Otherwise the floating-point constant will be first
 converted into a reduced-precision (e.g., 53-bit) binary number before
@@ -1235,11 +1237,12 @@ See @code{mpfr_set_str}.
 @cindex Conversion functions
 @section Conversion Functions
 
-@deftypefun double mpfr_get_d (mpfr_t @var{op}, mpfr_rnd_t @var{rnd})
+@deftypefun float mpfr_get_binary32 (mpfr_t @var{op}, mpfr_rnd_t @var{rnd})
+@deftypefunx double mpfr_get_d (mpfr_t @var{op}, mpfr_rnd_t @var{rnd})
 @deftypefunx {long double} mpfr_get_ld (mpfr_t @var{op}, mpfr_rnd_t @var{rnd})
 @deftypefunx _Decimal64 mpfr_get_decimal64 (mpfr_t @var{op}, mpfr_rnd_t @var{rnd})
-Convert @var{op} to a @code{double} (respectively @code{_Decimal64} or
-@code{long double}), using the rounding mode @var{rnd}.
+Convert @var{op} to a @code{float} (respectively @code{double},
+@code{long double} or @code{_Decimal64}), using the rounding mode @var{rnd}.
 If @var{op} is NaN, some fixed NaN (either quiet or signaling) or the result
 of 0.0/0.0 is returned. If @var{op} is @pom{}Inf, an infinity of the same
 sign or the result of @pom{}1.0/0.0 is returned. If @var{op} is zero, these
diff --git a/scale2.c b/scale2.c
new file mode 100644
index 000000000..19c0ae7c8
--- /dev/null
+++ b/scale2.c
@@ -0,0 +1,88 @@
+/* mpfr_scale2 -- multiply a double float by 2^exp
+
+Copyright 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 Free Software Foundation, Inc.
+Contributed by the Arenaire and Cacao projects, INRIA.
+
+This file is part of the GNU MPFR Library.
+
+The GNU MPFR Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Lesser General Public License as published by
+the Free Software Foundation; either version 3 of the License, or (at your
+option) any later version.
+
+The GNU MPFR Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
+License for more details.
+
+You should have received a copy of the GNU Lesser General Public License
+along with the GNU MPFR Library; see the file COPYING.LESSER.  If not, see
+http://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc.,
+51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */
+
+/* Note: we could use the ldexp function, but since we want not to depend on
+   math.h, we write our own implementation. */
+
+/* multiplies 1/2 <= d <= 1 by 2^exp */
+static double
+mpfr_scale2 (double d, int exp)
+{
+#if _GMP_IEEE_FLOATS
+  {
+    union ieee_double_extract x;
+
+    if (MPFR_UNLIKELY (d == 1.0))
+      {
+        d = 0.5;
+        exp ++;
+      }
+
+    /* now 1/2 <= d < 1 */
+
+    /* infinities and zeroes have already been checked */
+    MPFR_ASSERTD (-1073 <= exp && exp <= 1025);
+
+    x.d = d;
+    if (MPFR_UNLIKELY (exp < -1021)) /* subnormal case */
+      {
+        x.s.exp += exp + 52;
+        x.d *= DBL_EPSILON;
+      }
+    else /* normalized case */
+      {
+        x.s.exp += exp;
+      }
+    return x.d;
+  }
+#else /* _GMP_IEEE_FLOATS */
+  {
+    double factor;
+
+    /* An overflow may occurs (example: 0.5*2^1024) */
+    if (d < 1.0)
+      {
+        d += d;
+        exp--;
+      }
+    /* Now 1.0 <= d < 2.0 */
+
+    if (exp < 0)
+      {
+        factor = 0.5;
+        exp = -exp;
+      }
+    else
+      {
+        factor = 2.0;
+      }
+    while (exp != 0)
+      {
+        if ((exp & 1) != 0)
+          d *= factor;
+        exp >>= 1;
+        factor *= factor;
+      }
+    return d;
+  }
+#endif
+}
diff --git a/tests/Makefile.am b/tests/Makefile.am
index 212496965..dcf19e6d8 100644
--- a/tests/Makefile.am
+++ b/tests/Makefile.am
@@ -20,10 +20,10 @@ check_PROGRAMS = tversion tinternals tinits tisqrt tsgn		\
      tconst_log2 tconst_pi tcopysign tcos tcosh tcot tcoth	\
      tcsc tcsch td_div td_sub tdigamma tdim tdiv tdiv_d tdiv_ui \
      teint teq terf texp texp10 texp2 texpm1 tfactorial tfits   \
-     tfma tfmod tfms tfprintf tfrac tgamma tget_d tget_d_2exp	\
-     tget_f tget_ld_2exp tget_set_d64 tget_sj tget_str		\
-     tget_z tgmpop thyperbolic thypot tinp_str tj0 tj1 tjn	\
-     tl2b tlgamma tli2 tlngamma tlog tlog10 tlog1p tlog2	\
+     tfma tfmod tfms tfprintf tfrac tgamma tget_binary32 tget_d \
+     tget_d_2exp tget_f tget_ld_2exp tget_set_d64 tget_sj       \
+     tget_str tget_z tgmpop thyperbolic thypot tinp_str tj0 tj1 \
+     tjn tl2b tlgamma tli2 tlngamma tlog tlog10 tlog1p tlog2	\
      tmin_prec tminmax tmodf tmul tmul_2exp tmul_d tmul_ui	\
      tnext tout_str toutimpl tpow tpow3 tpow_all tpow_z		\
      tprintf trandom trec_sqrt tremquo trint troot		\
diff --git a/tests/tget_binary32.c b/tests/tget_binary32.c
new file mode 100644
index 000000000..045f2c9ab
--- /dev/null
+++ b/tests/tget_binary32.c
@@ -0,0 +1,115 @@
+/* Test file for mpfr_get_binary32 and mpfr_set_binary32
+
+Copyright 2009 Free Software Foundation, Inc.
+Contributed by the Arenaire and Cacao projects, INRIA.
+
+This file is part of the GNU MPFR Library.
+
+The GNU MPFR Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Lesser General Public License as published by
+the Free Software Foundation; either version 3 of the License, or (at your
+option) any later version.
+
+The GNU MPFR Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
+License for more details.
+
+You should have received a copy of the GNU Lesser General Public License
+along with the GNU MPFR Library; see the file COPYING.LESSER.  If not, see
+http://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc.,
+51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */
+
+#include <stdlib.h>
+#include "mpfr-test.h"
+
+int
+main (void)
+{
+  mpfr_t x;
+  float f;
+
+  tests_start_mpfr ();
+
+  mpfr_init2 (x, 24);
+
+  mpfr_set_nan (x);
+  f = mpfr_get_binary32 (x, MPFR_RNDN);
+  if (f == f)
+    {
+      printf ("Error for mpfr_get_binary32(NaN)\n");
+      exit (1);
+    }
+  mpfr_set_binary32 (x, f, MPFR_RNDN);
+  if (mpfr_nan_p (x) == 0)
+    {
+      printf ("Error for mpfr_set_binary32(NaN)\n");
+      exit (1);
+    }
+  
+  mpfr_set_inf (x, 1);
+  f = mpfr_get_binary32 (x, MPFR_RNDN);
+  mpfr_set_binary32 (x, f, MPFR_RNDN);
+  if (mpfr_inf_p (x) == 0 || mpfr_sgn (x) < 0)
+    {
+      printf ("Error for mpfr_set_binary32(mpfr_get_binary32(+Inf))n");
+      exit (1);
+    }
+  
+  mpfr_set_inf (x, -1);
+  f = mpfr_get_binary32 (x, MPFR_RNDN);
+  mpfr_set_binary32 (x, f, MPFR_RNDN);
+  if (mpfr_inf_p (x) == 0 || mpfr_sgn (x) > 0)
+    {
+      printf ("Error for mpfr_set_binary32(mpfr_get_binary32(-Inf))n");
+      exit (1);
+    }
+  
+  mpfr_set_ui (x, 0, MPFR_RNDN);
+  f = mpfr_get_binary32 (x, MPFR_RNDN);
+  mpfr_set_binary32 (x, f, MPFR_RNDN);
+  if (mpfr_zero_p (x) == 0 || MPFR_SIGN (x) < 0)
+    {
+      printf ("Error for mpfr_set_binary32(mpfr_get_binary32(+0))n");
+      exit (1);
+    }
+  
+  mpfr_set_ui (x, 0, MPFR_RNDN);
+  mpfr_neg (x, x, MPFR_RNDN);
+  f = mpfr_get_binary32 (x, MPFR_RNDN);
+  mpfr_set_binary32 (x, f, MPFR_RNDN);
+  if (mpfr_zero_p (x) == 0 || MPFR_SIGN (x) > 0)
+    {
+      printf ("Error for mpfr_set_binary32(mpfr_get_binary32(-0))n");
+      exit (1);
+    }
+  
+  mpfr_set_ui (x, 17, MPFR_RNDN);
+  f = mpfr_get_binary32 (x, MPFR_RNDN);
+  mpfr_set_binary32 (x, f, MPFR_RNDN);
+  if (mpfr_cmp_ui (x, 17) != 0)
+    {
+      printf ("Error for mpfr_set_binary32(mpfr_get_binary32(17))n");
+      printf ("expected 17\n");
+      printf ("got      ");
+      mpfr_dump (x);
+      exit (1);
+    }
+  
+  mpfr_set_si (x, -42, MPFR_RNDN);
+  f = mpfr_get_binary32 (x, MPFR_RNDN);
+  mpfr_set_binary32 (x, f, MPFR_RNDN);
+  if (mpfr_cmp_si (x, -42) != 0)
+    {
+      printf ("Error for mpfr_set_binary32(mpfr_get_binary32(-42))n");
+      printf ("expected -42\n");
+      printf ("got      ");
+      mpfr_dump (x);
+      exit (1);
+    }
+  
+  mpfr_clear (x);
+  
+  tests_end_mpfr ();
+  return 0;
+}
diff --git a/tests/tget_d.c b/tests/tget_d.c
index 9d0af0cef..be5f6160a 100644
--- a/tests/tget_d.c
+++ b/tests/tget_d.c
@@ -65,6 +65,26 @@ check_denorms (void)
         }
     }
 
+  mpfr_set_str_binary (x, "1e-1074");
+  dd = mpfr_get_d (x, MPFR_RNDA);
+  d2 = 0x1p-1074;
+  if (dd != d2)
+    {
+      printf ("Error for x=1e-1074, RNDA\n");
+      exit (1);
+    }
+
+  mpfr_set_str_binary (x, "1e-1075");
+  dd = mpfr_get_d (x, MPFR_RNDA);
+  d2 = 0x1p-1074;
+  if (dd != d2)
+    {
+      printf ("Error for x=1e-1075, RNDA\n");
+      printf ("expected %a\n", d2);
+      printf ("got      %a\n", dd);
+      exit (1);
+    }
+
   mpfr_clear (x);
   return fail;
 }