186d7f5d3SJohn Marino /* mpz_congruent_ui_p -- test congruence of mpz and ulong.
286d7f5d3SJohn Marino
386d7f5d3SJohn Marino Copyright 2000, 2001, 2002 Free Software Foundation, Inc.
486d7f5d3SJohn Marino
586d7f5d3SJohn Marino This file is part of the GNU MP Library.
686d7f5d3SJohn Marino
786d7f5d3SJohn Marino The GNU MP Library is free software; you can redistribute it and/or modify
886d7f5d3SJohn Marino it under the terms of the GNU Lesser General Public License as published by
986d7f5d3SJohn Marino the Free Software Foundation; either version 3 of the License, or (at your
1086d7f5d3SJohn Marino option) any later version.
1186d7f5d3SJohn Marino
1286d7f5d3SJohn Marino The GNU MP Library is distributed in the hope that it will be useful, but
1386d7f5d3SJohn Marino WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
1486d7f5d3SJohn Marino or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
1586d7f5d3SJohn Marino License for more details.
1686d7f5d3SJohn Marino
1786d7f5d3SJohn Marino You should have received a copy of the GNU Lesser General Public License
1886d7f5d3SJohn Marino along with the GNU MP Library. If not, see http://www.gnu.org/licenses/. */
1986d7f5d3SJohn Marino
2086d7f5d3SJohn Marino #include "gmp.h"
2186d7f5d3SJohn Marino #include "gmp-impl.h"
2286d7f5d3SJohn Marino #include "longlong.h"
2386d7f5d3SJohn Marino
2486d7f5d3SJohn Marino
2586d7f5d3SJohn Marino /* There's some explicit checks for c<d since it seems reasonably likely an
2686d7f5d3SJohn Marino application might use that in a test.
2786d7f5d3SJohn Marino
2886d7f5d3SJohn Marino Hopefully the compiler can generate something good for r==(c%d), though
2986d7f5d3SJohn Marino if modexact is being used exclusively then that's not reached. */
3086d7f5d3SJohn Marino
3186d7f5d3SJohn Marino int
mpz_congruent_ui_p(mpz_srcptr a,unsigned long cu,unsigned long du)3286d7f5d3SJohn Marino mpz_congruent_ui_p (mpz_srcptr a, unsigned long cu, unsigned long du)
3386d7f5d3SJohn Marino {
3486d7f5d3SJohn Marino mp_srcptr ap;
3586d7f5d3SJohn Marino mp_size_t asize;
3686d7f5d3SJohn Marino mp_limb_t c, d, r;
3786d7f5d3SJohn Marino
3886d7f5d3SJohn Marino if (UNLIKELY (du == 0))
3986d7f5d3SJohn Marino return (mpz_cmp_ui (a, cu) == 0);
4086d7f5d3SJohn Marino
4186d7f5d3SJohn Marino asize = SIZ(a);
4286d7f5d3SJohn Marino if (asize == 0)
4386d7f5d3SJohn Marino {
4486d7f5d3SJohn Marino if (cu < du)
4586d7f5d3SJohn Marino return cu == 0;
4686d7f5d3SJohn Marino else
4786d7f5d3SJohn Marino return (cu % du) == 0;
4886d7f5d3SJohn Marino }
4986d7f5d3SJohn Marino
5086d7f5d3SJohn Marino /* For nails don't try to be clever if c or d is bigger than a limb, just
5186d7f5d3SJohn Marino fake up some mpz_t's and go to the main mpz_congruent_p. */
5286d7f5d3SJohn Marino if (du > GMP_NUMB_MAX || cu > GMP_NUMB_MAX)
5386d7f5d3SJohn Marino {
5486d7f5d3SJohn Marino mp_limb_t climbs[2], dlimbs[2];
5586d7f5d3SJohn Marino mpz_t cz, dz;
5686d7f5d3SJohn Marino
5786d7f5d3SJohn Marino ALLOC(cz) = 2;
5886d7f5d3SJohn Marino PTR(cz) = climbs;
5986d7f5d3SJohn Marino ALLOC(dz) = 2;
6086d7f5d3SJohn Marino PTR(dz) = dlimbs;
6186d7f5d3SJohn Marino
6286d7f5d3SJohn Marino mpz_set_ui (cz, cu);
6386d7f5d3SJohn Marino mpz_set_ui (dz, du);
6486d7f5d3SJohn Marino return mpz_congruent_p (a, cz, dz);
6586d7f5d3SJohn Marino }
6686d7f5d3SJohn Marino
6786d7f5d3SJohn Marino /* NEG_MOD works on limbs, so convert ulong to limb */
6886d7f5d3SJohn Marino c = cu;
6986d7f5d3SJohn Marino d = du;
7086d7f5d3SJohn Marino
7186d7f5d3SJohn Marino if (asize < 0)
7286d7f5d3SJohn Marino {
7386d7f5d3SJohn Marino asize = -asize;
7486d7f5d3SJohn Marino NEG_MOD (c, c, d);
7586d7f5d3SJohn Marino }
7686d7f5d3SJohn Marino
7786d7f5d3SJohn Marino ap = PTR (a);
7886d7f5d3SJohn Marino
7986d7f5d3SJohn Marino if (ABOVE_THRESHOLD (asize, BMOD_1_TO_MOD_1_THRESHOLD))
8086d7f5d3SJohn Marino {
8186d7f5d3SJohn Marino r = mpn_mod_1 (ap, asize, d);
8286d7f5d3SJohn Marino if (c < d)
8386d7f5d3SJohn Marino return r == c;
8486d7f5d3SJohn Marino else
8586d7f5d3SJohn Marino return r == (c % d);
8686d7f5d3SJohn Marino }
8786d7f5d3SJohn Marino
8886d7f5d3SJohn Marino if ((d & 1) == 0)
8986d7f5d3SJohn Marino {
9086d7f5d3SJohn Marino /* Strip low zero bits to get odd d required by modexact. If
9186d7f5d3SJohn Marino d==e*2^n then a==c mod d if and only if both a==c mod 2^n
9286d7f5d3SJohn Marino and a==c mod e. */
9386d7f5d3SJohn Marino
9486d7f5d3SJohn Marino unsigned twos;
9586d7f5d3SJohn Marino
9686d7f5d3SJohn Marino if ((ap[0]-c) & LOW_ZEROS_MASK (d))
9786d7f5d3SJohn Marino return 0;
9886d7f5d3SJohn Marino
9986d7f5d3SJohn Marino count_trailing_zeros (twos, d);
10086d7f5d3SJohn Marino d >>= twos;
10186d7f5d3SJohn Marino }
10286d7f5d3SJohn Marino
10386d7f5d3SJohn Marino r = mpn_modexact_1c_odd (ap, asize, d, c);
10486d7f5d3SJohn Marino return r == 0 || r == d;
10586d7f5d3SJohn Marino }
106