1/* Test file for mpfr_set_ld and mpfr_get_ld.
2
3Copyright 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011 Free Software Foundation, Inc.
4Contributed by the Arenaire and Cacao projects, INRIA.
5
6This file is part of the GNU MPFR Library.
7
8The GNU MPFR Library is free software; you can redistribute it and/or modify
9it under the terms of the GNU Lesser General Public License as published by
10the Free Software Foundation; either version 3 of the License, or (at your
11option) any later version.
12
13The GNU MPFR Library is distributed in the hope that it will be useful, but
14WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
15or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
16License for more details.
17
18You should have received a copy of the GNU Lesser General Public License
19along with the GNU MPFR Library; see the file COPYING.LESSER.  If not, see
20http://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc.,
2151 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */
22
23#include <stdio.h>
24#include <stdlib.h>
25#include <float.h>
26#include <limits.h>
27#ifdef WITH_FPU_CONTROL
28#include <fpu_control.h>
29#endif
30
31#include "mpfr-test.h"
32
33static void
34check_gcc33_bug (void)
35{
36  volatile long double x;
37  x = (long double) 9007199254740992.0 + 1.0;
38  if (x != 0.0)
39    return;  /* OK */
40  printf
41    ("Detected optimization bug of gcc 3.3 on Alpha concerning long double\n"
42     "comparisons; set_ld tests might fail (set_ld won't work correctly).\n"
43     "See http://gcc.gnu.org/ml/gcc-bugs/2003-10/msg00853.html for more\n"
44     "information.\n");
45}
46
47static int
48Isnan_ld (long double d)
49{
50  double e = (double) d;
51  if (DOUBLE_ISNAN (e))
52    return 1;
53  LONGDOUBLE_NAN_ACTION (d, goto yes);
54  return 0;
55 yes:
56  return 1;
57}
58
59/* checks that a long double converted to a mpfr (with precision >=113),
60   then converted back to a long double gives the initial value,
61   or in other words mpfr_get_ld(mpfr_set_ld(d)) = d.
62*/
63static void
64check_set_get (long double d, mpfr_t x)
65{
66  int r;
67  long double e;
68  int inex;
69
70  for (r = 0; r < MPFR_RND_MAX; r++)
71    {
72      inex = mpfr_set_ld (x, d, (mpfr_rnd_t) r);
73      if (inex != 0)
74        {
75          printf ("Error: mpfr_set_ld should be exact\n");
76          printf ("d=%1.30Le inex=%d\n", d, inex);
77          printf ("emin=%ld emax=%ld\n", mpfr_get_emin (), mpfr_get_emax ());
78          mpfr_dump (x);
79          exit (1);
80        }
81      e = mpfr_get_ld (x, (mpfr_rnd_t) r);
82      if ((Isnan_ld(d) && ! Isnan_ld(e)) ||
83          (Isnan_ld(e) && ! Isnan_ld(d)) ||
84          (e != d && !(Isnan_ld(e) && Isnan_ld(d))))
85        {
86          printf ("Error: mpfr_get_ld o mpfr_set_ld <> Id\n");
87          printf ("  r=%d\n", r);
88          printf ("  d=%1.30Le get_ld(set_ld(d))=%1.30Le\n", d, e);
89          ld_trace ("  d", d);
90          printf ("  x="); mpfr_out_str (NULL, 16, 0, x, MPFR_RNDN);
91          printf ("\n");
92          ld_trace ("  e", e);
93#ifdef MPFR_NANISNAN
94          if (Isnan_ld(d) || Isnan_ld(e))
95            printf ("The reason is that NAN == NAN. Please look at the "
96                    "configure output\nand Section \"In case of problem\" "
97                    "of the INSTALL file.\n");
98#endif
99          exit (1);
100        }
101    }
102}
103
104static void
105test_small (void)
106{
107  mpfr_t x, y, z;
108  long double d;
109
110  mpfr_init2 (x, 64);
111  mpfr_init2 (y, 64);
112  mpfr_init2 (z, 64);
113
114  /* x = 11906603631607553907/2^(16381+64) */
115  mpfr_set_str (x, "0.1010010100111100110000001110101101000111010110000001111101110011E-16381", 2, MPFR_RNDN);
116  d = mpfr_get_ld (x, MPFR_RNDN);  /* infinite loop? */
117  mpfr_set_ld (y, d, MPFR_RNDN);
118  mpfr_sub (z, x, y, MPFR_RNDN);
119  mpfr_abs (z, z, MPFR_RNDN);
120  mpfr_clear_erangeflag ();
121  /* If long double = double, d should be equal to 0;
122     in this case, everything is OK. */
123  if (d != 0 && (mpfr_cmp_str (z, "1E-16434", 2, MPFR_RNDN) > 0 ||
124                 mpfr_erangeflag_p ()))
125    {
126      printf ("Error with x = ");
127      mpfr_out_str (NULL, 10, 21, x, MPFR_RNDN);
128      printf (" = ");
129      mpfr_out_str (NULL, 16, 0, x, MPFR_RNDN);
130      printf ("\n        -> d = %.21Lg", d);
131      printf ("\n        -> y = ");
132      mpfr_out_str (NULL, 10, 21, y, MPFR_RNDN);
133      printf (" = ");
134      mpfr_out_str (NULL, 16, 0, y, MPFR_RNDN);
135      printf ("\n        -> |x-y| = ");
136      mpfr_out_str (NULL, 16, 0, z, MPFR_RNDN);
137      printf ("\n");
138      exit (1);
139    }
140
141  mpfr_clear (x);
142  mpfr_clear (y);
143  mpfr_clear (z);
144}
145
146static void
147test_fixed_bugs (void)
148{
149  mpfr_t x;
150  long double l, m;
151
152  /* bug found by Steve Kargl (2009-03-14) */
153  mpfr_init2 (x, 64);
154  mpfr_set_ui_2exp (x, 1, -16447, MPFR_RNDN);
155  mpfr_get_ld (x, MPFR_RNDN);  /* an assertion failed in init2.c:50 */
156
157  /* bug reported by Jakub Jelinek (2010-10-17)
158     https://gforge.inria.fr/tracker/?func=detail&aid=11300 */
159  mpfr_set_prec (x, MPFR_LDBL_MANT_DIG);
160  /* l = 0x1.23456789abcdef0123456789abcdp-914L; */
161  l = 8.215640181713713164092636634579e-276;
162  mpfr_set_ld (x, l, MPFR_RNDN);
163  m = mpfr_get_ld (x, MPFR_RNDN);
164  if (m != l)
165    {
166      printf ("Error in get_ld o set_ld for l=%Le\n", l);
167      printf ("Got m=%Le instead of l\n", m);
168      exit (1);
169    }
170
171  /* another similar test which failed with extended double precision and the
172     generic code for mpfr_set_ld */
173  /* l = 0x1.23456789abcdef0123456789abcdp-968L; */
174  l = 4.560596445887084662336528403703e-292;
175  mpfr_set_ld (x, l, MPFR_RNDN);
176  m = mpfr_get_ld (x, MPFR_RNDN);
177  if (m != l)
178    {
179      printf ("Error in get_ld o set_ld for l=%Le\n", l);
180      printf ("Got m=%Le instead of l\n", m);
181      exit (1);
182    }
183
184  mpfr_clear (x);
185}
186
187int
188main (int argc, char *argv[])
189{
190  long double d, e;
191  mpfr_t x;
192  int i;
193  mpfr_exp_t emax;
194#ifdef WITH_FPU_CONTROL
195  fpu_control_t cw;
196
197  if (argc > 1)
198    {
199      cw = strtol(argv[1], NULL, 0);
200      printf ("FPU control word: 0x%x\n", (unsigned int) cw);
201      _FPU_SETCW (cw);
202    }
203#endif
204
205  check_gcc33_bug ();
206  test_fixed_bugs ();
207
208  tests_start_mpfr ();
209  mpfr_test_init ();
210
211  mpfr_init2 (x, MPFR_LDBL_MANT_DIG);
212
213  /* check NaN */
214  mpfr_set_nan (x);
215  d = mpfr_get_ld (x, MPFR_RNDN);
216  check_set_get (d, x);
217
218  /* check +0.0 and -0.0 */
219  d = 0.0;
220  check_set_get (d, x);
221  d = DBL_NEG_ZERO;
222  check_set_get (d, x);
223
224  /* check that the sign of -0.0 is set */
225  mpfr_set_ld (x, DBL_NEG_ZERO, MPFR_RNDN);
226  if (MPFR_SIGN(x) > 0)
227    {
228      printf ("Error: sign of -0.0 is not set correctly\n");
229#ifdef _GMP_IEEE_FLOATS
230      exit (1);
231      /* Non IEEE doesn't support negative zero yet */
232#endif
233    }
234
235  /* check +Inf */
236  mpfr_set_inf (x, 1);
237  d = mpfr_get_ld (x, MPFR_RNDN);
238  check_set_get (d, x);
239
240  /* check -Inf */
241  mpfr_set_inf (x, -1);
242  d = mpfr_get_ld (x, MPFR_RNDN);
243  check_set_get (d, x);
244
245  /* check the largest power of two */
246  d = 1.0; while (d < LDBL_MAX / 2.0) d += d;
247  check_set_get (d, x);
248  check_set_get (-d, x);
249
250  /* check largest long double */
251  d = LDBL_MAX;
252  check_set_get (d, x);
253  check_set_get (-d, x);
254
255  /* check the smallest power of two */
256  d = 1.0;
257  while ((e = d / 2.0) != (long double) 0.0 && e != d)
258    d = e;
259  check_set_get (d, x);
260  check_set_get (-d, x);
261
262  /* check largest 2^(2^k) that is representable as a long double */
263  d = (LDBL_MAX / 2) + (LDBL_MAX / 4 * LDBL_EPSILON);
264  check_set_get (d, x);
265
266  /* check that 2^i, 2^i+1 and 2^i-1 are correctly converted */
267  d = 1.0;
268  for (i = 1; i < MPFR_LDBL_MANT_DIG; i++)
269    {
270      d = 2.0 * d; /* d = 2^i */
271      check_set_get (d, x);
272      check_set_get (d + 1.0, x);
273      check_set_get (d - 1.0, x);
274    }
275
276  for (i = 0; i < 10000; i++)
277    {
278      mpfr_urandomb (x, RANDS);
279      d = mpfr_get_ld (x, MPFR_RNDN);
280      check_set_get (d, x);
281    }
282
283  /* check with reduced emax to exercise overflow */
284  emax = mpfr_get_emax ();
285  mpfr_set_prec (x, 2);
286  set_emax (1);
287  mpfr_set_ld (x, (long double) 2.0, MPFR_RNDN);
288  MPFR_ASSERTN(mpfr_inf_p (x) && mpfr_sgn (x) > 0);
289  for (d = (long double) 2.0, i = 0; i < 13; i++, d *= d);
290  /* now d = 2^8192 */
291  mpfr_set_ld (x, d, MPFR_RNDN);
292  MPFR_ASSERTN(mpfr_inf_p (x) && mpfr_sgn (x) > 0);
293  set_emax (emax);
294
295  mpfr_clear (x);
296
297  test_small ();
298
299  tests_end_mpfr ();
300
301  return 0;
302}
303