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 } 295