--- /dev/null
+/* _mpz_realloc -- make the mpz_t have NEW_SIZE digits allocated.
+
+Copyright (C) 1991, 1993, 1994, 1995 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+void *
+#if __STDC__
+_mpz_realloc (mpz_ptr m, mp_size_t new_size)
+#else
+_mpz_realloc (m, new_size)
+ mpz_ptr m;
+ mp_size_t new_size;
+#endif
+{
+ /* Never allocate zero space. */
+ if (new_size == 0)
+ new_size = 1;
+
+ m->_mp_d = (mp_ptr) (*_mp_reallocate_func) (m->_mp_d,
+ m->_mp_alloc * BYTES_PER_MP_LIMB,
+ new_size * BYTES_PER_MP_LIMB);
+ m->_mp_alloc = new_size;
+
+#if 0
+ /* This might break some code that reads the size field after
+ reallocation, in the case the reallocated destination and a
+ source argument are identical. */
+ if (ABS (m->_mp_size) > new_size)
+ m->_mp_size = 0;
+#endif
+
+ return (void *) m->_mp_d;
+}
--- /dev/null
+/* mpz_scan0(op, startbit) -- Scan for the next set bit, starting at startbit.
+
+Copyright (C) 1996 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+unsigned long int
+#if __STDC__
+mpz_scan0 (mpz_srcptr u, unsigned long int starting_bit)
+#else
+mpz_scan0 (u, starting_bit)
+ mpz_srcptr u;
+ unsigned long int starting_bit;
+#endif
+{
+ return mpn_scan0 (u->_mp_d, starting_bit);
+}
--- /dev/null
+/* mpz_scan1(op, startbit) -- Scan for the next set bit, starting at startbit.
+
+Copyright (C) 1996 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+unsigned long int
+#if __STDC__
+mpz_scan1 (mpz_srcptr u, unsigned long int starting_bit)
+#else
+mpz_scan1 (u, starting_bit)
+ mpz_srcptr u;
+ unsigned long int starting_bit;
+#endif
+{
+ return mpn_scan1 (u->_mp_d, starting_bit);
+}
--- /dev/null
+/* mpz_set (dest_integer, src_integer) -- Assign DEST_INTEGER from SRC_INTEGER.
+
+Copyright (C) 1991, 1993, 1994, 1995 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+void
+#if __STDC__
+mpz_set (mpz_ptr w, mpz_srcptr u)
+#else
+mpz_set (w, u)
+ mpz_ptr w;
+ mpz_srcptr u;
+#endif
+{
+ mp_ptr wp, up;
+ mp_size_t usize, size;
+
+ usize = u->_mp_size;
+ size = ABS (usize);
+
+ if (w->_mp_alloc < size)
+ _mpz_realloc (w, size);
+
+ wp = w->_mp_d;
+ up = u->_mp_d;
+
+ MPN_COPY (wp, up, size);
+ w->_mp_size = usize;
+}
--- /dev/null
+/* mpz_set_d(integer, val) -- Assign INTEGER with a double value VAL.
+
+Copyright (C) 1995 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+void
+#if __STDC__
+mpz_set_d (mpz_ptr r, double d)
+#else
+mpz_set_d (r, d)
+ mpz_ptr r;
+ double d;
+#endif
+{
+ int negative;
+ mp_size_t size;
+ mp_limb_t tp[3];
+ mp_ptr rp;
+
+ negative = d < 0;
+ d = ABS (d);
+
+ /* Handle small arguments quickly. */
+ if (d < MP_BASE_AS_DOUBLE)
+ {
+ mp_limb_t tmp;
+ tmp = d;
+ PTR(r)[0] = tmp;
+ SIZ(r) = negative ? -(tmp != 0) : (tmp != 0);
+ return;
+ }
+
+ size = __gmp_extract_double (tp, d);
+
+ if (ALLOC(r) < size)
+ _mpz_realloc (r, size);
+
+ rp = PTR (r);
+
+#if BITS_PER_MP_LIMB == 32
+ switch (size)
+ {
+ default:
+ MPN_ZERO (rp, size - 3);
+ rp += size - 3;
+ case 3:
+ rp[2] = tp[2];
+ rp[1] = tp[1];
+ rp[0] = tp[0];
+ break;
+ case 2:
+ rp[1] = tp[2];
+ rp[0] = tp[1];
+ break;
+ case 1:
+ abort ();
+ }
+#else
+ switch (size)
+ {
+ default:
+ MPN_ZERO (rp, size - 2);
+ rp += size - 2;
+ case 2:
+ rp[1] = tp[1];
+ rp[0] = tp[0];
+ break;
+ case 1:
+ abort ();
+ }
+#endif
+
+ SIZ(r) = negative ? -size : size;
+}
--- /dev/null
+/* mpz_set_f (dest_integer, src_float) -- Assign DEST_INTEGER from SRC_FLOAT.
+
+Copyright (C) 1996 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+void
+#if __STDC__
+mpz_set_f (mpz_ptr w, mpf_srcptr u)
+#else
+mpz_set_f (w, u)
+ mpz_ptr w;
+ mpf_srcptr u;
+#endif
+{
+ mp_ptr wp, up;
+ mp_size_t usize, size;
+ mp_exp_t exp;
+
+ usize = SIZ (u);
+ size = ABS (usize);
+ exp = EXP (u);
+
+ if (w->_mp_alloc < exp)
+ _mpz_realloc (w, exp);
+
+ wp = w->_mp_d;
+ up = u->_mp_d;
+
+ if (exp <= 0)
+ {
+ SIZ (w) = 0;
+ return;
+ }
+ if (exp < size)
+ {
+ MPN_COPY (wp, up + size - exp, exp);
+ }
+ else
+ {
+ MPN_ZERO (wp, exp - size);
+ MPN_COPY (wp + exp - size, up, size);
+ }
+
+ w->_mp_size = usize >= 0 ? exp : -exp;
+}
--- /dev/null
+/* mpz_set_q (dest_integer, src_rational) -- Assign DEST_INTEGER from
+ SRC_rational.
+
+Copyright (C) 1996 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+void
+#if __STDC__
+mpz_set_q (mpz_ptr w, mpq_srcptr u)
+#else
+mpz_set_q (w, u)
+ mpz_ptr w;
+ mpq_srcptr u;
+#endif
+{
+ mpz_tdiv_q (w, mpq_numref (u), mpq_denref (u));
+}
--- /dev/null
+/* mpz_set_si(integer, val) -- Assign INTEGER with a small value VAL.
+
+Copyright (C) 1991, 1993, 1994, 1995 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+void
+#if __STDC__
+mpz_set_si (mpz_ptr dest, signed long int val)
+#else
+mpz_set_si (dest, val)
+ mpz_ptr dest;
+ signed long int val;
+#endif
+{
+ /* We don't check if the allocation is enough, since the rest of the
+ package ensures it's at least 1, which is what we need here. */
+ if (val > 0)
+ {
+ dest->_mp_d[0] = val;
+ dest->_mp_size = 1;
+ }
+ else if (val < 0)
+ {
+ dest->_mp_d[0] = -val;
+ dest->_mp_size = -1;
+ }
+ else
+ dest->_mp_size = 0;
+}
--- /dev/null
+/* mpz_set_str(mp_dest, string, base) -- Convert the \0-terminated
+ string STRING in base BASE to multiple precision integer in
+ MP_DEST. Allow white space in the string. If BASE == 0 determine
+ the base in the C standard way, i.e. 0xhh...h means base 16,
+ 0oo...o means base 8, otherwise assume base 10.
+
+Copyright (C) 1991, 1993, 1994, Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include <ctype.h>
+#include "gmp.h"
+#include "gmp-impl.h"
+#include "longlong.h"
+
+static int
+digit_value_in_base (c, base)
+ int c;
+ int base;
+{
+ int digit;
+
+ if (isdigit (c))
+ digit = c - '0';
+ else if (islower (c))
+ digit = c - 'a' + 10;
+ else if (isupper (c))
+ digit = c - 'A' + 10;
+ else
+ return -1;
+
+ if (digit < base)
+ return digit;
+ return -1;
+}
+
+int
+#if __STDC__
+mpz_set_str (mpz_ptr x, const char *str, int base)
+#else
+mpz_set_str (x, str, base)
+ mpz_ptr x;
+ const char *str;
+ int base;
+#endif
+{
+ size_t str_size;
+ char *s, *begs;
+ size_t i;
+ mp_size_t xsize;
+ int c;
+ int negative;
+ TMP_DECL (marker);
+
+ /* Skip whitespace. */
+ do
+ c = *str++;
+ while (isspace (c));
+
+ negative = 0;
+ if (c == '-')
+ {
+ negative = 1;
+ c = *str++;
+ }
+
+ if (digit_value_in_base (c, base == 0 ? 10 : base) < 0)
+ return -1; /* error if no digits */
+
+ /* If BASE is 0, try to find out the base by looking at the initial
+ characters. */
+ if (base == 0)
+ {
+ base = 10;
+ if (c == '0')
+ {
+ base = 8;
+ c = *str++;
+ if (c == 'x' || c == 'X')
+ {
+ base = 16;
+ c = *str++;
+ }
+ }
+ }
+
+ TMP_MARK (marker);
+ str_size = strlen (str - 1);
+ s = begs = (char *) TMP_ALLOC (str_size + 1);
+
+ for (i = 0; i < str_size; i++)
+ {
+ if (!isspace (c))
+ {
+ int dig = digit_value_in_base (c, base);
+ if (dig < 0)
+ {
+ TMP_FREE (marker);
+ return -1;
+ }
+ *s++ = dig;
+ }
+ c = *str++;
+ }
+
+ str_size = s - begs;
+
+ xsize = str_size / __mp_bases[base].chars_per_limb + 1;
+ if (x->_mp_alloc < xsize)
+ _mpz_realloc (x, xsize);
+
+ xsize = mpn_set_str (x->_mp_d, (unsigned char *) begs, str_size, base);
+ x->_mp_size = negative ? -xsize : xsize;
+
+ TMP_FREE (marker);
+ return 0;
+}
--- /dev/null
+/* mpz_set_ui(integer, val) -- Assign INTEGER with a small value VAL.
+
+Copyright (C) 1991, 1993, 1994, 1995 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+void
+#if __STDC__
+mpz_set_ui (mpz_ptr dest, unsigned long int val)
+#else
+mpz_set_ui (dest, val)
+ mpz_ptr dest;
+ unsigned long int val;
+#endif
+{
+ /* We don't check if the allocation is enough, since the rest of the
+ package ensures it's at least 1, which is what we need here. */
+ if (val > 0)
+ {
+ dest->_mp_d[0] = val;
+ dest->_mp_size = 1;
+ }
+ else
+ dest->_mp_size = 0;
+}
--- /dev/null
+/* mpz_setbit -- set a specified bit.
+
+Copyright (C) 1991, 1993, 1994, 1995 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+void
+#if __STDC__
+mpz_setbit (mpz_ptr d, unsigned long int bit_index)
+#else
+mpz_setbit (d, bit_index)
+ mpz_ptr d;
+ unsigned long int bit_index;
+#endif
+{
+ mp_size_t dsize = d->_mp_size;
+ mp_ptr dp = d->_mp_d;
+ mp_size_t limb_index;
+
+ limb_index = bit_index / BITS_PER_MP_LIMB;
+ if (dsize >= 0)
+ {
+ if (limb_index < dsize)
+ {
+ dp[limb_index] |= (mp_limb_t) 1 << (bit_index % BITS_PER_MP_LIMB);
+ d->_mp_size = dsize;
+ }
+ else
+ {
+ /* Ugh. The bit should be set outside of the end of the
+ number. We have to increase the size of the number. */
+ if (d->_mp_alloc < limb_index + 1)
+ {
+ _mpz_realloc (d, limb_index + 1);
+ dp = d->_mp_d;
+ }
+ MPN_ZERO (dp + dsize, limb_index - dsize);
+ dp[limb_index] = (mp_limb_t) 1 << (bit_index % BITS_PER_MP_LIMB);
+ d->_mp_size = limb_index + 1;
+ }
+ }
+ else
+ {
+ mp_size_t zero_bound;
+
+ /* Simulate two's complement arithmetic, i.e. simulate
+ 1. Set OP = ~(OP - 1) [with infinitely many leading ones].
+ 2. set the bit.
+ 3. Set OP = ~OP + 1. */
+
+ dsize = -dsize;
+
+ /* No upper bound on this loop, we're sure there's a non-zero limb
+ sooner ot later. */
+ for (zero_bound = 0; ; zero_bound++)
+ if (dp[zero_bound] != 0)
+ break;
+
+ if (limb_index > zero_bound)
+ {
+ if (limb_index < dsize)
+ dp[limb_index] &= ~((mp_limb_t) 1 << (bit_index % BITS_PER_MP_LIMB));
+ else
+ ;
+ }
+ else if (limb_index == zero_bound)
+ {
+ dp[limb_index] = ((dp[limb_index] - 1)
+ & ~((mp_limb_t) 1 << (bit_index % BITS_PER_MP_LIMB))) + 1;
+ if (dp[limb_index] == 0)
+ {
+ mp_size_t i;
+ for (i = limb_index + 1; i < dsize; i++)
+ {
+ dp[i] += 1;
+ if (dp[i] != 0)
+ goto fin;
+ }
+ /* We got carry all way out beyond the end of D. Increase
+ its size (and allocation if necessary). */
+ dsize++;
+ if (d->_mp_alloc < dsize)
+ {
+ _mpz_realloc (d, dsize);
+ dp = d->_mp_d;
+ }
+ dp[i] = 1;
+ d->_mp_size = -dsize;
+ fin:;
+ }
+ }
+ else
+ ;
+ }
+}
--- /dev/null
+/* mpz_size(x) -- return the number of lims currently used by the
+ value of integer X.
+
+Copyright (C) 1991, 1993, 1994, 1995 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+size_t
+#if __STDC__
+mpz_size (mpz_srcptr x)
+#else
+mpz_size (x)
+ mpz_srcptr x;
+#endif
+{
+ return ABS (x->_mp_size);
+}
--- /dev/null
+/* mpz_sizeinbase(x, base) -- return an approximation to the number of
+ character the integer X would have printed in base BASE. The
+ approximation is never too small.
+
+Copyright (C) 1991, 1993, 1994, 1995 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+#include "longlong.h"
+
+size_t
+#if __STDC__
+mpz_sizeinbase (mpz_srcptr x, int base)
+#else
+mpz_sizeinbase (x, base)
+ mpz_srcptr x;
+ int base;
+#endif
+{
+ mp_size_t size = ABS (x->_mp_size);
+ int lb_base, cnt;
+ size_t totbits;
+
+ /* Special case for X == 0. */
+ if (size == 0)
+ return 1;
+
+ /* Calculate the total number of significant bits of X. */
+ count_leading_zeros (cnt, x->_mp_d[size - 1]);
+ totbits = size * BITS_PER_MP_LIMB - cnt;
+
+ if ((base & (base - 1)) == 0)
+ {
+ /* Special case for powers of 2, giving exact result. */
+
+ count_leading_zeros (lb_base, base);
+ lb_base = BITS_PER_MP_LIMB - lb_base - 1;
+
+ return (totbits + lb_base - 1) / lb_base;
+ }
+ else
+ return (size_t) (totbits * __mp_bases[base].chars_per_bit_exactly) + 1;
+}
--- /dev/null
+/* mpz_sqrt(root, u) -- Set ROOT to floor(sqrt(U)).
+
+Copyright (C) 1991, 1993, 1994, 1996 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+void
+#if __STDC__
+mpz_sqrt (mpz_ptr root, mpz_srcptr op)
+#else
+mpz_sqrt (root, op)
+ mpz_ptr root;
+ mpz_srcptr op;
+#endif
+{
+ mp_size_t op_size, root_size;
+ mp_ptr root_ptr, op_ptr;
+ mp_ptr free_me = NULL;
+ mp_size_t free_me_size;
+ TMP_DECL (marker);
+
+ TMP_MARK (marker);
+ op_size = op->_mp_size;
+ if (op_size < 0)
+ op_size = 1 / (op_size > 0); /* Divide by zero for negative OP. */
+
+ /* The size of the root is accurate after this simple calculation. */
+ root_size = (op_size + 1) / 2;
+
+ root_ptr = root->_mp_d;
+ op_ptr = op->_mp_d;
+
+ if (root->_mp_alloc < root_size)
+ {
+ if (root_ptr == op_ptr)
+ {
+ free_me = root_ptr;
+ free_me_size = root->_mp_alloc;
+ }
+ else
+ (*_mp_free_func) (root_ptr, root->_mp_alloc * BYTES_PER_MP_LIMB);
+
+ root->_mp_alloc = root_size;
+ root_ptr = (mp_ptr) (*_mp_allocate_func) (root_size * BYTES_PER_MP_LIMB);
+ root->_mp_d = root_ptr;
+ }
+ else
+ {
+ /* Make OP not overlap with ROOT. */
+ if (root_ptr == op_ptr)
+ {
+ /* ROOT and OP are identical. Allocate temporary space for OP. */
+ op_ptr = (mp_ptr) TMP_ALLOC (op_size * BYTES_PER_MP_LIMB);
+ /* Copy to the temporary space. Hack: Avoid temporary variable
+ by using ROOT_PTR. */
+ MPN_COPY (op_ptr, root_ptr, op_size);
+ }
+ }
+
+ mpn_sqrtrem (root_ptr, NULL, op_ptr, op_size);
+
+ root->_mp_size = root_size;
+
+ if (free_me != NULL)
+ (*_mp_free_func) (free_me, free_me_size * BYTES_PER_MP_LIMB);
+ TMP_FREE (marker);
+}
--- /dev/null
+/* mpz_sqrtrem(root,rem,x) -- Set ROOT to floor(sqrt(X)) and REM
+ to the remainder, i.e. X - ROOT**2.
+
+Copyright (C) 1991, 1993, 1994, 1996 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+#ifndef BERKELEY_MP
+void
+#if __STDC__
+mpz_sqrtrem (mpz_ptr root, mpz_ptr rem, mpz_srcptr op)
+#else
+mpz_sqrtrem (root, rem, op)
+ mpz_ptr root;
+ mpz_ptr rem;
+ mpz_srcptr op;
+#endif
+#else /* BERKELEY_MP */
+void
+#if __STDC__
+msqrt (mpz_srcptr op, mpz_ptr root, mpz_ptr rem)
+#else
+msqrt (op, root, rem)
+ mpz_srcptr op;
+ mpz_ptr root;
+ mpz_ptr rem;
+#endif
+#endif /* BERKELEY_MP */
+{
+ mp_size_t op_size, root_size, rem_size;
+ mp_ptr root_ptr, op_ptr;
+ mp_ptr free_me = NULL;
+ mp_size_t free_me_size;
+ TMP_DECL (marker);
+
+ TMP_MARK (marker);
+ op_size = op->_mp_size;
+ if (op_size < 0)
+ op_size = 1 / (op_size > 0); /* Divide by zero for negative OP. */
+
+ if (rem->_mp_alloc < op_size)
+ _mpz_realloc (rem, op_size);
+
+ /* The size of the root is accurate after this simple calculation. */
+ root_size = (op_size + 1) / 2;
+
+ root_ptr = root->_mp_d;
+ op_ptr = op->_mp_d;
+
+ if (root->_mp_alloc < root_size)
+ {
+ if (root_ptr == op_ptr)
+ {
+ free_me = root_ptr;
+ free_me_size = root->_mp_alloc;
+ }
+ else
+ (*_mp_free_func) (root_ptr, root->_mp_alloc * BYTES_PER_MP_LIMB);
+
+ root->_mp_alloc = root_size;
+ root_ptr = (mp_ptr) (*_mp_allocate_func) (root_size * BYTES_PER_MP_LIMB);
+ root->_mp_d = root_ptr;
+ }
+ else
+ {
+ /* Make OP not overlap with ROOT. */
+ if (root_ptr == op_ptr)
+ {
+ /* ROOT and OP are identical. Allocate temporary space for OP. */
+ op_ptr = (mp_ptr) TMP_ALLOC (op_size * BYTES_PER_MP_LIMB);
+ /* Copy to the temporary space. Hack: Avoid temporary variable
+ by using ROOT_PTR. */
+ MPN_COPY (op_ptr, root_ptr, op_size);
+ }
+ }
+
+ rem_size = mpn_sqrtrem (root_ptr, rem->_mp_d, op_ptr, op_size);
+
+ root->_mp_size = root_size;
+
+ /* Write remainder size last, to enable us to define this function to
+ give only the square root remainder, if the user calls if with
+ ROOT == REM. */
+ rem->_mp_size = rem_size;
+
+ if (free_me != NULL)
+ (*_mp_free_func) (free_me, free_me_size * BYTES_PER_MP_LIMB);
+ TMP_FREE (marker);
+}
--- /dev/null
+/* mpz_sub -- Subtract two integers.
+
+Copyright (C) 1991, 1993, 1994, 1996 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+#ifndef BERKELEY_MP
+void
+#if __STDC__
+mpz_sub (mpz_ptr w, mpz_srcptr u, mpz_srcptr v)
+#else
+mpz_sub (w, u, v)
+ mpz_ptr w;
+ mpz_srcptr u;
+ mpz_srcptr v;
+#endif
+#else /* BERKELEY_MP */
+void
+#if __STDC__
+msub (mpz_srcptr u, mpz_srcptr v, mpz_ptr w)
+#else
+msub (u, v, w)
+ mpz_srcptr u;
+ mpz_srcptr v;
+ mpz_ptr w;
+#endif
+#endif /* BERKELEY_MP */
+{
+ mp_srcptr up, vp;
+ mp_ptr wp;
+ mp_size_t usize, vsize, wsize;
+ mp_size_t abs_usize;
+ mp_size_t abs_vsize;
+
+ usize = u->_mp_size;
+ vsize = -v->_mp_size; /* The "-" makes the difference from mpz_add */
+ abs_usize = ABS (usize);
+ abs_vsize = ABS (vsize);
+
+ if (abs_usize < abs_vsize)
+ {
+ /* Swap U and V. */
+ {const __mpz_struct *t = u; u = v; v = t;}
+ {mp_size_t t = usize; usize = vsize; vsize = t;}
+ {mp_size_t t = abs_usize; abs_usize = abs_vsize; abs_vsize = t;}
+ }
+
+ /* True: ABS_USIZE >= ABS_VSIZE. */
+
+ /* If not space for w (and possible carry), increase space. */
+ wsize = abs_usize + 1;
+ if (w->_mp_alloc < wsize)
+ _mpz_realloc (w, wsize);
+
+ /* These must be after realloc (u or v may be the same as w). */
+ up = u->_mp_d;
+ vp = v->_mp_d;
+ wp = w->_mp_d;
+
+ if ((usize ^ vsize) < 0)
+ {
+ /* U and V have different sign. Need to compare them to determine
+ which operand to subtract from which. */
+
+ /* This test is right since ABS_USIZE >= ABS_VSIZE. */
+ if (abs_usize != abs_vsize)
+ {
+ mpn_sub (wp, up, abs_usize, vp, abs_vsize);
+ wsize = abs_usize;
+ MPN_NORMALIZE (wp, wsize);
+ if (usize < 0)
+ wsize = -wsize;
+ }
+ else if (mpn_cmp (up, vp, abs_usize) < 0)
+ {
+ mpn_sub_n (wp, vp, up, abs_usize);
+ wsize = abs_usize;
+ MPN_NORMALIZE (wp, wsize);
+ if (usize >= 0)
+ wsize = -wsize;
+ }
+ else
+ {
+ mpn_sub_n (wp, up, vp, abs_usize);
+ wsize = abs_usize;
+ MPN_NORMALIZE (wp, wsize);
+ if (usize < 0)
+ wsize = -wsize;
+ }
+ }
+ else
+ {
+ /* U and V have same sign. Add them. */
+ mp_limb_t cy_limb = mpn_add (wp, up, abs_usize, vp, abs_vsize);
+ wp[abs_usize] = cy_limb;
+ wsize = abs_usize + cy_limb;
+ if (usize < 0)
+ wsize = -wsize;
+ }
+
+ w->_mp_size = wsize;
+}
--- /dev/null
+/* mpz_sub_ui -- Subtract an unsigned one-word integer from an MP_INT.
+
+Copyright (C) 1991, 1993, 1994, 1996 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+void
+#if __STDC__
+mpz_sub_ui (mpz_ptr w, mpz_srcptr u, unsigned long int v)
+#else
+mpz_sub_ui (w, u, v)
+ mpz_ptr w;
+ mpz_srcptr u;
+ unsigned long int v;
+#endif
+{
+ mp_srcptr up;
+ mp_ptr wp;
+ mp_size_t usize, wsize;
+ mp_size_t abs_usize;
+
+ usize = u->_mp_size;
+ abs_usize = ABS (usize);
+
+ /* If not space for W (and possible carry), increase space. */
+ wsize = abs_usize + 1;
+ if (w->_mp_alloc < wsize)
+ _mpz_realloc (w, wsize);
+
+ /* These must be after realloc (U may be the same as W). */
+ up = u->_mp_d;
+ wp = w->_mp_d;
+
+ if (abs_usize == 0)
+ {
+ wp[0] = v;
+ w->_mp_size = -(v != 0);
+ return;
+ }
+
+ if (usize < 0)
+ {
+ mp_limb_t cy;
+ cy = mpn_add_1 (wp, up, abs_usize, v);
+ wp[abs_usize] = cy;
+ wsize = -(abs_usize + cy);
+ }
+ else
+ {
+ /* The signs are different. Need exact comparison to determine
+ which operand to subtract from which. */
+ if (abs_usize == 1 && up[0] < v)
+ {
+ wp[0] = v - up[0];
+ wsize = -1;
+ }
+ else
+ {
+ mpn_sub_1 (wp, up, abs_usize, v);
+ /* Size can decrease with at most one limb. */
+ wsize = abs_usize - (wp[abs_usize - 1] == 0);
+ }
+ }
+
+ w->_mp_size = wsize;
+}
--- /dev/null
+/* mpz_tdiv_q -- divide two integers and produce a quotient.
+
+Copyright (C) 1991, 1993, 1994, 1996 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+#include "longlong.h"
+
+void
+#if __STDC__
+mpz_tdiv_q (mpz_ptr quot, mpz_srcptr num, mpz_srcptr den)
+#else
+mpz_tdiv_q (quot, num, den)
+ mpz_ptr quot;
+ mpz_srcptr num;
+ mpz_srcptr den;
+#endif
+{
+ mp_srcptr np, dp;
+ mp_ptr qp, rp;
+ mp_size_t nsize = num->_mp_size;
+ mp_size_t dsize = den->_mp_size;
+ mp_size_t qsize, rsize;
+ mp_size_t sign_quotient = nsize ^ dsize;
+ unsigned normalization_steps;
+ mp_limb_t q_limb;
+ TMP_DECL (marker);
+
+ nsize = ABS (nsize);
+ dsize = ABS (dsize);
+
+ /* Ensure space is enough for quotient. */
+
+ qsize = nsize - dsize + 1; /* qsize cannot be bigger than this. */
+ if (qsize <= 0)
+ {
+ quot->_mp_size = 0;
+ return;
+ }
+
+ if (quot->_mp_alloc < qsize)
+ _mpz_realloc (quot, qsize);
+
+ qp = quot->_mp_d;
+ np = num->_mp_d;
+ dp = den->_mp_d;
+
+ /* Optimize division by a single-limb divisor. */
+ if (dsize == 1)
+ {
+ mpn_divmod_1 (qp, np, nsize, dp[0]);
+ qsize -= qp[qsize - 1] == 0;
+ quot->_mp_size = sign_quotient >= 0 ? qsize : -qsize;
+ return;
+ }
+
+ TMP_MARK (marker);
+
+ rp = (mp_ptr) TMP_ALLOC ((nsize + 1) * BYTES_PER_MP_LIMB);
+
+ count_leading_zeros (normalization_steps, dp[dsize - 1]);
+
+ /* Normalize the denominator, i.e. make its most significant bit set by
+ shifting it NORMALIZATION_STEPS bits to the left. Also shift the
+ numerator the same number of steps (to keep the quotient the same!). */
+ if (normalization_steps != 0)
+ {
+ mp_ptr tp;
+ mp_limb_t nlimb;
+
+ /* Shift up the denominator setting the most significant bit of
+ the most significant word. Use temporary storage not to clobber
+ the original contents of the denominator. */
+ tp = (mp_ptr) TMP_ALLOC (dsize * BYTES_PER_MP_LIMB);
+ mpn_lshift (tp, dp, dsize, normalization_steps);
+ dp = tp;
+
+ /* Shift up the numerator, possibly introducing a new most
+ significant word. Move the shifted numerator in the remainder
+ meanwhile. */
+ nlimb = mpn_lshift (rp, np, nsize, normalization_steps);
+ if (nlimb != 0)
+ {
+ rp[nsize] = nlimb;
+ rsize = nsize + 1;
+ }
+ else
+ rsize = nsize;
+ }
+ else
+ {
+ /* The denominator is already normalized, as required. Copy it to
+ temporary space if it overlaps with the quotient. */
+ if (dp == qp)
+ {
+ dp = (mp_ptr) TMP_ALLOC (dsize * BYTES_PER_MP_LIMB);
+ MPN_COPY ((mp_ptr) dp, qp, dsize);
+ }
+
+ /* Move the numerator to the remainder. */
+ MPN_COPY (rp, np, nsize);
+ rsize = nsize;
+ }
+
+ q_limb = mpn_divmod (qp, rp, rsize, dp, dsize);
+
+ qsize = rsize - dsize;
+ if (q_limb)
+ {
+ qp[qsize] = q_limb;
+ qsize += 1;
+ }
+
+ quot->_mp_size = sign_quotient >= 0 ? qsize : -qsize;
+ TMP_FREE (marker);
+}
--- /dev/null
+/* mpz_tdiv_q_2exp -- Divide an integer by 2**CNT. Round the quotient
+ towards -infinity.
+
+Copyright (C) 1991, 1993, 1994, 1996 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+void
+#if __STDC__
+mpz_tdiv_q_2exp (mpz_ptr w, mpz_srcptr u, unsigned long int cnt)
+#else
+mpz_tdiv_q_2exp (w, u, cnt)
+ mpz_ptr w;
+ mpz_srcptr u;
+ unsigned long int cnt;
+#endif
+{
+ mp_size_t usize, wsize;
+ mp_size_t limb_cnt;
+
+ usize = u->_mp_size;
+ limb_cnt = cnt / BITS_PER_MP_LIMB;
+ wsize = ABS (usize) - limb_cnt;
+ if (wsize <= 0)
+ w->_mp_size = 0;
+ else
+ {
+ mp_ptr wp;
+ mp_srcptr up;
+
+ if (w->_mp_alloc < wsize)
+ _mpz_realloc (w, wsize);
+
+ wp = w->_mp_d;
+ up = u->_mp_d;
+
+ cnt %= BITS_PER_MP_LIMB;
+ if (cnt != 0)
+ {
+ mpn_rshift (wp, up + limb_cnt, wsize, cnt);
+ wsize -= wp[wsize - 1] == 0;
+ }
+ else
+ {
+ MPN_COPY_INCR (wp, up + limb_cnt, wsize);
+ }
+
+ w->_mp_size = usize >= 0 ? wsize : -wsize;
+ }
+}
--- /dev/null
+/* mpz_tdiv_q_ui(quot, dividend, divisor_limb)
+ -- Divide DIVIDEND by DIVISOR_LIMB and store the result in QUOT.
+
+Copyright (C) 1991, 1993, 1994, 1996 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+void
+#if __STDC__
+mpz_tdiv_q_ui (mpz_ptr quot, mpz_srcptr dividend, unsigned long int divisor)
+#else
+mpz_tdiv_q_ui (quot, dividend, divisor)
+ mpz_ptr quot;
+ mpz_srcptr dividend;
+ unsigned long int divisor;
+#endif
+{
+ mp_size_t dividend_size;
+ mp_size_t size;
+ mp_ptr quot_ptr;
+
+ dividend_size = dividend->_mp_size;
+ size = ABS (dividend_size);
+
+ if (size == 0)
+ {
+ quot->_mp_size = 0;
+ return;
+ }
+
+ /* No need for temporary allocation and copying if QUOT == DIVIDEND as
+ the divisor is just one limb, and thus no intermediate remainders
+ need to be stored. */
+
+ if (quot->_mp_alloc < size)
+ _mpz_realloc (quot, size);
+
+ quot_ptr = quot->_mp_d;
+
+ mpn_divmod_1 (quot_ptr, dividend->_mp_d, size, (mp_limb_t) divisor);
+
+ /* The quotient is SIZE limbs, but the most significant might be zero. */
+ size -= quot_ptr[size - 1] == 0;
+ quot->_mp_size = dividend_size >= 0 ? size : -size;
+}
--- /dev/null
+/* mpz_tdiv_qr(quot,rem,dividend,divisor) -- Set QUOT to DIVIDEND/DIVISOR,
+ and REM to DIVIDEND mod DIVISOR.
+
+Copyright (C) 1991, 1993, 1994 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+#include "longlong.h"
+
+void
+#if __STDC__
+mpz_tdiv_qr (mpz_ptr quot, mpz_ptr rem, mpz_srcptr num, mpz_srcptr den)
+#else
+mpz_tdiv_qr (quot, rem, num, den)
+ mpz_ptr quot;
+ mpz_ptr rem;
+ mpz_srcptr num;
+ mpz_srcptr den;
+#endif
+
+#define COMPUTE_QUOTIENT
+#include "dmincl.c"
--- /dev/null
+/* mpz_tdiv_qr_ui(quot,rem,dividend,short_divisor) --
+ Set QUOT to DIVIDEND / SHORT_DIVISOR
+ and REM to DIVIDEND mod SHORT_DIVISOR.
+
+Copyright (C) 1991, 1993, 1994, 1996 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+void
+#if __STDC__
+mpz_tdiv_qr_ui (mpz_ptr quot, mpz_ptr rem, mpz_srcptr dividend, unsigned long int divisor)
+#else
+mpz_tdiv_qr_ui (quot, rem, dividend, divisor)
+ mpz_ptr quot;
+ mpz_ptr rem;
+ mpz_srcptr dividend;
+ unsigned long int divisor;
+#endif
+{
+ mp_size_t dividend_size;
+ mp_size_t size;
+ mp_ptr quot_ptr;
+ mp_limb_t remainder_limb;
+
+ dividend_size = dividend->_mp_size;
+ size = ABS (dividend_size);
+
+ if (size == 0)
+ {
+ quot->_mp_size = 0;
+ rem->_mp_size = 0;
+ return;
+ }
+
+ /* No need for temporary allocation and copying if QUOT == DIVIDEND as
+ the divisor is just one limb, and thus no intermediate remainders
+ need to be stored. */
+
+ if (quot->_mp_alloc < size)
+ _mpz_realloc (quot, size);
+
+ quot_ptr = quot->_mp_d;
+
+ remainder_limb = mpn_divmod_1 (quot_ptr, dividend->_mp_d, size,
+ (mp_limb_t) divisor);
+
+ if (remainder_limb == 0)
+ rem->_mp_size = 0;
+ else
+ {
+ /* Store the single-limb remainder. We don't check if there's space
+ for just one limb, since no function ever makes zero space. */
+ rem->_mp_size = dividend_size >= 0 ? 1 : -1;
+ rem->_mp_d[0] = remainder_limb;
+ }
+
+ /* The quotient is SIZE limbs, but the most significant might be zero. */
+ size -= quot_ptr[size - 1] == 0;
+ quot->_mp_size = dividend_size >= 0 ? size : -size;
+}
--- /dev/null
+/* mpz_tdiv_r(rem, dividend, divisor) -- Set REM to DIVIDEND mod DIVISOR.
+
+Copyright (C) 1991, 1993, 1994 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+#include "longlong.h"
+
+void
+#if __STDC__
+mpz_tdiv_r (mpz_ptr rem, mpz_srcptr num, mpz_srcptr den)
+#else
+mpz_tdiv_r (rem, num, den)
+ mpz_ptr rem;
+ mpz_srcptr num;
+ mpz_srcptr den;
+#endif
+
+#undef COMPUTE_QUOTIENT
+#include "dmincl.c"
--- /dev/null
+/* mpz_tdiv_r_2exp -- Divide a integer by 2**CNT and produce a remainder.
+
+Copyright (C) 1991, 1993, 1994, 1995 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+void
+#if __STDC__
+mpz_tdiv_r_2exp (mpz_ptr res, mpz_srcptr in, unsigned long int cnt)
+#else
+mpz_tdiv_r_2exp (res, in, cnt)
+ mpz_ptr res;
+ mpz_srcptr in;
+ unsigned long int cnt;
+#endif
+{
+ mp_size_t in_size = ABS (in->_mp_size);
+ mp_size_t res_size;
+ mp_size_t limb_cnt = cnt / BITS_PER_MP_LIMB;
+ mp_srcptr in_ptr = in->_mp_d;
+
+ if (in_size > limb_cnt)
+ {
+ /* The input operand is (probably) greater than 2**CNT. */
+ mp_limb_t x;
+
+ x = in_ptr[limb_cnt] & (((mp_limb_t) 1 << cnt % BITS_PER_MP_LIMB) - 1);
+ if (x != 0)
+ {
+ res_size = limb_cnt + 1;
+ if (res->_mp_alloc < res_size)
+ _mpz_realloc (res, res_size);
+
+ res->_mp_d[limb_cnt] = x;
+ }
+ else
+ {
+ res_size = limb_cnt;
+ MPN_NORMALIZE (in_ptr, res_size);
+
+ if (res->_mp_alloc < res_size)
+ _mpz_realloc (res, res_size);
+
+ limb_cnt = res_size;
+ }
+ }
+ else
+ {
+ /* The input operand is smaller than 2**CNT. We perform a no-op,
+ apart from that we might need to copy IN to RES. */
+ res_size = in_size;
+ if (res->_mp_alloc < res_size)
+ _mpz_realloc (res, res_size);
+
+ limb_cnt = res_size;
+ }
+
+ if (res != in)
+ MPN_COPY (res->_mp_d, in->_mp_d, limb_cnt);
+ res->_mp_size = in->_mp_size >= 0 ? res_size : -res_size;
+}
--- /dev/null
+/* mpz_tdiv_r_ui(rem, dividend, divisor_limb)
+ -- Set REM to DIVDEND mod DIVISOR_LIMB.
+
+Copyright (C) 1991, 1993, 1994, 1996 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+
+void
+#if __STDC__
+mpz_tdiv_r_ui (mpz_ptr rem, mpz_srcptr dividend, unsigned long int divisor)
+#else
+mpz_tdiv_r_ui (rem, dividend, divisor)
+ mpz_ptr rem;
+ mpz_srcptr dividend;
+ unsigned long int divisor;
+#endif
+{
+ mp_size_t dividend_size;
+ mp_size_t size;
+ mp_limb_t remainder_limb;
+
+ dividend_size = dividend->_mp_size;
+ size = ABS (dividend_size);
+
+ if (size == 0)
+ {
+ rem->_mp_size = 0;
+ return;
+ }
+
+ /* No need for temporary allocation and copying if QUOT == DIVIDEND as
+ the divisor is just one limb, and thus no intermediate remainders
+ need to be stored. */
+
+ remainder_limb = mpn_mod_1 (dividend->_mp_d, size, (mp_limb_t) divisor);
+
+ if (remainder_limb == 0)
+ rem->_mp_size = 0;
+ else
+ {
+ /* Store the single-limb remainder. We don't check if there's space
+ for just one limb, since no function ever makes zero space. */
+ rem->_mp_size = dividend_size >= 0 ? 1 : -1;
+ rem->_mp_d[0] = remainder_limb;
+ }
+}
--- /dev/null
+/* mpz_ui_pow_ui(res, base, exp) -- Set RES to BASE**EXP.
+
+Copyright (C) 1991, 1993, 1994, 1996 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the GNU MP Library; see the file COPYING.LIB. If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
+#include "gmp.h"
+#include "gmp-impl.h"
+#include "longlong.h"
+
+void
+#if __STDC__
+mpz_ui_pow_ui (mpz_ptr r, unsigned long int b, unsigned long int e)
+#else
+mpz_ui_pow_ui (r, b, e)
+ mpz_ptr r;
+ unsigned long int b;
+ unsigned long int e;
+#endif
+{
+ mp_ptr rp, tp, xp;
+ mp_size_t rsize;
+ int cnt, i;
+ mp_limb_t blimb = b;
+ TMP_DECL (marker);
+
+ /* Single out cases that give result == 0 or 1. These tests are here
+ to simplify the general code below, not to optimize. */
+ if (e == 0)
+ {
+ r->_mp_d[0] = 1;
+ r->_mp_size = 1;
+ return;
+ }
+ if (blimb == 0)
+ {
+ r->_mp_size = 0;
+ return;
+ }
+
+ if (blimb < 0x100)
+ {
+ /* Estimate space requirements accurately. Using the code from the
+ `else' path would over-estimate space requirements wildly. */
+ float lb = __mp_bases[blimb].chars_per_bit_exactly;
+ rsize = 2 + ((mp_size_t) (e / lb) / BITS_PER_MP_LIMB);
+ }
+ else
+ {
+ /* Over-estimate space requirements somewhat. */
+ count_leading_zeros (cnt, blimb);
+ rsize = e - cnt * e / BITS_PER_MP_LIMB + 1;
+ }
+
+ TMP_MARK (marker);
+
+ /* The two areas are used to alternatingly hold the input and recieve the
+ product for mpn_mul. (This scheme is used to fulfill the requirements
+ of mpn_mul; that the product space may not be the same as any of the
+ input operands.) */
+ rp = (mp_ptr) TMP_ALLOC (rsize * BYTES_PER_MP_LIMB);
+ tp = (mp_ptr) TMP_ALLOC (rsize * BYTES_PER_MP_LIMB);
+
+ rp[0] = blimb;
+ rsize = 1;
+ count_leading_zeros (cnt, e);
+
+ for (i = BITS_PER_MP_LIMB - cnt - 2; i >= 0; i--)
+ {
+ mpn_mul_n (tp, rp, rp, rsize);
+ rsize = 2 * rsize;
+ rsize -= tp[rsize - 1] == 0;
+ xp = tp; tp = rp; rp = xp;
+
+ if ((e & ((mp_limb_t) 1 << i)) != 0)
+ {
+ mp_limb_t cy;
+ cy = mpn_mul_1 (tp, rp, rsize, blimb);
+ if (cy != 0)
+ {
+ tp[rsize] = cy;
+ rsize++;
+ }
+ xp = tp; tp = rp; rp = xp;
+ }
+ }
+
+ /* Now then we know the exact space requirements, reallocate if
+ necessary. */
+ if (r->_mp_alloc < rsize)
+ _mpz_realloc (r, rsize);
+
+ MPN_COPY (r->_mp_d, rp, rsize);
+ r->_mp_size = rsize;
+ TMP_FREE (marker);
+}
projects / ghc-hetmet.git / commitdiff