Re: How to do fast integer division

[Max Alekseyev <maxale@gmail.com>]

I would try to use X[i,3]\X[j,3] instead of truncate(X[i,3]/X[j,3]).
Max

On Fri, Aug 24, 2012 at 3:31 PM, Dirk Laurie <dirk.laurie@gmail.com> wrote:
> I have an application that does lots of integer divisions with
> large operands but small quotients.  As a model for that,
> I handcoded the extended Euclidean algorithm as follows:
>
> {mbezout(a,b)= my(X); X=[1,0,a;0,1,b]; if(a<b, i=2;j=1, j=2;i=1);
>   while(X[j,3]>0, X[i,] -= truncate(X[i,3]/X[j,3])*X[j,];
>     i=j; j=3-i);
>   X[i,]
> }
>
> This code is pretty standard, so I was astounded to see how
> slow it is.  Timing on 10000-digit numbers:
>
> ? b1=bezout(a1,a2);
> time = 16 ms.
>
> ? b3=mybezout(a1,a2);
> time = 38,622 ms.
>
> That's a factor of over 2000 slower, not to be explained
> by overhead.
>
> Suspecting that "truncate(X[i,3]/X[j,3])" hides a multitude
> of sins (starting with reduction to lowest terms), I changed
> it to "truncate(X[i,3]*1./X[j,3])" and the time dropped down
> to 200ms.  There is a tiny chance that this is not exact, so
> I did "divrem(X[i,3],X[j,3])[1]" instead and got 221ms.
>
> Is that the closest one can get to a fast integer quotient
> in Pari/GP?
>