Fibonacci数 の一般項

一般項

F_{n} = \frac{1}{\sqrt{5}} ( ( \frac{1 + \sqrt{5}}{2} )^n - ( \frac{1 - \sqrt{5}}{2} )^n )

が、Fibonacci数を求めるプログラムにどの程度使えるか検証してみる。
無理数絡みなので、どのあたりで誤差が出るのかが興味の対象。

#include <iostream>
#include <cmath>
using namespace std;

const double R5 = sqrt(5.0);

//typedef int num_t;
typedef long long num_t;

num_t fib(int n) {
  double f = (pow((1+R5)/2, n) - pow((1-R5)/2, n)) / R5;
  return llround(f); // 四捨五入
}

void fibtest(num_t a, num_t b, int n) {
  num_t c = a + b;
  if (c < b) {
    cout << "overflow:  n = " << n << endl;
  } else if (c != fib(n)) {
    cout << "wrong answer:  n = " << n << endl
         << fib(n) << endl << c << endl;
  } else {
    fibtest(b, c, n+1);
  }
}

int main() {
  fibtest(0, 1, 2);
}

結果

# num_t = int
overflow: n = 47


# num_t = long long
wrong answer: n = 71
308061521170130
308061521170129

32bitの範囲なら大丈夫なものの、64bit程度の大きさでもうアウト。
ちなみに四捨五入ではなく単に切り捨てにすると、n = 72 でアウトという結果になりました。
あと割り算を逆数の掛け算にするとか試してみましたがあまり効果無いみたいです。