tset_ld.c revision 1.1.1.3
1/* Test file for mpfr_set_ld and mpfr_get_ld. 2 3Copyright 2002-2016 Free Software Foundation, Inc. 4Contributed by the AriC and Caramba 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 /* Do not convert d to double as this can give an overflow, which 51 may confuse compilers without IEEE 754 support (such as clang 52 -fsanitize=undefined), or trigger a trap if enabled. 53 The DOUBLE_ISNAN macro should work fine on long double. */ 54 if (DOUBLE_ISNAN (d)) 55 return 1; 56 LONGDOUBLE_NAN_ACTION (d, goto yes); 57 return 0; 58 yes: 59 return 1; 60} 61 62/* checks that a long double converted to a mpfr (with precision >=113), 63 then converted back to a long double gives the initial value, 64 or in other words mpfr_get_ld(mpfr_set_ld(d)) = d. 65*/ 66static void 67check_set_get (long double d, mpfr_t x) 68{ 69 int r; 70 long double e; 71 int inex; 72 73 for (r = 0; r < MPFR_RND_MAX; r++) 74 { 75 inex = mpfr_set_ld (x, d, (mpfr_rnd_t) r); 76 if (inex != 0) 77 { 78 mpfr_exp_t emin, emax; 79 emin = mpfr_get_emin (); 80 emax = mpfr_get_emax (); 81 printf ("Error: mpfr_set_ld should be exact\n"); 82 printf ("d=%1.30Le inex=%d\n", d, inex); 83 if (emin >= LONG_MIN) 84 printf ("emin=%ld\n", (long) emin); 85 if (emax <= LONG_MAX) 86 printf ("emax=%ld\n", (long) emax); 87 mpfr_dump (x); 88 exit (1); 89 } 90 e = mpfr_get_ld (x, (mpfr_rnd_t) r); 91 if ((Isnan_ld(d) && ! Isnan_ld(e)) || 92 (Isnan_ld(e) && ! Isnan_ld(d)) || 93 (e != d && !(Isnan_ld(e) && Isnan_ld(d)))) 94 { 95 printf ("Error: mpfr_get_ld o mpfr_set_ld <> Id\n"); 96 printf (" r=%d\n", r); 97 printf (" d=%1.30Le get_ld(set_ld(d))=%1.30Le\n", d, e); 98 ld_trace (" d", d); 99 printf (" x="); mpfr_out_str (NULL, 16, 0, x, MPFR_RNDN); 100 printf ("\n"); 101 ld_trace (" e", e); 102#ifdef MPFR_NANISNAN 103 if (Isnan_ld(d) || Isnan_ld(e)) 104 printf ("The reason is that NAN == NAN. Please look at the " 105 "configure output\nand Section \"In case of problem\" " 106 "of the INSTALL file.\n"); 107#endif 108 exit (1); 109 } 110 } 111} 112 113static void 114test_small (void) 115{ 116 mpfr_t x, y, z; 117 long double d; 118 119 mpfr_init2 (x, 64); 120 mpfr_init2 (y, 64); 121 mpfr_init2 (z, 64); 122 123 /* x = 11906603631607553907/2^(16381+64) */ 124 mpfr_set_str (x, "0.1010010100111100110000001110101101000111010110000001111101110011E-16381", 2, MPFR_RNDN); 125 d = mpfr_get_ld (x, MPFR_RNDN); /* infinite loop? */ 126 mpfr_set_ld (y, d, MPFR_RNDN); 127 mpfr_sub (z, x, y, MPFR_RNDN); 128 mpfr_abs (z, z, MPFR_RNDN); 129 mpfr_clear_erangeflag (); 130 /* If long double = double, d should be equal to 0; 131 in this case, everything is OK. */ 132 if (d != 0 && (mpfr_cmp_str (z, "1E-16434", 2, MPFR_RNDN) > 0 || 133 mpfr_erangeflag_p ())) 134 { 135 printf ("Error with x = "); 136 mpfr_out_str (NULL, 10, 21, x, MPFR_RNDN); 137 printf (" = "); 138 mpfr_out_str (NULL, 16, 0, x, MPFR_RNDN); 139 printf ("\n -> d = %.21Lg", d); 140 printf ("\n -> y = "); 141 mpfr_out_str (NULL, 10, 21, y, MPFR_RNDN); 142 printf (" = "); 143 mpfr_out_str (NULL, 16, 0, y, MPFR_RNDN); 144 printf ("\n -> |x-y| = "); 145 mpfr_out_str (NULL, 16, 0, z, MPFR_RNDN); 146 printf ("\n"); 147 exit (1); 148 } 149 150 mpfr_clear (x); 151 mpfr_clear (y); 152 mpfr_clear (z); 153} 154 155static void 156test_fixed_bugs (void) 157{ 158 mpfr_t x; 159 long double l, m; 160 161 /* bug found by Steve Kargl (2009-03-14) */ 162 mpfr_init2 (x, 64); 163 mpfr_set_ui_2exp (x, 1, -16447, MPFR_RNDN); 164 mpfr_get_ld (x, MPFR_RNDN); /* an assertion failed in init2.c:50 */ 165 166 /* bug reported by Jakub Jelinek (2010-10-17) 167 https://gforge.inria.fr/tracker/?func=detail&aid=11300 */ 168 mpfr_set_prec (x, MPFR_LDBL_MANT_DIG); 169 /* l = 0x1.23456789abcdef0123456789abcdp-914L; */ 170 l = 8.215640181713713164092636634579e-276; 171 mpfr_set_ld (x, l, MPFR_RNDN); 172 m = mpfr_get_ld (x, MPFR_RNDN); 173 if (m != l) 174 { 175 printf ("Error in get_ld o set_ld for l=%Le\n", l); 176 printf ("Got m=%Le instead of l\n", m); 177 exit (1); 178 } 179 180 /* another similar test which failed with extended double precision and the 181 generic code for mpfr_set_ld */ 182 /* l = 0x1.23456789abcdef0123456789abcdp-968L; */ 183 l = 4.560596445887084662336528403703e-292; 184 mpfr_set_ld (x, l, MPFR_RNDN); 185 m = mpfr_get_ld (x, MPFR_RNDN); 186 if (m != l) 187 { 188 printf ("Error in get_ld o set_ld for l=%Le\n", l); 189 printf ("Got m=%Le instead of l\n", m); 190 exit (1); 191 } 192 193 mpfr_clear (x); 194} 195 196/* bug reported by Walter Mascarenhas 197 https://sympa.inria.fr/sympa/arc/mpfr/2016-09/msg00005.html */ 198static void 199bug_20160907 (void) 200{ 201#if HAVE_LDOUBLE_IEEE_EXT_LITTLE 202 long double dn, ld; 203 mpfr_t mp; 204 long e; 205 mpfr_long_double_t x; 206 207 /* the following is the encoding of the smallest subnormal number 208 for HAVE_LDOUBLE_IEEE_EXT_LITTLE */ 209 x.s.manl = 1; 210 x.s.manh = 0; 211 x.s.expl = 0; 212 x.s.exph = 0; 213 x.s.sign= 0; 214 dn = x.ld; 215 e = -16445; 216 /* dn=2^e is now the smallest subnormal. */ 217 218 mpfr_init2 (mp, 64); 219 mpfr_set_ui_2exp (mp, 1, e - 1, MPFR_RNDN); 220 ld = mpfr_get_ld (mp, MPFR_RNDU); 221 /* since mp = 2^(e-1) and ld is rounded upwards, we should have 222 ld = 2^e */ 223 if (ld != dn) 224 { 225 printf ("Error, ld = %Le <> dn = %Le\n", ld, dn); 226 printf ("mp="); 227 mpfr_out_str (stdout, 10, 0, mp, MPFR_RNDN); 228 printf ("\n"); 229 exit (1); 230 } 231 232 /* check a few more numbers */ 233 for (e = -16446; e <= -16381; e++) 234 { 235 mpfr_set_ui_2exp (mp, 1, e, MPFR_RNDN); 236 ld = mpfr_get_ld (mp, MPFR_RNDU); 237 mpfr_set_ld (mp, ld, MPFR_RNDU); 238 /* mp is 2^e rounded up, thus should be >= 2^e */ 239 MPFR_ASSERTN(mpfr_cmp_ui_2exp (mp, 1, e) >= 0); 240 241 mpfr_set_ui_2exp (mp, 1, e, MPFR_RNDN); 242 ld = mpfr_get_ld (mp, MPFR_RNDD); 243 mpfr_set_ld (mp, ld, MPFR_RNDD); 244 /* mp is 2^e rounded down, thus should be <= 2^e */ 245 if (mpfr_cmp_ui_2exp (mp, 3, e) > 0) 246 { 247 printf ("Error, expected value <= 2^%ld\n", e); 248 printf ("got "); mpfr_dump (mp); 249 exit (1); 250 } 251 } 252 253 mpfr_clear (mp); 254#endif 255} 256 257int 258main (int argc, char *argv[]) 259{ 260 long double d, e; 261 mpfr_t x; 262 int i; 263 mpfr_exp_t emax; 264#ifdef WITH_FPU_CONTROL 265 fpu_control_t cw; 266 267 if (argc > 1) 268 { 269 cw = strtol(argv[1], NULL, 0); 270 printf ("FPU control word: 0x%x\n", (unsigned int) cw); 271 _FPU_SETCW (cw); 272 } 273#endif 274 275 tests_start_mpfr (); 276 mpfr_test_init (); 277 278 check_gcc33_bug (); 279 test_fixed_bugs (); 280 281 mpfr_init2 (x, MPFR_LDBL_MANT_DIG); 282 283#if !defined(MPFR_ERRDIVZERO) 284 /* check NaN */ 285 mpfr_set_nan (x); 286 d = mpfr_get_ld (x, MPFR_RNDN); 287 check_set_get (d, x); 288#endif 289 290 /* check +0.0 and -0.0 */ 291 d = 0.0; 292 check_set_get (d, x); 293 d = DBL_NEG_ZERO; 294 check_set_get (d, x); 295 296 /* check that the sign of -0.0 is set */ 297 mpfr_set_ld (x, DBL_NEG_ZERO, MPFR_RNDN); 298 if (MPFR_SIGN(x) > 0) 299 { 300 printf ("Error: sign of -0.0 is not set correctly\n"); 301#if _GMP_IEEE_FLOATS 302 exit (1); 303 /* Non IEEE doesn't support negative zero yet */ 304#endif 305 } 306 307#if !defined(MPFR_ERRDIVZERO) 308 /* check +Inf */ 309 mpfr_set_inf (x, 1); 310 d = mpfr_get_ld (x, MPFR_RNDN); 311 check_set_get (d, x); 312 313 /* check -Inf */ 314 mpfr_set_inf (x, -1); 315 d = mpfr_get_ld (x, MPFR_RNDN); 316 check_set_get (d, x); 317#endif 318 319 /* check the largest power of two */ 320 d = 1.0; while (d < LDBL_MAX / 2.0) d += d; 321 check_set_get (d, x); 322 check_set_get (-d, x); 323 324 /* check largest long double */ 325 d = LDBL_MAX; 326 check_set_get (d, x); 327 check_set_get (-d, x); 328 329 /* check the smallest power of two */ 330 d = 1.0; 331 while ((e = d / 2.0) != (long double) 0.0 && e != d) 332 d = e; 333 check_set_get (d, x); 334 check_set_get (-d, x); 335 336 /* check largest 2^(2^k) that is representable as a long double */ 337 d = (LDBL_MAX / 2) + (LDBL_MAX / 4 * LDBL_EPSILON); 338 check_set_get (d, x); 339 340 /* check that 2^i, 2^i+1 and 2^i-1 are correctly converted */ 341 d = 1.0; 342 for (i = 1; i < MPFR_LDBL_MANT_DIG; i++) 343 { 344 d = 2.0 * d; /* d = 2^i */ 345 check_set_get (d, x); 346 check_set_get (d + 1.0, x); 347 check_set_get (d - 1.0, x); 348 } 349 350 for (i = 0; i < 10000; i++) 351 { 352 mpfr_urandomb (x, RANDS); 353 d = mpfr_get_ld (x, MPFR_RNDN); 354 check_set_get (d, x); 355 } 356 357 /* check with reduced emax to exercise overflow */ 358 emax = mpfr_get_emax (); 359 mpfr_set_prec (x, 2); 360 set_emax (1); 361 mpfr_set_ld (x, (long double) 2.0, MPFR_RNDN); 362 MPFR_ASSERTN(mpfr_inf_p (x) && mpfr_sgn (x) > 0); 363 for (d = (long double) 2.0, i = 0; i < 13; i++, d *= d); 364 /* now d = 2^8192 */ 365 mpfr_set_ld (x, d, MPFR_RNDN); 366 MPFR_ASSERTN(mpfr_inf_p (x) && mpfr_sgn (x) > 0); 367 set_emax (emax); 368 369 mpfr_clear (x); 370 371 test_small (); 372 373 bug_20160907 (); 374 375 tests_end_mpfr (); 376 377 return 0; 378} 379