dead code in div_q.c?

Niels Möller nisse at lysator.liu.se
Thu Apr 26 20:12:09 UTC 2018


paul zimmermann <Paul.Zimmermann at inria.fr> writes:

> In the second case (lines 267-276) this is less obvious since you use an
> approximate division (divappr), and I don't know what is the maximal allowed
> difference between divappr and div. If divappr always return a quotient
> less or equal to the true quotient, then qh cannot be zero in this case.
> But if qh can be larger than the true quotient, even by one bit, then it
> is possible that qh is 1, in the case {new_np, new_nn} = 0111...111 and
> {new_dp, dn} = 1000...000.

You're right again... I thought we had more margin.

According to comments, quotient from sb and dc divappr is correct or
one too large, while quotient from mu divappr can be at most 4 too
large.

And the limit case, N = B^n/2 - 1, D = B^k/2 is 

   floor ((B^n/2 - 1) / (B^k/2)) = floor(B^{n-k} - B^{-k}/2) 
                                 = B^{n-k} - 1 

So even we only a single unit of error, we could get qh == 1.

I wonder what actual error can be produced by the divappr functions in
the corner cases that error could spill over into qh?

In the cy != 0 case, it looks like we we make a balanced divappr call,
(2qn+2) / (qn+1). Lets put k = qn+1. For divappr we then have

  N <= B^{2k}/2 - 1, 
  D >= B^k/2, lets say D = B/2 + e, e >= 0

To have any possibility to get qh > 0 with from the approximation error,
we must have Q >= B^k - 4, where Q denotes the correct quotient. So
assume that is the case. Now,

  D <= N/Q <= (B^{2k}/2 - 1) / (B^k - 4)

which implies

  e <= N/Q - B^k/2 <= (2 B^k - 1) / (B^k - 4) < 3

So we need to consider only D = B^k/2 + {0, 1, 2}. Then, any inverse not
taking the very least significant D limb into account (including mu inverse)
must be *all* ones, right?

It would be good to either prove that qh > 0 can't happen, or add a test
that will exercise the code handling that case.

> beware that in the case cy <> 0, one divides {new_np, nn+1} by {new_dp, dn}
> thus filling {qp, nn-dn+1}. 

Thanks for pointing this out. Now I understand better.

Regards,
/Niels

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


More information about the gmp-devel mailing list