1/* Test mpz_[cft]div_[qr]_2exp.
2
3Copyright 2001 Free Software Foundation, Inc.
4
5This file is part of the GNU MP Library.
6
7The GNU MP Library is free software; you can redistribute it and/or modify
8it under the terms of the GNU Lesser General Public License as published by
9the Free Software Foundation; either version 3 of the License, or (at your
10option) any later version.
11
12The GNU MP Library is distributed in the hope that it will be useful, but
13WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
14or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
15License for more details.
16
17You should have received a copy of the GNU Lesser General Public License
18along with the GNU MP Library.  If not, see http://www.gnu.org/licenses/.  */
19
20#include <stdio.h>
21#include <stdlib.h>
22
23#include "gmp.h"
24#include "gmp-impl.h"
25#include "tests.h"
26
27
28/* If the remainder is in the correct range and q*d+r is correct, then q
29   must have rounded correctly.  */
30
31void
32check_one (mpz_srcptr a, unsigned long d)
33{
34  mpz_t  q, r, p, d2exp;
35  int    inplace;
36
37  mpz_init (d2exp);
38  mpz_init (q);
39  mpz_init (r);
40  mpz_init (p);
41
42  mpz_set_ui (d2exp, 1L);
43  mpz_mul_2exp (d2exp, d2exp, d);
44
45#define INPLACE(fun,dst,src,d)  \
46  if (inplace)                  \
47    {                           \
48      mpz_set (dst, src);       \
49      fun (dst, dst, d);        \
50    }                           \
51  else                          \
52    fun (dst, src, d);
53
54  for (inplace = 0; inplace <= 1; inplace++)
55    {
56      INPLACE (mpz_fdiv_q_2exp, q, a, d);
57      INPLACE (mpz_fdiv_r_2exp, r, a, d);
58
59      mpz_mul_2exp (p, q, d);
60      mpz_add (p, p, r);
61      if (mpz_sgn (r) < 0 || mpz_cmp (r, d2exp) >= 0)
62	{
63	  printf ("mpz_fdiv_r_2exp result out of range\n");
64	  goto error;
65	}
66      if (mpz_cmp (p, a) != 0)
67	{
68	  printf ("mpz_fdiv_[qr]_2exp doesn't multiply back\n");
69	  goto error;
70	}
71
72
73      INPLACE (mpz_cdiv_q_2exp, q, a, d);
74      INPLACE (mpz_cdiv_r_2exp, r, a, d);
75
76      mpz_mul_2exp (p, q, d);
77      mpz_add (p, p, r);
78      if (mpz_sgn (r) > 0 || mpz_cmpabs (r, d2exp) >= 0)
79	{
80	  printf ("mpz_cdiv_r_2exp result out of range\n");
81	  goto error;
82	}
83      if (mpz_cmp (p, a) != 0)
84	{
85	  printf ("mpz_cdiv_[qr]_2exp doesn't multiply back\n");
86	  goto error;
87	}
88
89
90      INPLACE (mpz_tdiv_q_2exp, q, a, d);
91      INPLACE (mpz_tdiv_r_2exp, r, a, d);
92
93      mpz_mul_2exp (p, q, d);
94      mpz_add (p, p, r);
95      if (mpz_sgn (r) != 0 && mpz_sgn (r) != mpz_sgn (a))
96	{
97	  printf ("mpz_tdiv_r_2exp result wrong sign\n");
98	  goto error;
99	}
100      if (mpz_cmpabs (r, d2exp) >= 0)
101	{
102	  printf ("mpz_tdiv_r_2exp result out of range\n");
103	  goto error;
104	}
105      if (mpz_cmp (p, a) != 0)
106	{
107	  printf ("mpz_tdiv_[qr]_2exp doesn't multiply back\n");
108	  goto error;
109	}
110    }
111
112  mpz_clear (d2exp);
113  mpz_clear (q);
114  mpz_clear (r);
115  mpz_clear (p);
116  return;
117
118
119 error:
120  mpz_trace ("a", a);
121  printf    ("d=%lu\n", d);
122  mpz_trace ("q", q);
123  mpz_trace ("r", r);
124  mpz_trace ("p", p);
125
126  mp_trace_base = -16;
127  mpz_trace ("a", a);
128  printf    ("d=0x%lX\n", d);
129  mpz_trace ("q", q);
130  mpz_trace ("r", r);
131  mpz_trace ("p", p);
132
133  abort ();
134}
135
136
137void
138check_all (mpz_ptr a, unsigned long d)
139{
140  check_one (a, d);
141  mpz_neg (a, a);
142  check_one (a, d);
143}
144
145
146void
147check_various (void)
148{
149  static const unsigned long  table[] = {
150    0, 1, 2, 3, 4, 5,
151    GMP_NUMB_BITS-1, GMP_NUMB_BITS, GMP_NUMB_BITS+1,
152    2*GMP_NUMB_BITS-1, 2*GMP_NUMB_BITS, 2*GMP_NUMB_BITS+1,
153    3*GMP_NUMB_BITS-1, 3*GMP_NUMB_BITS, 3*GMP_NUMB_BITS+1,
154    4*GMP_NUMB_BITS-1, 4*GMP_NUMB_BITS, 4*GMP_NUMB_BITS+1
155  };
156
157  int            i, j;
158  unsigned long  n, d;
159  mpz_t          a;
160
161  mpz_init (a);
162
163  /* a==0, and various d */
164  mpz_set_ui (a, 0L);
165  for (i = 0; i < numberof (table); i++)
166    check_one (a, table[i]);
167
168  /* a==2^n, and various d */
169  for (i = 0; i < numberof (table); i++)
170    {
171      n = table[i];
172      mpz_set_ui (a, 1L);
173      mpz_mul_2exp (a, a, n);
174
175      for (j = 0; j < numberof (table); j++)
176	{
177	  d = table[j];
178	  check_all (a, d);
179	}
180    }
181
182  mpz_clear (a);
183}
184
185
186void
187check_random (int argc, char *argv[])
188{
189  gmp_randstate_ptr  rands = RANDS;
190  int            reps = 100;
191  mpz_t          a;
192  unsigned long  d;
193  int            i;
194
195  if (argc == 2)
196    reps = atoi (argv[1]);
197
198  mpz_init (a);
199
200  for (i = 0; i < reps; i++)
201    {
202      /* exponentially within 2 to 257 bits */
203      mpz_erandomb (a, rands, urandom () % 8 + 2);
204
205      d = urandom () % 256;
206
207      check_all (a, d);
208    }
209
210  mpz_clear (a);
211}
212
213
214int
215main (int argc, char *argv[])
216{
217  tests_start ();
218
219  check_various ();
220  check_random (argc, argv);
221
222  tests_end ();
223  exit (0);
224}
225