xref: /netbsd-src/external/lgpl3/mpfr/dist/tests/tset_ld.c (revision 0953dc8744b62dfdecb2f203329e730593755659)
1 /* Test file for mpfr_set_ld and mpfr_get_ld.
2 
3 Copyright 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011 Free Software Foundation, Inc.
4 Contributed by the Arenaire and Cacao projects, INRIA.
5 
6 This file is part of the GNU MPFR Library.
7 
8 The GNU MPFR Library is free software; you can redistribute it and/or modify
9 it under the terms of the GNU Lesser General Public License as published by
10 the Free Software Foundation; either version 3 of the License, or (at your
11 option) any later version.
12 
13 The GNU MPFR Library is distributed in the hope that it will be useful, but
14 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
15 or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
16 License for more details.
17 
18 You should have received a copy of the GNU Lesser General Public License
19 along with the GNU MPFR Library; see the file COPYING.LESSER.  If not, see
20 http://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc.,
21 51 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 
33 static void
34 check_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 
47 static int
48 Isnan_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 */
63 static void
64 check_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 
104 static void
105 test_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 
146 static void
147 test_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 
187 int
188 main (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