1 /* Test file for mpfr_mul. 2 3 Copyright 1999, 2001-2018 Free Software Foundation, Inc. 4 Contributed by the AriC and Caramba 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 "mpfr-test.h" 24 25 #ifdef CHECK_EXTERNAL 26 static int 27 test_mul (mpfr_ptr a, mpfr_srcptr b, mpfr_srcptr c, mpfr_rnd_t rnd_mode) 28 { 29 int res; 30 int ok = rnd_mode == MPFR_RNDN && mpfr_number_p (b) && mpfr_number_p (c); 31 if (ok) 32 { 33 mpfr_print_raw (b); 34 printf (" "); 35 mpfr_print_raw (c); 36 } 37 res = mpfr_mul (a, b, c, rnd_mode); 38 if (ok) 39 { 40 printf (" "); 41 mpfr_print_raw (a); 42 printf ("\n"); 43 } 44 return res; 45 } 46 #else 47 #define test_mul mpfr_mul 48 #endif 49 50 /* checks that xs * ys gives the expected result res */ 51 static void 52 check (const char *xs, const char *ys, mpfr_rnd_t rnd_mode, 53 unsigned int px, unsigned int py, unsigned int pz, const char *res) 54 { 55 mpfr_t xx, yy, zz; 56 57 mpfr_init2 (xx, px); 58 mpfr_init2 (yy, py); 59 mpfr_init2 (zz, pz); 60 mpfr_set_str1 (xx, xs); 61 mpfr_set_str1 (yy, ys); 62 test_mul(zz, xx, yy, rnd_mode); 63 if (mpfr_cmp_str1 (zz, res) ) 64 { 65 printf ("(1) mpfr_mul failed for x=%s y=%s with rnd=%s\n", 66 xs, ys, mpfr_print_rnd_mode (rnd_mode)); 67 printf ("correct is %s, mpfr_mul gives ", res); 68 mpfr_out_str (stdout, 10, 0, zz, MPFR_RNDN); 69 putchar ('\n'); 70 exit (1); 71 } 72 mpfr_clear (xx); 73 mpfr_clear (yy); 74 mpfr_clear (zz); 75 } 76 77 static void 78 check53 (const char *xs, const char *ys, mpfr_rnd_t rnd_mode, const char *zs) 79 { 80 mpfr_t xx, yy, zz; 81 82 mpfr_inits2 (53, xx, yy, zz, (mpfr_ptr) 0); 83 mpfr_set_str1 (xx, xs); 84 mpfr_set_str1 (yy, ys); 85 test_mul (zz, xx, yy, rnd_mode); 86 if (mpfr_cmp_str1 (zz, zs) ) 87 { 88 printf ("(2) mpfr_mul failed for x=%s y=%s with rnd=%s\n", 89 xs, ys, mpfr_print_rnd_mode(rnd_mode)); 90 printf ("correct result is %s,\n mpfr_mul gives ", zs); 91 mpfr_out_str (stdout, 10, 0, zz, MPFR_RNDN); 92 putchar ('\n'); 93 exit (1); 94 } 95 mpfr_clears (xx, yy, zz, (mpfr_ptr) 0); 96 } 97 98 /* checks that x*y gives the right result with 24 bits of precision */ 99 static void 100 check24 (const char *xs, const char *ys, mpfr_rnd_t rnd_mode, const char *zs) 101 { 102 mpfr_t xx, yy, zz; 103 104 mpfr_inits2 (24, xx, yy, zz, (mpfr_ptr) 0); 105 mpfr_set_str1 (xx, xs); 106 mpfr_set_str1 (yy, ys); 107 test_mul (zz, xx, yy, rnd_mode); 108 if (mpfr_cmp_str1 (zz, zs) ) 109 { 110 printf ("(3) mpfr_mul failed for x=%s y=%s with " 111 "rnd=%s\n", xs, ys, mpfr_print_rnd_mode(rnd_mode)); 112 printf ("correct result is gives %s, mpfr_mul gives ", zs); 113 mpfr_out_str (stdout, 10, 0, zz, MPFR_RNDN); 114 putchar ('\n'); 115 exit (1); 116 } 117 mpfr_clears (xx, yy, zz, (mpfr_ptr) 0); 118 } 119 120 /* the following examples come from the paper "Number-theoretic Test 121 Generation for Directed Rounding" from Michael Parks, Table 1 */ 122 static void 123 check_float (void) 124 { 125 check24("8388609.0", "8388609.0", MPFR_RNDN, "70368760954880.0"); 126 check24("16777213.0", "8388609.0", MPFR_RNDN, "140737479966720.0"); 127 check24("8388611.0", "8388609.0", MPFR_RNDN, "70368777732096.0"); 128 check24("12582911.0", "8388610.0", MPFR_RNDN, "105553133043712.0"); 129 check24("12582914.0", "8388610.0", MPFR_RNDN, "105553158209536.0"); 130 check24("13981013.0", "8388611.0", MPFR_RNDN, "117281279442944.0"); 131 check24("11184811.0", "8388611.0", MPFR_RNDN, "93825028587520.0"); 132 check24("11184810.0", "8388611.0", MPFR_RNDN, "93825020198912.0"); 133 check24("13981014.0", "8388611.0", MPFR_RNDN, "117281287831552.0"); 134 135 check24("8388609.0", "8388609.0", MPFR_RNDZ, "70368760954880.0"); 136 check24("16777213.0", "8388609.0", MPFR_RNDZ, "140737471578112.0"); 137 check24("8388611.0", "8388609.0", MPFR_RNDZ, "70368777732096.0"); 138 check24("12582911.0", "8388610.0", MPFR_RNDZ, "105553124655104.0"); 139 check24("12582914.0", "8388610.0", MPFR_RNDZ, "105553158209536.0"); 140 check24("13981013.0", "8388611.0", MPFR_RNDZ, "117281271054336.0"); 141 check24("11184811.0", "8388611.0", MPFR_RNDZ, "93825028587520.0"); 142 check24("11184810.0", "8388611.0", MPFR_RNDZ, "93825011810304.0"); 143 check24("13981014.0", "8388611.0", MPFR_RNDZ, "117281287831552.0"); 144 145 check24("8388609.0", "8388609.0", MPFR_RNDU, "70368769343488.0"); 146 check24("16777213.0", "8388609.0", MPFR_RNDU, "140737479966720.0"); 147 check24("8388611.0", "8388609.0", MPFR_RNDU, "70368786120704.0"); 148 check24("12582911.0", "8388610.0", MPFR_RNDU, "105553133043712.0"); 149 check24("12582914.0", "8388610.0", MPFR_RNDU, "105553166598144.0"); 150 check24("13981013.0", "8388611.0", MPFR_RNDU, "117281279442944.0"); 151 check24("11184811.0", "8388611.0", MPFR_RNDU, "93825036976128.0"); 152 check24("11184810.0", "8388611.0", MPFR_RNDU, "93825020198912.0"); 153 check24("13981014.0", "8388611.0", MPFR_RNDU, "117281296220160.0"); 154 155 check24("8388609.0", "8388609.0", MPFR_RNDD, "70368760954880.0"); 156 check24("16777213.0", "8388609.0", MPFR_RNDD, "140737471578112.0"); 157 check24("8388611.0", "8388609.0", MPFR_RNDD, "70368777732096.0"); 158 check24("12582911.0", "8388610.0", MPFR_RNDD, "105553124655104.0"); 159 check24("12582914.0", "8388610.0", MPFR_RNDD, "105553158209536.0"); 160 check24("13981013.0", "8388611.0", MPFR_RNDD, "117281271054336.0"); 161 check24("11184811.0", "8388611.0", MPFR_RNDD, "93825028587520.0"); 162 check24("11184810.0", "8388611.0", MPFR_RNDD, "93825011810304.0"); 163 check24("13981014.0", "8388611.0", MPFR_RNDD, "117281287831552.0"); 164 } 165 166 /* check sign of result */ 167 static void 168 check_sign (void) 169 { 170 mpfr_t a, b; 171 172 mpfr_init2 (a, 53); 173 mpfr_init2 (b, 53); 174 mpfr_set_si (a, -1, MPFR_RNDN); 175 mpfr_set_ui (b, 2, MPFR_RNDN); 176 test_mul(a, b, b, MPFR_RNDN); 177 if (mpfr_cmp_ui (a, 4) ) 178 { 179 printf ("2.0*2.0 gives \n"); 180 mpfr_out_str (stdout, 10, 0, a, MPFR_RNDN); 181 putchar ('\n'); 182 exit (1); 183 } 184 mpfr_clear(a); mpfr_clear(b); 185 } 186 187 /* checks that the inexact return value is correct */ 188 static void 189 check_exact (void) 190 { 191 mpfr_t a, b, c, d; 192 mpfr_prec_t prec; 193 int i, inexact; 194 mpfr_rnd_t rnd; 195 196 mpfr_init (a); 197 mpfr_init (b); 198 mpfr_init (c); 199 mpfr_init (d); 200 201 mpfr_set_prec (a, 17); 202 mpfr_set_prec (b, 17); 203 mpfr_set_prec (c, 32); 204 mpfr_set_str_binary (a, "1.1000111011000100e-1"); 205 mpfr_set_str_binary (b, "1.0010001111100111e-1"); 206 if (test_mul (c, a, b, MPFR_RNDZ)) 207 { 208 printf ("wrong return value (1)\n"); 209 exit (1); 210 } 211 212 for (prec = MPFR_PREC_MIN; prec < 100; prec++) 213 { 214 mpfr_set_prec (a, prec); 215 mpfr_set_prec (b, prec); 216 /* for prec=1, ensure PREC(c) >= 1 */ 217 mpfr_set_prec (c, 2 * prec - 2 + (prec == 1)); 218 mpfr_set_prec (d, 2 * prec); 219 for (i = 0; i < 1000; i++) 220 { 221 mpfr_urandomb (a, RANDS); 222 mpfr_urandomb (b, RANDS); 223 rnd = RND_RAND (); 224 inexact = test_mul (c, a, b, rnd); 225 if (test_mul (d, a, b, rnd)) /* should be always exact */ 226 { 227 printf ("unexpected inexact return value\n"); 228 exit (1); 229 } 230 if ((inexact == 0) && mpfr_cmp (c, d) && rnd != MPFR_RNDF) 231 { 232 printf ("rnd=%s: inexact=0 but results differ\n", 233 mpfr_print_rnd_mode (rnd)); 234 printf ("a="); 235 mpfr_out_str (stdout, 2, 0, a, rnd); 236 printf ("\nb="); 237 mpfr_out_str (stdout, 2, 0, b, rnd); 238 printf ("\nc="); 239 mpfr_out_str (stdout, 2, 0, c, rnd); 240 printf ("\nd="); 241 mpfr_out_str (stdout, 2, 0, d, rnd); 242 printf ("\n"); 243 exit (1); 244 } 245 else if (inexact && (mpfr_cmp (c, d) == 0) && rnd != MPFR_RNDF) 246 { 247 printf ("inexact!=0 but results agree\n"); 248 printf ("prec=%u rnd=%s a=", (unsigned int) prec, 249 mpfr_print_rnd_mode (rnd)); 250 mpfr_out_str (stdout, 2, 0, a, rnd); 251 printf ("\nb="); 252 mpfr_out_str (stdout, 2, 0, b, rnd); 253 printf ("\nc="); 254 mpfr_out_str (stdout, 2, 0, c, rnd); 255 printf ("\nd="); 256 mpfr_out_str (stdout, 2, 0, d, rnd); 257 printf ("\n"); 258 exit (1); 259 } 260 } 261 } 262 263 mpfr_clear (a); 264 mpfr_clear (b); 265 mpfr_clear (c); 266 mpfr_clear (d); 267 } 268 269 static void 270 check_max(void) 271 { 272 mpfr_t xx, yy, zz; 273 mpfr_exp_t emin; 274 275 mpfr_init2(xx, 4); 276 mpfr_init2(yy, 4); 277 mpfr_init2(zz, 4); 278 mpfr_set_str1 (xx, "0.68750"); 279 mpfr_mul_2si(xx, xx, MPFR_EMAX_DEFAULT/2, MPFR_RNDN); 280 mpfr_set_str1 (yy, "0.68750"); 281 mpfr_mul_2si(yy, yy, MPFR_EMAX_DEFAULT - MPFR_EMAX_DEFAULT/2 + 1, MPFR_RNDN); 282 mpfr_clear_flags(); 283 test_mul(zz, xx, yy, MPFR_RNDU); 284 if (!(mpfr_overflow_p() && MPFR_IS_INF(zz))) 285 { 286 printf("check_max failed (should be an overflow)\n"); 287 exit(1); 288 } 289 290 mpfr_clear_flags(); 291 test_mul(zz, xx, yy, MPFR_RNDD); 292 if (mpfr_overflow_p() || MPFR_IS_INF(zz)) 293 { 294 printf("check_max failed (should NOT be an overflow)\n"); 295 exit(1); 296 } 297 mpfr_set_str1 (xx, "0.93750"); 298 mpfr_mul_2si(xx, xx, MPFR_EMAX_DEFAULT, MPFR_RNDN); 299 if (!(MPFR_IS_FP(xx) && MPFR_IS_FP(zz))) 300 { 301 printf("check_max failed (internal error)\n"); 302 exit(1); 303 } 304 if (mpfr_cmp(xx, zz) != 0) 305 { 306 printf("check_max failed: got "); 307 mpfr_out_str(stdout, 2, 0, zz, MPFR_RNDZ); 308 printf(" instead of "); 309 mpfr_out_str(stdout, 2, 0, xx, MPFR_RNDZ); 310 printf("\n"); 311 exit(1); 312 } 313 314 /* check underflow */ 315 emin = mpfr_get_emin (); 316 set_emin (0); 317 mpfr_set_str_binary (xx, "0.1E0"); 318 mpfr_set_str_binary (yy, "0.1E0"); 319 test_mul (zz, xx, yy, MPFR_RNDN); 320 /* exact result is 0.1E-1, which should round to 0 */ 321 MPFR_ASSERTN(mpfr_cmp_ui (zz, 0) == 0 && MPFR_IS_POS(zz)); 322 set_emin (emin); 323 324 /* coverage test for mpfr_powerof2_raw */ 325 emin = mpfr_get_emin (); 326 set_emin (0); 327 mpfr_set_prec (xx, mp_bits_per_limb + 1); 328 mpfr_set_str_binary (xx, "0.1E0"); 329 mpfr_nextabove (xx); 330 mpfr_set_str_binary (yy, "0.1E0"); 331 test_mul (zz, xx, yy, MPFR_RNDN); 332 /* exact result is just above 0.1E-1, which should round to minfloat */ 333 MPFR_ASSERTN(mpfr_cmp (zz, yy) == 0); 334 set_emin (emin); 335 336 mpfr_clear(xx); 337 mpfr_clear(yy); 338 mpfr_clear(zz); 339 } 340 341 static void 342 check_min(void) 343 { 344 mpfr_t xx, yy, zz; 345 346 mpfr_init2(xx, 4); 347 mpfr_init2(yy, 4); 348 mpfr_init2(zz, 3); 349 mpfr_set_str1(xx, "0.9375"); 350 mpfr_mul_2si(xx, xx, MPFR_EMIN_DEFAULT/2, MPFR_RNDN); 351 mpfr_set_str1(yy, "0.9375"); 352 mpfr_mul_2si(yy, yy, MPFR_EMIN_DEFAULT - MPFR_EMIN_DEFAULT/2 - 1, MPFR_RNDN); 353 test_mul(zz, xx, yy, MPFR_RNDD); 354 if (MPFR_NOTZERO (zz)) 355 { 356 printf("check_min failed: got "); 357 mpfr_out_str(stdout, 2, 0, zz, MPFR_RNDZ); 358 printf(" instead of 0\n"); 359 exit(1); 360 } 361 362 test_mul(zz, xx, yy, MPFR_RNDU); 363 mpfr_set_str1 (xx, "0.5"); 364 mpfr_mul_2si(xx, xx, MPFR_EMIN_DEFAULT, MPFR_RNDN); 365 if (mpfr_sgn(xx) <= 0) 366 { 367 printf("check_min failed (internal error)\n"); 368 exit(1); 369 } 370 if (mpfr_cmp(xx, zz) != 0) 371 { 372 printf("check_min failed: got "); 373 mpfr_out_str(stdout, 2, 0, zz, MPFR_RNDZ); 374 printf(" instead of "); 375 mpfr_out_str(stdout, 2, 0, xx, MPFR_RNDZ); 376 printf("\n"); 377 exit(1); 378 } 379 380 mpfr_clear(xx); 381 mpfr_clear(yy); 382 mpfr_clear(zz); 383 } 384 385 static void 386 check_nans (void) 387 { 388 mpfr_t p, x, y; 389 390 mpfr_init2 (x, 123L); 391 mpfr_init2 (y, 123L); 392 mpfr_init2 (p, 123L); 393 394 /* nan * 0 == nan */ 395 mpfr_set_nan (x); 396 mpfr_set_ui (y, 0L, MPFR_RNDN); 397 test_mul (p, x, y, MPFR_RNDN); 398 MPFR_ASSERTN (mpfr_nan_p (p)); 399 400 /* 1 * nan == nan */ 401 mpfr_set_ui (x, 1L, MPFR_RNDN); 402 mpfr_set_nan (y); 403 test_mul (p, x, y, MPFR_RNDN); 404 MPFR_ASSERTN (mpfr_nan_p (p)); 405 406 /* 0 * +inf == nan */ 407 mpfr_set_ui (x, 0L, MPFR_RNDN); 408 mpfr_set_nan (y); 409 test_mul (p, x, y, MPFR_RNDN); 410 MPFR_ASSERTN (mpfr_nan_p (p)); 411 412 /* +1 * +inf == +inf */ 413 mpfr_set_ui (x, 1L, MPFR_RNDN); 414 mpfr_set_inf (y, 1); 415 test_mul (p, x, y, MPFR_RNDN); 416 MPFR_ASSERTN (mpfr_inf_p (p)); 417 MPFR_ASSERTN (mpfr_sgn (p) > 0); 418 419 /* -1 * +inf == -inf */ 420 mpfr_set_si (x, -1L, MPFR_RNDN); 421 mpfr_set_inf (y, 1); 422 test_mul (p, x, y, MPFR_RNDN); 423 MPFR_ASSERTN (mpfr_inf_p (p)); 424 MPFR_ASSERTN (mpfr_sgn (p) < 0); 425 426 mpfr_clear (x); 427 mpfr_clear (y); 428 mpfr_clear (p); 429 } 430 431 #define BUFSIZE 1552 432 433 static void 434 get_string (char *s, FILE *fp) 435 { 436 int c, n = BUFSIZE; 437 438 while ((c = getc (fp)) != '\n') 439 { 440 if (c == EOF) 441 { 442 printf ("Error in get_string: end of file\n"); 443 exit (1); 444 } 445 *(unsigned char *)s++ = c; 446 if (--n == 0) 447 { 448 printf ("Error in get_string: buffer is too small\n"); 449 exit (1); 450 } 451 } 452 *s = '\0'; 453 } 454 455 static void 456 check_regression (void) 457 { 458 mpfr_t x, y, z; 459 int i; 460 FILE *fp; 461 char s[BUFSIZE]; 462 463 mpfr_inits2 (6177, x, y, z, (mpfr_ptr) 0); 464 /* we read long strings from a file since ISO C90 does not support strings of 465 length > 509 */ 466 fp = src_fopen ("tmul.dat", "r"); 467 if (fp == NULL) 468 { 469 fprintf (stderr, "Error, cannot open tmul.dat in srcdir\n"); 470 exit (1); 471 } 472 get_string (s, fp); 473 mpfr_set_str (y, s, 16, MPFR_RNDN); 474 get_string (s, fp); 475 mpfr_set_str (z, s, 16, MPFR_RNDN); 476 i = mpfr_mul (x, y, z, MPFR_RNDN); 477 get_string (s, fp); 478 if (mpfr_cmp_str (x, s, 16, MPFR_RNDN) != 0 || i != -1) 479 { 480 printf ("Regression test 1 failed (i=%d, expected -1)\nx=", i); 481 mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN); putchar ('\n'); 482 exit (1); 483 } 484 fclose (fp); 485 486 mpfr_set_prec (x, 606); 487 mpfr_set_prec (y, 606); 488 mpfr_set_prec (z, 606); 489 490 mpfr_set_str (y, "-f.ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff92daefc3f8052ca9f58736564d9e93e62d324@-1", 16, MPFR_RNDN); 491 mpfr_set_str (z, "-f.ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff92daefc3f8052ca9f58736564d9e93e62d324@-1", 16, MPFR_RNDN); 492 i = mpfr_mul (x, y, z, MPFR_RNDU); 493 mpfr_set_str (y, "f.ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff25b5df87f00a5953eb0e6cac9b3d27cc5a64c@-1", 16, MPFR_RNDN); 494 if (mpfr_cmp (x, y) || i <= 0) 495 { 496 printf ("Regression test (2) failed! (i=%d - Expected 1)\n", i); 497 mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN); putchar ('\n'); 498 exit (1); 499 } 500 501 mpfr_set_prec (x, 184); 502 mpfr_set_prec (y, 92); 503 mpfr_set_prec (z, 1023); 504 505 mpfr_set_str (y, "6.9b8c8498882770d8038c3b0@-1", 16, MPFR_RNDN); 506 mpfr_set_str (z, "7.44e24b986e7fb296f1e936ce749fec3504cbf0d5ba769466b1c9f1578115efd5d29b4c79271191a920a99280c714d3a657ad6e3afbab77ffce9d697e9bb9110e26d676069afcea8b69f1d1541f2365042d80a97c21dcccd8ace4f1bb58b49922003e738e6f37bb82ef653cb2e87f763974e6ae50ae54e7724c38b80653e3289@255", 16, MPFR_RNDN); 507 i = mpfr_mul (x, y, z, MPFR_RNDU); 508 mpfr_set_prec (y, 184); 509 mpfr_set_str (y, "3.0080038f2ac5054e3e71ccbb95f76aaab2221715025a28@255", 510 16, MPFR_RNDN); 511 if (mpfr_cmp (x, y) || i <= 0) 512 { 513 printf ("Regression test (4) failed! (i=%d - expected 1)\n", i); 514 printf ("Ref: 3.0080038f2ac5054e3e71ccbb95f76aaab2221715025a28@255\n" 515 "Got: "); 516 mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN); 517 printf ("\n"); 518 exit (1); 519 } 520 521 mpfr_set_prec (x, 908); 522 mpfr_set_prec (y, 908); 523 mpfr_set_prec (z, 908); 524 mpfr_set_str (y, "-f.fffffffffffffffffffffffffffffffffffffffffffffffffffffff" 525 "fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff" 526 "ffffffffffffffffffffffffffffffffffffffffffffffffffffff99be91f83ec6f0ed28a3d42" 527 "e6e9a327230345ea6@-1", 16, MPFR_RNDN); 528 mpfr_set_str (z, "-f.fffffffffffffffffffffffffffffffffffffffffffffffffffffff" 529 "fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff" 530 "ffffffffffffffffffffffffffffffffffffffffffffffffffffff99be91f83ec6f0ed28a3d42" 531 "e6e9a327230345ea6@-1", 16, MPFR_RNDN); 532 i = mpfr_mul (x, y, z, MPFR_RNDU); 533 mpfr_set_str (y, "f.ffffffffffffffffffffffffffffffffffffffffffffffffffffffff" 534 "fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff" 535 "fffffffffffffffffffffffffffffffffffffffffffffffffffff337d23f07d8de1da5147a85c" 536 "dd3464e46068bd4d@-1", 16, MPFR_RNDN); 537 if (mpfr_cmp (x, y) || i <= 0) 538 { 539 printf ("Regression test (5) failed! (i=%d - expected 1)\n", i); 540 mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN); 541 printf ("\n"); 542 exit (1); 543 } 544 545 546 mpfr_set_prec (x, 50); 547 mpfr_set_prec (y, 40); 548 mpfr_set_prec (z, 53); 549 mpfr_set_str (y, "4.1ffffffff8", 16, MPFR_RNDN); 550 mpfr_set_str (z, "4.2000000ffe0000@-4", 16, MPFR_RNDN); 551 i = mpfr_mul (x, y, z, MPFR_RNDN); 552 if (mpfr_cmp_str (x, "1.104000041d6c0@-3", 16, MPFR_RNDN) != 0 553 || i <= 0) 554 { 555 printf ("Regression test (6) failed! (i=%d - expected 1)\nx=", i); 556 mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN); 557 printf ("\nMore prec="); 558 mpfr_set_prec (x, 93); 559 mpfr_mul (x, y, z, MPFR_RNDN); 560 mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN); 561 printf ("\n"); 562 exit (1); 563 } 564 565 mpfr_set_prec (x, 439); 566 mpfr_set_prec (y, 393); 567 mpfr_set_str (y, "-1.921fb54442d18469898cc51701b839a252049c1114cf98e804177d" 568 "4c76273644a29410f31c6809bbdf2a33679a748636600", 569 16, MPFR_RNDN); 570 i = mpfr_mul (x, y, y, MPFR_RNDU); 571 if (mpfr_cmp_str (x, "2.77a79937c8bbcb495b89b36602306b1c2159a8ff834288a19a08" 572 "84094f1cda3dc426da61174c4544a173de83c2500f8bfea2e0569e3698", 573 16, MPFR_RNDN) != 0 574 || i <= 0) 575 { 576 printf ("Regression test (7) failed! (i=%d - expected 1)\nx=", i); 577 mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN); 578 printf ("\n"); 579 exit (1); 580 } 581 582 mpfr_set_prec (x, 1023); 583 mpfr_set_prec (y, 1023); 584 mpfr_set_prec (z, 511); 585 mpfr_set_ui (x, 17, MPFR_RNDN); 586 mpfr_set_ui (y, 42, MPFR_RNDN); 587 i = mpfr_mul (z, x, y, MPFR_RNDN); 588 if (mpfr_cmp_ui (z, 17*42) != 0 || i != 0) 589 { 590 printf ("Regression test (8) failed! (i=%d - expected 0)\nz=", i); 591 mpfr_out_str (stdout, 16, 0, z, MPFR_RNDN); 592 printf ("\n"); 593 exit (1); 594 } 595 596 mpfr_clears (x, y, z, (mpfr_ptr) 0); 597 } 598 599 #define TEST_FUNCTION test_mul 600 #define TWO_ARGS 601 #define RAND_FUNCTION(x) mpfr_random2(x, MPFR_LIMB_SIZE (x), randlimb () % 100, RANDS) 602 #include "tgeneric.c" 603 604 /* multiplies x by 53-bit approximation of Pi */ 605 static int 606 mpfr_mulpi (mpfr_ptr y, mpfr_srcptr x, mpfr_rnd_t r) 607 { 608 mpfr_t z; 609 int inex; 610 611 mpfr_init2 (z, 53); 612 mpfr_set_str_binary (z, "11.001001000011111101101010100010001000010110100011"); 613 inex = mpfr_mul (y, x, z, r); 614 mpfr_clear (z); 615 return inex; 616 } 617 618 static void 619 valgrind20110503 (void) 620 { 621 mpfr_t a, b, c; 622 623 mpfr_init2 (a, 2); 624 mpfr_init2 (b, 2005); 625 mpfr_init2 (c, 2); 626 627 mpfr_set_ui (b, 5, MPFR_RNDN); 628 mpfr_nextabove (b); 629 mpfr_set_ui (c, 1, MPFR_RNDN); 630 mpfr_mul (a, b, c, MPFR_RNDZ); 631 /* After the call to mpfr_mulhigh_n, valgrind complains: 632 Conditional jump or move depends on uninitialised value(s) */ 633 634 mpfr_clears (a, b, c, (mpfr_ptr) 0); 635 } 636 637 static void 638 testall_rndf (mpfr_prec_t pmax) 639 { 640 mpfr_t a, b, c, d; 641 mpfr_prec_t pa, pb, pc; 642 643 for (pa = MPFR_PREC_MIN; pa <= pmax; pa++) 644 { 645 mpfr_init2 (a, pa); 646 mpfr_init2 (d, pa); 647 for (pb = MPFR_PREC_MIN; pb <= pmax; pb++) 648 { 649 mpfr_init2 (b, pb); 650 mpfr_set_ui (b, 1, MPFR_RNDN); 651 while (mpfr_cmp_ui (b, 2) < 0) 652 { 653 for (pc = MPFR_PREC_MIN; pc <= pmax; pc++) 654 { 655 mpfr_init2 (c, pc); 656 mpfr_set_ui (c, 1, MPFR_RNDN); 657 while (mpfr_cmp_ui (c, 2) < 0) 658 { 659 mpfr_mul (a, b, c, MPFR_RNDF); 660 mpfr_mul (d, b, c, MPFR_RNDD); 661 if (!mpfr_equal_p (a, d)) 662 { 663 mpfr_mul (d, b, c, MPFR_RNDU); 664 if (!mpfr_equal_p (a, d)) 665 { 666 printf ("Error: mpfr_mul(a,b,c,RNDF) does not " 667 "match RNDD/RNDU\n"); 668 printf ("b="); mpfr_dump (b); 669 printf ("c="); mpfr_dump (c); 670 printf ("a="); mpfr_dump (a); 671 exit (1); 672 } 673 } 674 mpfr_nextabove (c); 675 } 676 mpfr_clear (c); 677 } 678 mpfr_nextabove (b); 679 } 680 mpfr_clear (b); 681 } 682 mpfr_clear (a); 683 mpfr_clear (d); 684 } 685 } 686 687 /* Check underflow flag corresponds to *after* rounding. 688 * 689 * More precisely, we want to test mpfr_mul on inputs b and c such that 690 * EXP(b*c) < emin but EXP(round(b*c, p, rnd)) = emin. Thus an underflow 691 * must not be generated. 692 */ 693 static void 694 test_underflow (mpfr_prec_t pmax) 695 { 696 mpfr_exp_t emin; 697 mpfr_prec_t p; 698 mpfr_t a0, a, b, c; 699 int inex; 700 701 mpfr_init2 (a0, MPFR_PREC_MIN); 702 emin = mpfr_get_emin (); 703 mpfr_setmin (a0, emin); /* 0.5 * 2^emin */ 704 705 /* for RNDN, we want b*c < 0.5 * 2^emin but RNDN(b*c, p) = 0.5 * 2^emin, 706 thus b*c >= (0.5 - 1/4 * ulp_p(0.5)) * 2^emin */ 707 for (p = MPFR_PREC_MIN; p <= pmax; p++) 708 { 709 mpfr_init2 (a, p + 1); 710 mpfr_init2 (b, p + 10); 711 mpfr_init2 (c, p + 10); 712 do mpfr_urandomb (b, RANDS); while (MPFR_IS_ZERO (b)); 713 inex = mpfr_set_ui_2exp (a, 1, -1, MPFR_RNDZ); /* a = 0.5 */ 714 MPFR_ASSERTN (inex == 0); 715 mpfr_nextbelow (a); /* 0.5 - 1/2*ulp_{p+1}(0.5) = 0.5 - 1/4*ulp_p(0.5) */ 716 inex = mpfr_div (c, a, b, MPFR_RNDU); 717 /* 0.5 - 1/4 * ulp_p(0.5) = a <= b*c < 0.5 */ 718 mpfr_mul_2si (b, b, emin / 2, MPFR_RNDZ); 719 mpfr_mul_2si (c, c, (emin - 1) / 2, MPFR_RNDZ); 720 /* now (0.5 - 1/4 * ulp_p(0.5)) * 2^emin <= b*c < 0.5 * 2^emin, 721 thus b*c should be rounded to 0.5 * 2^emin */ 722 mpfr_set_prec (a, p); 723 mpfr_clear_underflow (); 724 mpfr_mul (a, b, c, MPFR_RNDN); 725 if (mpfr_underflow_p () || ! mpfr_equal_p (a, a0)) 726 { 727 printf ("Error, b*c incorrect or underflow flag incorrectly set" 728 " for emin=%" MPFR_EXP_FSPEC "d, rnd=%s\n", 729 (mpfr_eexp_t) emin, mpfr_print_rnd_mode (MPFR_RNDN)); 730 printf ("b="); mpfr_dump (b); 731 printf ("c="); mpfr_dump (c); 732 printf ("a="); mpfr_dump (a); 733 mpfr_set_prec (a, mpfr_get_prec (b) + mpfr_get_prec (c)); 734 mpfr_mul_2exp (b, b, 1, MPFR_RNDN); 735 inex = mpfr_mul (a, b, c, MPFR_RNDN); 736 MPFR_ASSERTN (inex == 0); 737 printf ("Exact 2*a="); mpfr_dump (a); 738 exit (1); 739 } 740 mpfr_clear (a); 741 mpfr_clear (b); 742 mpfr_clear (c); 743 } 744 745 /* for RNDU, we want b*c < 0.5*2^emin but RNDU(b*c, p) = 0.5*2^emin thus 746 b*c > (0.5 - 1/2 * ulp_p(0.5)) * 2^emin */ 747 for (p = MPFR_PREC_MIN; p <= pmax; p++) 748 { 749 mpfr_init2 (a, p); 750 mpfr_init2 (b, p + 10); 751 mpfr_init2 (c, p + 10); 752 do mpfr_urandomb (b, RANDS); while (MPFR_IS_ZERO (b)); 753 inex = mpfr_set_ui_2exp (a, 1, -1, MPFR_RNDZ); /* a = 0.5 */ 754 MPFR_ASSERTN (inex == 0); 755 mpfr_nextbelow (a); /* 0.5 - 1/2 * ulp_p(0.5) */ 756 inex = mpfr_div (c, a, b, MPFR_RNDU); 757 /* 0.5 - 1/2 * ulp_p(0.5) <= b*c < 0.5 */ 758 mpfr_mul_2si (b, b, emin / 2, MPFR_RNDZ); 759 mpfr_mul_2si (c, c, (emin - 1) / 2, MPFR_RNDZ); 760 if (inex == 0) 761 mpfr_nextabove (c); /* ensures b*c > (0.5 - 1/2 * ulp_p(0.5)) * 2^emin. 762 Warning: for p=1, 0.5 - 1/2 * ulp_p(0.5) 763 = 0.25, thus b*c > 2^(emin-2), which should 764 also be rounded up with p=1 to 0.5 * 2^emin 765 with an unbounded exponent range. */ 766 mpfr_clear_underflow (); 767 mpfr_mul (a, b, c, MPFR_RNDU); 768 if (mpfr_underflow_p () || ! mpfr_equal_p (a, a0)) 769 { 770 printf ("Error, b*c incorrect or underflow flag incorrectly set" 771 " for emin=%" MPFR_EXP_FSPEC "d, rnd=%s\n", 772 (mpfr_eexp_t) emin, mpfr_print_rnd_mode (MPFR_RNDU)); 773 printf ("b="); mpfr_dump (b); 774 printf ("c="); mpfr_dump (c); 775 printf ("a="); mpfr_dump (a); 776 mpfr_set_prec (a, mpfr_get_prec (b) + mpfr_get_prec (c)); 777 mpfr_mul_2exp (b, b, 1, MPFR_RNDN); 778 inex = mpfr_mul (a, b, c, MPFR_RNDN); 779 MPFR_ASSERTN (inex == 0); 780 printf ("Exact 2*a="); mpfr_dump (a); 781 exit (1); 782 } 783 mpfr_clear (a); 784 mpfr_clear (b); 785 mpfr_clear (c); 786 } 787 788 mpfr_clear (a0); 789 } 790 791 /* checks special case where no underflow should occur */ 792 static void 793 bug20161209 (void) 794 { 795 mpfr_exp_t emin; 796 mpfr_t x, y, z; 797 798 emin = mpfr_get_emin (); 799 set_emin (-1); 800 801 /* test for mpfr_mul_1 for 64-bit limb, mpfr_mul_2 for 32-bit limb */ 802 mpfr_init2 (x, 53); 803 mpfr_init2 (y, 53); 804 mpfr_init2 (z, 53); 805 mpfr_set_str_binary (x, "0.101000001E-1"); /* x = 321/2^10 */ 806 mpfr_set_str_binary (y, "0.110011000010100101111000011011000111011000001E-1"); 807 /* y = 28059810762433/2^46 */ 808 /* x * y = (2^53+1)/2^56 = 0.001...000[1]000..., and should round to 0.25 */ 809 mpfr_mul (z, x, y, MPFR_RNDN); 810 MPFR_ASSERTN(mpfr_cmp_ui_2exp (z, 1, -2) == 0); 811 812 /* test for mpfr_mul_2 for 64-bit limb */ 813 mpfr_set_prec (x, 65); 814 mpfr_set_prec (y, 65); 815 mpfr_set_prec (z, 65); 816 mpfr_set_str_binary (x, "0.101101000010010110100001E-1"); /* 11806113/2^25 */ 817 mpfr_set_str_binary (y, "0.101101011110010101011010001111111001100001E-1"); 818 /* y = 3124947910241/2^43 */ 819 /* x * y = (2^65+1)/2^68 = 0.001...000[1]000..., and should round to 0.25 */ 820 mpfr_mul (z, x, y, MPFR_RNDN); 821 MPFR_ASSERTN(mpfr_cmp_ui_2exp (z, 1, -2) == 0); 822 823 /* test for the generic code */ 824 mpfr_set_prec (x, 54); 825 mpfr_set_prec (y, 55); 826 mpfr_set_prec (z, 53); 827 mpfr_set_str_binary (x, "0.101000001E-1"); 828 mpfr_set_str_binary (y, "0.110011000010100101111000011011000111011000001E-1"); 829 mpfr_mul (z, x, y, MPFR_RNDN); 830 MPFR_ASSERTN(mpfr_cmp_ui_2exp (z, 1, -2) == 0); 831 832 /* another test for the generic code */ 833 mpfr_set_prec (x, 66); 834 mpfr_set_prec (y, 67); 835 mpfr_set_prec (z, 65); 836 mpfr_set_str_binary (x, "0.101101000010010110100001E-1"); 837 mpfr_set_str_binary (y, "0.101101011110010101011010001111111001100001E-1"); 838 mpfr_mul (z, x, y, MPFR_RNDN); 839 MPFR_ASSERTN(mpfr_cmp_ui_2exp (z, 1, -2) == 0); 840 841 mpfr_clear (x); 842 mpfr_clear (y); 843 mpfr_clear (z); 844 set_emin (emin); 845 } 846 847 /* test for case in mpfr_mul_1() where: 848 ax = __gmpfr_emin - 1 849 ap[0] == ~mask 850 rnd_mode = MPFR_RNDZ. 851 Whatever the values of rb and sb, we should round to zero (underflow). */ 852 static void 853 bug20161209a (void) 854 { 855 mpfr_exp_t emin; 856 mpfr_t x, y, z; 857 858 emin = mpfr_get_emin (); 859 set_emin (-1); 860 861 mpfr_init2 (x, 53); 862 mpfr_init2 (y, 53); 863 mpfr_init2 (z, 53); 864 865 /* case rb = sb = 0 */ 866 mpfr_set_str_binary (x, "0.11010010100110000110110011111E-1"); 867 mpfr_set_str_binary (y, "0.1001101110011000110100001"); 868 /* x = 441650591/2^30, y = 20394401/2^25 */ 869 /* x * y = (2^53-1)/2^55 = 0.00111...111[0]000..., and should round to 0 */ 870 mpfr_mul (z, x, y, MPFR_RNDZ); 871 MPFR_ASSERTN(mpfr_zero_p (z)); 872 873 /* case rb = 1, sb = 0 */ 874 mpfr_set_str_binary (x, "0.111111111000000000000000000111111111E-1"); 875 mpfr_set_str_binary (y, "0.1000000001000000001"); 876 /* x = 68585259519/2^37, y = 262657/2^19 */ 877 /* x * y = (2^54-1)/2^56 = 0.00111...111[1]000..., and should round to 0 */ 878 mpfr_mul (z, x, y, MPFR_RNDZ); 879 MPFR_ASSERTN(mpfr_zero_p (z)); 880 881 /* case rb = 0, sb = 1 */ 882 mpfr_set_str_binary (x, "0.110010011001011110001100100001000001E-1"); 883 mpfr_set_str_binary (y, "0.10100010100010111101"); 884 /* x = 541144371852^37, y = 665789/2^20 */ 885 /* x * y = (2^55-3)/2^57 = 0.00111...111[0]100..., and should round to 0 */ 886 mpfr_mul (z, x, y, MPFR_RNDZ); 887 MPFR_ASSERTN(mpfr_zero_p (z)); 888 889 /* case rb = sb = 1 */ 890 mpfr_set_str_binary (x, "0.10100110001001001010001111110010100111E-1"); 891 mpfr_set_str_binary (y, "0.110001010011101001"); 892 /* x = 178394823847/2^39, y = 201961/2^18 */ 893 /* x * y = (2^55-1)/2^57 = 0.00111...111[1]100..., and should round to 0 */ 894 mpfr_mul (z, x, y, MPFR_RNDZ); 895 MPFR_ASSERTN(mpfr_zero_p (z)); 896 897 /* similar test for mpfr_mul_2 (we only check rb = sb = 1 here) */ 898 mpfr_set_prec (x, 65); 899 mpfr_set_prec (y, 65); 900 mpfr_set_prec (z, 65); 901 /* 2^67-1 = 193707721 * 761838257287 */ 902 mpfr_set_str_binary (x, "0.1011100010111011111011001001E-1"); 903 mpfr_set_str_binary (y, "0.1011000101100001000110010100010010000111"); 904 mpfr_mul (z, x, y, MPFR_RNDZ); 905 MPFR_ASSERTN(mpfr_zero_p (z)); 906 907 mpfr_clear (x); 908 mpfr_clear (y); 909 mpfr_clear (z); 910 set_emin (emin); 911 } 912 913 /* bug for RNDF */ 914 static void 915 bug20170602 (void) 916 { 917 mpfr_t x, u, y, yd, yu; 918 919 mpfr_init2 (x, 493); 920 mpfr_init2 (u, 493); 921 mpfr_init2 (y, 503); 922 mpfr_init2 (yd, 503); 923 mpfr_init2 (yu, 503); 924 mpfr_set_str_binary (x, "0.1111100000000000001111111110000000001111111111111000000000000000000011111111111111111111111000000000000000000001111111111111111111111111111111111111111000000000011111111111111111111000000000011111111111111000000000000001110000000000000000000000000000000000000000011111111111110011111111111100000000000000011111111111111111110000000011111111111111111110011111111111110000000000001111111111111111000000000000000000000000000000000000111111111111111111111111111111011111111111111111111111111111100E44"); 925 mpfr_set_str_binary (u, "0.1110000000000000001111111111111111111111111111111111111000000000111111111111111111111111111111000000000000000000001111111000000000000000011111111111111111111111111111111111111111111111111111111000000000000000011111111111111000000011111111111111111110000000000000001111111111111111111111111111111111111110000000000001111111111111111111111111111111111111000000000000000000000000000000000001111111111111000000000000000000001111111111100000000000000011111111111111111111111111111111111111111111111E35"); 926 mpfr_mul (y, x, u, MPFR_RNDF); 927 mpfr_mul (yd, x, u, MPFR_RNDD); 928 mpfr_mul (yu, x, u, MPFR_RNDU); 929 if (mpfr_cmp (y, yd) != 0 && mpfr_cmp (y, yu) != 0) 930 { 931 printf ("RNDF is neither RNDD nor RNDU\n"); 932 printf ("x"); mpfr_dump (x); 933 printf ("u"); mpfr_dump (u); 934 printf ("y(RNDF)="); mpfr_dump (y); 935 printf ("y(RNDD)="); mpfr_dump (yd); 936 printf ("y(RNDU)="); mpfr_dump (yu); 937 exit (1); 938 } 939 mpfr_clear (x); 940 mpfr_clear (u); 941 mpfr_clear (y); 942 mpfr_clear (yd); 943 mpfr_clear (yu); 944 } 945 946 /* Test for 1 to 3 limbs. */ 947 static void 948 small_prec (void) 949 { 950 mpfr_exp_t emin, emax; 951 mpfr_t x, y, z1, z2, zz; 952 int xq, yq, zq; 953 mpfr_rnd_t rnd; 954 mpfr_flags_t flags1, flags2; 955 int inex1, inex2; 956 int i, j, r, s, ediff; 957 958 emin = mpfr_get_emin (); 959 emax = mpfr_get_emax (); 960 961 /* The mpfr_mul implementation doesn't extend the exponent range, 962 so that it is OK to reduce it here for the test to make sure that 963 mpfr_mul_2si can be used. */ 964 set_emin (-1000); 965 set_emax (1000); 966 967 mpfr_inits2 (3 * GMP_NUMB_BITS, x, y, z1, z2, (mpfr_ptr) 0); 968 mpfr_init2 (zz, 6 * GMP_NUMB_BITS); 969 for (i = 0; i < 3; i++) 970 for (j = 0; j < 10000; j++) 971 { 972 xq = i * GMP_NUMB_BITS + 1 + randlimb () % GMP_NUMB_BITS; 973 mpfr_set_prec (x, xq); 974 yq = i * GMP_NUMB_BITS + 1 + randlimb () % GMP_NUMB_BITS; 975 mpfr_set_prec (y, yq); 976 zq = i * GMP_NUMB_BITS + 1 + randlimb () % (GMP_NUMB_BITS-1); 977 mpfr_set_prec (z1, zq); 978 mpfr_set_prec (z2, zq); 979 s = r = randlimb () & 0x7f; 980 do mpfr_urandomb (x, RANDS); while (mpfr_zero_p (x)); 981 if (s & 1) 982 mpfr_neg (x, x, MPFR_RNDN); 983 s >>= 1; 984 if (s & 1) 985 { 986 do mpfr_urandomb (y, RANDS); while (mpfr_zero_p (y)); 987 } 988 else 989 { 990 mpfr_ui_div (y, 1, x, MPFR_RNDN); 991 mpfr_set_exp (y, 0); 992 } 993 s >>= 1; 994 if (s & 1) 995 mpfr_neg (y, y, MPFR_RNDN); 996 s >>= 1; 997 rnd = RND_RAND_NO_RNDF (); 998 inex1 = mpfr_mul (zz, x, y, MPFR_RNDN); 999 MPFR_ASSERTN (inex1 == 0); 1000 if (s == 0) 1001 { 1002 ediff = __gmpfr_emin - MPFR_EXP (x); 1003 mpfr_set_exp (x, __gmpfr_emin); 1004 } 1005 else if (s == 1) 1006 { 1007 ediff = __gmpfr_emax - MPFR_EXP (x) + 1; 1008 mpfr_set_exp (x, __gmpfr_emax); 1009 mpfr_mul_2ui (y, y, 1, MPFR_RNDN); 1010 } 1011 else 1012 ediff = 0; 1013 mpfr_clear_flags (); 1014 inex1 = mpfr_mul_2si (z1, zz, ediff, rnd); 1015 flags1 = __gmpfr_flags; 1016 mpfr_clear_flags (); 1017 inex2 = mpfr_mul (z2, x, y, rnd); 1018 flags2 = __gmpfr_flags; 1019 if (!(mpfr_equal_p (z1, z2) && 1020 SAME_SIGN (inex1, inex2) && 1021 flags1 == flags2)) 1022 { 1023 printf ("Error in small_prec on i = %d, j = %d\n", i, j); 1024 printf ("r = 0x%x, xq = %d, yq = %d, zq = %d, rnd = %s\n", 1025 r, xq, yq, zq, mpfr_print_rnd_mode (rnd)); 1026 printf ("x = "); 1027 mpfr_dump (x); 1028 printf ("y = "); 1029 mpfr_dump (y); 1030 printf ("ediff = %d\n", ediff); 1031 printf ("zz = "); 1032 mpfr_dump (zz); 1033 printf ("Expected "); 1034 mpfr_dump (z1); 1035 printf ("with inex = %d and flags =", inex1); 1036 flags_out (flags1); 1037 printf ("Got "); 1038 mpfr_dump (z2); 1039 printf ("with inex = %d and flags =", inex2); 1040 flags_out (flags2); 1041 exit (1); 1042 } 1043 } 1044 mpfr_clears (x, y, z1, z2, zz, (mpfr_ptr) 0); 1045 1046 set_emin (emin); 1047 set_emax (emax); 1048 } 1049 1050 /* check ax < __gmpfr_emin with rnd_mode == MPFR_RNDN, rb = 0 and sb <> 0 */ 1051 static void 1052 test_underflow2 (void) 1053 { 1054 mpfr_t x, y, z; 1055 mpfr_exp_t emin; 1056 1057 emin = mpfr_get_emin (); 1058 mpfr_set_emin (0); 1059 1060 mpfr_init2 (x, 24); 1061 mpfr_init2 (y, 24); 1062 mpfr_init2 (z, 24); 1063 1064 mpfr_set_ui_2exp (x, 12913, -14, MPFR_RNDN); 1065 mpfr_set_ui_2exp (y, 5197, -13, MPFR_RNDN); 1066 mpfr_clear_underflow (); 1067 /* x*y = 0.0111111111111111111111111[01] thus underflow */ 1068 mpfr_mul (z, y, x, MPFR_RNDN); 1069 MPFR_ASSERTN(mpfr_cmp_ui_2exp (z, 1, -1) == 0); 1070 MPFR_ASSERTN(mpfr_underflow_p ()); 1071 1072 mpfr_set_prec (x, 65); 1073 mpfr_set_prec (y, 65); 1074 mpfr_set_prec (z, 65); 1075 mpfr_set_str_binary (x, "0.10011101000110111011111100101001111"); 1076 mpfr_set_str_binary (y, "0.110100001001000111000011011110011"); 1077 mpfr_clear_underflow (); 1078 /* x*y = 0.011111111111111111111111111111111111111111111111111111111111111111[01] thus underflow */ 1079 mpfr_mul (z, y, x, MPFR_RNDN); 1080 MPFR_ASSERTN(mpfr_cmp_ui_2exp (z, 1, -1) == 0); 1081 MPFR_ASSERTN(mpfr_underflow_p ()); 1082 1083 mpfr_clear (y); 1084 mpfr_clear (x); 1085 mpfr_clear (z); 1086 1087 mpfr_set_emin (emin); 1088 } 1089 1090 int 1091 main (int argc, char *argv[]) 1092 { 1093 tests_start_mpfr (); 1094 1095 testall_rndf (9); 1096 check_nans (); 1097 check_exact (); 1098 check_float (); 1099 1100 check53("6.9314718055994530941514e-1", "0.0", MPFR_RNDZ, "0.0"); 1101 check53("0.0", "6.9314718055994530941514e-1", MPFR_RNDZ, "0.0"); 1102 check_sign(); 1103 check53("-4.165000000e4", "-0.00004801920768307322868063274915", MPFR_RNDN, 1104 "2.0"); 1105 check53("2.71331408349172961467e-08", "-6.72658901114033715233e-165", 1106 MPFR_RNDZ, "-1.8251348697787782844e-172"); 1107 check53("2.71331408349172961467e-08", "-6.72658901114033715233e-165", 1108 MPFR_RNDA, "-1.8251348697787786e-172"); 1109 check53("0.31869277231188065", "0.88642843322303122", MPFR_RNDZ, 1110 "2.8249833483992453642e-1"); 1111 check("8.47622108205396074254e-01", "3.24039313247872939883e-01", MPFR_RNDU, 1112 28, 45, 2, "0.375"); 1113 check("8.47622108205396074254e-01", "3.24039313247872939883e-01", MPFR_RNDA, 1114 28, 45, 2, "0.375"); 1115 check("2.63978122803639081440e-01", "6.8378615379333496093e-1", MPFR_RNDN, 1116 34, 23, 31, "0.180504585267044603"); 1117 check("1.0", "0.11835170935876249132", MPFR_RNDU, 6, 41, 36, 1118 "0.1183517093595583"); 1119 check53("67108865.0", "134217729.0", MPFR_RNDN, "9.007199456067584e15"); 1120 check("1.37399642157394197284e-01", "2.28877275604219221350e-01", MPFR_RNDN, 1121 49, 15, 32, "0.0314472340833162888"); 1122 check("4.03160720978664954828e-01", "5.854828e-1" 1123 /*"5.85483042917246621073e-01"*/, MPFR_RNDZ, 1124 51, 22, 32, "0.2360436821472831"); 1125 check("3.90798504668055102229e-14", "9.85394674650308388664e-04", MPFR_RNDN, 1126 46, 22, 12, "0.385027296503914762e-16"); 1127 check("4.58687081072827851358e-01", "2.20543551472118792844e-01", MPFR_RNDN, 1128 49, 3, 2, "0.09375"); 1129 check_max(); 1130 check_min(); 1131 small_prec (); 1132 1133 check_regression (); 1134 test_generic (MPFR_PREC_MIN, 500, 100); 1135 1136 data_check ("data/mulpi", mpfr_mulpi, "mpfr_mulpi"); 1137 1138 valgrind20110503 (); 1139 test_underflow (128); 1140 bug20161209 (); 1141 bug20161209a (); 1142 bug20170602 (); 1143 test_underflow2 (); 1144 1145 tests_end_mpfr (); 1146 return 0; 1147 } 1148