1 /* Test file for mpfr_sin_cos. 2 3 Copyright 2000-2004, 2006-2016 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 <stdio.h> 24 #include <stdlib.h> 25 26 #include "mpfr-test.h" 27 28 static void 29 large_test (char *X, int prec, int N) 30 { 31 int i; 32 mpfr_t x, s, c; 33 34 mpfr_init2 (x, prec); 35 mpfr_init2 (s, prec); 36 mpfr_init2 (c, prec); 37 mpfr_set_str (x, X, 10, MPFR_RNDN); 38 39 for (i = 0; i < N; i++) 40 mpfr_sin_cos (s, c, x, MPFR_RNDN); 41 42 mpfr_clear (x); 43 mpfr_clear (s); 44 mpfr_clear (c); 45 } 46 47 static void 48 check53 (const char *xs, const char *sin_xs, const char *cos_xs, 49 mpfr_rnd_t rnd_mode) 50 { 51 mpfr_t xx, s, c; 52 53 mpfr_inits2 (53, xx, s, c, (mpfr_ptr) 0); 54 mpfr_set_str1 (xx, xs); /* should be exact */ 55 mpfr_sin_cos (s, c, xx, rnd_mode); 56 if (mpfr_cmp_str1 (s, sin_xs)) 57 { 58 printf ("mpfr_sin_cos failed for x=%s, rnd=%s\n", 59 xs, mpfr_print_rnd_mode (rnd_mode)); 60 printf ("mpfr_sin_cos gives sin(x)="); 61 mpfr_out_str(stdout, 10, 0, s, MPFR_RNDN); 62 printf(", expected %s\n", sin_xs); 63 exit (1); 64 } 65 if (mpfr_cmp_str1 (c, cos_xs)) 66 { 67 printf ("mpfr_sin_cos failed for x=%s, rnd=%s\n", 68 xs, mpfr_print_rnd_mode (rnd_mode)); 69 printf ("mpfr_sin_cos gives cos(x)="); 70 mpfr_out_str(stdout, 10, 0, c, MPFR_RNDN); 71 printf(", expected %s\n", cos_xs); 72 exit (1); 73 } 74 mpfr_clears (xx, s, c, (mpfr_ptr) 0); 75 } 76 77 static void 78 check53sin (const char *xs, const char *sin_xs, mpfr_rnd_t rnd_mode) 79 { 80 mpfr_t xx, s, c; 81 82 mpfr_inits2 (53, xx, s, c, (mpfr_ptr) 0); 83 mpfr_set_str1 (xx, xs); /* should be exact */ 84 mpfr_sin_cos (s, c, xx, rnd_mode); 85 if (mpfr_cmp_str1 (s, sin_xs)) 86 { 87 printf ("mpfr_sin_cos failed for x=%s, rnd=%s\n", 88 xs, mpfr_print_rnd_mode (rnd_mode)); 89 printf ("mpfr_sin_cos gives sin(x)="); 90 mpfr_out_str(stdout, 10, 0, s, MPFR_RNDN); 91 printf(", expected %s\n", sin_xs); 92 exit (1); 93 } 94 mpfr_clears (xx, s, c, (mpfr_ptr) 0); 95 } 96 97 static void 98 check53cos (const char *xs, const char *cos_xs, mpfr_rnd_t rnd_mode) 99 { 100 mpfr_t xx, c, s; 101 102 mpfr_inits2 (53, xx, s, c, (mpfr_ptr) 0); 103 mpfr_set_str1 (xx, xs); /* should be exact */ 104 mpfr_sin_cos (s, c, xx, rnd_mode); 105 if (mpfr_cmp_str1 (c, cos_xs)) 106 { 107 printf ("mpfr_sin_cos failed for x=%s, rnd=%s\n", 108 xs, mpfr_print_rnd_mode (rnd_mode)); 109 printf ("mpfr_sin_cos gives cos(x)="); 110 mpfr_out_str(stdout, 10, 0, c, MPFR_RNDN); 111 printf(", expected %s\n", cos_xs); 112 exit (1); 113 } 114 mpfr_clears (xx, s, c, (mpfr_ptr) 0); 115 } 116 117 static void 118 check_nans (void) 119 { 120 mpfr_t x, s, c; 121 122 mpfr_init2 (x, 123L); 123 mpfr_init2 (s, 123L); 124 mpfr_init2 (c, 123L); 125 126 /* sin(NaN)==NaN, cos(NaN)==NaN */ 127 mpfr_set_nan (x); 128 mpfr_sin_cos (s, c, x, MPFR_RNDN); 129 MPFR_ASSERTN (mpfr_nan_p (s)); 130 MPFR_ASSERTN (mpfr_nan_p (c)); 131 132 /* sin(+Inf)==NaN, cos(+Inf)==NaN */ 133 mpfr_set_inf (x, 1); 134 mpfr_sin_cos (s, c, x, MPFR_RNDN); 135 MPFR_ASSERTN (mpfr_nan_p (s)); 136 MPFR_ASSERTN (mpfr_nan_p (c)); 137 138 /* sin(-Inf)==NaN, cos(-Inf)==NaN */ 139 mpfr_set_inf (x, -1); 140 mpfr_sin_cos (s, c, x, MPFR_RNDN); 141 MPFR_ASSERTN (mpfr_nan_p (s)); 142 MPFR_ASSERTN (mpfr_nan_p (c)); 143 144 /* check zero */ 145 mpfr_set_ui (x, 0, MPFR_RNDN); 146 mpfr_sin_cos (s, c, x, MPFR_RNDN); 147 MPFR_ASSERTN (mpfr_cmp_ui (s, 0) == 0 && MPFR_IS_POS (s)); 148 MPFR_ASSERTN (mpfr_cmp_ui (c, 1) == 0); 149 mpfr_neg (x, x, MPFR_RNDN); 150 mpfr_sin_cos (s, c, x, MPFR_RNDN); 151 MPFR_ASSERTN (mpfr_cmp_ui (s, 0) == 0 && MPFR_IS_NEG (s)); 152 MPFR_ASSERTN (mpfr_cmp_ui (c, 1) == 0); 153 154 /* coverage test */ 155 mpfr_set_prec (x, 2); 156 mpfr_set_ui (x, 4, MPFR_RNDN); 157 mpfr_set_prec (s, 2); 158 mpfr_set_prec (c, 2); 159 mpfr_sin_cos (s, c, x, MPFR_RNDN); 160 MPFR_ASSERTN (mpfr_cmp_si_2exp (s, -3, -2) == 0); 161 MPFR_ASSERTN (mpfr_cmp_si_2exp (c, -3, -2) == 0); 162 163 mpfr_clear (x); 164 mpfr_clear (s); 165 mpfr_clear (c); 166 } 167 168 static void 169 overflowed_sin_cos0 (void) 170 { 171 mpfr_t x, y, z; 172 int emax, inex, rnd, err = 0; 173 mpfr_exp_t old_emax; 174 175 old_emax = mpfr_get_emax (); 176 177 mpfr_init2 (x, 8); 178 mpfr_init2 (y, 8); 179 mpfr_init2 (z, 8); 180 181 for (emax = -1; emax <= 0; emax++) 182 { 183 mpfr_set_ui_2exp (z, 1, emax, MPFR_RNDN); 184 mpfr_nextbelow (z); 185 set_emax (emax); /* 1 is not representable. */ 186 /* and if emax < 0, 1 - eps is not representable either. */ 187 RND_LOOP (rnd) 188 { 189 mpfr_set_si (x, 0, MPFR_RNDN); 190 mpfr_neg (x, x, MPFR_RNDN); 191 mpfr_clear_flags (); 192 inex = mpfr_sin_cos (x, y, x, (mpfr_rnd_t) rnd); 193 if (! mpfr_overflow_p ()) 194 { 195 printf ("Error in overflowed_sin_cos0 (rnd = %s):\n" 196 " The overflow flag is not set.\n", 197 mpfr_print_rnd_mode ((mpfr_rnd_t) rnd)); 198 err = 1; 199 } 200 if (! (mpfr_zero_p (x) && MPFR_SIGN (x) < 0)) 201 { 202 printf ("Error in overflowed_sin_cos0 (rnd = %s):\n" 203 " Got sin = ", mpfr_print_rnd_mode ((mpfr_rnd_t) rnd)); 204 mpfr_print_binary (x); 205 printf (" instead of -0.\n"); 206 err = 1; 207 } 208 if (rnd == MPFR_RNDZ || rnd == MPFR_RNDD) 209 { 210 if (inex == 0) 211 { 212 printf ("Error in overflowed_sin_cos0 (rnd = %s):\n" 213 " The inexact value must be non-zero.\n", 214 mpfr_print_rnd_mode ((mpfr_rnd_t) rnd)); 215 err = 1; 216 } 217 if (! mpfr_equal_p (y, z)) 218 { 219 printf ("Error in overflowed_sin_cos0 (rnd = %s):\n" 220 " Got cos = ", 221 mpfr_print_rnd_mode ((mpfr_rnd_t) rnd)); 222 mpfr_print_binary (y); 223 printf (" instead of 0.11111111E%d.\n", emax); 224 err = 1; 225 } 226 } 227 else 228 { 229 if (inex == 0) 230 { 231 printf ("Error in overflowed_sin_cos0 (rnd = %s):\n" 232 " The inexact value must be non-zero.\n", 233 mpfr_print_rnd_mode ((mpfr_rnd_t) rnd)); 234 err = 1; 235 } 236 if (! (mpfr_inf_p (y) && MPFR_SIGN (y) > 0)) 237 { 238 printf ("Error in overflowed_sin_cos0 (rnd = %s):\n" 239 " Got cos = ", 240 mpfr_print_rnd_mode ((mpfr_rnd_t) rnd)); 241 mpfr_print_binary (y); 242 printf (" instead of +Inf.\n"); 243 err = 1; 244 } 245 } 246 } 247 set_emax (old_emax); 248 } 249 250 if (err) 251 exit (1); 252 mpfr_clear (x); 253 mpfr_clear (y); 254 mpfr_clear (z); 255 } 256 257 static void 258 tiny (void) 259 { 260 mpfr_t x, s, c; 261 int i, inex; 262 263 mpfr_inits2 (64, x, s, c, (mpfr_ptr) 0); 264 265 for (i = -1; i <= 1; i += 2) 266 { 267 mpfr_set_si (x, i, MPFR_RNDN); 268 mpfr_set_exp (x, mpfr_get_emin ()); 269 inex = mpfr_sin_cos (s, c, x, MPFR_RNDN); 270 MPFR_ASSERTN (inex != 0); 271 MPFR_ASSERTN (mpfr_equal_p (s, x)); 272 MPFR_ASSERTN (!mpfr_nan_p (c) && mpfr_cmp_ui (c, 1) == 0); 273 } 274 275 mpfr_clears (x, s, c, (mpfr_ptr) 0); 276 } 277 278 /* bug found in nightly tests */ 279 static void 280 test20071214 (void) 281 { 282 mpfr_t a, b; 283 int inex; 284 285 mpfr_init2 (a, 4); 286 mpfr_init2 (b, 4); 287 288 mpfr_set_ui_2exp (a, 3, -4, MPFR_RNDN); 289 inex = mpfr_sin_cos (a, b, a, MPFR_RNDD); 290 MPFR_ASSERTN(mpfr_cmp_ui_2exp (a, 11, -6) == 0); 291 MPFR_ASSERTN(mpfr_cmp_ui_2exp (b, 15, -4) == 0); 292 MPFR_ASSERTN(inex == 10); 293 294 mpfr_set_ui_2exp (a, 3, -4, MPFR_RNDN); 295 inex = mpfr_sin_cos (a, b, a, MPFR_RNDU); 296 MPFR_ASSERTN(mpfr_cmp_ui_2exp (a, 3, -4) == 0); 297 MPFR_ASSERTN(mpfr_cmp_ui (b, 1) == 0); 298 MPFR_ASSERTN(inex == 5); 299 300 mpfr_set_ui_2exp (a, 3, -4, MPFR_RNDN); 301 inex = mpfr_sin_cos (a, b, a, MPFR_RNDN); 302 MPFR_ASSERTN(mpfr_cmp_ui_2exp (a, 3, -4) == 0); 303 MPFR_ASSERTN(mpfr_cmp_ui (b, 1) == 0); 304 MPFR_ASSERTN(inex == 5); 305 306 mpfr_clear (a); 307 mpfr_clear (b); 308 } 309 310 /* check that mpfr_sin_cos and test_mpfr_sincos_fast agree */ 311 static void 312 test_mpfr_sincos_fast (void) 313 { 314 mpfr_t x, y, z, yref, zref, h; 315 mpfr_prec_t p = 1000; 316 int i, inex, inexref; 317 mpfr_rnd_t r; 318 319 mpfr_init2 (x, p); 320 mpfr_init2 (y, p); 321 mpfr_init2 (z, p); 322 mpfr_init2 (yref, p); 323 mpfr_init2 (zref, p); 324 mpfr_init2 (h, p); 325 mpfr_set_ui (x, 0, MPFR_RNDN); 326 /* we generate a random value x, compute sin(x) and cos(x) with both 327 mpfr_sin_cos and mpfr_sincos_fast, and check the values and the flags 328 agree */ 329 for (i = 0; i < 100; i++) 330 { 331 mpfr_urandomb (h, RANDS); 332 mpfr_add (x, x, h, MPFR_RNDN); 333 r = RND_RAND (); 334 inexref = mpfr_sin_cos (yref, zref, x, r); 335 inex = mpfr_sincos_fast (y, z, x, r); 336 if (mpfr_cmp (y, yref)) 337 { 338 printf ("mpfr_sin_cos and mpfr_sincos_fast disagree\n"); 339 printf ("x="); mpfr_dump (x); 340 printf ("rnd=%s\n", mpfr_print_rnd_mode (r)); 341 printf ("yref="); mpfr_dump (yref); 342 printf ("y="); mpfr_dump (y); 343 exit (1); 344 } 345 if (mpfr_cmp (z, zref)) 346 { 347 printf ("mpfr_sin_cos and mpfr_sincos_fast disagree\n"); 348 printf ("x="); mpfr_dump (x); 349 printf ("rnd=%s\n", mpfr_print_rnd_mode (r)); 350 printf ("zref="); mpfr_dump (zref); 351 printf ("z="); mpfr_dump (z); 352 exit (1); 353 } 354 if (inex != inexref) 355 { 356 printf ("mpfr_sin_cos and mpfr_sincos_fast disagree\n"); 357 printf ("x="); mpfr_dump (x); 358 printf ("rnd=%s\n", mpfr_print_rnd_mode (r)); 359 printf ("inexref=%d inex=%d\n", inexref, inex); 360 exit (1); 361 } 362 } 363 mpfr_clear (x); 364 mpfr_clear (y); 365 mpfr_clear (z); 366 mpfr_clear (yref); 367 mpfr_clear (zref); 368 mpfr_clear (h); 369 } 370 371 static void 372 bug20091007 (void) 373 { 374 mpfr_t x, y, z, yref, zref; 375 mpfr_prec_t p = 1000; 376 int inex, inexref; 377 mpfr_rnd_t r = MPFR_RNDZ; 378 379 mpfr_init2 (x, p); 380 mpfr_init2 (y, p); 381 mpfr_init2 (z, p); 382 mpfr_init2 (yref, p); 383 mpfr_init2 (zref, p); 384 385 mpfr_set_str (x, "1.9ecdc22ba77a5ab2560f7e84289e2a328906f47377ea3fd4c82d1bb2f13ee05c032cffc1933eadab7b0a5498e03e3bd0508968e59c25829d97a0b54f20cd4662c8dfffa54e714de41fc8ee3e0e0b244d110a194db05b70022b7d77f88955d415b09f17dd404576098dc51a583a3e49c35839551646e880c7eb790a01a4@1", 16, MPFR_RNDN); 386 inexref = mpfr_sin_cos (yref, zref, x, r); 387 inex = mpfr_sincos_fast (y, z, x, r); 388 389 if (mpfr_cmp (y, yref)) 390 { 391 printf ("mpfr_sin_cos and mpfr_sincos_fast disagree (bug20091007)\n"); 392 printf ("yref="); mpfr_dump (yref); 393 printf ("y="); mpfr_dump (y); 394 exit (1); 395 } 396 if (mpfr_cmp (z, zref)) 397 { 398 printf ("mpfr_sin_cos and mpfr_sincos_fast disagree (bug20091007)\n"); 399 printf ("zref="); mpfr_dump (zref); 400 printf ("z="); mpfr_dump (z); 401 exit (1); 402 } 403 if (inex != inexref) 404 { 405 printf ("mpfr_sin_cos and mpfr_sincos_fast disagree (bug20091007)\n"); 406 printf ("inexref=%d inex=%d\n", inexref, inex); 407 exit (1); 408 } 409 410 mpfr_clear (x); 411 mpfr_clear (y); 412 mpfr_clear (z); 413 mpfr_clear (yref); 414 mpfr_clear (zref); 415 } 416 417 /* Note: with the sin_cos.c code before r6507, the disagreement occurs 418 only on the return ("inexact") value, which is new in r6444. */ 419 static void 420 bug20091008 (void) 421 { 422 mpfr_t x, y, z, yref, zref; 423 mpfr_prec_t p = 1000; 424 int inex, inexref; 425 mpfr_rnd_t r = MPFR_RNDN; 426 427 mpfr_init2 (x, p); 428 mpfr_init2 (y, p); 429 mpfr_init2 (z, p); 430 mpfr_init2 (yref, p); 431 mpfr_init2 (zref, p); 432 433 mpfr_set_str (x, "c.91813724e28ef6a711d33e6505984699daef7fe93636c1ed5d0168bc96989cc6802f7f9e405c902ec62fb90cd39c9d21084c8ad8b5af4c4aa87bf402e2e4a78e6fe1ffeb6dbbbdbbc2983c196c518966ccc1e094ed39ee77984ef2428069d65de37928e75247edbe7007245e682616b5ebbf05f2fdefc74ad192024f10", 16, MPFR_RNDN); 434 inexref = mpfr_sin_cos (yref, zref, x, r); 435 inex = mpfr_sincos_fast (y, z, x, r); 436 437 if (mpfr_cmp (y, yref)) 438 { 439 printf ("mpfr_sin_cos and mpfr_sincos_fast disagree (bug20091008)\n"); 440 printf ("yref="); mpfr_dump (yref); 441 printf ("y="); mpfr_dump (y); 442 exit (1); 443 } 444 if (mpfr_cmp (z, zref)) 445 { 446 printf ("mpfr_sin_cos and mpfr_sincos_fast disagree (bug20091008)\n"); 447 printf ("zref="); mpfr_dump (zref); 448 printf ("z="); mpfr_dump (z); 449 exit (1); 450 } 451 /* sin(x) is rounded up, cos(x) is rounded up too, thus we should get 5 452 for the return value */ 453 if (inex != inexref) 454 { 455 printf ("mpfr_sin_cos and mpfr_sincos_fast disagree (bug20091008)\n"); 456 printf ("inexref=%d inex=%d (5 expected)\n", inexref, inex); 457 exit (1); 458 } 459 460 mpfr_clear (x); 461 mpfr_clear (y); 462 mpfr_clear (z); 463 mpfr_clear (yref); 464 mpfr_clear (zref); 465 } 466 467 static void 468 bug20091013 (void) 469 { 470 mpfr_t x, y, z, yref, zref; 471 mpfr_prec_t p = 1000; 472 int inex, inexref; 473 mpfr_rnd_t r = MPFR_RNDN; 474 475 mpfr_init2 (x, p); 476 mpfr_init2 (y, p); 477 mpfr_init2 (z, p); 478 mpfr_init2 (yref, p); 479 mpfr_init2 (zref, p); 480 481 mpfr_set_str (x, "3.240ff3fdcb1ee7cd667b96287593ae24e20fb63ed7c2d5bf4bd0f2cc5509283b04e7628e66382605f14ed5967cef15296041539a1bdaa626c777c7fbb6f2068414759b78cee14f37848689b3a170f583656be4e0837f464d8210556a3a822d4ecfdd59f4e0d5fdb76bf7e15b8a57234e2160b98e14c17bbdf27c4643b8@1", 16, MPFR_RNDN); 482 inexref = mpfr_sin_cos (yref, zref, x, r); 483 inex = mpfr_sincos_fast (y, z, x, r); 484 485 if (mpfr_cmp (y, yref)) 486 { 487 printf ("mpfr_sin_cos and mpfr_sincos_fast disagree (bug20091013)\n"); 488 printf ("yref="); mpfr_dump (yref); 489 printf ("y="); mpfr_dump (y); 490 exit (1); 491 } 492 if (mpfr_cmp (z, zref)) 493 { 494 printf ("mpfr_sin_cos and mpfr_sincos_fast disagree (bug20091013)\n"); 495 printf ("zref="); mpfr_dump (zref); 496 printf ("z="); mpfr_dump (z); 497 exit (1); 498 } 499 /* sin(x) is rounded down and cos(x) is rounded down, thus we should get 500 2+4*2 = 10 as return value */ 501 if (inex != inexref) 502 { 503 printf ("mpfr_sin_cos and mpfr_sincos_fast disagree (bug20091013)\n"); 504 printf ("inexref=%d inex=%d (10 expected)\n", inexref, inex); 505 exit (1); 506 } 507 508 mpfr_clear (x); 509 mpfr_clear (y); 510 mpfr_clear (z); 511 mpfr_clear (yref); 512 mpfr_clear (zref); 513 } 514 515 /* Bug reported by Laurent Fousse for the 2.4 branch. 516 No problem in the trunk. 517 https://sympa.inria.fr/sympa/arc/mpfr/2009-11/msg00044.html */ 518 static void 519 bug20091122 (void) 520 { 521 mpfr_t x, y, z, yref, zref; 522 mpfr_prec_t p = 3; 523 mpfr_rnd_t r = MPFR_RNDN; 524 525 mpfr_init2 (x, 5); 526 mpfr_init2 (y, p); 527 mpfr_init2 (z, p); 528 mpfr_init2 (yref, p); 529 mpfr_init2 (zref, p); 530 531 mpfr_set_str (x, "0.11111E49", 2, MPFR_RNDN); 532 mpfr_sin_cos (yref, zref, x, r); 533 534 mpfr_sin (y, x, r); 535 mpfr_cos (z, x, r); 536 537 if (! mpfr_equal_p (y, yref)) 538 { 539 printf ("mpfr_sin_cos and mpfr_sin disagree (bug20091122)\n"); 540 printf ("yref = "); mpfr_dump (yref); 541 printf ("y = "); mpfr_dump (y); 542 exit (1); 543 } 544 if (! mpfr_equal_p (z, zref)) 545 { 546 printf ("mpfr_sin_cos and mpfr_cos disagree (bug20091122)\n"); 547 printf ("zref = "); mpfr_dump (zref); 548 printf ("z = "); mpfr_dump (z); 549 exit (1); 550 } 551 552 mpfr_clear (x); 553 mpfr_clear (y); 554 mpfr_clear (z); 555 mpfr_clear (yref); 556 mpfr_clear (zref); 557 } 558 559 static void 560 consistency (void) 561 { 562 mpfr_t x, s1, s2, c1, c2; 563 mpfr_exp_t emin, emax; 564 mpfr_rnd_t rnd; 565 unsigned int flags_sin, flags_cos, flags, flags_before, flags_ref; 566 int inex_sin, is, inex_cos, ic, inex, inex_ref; 567 int i; 568 569 emin = mpfr_get_emin (); 570 emax = mpfr_get_emax (); 571 572 for (i = 0; i <= 10000; i++) 573 { 574 mpfr_init2 (x, MPFR_PREC_MIN + (randlimb () % 8)); 575 mpfr_inits2 (MPFR_PREC_MIN + (randlimb () % 8), s1, s2, c1, c2, 576 (mpfr_ptr) 0); 577 if (i < 8 * MPFR_RND_MAX) 578 { 579 int j = i / MPFR_RND_MAX; 580 if (j & 1) 581 mpfr_set_emin (MPFR_EMIN_MIN); 582 mpfr_set_si (x, (j & 2) ? 1 : -1, MPFR_RNDN); 583 mpfr_set_exp (x, mpfr_get_emin ()); 584 rnd = (mpfr_rnd_t) (i % MPFR_RND_MAX); 585 flags_before = 0; 586 if (j & 4) 587 mpfr_set_emax (-17); 588 } 589 else 590 { 591 tests_default_random (x, 256, -5, 50, 0); 592 rnd = RND_RAND (); 593 flags_before = (randlimb () & 1) ? 594 (unsigned int) (MPFR_FLAGS_ALL ^ MPFR_FLAGS_ERANGE) : 595 (unsigned int) 0; 596 } 597 __gmpfr_flags = flags_before; 598 inex_sin = mpfr_sin (s1, x, rnd); 599 is = inex_sin < 0 ? 2 : inex_sin > 0 ? 1 : 0; 600 flags_sin = __gmpfr_flags; 601 __gmpfr_flags = flags_before; 602 inex_cos = mpfr_cos (c1, x, rnd); 603 ic = inex_cos < 0 ? 2 : inex_cos > 0 ? 1 : 0; 604 flags_cos = __gmpfr_flags; 605 __gmpfr_flags = flags_before; 606 inex = mpfr_sin_cos (s2, c2, x, rnd); 607 flags = __gmpfr_flags; 608 inex_ref = is + 4 * ic; 609 flags_ref = flags_sin | flags_cos; 610 if (!(mpfr_equal_p (s1, s2) && mpfr_equal_p (c1, c2)) || 611 inex != inex_ref || flags != flags_ref) 612 { 613 printf ("mpfr_sin_cos and mpfr_sin/mpfr_cos disagree on %s," 614 " i = %d\nx = ", mpfr_print_rnd_mode (rnd), i); 615 mpfr_dump (x); 616 printf ("s1 = "); 617 mpfr_dump (s1); 618 printf ("s2 = "); 619 mpfr_dump (s2); 620 printf ("c1 = "); 621 mpfr_dump (c1); 622 printf ("c2 = "); 623 mpfr_dump (c2); 624 printf ("inex_sin = %d (s = %d), inex_cos = %d (c = %d), " 625 "inex = %d (expected %d)\n", 626 inex_sin, is, inex_cos, ic, inex, inex_ref); 627 printf ("flags_sin = 0x%x, flags_cos = 0x%x, " 628 "flags = 0x%x (expected 0x%x)\n", 629 flags_sin, flags_cos, flags, flags_ref); 630 exit (1); 631 } 632 mpfr_clears (x, s1, s2, c1, c2, (mpfr_ptr) 0); 633 mpfr_set_emin (emin); 634 mpfr_set_emax (emax); 635 } 636 } 637 638 static void 639 coverage_01032011 (void) 640 { 641 mpfr_t val, cval, sval, svalf; 642 int status_f, status; 643 644 mpfr_init2 (val, MPFR_PREC_MIN); 645 mpfr_init2 (cval, MPFR_PREC_MIN); 646 mpfr_init2 (sval, MPFR_PREC_MIN); 647 mpfr_init2 (svalf, MPFR_PREC_MIN); 648 649 mpfr_set_str1 (val, "-0.7"); 650 651 status_f = mpfr_sincos_fast (svalf, NULL, val, MPFR_RNDN); 652 status = mpfr_sin_cos (sval, cval, val, MPFR_RNDN); 653 if (! mpfr_equal_p (svalf, sval) || SIGN (status_f) != SIGN (status)) 654 { 655 printf ("mpfr_sincos_fast differ from mpfr_sin_cos result:\n" 656 " sin fast is "); 657 mpfr_dump (svalf); 658 printf (" sin is "); 659 mpfr_dump (sval); 660 printf ("status_f = %d, status = %d\n", status_f, status); 661 exit (1); 662 } 663 664 mpfr_clears(val, cval, sval, svalf, (mpfr_ptr) 0); 665 } 666 667 /* tsin_cos prec [N] performs N tests with prec bits */ 668 int 669 main (int argc, char *argv[]) 670 { 671 tests_start_mpfr (); 672 673 if (argc > 1) 674 { 675 if (argc != 3 && argc != 4) 676 { 677 fprintf (stderr, "Usage: tsin_cos x prec [n]\n"); 678 exit (1); 679 } 680 large_test (argv[1], atoi (argv[2]), (argc > 3) ? atoi (argv[3]) : 1); 681 goto end; 682 } 683 684 bug20091013 (); 685 bug20091008 (); 686 bug20091007 (); 687 bug20091122 (); 688 consistency (); 689 690 test_mpfr_sincos_fast (); 691 692 check_nans (); 693 694 /* worst case from PhD thesis of Vincent Lefe`vre: x=8980155785351021/2^54 */ 695 check53 ("4.984987858808754279e-1", "4.781075595393330379e-1", 696 "8.783012931285841817e-1", MPFR_RNDN); 697 check53 ("4.984987858808754279e-1", "4.781075595393329824e-1", 698 "8.783012931285840707e-1", MPFR_RNDD); 699 check53 ("4.984987858808754279e-1", "4.781075595393329824e-1", 700 "8.783012931285840707e-1", MPFR_RNDZ); 701 check53 ("4.984987858808754279e-1", "4.781075595393330379e-1", 702 "8.783012931285841817e-1", MPFR_RNDU); 703 check53 ("1.00031274099908640274", "8.416399183372403892e-1", 704 "0.540039116973283217504", MPFR_RNDN); 705 check53 ("1.00229256850978698523", "8.427074524447979442e-1", 706 "0.538371757797526551137", MPFR_RNDZ); 707 check53 ("1.00288304857059840103", "8.430252033025980029e-1", 708 "0.537874062022526966409", MPFR_RNDZ); 709 check53 ("1.00591265847407274059", "8.446508805292128885e-1", 710 "0.53531755997839769456", MPFR_RNDN); 711 712 /* check one argument only */ 713 check53sin ("1.00591265847407274059", "8.446508805292128885e-1", MPFR_RNDN); 714 check53cos ("1.00591265847407274059", "0.53531755997839769456", MPFR_RNDN); 715 716 overflowed_sin_cos0 (); 717 tiny (); 718 test20071214 (); 719 720 coverage_01032011 (); 721 722 end: 723 tests_end_mpfr (); 724 return 0; 725 } 726