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