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