Annotation of OpenXM_contrib/gmp/mpn/generic/divrem_2.c, Revision 1.1.1.2
1.1 maekawa 1: /* mpn_divrem_2 -- Divide natural numbers, producing both remainder and
2: quotient. The divisor is two limbs.
3:
4: THIS FILE CONTAINS INTERNAL FUNCTIONS WITH MUTABLE INTERFACES. IT IS
5: ONLY SAFE TO REACH THEM THROUGH DOCUMENTED INTERFACES. IN FACT, IT IS
6: ALMOST GUARANTEED THAT THEY'LL CHANGE OR DISAPPEAR IN A FUTURE GNU MP
7: RELEASE.
8:
9:
1.1.1.2 ! ohara 10: Copyright 1993, 1994, 1995, 1996, 1999, 2000, 2001, 2002 Free Software
! 11: Foundation, Inc.
1.1 maekawa 12:
13: This file is part of the GNU MP Library.
14:
15: The GNU MP Library is free software; you can redistribute it and/or modify
16: it under the terms of the GNU Lesser General Public License as published by
17: the Free Software Foundation; either version 2.1 of the License, or (at your
18: option) any later version.
19:
20: The GNU MP Library is distributed in the hope that it will be useful, but
21: WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
22: or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
23: License for more details.
24:
25: You should have received a copy of the GNU Lesser General Public License
26: along with the GNU MP Library; see the file COPYING.LIB. If not, write to
27: the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
28: MA 02111-1307, USA. */
29:
30: #include "gmp.h"
31: #include "gmp-impl.h"
32: #include "longlong.h"
33:
1.1.1.2 ! ohara 34:
! 35: /* The size where udiv_qrnnd_preinv should be used rather than udiv_qrnnd,
! 36: meaning the quotient size where that should happen, the quotient size
! 37: being how many udiv divisions will be done.
! 38:
! 39: The default is to use preinv always, CPUs where this doesn't suit have
! 40: tuned thresholds. Note in particular that preinv should certainly be
! 41: used if that's the only division available (USE_PREINV_ALWAYS). */
! 42:
! 43: #ifndef DIVREM_2_THRESHOLD
! 44: #define DIVREM_2_THRESHOLD 0
! 45: #endif
! 46:
! 47:
1.1 maekawa 48: /* Divide num (NP/NSIZE) by den (DP/2) and write
49: the NSIZE-2 least significant quotient limbs at QP
50: and the 2 long remainder at NP. If QEXTRA_LIMBS is
51: non-zero, generate that many fraction bits and append them after the
52: other quotient limbs.
53: Return the most significant limb of the quotient, this is always 0 or 1.
54:
55: Preconditions:
56: 0. NSIZE >= 2.
57: 1. The most significant bit of the divisor must be set.
58: 2. QP must either not overlap with the input operands at all, or
59: QP + 2 >= NP must hold true. (This means that it's
60: possible to put the quotient in the high part of NUM, right after the
61: remainder in NUM.
62: 3. NSIZE >= 2, even if QEXTRA_LIMBS is non-zero. */
63:
64: mp_limb_t
65: mpn_divrem_2 (mp_ptr qp, mp_size_t qxn,
1.1.1.2 ! ohara 66: mp_ptr np, mp_size_t nn,
1.1 maekawa 67: mp_srcptr dp)
68: {
69: mp_limb_t most_significant_q_limb = 0;
70: mp_size_t i;
71: mp_limb_t n1, n0, n2;
72: mp_limb_t d1, d0;
73: mp_limb_t d1inv;
1.1.1.2 ! ohara 74: int use_preinv;
! 75:
! 76: ASSERT (nn >= 2);
! 77: ASSERT (qxn >= 0);
! 78: ASSERT (dp[1] & GMP_NUMB_HIGHBIT);
! 79: ASSERT (! MPN_OVERLAP_P (qp, nn-2+qxn, np, nn) || qp+2 >= np);
! 80: ASSERT_MPN (np, nn);
! 81: ASSERT_MPN (dp, 2);
1.1 maekawa 82:
1.1.1.2 ! ohara 83: np += nn - 2;
1.1 maekawa 84: d1 = dp[1];
85: d0 = dp[0];
86: n1 = np[1];
87: n0 = np[0];
88:
89: if (n1 >= d1 && (n1 > d1 || n0 >= d0))
90: {
1.1.1.2 ! ohara 91: #if GMP_NAIL_BITS == 0
1.1 maekawa 92: sub_ddmmss (n1, n0, n1, n0, d1, d0);
1.1.1.2 ! ohara 93: #else
! 94: n0 = n0 - d0;
! 95: n1 = n1 - d1 - (n0 >> GMP_LIMB_BITS - 1);
! 96: n0 &= GMP_NUMB_MASK;
! 97: #endif
1.1 maekawa 98: most_significant_q_limb = 1;
99: }
100:
1.1.1.2 ! ohara 101: use_preinv = ABOVE_THRESHOLD (qxn + nn - 2, DIVREM_2_THRESHOLD);
! 102: if (use_preinv)
! 103: invert_limb (d1inv, d1);
1.1 maekawa 104:
1.1.1.2 ! ohara 105: for (i = qxn + nn - 2 - 1; i >= 0; i--)
1.1 maekawa 106: {
107: mp_limb_t q;
108: mp_limb_t r;
109:
110: if (i >= qxn)
111: np--;
112: else
113: np[0] = 0;
114:
115: if (n1 == d1)
116: {
1.1.1.2 ! ohara 117: /* Q should be either 111..111 or 111..110. Need special handling
1.1 maekawa 118: of this rare case as normal division would give overflow. */
1.1.1.2 ! ohara 119: q = GMP_NUMB_MASK;
1.1 maekawa 120:
1.1.1.2 ! ohara 121: r = (n0 + d1) & GMP_NUMB_MASK;
1.1 maekawa 122: if (r < d1) /* Carry in the addition? */
123: {
1.1.1.2 ! ohara 124: #if GMP_NAIL_BITS == 0
1.1 maekawa 125: add_ssaaaa (n1, n0, r - d0, np[0], 0, d0);
1.1.1.2 ! ohara 126: #else
! 127: n0 = np[0] + d0;
! 128: n1 = (r - d0 + (n0 >> GMP_NUMB_BITS)) & GMP_NUMB_MASK;
! 129: n0 &= GMP_NUMB_MASK;
! 130: #endif
1.1 maekawa 131: qp[i] = q;
132: continue;
133: }
134: n1 = d0 - (d0 != 0);
1.1.1.2 ! ohara 135: n0 = -d0 & GMP_NUMB_MASK;
1.1 maekawa 136: }
137: else
138: {
1.1.1.2 ! ohara 139: if (use_preinv)
1.1 maekawa 140: udiv_qrnnd_preinv (q, r, n1, n0, d1, d1inv);
141: else
1.1.1.2 ! ohara 142: udiv_qrnnd (q, r, n1, n0 << GMP_NAIL_BITS, d1 << GMP_NAIL_BITS);
! 143: r >>= GMP_NAIL_BITS;
! 144: umul_ppmm (n1, n0, d0, q << GMP_NAIL_BITS);
! 145: n0 >>= GMP_NAIL_BITS;
1.1 maekawa 146: }
147:
148: n2 = np[0];
149:
150: q_test:
151: if (n1 > r || (n1 == r && n0 > n2))
152: {
153: /* The estimated Q was too large. */
154: q--;
155:
1.1.1.2 ! ohara 156: #if GMP_NAIL_BITS == 0
1.1 maekawa 157: sub_ddmmss (n1, n0, n1, n0, 0, d0);
1.1.1.2 ! ohara 158: #else
! 159: n0 = n0 - d0;
! 160: n1 = n1 - (n0 >> GMP_LIMB_BITS - 1);
! 161: n0 &= GMP_NUMB_MASK;
! 162: #endif
1.1 maekawa 163: r += d1;
164: if (r >= d1) /* If not carry, test Q again. */
165: goto q_test;
166: }
167:
168: qp[i] = q;
1.1.1.2 ! ohara 169: #if GMP_NAIL_BITS == 0
1.1 maekawa 170: sub_ddmmss (n1, n0, r, n2, n1, n0);
1.1.1.2 ! ohara 171: #else
! 172: n0 = n2 - n0;
! 173: n1 = r - n1 - (n0 >> GMP_LIMB_BITS - 1);
! 174: n0 &= GMP_NUMB_MASK;
! 175: #endif
1.1 maekawa 176: }
177: np[1] = n1;
178: np[0] = n0;
179:
180: return most_significant_q_limb;
181: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to divrem_2.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / gmp / mpn / generic