1 /* Copyright 2006, 2007, 2009, 2010, 2017 Free Software Foundation, Inc. 2 3 This file is part of the GNU MP Library test suite. 4 5 The GNU MP Library test suite is free software; you can redistribute it 6 and/or modify it under the terms of the GNU General Public License as 7 published by the Free Software Foundation; either version 3 of the License, 8 or (at your option) any later version. 9 10 The GNU MP Library test suite is distributed in the hope that it will be 11 useful, but WITHOUT ANY WARRANTY; without even the implied warranty of 12 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General 13 Public License for more details. 14 15 You should have received a copy of the GNU General Public License along with 16 the GNU MP Library test suite. If not, see https://www.gnu.org/licenses/. */ 17 18 19 #include <stdlib.h> /* for strtol */ 20 #include <stdio.h> /* for printf */ 21 22 #include "gmp-impl.h" 23 #include "longlong.h" 24 #include "tests/tests.h" 25 26 27 static void 28 dumpy (mp_srcptr p, mp_size_t n) 29 { 30 mp_size_t i; 31 if (n > 20) 32 { 33 for (i = n - 1; i >= n - 4; i--) 34 { 35 printf ("%0*lx", (int) (2 * sizeof (mp_limb_t)), p[i]); 36 printf (" "); 37 } 38 printf ("... "); 39 for (i = 3; i >= 0; i--) 40 { 41 printf ("%0*lx", (int) (2 * sizeof (mp_limb_t)), p[i]); 42 printf (i == 0 ? "" : " "); 43 } 44 } 45 else 46 { 47 for (i = n - 1; i >= 0; i--) 48 { 49 printf ("%0*lx", (int) (2 * sizeof (mp_limb_t)), p[i]); 50 printf (i == 0 ? "" : " "); 51 } 52 } 53 puts (""); 54 } 55 56 static unsigned long test; 57 58 void 59 check_one (mp_ptr qp, mp_srcptr rp, mp_limb_t rh, 60 mp_srcptr np, mp_size_t nn, mp_srcptr dp, mp_size_t dn, const char *fname) 61 { 62 mp_size_t qn; 63 mp_ptr tp; 64 mp_limb_t cy = 4711; /* silence warnings */ 65 TMP_DECL; 66 67 qn = nn - dn; 68 69 if (qn == 0) 70 return; 71 72 TMP_MARK; 73 74 tp = TMP_ALLOC_LIMBS (nn + 1); 75 76 if (dn >= qn) 77 mpn_mul (tp, dp, dn, qp, qn); 78 else 79 mpn_mul (tp, qp, qn, dp, dn); 80 81 cy = mpn_add_n (tp, tp, np, nn); 82 83 if (! mpn_zero_p (tp, qn) 84 || (rp != NULL && (cy != rh || mpn_cmp (tp + qn, rp, dn) != 0))) 85 { 86 printf ("\r*******************************************************************************\n"); 87 printf ("%s inconsistent in test %lu\n", fname, test); 88 printf ("N= "); dumpy (np, nn); 89 printf ("D= "); dumpy (dp, dn); 90 printf ("Q= "); dumpy (qp, qn); 91 if (rp != NULL) 92 { 93 printf ("R= "); dumpy (rp, dn); 94 printf ("Rb= %d, Cy=%d\n", (int) cy, (int) rh); 95 } 96 printf ("T= "); dumpy (tp, nn); 97 printf ("nn = %ld, dn = %ld, qn = %ld", nn, dn, qn); 98 printf ("\n*******************************************************************************\n"); 99 abort (); 100 } 101 102 TMP_FREE; 103 } 104 105 106 /* These are *bit* sizes. */ 107 #define SIZE_LOG 16 108 #define MAX_DN (1L << SIZE_LOG) 109 #define MAX_NN (1L << (SIZE_LOG + 1)) 110 111 #define COUNT 500 112 113 mp_limb_t 114 random_word (gmp_randstate_ptr rs) 115 { 116 mpz_t x; 117 mp_limb_t r; 118 TMP_DECL; 119 TMP_MARK; 120 121 MPZ_TMP_INIT (x, 2); 122 mpz_urandomb (x, rs, 32); 123 r = mpz_get_ui (x); 124 TMP_FREE; 125 return r; 126 } 127 128 int 129 main (int argc, char **argv) 130 { 131 gmp_randstate_ptr rands; 132 unsigned long maxnbits, maxdbits, nbits, dbits; 133 mpz_t n, d, tz; 134 mp_size_t maxnn, maxdn, nn, dn, clearn, i; 135 mp_ptr np, dp, qp, rp; 136 mp_limb_t rh; 137 mp_limb_t t; 138 mp_limb_t dinv; 139 int count = COUNT; 140 mp_ptr scratch; 141 mp_limb_t ran; 142 mp_size_t alloc, itch; 143 mp_limb_t rran0, rran1, qran0, qran1; 144 TMP_DECL; 145 146 TESTS_REPS (count, argv, argc); 147 148 maxdbits = MAX_DN; 149 maxnbits = MAX_NN; 150 151 tests_start (); 152 rands = RANDS; 153 154 mpz_init (n); 155 mpz_init (d); 156 mpz_init (tz); 157 158 maxnn = maxnbits / GMP_NUMB_BITS + 1; 159 maxdn = maxdbits / GMP_NUMB_BITS + 1; 160 161 TMP_MARK; 162 163 qp = TMP_ALLOC_LIMBS (maxnn + 2) + 1; 164 rp = TMP_ALLOC_LIMBS (maxnn + 2) + 1; 165 166 alloc = 1; 167 scratch = __GMP_ALLOCATE_FUNC_LIMBS (alloc); 168 169 for (test = 0; test < count;) 170 { 171 nbits = random_word (rands) % (maxnbits - GMP_NUMB_BITS) + 2 * GMP_NUMB_BITS; 172 if (maxdbits > nbits) 173 dbits = random_word (rands) % nbits + 1; 174 else 175 dbits = random_word (rands) % maxdbits + 1; 176 177 #if RAND_UNIFORM 178 #define RANDFUNC mpz_urandomb 179 #else 180 #define RANDFUNC mpz_rrandomb 181 #endif 182 183 do 184 { 185 RANDFUNC (n, rands, nbits); 186 do 187 { 188 RANDFUNC (d, rands, dbits); 189 } 190 while (mpz_sgn (d) == 0); 191 192 np = PTR (n); 193 dp = PTR (d); 194 nn = SIZ (n); 195 dn = SIZ (d); 196 } 197 while (nn < dn); 198 199 dp[0] |= 1; 200 201 mpz_urandomb (tz, rands, 32); 202 t = mpz_get_ui (tz); 203 204 if (t % 17 == 0) 205 dp[0] = GMP_NUMB_MAX; 206 207 switch ((int) t % 16) 208 { 209 case 0: 210 clearn = random_word (rands) % nn; 211 for (i = 0; i <= clearn; i++) 212 np[i] = 0; 213 break; 214 case 1: 215 mpn_sub_1 (np + nn - dn, dp, dn, random_word (rands)); 216 break; 217 case 2: 218 mpn_add_1 (np + nn - dn, dp, dn, random_word (rands)); 219 break; 220 } 221 222 test++; 223 224 binvert_limb (dinv, dp[0]); 225 226 rran0 = random_word (rands); 227 rran1 = random_word (rands); 228 qran0 = random_word (rands); 229 qran1 = random_word (rands); 230 231 qp[-1] = qran0; 232 qp[nn - dn + 1] = qran1; 233 rp[-1] = rran0; 234 235 ran = random_word (rands); 236 237 if ((double) (nn - dn) * dn < 1e5) 238 { 239 if (nn > dn) 240 { 241 /* Test mpn_sbpi1_bdiv_qr */ 242 MPN_ZERO (qp, nn - dn); 243 MPN_ZERO (rp, dn); 244 MPN_COPY (rp, np, nn); 245 rh = mpn_sbpi1_bdiv_qr (qp, rp, nn, dp, dn, -dinv); 246 ASSERT_ALWAYS (qp[-1] == qran0); ASSERT_ALWAYS (qp[nn - dn + 1] == qran1); 247 ASSERT_ALWAYS (rp[-1] == rran0); 248 check_one (qp, rp + nn - dn, rh, np, nn, dp, dn, "mpn_sbpi1_bdiv_qr"); 249 250 /* Test mpn_sbpi1_bdiv_q */ 251 MPN_COPY (rp, np, nn); 252 MPN_ZERO (qp, nn - dn); 253 mpn_sbpi1_bdiv_q (qp, rp, nn - dn, dp, MIN(dn,nn-dn), -dinv); 254 ASSERT_ALWAYS (qp[-1] == qran0); ASSERT_ALWAYS (qp[nn - dn + 1] == qran1); 255 ASSERT_ALWAYS (rp[-1] == rran0); 256 check_one (qp, NULL, 0, np, nn, dp, dn, "mpn_sbpi1_bdiv_q"); 257 258 /* Test mpn_sbpi1_bdiv_r; we use mpn_sbpi1_bdiv_q's quotient. */ 259 MPN_COPY (rp, np, nn); 260 mpn_sbpi1_bdiv_r (rp, nn, dp, dn, -dinv); 261 ASSERT_ALWAYS (rp[-1] == rran0); 262 check_one (qp, NULL, 0, np, nn, dp, dn, "mpn_sbpi1_bdiv_r"); 263 } 264 } 265 266 if (dn >= 4 && nn - dn >= 2) 267 { 268 /* Test mpn_dcpi1_bdiv_qr */ 269 MPN_COPY (rp, np, nn); 270 MPN_ZERO (qp, nn - dn); 271 rh = mpn_dcpi1_bdiv_qr (qp, rp, nn, dp, dn, -dinv); 272 ASSERT_ALWAYS (qp[-1] == qran0); ASSERT_ALWAYS (qp[nn - dn + 1] == qran1); 273 ASSERT_ALWAYS (rp[-1] == rran0); 274 check_one (qp, rp + nn - dn, rh, np, nn, dp, dn, "mpn_dcpi1_bdiv_qr"); 275 } 276 277 if (dn >= 4 && nn - dn >= 2) 278 { 279 /* Test mpn_dcpi1_bdiv_q */ 280 MPN_COPY (rp, np, nn); 281 MPN_ZERO (qp, nn - dn); 282 mpn_dcpi1_bdiv_q (qp, rp, nn - dn, dp, MIN(dn,nn-dn), -dinv); 283 ASSERT_ALWAYS (qp[-1] == qran0); ASSERT_ALWAYS (qp[nn - dn + 1] == qran1); 284 ASSERT_ALWAYS (rp[-1] == rran0); 285 check_one (qp, NULL, 0, np, nn, dp, dn, "mpn_dcpi1_bdiv_q"); 286 } 287 288 if (nn > dn) 289 { 290 /* Test mpn_bdiv_qr */ 291 itch = mpn_bdiv_qr_itch (nn, dn); 292 if (itch + 1 > alloc) 293 { 294 scratch = __GMP_REALLOCATE_FUNC_LIMBS (scratch, alloc, itch + 1); 295 alloc = itch + 1; 296 } 297 scratch[itch] = ran; 298 MPN_ZERO (qp, nn - dn); 299 MPN_ZERO (rp, dn); 300 rp[dn] = rran1; 301 rh = mpn_bdiv_qr (qp, rp, np, nn, dp, dn, scratch); 302 ASSERT_ALWAYS (ran == scratch[itch]); 303 ASSERT_ALWAYS (qp[-1] == qran0); ASSERT_ALWAYS (qp[nn - dn + 1] == qran1); 304 ASSERT_ALWAYS (rp[-1] == rran0); ASSERT_ALWAYS (rp[dn] == rran1); 305 306 check_one (qp, rp, rh, np, nn, dp, dn, "mpn_bdiv_qr"); 307 } 308 309 if (nn - dn < 2 || dn < 2) 310 continue; 311 312 /* Test mpn_mu_bdiv_qr */ 313 itch = mpn_mu_bdiv_qr_itch (nn, dn); 314 if (itch + 1 > alloc) 315 { 316 scratch = __GMP_REALLOCATE_FUNC_LIMBS (scratch, alloc, itch + 1); 317 alloc = itch + 1; 318 } 319 scratch[itch] = ran; 320 MPN_ZERO (qp, nn - dn); 321 MPN_ZERO (rp, dn); 322 rp[dn] = rran1; 323 rh = mpn_mu_bdiv_qr (qp, rp, np, nn, dp, dn, scratch); 324 ASSERT_ALWAYS (ran == scratch[itch]); 325 ASSERT_ALWAYS (qp[-1] == qran0); ASSERT_ALWAYS (qp[nn - dn + 1] == qran1); 326 ASSERT_ALWAYS (rp[-1] == rran0); ASSERT_ALWAYS (rp[dn] == rran1); 327 check_one (qp, rp, rh, np, nn, dp, dn, "mpn_mu_bdiv_qr"); 328 329 /* Test mpn_mu_bdiv_q */ 330 itch = mpn_mu_bdiv_q_itch (nn, dn); 331 if (itch + 1 > alloc) 332 { 333 scratch = __GMP_REALLOCATE_FUNC_LIMBS (scratch, alloc, itch + 1); 334 alloc = itch + 1; 335 } 336 scratch[itch] = ran; 337 MPN_ZERO (qp, nn - dn + 1); 338 mpn_mu_bdiv_q (qp, np, nn - dn, dp, dn, scratch); 339 ASSERT_ALWAYS (ran == scratch[itch]); 340 ASSERT_ALWAYS (qp[-1] == qran0); ASSERT_ALWAYS (qp[nn - dn + 1] == qran1); 341 check_one (qp, NULL, 0, np, nn, dp, dn, "mpn_mu_bdiv_q"); 342 } 343 344 __GMP_FREE_FUNC_LIMBS (scratch, alloc); 345 346 TMP_FREE; 347 348 mpz_clear (n); 349 mpz_clear (d); 350 mpz_clear (tz); 351 352 tests_end (); 353 return 0; 354 } 355