bdiv vs redc

[Niels Möller]

Torbjorn Granlund <tg at gmplib.org> writes:

> Could you perhaps write a new proposed sbpi1_bdiv_qr with the result
> normalisation you suggest, using my proposed sbpi1_bdiv_r style?  That
> would allow for a single outer loop, unlike the current sbpi1_bdiv_qr.

Below. I'm also attaching a larger patch which also updates
sbpi1_bdiv_q, the dc code, the t-bdiv testcode, and some of the bdiv
callers (divexact and remove). I haven't touched the mu_bdiv code. I get
lots of failures, I guess some but not all are caused by mu_bdiv
following the old rather than the new convention (I tried setting the
MU_BDIV thresholds to infinity, and I still have lots of failures.).

Regards,
/Niels

/* mpn_sbpi1_bdiv_qr -- schoolbook Hensel division with precomputed inverse,
   returning quotient and remainder.

   Contributed to the GNU project by Niels Möller.

   THE FUNCTIONS IN THIS FILE ARE INTERNAL FUNCTIONS WITH MUTABLE INTERFACES.
   IT IS ONLY SAFE TO REACH THEM THROUGH DOCUMENTED INTERFACES.  IN FACT, IT IS
   ALMOST GUARANTEED THAT THEY'LL CHANGE OR DISAPPEAR IN A FUTURE GMP RELEASE.

Copyright 2006, 2009, 2011 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 Lesser General Public License as published by
the Free Software Foundation; either version 3 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 Lesser General Public
License for more details.

You should have received a copy of the GNU Lesser General Public License
along with the GNU MP Library.  If not, see http://www.gnu.org/licenses/.  */

#include "gmp.h"
#include "gmp-impl.h"


/* Computes a binary quotient of size qn = un - dn.
   Output:

      Q = -U * D^{-1} mod B^qn,

      R = (U + Q * D) * B^(-qn)

   Stores the dn least significant limbs of R at {up + un - dn, dn},
   and returns the carry from the addition N + Q*D.

   D must be odd. dinv is (-D)^-1 mod B. */

mp_limb_t
mpn_sbpi1_bdiv_qr (mp_ptr qp,
		   mp_ptr up, mp_size_t un,
		   mp_srcptr dp, mp_size_t dn, mp_limb_t dinv)
{
  mp_size_t i;
  mp_limb_t cy;

  ASSERT (dn > 0);
  ASSERT (un > dn);
  ASSERT ((dp[0] & 1) != 0);
  ASSERT (up == qp || !MPN_OVERLAP_P (up, un, qp, un - dn));

  for (i = un - dn, cy = 0; i != 0; i--)
    {
      mp_limb_t q = dinv * up[0];
      mp_limb_t hi = mpn_addmul_1 (up, dp, dn, q);
      *qp++ = q;
      
      hi += cy;
      cy = hi < cy;
      hi += up[dn];
      cy += hi < up[dn];
      up[dn] = hi;
      up++;
    }

  return cy;
}

-------------- next part --------------
A non-text attachment was scrubbed...
Name: bdiv.patch
Type: text/x-patch
Size: 12157 bytes
Desc: not available
URL: <http://gmplib.org/list-archives/gmp-devel/attachments/20120703/a4c261fa/attachment.bin>
-------------- next part --------------

-- 
Niels Möller. PGP-encrypted email is preferred. Keyid C0B98E26.
Internet email is subject to wholesale government surveillance.