1/* Test file for mpfr_add and mpfr_sub.
2
3Copyright 1999, 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#define NUM 30000
24
25#include <stdio.h>
26#include <stdlib.h>
27#include <float.h>
28
29#include "mpfr-test.h"
30
31/* If the precisions are the same, we want to test both mpfr_add1sp
32   and mpfr_add1. */
33
34static int usesp;
35
36static int
37test_add (mpfr_ptr a, mpfr_srcptr b, mpfr_srcptr c, mpfr_rnd_t rnd_mode)
38{
39  int res;
40#ifdef CHECK_EXTERNAL
41  int ok = rnd_mode == MPFR_RNDN && mpfr_number_p (b) && mpfr_number_p (c);
42  if (ok)
43    {
44      mpfr_print_raw (b);
45      printf (" ");
46      mpfr_print_raw (c);
47    }
48#endif
49  if (usesp || MPFR_ARE_SINGULAR(b,c) || MPFR_SIGN(b) != MPFR_SIGN(c))
50    res = mpfr_add (a, b, c, rnd_mode);
51  else
52    {
53      if (MPFR_GET_EXP(b) < MPFR_GET_EXP(c))
54        res = mpfr_add1(a, c, b, rnd_mode);
55      else
56        res = mpfr_add1(a, b, c, rnd_mode);
57    }
58#ifdef CHECK_EXTERNAL
59  if (ok)
60    {
61      printf (" ");
62      mpfr_print_raw (a);
63      printf ("\n");
64    }
65#endif
66  return res;
67}
68
69/* checks that xs+ys gives the expected result zs */
70static void
71check (const char *xs, const char *ys, mpfr_rnd_t rnd_mode,
72        unsigned int px, unsigned int py, unsigned int pz, const char *zs)
73{
74  mpfr_t xx,yy,zz;
75
76  mpfr_init2 (xx, px);
77  mpfr_init2 (yy, py);
78  mpfr_init2 (zz, pz);
79
80  mpfr_set_str1 (xx, xs);
81  mpfr_set_str1 (yy, ys);
82  test_add (zz, xx, yy, rnd_mode);
83  if (mpfr_cmp_str1 (zz, zs) )
84    {
85      printf ("expected sum is %s, got ", zs);
86      mpfr_out_str(stdout, 10, 0, zz, MPFR_RNDN);
87      printf ("mpfr_add failed for x=%s y=%s with rnd_mode=%s\n",
88              xs, ys, mpfr_print_rnd_mode (rnd_mode));
89      exit (1);
90    }
91  mpfr_clears (xx, yy, zz, (mpfr_ptr) 0);
92}
93
94static void
95check2b (const char *xs, int px,
96         const char *ys, int py,
97         const char *rs, int pz,
98         mpfr_rnd_t rnd_mode)
99{
100  mpfr_t xx, yy, zz;
101
102  mpfr_init2 (xx,px);
103  mpfr_init2 (yy,py);
104  mpfr_init2 (zz,pz);
105  mpfr_set_str_binary (xx, xs);
106  mpfr_set_str_binary (yy, ys);
107  test_add (zz, xx, yy, rnd_mode);
108  if (mpfr_cmp_str (zz, rs, 2, MPFR_RNDN))
109    {
110      printf ("(2) x=%s,%d y=%s,%d pz=%d,rnd=%s\n",
111              xs, px, ys, py, pz, mpfr_print_rnd_mode (rnd_mode));
112      printf ("got        "); mpfr_print_binary(zz); puts ("");
113      mpfr_set_str(zz, rs, 2, MPFR_RNDN);
114      printf ("instead of "); mpfr_print_binary(zz); puts ("");
115      exit (1);
116    }
117  mpfr_clear(xx); mpfr_clear(yy); mpfr_clear(zz);
118}
119
120static void
121check64 (void)
122{
123  mpfr_t x, t, u;
124
125  mpfr_init (x);
126  mpfr_init (t);
127  mpfr_init (u);
128
129  mpfr_set_prec (x, 29);
130  mpfr_set_str_binary (x, "1.1101001000101111011010010110e-3");
131  mpfr_set_prec (t, 58);
132  mpfr_set_str_binary (t, "0.11100010011111001001100110010111110110011000000100101E-1");
133  mpfr_set_prec (u, 29);
134  test_add (u, x, t, MPFR_RNDD);
135  mpfr_set_str_binary (t, "1.0101011100001000011100111110e-1");
136  if (mpfr_cmp (u, t))
137    {
138      printf ("mpfr_add(u, x, t) failed for prec(x)=29, prec(t)=58\n");
139      printf ("expected "); mpfr_out_str (stdout, 2, 29, t, MPFR_RNDN);
140      puts ("");
141      printf ("got      "); mpfr_out_str (stdout, 2, 29, u, MPFR_RNDN);
142      puts ("");
143      exit(1);
144    }
145
146  mpfr_set_prec (x, 4);
147  mpfr_set_str_binary (x, "-1.0E-2");
148  mpfr_set_prec (t, 2);
149  mpfr_set_str_binary (t, "-1.1e-2");
150  mpfr_set_prec (u, 2);
151  test_add (u, x, t, MPFR_RNDN);
152  if (MPFR_MANT(u)[0] << 2)
153    {
154      printf ("result not normalized for prec=2\n");
155      mpfr_print_binary (u); puts ("");
156      exit (1);
157    }
158  mpfr_set_str_binary (t, "-1.0e-1");
159  if (mpfr_cmp (u, t))
160    {
161      printf ("mpfr_add(u, x, t) failed for prec(x)=4, prec(t)=2\n");
162      printf ("expected -1.0e-1\n");
163      printf ("got      "); mpfr_out_str (stdout, 2, 4, u, MPFR_RNDN);
164      puts ("");
165      exit (1);
166    }
167
168  mpfr_set_prec (x, 8);
169  mpfr_set_str_binary (x, "-0.10011010"); /* -77/128 */
170  mpfr_set_prec (t, 4);
171  mpfr_set_str_binary (t, "-1.110e-5"); /* -7/128 */
172  mpfr_set_prec (u, 4);
173  test_add (u, x, t, MPFR_RNDN); /* should give -5/8 */
174  mpfr_set_str_binary (t, "-1.010e-1");
175  if (mpfr_cmp (u, t)) {
176    printf ("mpfr_add(u, x, t) failed for prec(x)=8, prec(t)=4\n");
177    printf ("expected -1.010e-1\n");
178    printf ("got      "); mpfr_out_str (stdout, 2, 4, u, MPFR_RNDN);
179    puts ("");
180    exit (1);
181  }
182
183  mpfr_set_prec (x, 112); mpfr_set_prec (t, 98); mpfr_set_prec (u, 54);
184  mpfr_set_str_binary (x, "-0.11111100100000000011000011100000101101010001000111E-401");
185  mpfr_set_str_binary (t, "0.10110000100100000101101100011111111011101000111000101E-464");
186  test_add (u, x, t, MPFR_RNDN);
187  if (mpfr_cmp (u, x))
188    {
189      printf ("mpfr_add(u, x, t) failed for prec(x)=112, prec(t)=98\n");
190      exit (1);
191    }
192
193  mpfr_set_prec (x, 92); mpfr_set_prec (t, 86); mpfr_set_prec (u, 53);
194  mpfr_set_str (x, "-5.03525136761487735093e-74", 10, MPFR_RNDN);
195  mpfr_set_str (t, "8.51539046314262304109e-91", 10, MPFR_RNDN);
196  test_add (u, x, t, MPFR_RNDN);
197  if (mpfr_cmp_str1 (u, "-5.0352513676148773509283672e-74") )
198    {
199      printf ("mpfr_add(u, x, t) failed for prec(x)=92, prec(t)=86\n");
200      exit (1);
201    }
202
203  mpfr_set_prec(x, 53); mpfr_set_prec(t, 76); mpfr_set_prec(u, 76);
204  mpfr_set_str_binary(x, "-0.10010010001001011011110000000000001010011011011110001E-32");
205  mpfr_set_str_binary(t, "-0.1011000101110010000101111111011111010001110011110111100110101011110010011111");
206  mpfr_sub(u, x, t, MPFR_RNDU);
207  mpfr_set_str_binary(t, "0.1011000101110010000101111111011100111111101010011011110110101011101000000100");
208  if (mpfr_cmp(u,t))
209    {
210      printf ("expect "); mpfr_print_binary(t); puts ("");
211      printf ("mpfr_add failed for precisions 53-76\n");
212      exit (1);
213    }
214  mpfr_set_prec(x, 53); mpfr_set_prec(t, 108); mpfr_set_prec(u, 108);
215  mpfr_set_str_binary(x, "-0.10010010001001011011110000000000001010011011011110001E-32");
216  mpfr_set_str_binary(t, "-0.101100010111001000010111111101111101000111001111011110011010101111001001111000111011001110011000000000111111");
217  mpfr_sub(u, x, t, MPFR_RNDU);
218  mpfr_set_str_binary(t, "0.101100010111001000010111111101110011111110101001101111011010101110100000001011000010101110011000000000111111");
219  if (mpfr_cmp(u,t))
220    {
221      printf ("expect "); mpfr_print_binary(t); puts ("");
222      printf ("mpfr_add failed for precisions 53-108\n");
223      exit (1);
224    }
225  mpfr_set_prec(x, 97); mpfr_set_prec(t, 97); mpfr_set_prec(u, 97);
226  mpfr_set_str_binary(x, "0.1111101100001000000001011000110111101000001011111000100001000101010100011111110010000000000000000E-39");
227  mpfr_set_ui(t, 1, MPFR_RNDN);
228  test_add (u, x, t, MPFR_RNDN);
229  mpfr_set_str_binary(x, "0.1000000000000000000000000000000000000000111110110000100000000101100011011110100000101111100010001E1");
230  if (mpfr_cmp(u,x))
231    {
232      printf ("mpfr_add failed for precision 97\n");
233      exit (1);
234    }
235  mpfr_set_prec(x, 128); mpfr_set_prec(t, 128); mpfr_set_prec(u, 128);
236  mpfr_set_str_binary(x, "0.10101011111001001010111011001000101100111101000000111111111011010100001100011101010001010111111101111010100110111111100101100010E-4");
237  mpfr_set(t, x, MPFR_RNDN);
238  mpfr_sub(u, x, t, MPFR_RNDN);
239  mpfr_set_prec(x, 96); mpfr_set_prec(t, 96); mpfr_set_prec(u, 96);
240  mpfr_set_str_binary(x, "0.111000000001110100111100110101101001001010010011010011100111100011010100011001010011011011000010E-4");
241  mpfr_set(t, x, MPFR_RNDN);
242  mpfr_sub(u, x, t, MPFR_RNDN);
243  mpfr_set_prec(x, 85); mpfr_set_prec(t, 85); mpfr_set_prec(u, 85);
244  mpfr_set_str_binary(x, "0.1111101110100110110110100010101011101001100010100011110110110010010011101100101111100E-4");
245  mpfr_set_str_binary(t, "0.1111101110100110110110100010101001001000011000111000011101100101110100001110101010110E-4");
246  mpfr_sub(u, x, t, MPFR_RNDU);
247  mpfr_sub(x, x, t, MPFR_RNDU);
248  if (mpfr_cmp(x, u) != 0)
249    {
250      printf ("Error in mpfr_sub: u=x-t and x=x-t give different results\n");
251      exit (1);
252    }
253  if ((MPFR_MANT(u)[(MPFR_PREC(u)-1)/mp_bits_per_limb] &
254       ((mp_limb_t)1<<(mp_bits_per_limb-1)))==0)
255    {
256      printf ("Error in mpfr_sub: result is not msb-normalized (1)\n");
257      exit (1);
258    }
259  mpfr_set_prec(x, 65); mpfr_set_prec(t, 65); mpfr_set_prec(u, 65);
260  mpfr_set_str_binary(x, "0.10011010101000110101010000000011001001001110001011101011111011101E623");
261  mpfr_set_str_binary(t, "0.10011010101000110101010000000011001001001110001011101011111011100E623");
262  mpfr_sub(u, x, t, MPFR_RNDU);
263  if (mpfr_cmp_ui_2exp(u, 1, 558))
264    { /* 2^558 */
265      printf ("Error (1) in mpfr_sub\n");
266      exit (1);
267    }
268
269  mpfr_set_prec(x, 64); mpfr_set_prec(t, 64); mpfr_set_prec(u, 64);
270  mpfr_set_str_binary(x, "0.1000011110101111011110111111000011101011101111101101101100000100E-220");
271  mpfr_set_str_binary(t, "0.1000011110101111011110111111000011101011101111101101010011111101E-220");
272  test_add (u, x, t, MPFR_RNDU);
273  if ((MPFR_MANT(u)[0] & 1) != 1)
274    {
275      printf ("error in mpfr_add with rnd_mode=MPFR_RNDU\n");
276      printf ("b=  "); mpfr_print_binary(x); puts ("");
277      printf ("c=  "); mpfr_print_binary(t); puts ("");
278      printf ("b+c="); mpfr_print_binary(u); puts ("");
279      exit (1);
280    }
281
282  /* bug found by Norbert Mueller, 14 Sep 2000 */
283  mpfr_set_prec(x, 56); mpfr_set_prec(t, 83); mpfr_set_prec(u, 10);
284  mpfr_set_str_binary(x, "0.10001001011011001111101100110100000101111010010111010111E-7");
285  mpfr_set_str_binary(t, "0.10001001011011001111101100110100000101111010010111010111000000000111110110110000100E-7");
286  mpfr_sub(u, x, t, MPFR_RNDU);
287
288  /* array bound write found by Norbert Mueller, 26 Sep 2000 */
289  mpfr_set_prec(x, 109); mpfr_set_prec(t, 153); mpfr_set_prec(u, 95);
290  mpfr_set_str_binary(x,"0.1001010000101011101100111000110001111111111111111111111111111111111111111111111111111111111111100000000000000E33");
291  mpfr_set_str_binary(t,"-0.100101000010101110110011100011000000000000000000000000000000000000000000000000000000000000000000000000000000000000000011100101101000000100100001100110111E33");
292  test_add (u, x, t, MPFR_RNDN);
293
294  /* array bound writes found by Norbert Mueller, 27 Sep 2000 */
295  mpfr_set_prec(x, 106); mpfr_set_prec(t, 53); mpfr_set_prec(u, 23);
296  mpfr_set_str_binary(x, "-0.1000011110101111111001010001000100001011000000000000000000000000000000000000000000000000000000000000000000E-59");
297  mpfr_set_str_binary(t, "-0.10000111101011111110010100010001101100011100110100000E-59");
298  mpfr_sub(u, x, t, MPFR_RNDN);
299  mpfr_set_prec(x, 177); mpfr_set_prec(t, 217); mpfr_set_prec(u, 160);
300  mpfr_set_str_binary(x, "-0.111010001011010000111001001010010000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000E35");
301  mpfr_set_str_binary(t, "0.1110100010110100001110010010100100000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000111011010011100001111001E35");
302  test_add (u, x, t, MPFR_RNDN);
303  mpfr_set_prec(x, 214); mpfr_set_prec(t, 278); mpfr_set_prec(u, 207);
304  mpfr_set_str_binary(x, "0.1000100110100110101101101101000000010000100111000001001110001000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000E66");
305  mpfr_set_str_binary(t, "-0.10001001101001101011011011010000000100001001110000010011100010000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000001111011111001001100011E66");
306  test_add (u, x, t, MPFR_RNDN);
307  mpfr_set_prec(x, 32); mpfr_set_prec(t, 247); mpfr_set_prec(u, 223);
308  mpfr_set_str_binary(x, "0.10000000000000000000000000000000E1");
309  mpfr_set_str_binary(t, "0.1111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111100000110001110100000100011110000101110110011101110100110110111111011010111100100000000000000000000000000E0");
310  mpfr_sub(u, x, t, MPFR_RNDN);
311  if ((MPFR_MANT(u)[(MPFR_PREC(u)-1)/mp_bits_per_limb] &
312       ((mp_limb_t)1<<(mp_bits_per_limb-1)))==0)
313    {
314      printf ("Error in mpfr_sub: result is not msb-normalized (2)\n");
315      exit (1);
316    }
317
318  /* bug found by Nathalie Revol, 21 March 2001 */
319  mpfr_set_prec (x, 65);
320  mpfr_set_prec (t, 65);
321  mpfr_set_prec (u, 65);
322  mpfr_set_str_binary (x, "0.11100100101101001100111011111111110001101001000011101001001010010E-35");
323  mpfr_set_str_binary (t, "0.10000000000000000000000000000000000001110010010110100110011110000E1");
324  mpfr_sub (u, t, x, MPFR_RNDU);
325  if ((MPFR_MANT(u)[(MPFR_PREC(u)-1)/mp_bits_per_limb] &
326       ((mp_limb_t)1<<(mp_bits_per_limb-1)))==0)
327    {
328      printf ("Error in mpfr_sub: result is not msb-normalized (3)\n");
329      exit (1);
330    }
331
332  /* bug found by Fabrice Rouillier, 27 Mar 2001 */
333  mpfr_set_prec (x, 107);
334  mpfr_set_prec (t, 107);
335  mpfr_set_prec (u, 107);
336  mpfr_set_str_binary (x, "0.10111001001111010010001000000010111111011011011101000001001000101000000000000000000000000000000000000000000E315");
337  mpfr_set_str_binary (t, "0.10000000000000000000000000000000000101110100100101110110000001100101011111001000011101111100100100111011000E350");
338  mpfr_sub (u, x, t, MPFR_RNDU);
339  if ((MPFR_MANT(u)[(MPFR_PREC(u)-1)/mp_bits_per_limb] &
340       ((mp_limb_t)1<<(mp_bits_per_limb-1)))==0)
341    {
342      printf ("Error in mpfr_sub: result is not msb-normalized (4)\n");
343      exit (1);
344    }
345
346  /* checks that NaN flag is correctly reset */
347  mpfr_set_ui (t, 1, MPFR_RNDN);
348  mpfr_set_ui (u, 1, MPFR_RNDN);
349  mpfr_set_nan (x);
350  test_add (x, t, u, MPFR_RNDN);
351  if (mpfr_cmp_ui (x, 2))
352    {
353      printf ("Error in mpfr_add: 1+1 gives ");
354      mpfr_out_str(stdout, 10, 0, x, MPFR_RNDN);
355      exit (1);
356    }
357
358  mpfr_clear(x); mpfr_clear(t); mpfr_clear(u);
359}
360
361/* check case when c does not overlap with a, but both b and c count
362   for rounding */
363static void
364check_case_1b (void)
365{
366  mpfr_t a, b, c;
367  unsigned int prec_a, prec_b, prec_c, dif;
368
369  mpfr_init (a);
370  mpfr_init (b);
371  mpfr_init (c);
372
373    {
374      prec_a = MPFR_PREC_MIN + (randlimb () % 63);
375      mpfr_set_prec (a, prec_a);
376      for (prec_b = prec_a + 2; prec_b <= 64; prec_b++)
377        {
378          dif = prec_b - prec_a;
379          mpfr_set_prec (b, prec_b);
380          /* b = 1 - 2^(-prec_a) + 2^(-prec_b) */
381          mpfr_set_ui (b, 1, MPFR_RNDN);
382          mpfr_div_2exp (b, b, dif, MPFR_RNDN);
383          mpfr_sub_ui (b, b, 1, MPFR_RNDN);
384          mpfr_div_2exp (b, b, prec_a, MPFR_RNDN);
385          mpfr_add_ui (b, b, 1, MPFR_RNDN);
386          for (prec_c = dif; prec_c <= 64; prec_c++)
387            {
388              /* c = 2^(-prec_a) - 2^(-prec_b) */
389              mpfr_set_prec (c, prec_c);
390              mpfr_set_si (c, -1, MPFR_RNDN);
391              mpfr_div_2exp (c, c, dif, MPFR_RNDN);
392              mpfr_add_ui (c, c, 1, MPFR_RNDN);
393              mpfr_div_2exp (c, c, prec_a, MPFR_RNDN);
394              test_add (a, b, c, MPFR_RNDN);
395              if (mpfr_cmp_ui (a, 1) != 0)
396                {
397                  printf ("case (1b) failed for prec_a=%u, prec_b=%u,"
398                          " prec_c=%u\n", prec_a, prec_b, prec_c);
399                  printf ("b="); mpfr_print_binary(b); puts ("");
400                  printf ("c="); mpfr_print_binary(c); puts ("");
401                  printf ("a="); mpfr_print_binary(a); puts ("");
402                  exit (1);
403                }
404            }
405        }
406    }
407
408  mpfr_clear (a);
409  mpfr_clear (b);
410  mpfr_clear (c);
411}
412
413/* check case when c overlaps with a */
414static void
415check_case_2 (void)
416{
417  mpfr_t a, b, c, d;
418
419  mpfr_init2 (a, 300);
420  mpfr_init2 (b, 800);
421  mpfr_init2 (c, 500);
422  mpfr_init2 (d, 800);
423
424  mpfr_set_str_binary(a, "1E110");  /* a = 2^110 */
425  mpfr_set_str_binary(b, "1E900");  /* b = 2^900 */
426  mpfr_set_str_binary(c, "1E500");  /* c = 2^500 */
427  test_add (c, c, a, MPFR_RNDZ);   /* c = 2^500 + 2^110 */
428  mpfr_sub (d, b, c, MPFR_RNDZ);   /* d = 2^900 - 2^500 - 2^110 */
429  test_add (b, b, c, MPFR_RNDZ);   /* b = 2^900 + 2^500 + 2^110 */
430  test_add (a, b, d, MPFR_RNDZ);   /* a = 2^901 */
431  if (mpfr_cmp_ui_2exp (a, 1, 901))
432    {
433      printf ("b + d fails for b=2^900+2^500+2^110, d=2^900-2^500-2^110\n");
434      printf ("expected 1.0e901, got ");
435      mpfr_out_str (stdout, 2, 0, a, MPFR_RNDN);
436      printf ("\n");
437      exit (1);
438    }
439
440  mpfr_clear (a);
441  mpfr_clear (b);
442  mpfr_clear (c);
443  mpfr_clear (d);
444}
445
446/* checks when source and destination are equal */
447static void
448check_same (void)
449{
450  mpfr_t x;
451
452  mpfr_init(x); mpfr_set_ui(x, 1, MPFR_RNDZ);
453  test_add (x, x, x, MPFR_RNDZ);
454  if (mpfr_cmp_ui (x, 2))
455    {
456      printf ("Error when all 3 operands are equal\n");
457      exit (1);
458    }
459  mpfr_clear(x);
460}
461
462#define check53(x, y, r, z) check(x, y, r, 53, 53, 53, z)
463
464#define MAX_PREC 256
465
466static void
467check_inexact (void)
468{
469  mpfr_t x, y, z, u;
470  mpfr_prec_t px, py, pu, pz;
471  int inexact, cmp;
472  mpfr_rnd_t rnd;
473
474  mpfr_init (x);
475  mpfr_init (y);
476  mpfr_init (z);
477  mpfr_init (u);
478
479  mpfr_set_prec (x, 2);
480  mpfr_set_str_binary (x, "0.1E-4");
481  mpfr_set_prec (u, 33);
482  mpfr_set_str_binary (u, "0.101110100101101100000000111100000E-1");
483  mpfr_set_prec (y, 31);
484  if ((inexact = test_add (y, x, u, MPFR_RNDN)))
485    {
486      printf ("Wrong inexact flag (2): expected 0, got %d\n", inexact);
487      exit (1);
488    }
489
490  mpfr_set_prec (x, 2);
491  mpfr_set_str_binary (x, "0.1E-4");
492  mpfr_set_prec (u, 33);
493  mpfr_set_str_binary (u, "0.101110100101101100000000111100000E-1");
494  mpfr_set_prec (y, 28);
495  if ((inexact = test_add (y, x, u, MPFR_RNDN)))
496    {
497      printf ("Wrong inexact flag (1): expected 0, got %d\n", inexact);
498      exit (1);
499    }
500
501  for (px=2; px<MAX_PREC; px++)
502    {
503      mpfr_set_prec (x, px);
504      do
505        {
506          mpfr_urandomb (x, RANDS);
507        }
508      while (mpfr_cmp_ui (x, 0) == 0);
509      for (pu=2; pu<MAX_PREC; pu++)
510        {
511          mpfr_set_prec (u, pu);
512          do
513            {
514              mpfr_urandomb (u, RANDS);
515            }
516          while (mpfr_cmp_ui (u, 0) == 0);
517          {
518              py = MPFR_PREC_MIN + (randlimb () % (MAX_PREC - 1));
519              mpfr_set_prec (y, py);
520              pz =  (mpfr_cmpabs (x, u) >= 0) ? MPFR_EXP(x) - MPFR_EXP(u)
521                : MPFR_EXP(u) - MPFR_EXP(x);
522              /* x + u is exactly representable with precision
523                 abs(EXP(x)-EXP(u)) + max(prec(x), prec(u)) + 1 */
524              pz = pz + MAX(MPFR_PREC(x), MPFR_PREC(u)) + 1;
525              mpfr_set_prec (z, pz);
526              rnd = RND_RAND ();
527              if (test_add (z, x, u, rnd))
528                {
529                  printf ("z <- x + u should be exact\n");
530                  printf ("x="); mpfr_print_binary (x); puts ("");
531                  printf ("u="); mpfr_print_binary (u); puts ("");
532                  printf ("z="); mpfr_print_binary (z); puts ("");
533                  exit (1);
534                }
535                {
536                  rnd = RND_RAND ();
537                  inexact = test_add (y, x, u, rnd);
538                  cmp = mpfr_cmp (y, z);
539                  if (((inexact == 0) && (cmp != 0)) ||
540                      ((inexact > 0) && (cmp <= 0)) ||
541                      ((inexact < 0) && (cmp >= 0)))
542                    {
543                      printf ("Wrong inexact flag for rnd=%s\n",
544                              mpfr_print_rnd_mode(rnd));
545                      printf ("expected %d, got %d\n", cmp, inexact);
546                      printf ("x="); mpfr_print_binary (x); puts ("");
547                      printf ("u="); mpfr_print_binary (u); puts ("");
548                      printf ("y=  "); mpfr_print_binary (y); puts ("");
549                      printf ("x+u="); mpfr_print_binary (z); puts ("");
550                      exit (1);
551                    }
552                }
553            }
554        }
555    }
556
557  mpfr_clear (x);
558  mpfr_clear (y);
559  mpfr_clear (z);
560  mpfr_clear (u);
561}
562
563static void
564check_nans (void)
565{
566  mpfr_t  s, x, y;
567
568  mpfr_init2 (x, 8L);
569  mpfr_init2 (y, 8L);
570  mpfr_init2 (s, 8L);
571
572  /* +inf + -inf == nan */
573  mpfr_set_inf (x, 1);
574  mpfr_set_inf (y, -1);
575  test_add (s, x, y, MPFR_RNDN);
576  MPFR_ASSERTN (mpfr_nan_p (s));
577
578  /* +inf + 1 == +inf */
579  mpfr_set_inf (x, 1);
580  mpfr_set_ui (y, 1L, MPFR_RNDN);
581  test_add (s, x, y, MPFR_RNDN);
582  MPFR_ASSERTN (mpfr_inf_p (s));
583  MPFR_ASSERTN (mpfr_sgn (s) > 0);
584
585  /* -inf + 1 == -inf */
586  mpfr_set_inf (x, -1);
587  mpfr_set_ui (y, 1L, MPFR_RNDN);
588  test_add (s, x, y, MPFR_RNDN);
589  MPFR_ASSERTN (mpfr_inf_p (s));
590  MPFR_ASSERTN (mpfr_sgn (s) < 0);
591
592  /* 1 + +inf == +inf */
593  mpfr_set_ui (x, 1L, MPFR_RNDN);
594  mpfr_set_inf (y, 1);
595  test_add (s, x, y, MPFR_RNDN);
596  MPFR_ASSERTN (mpfr_inf_p (s));
597  MPFR_ASSERTN (mpfr_sgn (s) > 0);
598
599  /* 1 + -inf == -inf */
600  mpfr_set_ui (x, 1L, MPFR_RNDN);
601  mpfr_set_inf (y, -1);
602  test_add (s, x, y, MPFR_RNDN);
603  MPFR_ASSERTN (mpfr_inf_p (s));
604  MPFR_ASSERTN (mpfr_sgn (s) < 0);
605
606  mpfr_clear (x);
607  mpfr_clear (y);
608  mpfr_clear (s);
609}
610
611static void
612check_alloc (void)
613{
614  mpfr_t a;
615
616  mpfr_init2 (a, 10000);
617  mpfr_set_prec (a, 53);
618  mpfr_set_ui (a, 15236, MPFR_RNDN);
619  test_add (a, a, a, MPFR_RNDN);
620  mpfr_mul (a, a, a, MPFR_RNDN);
621  mpfr_div (a, a, a, MPFR_RNDN);
622  mpfr_sub (a, a, a, MPFR_RNDN);
623  mpfr_clear (a);
624}
625
626static void
627check_overflow (void)
628{
629  mpfr_t a, b, c;
630  mpfr_prec_t prec_a;
631  int r;
632
633  mpfr_init2 (a, 256);
634  mpfr_init2 (b, 256);
635  mpfr_init2 (c, 256);
636
637  mpfr_set_ui (b, 1, MPFR_RNDN);
638  mpfr_setmax (b, mpfr_get_emax ());
639  mpfr_set_ui (c, 1, MPFR_RNDN);
640  mpfr_set_exp (c, mpfr_get_emax () - 192);
641  RND_LOOP(r)
642    for (prec_a = 128; prec_a < 512; prec_a += 64)
643      {
644        mpfr_set_prec (a, prec_a);
645        mpfr_clear_overflow ();
646        test_add (a, b, c, (mpfr_rnd_t) r);
647        if (!mpfr_overflow_p ())
648          {
649            printf ("No overflow in check_overflow\n");
650            exit (1);
651          }
652      }
653
654  mpfr_set_exp (c, mpfr_get_emax () - 512);
655  mpfr_set_prec (a, 256);
656  mpfr_clear_overflow ();
657  test_add (a, b, c, MPFR_RNDU);
658  if (!mpfr_overflow_p ())
659    {
660      printf ("No overflow in check_overflow\n");
661      exit (1);
662    }
663
664  mpfr_clear (a);
665  mpfr_clear (b);
666  mpfr_clear (c);
667}
668
669static void
670check_1111 (void)
671{
672  mpfr_t one;
673  long n;
674
675  mpfr_init2 (one, MPFR_PREC_MIN);
676  mpfr_set_ui (one, 1, MPFR_RNDN);
677  for (n = 0; n < NUM; n++)
678    {
679      mpfr_prec_t prec_a, prec_b, prec_c;
680      mpfr_exp_t tb=0, tc, diff;
681      mpfr_t a, b, c, s;
682      int m = 512;
683      int sb, sc;
684      int inex_a, inex_s;
685      mpfr_rnd_t rnd_mode;
686
687      prec_a = MPFR_PREC_MIN + (randlimb () % m);
688      prec_b = MPFR_PREC_MIN + (randlimb () % m);
689      prec_c = MPFR_PREC_MIN + (randlimb () % m);
690      mpfr_init2 (a, prec_a);
691      mpfr_init2 (b, prec_b);
692      mpfr_init2 (c, prec_c);
693      sb = randlimb () % 3;
694      if (sb != 0)
695        {
696          tb = 1 + (randlimb () % (prec_b - (sb != 2)));
697          mpfr_div_2ui (b, one, tb, MPFR_RNDN);
698          if (sb == 2)
699            mpfr_neg (b, b, MPFR_RNDN);
700          test_add (b, b, one, MPFR_RNDN);
701        }
702      else
703        mpfr_set (b, one, MPFR_RNDN);
704      tc = 1 + (randlimb () % (prec_c - 1));
705      mpfr_div_2ui (c, one, tc, MPFR_RNDN);
706      sc = randlimb () % 2;
707      if (sc)
708        mpfr_neg (c, c, MPFR_RNDN);
709      test_add (c, c, one, MPFR_RNDN);
710      diff = (randlimb () % (2*m)) - m;
711      mpfr_mul_2si (c, c, diff, MPFR_RNDN);
712      rnd_mode = RND_RAND ();
713      inex_a = test_add (a, b, c, rnd_mode);
714      mpfr_init2 (s, MPFR_PREC_MIN + 2*m);
715      inex_s = test_add (s, b, c, MPFR_RNDN); /* exact */
716      if (inex_s)
717        {
718          printf ("check_1111: result should have been exact.\n");
719          exit (1);
720        }
721      inex_s = mpfr_prec_round (s, prec_a, rnd_mode);
722      if ((inex_a < 0 && inex_s >= 0) ||
723          (inex_a == 0 && inex_s != 0) ||
724          (inex_a > 0 && inex_s <= 0) ||
725          !mpfr_equal_p (a, s))
726        {
727          printf ("check_1111: results are different.\n");
728          printf ("prec_a = %d, prec_b = %d, prec_c = %d\n",
729                  (int) prec_a, (int) prec_b, (int) prec_c);
730          printf ("tb = %d, tc = %d, diff = %d, rnd = %s\n",
731                  (int) tb, (int) tc, (int) diff,
732                  mpfr_print_rnd_mode (rnd_mode));
733          printf ("sb = %d, sc = %d\n", sb, sc);
734          printf ("a = "); mpfr_print_binary (a); puts ("");
735          printf ("s = "); mpfr_print_binary (s); puts ("");
736          printf ("inex_a = %d, inex_s = %d\n", inex_a, inex_s);
737          exit (1);
738        }
739      mpfr_clear (a);
740      mpfr_clear (b);
741      mpfr_clear (c);
742      mpfr_clear (s);
743    }
744  mpfr_clear (one);
745}
746
747static void
748check_1minuseps (void)
749{
750  static mpfr_prec_t prec_a[] = {
751    MPFR_PREC_MIN, 30, 31, 32, 33, 62, 63, 64, 65, 126, 127, 128, 129
752  };
753  static int supp_b[] = {
754    0, 1, 2, 3, 4, 29, 30, 31, 32, 33, 34, 35, 61, 62, 63, 64, 65, 66, 67
755  };
756  mpfr_t a, b, c;
757  unsigned int ia, ib, ic;
758
759  mpfr_init2 (c, MPFR_PREC_MIN);
760
761  for (ia = 0; ia < numberof (prec_a); ia++)
762    for (ib = 0; ib < numberof(supp_b); ib++)
763      {
764        mpfr_prec_t prec_b;
765        int rnd_mode;
766
767        prec_b = prec_a[ia] + supp_b[ib];
768
769        mpfr_init2 (a, prec_a[ia]);
770        mpfr_init2 (b, prec_b);
771
772        mpfr_set_ui (c, 1, MPFR_RNDN);
773        mpfr_div_ui (b, c, prec_a[ia], MPFR_RNDN);
774        mpfr_sub (b, c, b, MPFR_RNDN);  /* b = 1 - 2^(-prec_a) */
775
776        for (ic = 0; ic < numberof(supp_b); ic++)
777          for (rnd_mode = 0; rnd_mode < MPFR_RND_MAX; rnd_mode++)
778            {
779              mpfr_t s;
780              int inex_a, inex_s;
781
782              mpfr_set_ui (c, 1, MPFR_RNDN);
783              mpfr_div_ui (c, c, prec_a[ia] + supp_b[ic], MPFR_RNDN);
784              inex_a = test_add (a, b, c, (mpfr_rnd_t) rnd_mode);
785              mpfr_init2 (s, 256);
786              inex_s = test_add (s, b, c, MPFR_RNDN); /* exact */
787              if (inex_s)
788                {
789                  printf ("check_1minuseps: result should have been exact "
790                          "(ia = %u, ib = %u, ic = %u)\n", ia, ib, ic);
791                  exit (1);
792                }
793              inex_s = mpfr_prec_round (s, prec_a[ia], (mpfr_rnd_t) rnd_mode);
794              if ((inex_a < 0 && inex_s >= 0) ||
795                  (inex_a == 0 && inex_s != 0) ||
796                  (inex_a > 0 && inex_s <= 0) ||
797                  !mpfr_equal_p (a, s))
798                {
799                  printf ("check_1minuseps: results are different.\n");
800                  printf ("ia = %u, ib = %u, ic = %u\n", ia, ib, ic);
801                  exit (1);
802                }
803              mpfr_clear (s);
804            }
805
806        mpfr_clear (a);
807        mpfr_clear (b);
808      }
809
810  mpfr_clear (c);
811}
812
813/* Test case bk == 0 in add1.c (b has entirely been read and
814   c hasn't been taken into account). */
815static void
816coverage_bk_eq_0 (void)
817{
818  mpfr_t a, b, c;
819  int inex;
820
821  mpfr_init2 (a, GMP_NUMB_BITS);
822  mpfr_init2 (b, 2 * GMP_NUMB_BITS);
823  mpfr_init2 (c, GMP_NUMB_BITS);
824
825  mpfr_set_ui_2exp (b, 1, 2 * GMP_NUMB_BITS, MPFR_RNDN);
826  mpfr_sub_ui (b, b, 1, MPFR_RNDN);
827  /* b = 111...111 (in base 2) where the 1's fit 2 whole limbs */
828
829  mpfr_set_ui_2exp (c, 1, -1, MPFR_RNDN);  /* c = 1/2 */
830
831  inex = mpfr_add (a, b, c, MPFR_RNDU);
832  mpfr_set_ui_2exp (c, 1, 2 * GMP_NUMB_BITS, MPFR_RNDN);
833  if (! mpfr_equal_p (a, c))
834    {
835      printf ("Error in coverage_bk_eq_0\n");
836      printf ("Expected ");
837      mpfr_dump (c);
838      printf ("Got      ");
839      mpfr_dump (a);
840      exit (1);
841    }
842  MPFR_ASSERTN (inex > 0);
843
844  mpfr_clear (a);
845  mpfr_clear (b);
846  mpfr_clear (c);
847}
848
849static void
850tests (void)
851{
852  check_alloc ();
853  check_nans ();
854  check_inexact ();
855  check_case_1b ();
856  check_case_2 ();
857  check64();
858  coverage_bk_eq_0 ();
859
860  check("293607738.0", "1.9967571564050541e-5", MPFR_RNDU, 64, 53, 53,
861        "2.9360773800002003e8");
862  check("880524.0", "-2.0769715792901673e-5", MPFR_RNDN, 64, 53, 53,
863        "8.8052399997923023e5");
864  check("1196426492.0", "-1.4218093058435347e-3", MPFR_RNDN, 64, 53, 53,
865        "1.1964264919985781e9");
866  check("982013018.0", "-8.941829477291838e-7", MPFR_RNDN, 64, 53, 53,
867        "9.8201301799999905e8");
868  check("1092583421.0", "1.0880649218158844e9", MPFR_RNDN, 64, 53, 53,
869        "2.1806483428158846e9");
870  check("1.8476886419022969e-6", "961494401.0", MPFR_RNDN, 53, 64, 53,
871        "9.6149440100000179e8");
872  check("-2.3222118418069868e5", "1229318102.0", MPFR_RNDN, 53, 64, 53,
873        "1.2290858808158193e9");
874  check("-3.0399171300395734e-6", "874924868.0", MPFR_RNDN, 53, 64, 53,
875        "8.749248679999969e8");
876  check("9.064246624706179e1", "663787413.0", MPFR_RNDN, 53, 64, 53,
877        "6.6378750364246619e8");
878  check("-1.0954322421551264e2", "281806592.0", MPFR_RNDD, 53, 64, 53,
879        "2.8180648245677572e8");
880  check("5.9836930386056659e-8", "1016217213.0", MPFR_RNDN, 53, 64, 53,
881        "1.0162172130000001e9");
882  check("-1.2772161928500301e-7", "1237734238.0", MPFR_RNDN, 53, 64, 53,
883        "1.2377342379999998e9");
884  check("-4.567291988483277e8", "1262857194.0", MPFR_RNDN, 53, 64, 53,
885        "8.0612799515167236e8");
886  check("4.7719471752925262e7", "196089880.0", MPFR_RNDN, 53, 53, 53,
887        "2.4380935175292528e8");
888  check("4.7719471752925262e7", "196089880.0", MPFR_RNDN, 53, 64, 53,
889        "2.4380935175292528e8");
890  check("-1.716113812768534e-140", "1271212614.0", MPFR_RNDZ, 53, 64, 53,
891        "1.2712126139999998e9");
892  check("-1.2927455200185474e-50", "1675676122.0", MPFR_RNDD, 53, 64, 53,
893        "1.6756761219999998e9");
894
895  check53("1.22191250737771397120e+20", "948002822.0", MPFR_RNDN,
896          "122191250738719408128.0");
897  check53("9966027674114492.0", "1780341389094537.0", MPFR_RNDN,
898          "11746369063209028.0");
899  check53("2.99280481918991653800e+272", "5.34637717585790933424e+271",
900          MPFR_RNDN, "3.5274425367757071711e272");
901  check_same();
902  check53("6.14384195492641560499e-02", "-6.14384195401037683237e-02",
903          MPFR_RNDU, "9.1603877261370314499e-12");
904  check53("1.16809465359248765399e+196", "7.92883212101990665259e+196",
905          MPFR_RNDU, "9.0969267746123943065e196");
906  check53("3.14553393112021279444e-67", "3.14553401015952024126e-67", MPFR_RNDU,
907          "6.2910679412797336946e-67");
908
909  check53("5.43885304644369509058e+185","-1.87427265794105342763e-57",MPFR_RNDN,
910          "5.4388530464436950905e185");
911  check53("5.43885304644369509058e+185","-1.87427265794105342763e-57",MPFR_RNDZ,
912          "5.4388530464436944867e185");
913  check53("5.43885304644369509058e+185","-1.87427265794105342763e-57",MPFR_RNDU,
914          "5.4388530464436950905e185");
915  check53("5.43885304644369509058e+185","-1.87427265794105342763e-57",MPFR_RNDD,
916          "5.4388530464436944867e185");
917
918  check2b("1.001010101110011000000010100101110010111001010000001000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000e358",187,
919          "-1.11100111001101100010001111111110101101110001000000000000000000000000000000000000000000e160",87,
920          "1.001010101110011000000010100101110010111001010000000111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111e358",178,
921          MPFR_RNDD);
922  check2b("-1.111100100011100111010101010101001010100100000111001000000000000000000e481",70,
923          "1.1111000110100011110101111110110010010000000110101000000000000000e481",65,
924          "-1.001010111111101011010000001100011101100101000000000000000000e472",61,
925          MPFR_RNDD);
926  check2b("1.0100010111010000100101000000111110011100011001011010000000000000000000000000000000e516",83,
927          "-1.1001111000100001011100000001001100110011110010111111000000e541",59,
928          "-1.1001111000100001011011110111000001001011100000011110100000110001110011010011000000000000000000000000000000000000000000000000e541",125,
929          MPFR_RNDZ);
930  check2b("-1.0010111100000100110001011011010000000011000111101000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000e261",155,
931          "-1.00111110100011e239",15,
932          "-1.00101111000001001100101010101110001100110001111010000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000e261",159,
933          MPFR_RNDD);
934  check2b("-1.110111000011111011000000001001111101101001010100111000000000000000000000000e880",76,
935          "-1.1010010e-634",8,
936          "-1.11011100001111101100000000100111110110100101010011100000000000000000000000e880",75,
937          MPFR_RNDZ);
938  check2b("1.00100100110110101001010010101111000001011100100101010000000000000000000000000000e-530",81,
939          "-1.101101111100000111000011001010110011001011101001110100000e-908",58,
940          "1.00100100110110101001010010101111000001011100100101010e-530",54,
941          MPFR_RNDN);
942  check2b("1.0101100010010111101000000001000010010010011000111011000000000000000000000000000000000000000000000000000000000000000000e374",119,
943          "1.11100101100101e358",15,
944          "1.01011000100110011000010110100100100100100110001110110000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000e374",150,
945          MPFR_RNDZ);
946  check2b("-1.10011001000010100000010100100110110010011111101111000000000000000000000000000000000000000000000000000000000000000000e-172",117,
947          "1.111011100000101010110000100100110100100001001000011100000000e-173",61,
948          "-1.0100010000001001010110011011101001001011101011110001000000000000000e-173",68,
949          MPFR_RNDZ);
950  check2b("-1.011110000111101011100001100110100011100101000011011000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000e-189",175,
951          "1.1e631",2,
952          "1.011111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111e631",115,
953          MPFR_RNDZ);
954  check2b("-1.101011101001011101100011001000001100010100001101011000000000000000000000000000000000000000000e-449",94,
955          "-1.01111101010111000011000110011101000111001100110111100000000000000e-429",66,
956          "-1.01111101010111000100110010000110100100101111111111101100010100001101011000000000000000000000000000000000000000e-429",111,
957          MPFR_RNDU);
958  check2b("-1.101011101001011101100011001000001100010100001101011000000000000000000000000000000000000000000e-449",94,
959          "-1.01111101010111000011000110011101000111001100110111100000000000000e-429",66,
960          "-1.01111101010111000100110010000110100100101111111111101100010100001101011000000000000000000000000000000000000000e-429",111,
961          MPFR_RNDD);
962  check2b("-1.1001000011101000110000111110010100100101110101111100000000000000000000000000000000000000000000000000000000e-72",107,
963          "-1.001100011101100100010101101010101011010010111111010000000000000000000000000000e521",79,
964          "-1.00110001110110010001010110101010101101001011111101000000000000000000000000000000000000000000000001e521",99,
965          MPFR_RNDD);
966  check2b("-1.01010001111000000101010100100100110101011011100001110000000000e498",63,
967          "1.010000011010101111000100111100011100010101011110010100000000000e243",64,
968          "-1.010100011110000001010101001001001101010110111000011100000000000e498",64,
969          MPFR_RNDN);
970  check2b("1.00101100010101000011010000011000111111011110010111000000000000000000000000000000000000000000000000000000000e178",108,
971          "-1.10101101010101000110011011111001001101111111110000100000000e160",60,
972          "1.00101100010100111100100011000011111001000010011101110010000000001111100000000000000000000000000000000000e178",105,
973          MPFR_RNDN);
974  check2b("1.00110011010100111110011010110100111101110101100100110000000000000000000000000000000000000000000000e559",99,
975          "-1.011010110100111011100110100110011100000000111010011000000000000000e559",67,
976          "-1.101111111101011111111111001001100100011100001001100000000000000000000000000000000000000000000e556",94,
977          MPFR_RNDU);
978  check2b("-1.100000111100101001100111011100011011000001101001111100000000000000000000000000e843",79,
979          "-1.1101101010110000001001000100001100110011000110110111000000000000000000000000000000000000000000e414",95,
980          "-1.1000001111001010011001110111000110110000011010100000e843",53,
981          MPFR_RNDD);
982  check2b("-1.110110010110100010100011000110111001010000010111110000000000e-415",61,
983          "-1.0000100101100001111100110011111111110100011101101011000000000000000000e751",71,
984          "-1.00001001011000011111001100111111111101000111011010110e751",54,
985          MPFR_RNDN);
986  check2b("-1.1011011011110001001101010101001000010100010110111101000000000000000000000e258",74,
987          "-1.00011100010110110101001011000100100000100010101000010000000000000000000000000000000000000000000000e268",99,
988          "-1.0001110011001001000011110001000111010110101011110010011011110100000000000000000000000000000000000000e268",101,
989          MPFR_RNDD);
990  check2b("-1.1011101010011101011000000100100110101101101110000001000000000e629",62,
991          "1.111111100000011100100011100000011101100110111110111000000000000000000000000000000000000000000e525",94,
992          "-1.101110101001110101100000010010011010110110111000000011111111111111111111111111111111111111111111111111101e629",106,
993          MPFR_RNDD);
994  check2b("1.111001000010001100010000001100000110001011110111011000000000000000000000000000000000000e152",88,
995          "1.111110111001100100000100111111010111000100111111001000000000000000e152",67,
996          "1.1110111111011110000010101001011011101010000110110100e153",53,
997          MPFR_RNDN);
998  check2b("1.000001100011110010110000110100001010101101111011110100e696",55,
999          "-1.1011001111011100100001011110100101010101110111010101000000000000000000000000000000000000000000000000000000000000e730",113,
1000          "-1.1011001111011100100001011110100100010100010011100010e730",53,
1001          MPFR_RNDN);
1002  check2b("-1.11010111100001001111000001110101010010001111111001100000000000000000000000000000000000000000000000000000000000e530",111,
1003          "1.01110100010010000000010110111101011101000001111101100000000000000000000000000000000000000000000000e530",99,
1004          "-1.1000110011110011101010101101111101010011011111000000000000000e528",62,
1005          MPFR_RNDD);
1006  check2b("-1.0001100010010100111101101011101000100100010011100011000000000000000000000000000000000000000000000000000000000e733",110,
1007          "-1.001000000111110010100101010100110111001111011011001000000000000000000000000000000000000000000000000000000000e710",109,
1008          "-1.000110001001010011111000111110110001110110011000110110e733",55,
1009          MPFR_RNDN);
1010  check2b("-1.1101011110000100111100000111010101001000111111100110000000000000000000000e530",74,
1011          "1.01110100010010000000010110111101011101000001111101100000000000000000000000000000000000000000000000000000000000e530",111,
1012          "-1.10001100111100111010101011011111010100110111110000000000000000000000000000e528",75,
1013          MPFR_RNDU);
1014  check2b("1.00110011010100111110011010110100111101110101100100110000000000000000000000000000000000000000000000e559",99,
1015          "-1.011010110100111011100110100110011100000000111010011000000000000000e559",67,
1016          "-1.101111111101011111111111001001100100011100001001100000000000000000000000000000000000000000000e556",94,
1017          MPFR_RNDU);
1018  check2b("-1.100101111110110000000110111111011010011101101111100100000000000000e-624",67,
1019          "1.10111010101110100000010110101000000000010011100000100000000e-587",60,
1020          "1.1011101010111010000001011010011111110100011110001011111111001000000100101100010010000011100000000000000000000e-587",110,
1021          MPFR_RNDU);
1022  check2b("-1.10011001000010100000010100100110110010011111101111000000000000000000000000000000000000000000000000000000000000000000e-172",117,
1023          "1.111011100000101010110000100100110100100001001000011100000000e-173",61,
1024          "-1.0100010000001001010110011011101001001011101011110001000000000000000e-173",68,
1025          MPFR_RNDZ);
1026  check2b("1.1000111000110010101001010011010011101100010110001001000000000000000000000000000000000000000000000000e167",101,
1027          "1.0011110010000110000000101100100111000001110110110000000000000000000000000e167",74,
1028          "1.01100101010111000101001111111111010101110001100111001000000000000000000000000000000000000000000000000000e168",105,
1029          MPFR_RNDZ);
1030  check2b("1.100101111111110010100101110111100001110000100001010000000000000000000000000000000000000000000000e808",97,
1031          "-1.1110011001100000100000111111110000110010100111001011000000000000000000000000000000e807",83,
1032          "1.01001001100110001100011111000000000001011010010111010000000000000000000000000000000000000000000e807",96,
1033          MPFR_RNDN);
1034  check2b("1e128",128,
1035          "1e0",128,
1036          "100000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000001e0",256,
1037          MPFR_RNDN);
1038
1039  /* Checking double precision (53 bits) */
1040  check53("-8.22183238641455905806e-19", "7.42227178769761587878e-19",MPFR_RNDD,
1041          "-7.9956059871694317927e-20");
1042  check53("5.82106394662028628236e+234","-5.21514064202368477230e+89",MPFR_RNDD,
1043          "5.8210639466202855763e234");
1044  check53("5.72931679569871602371e+122","-5.72886070363264321230e+122",
1045          MPFR_RNDN, "4.5609206607281141508e118");
1046  check53("-5.09937369394650450820e+238", "2.70203299854862982387e+250",
1047          MPFR_RNDD, "2.7020329985435301323e250");
1048  check53("-2.96695924472363684394e+27", "1.22842938251111500000e+16",MPFR_RNDD,
1049          "-2.96695924471135255027e27");
1050  check53("1.74693641655743793422e-227", "-7.71776956366861843469e-229",
1051          MPFR_RNDN, "1.669758720920751867e-227");
1052  /*  x = -7883040437021647.0; for (i=0; i<468; i++) x = x / 2.0;*/
1053  check53("-1.03432206392780011159e-125", "1.30127034799251347548e-133",
1054          MPFR_RNDN,
1055          "-1.0343220509150965661100887242027378881805094180354e-125");
1056  check53("1.05824655795525779205e+71", "-1.06022698059744327881e+71",MPFR_RNDZ,
1057          "-1.9804226421854867632e68");
1058  check53("-5.84204911040921732219e+240", "7.26658169050749590763e+240",
1059          MPFR_RNDD, "1.4245325800982785854e240");
1060  check53("1.00944884131046636376e+221","2.33809162651471520268e+215",MPFR_RNDN,
1061          "1.0094511794020929787e221");
1062  /*x = 7045852550057985.0; for (i=0; i<986; i++) x = x / 2.0;*/
1063  check53("4.29232078932667367325e-278",
1064          "1.0773525047389793833221116707010783793203080117586e-281"
1065          , MPFR_RNDU, "4.2933981418314132787e-278");
1066  check53("5.27584773801377058681e-80", "8.91207657803547196421e-91", MPFR_RNDN,
1067          "5.2758477381028917269e-80");
1068  check53("2.99280481918991653800e+272", "5.34637717585790933424e+271",
1069          MPFR_RNDN, "3.5274425367757071711e272");
1070  check53("4.67302514390488041733e-184", "2.18321376145645689945e-190",
1071          MPFR_RNDN, "4.6730273271186420541e-184");
1072  check53("5.57294120336300389254e+71", "2.60596167942024924040e+65", MPFR_RNDZ,
1073          "5.5729438093246831053e71");
1074  check53("6.6052588496951015469e24", "4938448004894539.0", MPFR_RNDU,
1075          "6.6052588546335505068e24");
1076  check53("1.23056185051606761523e-190", "1.64589756643433857138e-181",
1077          MPFR_RNDU, "1.6458975676649006598e-181");
1078  check53("2.93231171510175981584e-280", "3.26266919161341483877e-273",
1079          MPFR_RNDU, "3.2626694848445867288e-273");
1080  check53("5.76707395945001907217e-58", "4.74752971449827687074e-51", MPFR_RNDD,
1081          "4.747530291205672325e-51");
1082  check53("277363943109.0", "11.0", MPFR_RNDN, "277363943120.0");
1083  check53("1.44791789689198883921e-140", "-1.90982880222349071284e-121",
1084          MPFR_RNDN, "-1.90982880222349071e-121");
1085
1086
1087  /* tests for particular cases (Vincent Lefevre, 22 Aug 2001) */
1088  check53("9007199254740992.0", "1.0", MPFR_RNDN, "9007199254740992.0");
1089  check53("9007199254740994.0", "1.0", MPFR_RNDN, "9007199254740996.0");
1090  check53("9007199254740992.0", "-1.0", MPFR_RNDN, "9007199254740991.0");
1091  check53("9007199254740994.0", "-1.0", MPFR_RNDN, "9007199254740992.0");
1092  check53("9007199254740996.0", "-1.0", MPFR_RNDN, "9007199254740996.0");
1093
1094  check_overflow ();
1095  check_1111 ();
1096  check_1minuseps ();
1097}
1098
1099#define TEST_FUNCTION test_add
1100#define TWO_ARGS
1101#define RAND_FUNCTION(x) mpfr_random2(x, MPFR_LIMB_SIZE (x), randlimb () % 100, RANDS)
1102#include "tgeneric.c"
1103
1104int
1105main (int argc, char *argv[])
1106{
1107  tests_start_mpfr ();
1108
1109  usesp = 0;
1110  tests ();
1111
1112#ifndef CHECK_EXTERNAL /* no need to check twice */
1113  usesp = 1;
1114  tests ();
1115#endif
1116  test_generic (2, 1000, 100);
1117
1118  tests_end_mpfr ();
1119  return 0;
1120}
1121