divappr interface

[paul zimmermann]

       Niels,

> From: nisse at lysator.liu.se (Niels Möller)
> Date: Fri, 27 Apr 2018 22:28:39 +0200
> 
> nisse at lysator.liu.se (Niels Möller) writes:
> 
> > Once we reach qn = dn - 1, keep looping to produce quotient limbs, but
> > also discard one limb of dp in each interation, until we in the final
> > iteration have qn = 1, qn+2 = 3 numerator limbs, and qn+1 = 2 divisor
> > limbs, i.e., a single udiv_qr_3by2 (where we might consider skipping the
> > adjustment steps). I haven't done the error analysis, but I would expect
> > that errors are similar to the current code.
> 
> In the current mpn_sbpi1_divappr_q, there's a curious flag (or mask)
> variable used in this loop. It's initially all ones, cleared if we ever
> fail to cancel the top limb of the current partial remainder. When the
> flag is cleared, remaining quotient limbs are set to GMP_NUMB_MAX.
> 
> Torbjörn, Paul, do you remember the analysis behind this? (Code is since
> 2009...).
> 
> I would guess it might happen that even if higher quotient limbs are all
> correct, we might get non-zero high limb in
> 
>   partial remainder - GMP_NUMB_MAX * truncated D 
> 
> If we set the rest of the quotient limbs to GMP_NUMB_MAX, then the
> quotient should be large enough thanks to the low end divisor limbs
> which we ignored in the truncation.
> 
> Regards,
> /Niels

I don't believe I have written that code. Anyway I don't quite understand
the code from sbpi1_divappr_q.c (say in GMP 6.1.2):

* when flag becomes 0, the condition UNLIKELY (n1 >= (d1 & flag)) is always
  true, thus we always take that branch. Is that wanted?

* when n1 != cy on line 129, I guess that q = GMP_NUMB_MASK was too large,
  thus we should decrease q and add back the divisor. But the code in lines
  133-134 is never executed (https://gmplib.org/devel/lcov/shell/var/tmp/lcov/gmp/mpn/sbpi1_divappr_q.c.gcov.html)

* idem for the code in lines 177-178

Paul