Annotation of OpenXM_contrib/gmp/mpz/divexact.c, Revision 1.1.1.3
1.1 maekawa 1: /* mpz_divexact -- finds quotient when known that quot * den == num && den != 0.
2:
1.1.1.3 ! ohara 3: Copyright 1991, 1993, 1994, 1995, 1996, 1997, 1998, 2000, 2001, 2002 Free
! 4: Software Foundation, Inc.
1.1 maekawa 5:
6: This file is part of the GNU MP Library.
7:
8: The GNU MP Library is free software; you can redistribute it and/or modify
1.1.1.2 maekawa 9: it under the terms of the GNU Lesser General Public License as published by
10: the Free Software Foundation; either version 2.1 of the License, or (at your
1.1 maekawa 11: option) any later version.
12:
13: The GNU MP Library is distributed in the hope that it will be useful, but
14: WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
1.1.1.2 maekawa 15: or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
1.1 maekawa 16: License for more details.
17:
1.1.1.2 maekawa 18: You should have received a copy of the GNU Lesser General Public License
1.1 maekawa 19: along with the GNU MP Library; see the file COPYING.LIB. If not, write to
20: the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
21: MA 02111-1307, USA. */
22:
23: /* Ken Weber (kweber@mat.ufrgs.br, kweber@mcs.kent.edu)
24:
25: Funding for this work has been partially provided by Conselho Nacional
26: de Desenvolvimento Cienti'fico e Tecnolo'gico (CNPq) do Brazil, Grant
27: 301314194-2, and was done while I was a visiting reseacher in the Instituto
28: de Matema'tica at Universidade Federal do Rio Grande do Sul (UFRGS).
29:
30: References:
1.1.1.2 maekawa 31: T. Jebelean, An algorithm for exact division, Journal of Symbolic
32: Computation, v. 15, 1993, pp. 169-180. */
1.1 maekawa 33:
34: #include "gmp.h"
35: #include "gmp-impl.h"
36: #include "longlong.h"
37:
38: void
39: mpz_divexact (mpz_ptr quot, mpz_srcptr num, mpz_srcptr den)
40: {
41: mp_ptr qp, tp;
42: mp_size_t qsize, tsize;
1.1.1.2 maekawa 43: mp_srcptr np, dp;
44: mp_size_t nsize, dsize;
1.1 maekawa 45: TMP_DECL (marker);
46:
1.1.1.2 maekawa 47: nsize = ABS (num->_mp_size);
48: dsize = ABS (den->_mp_size);
49:
50: qsize = nsize - dsize + 1;
51: if (quot->_mp_alloc < qsize)
52: _mpz_realloc (quot, qsize);
53:
54: np = num->_mp_d;
55: dp = den->_mp_d;
56: qp = quot->_mp_d;
1.1 maekawa 57:
58: if (nsize == 0)
59: {
1.1.1.2 maekawa 60: if (dsize == 0)
61: DIVIDE_BY_ZERO;
1.1 maekawa 62: quot->_mp_size = 0;
63: return;
64: }
65:
1.1.1.2 maekawa 66: if (dsize <= 1)
67: {
68: if (dsize == 1)
69: {
1.1.1.3 ! ohara 70: MPN_DIVREM_OR_DIVEXACT_1 (qp, np, nsize, dp[0]);
1.1.1.2 maekawa 71: qsize -= qp[qsize - 1] == 0;
72: quot->_mp_size = (num->_mp_size ^ den->_mp_size) >= 0 ? qsize : -qsize;
73: return;
74: }
75:
76: /* Generate divide-by-zero error since dsize == 0. */
77: DIVIDE_BY_ZERO;
78: }
1.1 maekawa 79:
80: TMP_MARK (marker);
81:
82: /* QUOT <-- NUM/2^r, T <-- DEN/2^r where = r number of twos in DEN. */
83: while (dp[0] == 0)
84: np += 1, nsize -= 1, dp += 1, dsize -= 1;
85: tsize = MIN (qsize, dsize);
1.1.1.2 maekawa 86: if ((dp[0] & 1) != 0)
1.1 maekawa 87: {
1.1.1.2 maekawa 88: if (quot == den) /* QUOT and DEN overlap. */
1.1 maekawa 89: {
1.1.1.2 maekawa 90: tp = (mp_ptr) TMP_ALLOC (tsize * BYTES_PER_MP_LIMB);
1.1 maekawa 91: MPN_COPY (tp, dp, tsize);
92: }
93: else
94: tp = (mp_ptr) dp;
1.1.1.2 maekawa 95: if (qp != np)
96: MPN_COPY_INCR (qp, np, qsize);
1.1 maekawa 97: }
98: else
99: {
1.1.1.2 maekawa 100: unsigned int r;
101: tp = (mp_ptr) TMP_ALLOC (tsize * BYTES_PER_MP_LIMB);
1.1 maekawa 102: count_trailing_zeros (r, dp[0]);
103: mpn_rshift (tp, dp, tsize, r);
104: if (dsize > tsize)
1.1.1.3 ! ohara 105: tp[tsize - 1] |= (dp[tsize] << (GMP_NUMB_BITS - r)) & GMP_NUMB_MASK;
1.1 maekawa 106: mpn_rshift (qp, np, qsize, r);
107: if (nsize > qsize)
1.1.1.3 ! ohara 108: qp[qsize - 1] |= (np[qsize] << (GMP_NUMB_BITS - r)) & GMP_NUMB_MASK;
1.1 maekawa 109: }
110:
111: /* Now QUOT <-- QUOT/T. */
1.1.1.3 ! ohara 112: mpn_bdivmod (qp, qp, qsize, tp, tsize, qsize * GMP_NUMB_BITS);
1.1 maekawa 113: MPN_NORMALIZE (qp, qsize);
114:
1.1.1.2 maekawa 115: quot->_mp_size = (num->_mp_size ^ den->_mp_size) >= 0 ? qsize : -qsize;
1.1 maekawa 116:
117: TMP_FREE (marker);
118: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>