1 /* mpz_congruent_p -- test congruence of two mpz's.
2
3 Copyright 2001, 2002, 2005 Free Software Foundation, Inc.
4
5 This file is part of the GNU MP Library.
6
7 The GNU MP Library is free software; you can redistribute it and/or modify
8 it under the terms of the GNU Lesser General Public License as published by
9 the Free Software Foundation; either version 3 of the License, or (at your
10 option) any later version.
11
12 The GNU MP Library is distributed in the hope that it will be useful, but
13 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
14 or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
15 License for more details.
16
17 You should have received a copy of the GNU Lesser General Public License
18 along with the GNU MP Library. If not, see http://www.gnu.org/licenses/. */
19
20 #include "gmp.h"
21 #include "gmp-impl.h"
22 #include "longlong.h"
23
24
25 /* For big divisors this code is only very slightly better than the user
26 doing a combination of mpz_sub and mpz_tdiv_r, but it's quite convenient,
27 and perhaps in the future can be improved, in similar ways to
28 mpn_divisible_p perhaps.
29
30 The csize==1 / dsize==1 special case makes mpz_congruent_p as good as
31 mpz_congruent_ui_p on relevant operands, though such a combination
32 probably doesn't occur often.
33
34 Alternatives:
35
36 If c<d then it'd work to just form a%d and compare a and c (either as
37 a==c or a+c==d depending on the signs), but the saving from avoiding the
38 abs(a-c) calculation would be small compared to the division.
39
40 Similarly if both a<d and c<d then it would work to just compare a and c
41 (a==c or a+c==d), but this isn't considered a particularly important case
42 and so isn't done for the moment.
43
44 Low zero limbs on d could be stripped and the corresponding limbs of a
45 and c tested and skipped, but doing so would introduce a borrow when a
46 and c differ in sign and have non-zero skipped limbs. It doesn't seem
47 worth the complications to do this, since low zero limbs on d should
48 occur only rarely. */
49
50 int
mpz_congruent_p(mpz_srcptr a,mpz_srcptr c,mpz_srcptr d)51 mpz_congruent_p (mpz_srcptr a, mpz_srcptr c, mpz_srcptr d)
52 {
53 mp_size_t asize, csize, dsize, sign;
54 mp_srcptr ap, cp, dp;
55 mp_ptr xp;
56 mp_limb_t alow, clow, dlow, dmask, r;
57 int result;
58 TMP_DECL;
59
60 dsize = SIZ(d);
61 if (UNLIKELY (dsize == 0))
62 return (mpz_cmp (a, c) == 0);
63
64 dsize = ABS(dsize);
65 dp = PTR(d);
66
67 if (ABSIZ(a) < ABSIZ(c))
68 MPZ_SRCPTR_SWAP (a, c);
69
70 asize = SIZ(a);
71 csize = SIZ(c);
72 sign = (asize ^ csize);
73
74 asize = ABS(asize);
75 ap = PTR(a);
76
77 if (csize == 0)
78 return mpn_divisible_p (ap, asize, dp, dsize);
79
80 csize = ABS(csize);
81 cp = PTR(c);
82
83 alow = ap[0];
84 clow = cp[0];
85 dlow = dp[0];
86
87 /* Check a==c mod low zero bits of dlow. This might catch a few cases of
88 a!=c quickly, and it helps the csize==1 special cases below. */
89 dmask = LOW_ZEROS_MASK (dlow) & GMP_NUMB_MASK;
90 alow = (sign >= 0 ? alow : -alow);
91 if (((alow-clow) & dmask) != 0)
92 return 0;
93
94 if (csize == 1)
95 {
96 if (dsize == 1)
97 {
98 cong_1:
99 if (sign < 0)
100 NEG_MOD (clow, clow, dlow);
101
102 if (ABOVE_THRESHOLD (asize, BMOD_1_TO_MOD_1_THRESHOLD))
103 {
104 r = mpn_mod_1 (ap, asize, dlow);
105 if (clow < dlow)
106 return r == clow;
107 else
108 return r == (clow % dlow);
109 }
110
111 if ((dlow & 1) == 0)
112 {
113 /* Strip low zero bits to get odd d required by modexact. If
114 d==e*2^n then a==c mod d if and only if both a==c mod e and
115 a==c mod 2^n, the latter having been done above. */
116 unsigned twos;
117 count_trailing_zeros (twos, dlow);
118 dlow >>= twos;
119 }
120
121 r = mpn_modexact_1c_odd (ap, asize, dlow, clow);
122 return r == 0 || r == dlow;
123 }
124
125 /* dlow==0 is avoided since we don't want to bother handling extra low
126 zero bits if dsecond is even (would involve borrow if a,c differ in
127 sign and alow,clow!=0). */
128 if (dsize == 2 && dlow != 0)
129 {
130 mp_limb_t dsecond = dp[1];
131
132 if (dsecond <= dmask)
133 {
134 unsigned twos;
135 count_trailing_zeros (twos, dlow);
136 dlow = (dlow >> twos) | (dsecond << (GMP_NUMB_BITS-twos));
137 ASSERT_LIMB (dlow);
138
139 /* dlow will be odd here, so the test for it even under cong_1
140 is unnecessary, but the rest of that code is wanted. */
141 goto cong_1;
142 }
143 }
144 }
145
146 TMP_MARK;
147 xp = TMP_ALLOC_LIMBS (asize+1);
148
149 /* calculate abs(a-c) */
150 if (sign >= 0)
151 {
152 /* same signs, subtract */
153 if (asize > csize || mpn_cmp (ap, cp, asize) >= 0)
154 ASSERT_NOCARRY (mpn_sub (xp, ap, asize, cp, csize));
155 else
156 ASSERT_NOCARRY (mpn_sub_n (xp, cp, ap, asize));
157 MPN_NORMALIZE (xp, asize);
158 }
159 else
160 {
161 /* different signs, add */
162 mp_limb_t carry;
163 carry = mpn_add (xp, ap, asize, cp, csize);
164 xp[asize] = carry;
165 asize += (carry != 0);
166 }
167
168 result = mpn_divisible_p (xp, asize, dp, dsize);
169
170 TMP_FREE;
171 return result;
172 }
173