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