xref: /netbsd-src/external/lgpl3/mpfr/dist/tests/tui_sub.c (revision 87d689fb734c654d2486f87f7be32f1b53ecdbec)
1 /* Test file for mpfr_ui_sub.
2 
3 Copyright 2000-2016 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 http://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 <stdio.h>
24 #include <stdlib.h>
25 #include <float.h>
26 
27 #include "mpfr-test.h"
28 
29 static void
30 special (void)
31 {
32   mpfr_t x, y, res;
33   int inexact;
34 
35   mpfr_init (x);
36   mpfr_init (y);
37   mpfr_init (res);
38 
39   mpfr_set_prec (x, 24);
40   mpfr_set_prec (y, 24);
41   mpfr_set_str_binary (y, "0.111100110001011010111");
42   inexact = mpfr_ui_sub (x, 1, y, MPFR_RNDN);
43   if (inexact)
44     {
45       printf ("Wrong inexact flag: got %d, expected 0\n", inexact);
46       exit (1);
47     }
48 
49   mpfr_set_prec (x, 24);
50   mpfr_set_prec (y, 24);
51   mpfr_set_str_binary (y, "0.111100110001011010111");
52   if ((inexact = mpfr_ui_sub (x, 38181761, y, MPFR_RNDN)) >= 0)
53     {
54       printf ("Wrong inexact flag: got %d, expected -1\n", inexact);
55       exit (1);
56     }
57 
58   mpfr_set_prec (x, 63);
59   mpfr_set_prec (y, 63);
60   mpfr_set_str_binary (y, "0.111110010010100100110101101010001001100101110001000101110111111E-1");
61   if ((inexact = mpfr_ui_sub (x, 1541116494, y, MPFR_RNDN)) <= 0)
62     {
63       printf ("Wrong inexact flag: got %d, expected +1\n", inexact);
64       exit (1);
65     }
66 
67   mpfr_set_prec (x, 32);
68   mpfr_set_prec (y, 32);
69   mpfr_set_str_binary (y, "0.11011000110111010001011100011100E-1");
70   if ((inexact = mpfr_ui_sub (x, 2000375416, y, MPFR_RNDN)) >= 0)
71     {
72       printf ("Wrong inexact flag: got %d, expected -1\n", inexact);
73       exit (1);
74     }
75 
76   mpfr_set_prec (x, 24);
77   mpfr_set_prec (y, 24);
78   mpfr_set_str_binary (y, "0.110011011001010011110111E-2");
79   if ((inexact = mpfr_ui_sub (x, 927694848, y, MPFR_RNDN)) <= 0)
80     {
81       printf ("Wrong inexact flag: got %d, expected +1\n", inexact);
82       exit (1);
83     }
84 
85   /* bug found by Mathieu Dutour, 12 Apr 2001 */
86   mpfr_set_prec (x, 5);
87   mpfr_set_prec (y, 5);
88   mpfr_set_prec (res, 5);
89   mpfr_set_str_binary (x, "1e-12");
90 
91   mpfr_ui_sub (y, 1, x, MPFR_RNDD);
92   mpfr_set_str_binary (res, "0.11111");
93   if (mpfr_cmp (y, res))
94     {
95       printf ("Error in mpfr_ui_sub (y, 1, x, MPFR_RNDD) for x=2^(-12)\nexpected 1.1111e-1, got ");
96       mpfr_out_str (stdout, 2, 0, y, MPFR_RNDN);
97       printf ("\n");
98       exit (1);
99     }
100 
101   mpfr_ui_sub (y, 1, x, MPFR_RNDU);
102   mpfr_set_str_binary (res, "1.0");
103   if (mpfr_cmp (y, res))
104     {
105       printf ("Error in mpfr_ui_sub (y, 1, x, MPFR_RNDU) for x=2^(-12)\n"
106               "expected 1.0, got ");
107       mpfr_out_str (stdout, 2, 0, y, MPFR_RNDN);
108       printf ("\n");
109       exit (1);
110     }
111 
112   mpfr_ui_sub (y, 1, x, MPFR_RNDN);
113   mpfr_set_str_binary (res, "1.0");
114   if (mpfr_cmp (y, res))
115     {
116       printf ("Error in mpfr_ui_sub (y, 1, x, MPFR_RNDN) for x=2^(-12)\n"
117               "expected 1.0, got ");
118       mpfr_out_str (stdout, 2, 0, y, MPFR_RNDN);
119       printf ("\n");
120       exit (1);
121     }
122 
123   mpfr_set_prec (x, 10);
124   mpfr_set_prec (y, 10);
125   mpfr_urandomb (x, RANDS);
126   mpfr_ui_sub (y, 0, x, MPFR_RNDN);
127   if (MPFR_IS_ZERO(x))
128     MPFR_ASSERTN(MPFR_IS_ZERO(y));
129   else
130     MPFR_ASSERTN(mpfr_cmpabs (x, y) == 0 && mpfr_sgn (x) != mpfr_sgn (y));
131 
132   mpfr_set_prec (x, 73);
133   mpfr_set_str_binary (x, "0.1101111010101011011011100011010000000101110001011111001011011000101111101E-99");
134   mpfr_ui_sub (x, 1, x, MPFR_RNDZ);
135   mpfr_nextabove (x);
136   MPFR_ASSERTN(mpfr_cmp_ui (x, 1) == 0);
137 
138   mpfr_clear (x);
139   mpfr_clear (y);
140   mpfr_clear (res);
141 }
142 
143 /* checks that (y-x) gives the right results with 53 bits of precision */
144 static void
145 check (unsigned long y, const char *xs, mpfr_rnd_t rnd_mode, const char *zs)
146 {
147   mpfr_t xx, zz;
148 
149   mpfr_inits2 (53, xx, zz, (mpfr_ptr) 0);
150   mpfr_set_str1 (xx, xs);
151   mpfr_ui_sub (zz, y, xx, rnd_mode);
152   if (mpfr_cmp_str1 (zz, zs) )
153     {
154       printf ("expected difference is %s, got\n",zs);
155       mpfr_out_str(stdout, 10, 0, zz, MPFR_RNDN);
156       printf ("mpfr_ui_sub failed for y=%lu x=%s with rnd_mode=%s\n",
157               y, xs, mpfr_print_rnd_mode (rnd_mode));
158       exit (1);
159     }
160   mpfr_clears (xx, zz, (mpfr_ptr) 0);
161 }
162 
163 /* if u = o(x-y), v = o(u-x), w = o(v+y), then x-y = u-w */
164 static void
165 check_two_sum (mpfr_prec_t p)
166 {
167   unsigned int x;
168   mpfr_t y, u, v, w;
169   mpfr_rnd_t rnd;
170   int inexact, cmp;
171 
172   mpfr_inits2 (p, y, u, v, w, (mpfr_ptr) 0);
173   do
174     {
175       x = randlimb ();
176     }
177   while (x < 1);
178   mpfr_urandomb (y, RANDS);
179   rnd = MPFR_RNDN;
180   inexact = mpfr_ui_sub (u, x, y, rnd);
181   mpfr_sub_ui (v, u, x, rnd);
182   mpfr_add (w, v, y, rnd);
183   cmp = mpfr_cmp_ui (w, 0);
184   /* as u = (x-y) + w, we should have inexact and w of same sign */
185   if (! SAME_SIGN (inexact, cmp))
186     {
187       printf ("Wrong inexact flag for prec=%u, rnd=%s\n",
188               (unsigned int) p, mpfr_print_rnd_mode (rnd));
189       printf ("x = %u\n", x);
190       printf ("y = "); mpfr_dump (y);
191       printf ("u = "); mpfr_dump (u);
192       printf ("v = "); mpfr_dump (v);
193       printf ("w = "); mpfr_dump (w);
194       printf ("inexact = %d\n", inexact);
195       exit (1);
196     }
197   mpfr_clears (y, u, v, w, (mpfr_ptr) 0);
198 }
199 
200 static void
201 check_nans (void)
202 {
203   mpfr_t  x, y;
204 
205   mpfr_init2 (x, 123L);
206   mpfr_init2 (y, 123L);
207 
208   /* 1 - nan == nan */
209   mpfr_set_nan (x);
210   mpfr_ui_sub (y, 1L, x, MPFR_RNDN);
211   MPFR_ASSERTN (mpfr_nan_p (y));
212 
213   /* 1 - +inf == -inf */
214   mpfr_set_inf (x, 1);
215   mpfr_ui_sub (y, 1L, x, MPFR_RNDN);
216   MPFR_ASSERTN (mpfr_inf_p (y));
217   MPFR_ASSERTN (mpfr_sgn (y) < 0);
218 
219   /* 1 - -inf == +inf */
220   mpfr_set_inf (x, -1);
221   mpfr_ui_sub (y, 1L, x, MPFR_RNDN);
222   MPFR_ASSERTN (mpfr_inf_p (y));
223   MPFR_ASSERTN (mpfr_sgn (y) > 0);
224 
225   mpfr_clear (x);
226   mpfr_clear (y);
227 }
228 
229 /* Check mpfr_ui_sub with u = 0 (unsigned). */
230 static void check_neg (void)
231 {
232   mpfr_t x, yneg, ysub;
233   int i, s;
234   int r;
235 
236   mpfr_init2 (x, 64);
237   mpfr_init2 (yneg, 32);
238   mpfr_init2 (ysub, 32);
239 
240   for (i = 0; i <= 25; i++)
241     {
242       mpfr_sqrt_ui (x, i, MPFR_RNDN);
243       for (s = 0; s <= 1; s++)
244         {
245           RND_LOOP (r)
246             {
247               int tneg, tsub;
248 
249               tneg = mpfr_neg (yneg, x, (mpfr_rnd_t) r);
250               tsub = mpfr_ui_sub (ysub, 0, x, (mpfr_rnd_t) r);
251               MPFR_ASSERTN (mpfr_equal_p (yneg, ysub));
252               MPFR_ASSERTN (!(MPFR_IS_POS (yneg) ^ MPFR_IS_POS (ysub)));
253               MPFR_ASSERTN (tneg == tsub);
254             }
255           mpfr_neg (x, x, MPFR_RNDN);
256         }
257     }
258 
259   mpfr_clear (x);
260   mpfr_clear (yneg);
261   mpfr_clear (ysub);
262 }
263 
264 static void
265 check_overflow (void)
266 {
267   mpfr_exp_t emin, emax;
268   mpfr_t x, y1, y2;
269   int inex1, inex2, rnd_mode;
270   unsigned int flags1, flags2;
271 
272   emin = mpfr_get_emin ();
273   emax = mpfr_get_emax ();
274   set_emin (MPFR_EMIN_MIN);
275   set_emax (MPFR_EMAX_MAX);
276 
277   mpfr_inits2 (32, x, y1, y2, (mpfr_ptr) 0);
278   mpfr_setmax (x, MPFR_EMAX_MAX);
279   mpfr_neg (x, x, MPFR_RNDN);
280   RND_LOOP (rnd_mode)
281   {
282     if (rnd_mode == MPFR_RNDU || rnd_mode == MPFR_RNDA)
283       {
284         inex1 = mpfr_overflow (y1, (mpfr_rnd_t) rnd_mode, 1);
285         flags1 = MPFR_FLAGS_OVERFLOW | MPFR_FLAGS_INEXACT;
286       }
287     else
288       {
289         mpfr_neg (y1, x, MPFR_RNDN);
290         inex1 = -1;
291         flags1 = MPFR_FLAGS_INEXACT;
292       }
293     mpfr_clear_flags ();
294     inex2 = mpfr_ui_sub (y2, 1, x, (mpfr_rnd_t) rnd_mode);
295     flags2 = __gmpfr_flags;
296     if (!(mpfr_equal_p (y1, y2) &&
297           SAME_SIGN (inex1, inex2) &&
298           flags1 == flags2))
299       {
300         printf ("Error in check_overflow for %s\n",
301                 mpfr_print_rnd_mode ((mpfr_rnd_t) rnd_mode));
302         printf ("Expected ");
303         mpfr_dump (y1);
304         printf ("  with inex = %d, flags =", inex1);
305         flags_out (flags1);
306         printf ("Got      ");
307         mpfr_dump (y2);
308         printf ("  with inex = %d, flags =", inex2);
309         flags_out (flags2);
310         exit (1);
311       }
312   }
313   mpfr_clears (x, y1, y2, (mpfr_ptr) 0);
314 
315   set_emin (emin);
316   set_emax (emax);
317 }
318 
319 #define TEST_FUNCTION mpfr_ui_sub
320 #define ULONG_ARG1
321 #include "tgeneric.c"
322 
323 int
324 main (int argc, char *argv[])
325 {
326   mpfr_prec_t p;
327   unsigned k;
328 
329   tests_start_mpfr ();
330 
331   check_nans ();
332 
333   special ();
334   for (p=2; p<100; p++)
335     for (k=0; k<100; k++)
336       check_two_sum (p);
337 
338   check(1196426492, "1.4218093058435347e-3", MPFR_RNDN,
339         "1.1964264919985781e9");
340   check(1092583421, "-1.0880649218158844e9", MPFR_RNDN,
341         "2.1806483428158845901e9");
342   check(948002822, "1.22191250737771397120e+20", MPFR_RNDN,
343         "-1.2219125073682338611e20");
344   check(832100416, "4.68311314939691330000e-215", MPFR_RNDD,
345         "8.3210041599999988079e8");
346   check(1976245324, "1.25296395864546893357e+232", MPFR_RNDZ,
347         "-1.2529639586454686577e232");
348   check(2128997392, "-1.08496826129284207724e+187", MPFR_RNDU,
349         "1.0849682612928422704e187");
350   check(293607738, "-1.9967571564050541e-5", MPFR_RNDU,
351         "2.9360773800002003e8");
352   check(354270183, "2.9469161763489528e3", MPFR_RNDN,
353         "3.5426723608382362e8");
354   check_overflow ();
355 
356   check_neg ();
357 
358   test_generic (2, 1000, 100);
359 
360   tests_end_mpfr ();
361   return 0;
362 }
363