1/* Test file for mpfr_sin_cos. 2 3Copyright 2000, 2001, 2002, 2003, 2004, 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 26#include "mpfr-test.h" 27 28static void 29large_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 47static void 48check53 (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 77static void 78check53sin (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 97static void 98check53cos (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 117static void 118check_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 168static void 169overflowed_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 257static void 258tiny (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 */ 279static void 280test20071214 (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 */ 311static void 312test_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 371static void 372bug20091007 (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. */ 419static void 420bug20091008 (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 467static void 468bug20091013 (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 http://websympa.loria.fr/wwsympa/arc/mpfr/2009-11/msg00044.html */ 518static void 519bug20091122 (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 559static void 560consistency (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); 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/* tsin_cos prec [N] performs N tests with prec bits */ 639int 640main (int argc, char *argv[]) 641{ 642 tests_start_mpfr (); 643 644 if (argc > 1) 645 { 646 if (argc != 3 && argc != 4) 647 { 648 fprintf (stderr, "Usage: tsin_cos x prec [n]\n"); 649 exit (1); 650 } 651 large_test (argv[1], atoi (argv[2]), (argc > 3) ? atoi (argv[3]) : 1); 652 goto end; 653 } 654 655 bug20091013 (); 656 bug20091008 (); 657 bug20091007 (); 658 bug20091122 (); 659 consistency (); 660 661 test_mpfr_sincos_fast (); 662 663 check_nans (); 664 665 /* worst case from PhD thesis of Vincent Lefe`vre: x=8980155785351021/2^54 */ 666 check53 ("4.984987858808754279e-1", "4.781075595393330379e-1", 667 "8.783012931285841817e-1", MPFR_RNDN); 668 check53 ("4.984987858808754279e-1", "4.781075595393329824e-1", 669 "8.783012931285840707e-1", MPFR_RNDD); 670 check53 ("4.984987858808754279e-1", "4.781075595393329824e-1", 671 "8.783012931285840707e-1", MPFR_RNDZ); 672 check53 ("4.984987858808754279e-1", "4.781075595393330379e-1", 673 "8.783012931285841817e-1", MPFR_RNDU); 674 check53 ("1.00031274099908640274", "8.416399183372403892e-1", 675 "0.540039116973283217504", MPFR_RNDN); 676 check53 ("1.00229256850978698523", "8.427074524447979442e-1", 677 "0.538371757797526551137", MPFR_RNDZ); 678 check53 ("1.00288304857059840103", "8.430252033025980029e-1", 679 "0.537874062022526966409", MPFR_RNDZ); 680 check53 ("1.00591265847407274059", "8.446508805292128885e-1", 681 "0.53531755997839769456", MPFR_RNDN); 682 683 /* check one argument only */ 684 check53sin ("1.00591265847407274059", "8.446508805292128885e-1", MPFR_RNDN); 685 check53cos ("1.00591265847407274059", "0.53531755997839769456", MPFR_RNDN); 686 687 overflowed_sin_cos0 (); 688 tiny (); 689 test20071214 (); 690 691 end: 692 tests_end_mpfr (); 693 return 0; 694} 695