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 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