mpz_rootrem (mpz_t x, mpz_t r, const mpz_t y, unsigned long z)
{
int sgn;
+ mp_bitcnt_t bc;
mpz_t t, u;
sgn = y->_mp_size < 0;
mpz_init (u);
mpz_init (t);
- mpz_setbit (t, mpz_sizeinbase (y, 2) / z + 1);
+ bc = (mpz_sizeinbase (y, 2) - 1) / z + 1;
+ mpz_setbit (t, bc);
if (z == 2) /* simplify sqrt loop: z-1 == 1 */
do {
mpz_init (V);
/* n-(D/n) = n+1 = d*2^{b0}, with d = (n>>b0) | 1 */
- b0 = mpz_scan0 (n, 0);
+ b0 = mpn_common_scan (~ n->_mp_d[0], 0, n->_mp_d, n->_mp_size, GMP_LIMB_MAX);
+ /* b0 = mpz_scan0 (n, 0); */
/* D= P^2 - 4Q; P = 1; Q = (1-D)/4 */
Q = (D & 2) ? (long) (D >> 2) + 1 : -(long) (D >> 2);
mpz_powm_ui (y, y, 2, n);
if (mpz_cmp (y, nm1) == 0)
return 1;
- /* y == 1 means that the previous y was a non-trivial square root
- of 1 (mod n). y == 0 means that n is a power of the base.
- In either case, n is not prime. */
- if (mpz_cmp_ui (y, 1) <= 0)
- return 0;
}
return 0;
}