1/* Test file for mpfr_ui_sub.
2
3Copyright 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011 Free Software Foundation, Inc.
4Contributed by the Arenaire and Cacao projects, INRIA.
5
6This file is part of the GNU MPFR Library.
7
8The GNU MPFR Library is free software; you can redistribute it and/or modify
9it under the terms of the GNU Lesser General Public License as published by
10the Free Software Foundation; either version 3 of the License, or (at your
11option) any later version.
12
13The GNU MPFR Library is distributed in the hope that it will be useful, but
14WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
15or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
16License for more details.
17
18You should have received a copy of the GNU Lesser General Public License
19along with the GNU MPFR Library; see the file COPYING.LESSER.  If not, see
20http://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc.,
2151 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
29static void
30special (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 */
144static void
145check (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 */
164static void
165check_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;
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  /* as u = (x-y) + w, we should have inexact and w of same sign */
184  if (((inexact == 0) && mpfr_cmp_ui (w, 0)) ||
185      ((inexact > 0) && (mpfr_cmp_ui (w, 0) <= 0)) ||
186      ((inexact < 0) && (mpfr_cmp_ui (w, 0) >= 0)))
187    {
188      printf ("Wrong inexact flag for prec=%u, rnd=%s\n",
189              (unsigned int) p, mpfr_print_rnd_mode (rnd));
190      printf ("x=%u\n", x);
191      printf ("y="); mpfr_print_binary(y); puts ("");
192      printf ("u="); mpfr_print_binary(u); puts ("");
193      printf ("v="); mpfr_print_binary(v); puts ("");
194      printf ("w="); mpfr_print_binary(w); puts ("");
195      printf ("inexact = %d\n", inexact);
196      exit (1);
197    }
198  mpfr_clears (y, u, v, w, (mpfr_ptr) 0);
199}
200
201static void
202check_nans (void)
203{
204  mpfr_t  x, y;
205
206  mpfr_init2 (x, 123L);
207  mpfr_init2 (y, 123L);
208
209  /* 1 - nan == nan */
210  mpfr_set_nan (x);
211  mpfr_ui_sub (y, 1L, x, MPFR_RNDN);
212  MPFR_ASSERTN (mpfr_nan_p (y));
213
214  /* 1 - +inf == -inf */
215  mpfr_set_inf (x, 1);
216  mpfr_ui_sub (y, 1L, x, MPFR_RNDN);
217  MPFR_ASSERTN (mpfr_inf_p (y));
218  MPFR_ASSERTN (mpfr_sgn (y) < 0);
219
220  /* 1 - -inf == +inf */
221  mpfr_set_inf (x, -1);
222  mpfr_ui_sub (y, 1L, x, MPFR_RNDN);
223  MPFR_ASSERTN (mpfr_inf_p (y));
224  MPFR_ASSERTN (mpfr_sgn (y) > 0);
225
226  mpfr_clear (x);
227  mpfr_clear (y);
228}
229
230/* Check mpfr_ui_sub with u = 0 (unsigned). */
231static void check_neg (void)
232{
233  mpfr_t x, yneg, ysub;
234  int i, s;
235  int r;
236
237  mpfr_init2 (x, 64);
238  mpfr_init2 (yneg, 32);
239  mpfr_init2 (ysub, 32);
240
241  for (i = 0; i <= 25; i++)
242    {
243      mpfr_sqrt_ui (x, i, MPFR_RNDN);
244      for (s = 0; s <= 1; s++)
245        {
246          RND_LOOP (r)
247            {
248              int tneg, tsub;
249
250              tneg = mpfr_neg (yneg, x, (mpfr_rnd_t) r);
251              tsub = mpfr_ui_sub (ysub, 0, x, (mpfr_rnd_t) r);
252              MPFR_ASSERTN (mpfr_equal_p (yneg, ysub));
253              MPFR_ASSERTN (!(MPFR_IS_POS (yneg) ^ MPFR_IS_POS (ysub)));
254              MPFR_ASSERTN (tneg == tsub);
255            }
256          mpfr_neg (x, x, MPFR_RNDN);
257        }
258    }
259
260  mpfr_clear (x);
261  mpfr_clear (yneg);
262  mpfr_clear (ysub);
263}
264
265int
266main (int argc, char *argv[])
267{
268  mpfr_prec_t p;
269  unsigned k;
270
271  tests_start_mpfr ();
272
273  check_nans ();
274
275  special ();
276  for (p=2; p<100; p++)
277    for (k=0; k<100; k++)
278      check_two_sum (p);
279
280  check(1196426492, "1.4218093058435347e-3", MPFR_RNDN,
281        "1.1964264919985781e9");
282  check(1092583421, "-1.0880649218158844e9", MPFR_RNDN,
283        "2.1806483428158845901e9");
284  check(948002822, "1.22191250737771397120e+20", MPFR_RNDN,
285        "-1.2219125073682338611e20");
286  check(832100416, "4.68311314939691330000e-215", MPFR_RNDD,
287        "8.3210041599999988079e8");
288  check(1976245324, "1.25296395864546893357e+232", MPFR_RNDZ,
289        "-1.2529639586454686577e232");
290  check(2128997392, "-1.08496826129284207724e+187", MPFR_RNDU,
291        "1.0849682612928422704e187");
292  check(293607738, "-1.9967571564050541e-5", MPFR_RNDU,
293        "2.9360773800002003e8");
294  check(354270183, "2.9469161763489528e3", MPFR_RNDN,
295        "3.5426723608382362e8");
296
297  check_neg ();
298
299  tests_end_mpfr ();
300  return 0;
301}
302