On 20.04.2017 20:29, H. S. Teoh via Digitalmars-d wrote:
On Thu, Apr 20, 2017 at 02:51:12PM +0200, Timon Gehr via Digitalmars-d wrote:
On 20.04.2017 03:00, H. S. Teoh via Digitalmars-d wrote:
On Thu, Apr 20, 2017 at 02:01:20AM +0200, Timon Gehr via Digitalmars-d wrote:
[...]
Yes, there is in fact a beautifully simple way to do better. :)
...

Ahh, so *that's* what it's all about. I figured that's what I was
missing. :-D  Thanks, I'll include this in QRat soon.

Update: QRat now supports ^^. :-) Integral exponents only, of course.  I
also implemented negative exponents, because QRat supports division and
the same algorithm can be easily reused for that purpose.
...

Nice! :)

...

[...]
BTW, you are right that with std.bigint, computation using a linear
number of additions is actually faster for my example (100000th
Fibonacci number).  The asymptotic running time of the version with
pow on QRats is lower though, so there ought to be a crossover point.
(It is Θ(n^2) vs.  O(n^log₂(3)·log(n)). std.bigint does not implement
anything that is asymptotically faster than Karatsuba.)

Yeah, probably there is a crossover point. But it might be quite large.
I suppose one could make a graph of the running times for increasing n,
and either find the crossover point that way or extrapolate using the
known curve shapes.

Having said that, I haven't scrutinized the performance characteristics
of QRat too carefully just yet -- there is probably room for
optimization.


Gcd is the problem. The following code which implements a strategy based on matrix multiplication instead of QRat multiplication is significantly faster than naive linear computation:

BigInt fib(long n){
BigInt[2] a=[BigInt(0),BigInt(1)],b=[BigInt(1),BigInt(2)],c=[BigInt(1),BigInt(-1)];
     for(;n;n>>=1){
            foreach(i;1-n&1..2){
            auto d=a[i]*a[1];
            a[i]=a[i]*b[1]+c[i]*a[1];
            b[i]=b[i]*b[1]-d;
            c[i]=c[i]*c[1]-d;
        }
    }
    return a[0];
}

If I change the gcd computation in QRat (line 233) from

auto g = gcd(abs(a), abs(b), c);

to

auto g = gcd(abs(a), c, abs(b));

I get performance that is a lot closer to the matrix version and also beats the linear computation. (This is because if one of the operands is 1, gcd is cheap to compute.)


...

Yes, at most, except you don't need "+1". (For each radical ri, you
will at most need to pick a power between 0 to deg(ri)-1 to index into
the coefficients.)
[...]

The +1 is for the denominator, assuming integer coefficients.  Since
having 2^n rational coefficients is equivalent to having 2^n integer
coefficients (which are half the size of rational coefficients,
computer-representation-wise) + 1 denominator.
...

Ah, I see. Personally, I'm more in the one denominator per coefficient camp. :) I think having a designated ℚ type is cleaner, and it might even be more performant.

Reply via email to