フィボナッチ数列の計算法(2)

前回の宿題。
べき乗計算に指数の計算を使う(powを使う)場合、80ビット浮動小数点数を使うと、フィボナッチ数列を何項目まで正しく求められるだろうか。


D言語Intelの場合、realは80ビット浮動小数点数となる。これを用いた計算と、64ビット整数で前回の方法を用いた計算と比較する。fib1が実数計算、fib2が整数計算。

fib1(81) = 37889062373143906.081
fib2(81) = 37889062373143906
fib1(82) = 61305790721611591.134
fib2(82) = 61305790721611591
fib1(83) = 99194853094755497.226
fib2(83) = 99194853094755497
fib1(84) = 160500643816367088.370
fib2(84) = 160500643816367088
fib1(85) = 259695496911122585.600
fib2(85) = 259695496911122585
fib1(86) = 420196140727489673.930
fib2(86) = 420196140727489673
fib1(87) = 679891637638612259.620
fib2(87) = 679891637638612258
fib1(88) = 1100087778366101933.600
fib2(88) = 1100087778366101931
fib1(89) = 1779979416004714193.100
fib2(89) = 1779979416004714189
fib1(90) = 2880067194370816127.000
fib2(90) = 2880067194370816120
fib1(91) = 4660046610375530320.500
fib2(91) = 4660046610375530309
fib1(92) = 7540113804746346448.000
fib2(92) = 7540113804746346429

整数計算は92項まで計算できる。
実数計算は85項で間違える。




import std.cstream;
import std.math;

void main() {
for(int i = 81; i <= 92; i++) {
dout.writefln("fib1(%d) = %.3f", i, fib1(i));
dout.writefln("fib2(%d) = %d", i, fib2(i));
}
}

real fib1(int n) {
return (pow(((1 + sqrt(5.0L)) / 2), n)
- pow(((1 - sqrt(5.0L)) / 2), n)) / sqrt(5.0L);
}

ulong fib2(int n) {
cGal c = new cGal(1, 1);
return c.pow(n).getB;
}

class cGal {
private:
ulong a;
ulong b;
public:
this(ulong a, ulong b) {
this.a = a;
this.b = b;
}
ulong getA() { return a; }
ulong getB() { return b; }
cGal dup() { return new cGal(a, b); }
cGal opMul(cGal g) {
return new cGal((a * g.a + 5 * b * g.b) / 2,
(a * g.b + b * g.a) / 2);
}
cGal pow(int n) {
cGal x = this.dup;
int m = calcFigure(n);
for(int i = m - 2; i >= 0; i--) {
x = x * x;
if((n >> i) & 1)
x = x * this;
}

return x;
}
private:
int calcFigure(int n) {
int m;
for(m = 1; n >> m != 0; m++) {

}
return m;
}
}