xref: /dflybsd-src/contrib/gmp/mpn/generic/divrem.c (revision 86d7f5d305c6adaa56ff4582ece9859d73106103)
186d7f5d3SJohn Marino /* mpn_divrem -- Divide natural numbers, producing both remainder and
286d7f5d3SJohn Marino    quotient.  This is now just a middle layer for calling the new
386d7f5d3SJohn Marino    internal mpn_tdiv_qr.
486d7f5d3SJohn Marino 
586d7f5d3SJohn Marino Copyright 1993, 1994, 1995, 1996, 1997, 1999, 2000, 2001, 2002, 2005 Free
686d7f5d3SJohn Marino Software Foundation, Inc.
786d7f5d3SJohn Marino 
886d7f5d3SJohn Marino This file is part of the GNU MP Library.
986d7f5d3SJohn Marino 
1086d7f5d3SJohn Marino The GNU MP Library is free software; you can redistribute it and/or modify
1186d7f5d3SJohn Marino it under the terms of the GNU Lesser General Public License as published by
1286d7f5d3SJohn Marino the Free Software Foundation; either version 3 of the License, or (at your
1386d7f5d3SJohn Marino option) any later version.
1486d7f5d3SJohn Marino 
1586d7f5d3SJohn Marino The GNU MP Library is distributed in the hope that it will be useful, but
1686d7f5d3SJohn Marino WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
1786d7f5d3SJohn Marino or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
1886d7f5d3SJohn Marino License for more details.
1986d7f5d3SJohn Marino 
2086d7f5d3SJohn Marino You should have received a copy of the GNU Lesser General Public License
2186d7f5d3SJohn Marino along with the GNU MP Library.  If not, see http://www.gnu.org/licenses/.  */
2286d7f5d3SJohn Marino 
2386d7f5d3SJohn Marino #include "gmp.h"
2486d7f5d3SJohn Marino #include "gmp-impl.h"
2586d7f5d3SJohn Marino #include "longlong.h"
2686d7f5d3SJohn Marino 
2786d7f5d3SJohn Marino mp_limb_t
mpn_divrem(mp_ptr qp,mp_size_t qxn,mp_ptr np,mp_size_t nn,mp_srcptr dp,mp_size_t dn)2886d7f5d3SJohn Marino mpn_divrem (mp_ptr qp, mp_size_t qxn,
2986d7f5d3SJohn Marino 	    mp_ptr np, mp_size_t nn,
3086d7f5d3SJohn Marino 	    mp_srcptr dp, mp_size_t dn)
3186d7f5d3SJohn Marino {
3286d7f5d3SJohn Marino   ASSERT (qxn >= 0);
3386d7f5d3SJohn Marino   ASSERT (nn >= dn);
3486d7f5d3SJohn Marino   ASSERT (dn >= 1);
3586d7f5d3SJohn Marino   ASSERT (dp[dn-1] & GMP_NUMB_HIGHBIT);
3686d7f5d3SJohn Marino   ASSERT (! MPN_OVERLAP_P (np, nn, dp, dn));
3786d7f5d3SJohn Marino   ASSERT (! MPN_OVERLAP_P (qp, nn-dn+qxn, np, nn) || qp==np+dn+qxn);
3886d7f5d3SJohn Marino   ASSERT (! MPN_OVERLAP_P (qp, nn-dn+qxn, dp, dn));
3986d7f5d3SJohn Marino   ASSERT_MPN (np, nn);
4086d7f5d3SJohn Marino   ASSERT_MPN (dp, dn);
4186d7f5d3SJohn Marino 
4286d7f5d3SJohn Marino   if (dn == 1)
4386d7f5d3SJohn Marino     {
4486d7f5d3SJohn Marino       mp_limb_t ret;
4586d7f5d3SJohn Marino       mp_ptr q2p;
4686d7f5d3SJohn Marino       mp_size_t qn;
4786d7f5d3SJohn Marino       TMP_DECL;
4886d7f5d3SJohn Marino 
4986d7f5d3SJohn Marino       TMP_MARK;
5086d7f5d3SJohn Marino       q2p = TMP_ALLOC_LIMBS (nn + qxn);
5186d7f5d3SJohn Marino 
5286d7f5d3SJohn Marino       np[0] = mpn_divrem_1 (q2p, qxn, np, nn, dp[0]);
5386d7f5d3SJohn Marino       qn = nn + qxn - 1;
5486d7f5d3SJohn Marino       MPN_COPY (qp, q2p, qn);
5586d7f5d3SJohn Marino       ret = q2p[qn];
5686d7f5d3SJohn Marino 
5786d7f5d3SJohn Marino       TMP_FREE;
5886d7f5d3SJohn Marino       return ret;
5986d7f5d3SJohn Marino     }
6086d7f5d3SJohn Marino   else if (dn == 2)
6186d7f5d3SJohn Marino     {
6286d7f5d3SJohn Marino       return mpn_divrem_2 (qp, qxn, np, nn, dp);
6386d7f5d3SJohn Marino     }
6486d7f5d3SJohn Marino   else
6586d7f5d3SJohn Marino     {
6686d7f5d3SJohn Marino       mp_ptr rp, q2p;
6786d7f5d3SJohn Marino       mp_limb_t qhl;
6886d7f5d3SJohn Marino       mp_size_t qn;
6986d7f5d3SJohn Marino       TMP_DECL;
7086d7f5d3SJohn Marino 
7186d7f5d3SJohn Marino       TMP_MARK;
7286d7f5d3SJohn Marino       if (UNLIKELY (qxn != 0))
7386d7f5d3SJohn Marino 	{
7486d7f5d3SJohn Marino 	  mp_ptr n2p;
7586d7f5d3SJohn Marino 	  n2p = TMP_ALLOC_LIMBS (nn + qxn);
7686d7f5d3SJohn Marino 	  MPN_ZERO (n2p, qxn);
7786d7f5d3SJohn Marino 	  MPN_COPY (n2p + qxn, np, nn);
7886d7f5d3SJohn Marino 	  q2p = TMP_ALLOC_LIMBS (nn - dn + qxn + 1);
7986d7f5d3SJohn Marino 	  rp = TMP_ALLOC_LIMBS (dn);
8086d7f5d3SJohn Marino 	  mpn_tdiv_qr (q2p, rp, 0L, n2p, nn + qxn, dp, dn);
8186d7f5d3SJohn Marino 	  MPN_COPY (np, rp, dn);
8286d7f5d3SJohn Marino 	  qn = nn - dn + qxn;
8386d7f5d3SJohn Marino 	  MPN_COPY (qp, q2p, qn);
8486d7f5d3SJohn Marino 	  qhl = q2p[qn];
8586d7f5d3SJohn Marino 	}
8686d7f5d3SJohn Marino       else
8786d7f5d3SJohn Marino 	{
8886d7f5d3SJohn Marino 	  q2p = TMP_ALLOC_LIMBS (nn - dn + 1);
8986d7f5d3SJohn Marino 	  rp = TMP_ALLOC_LIMBS (dn);
9086d7f5d3SJohn Marino 	  mpn_tdiv_qr (q2p, rp, 0L, np, nn, dp, dn);
9186d7f5d3SJohn Marino 	  MPN_COPY (np, rp, dn);	/* overwrite np area with remainder */
9286d7f5d3SJohn Marino 	  qn = nn - dn;
9386d7f5d3SJohn Marino 	  MPN_COPY (qp, q2p, qn);
9486d7f5d3SJohn Marino 	  qhl = q2p[qn];
9586d7f5d3SJohn Marino 	}
9686d7f5d3SJohn Marino       TMP_FREE;
9786d7f5d3SJohn Marino       return qhl;
9886d7f5d3SJohn Marino     }
9986d7f5d3SJohn Marino }
100