xref: /netbsd-src/external/lgpl3/mpfr/dist/tests/tcos.c (revision ba125506a622fe649968631a56eba5d42ff57863)
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
test_cos(mpfr_ptr a,mpfr_srcptr b,mpfr_rnd_t rnd_mode)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
check53(const char * xs,const char * cos_xs,mpfr_rnd_t rnd_mode)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
check_nans(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
special_overflow(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
overflowed_cos0(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
bug20091030(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
main(int argc,char * argv[])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 }
383