1 /* Test file for mpfr_cos.
2
3 Copyright 2001-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 #ifdef CHECK_EXTERNAL
26 static int
27 test_cos (mpfr_ptr a, mpfr_srcptr b, mpfr_rnd_t rnd_mode)
28 {
29 int res;
30 int ok = rnd_mode == MPFR_RNDN && mpfr_number_p (b) && mpfr_get_prec (a)>=53;
31 if (ok)
32 {
33 mpfr_print_raw (b);
34 }
35 res = mpfr_cos (a, b, rnd_mode);
36 if (ok)
37 {
38 printf (" ");
39 mpfr_print_raw (a);
40 printf ("\n");
41 }
42 return res;
43 }
44 #else
45 #define test_cos mpfr_cos
46 #endif
47
48 static void
49 check53 (const char *xs, const char *cos_xs, mpfr_rnd_t rnd_mode)
50 {
51 mpfr_t xx, c;
52
53 mpfr_inits2 (53, xx, c, (mpfr_ptr) 0);
54 mpfr_set_str1 (xx, xs); /* should be exact */
55 test_cos (c, xx, rnd_mode);
56 if (mpfr_cmp_str1 (c, cos_xs))
57 {
58 printf ("mpfr_cos failed for x=%s, rnd=%s\n",
59 xs, mpfr_print_rnd_mode (rnd_mode));
60 printf ("mpfr_cos gives cos(x)=");
61 mpfr_out_str(stdout, 10, 0, c, MPFR_RNDN);
62 printf(", expected %s\n", cos_xs);
63 exit (1);
64 }
65 mpfr_clears (xx, c, (mpfr_ptr) 0);
66 }
67
68 #define TEST_FUNCTION test_cos
69 #define REDUCE_EMAX 262143 /* otherwise arg. reduction is too expensive */
70 #include "tgeneric.c"
71
72 static void
73 check_nans (void)
74 {
75 mpfr_t x, y;
76
77 mpfr_init2 (x, 123L);
78 mpfr_init2 (y, 123L);
79
80 mpfr_set_nan (x);
81 test_cos (y, x, MPFR_RNDN);
82 if (! mpfr_nan_p (y))
83 {
84 printf ("Error: cos(NaN) != NaN\n");
85 exit (1);
86 }
87
88 mpfr_set_inf (x, 1);
89 test_cos (y, x, MPFR_RNDN);
90 if (! mpfr_nan_p (y))
91 {
92 printf ("Error: cos(Inf) != NaN\n");
93 exit (1);
94 }
95
96 mpfr_set_inf (x, -1);
97 test_cos (y, x, MPFR_RNDN);
98 if (! mpfr_nan_p (y))
99 {
100 printf ("Error: cos(-Inf) != NaN\n");
101 exit (1);
102 }
103
104 /* cos(+/-0) = 1 */
105 mpfr_set_ui (x, 0, MPFR_RNDN);
106 test_cos (y, x, MPFR_RNDN);
107 if (mpfr_cmp_ui (y, 1))
108 {
109 printf ("Error: cos(+0) != 1\n");
110 exit (1);
111 }
112 mpfr_neg (x, x, MPFR_RNDN);
113 test_cos (y, x, MPFR_RNDN);
114 if (mpfr_cmp_ui (y, 1))
115 {
116 printf ("Error: cos(-0) != 1\n");
117 exit (1);
118 }
119
120 /* Compute ~Pi/2 to check */
121 /* FIXME: Too slow!
122 mpfr_set_prec (x, 20000);
123 mpfr_const_pi (x, MPFR_RNDD); mpfr_div_2ui (x, x, 1, MPFR_RNDN);
124 mpfr_set_prec (y, 24);
125 test_cos (y, x, MPFR_RNDN);
126 if (mpfr_cmp_str (y, "0.111001010110100011000001E-20000", 2, MPFR_RNDN))
127 {
128 printf("Error computing cos(~Pi/2)\n");
129 mpfr_dump (y);
130 exit (1);
131 } */
132
133 mpfr_clear (x);
134 mpfr_clear (y);
135 }
136
137 static void
138 special_overflow (void)
139 {
140 mpfr_t x, y;
141 mpfr_exp_t emin, emax;
142
143 emin = mpfr_get_emin ();
144 emax = mpfr_get_emax ();
145
146 mpfr_init2 (x, 24);
147 mpfr_init2 (y, 73);
148
149 /* Check special case: An overflow in const_pi could occurs! */
150 set_emin (-125);
151 set_emax (128);
152 mpfr_set_str_binary (x, "0.111101010110110011101101E6");
153 test_cos (y, x, MPFR_RNDZ);
154 set_emin (emin);
155 set_emax (emax);
156
157 mpfr_clear (x);
158 mpfr_clear (y);
159 }
160
161 static void
162 overflowed_cos0 (void)
163 {
164 mpfr_t x, y;
165 int emax, i, inex, rnd, err = 0;
166 mpfr_exp_t old_emax;
167
168 old_emax = mpfr_get_emax ();
169
170 mpfr_init2 (x, 8);
171 mpfr_init2 (y, 8);
172
173 for (emax = -1; emax <= 0; emax++)
174 {
175 mpfr_set_ui_2exp (y, 1, emax, MPFR_RNDN);
176 mpfr_nextbelow (y);
177 set_emax (emax); /* 1 is not representable. */
178 /* and if emax < 0, 1 - eps is not representable either. */
179 for (i = -1; i <= 1; i++)
180 RND_LOOP (rnd)
181 {
182 mpfr_set_si_2exp (x, i, -512 * ABS (i), MPFR_RNDN);
183 mpfr_clear_flags ();
184 inex = mpfr_cos (x, x, (mpfr_rnd_t) rnd);
185 if ((i == 0 || emax < 0 || rnd == MPFR_RNDN || rnd == MPFR_RNDU) &&
186 ! mpfr_overflow_p ())
187 {
188 printf ("Error in overflowed_cos0 (i = %d, rnd = %s):\n"
189 " The overflow flag is not set.\n",
190 i, mpfr_print_rnd_mode ((mpfr_rnd_t) rnd));
191 err = 1;
192 }
193 if (rnd == MPFR_RNDZ || rnd == MPFR_RNDD)
194 {
195 if (inex >= 0)
196 {
197 printf ("Error in overflowed_cos0 (i = %d, rnd = %s):\n"
198 " The inexact value must be negative.\n",
199 i, mpfr_print_rnd_mode ((mpfr_rnd_t) rnd));
200 err = 1;
201 }
202 if (! mpfr_equal_p (x, y))
203 {
204 printf ("Error in overflowed_cos0 (i = %d, rnd = %s):\n"
205 " Got ", i,
206 mpfr_print_rnd_mode ((mpfr_rnd_t) rnd));
207 mpfr_dump (x);
208 printf (" instead of 0.11111111E%d.\n", emax);
209 err = 1;
210 }
211 }
212 else if (rnd != MPFR_RNDF)
213 {
214 if (inex <= 0)
215 {
216 printf ("Error in overflowed_cos0 (i = %d, rnd = %s):\n"
217 " The inexact value must be positive.\n",
218 i, mpfr_print_rnd_mode ((mpfr_rnd_t) rnd));
219 err = 1;
220 }
221 if (! (mpfr_inf_p (x) && MPFR_IS_POS (x)))
222 {
223 printf ("Error in overflowed_cos0 (i = %d, rnd = %s):\n"
224 " Got ", i,
225 mpfr_print_rnd_mode ((mpfr_rnd_t) rnd));
226 mpfr_dump (x);
227 printf (" instead of +Inf.\n");
228 err = 1;
229 }
230 }
231 }
232 set_emax (old_emax);
233 }
234
235 if (err)
236 exit (1);
237 mpfr_clear (x);
238 mpfr_clear (y);
239 }
240
241 static void
242 bug20091030 (void)
243 {
244 mpfr_t x, y;
245
246 mpfr_init2 (x, 5);
247 mpfr_init2 (y, 2);
248 mpfr_set_str (x, "-0.11001E3", 2, MPFR_RNDN);
249 mpfr_cos (y, x, MPFR_RNDN);
250 mpfr_clear (x);
251 mpfr_clear (y);
252 }
253
254 int
255 main (int argc, char *argv[])
256 {
257 mpfr_t x, y;
258 int inex;
259
260 tests_start_mpfr ();
261
262 special_overflow ();
263 check_nans ();
264
265 mpfr_init (x);
266 mpfr_init (y);
267
268 mpfr_set_prec (x, 53);
269 mpfr_set_prec (y, 2);
270 mpfr_set_str (x, "9.81333845856942e-1", 10, MPFR_RNDN);
271 test_cos (y, x, MPFR_RNDN);
272
273 mpfr_set_prec (x, 30);
274 mpfr_set_prec (y, 30);
275 mpfr_set_str_binary (x, "1.00001010001101110010100010101e-1");
276 test_cos (y, x, MPFR_RNDU);
277 mpfr_set_str_binary (x, "1.10111100010101011110101010100e-1");
278 if (mpfr_cmp (y, x))
279 {
280 printf ("Error for prec=30, rnd=MPFR_RNDU\n");
281 printf ("expected "); mpfr_dump (x);
282 printf (" got "); mpfr_dump (y);
283 exit (1);
284 }
285
286 mpfr_set_prec (x, 59);
287 mpfr_set_prec (y, 59);
288 mpfr_set_str_binary (x, "1.01101011101111010011111110111111111011011101100111100011e-3");
289 test_cos (y, x, MPFR_RNDU);
290 mpfr_set_str_binary (x, "1.1111011111110010001001001011100111101110100010000010010011e-1");
291 if (mpfr_cmp (y, x))
292 {
293 printf ("Error for prec=59, rnd=MPFR_RNDU\n");
294 printf ("expected "); mpfr_dump (x);
295 printf (" got "); mpfr_dump (y);
296 exit (1);
297 }
298
299 mpfr_set_prec (x, 5);
300 mpfr_set_prec (y, 5);
301 mpfr_set_str_binary (x, "1.1100e-2");
302 test_cos (y, x, MPFR_RNDD);
303 mpfr_set_str_binary (x, "1.1100e-1");
304 if (mpfr_cmp (y, x))
305 {
306 printf ("Error for x=1.1100e-2, rnd=MPFR_RNDD\n");
307 printf ("expected 1.1100e-1, got "); mpfr_dump (y);
308 exit (1);
309 }
310
311 mpfr_set_prec (x, 32);
312 mpfr_set_prec (y, 32);
313
314 mpfr_set_str_binary (x, "0.10001000001001011000100001E-6");
315 mpfr_set_str_binary (y, "0.1111111111111101101111001100001");
316 test_cos (x, x, MPFR_RNDN);
317 if (mpfr_cmp (x, y))
318 {
319 printf ("Error for prec=32 (1)\n");
320 exit (1);
321 }
322
323 mpfr_set_str_binary (x, "-0.1101011110111100111010011001011E-1");
324 mpfr_set_str_binary (y, "0.11101001100110111011011010100011");
325 test_cos (x, x, MPFR_RNDN);
326 if (mpfr_cmp (x, y))
327 {
328 printf ("Error for prec=32 (2)\n");
329 exit (1);
330 }
331
332 /* huge argument reduction */
333 mpfr_set_str_binary (x, "0.10000010000001101011101111001011E40");
334 mpfr_set_str_binary (y, "0.10011000001111010000101011001011E-1");
335 test_cos (x, x, MPFR_RNDN);
336 if (mpfr_cmp (x, y))
337 {
338 printf ("Error for prec=32 (3)\n");
339 exit (1);
340 }
341
342 mpfr_set_prec (x, 3);
343 mpfr_set_prec (y, 3);
344 mpfr_set_str_binary (x, "0.110E60");
345 inex = mpfr_cos (y, x, MPFR_RNDD);
346 MPFR_ASSERTN(inex < 0);
347
348 /* worst case from PhD thesis of Vincent Lefe`vre: x=8980155785351021/2^54 */
349 check53 ("4.984987858808754279e-1", "8.783012931285841817e-1", MPFR_RNDN);
350 check53 ("4.984987858808754279e-1", "8.783012931285840707e-1", MPFR_RNDD);
351 check53 ("4.984987858808754279e-1", "8.783012931285840707e-1", MPFR_RNDZ);
352 check53 ("4.984987858808754279e-1", "8.783012931285841817e-1", MPFR_RNDU);
353 check53 ("1.00031274099908640274", "0.540039116973283217504", MPFR_RNDN);
354 check53 ("1.00229256850978698523", "0.538371757797526551137", MPFR_RNDZ);
355 check53 ("1.00288304857059840103", "0.537874062022526966409", MPFR_RNDZ);
356 check53 ("1.00591265847407274059", "0.53531755997839769456", MPFR_RNDN);
357
358 check53 ("1.00591265847407274059", "0.53531755997839769456", MPFR_RNDN);
359
360 overflowed_cos0 ();
361 test_generic (MPFR_PREC_MIN, 100, 15);
362
363 /* check inexact flag */
364 mpfr_set_prec (x, 3);
365 mpfr_set_prec (y, 13);
366 mpfr_set_str_binary (x, "-0.100E196");
367 inex = mpfr_cos (y, x, MPFR_RNDU);
368 mpfr_set_prec (x, 13);
369 mpfr_set_str_binary (x, "0.1111111100101");
370 MPFR_ASSERTN (inex > 0 && mpfr_equal_p (x, y));
371
372 mpfr_clear (x);
373 mpfr_clear (y);
374
375 bug20091030 ();
376
377 data_check ("data/cos", mpfr_cos, "mpfr_cos");
378 bad_cases (mpfr_cos, mpfr_acos, "mpfr_cos", 256, -40, 0, 4, 128, 800, 50);
379
380 tests_end_mpfr ();
381 return 0;
382 }