xref: /netbsd-src/external/lgpl3/gmp/dist/tests/mpf/t-sub.c (revision 72c7faa4dbb41dbb0238d6b4a109da0d4b236dd4)
1 /* Test mpf_sub.
2 
3 Copyright 1996, 2001, 2004, 2014 Free Software Foundation, Inc.
4 
5 This file is part of the GNU MP Library test suite.
6 
7 The GNU MP Library test suite is free software; you can redistribute it
8 and/or modify it under the terms of the GNU General Public License as
9 published by the Free Software Foundation; either version 3 of the License,
10 or (at your option) any later version.
11 
12 The GNU MP Library test suite is distributed in the hope that it will be
13 useful, but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General
15 Public License for more details.
16 
17 You should have received a copy of the GNU General Public License along with
18 the GNU MP Library test suite.  If not, see https://www.gnu.org/licenses/.  */
19 
20 #include <stdio.h>
21 #include <stdlib.h>
22 
23 #include "gmp-impl.h"
24 #include "tests.h"
25 
26 #ifndef SIZE
27 #define SIZE 16
28 #endif
29 
30 void
check_rand(int argc,char ** argv)31 check_rand (int argc, char **argv)
32 {
33   mp_size_t size;
34   mp_exp_t exp;
35   int reps = 20000;
36   int i;
37   mpf_t u, v, w, wref;
38   mp_size_t bprec = 100;
39   mpf_t rerr, max_rerr, limit_rerr;
40 
41   if (argc > 1)
42     {
43       reps = strtol (argv[1], 0, 0);
44       if (argc > 2)
45 	bprec = strtol (argv[2], 0, 0);
46     }
47 
48   mpf_set_default_prec (bprec);
49 
50   mpf_init_set_ui (limit_rerr, 1);
51   mpf_div_2exp (limit_rerr, limit_rerr, bprec);
52 #if VERBOSE
53   mpf_dump (limit_rerr);
54 #endif
55   mpf_init (rerr);
56   mpf_init_set_ui (max_rerr, 0);
57 
58   mpf_init (u);
59   mpf_init (v);
60   mpf_init (w);
61   mpf_init (wref);
62   for (i = 0; i < reps; i++)
63     {
64       size = urandom () % (2 * SIZE) - SIZE;
65       exp = urandom () % SIZE;
66       mpf_random2 (u, size, exp);
67 
68       size = urandom () % (2 * SIZE) - SIZE;
69       exp = urandom () % SIZE;
70       mpf_random2 (v, size, exp);
71 
72       if ((urandom () & 1) != 0)
73 	mpf_add_ui (u, v, 1);
74       else if ((urandom () & 1) != 0)
75 	mpf_sub_ui (u, v, 1);
76 
77       mpf_sub (w, u, v);
78       refmpf_sub (wref, u, v);
79 
80       mpf_reldiff (rerr, w, wref);
81       if (mpf_cmp (rerr, max_rerr) > 0)
82 	{
83 	  mpf_set (max_rerr, rerr);
84 #if VERBOSE
85 	  mpf_dump (max_rerr);
86 #endif
87 	  if (mpf_cmp (rerr, limit_rerr) > 0)
88 	    {
89 	      printf ("ERROR after %d tests\n", i);
90 	      printf ("   u = "); mpf_dump (u);
91 	      printf ("   v = "); mpf_dump (v);
92 	      printf ("wref = "); mpf_dump (wref);
93 	      printf ("   w = "); mpf_dump (w);
94 	      abort ();
95 	    }
96 	}
97     }
98 
99   mpf_clear (limit_rerr);
100   mpf_clear (rerr);
101   mpf_clear (max_rerr);
102 
103   mpf_clear (u);
104   mpf_clear (v);
105   mpf_clear (w);
106   mpf_clear (wref);
107 }
108 
109 #define W GMP_NUMB_MAX
110 
111 void
check_data(void)112 check_data (void)
113 {
114   static const struct {
115     struct {
116       int        exp, size;
117       mp_limb_t  d[10];
118     } x, y, want;
119 
120   } data[] = {
121     { { 123, 2, { 8, 9 } },             { 123, 1, { 9 } }, { 122, 1, { 8 } } },
122     { { 1, 1, { 9 } },                  { 1, 1, { 8 } },   { 1, 1, { 1 } } },
123     { { 1, 1, { 9 } },                 { 1, -1, { 6 } },   { 1, 1, { 15 } } },
124     { { 1, 2, { 8, 9 } },               { 1, 1, { 8 } },   { 1, 2, { 8, 1 } } },
125     { { 2, 2, { 8, 1 } },               { 1, 1, { 9 } },   { 1, 1, { W } } },
126     { { 2, 2, { 9, 8 } },               { 1, 1, { 9 } },   { 2, 1, { 8 } } },
127     { { 2, 1, { 1 } },                  { 1, 1, { 1 } },   { 1, 1, { W } } },
128     { { 2, 1, { 9 } },                  { 1, 1, { W } },   { 2, 2, { 1, 8 } } },
129 
130     { { 1, 2, { W, 8 } },             { 1, 1, { 9 } },   { 0, -1, { 1 } } },
131     { { 1, 2, { W, 7 } },             { 1, 1, { 9 } },   { 1, -2, { 1, 1 } } },
132     { { 1, 2, { 1, 8 } },             { 1, 1, { 9 } },   { 0, -1, { W } } },
133     { { 1, 2, { 1, 7 } },             { 1, 1, { 9 } },   { 1, -2, { W, 1 } } },
134     { { 1, 2, { 0, 8 } },             { 1, 1, { 9 } },   { 1, -1, { 1 } } },
135     { { 2, 3, { 5, 8, 1 } },          { 1, 1, { 9 } },   { 1, 2, { 5, W } } },
136     { { 3, 1, { 1 } },                { 1, 1, { 1 } },   { 2, 2, { W, W } } },
137     { { 1, 6, { W, W, W, W, W, 8 } }, { 1, 1, { 9 } },   { -4, -1, { 1 } } },
138     { { 5, 5, { W-6, W, W, W, W } },  { 6, 1, { 1 } },   { 1, -1, { 7 } } },
139 
140     /* f - f == 0, various sizes.
141        These exercise a past problem (gmp 4.1.3 and earlier) where the
142        result exponent was not zeroed on a zero result like this.  */
143     { { 0, 0 }, { 0, 0 }, { 0, 0 } },
144     { { 99, 3, { 0, 0, 1 } },       { 99, 1, { 1 } },             { 0, 0 } },
145     { { 99, 3, { 0, 123, 456 } },   { 99, 2, { 123, 456 } },      { 0, 0 } },
146     { { 99, 3, { 123, 456, 789 } }, { 99, 3, { 123, 456, 789 } }, { 0, 0 } },
147 
148     /* High limbs cancel, leaving just the low limbs of the longer operand.
149        This exercises a past problem (gmp 4.1.3 and earlier) where high zero
150        limbs on the remainder were not stripped before truncating to the
151        destination, causing loss of precision.  */
152     { { 123, 2, { 8, 9 } },             { 123, 1, { 9 } }, { 122, 1, { 8 } } },
153     { { 123, 3, { 8, 0, 9 } },          { 123, 1, { 9 } }, { 121, 1, { 8 } } },
154     { { 123, 4, { 8, 0, 0, 9 } },       { 123, 1, { 9 } }, { 120, 1, { 8 } } },
155     { { 123, 5, { 8, 0, 0, 0, 9 } },    { 123, 1, { 9 } }, { 119, 1, { 8 } } },
156     { { 123, 6, { 8, 0, 0, 0, 0, 9 } }, { 123, 1, { 9 } }, { 118, 1, { 8 } } },
157     /* { { 123, 6, { 8, 0, 0, 0, 0, 9 } }, { 123, 6, { 9, 0, 0, 0, 0, 8 } }, { 122, 5, { W, W, W, W, W } } }, */
158 
159   };
160 
161   mpf_t  x, y, got, want;
162   int  i, swap, fail;
163 
164   fail = 0;
165   mp_trace_base = 16;
166   mpf_init (got);
167 
168   for (i = 0; i < numberof (data); i++)
169     {
170       for (swap = 0; swap <= 7; swap++)
171         {
172           PTR(x) = (mp_ptr) data[i].x.d;
173           SIZ(x) = data[i].x.size;
174           EXP(x) = data[i].x.exp;
175           PREC(x) = numberof (data[i].x.d);
176           MPF_CHECK_FORMAT (x);
177 
178           PTR(y) = (mp_ptr) data[i].y.d;
179           SIZ(y) = data[i].y.size;
180           EXP(y) = data[i].y.exp;
181           PREC(y) = numberof (data[i].y.d);
182           MPF_CHECK_FORMAT (y);
183 
184           PTR(want) = (mp_ptr) data[i].want.d;
185           SIZ(want) = data[i].want.size;
186           EXP(want) = data[i].want.exp;
187           PREC(want) = numberof (data[i].want.d);
188           MPF_CHECK_FORMAT (want);
189 
190           if (swap & 4)
191             {
192               mpf_swap (want, y);
193             }
194 
195 	  if ((SIZ (x) ^ SIZ (y)) < 0)
196 	    continue; /* It's an addition, not a subtraction (TO BE REMOVED) */
197 
198           if (swap & 1)
199             {
200               mpf_swap (x, y);
201               SIZ(want) = - SIZ(want);
202             }
203 
204           if (swap & 2)
205             {
206               SIZ(want) = - SIZ(want);
207               SIZ(x) = - SIZ(x);
208               SIZ(y) = - SIZ(y);
209             }
210 
211           mpf_sub (got, x, y);
212 /*           MPF_CHECK_FORMAT (got); */
213 
214           if (! refmpf_validate ("mpf_sub", got, want))
215             {
216               printf ("check_data() wrong result at data[%d] (operands%s swapped)\n", i, swap ? "" : " not");
217               mpf_trace ("x   ", x);
218               mpf_trace ("y   ", y);
219               mpf_trace ("got ", got);
220               mpf_trace ("want", want);
221 	      fail = 1;
222             }
223 
224 	  if (SIZ (x) == 1 || SIZ (x) == 0 )
225 	    {
226 	      if (SIZ (y)) EXP (y) -= EXP (x) - (mp_exp_t) SIZ (x);
227 	      if (SIZ (want)) EXP (want) -= EXP (x) - (mp_exp_t) SIZ (x);
228 	      EXP (x) = (mp_exp_t) SIZ (x);
229 
230 	      if (mpf_fits_uint_p (x))
231 		{
232 		  mpf_ui_sub (got, mpf_get_ui (x), y);
233 
234 		  if (! refmpf_validate ("mpf_ui_sub", got, want))
235 		    {
236 		      printf ("check_data() wrong result at data[%d] (operands%s swapped)\n", i, swap ? "" : " not");
237 		      mpf_trace ("x   ", x);
238 		      mpf_trace ("y   ", y);
239 		      mpf_trace ("got ", got);
240 		      mpf_trace ("want", want);
241 		      fail = 1;
242 		    }
243 		}
244 	    }
245 
246 	  if (SIZ (y) == 1 || SIZ (y) == 0)
247 	    {
248 	      if (SIZ (x)) EXP (x) -= EXP (y) - (mp_exp_t) SIZ (y);
249 	      if (SIZ (want)) EXP (want) -= EXP (y) - (mp_exp_t) SIZ (y);
250 	      EXP (y) = (mp_exp_t) SIZ (y);
251 
252 	      if (mpf_fits_uint_p (x))
253 		{
254 		  mpf_sub_ui (got, x, mpf_get_ui (y));
255 
256 		  if (! refmpf_validate ("mpf_sub_ui", got, want))
257 		    {
258 		      printf ("check_data() wrong result at data[%d] (operands%s swapped)\n", i, swap ? "" : " not");
259 		      mpf_trace ("x   ", x);
260 		      mpf_trace ("y   ", y);
261 		      mpf_trace ("got ", got);
262 		      mpf_trace ("want", want);
263 		      fail = 1;
264 		    }
265 		}
266 	    }
267 
268         }
269     }
270 
271   mpf_clear (got);
272   if (fail)
273     abort ();
274 }
275 
276 
277 int
main(int argc,char ** argv)278 main (int argc, char **argv)
279 {
280   tests_start ();
281 
282   check_data ();
283   check_rand (argc, argv);
284 
285   tests_end ();
286   exit (0);
287 }
288