Annotation of OpenXM_contrib/gmp/mpn/generic/mod_1.c, Revision 1.1.1.3
1.1 maekawa 1: /* mpn_mod_1(dividend_ptr, dividend_size, divisor_limb) --
2: Divide (DIVIDEND_PTR,,DIVIDEND_SIZE) by DIVISOR_LIMB.
3: Return the single-limb remainder.
4: There are no constraints on the value of the divisor.
5:
1.1.1.3 ! ohara 6: Copyright 1991, 1993, 1994, 1999, 2000, 2002 Free Software Foundation, Inc.
1.1 maekawa 7:
8: This file is part of the GNU MP Library.
9:
10: The GNU MP Library is free software; you can redistribute it and/or modify
1.1.1.2 maekawa 11: it under the terms of the GNU Lesser General Public License as published by
12: the Free Software Foundation; either version 2.1 of the License, or (at your
1.1 maekawa 13: option) any later version.
14:
15: The GNU MP Library is distributed in the hope that it will be useful, but
16: WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
1.1.1.2 maekawa 17: or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
1.1 maekawa 18: License for more details.
19:
1.1.1.2 maekawa 20: You should have received a copy of the GNU Lesser General Public License
1.1 maekawa 21: along with the GNU MP Library; see the file COPYING.LIB. If not, write to
22: the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
23: MA 02111-1307, USA. */
24:
25: #include "gmp.h"
26: #include "gmp-impl.h"
27: #include "longlong.h"
28:
29:
1.1.1.3 ! ohara 30: /* The size where udiv_qrnnd_preinv should be used rather than udiv_qrnnd,
! 31: meaning the quotient size where that should happen, the quotient size
! 32: being how many udiv divisions will be done.
! 33:
! 34: The default is to use preinv always, CPUs where this doesn't suit have
! 35: tuned thresholds. Note in particular that preinv should certainly be
! 36: used if that's the only division available (USE_PREINV_ALWAYS). */
1.1 maekawa 37:
1.1.1.3 ! ohara 38: #ifndef MOD_1_NORM_THRESHOLD
! 39: #define MOD_1_NORM_THRESHOLD 0
! 40: #endif
! 41: #ifndef MOD_1_UNNORM_THRESHOLD
! 42: #define MOD_1_UNNORM_THRESHOLD 0
1.1 maekawa 43: #endif
44:
45:
1.1.1.3 ! ohara 46: /* The comments in mpn/generic/divrem_1.c apply here too.
1.1 maekawa 47:
1.1.1.3 ! ohara 48: As noted in the algorithms section of the manual, the shifts in the loop
! 49: for the unnorm case can be avoided by calculating r = a%(d*2^n), followed
! 50: by a final (r*2^n)%(d*2^n). In fact if it happens that a%(d*2^n) can
! 51: skip a division where (a*2^n)%(d*2^n) can't then there's the same number
! 52: of divide steps, though how often that happens depends on the assumed
! 53: distributions of dividend and divisor. In any case this idea is left to
! 54: CPU specific implementations to consider. */
1.1 maekawa 55:
1.1.1.3 ! ohara 56: mp_limb_t
! 57: mpn_mod_1 (mp_srcptr up, mp_size_t un, mp_limb_t d)
! 58: {
! 59: mp_size_t i;
! 60: mp_limb_t n1, n0, r;
! 61: mp_limb_t dummy;
! 62:
! 63: ASSERT (un >= 0);
! 64: ASSERT (d != 0);
! 65:
! 66: /* Botch: Should this be handled at all? Rely on callers?
! 67: But note un==0 is currently required by mpz/fdiv_r_ui.c and possibly
! 68: other places. */
! 69: if (un == 0)
! 70: return 0;
1.1 maekawa 71:
1.1.1.3 ! ohara 72: d <<= GMP_NAIL_BITS;
1.1 maekawa 73:
1.1.1.3 ! ohara 74: if ((d & GMP_LIMB_HIGHBIT) != 0)
! 75: {
! 76: /* High limb is initial remainder, possibly with one subtract of
! 77: d to get r<d. */
! 78: r = up[un - 1] << GMP_NAIL_BITS;
! 79: if (r >= d)
! 80: r -= d;
! 81: r >>= GMP_NAIL_BITS;
! 82: un--;
! 83: if (un == 0)
! 84: return r;
1.1 maekawa 85:
1.1.1.3 ! ohara 86: if (BELOW_THRESHOLD (un, MOD_1_NORM_THRESHOLD))
! 87: {
! 88: plain:
! 89: for (i = un - 1; i >= 0; i--)
1.1 maekawa 90: {
1.1.1.3 ! ohara 91: n0 = up[i] << GMP_NAIL_BITS;
! 92: udiv_qrnnd (dummy, r, r, n0, d);
! 93: r >>= GMP_NAIL_BITS;
1.1 maekawa 94: }
1.1.1.3 ! ohara 95: return r;
1.1 maekawa 96: }
97: else
98: {
1.1.1.3 ! ohara 99: mp_limb_t inv;
! 100: invert_limb (inv, d);
! 101: for (i = un - 1; i >= 0; i--)
1.1 maekawa 102: {
1.1.1.3 ! ohara 103: n0 = up[i] << GMP_NAIL_BITS;
! 104: udiv_qrnnd_preinv (dummy, r, r, n0, d, inv);
! 105: r >>= GMP_NAIL_BITS;
1.1 maekawa 106: }
107: return r;
108: }
109: }
110: else
111: {
1.1.1.3 ! ohara 112: int norm;
! 113:
! 114: /* Skip a division if high < divisor. Having the test here before
! 115: normalizing will still skip as often as possible. */
! 116: r = up[un - 1] << GMP_NAIL_BITS;
! 117: if (r < d)
1.1 maekawa 118: {
1.1.1.3 ! ohara 119: r >>= GMP_NAIL_BITS;
! 120: un--;
! 121: if (un == 0)
! 122: return r;
! 123: }
! 124: else
! 125: r = 0;
1.1 maekawa 126:
1.1.1.3 ! ohara 127: /* If udiv_qrnnd doesn't need a normalized divisor, can use the simple
! 128: code above. */
! 129: if (! UDIV_NEEDS_NORMALIZATION
! 130: && BELOW_THRESHOLD (un, MOD_1_UNNORM_THRESHOLD))
! 131: goto plain;
! 132:
! 133: count_leading_zeros (norm, d);
! 134: d <<= norm;
1.1 maekawa 135:
1.1.1.3 ! ohara 136: n1 = up[un - 1] << GMP_NAIL_BITS;
! 137: r = (r << norm) | (n1 >> (GMP_LIMB_BITS - norm));
1.1 maekawa 138:
1.1.1.3 ! ohara 139: if (UDIV_NEEDS_NORMALIZATION
! 140: && BELOW_THRESHOLD (un, MOD_1_UNNORM_THRESHOLD))
! 141: {
! 142: for (i = un - 2; i >= 0; i--)
! 143: {
! 144: n0 = up[i] << GMP_NAIL_BITS;
1.1 maekawa 145: udiv_qrnnd (dummy, r, r,
1.1.1.3 ! ohara 146: (n1 << norm) | (n0 >> (GMP_NUMB_BITS - norm)),
! 147: d);
! 148: r >>= GMP_NAIL_BITS;
! 149: n1 = n0;
1.1 maekawa 150: }
1.1.1.3 ! ohara 151: udiv_qrnnd (dummy, r, r, n1 << norm, d);
! 152: r >>= GMP_NAIL_BITS;
! 153: return r >> norm;
1.1 maekawa 154: }
155: else
156: {
1.1.1.3 ! ohara 157: mp_limb_t inv;
! 158: invert_limb (inv, d);
! 159:
! 160: for (i = un - 2; i >= 0; i--)
! 161: {
! 162: n0 = up[i] << GMP_NAIL_BITS;
! 163: udiv_qrnnd_preinv (dummy, r, r,
! 164: (n1 << norm) | (n0 >> (GMP_NUMB_BITS - norm)),
! 165: d, inv);
! 166: r >>= GMP_NAIL_BITS;
! 167: n1 = n0;
! 168: }
! 169: udiv_qrnnd_preinv (dummy, r, r, n1 << norm, d, inv);
! 170: r >>= GMP_NAIL_BITS;
! 171: return r >> norm;
1.1 maekawa 172: }
173: }
174: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>