1 /* Test file for round away.
2
3 Copyright 2000-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 #define DISP(s,t) \
26 do \
27 { \
28 printf (s); \
29 mpfr_out_str (stdout, 2, 0, t, MPFR_RNDN); \
30 } \
31 while (0)
32
33 #define DISP2(s,t) do { DISP(s,t); putchar ('\n'); } while (0)
34
35 #define SPECIAL_MAX 12
36
37 static void
38 set_special (mpfr_ptr x, unsigned int select)
39 {
40 MPFR_ASSERTN (select < SPECIAL_MAX);
41 switch (select)
42 {
43 case 0:
44 MPFR_SET_NAN (x);
45 break;
46 case 1:
47 MPFR_SET_INF (x);
48 MPFR_SET_POS (x);
49 break;
50 case 2:
51 MPFR_SET_INF (x);
52 MPFR_SET_NEG (x);
53 break;
54 case 3:
55 MPFR_SET_ZERO (x);
56 MPFR_SET_POS (x);
57 break;
58 case 4:
59 MPFR_SET_ZERO (x);
60 MPFR_SET_NEG (x);
61 break;
62 case 5:
63 mpfr_set_str_binary (x, "1");
64 break;
65 case 6:
66 mpfr_set_str_binary (x, "-1");
67 break;
68 case 7:
69 mpfr_set_str_binary (x, "1e-1");
70 break;
71 case 8:
72 mpfr_set_str_binary (x, "1e+1");
73 break;
74 case 9:
75 mpfr_const_pi (x, MPFR_RNDN);
76 break;
77 case 10:
78 mpfr_const_pi (x, MPFR_RNDN);
79 MPFR_SET_EXP (x, MPFR_GET_EXP (x)-1);
80 break;
81 default:
82 mpfr_urandomb (x, RANDS);
83 break;
84 }
85 }
86
87 /* same as mpfr_cmp, but returns 0 for both NaN's */
88 static int
89 mpfr_compare (mpfr_srcptr a, mpfr_srcptr b)
90 {
91 return (MPFR_IS_NAN(a)) ? !MPFR_IS_NAN(b) :
92 (MPFR_IS_NAN(b) || mpfr_cmp(a, b));
93 }
94
95 static void
96 test3 (int (*testfunc)(mpfr_ptr, mpfr_srcptr, mpfr_srcptr, mpfr_rnd_t),
97 const char *foo)
98 {
99 mpfr_t ref1, ref2, ref3;
100 mpfr_t res1;
101 mpfr_prec_t p1, p2, p3;
102 int i, inexa, inexd;
103 mpfr_rnd_t r;
104
105 p1 = (randlimb () % 200) + MPFR_PREC_MIN;
106 p2 = (randlimb () % 200) + MPFR_PREC_MIN;
107 p3 = (randlimb () % 200) + MPFR_PREC_MIN;
108
109 mpfr_init2 (ref1, p1);
110 mpfr_init2 (ref2, p2);
111 mpfr_init2 (ref3, p3);
112 mpfr_init2 (res1, p1);
113
114 /* for each variable, consider each of the following 6 possibilities:
115 NaN, +Infinity, -Infinity, +0, -0 or a random number */
116 for (i = 0; i < SPECIAL_MAX * SPECIAL_MAX; i++)
117 {
118 set_special (ref2, i%SPECIAL_MAX);
119 set_special (ref3, i/SPECIAL_MAX);
120
121 inexa = testfunc (res1, ref2, ref3, MPFR_RNDA);
122 r = MPFR_IS_POS (res1) ? MPFR_RNDU : MPFR_RNDD;
123 inexd = testfunc (ref1, ref2, ref3, r);
124
125 if (mpfr_compare (res1, ref1) || inexa != inexd)
126 {
127 printf ("Error with RNDA for %s with ", foo);
128 DISP("x=",ref2); DISP2(", y=",ref3);
129 printf ("inexa=%d inexd=%d\n", inexa, inexd);
130 printf ("expected "); mpfr_dump (ref1);
131 printf ("got "); mpfr_dump (res1);
132 exit (1);
133 }
134 }
135
136 mpfr_clear (ref1);
137 mpfr_clear (ref2);
138 mpfr_clear (ref3);
139 mpfr_clear (res1);
140 }
141
142 static void
143 test4 (int (*testfunc)(mpfr_ptr, mpfr_srcptr, mpfr_srcptr, mpfr_srcptr,
144 mpfr_rnd_t), const char *foo)
145 {
146 mpfr_t ref, op1, op2, op3;
147 mpfr_prec_t pout, p1, p2, p3;
148 mpfr_t res;
149 int i, j, k, inexa, inexd;
150 mpfr_rnd_t r;
151
152 pout = (randlimb () % 200) + MPFR_PREC_MIN;
153 p1 = (randlimb () % 200) + MPFR_PREC_MIN;
154 p2 = (randlimb () % 200) + MPFR_PREC_MIN;
155 p3 = (randlimb () % 200) + MPFR_PREC_MIN;
156
157 mpfr_init2 (ref, pout);
158 mpfr_init2 (res, pout);
159 mpfr_init2 (op1, p1);
160 mpfr_init2 (op2, p2);
161 mpfr_init2 (op3, p3);
162
163 /* for each variable, consider each of the following 6 possibilities:
164 NaN, +Infinity, -Infinity, +0, -0 or a random number */
165
166 for (i = 0; i < SPECIAL_MAX; i++)
167 {
168 set_special (op1, i);
169 for (j = 0; j < SPECIAL_MAX; j++)
170 {
171 set_special (op2, j);
172 for (k = 0; k < SPECIAL_MAX; k++)
173 {
174 set_special (op3, k);
175
176 inexa = testfunc (res, op1, op2, op3, MPFR_RNDA);
177 r = MPFR_IS_POS (res) ? MPFR_RNDU : MPFR_RNDD;
178 inexd = testfunc (ref, op1, op2, op3, r);
179
180 if (mpfr_compare (res, ref) || inexa != inexd)
181 {
182 printf ("Error with RNDA for %s with ", foo);
183 DISP("a=", op1); DISP(", b=", op2); DISP2(", c=", op3);
184 printf ("inexa=%d inexd=%d\n", inexa, inexd);
185 DISP("expected ", ref); DISP2(", got ", res);
186 exit (1);
187 }
188 }
189 }
190 }
191
192 mpfr_clear (ref);
193 mpfr_clear (op1);
194 mpfr_clear (op2);
195 mpfr_clear (op3);
196 mpfr_clear (res);
197 }
198
199 static void
200 test2ui (int (*testfunc)(mpfr_ptr, mpfr_srcptr, unsigned long int, mpfr_rnd_t),
201 const char *foo)
202 {
203 mpfr_t ref1, ref2;
204 unsigned int ref3;
205 mpfr_t res1;
206 mpfr_prec_t p1, p2;
207 int i, inexa, inexd;
208 mpfr_rnd_t r;
209
210 p1 = (randlimb () % 200) + MPFR_PREC_MIN;
211 p2 = (randlimb () % 200) + MPFR_PREC_MIN;
212
213 mpfr_init2 (ref1, p1);
214 mpfr_init2 (ref2, p2);
215 mpfr_init2 (res1, p1);
216
217 /* ref2 can be NaN, +Inf, -Inf, +0, -0 or any number
218 ref3 can be 0 or any number */
219 for (i = 0; i < SPECIAL_MAX * 2; i++)
220 {
221 set_special (ref2, i % SPECIAL_MAX);
222 ref3 = i / SPECIAL_MAX == 0 ? 0 : randlimb ();
223
224 inexa = testfunc (res1, ref2, ref3, MPFR_RNDA);
225 r = MPFR_IS_POS (res1) ? MPFR_RNDU : MPFR_RNDD;
226 inexd = testfunc (ref1, ref2, ref3, r);
227
228 if (mpfr_compare (res1, ref1) || inexa != inexd)
229 {
230 printf ("Error with RNDA for %s for c=%u\n", foo, ref3);
231 DISP2("a=",ref2);
232 printf ("inexa=%d inexd=%d\n", inexa, inexd);
233 printf ("expected "); mpfr_dump (ref1);
234 printf ("got "); mpfr_dump (res1);
235 exit (1);
236 }
237 }
238
239 mpfr_clear (ref1);
240 mpfr_clear (ref2);
241 mpfr_clear (res1);
242 }
243
244 static void
245 testui2 (int (*testfunc)(mpfr_ptr, unsigned long int, mpfr_srcptr, mpfr_rnd_t),
246 const char *foo)
247 {
248 mpfr_t ref1, ref3;
249 unsigned int ref2;
250 mpfr_t res1;
251 mpfr_prec_t p1, p3;
252 int i, inexa, inexd;
253 mpfr_rnd_t r;
254
255 p1 = (randlimb () % 200) + MPFR_PREC_MIN;
256 p3 = (randlimb () % 200) + MPFR_PREC_MIN;
257
258 mpfr_init2 (ref1, p1);
259 mpfr_init2 (ref3, p3);
260 mpfr_init2 (res1, p1);
261
262 for (i = 0; i < SPECIAL_MAX * 2; i++)
263 {
264 set_special (ref3, i % SPECIAL_MAX);
265 ref2 = i / SPECIAL_MAX == 0 ? 0 : randlimb ();
266
267 inexa = testfunc (res1, ref2, ref3, MPFR_RNDA);
268 r = MPFR_IS_POS (res1) ? MPFR_RNDU : MPFR_RNDD;
269 inexd = testfunc (ref1, ref2, ref3, r);
270
271 if (mpfr_compare (res1, ref1) || inexa != inexd)
272 {
273 printf ("Error with RNDA for %s for b=%u\n", foo, ref2);
274 DISP2("a=", ref3);
275 printf ("inexa=%d inexd=%d\n", inexa, inexd);
276 DISP("expected ", ref1); DISP2(", got ", res1);
277 exit (1);
278 }
279 }
280
281 mpfr_clear (ref1);
282 mpfr_clear (ref3);
283 mpfr_clear (res1);
284 }
285
286 /* foo(mpfr_ptr, mpfr_srcptr, mp_rndt) */
287 static void
288 test2 (int (*testfunc)(mpfr_ptr, mpfr_srcptr, mpfr_rnd_t), const char *foo)
289 {
290 mpfr_t ref1, ref2;
291 mpfr_t res1;
292 mpfr_prec_t p1, p2;
293 int i, inexa, inexd;
294 mpfr_rnd_t r;
295
296 p1 = (randlimb () % 200) + MPFR_PREC_MIN;
297 p2 = (randlimb () % 200) + MPFR_PREC_MIN;
298
299 mpfr_init2 (ref1, p1);
300 mpfr_init2 (ref2, p2);
301 mpfr_init2 (res1, p1);
302
303 for (i = 0; i < SPECIAL_MAX; i++)
304 {
305 set_special (ref2, i);
306
307 /* first round to away */
308 inexa = testfunc (res1, ref2, MPFR_RNDA);
309
310 r = MPFR_IS_POS (res1) ? MPFR_RNDU : MPFR_RNDD;
311 inexd = testfunc (ref1, ref2, r);
312 if (mpfr_compare (res1, ref1) || inexa != inexd)
313 {
314 printf ("Error with RNDA for %s with ", foo);
315 DISP2("x=", ref2);
316 printf ("inexa=%d inexd=%d\n", inexa, inexd);
317 DISP("expected ", ref1); DISP2(", got ", res1);
318 exit (1);
319 }
320 }
321
322 mpfr_clear (ref1);
323 mpfr_clear (ref2);
324 mpfr_clear (res1);
325 }
326
327 /* one operand, two results, like mpfr_sin_cos */
328 static void
329 test3a (int (*testfunc)(mpfr_ptr, mpfr_ptr, mpfr_srcptr, mpfr_rnd_t),
330 const char *foo)
331 {
332 mpfr_t ref1, ref2, ref3;
333 mpfr_t res1, res2;
334 mpfr_prec_t p1, p2, p3;
335 int i, inexa, inexd;
336 mpfr_rnd_t r;
337
338 p1 = (randlimb () % 200) + MPFR_PREC_MIN;
339 p2 = (randlimb () % 200) + MPFR_PREC_MIN;
340 p3 = (randlimb () % 200) + MPFR_PREC_MIN;
341
342 mpfr_init2 (ref1, p1);
343 mpfr_init2 (ref2, p2);
344 mpfr_init2 (ref3, p3);
345 mpfr_init2 (res1, p1);
346 mpfr_init2 (res2, p2);
347
348 for (i = 0; i < SPECIAL_MAX; i++)
349 {
350 set_special (ref3, i);
351
352 inexa = testfunc (res1, res2, ref3, MPFR_RNDA);
353
354 /* first check wrt the first operand */
355 r = MPFR_IS_POS (res1) ? MPFR_RNDU : MPFR_RNDD;
356 inexd = testfunc (ref1, ref2, ref3, r);
357 /* the low 2 bits of the inexact flag concern the 1st operand */
358 if (mpfr_compare (res1, ref1) || (inexa & 3) != (inexd & 3))
359 {
360 printf ("Error with RNDA for %s (1st operand)\n", foo);
361 DISP2("a=",ref3);
362 DISP("expected ", ref1); printf ("\n");
363 DISP("got ", res1); printf ("\n");
364 printf ("inexa=%d inexd=%d\n", inexa & 3, inexd & 3);
365 exit (1);
366 }
367
368 /* now check wrt the second operand */
369 r = MPFR_IS_POS (res2) ? MPFR_RNDU : MPFR_RNDD;
370 inexd = testfunc (ref1, ref2, ref3, r);
371 /* bits 2..3 of the inexact flag concern the 2nd operand */
372 if (mpfr_compare (res2, ref2) || (inexa >> 2) != (inexd >> 2))
373 {
374 printf ("Error with RNDA for %s (2nd operand)\n", foo);
375 DISP2("a=",ref3);
376 DISP("expected ", ref2); printf ("\n");
377 DISP("got ", res2); printf ("\n");
378 printf ("inexa=%d inexd=%d\n", inexa >> 2, inexd >> 2);
379 exit (1);
380 }
381
382 }
383
384 mpfr_clear (ref1);
385 mpfr_clear (ref2);
386 mpfr_clear (ref3);
387 mpfr_clear (res1);
388 mpfr_clear (res2);
389 }
390
391 int
392 main (void)
393 {
394 int N = 20;
395
396 tests_start_mpfr ();
397
398 while (N--)
399 {
400 /* no need to test mpfr_round, mpfr_ceil, mpfr_floor, mpfr_trunc since
401 they take no rounding mode */
402
403 test2ui (mpfr_add_ui, "mpfr_add_ui");
404 test2ui (mpfr_div_2exp, "mpfr_div_2exp");
405 test2ui (mpfr_div_2ui, "mpfr_div_2ui");
406 test2ui (mpfr_div_ui, "mpfr_div_ui");
407 test2ui (mpfr_mul_2exp, "mpfr_mul_2exp");
408 test2ui (mpfr_mul_2ui, "mpfr_mul_2ui");
409 test2ui (mpfr_mul_ui, "mpfr_mul_ui");
410 test2ui (mpfr_pow_ui, "mpfr_pow_ui");
411 test2ui (mpfr_sub_ui, "mpfr_sub_ui");
412
413 testui2 (mpfr_ui_div, "mpfr_ui_div");
414 testui2 (mpfr_ui_sub, "mpfr_ui_sub");
415 testui2 (mpfr_ui_pow, "mpfr_ui_pow");
416
417 test2 (mpfr_sqr, "mpfr_sqr");
418 test2 (mpfr_sqrt, "mpfr_sqrt");
419 test2 (mpfr_abs, "mpfr_abs");
420 test2 (mpfr_neg, "mpfr_neg");
421
422 test2 (mpfr_log, "mpfr_log");
423 test2 (mpfr_log2, "mpfr_log2");
424 test2 (mpfr_log10, "mpfr_log10");
425 test2 (mpfr_log1p, "mpfr_log1p");
426
427 test2 (mpfr_exp, "mpfr_exp");
428 test2 (mpfr_exp2, "mpfr_exp2");
429 test2 (mpfr_exp10, "mpfr_exp10");
430 test2 (mpfr_expm1, "mpfr_expm1");
431 test2 (mpfr_eint, "mpfr_eint");
432
433 test2 (mpfr_sinh, "mpfr_sinh");
434 test2 (mpfr_cosh, "mpfr_cosh");
435 test2 (mpfr_tanh, "mpfr_tanh");
436 test2 (mpfr_asinh, "mpfr_asinh");
437 test2 (mpfr_acosh, "mpfr_acosh");
438 test2 (mpfr_atanh, "mpfr_atanh");
439 test2 (mpfr_sech, "mpfr_sech");
440 test2 (mpfr_csch, "mpfr_csch");
441 test2 (mpfr_coth, "mpfr_coth");
442
443 test2 (mpfr_asin, "mpfr_asin");
444 test2 (mpfr_acos, "mpfr_acos");
445 test2 (mpfr_atan, "mpfr_atan");
446 test2 (mpfr_cos, "mpfr_cos");
447 test2 (mpfr_sin, "mpfr_sin");
448 test2 (mpfr_tan, "mpfr_tan");
449 test2 (mpfr_sec, "mpfr_sec");
450 test2 (mpfr_csc, "mpfr_csc");
451 test2 (mpfr_cot, "mpfr_cot");
452
453 test2 (mpfr_erf, "mpfr_erf");
454 test2 (mpfr_erfc, "mpfr_erfc");
455 test2 (mpfr_j0, "mpfr_j0");
456 test2 (mpfr_j1, "mpfr_j1");
457 test2 (mpfr_y0, "mpfr_y0");
458 test2 (mpfr_y1, "mpfr_y1");
459 test2 (mpfr_zeta, "mpfr_zeta");
460 test2 (mpfr_gamma, "mpfr_gamma");
461 test2 (mpfr_lngamma, "mpfr_lngamma");
462
463 test2 (mpfr_rint, "mpfr_rint");
464 test2 (mpfr_rint_ceil, "mpfr_rint_ceil");
465 test2 (mpfr_rint_floor, "mpfr_rint_floor");
466 test2 (mpfr_rint_round, "mpfr_rint_round");
467 test2 (mpfr_rint_trunc, "mpfr_rint_trunc");
468 test2 (mpfr_frac, "mpfr_frac");
469
470 test3 (mpfr_add, "mpfr_add");
471 test3 (mpfr_sub, "mpfr_sub");
472 test3 (mpfr_mul, "mpfr_mul");
473 test3 (mpfr_div, "mpfr_div");
474
475 test3 (mpfr_agm, "mpfr_agm");
476 test3 (mpfr_min, "mpfr_min");
477 test3 (mpfr_max, "mpfr_max");
478
479 /* we don't test reldiff since it does not guarantee correct rounding,
480 thus we can get different results with RNDA and RNDU or RNDD. */
481 test3 (mpfr_dim, "mpfr_dim");
482
483 test3 (mpfr_remainder, "mpfr_remainder");
484 test3 (mpfr_pow, "mpfr_pow");
485 test3 (mpfr_atan2, "mpfr_atan2");
486 test3 (mpfr_hypot, "mpfr_hypot");
487
488 test3a (mpfr_sin_cos, "mpfr_sin_cos");
489
490 test4 (mpfr_fma, "mpfr_fma");
491 test4 (mpfr_fms, "mpfr_fms");
492
493 test2 (mpfr_li2, "mpfr_li2");
494 test2 (mpfr_rec_sqrt, "mpfr_rec_sqrt");
495 test3 (mpfr_fmod, "mpfr_fmod");
496 test3a (mpfr_modf, "mpfr_modf");
497 test3a (mpfr_sinh_cosh, "mpfr_sinh_cosh");
498
499 #if MPFR_VERSION >= MPFR_VERSION_NUM(3,0,0)
500 test2 (mpfr_ai, "mpfr_ai");
501 test2 (mpfr_digamma, "mpfr_digamma");
502 #endif
503 }
504
505 tests_end_mpfr ();
506 return 0;
507 }