1 /* mpn_divisible_p -- mpn by mpn divisibility test
3 THE FUNCTIONS IN THIS FILE ARE FOR INTERNAL USE ONLY. THEY'RE ALMOST
4 CERTAIN TO BE SUBJECT TO INCOMPATIBLE CHANGES OR DISAPPEAR COMPLETELY IN
5 FUTURE GNU MP RELEASES.
7 Copyright 2001, 2002, 2005, 2009 Free Software Foundation, Inc.
9 This file is part of the GNU MP Library.
11 The GNU MP Library is free software; you can redistribute it and/or modify
12 it under the terms of the GNU Lesser General Public License as published by
13 the Free Software Foundation; either version 3 of the License, or (at your
14 option) any later version.
16 The GNU MP Library is distributed in the hope that it will be useful, but
17 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
18 or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
19 License for more details.
21 You should have received a copy of the GNU Lesser General Public License
22 along with the GNU MP Library. If not, see http://www.gnu.org/licenses/. */
29 /* Determine whether {ap,an} is divisible by {dp,dn}. Must have both
30 operands normalized, meaning high limbs non-zero, except that an==0 is
33 There usually won't be many low zero bits on d, but the checks for this
34 are fast and might pick up a few operand combinations, in particular they
35 might reduce d to fit the single-limb mod_1/modexact_1 code.
39 Getting the remainder limb by limb would make an early exit possible on
40 finding a non-zero. This would probably have to be bdivmod style so
41 there's no addback, but it would need a multi-precision inverse and so
42 might be slower than the plain method (on small sizes at least).
44 When d must be normalized (shifted to high bit set), it's possible to
45 just append a low zero limb to "a" rather than bit-shifting as
46 mpn_tdiv_qr does internally, so long as it's already been checked that a
47 has at least as many trailing zeros bits as d. Or equivalently, pass
48 qxn==1 to mpn_tdiv_qr, if/when it accepts that. */
51 mpn_divisible_p (mp_srcptr ap, mp_size_t an,
52 mp_srcptr dp, mp_size_t dn)
54 mp_limb_t alow, dlow, dmask;
62 ASSERT (an == 0 || ap[an-1] != 0);
64 ASSERT (dp[dn-1] != 0);
68 /* When a<d only a==0 is divisible.
69 Notice this test covers all cases of an==0. */
73 /* Strip low zero limbs from d, requiring a==0 on those. */
83 return 0; /* a has fewer low zero limbs than d, so not divisible */
85 /* a!=0 and d!=0 so won't get to n==0 */
86 an--; ASSERT (an >= 1);
87 dn--; ASSERT (dn >= 1);
92 /* a must have at least as many low zero bits as d */
93 dmask = LOW_ZEROS_MASK (dlow);
94 if ((alow & dmask) != 0)
99 if (ABOVE_THRESHOLD (an, BMOD_1_TO_MOD_1_THRESHOLD))
100 return mpn_mod_1 (ap, an, dlow) == 0;
102 count_trailing_zeros (twos, dlow);
104 return mpn_modexact_1_odd (ap, an, dlow) == 0;
109 mp_limb_t dsecond = dp[1];
110 if (dsecond <= dmask)
112 count_trailing_zeros (twos, dlow);
113 dlow = (dlow >> twos) | (dsecond << (GMP_NUMB_BITS-twos));
115 return MPN_MOD_OR_MODEXACT_1_ODD (ap, an, dlow) == 0;
119 /* Should we compute Q = A * D^(-1) mod B^k,
120 R = A - Q * D mod B^k
121 here, for some small values of k? Then check if R = 0 (mod B^k). */
123 /* We could also compute A' = A mod T and D' = D mod P, for some
124 P = 3 * 5 * 7 * 11 ..., and then check if any prime factor from P
125 dividing D' also divides A'. */
129 rp = TMP_ALLOC_LIMBS (an + 1);
130 qp = TMP_ALLOC_LIMBS (an - dn + 1); /* FIXME: Could we avoid this */
132 count_trailing_zeros (twos, dp[0]);
136 tp = TMP_ALLOC_LIMBS (dn);
137 ASSERT_NOCARRY (mpn_rshift (tp, dp, dn, twos));
140 ASSERT_NOCARRY (mpn_rshift (rp, ap, an, twos));
144 MPN_COPY (rp, ap, an);
146 if (rp[an - 1] >= dp[dn - 1])
157 ASSERT (an > dn); /* requirement of functions below */
159 if (BELOW_THRESHOLD (dn, DC_BDIV_QR_THRESHOLD) ||
160 BELOW_THRESHOLD (an - dn, DC_BDIV_QR_THRESHOLD))
162 binvert_limb (di, dp[0]);
163 mpn_sbpi1_bdiv_qr (qp, rp, an, dp, dn, -di);
166 else if (BELOW_THRESHOLD (dn, MU_BDIV_QR_THRESHOLD))
168 binvert_limb (di, dp[0]);
169 mpn_dcpi1_bdiv_qr (qp, rp, an, dp, dn, -di);
174 tp = TMP_ALLOC_LIMBS (mpn_mu_bdiv_qr_itch (an, dn));
175 mpn_mu_bdiv_qr (qp, rp, rp, an, dp, dn, tp);
178 /* test for {rp,dn} zero or non-zero */