mpn_invert implementing ApproximateReciprocal.

[Niels Möller]

bodrato at mail.dm.unipi.it writes:

> I'm not sure I fully understand the interface of divappr_q,

Me neither... Documentation could be improved. I think it uses 3/2
division, and hence the inverse should be computed using invert_pi1.

So

  mpn_tdiv_qr (rp, ip - rn, 0, xp, 2 * rn, dp - rn, rn);

should be replaced (assuming inputs are normalized!) by something like

  invert_pi1 (inv, dp[rn-1], dp[rn-2]);
  if (BELOW_THRESHOLD (rn, ...))
    mpn_sbpi1_divappr_q (ip - rn, xp, 2*rn, dp-rn, &inv, rn);
  else
    {
      tp = some scratch space;
      mpn_dcpi1_divappr_q_n (ip - rn, xp, dp-rn, rn, &inv, tp);
    }
  MPN_DECR_U (ip - rn, rn, 1);

> moreover the mpn_sbpi1_ flavour begins with ASSERT (dn > 2)... which means
> we need a basecase for dn=1 or 2

dn = 1 is invert_limb. The best way to handle dn = 2 is somewhat
unclear, one possible way to do it is to first invert the most
significant limb using invert_limb, and then use that for a naive 4/2
schoolbook division using udiv_qrnnd_preinv. (And then I'm not sure what
to do with normalization, if mpn_invert doesn't require dp[dn-1] to have
the most significant bit set). I don't think these cases are terribly
important for current applications of mpn_invert.

Regarding divappr (and for simplicity consider the balanced case, 2n/n),
I'm not sure why it needs inputs of size 2n resp. n, shouldn't it be
sufficient with 2n+1 and n, totally ignoring least significant limbs of
tha larger input? It matters for mpn_invert, since then it wouldn't have
to fill in (and allocate) all limbs of xp.

In principle, couldn't divappr u/v (u of size 2n, v of size n) be
implemented by computing an n+1 limb inverse of v, and then that and
the most significant n+1 limbs of u? Multiplication could use
mulhi_appr. I.e., both the low limbs of u, and all low limb products in
the multiplication, are ignored.

Regards,
/Niels

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