Niels, > From: [email protected] (Niels Möller) > Date: Fri, 27 Apr 2018 22:28:39 +0200 > > [email protected] (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 _______________________________________________ gmp-devel mailing list [email protected] https://gmplib.org/mailman/listinfo/gmp-devel
