rts/gmp/mpz/tdiv_qr.c

   1 /* mpz_tdiv_qr(quot,rem,dividend,divisor) -- Set QUOT to DIVIDEND/DIVISOR,
   2    and REM to DIVIDEND mod DIVISOR.
   3
   4 Copyright (C) 1991, 1993, 1994, 2000 Free Software Foundation, Inc.
   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
   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
  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
  15 or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
  16 License for more details.
  17
  18 You should have received a copy of the GNU Lesser General Public License
  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 #include "gmp.h"
  24 #include "gmp-impl.h"
  25 #include "longlong.h"
  26 #ifdef BERKELEY_MP
  27 #include "mp.h"
  28 #endif
  29
  30
  31 #ifndef BERKELEY_MP
  32
  33 void
  34 #if __STDC__
  35 mpz_tdiv_qr (mpz_ptr quot, mpz_ptr rem, mpz_srcptr num, mpz_srcptr den)
  36 #else
  37 mpz_tdiv_qr (quot, rem, num, den)
  38      mpz_ptr quot;
  39      mpz_ptr rem;
  40      mpz_srcptr num;
  41      mpz_srcptr den;
  42 #endif
  43
  44 #else /* BERKELEY_MP */
  45
  46 void
  47 #if __STDC__
  48 mdiv (mpz_srcptr num, mpz_srcptr den, mpz_ptr quot, mpz_ptr rem)
  49 #else
  50 mdiv (num, den, quot, rem)
  51      mpz_srcptr num;
  52      mpz_srcptr den;
  53      mpz_ptr    quot;
  54      mpz_ptr    rem;
  55 #endif
  56
  57 #endif /* BERKELEY_MP */
  58 {
  59   mp_size_t ql;
  60   mp_size_t ns, ds, nl, dl;
  61   mp_ptr np, dp, qp, rp;
  62   TMP_DECL (marker);
  63
  64   ns = SIZ (num);
  65   ds = SIZ (den);
  66   nl = ABS (ns);
  67   dl = ABS (ds);
  68   ql = nl - dl + 1;
  69
  70   if (dl == 0)
  71     DIVIDE_BY_ZERO;
  72
  73   MPZ_REALLOC (rem, dl);
  74
  75   if (ql <= 0)
  76     {
  77       if (num != rem)
  78         {
  79           mp_ptr np, rp;
  80           np = PTR (num);
  81           rp = PTR (rem);
  82           MPN_COPY (rp, np, nl);
  83           SIZ (rem) = SIZ (num);
  84         }
  85       /* This needs to follow the assignment to rem, in case the
  86          numerator and quotient are the same.  */
  87       SIZ (quot) = 0;
  88       return;
  89     }
  90
  91   MPZ_REALLOC (quot, ql);
  92
  93   TMP_MARK (marker);
  94   qp = PTR (quot);
  95   rp = PTR (rem);
  96   np = PTR (num);
  97   dp = PTR (den);
  98
  99   /* FIXME: We should think about how to handle the temporary allocation.
 100      Perhaps mpn_tdiv_qr should handle it, since it anyway often needs to
 101      allocate temp space.  */
 102
 103   /* Copy denominator to temporary space if it overlaps with the quotient
 104      or remainder.  */
 105   if (dp == rp || dp == qp)
 106     {
 107       mp_ptr tp;
 108       tp = (mp_ptr) TMP_ALLOC (dl * BYTES_PER_MP_LIMB);
 109       MPN_COPY (tp, dp, dl);
 110       dp = tp;
 111     }
 112   /* Copy numerator to temporary space if it overlaps with the quotient or
 113      remainder.  */
 114   if (np == rp || np == qp)
 115     {
 116       mp_ptr tp;
 117       tp = (mp_ptr) TMP_ALLOC (nl * BYTES_PER_MP_LIMB);
 118       MPN_COPY (tp, np, nl);
 119       np = tp;
 120     }
 121
 122   mpn_tdiv_qr (qp, rp, 0L, np, nl, dp, dl);
 123
 124   ql -=  qp[ql - 1] == 0;
 125   MPN_NORMALIZE (rp, dl);
 126
 127   SIZ (quot) = (ns ^ ds) >= 0 ? ql : -ql;
 128   SIZ (rem) = ns >= 0 ? dl : -dl;
 129   TMP_FREE (marker);
 130 }