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