1 /* mpz_congruent_2exp_p -- test congruence of mpz mod 2^n. 2 3 Copyright 2001, 2002 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 23 24 int 25 mpz_congruent_2exp_p (mpz_srcptr a, mpz_srcptr c, mp_bitcnt_t d) __GMP_NOTHROW 26 { 27 mp_size_t i, dlimbs; 28 unsigned dbits; 29 mp_ptr ap, cp; 30 mp_limb_t dmask, alimb, climb, sum; 31 mp_size_t asize_signed, csize_signed, asize, csize; 32 33 if (ABSIZ(a) < ABSIZ(c)) 34 MPZ_SRCPTR_SWAP (a, c); 35 36 dlimbs = d / GMP_NUMB_BITS; 37 dbits = d % GMP_NUMB_BITS; 38 dmask = (CNST_LIMB(1) << dbits) - 1; 39 40 ap = PTR(a); 41 cp = PTR(c); 42 43 asize_signed = SIZ(a); 44 asize = ABS(asize_signed); 45 46 csize_signed = SIZ(c); 47 csize = ABS(csize_signed); 48 49 if (csize_signed == 0) 50 goto a_zeros; 51 52 if ((asize_signed ^ csize_signed) >= 0) 53 { 54 /* same signs, direct comparison */ 55 56 /* a==c for limbs in common */ 57 if (mpn_cmp (ap, cp, MIN (csize, dlimbs)) != 0) 58 return 0; 59 60 /* if that's all of dlimbs, then a==c for remaining bits */ 61 if (csize > dlimbs) 62 return ((ap[dlimbs]-cp[dlimbs]) & dmask) == 0; 63 64 a_zeros: 65 /* a remains, need all zero bits */ 66 67 /* if d covers all of a and c, then must be exactly equal */ 68 if (asize <= dlimbs) 69 return asize == csize; 70 71 /* whole limbs zero */ 72 for (i = csize; i < dlimbs; i++) 73 if (ap[i] != 0) 74 return 0; 75 76 /* partial limb zero */ 77 return (ap[dlimbs] & dmask) == 0; 78 } 79 else 80 { 81 /* different signs, negated comparison */ 82 83 /* common low zero limbs, stopping at first non-zeros, which must 84 match twos complement */ 85 i = 0; 86 for (;;) 87 { 88 ASSERT (i < csize); /* always have a non-zero limb on c */ 89 alimb = ap[i]; 90 climb = cp[i]; 91 sum = (alimb + climb) & GMP_NUMB_MASK; 92 93 if (i >= dlimbs) 94 return (sum & dmask) == 0; 95 i++; 96 97 /* require both zero, or first non-zeros as twos-complements */ 98 if (sum != 0) 99 return 0; 100 101 if (alimb != 0) 102 break; 103 } 104 105 /* further limbs matching as ones-complement */ 106 for (;;) 107 { 108 if (i >= csize) 109 break; 110 111 alimb = ap[i]; 112 climb = cp[i]; 113 sum = (alimb + climb + 1) & GMP_NUMB_MASK; 114 115 if (i >= dlimbs) 116 return (sum & dmask) == 0; 117 118 if (sum != 0) 119 return 0; 120 121 i++; 122 } 123 124 /* no more c, so require all 1 bits in a */ 125 126 if (asize < dlimbs) 127 return 0; /* not enough a */ 128 129 /* whole limbs */ 130 for ( ; i < dlimbs; i++) 131 if (ap[i] != GMP_NUMB_MAX) 132 return 0; 133 134 /* if only whole limbs, no further fetches from a */ 135 if (dbits == 0) 136 return 1; 137 138 /* need enough a */ 139 if (asize == dlimbs) 140 return 0; 141 142 return ((ap[dlimbs]+1) & dmask) == 0; 143 } 144 } 145