1 /* Test mpq_get_d and mpq_set_d
2
3 Copyright 1991, 1993, 1994, 1996, 2000-2003, 2012, 2013 Free Software
4 Foundation, Inc.
5
6 This file is part of the GNU MP Library test suite.
7
8 The GNU MP Library test suite is free software; you can redistribute it
9 and/or modify it under the terms of the GNU General Public License as
10 published by the Free Software Foundation; either version 3 of the License,
11 or (at your option) any later version.
12
13 The GNU MP Library test suite is distributed in the hope that it will be
14 useful, but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
16 Public License for more details.
17
18 You should have received a copy of the GNU General Public License along with
19 the GNU MP Library test suite. If not, see https://www.gnu.org/licenses/. */
20
21 #include <stdio.h>
22 #include <stdlib.h>
23
24 #include "gmp-impl.h"
25 #include "tests.h"
26
27 #ifndef SIZE
28 #define SIZE 8
29 #endif
30
31 /* VAX D floats only have an 8 bit signed exponent, so anything 2^128 or
32 bigger will overflow, that being 4 limbs. */
33 #if defined (__vax) || defined (__vax__) && SIZE > 4
34 #undef SIZE
35 #define SIZE 4
36 #define EPSIZE 3
37 #else
38 #define EPSIZE SIZE
39 #endif
40
41 void dump (mpq_t);
42
43 void
44 check_monotonic (int argc, char **argv)
45 {
46 mpq_t a;
47 mp_size_t size;
48 int reps = 100;
49 int i, j;
50 double last_d, new_d;
51 mpq_t qlast_d, qnew_d;
52 mpq_t eps;
53
54 if (argc == 2)
55 reps = atoi (argv[1]);
56
57 /* The idea here is to test the monotonousness of mpq_get_d by adding
58 numbers to the numerator and denominator. */
59
60 mpq_init (a);
61 mpq_init (eps);
62 mpq_init (qlast_d);
63 mpq_init (qnew_d);
64
65 for (i = 0; i < reps; i++)
66 {
67 size = urandom () % SIZE - SIZE/2;
68 mpz_random2 (mpq_numref (a), size);
69 do
70 {
71 size = urandom () % SIZE - SIZE/2;
72 mpz_random2 (mpq_denref (a), size);
73 }
74 while (mpz_cmp_ui (mpq_denref (a), 0) == 0);
75
76 mpq_canonicalize (a);
77
78 last_d = mpq_get_d (a);
79 mpq_set_d (qlast_d, last_d);
80 for (j = 0; j < 10; j++)
81 {
82 size = urandom () % EPSIZE + 1;
83 mpz_random2 (mpq_numref (eps), size);
84 size = urandom () % EPSIZE + 1;
85 mpz_random2 (mpq_denref (eps), size);
86 mpq_canonicalize (eps);
87
88 mpq_add (a, a, eps);
89 mpq_canonicalize (a);
90 new_d = mpq_get_d (a);
91 if (last_d > new_d)
92 {
93 printf ("\nERROR (test %d/%d): bad mpq_get_d results\n", i, j);
94 printf ("last: %.16g\n", last_d);
95 printf (" new: %.16g\n", new_d); dump (a);
96 abort ();
97 }
98 mpq_set_d (qnew_d, new_d);
99 MPQ_CHECK_FORMAT (qnew_d);
100 if (mpq_cmp (qlast_d, qnew_d) > 0)
101 {
102 printf ("ERROR (test %d/%d): bad mpq_set_d results\n", i, j);
103 printf ("last: %.16g\n", last_d); dump (qlast_d);
104 printf (" new: %.16g\n", new_d); dump (qnew_d);
105 abort ();
106 }
107 last_d = new_d;
108 mpq_set (qlast_d, qnew_d);
109 }
110 }
111
112 mpq_clear (a);
113 mpq_clear (eps);
114 mpq_clear (qlast_d);
115 mpq_clear (qnew_d);
116 }
117
118 double
119 my_ldexp (double d, int e)
120 {
121 for (;;)
122 {
123 if (e > 0)
124 {
125 if (e >= 16)
126 {
127 d *= 65536.0;
128 e -= 16;
129 }
130 else
131 {
132 d *= 2.0;
133 e -= 1;
134 }
135 }
136 else if (e < 0)
137 {
138
139 if (e <= -16)
140 {
141 d /= 65536.0;
142 e += 16;
143 }
144 else
145 {
146 d /= 2.0;
147 e += 1;
148 }
149 }
150 else
151 return d;
152 }
153 }
154
155 #define MAXEXP 500
156
157 #if defined (__vax) || defined (__vax__)
158 #undef MAXEXP
159 #define MAXEXP 30
160 #endif
161
162 void
163 check_random (int argc, char **argv)
164 {
165 gmp_randstate_ptr rands = RANDS;
166
167 double d;
168 mpq_t q;
169 mpz_t a, t;
170 int exp;
171
172 int test, reps = 100000;
173
174 if (argc == 2)
175 reps = 100 * atoi (argv[1]);
176
177 mpq_init (q);
178 mpz_init (a);
179 mpz_init (t);
180
181 for (test = 0; test < reps; test++)
182 {
183 mpz_rrandomb (a, rands, 53);
184 mpz_urandomb (t, rands, 32);
185 exp = mpz_get_ui (t) % (2*MAXEXP) - MAXEXP;
186
187 d = my_ldexp (mpz_get_d (a), exp);
188 mpq_set_d (q, d);
189 /* Check that n/d = a * 2^exp, or
190 d*a 2^{exp} = n */
191 mpz_mul (t, a, mpq_denref (q));
192 if (exp > 0)
193 mpz_mul_2exp (t, t, exp);
194 else
195 {
196 if (!mpz_divisible_2exp_p (t, -exp))
197 goto fail;
198 mpz_div_2exp (t, t, -exp);
199 }
200 if (mpz_cmp (t, mpq_numref (q)) != 0)
201 {
202 fail:
203 printf ("ERROR (check_random test %d): bad mpq_set_d results\n", test);
204 printf ("%.16g\n", d);
205 gmp_printf ("%Qd\n", q);
206 abort ();
207 }
208 }
209 mpq_clear (q);
210 mpz_clear (t);
211 mpz_clear (a);
212 }
213
214 void
215 dump (mpq_t x)
216 {
217 mpz_out_str (stdout, 10, mpq_numref (x));
218 printf ("/");
219 mpz_out_str (stdout, 10, mpq_denref (x));
220 printf ("\n");
221 }
222
223 /* Check various values 2^n and 1/2^n. */
224 void
225 check_onebit (void)
226 {
227 static const long data[] = {
228 -3*GMP_NUMB_BITS-1, -3*GMP_NUMB_BITS, -3*GMP_NUMB_BITS+1,
229 -2*GMP_NUMB_BITS-1, -2*GMP_NUMB_BITS, -2*GMP_NUMB_BITS+1,
230 -GMP_NUMB_BITS-1, -GMP_NUMB_BITS, -GMP_NUMB_BITS+1,
231 -5, -2, -1, 0, 1, 2, 5,
232 GMP_NUMB_BITS-1, GMP_NUMB_BITS, GMP_NUMB_BITS+1,
233 2*GMP_NUMB_BITS-1, 2*GMP_NUMB_BITS, 2*GMP_NUMB_BITS+1,
234 3*GMP_NUMB_BITS-1, 3*GMP_NUMB_BITS, 3*GMP_NUMB_BITS+1,
235 };
236
237 int i, neg;
238 long exp, l;
239 mpq_t q;
240 double got, want;
241
242 mpq_init (q);
243
244 for (i = 0; i < numberof (data); i++)
245 {
246 exp = data[i];
247
248 mpq_set_ui (q, 1L, 1L);
249 if (exp >= 0)
250 mpq_mul_2exp (q, q, exp);
251 else
252 mpq_div_2exp (q, q, -exp);
253
254 want = 1.0;
255 for (l = 0; l < exp; l++)
256 want *= 2.0;
257 for (l = 0; l > exp; l--)
258 want /= 2.0;
259
260 for (neg = 0; neg <= 1; neg++)
261 {
262 if (neg)
263 {
264 mpq_neg (q, q);
265 want = -want;
266 }
267
268 got = mpq_get_d (q);
269
270 if (got != want)
271 {
272 printf ("mpq_get_d wrong on %s2**%ld\n", neg ? "-" : "", exp);
273 mpq_trace (" q ", q);
274 d_trace (" want ", want);
275 d_trace (" got ", got);
276 abort();
277 }
278 }
279 }
280 mpq_clear (q);
281 }
282
283 int
284 main (int argc, char **argv)
285 {
286 tests_start ();
287
288 check_onebit ();
289 check_monotonic (argc, argv);
290 check_random (argc, argv);
291
292 tests_end ();
293 exit (0);
294 }