バーゼル問題と数値計算

はじめに

この記事はBasel問題

\sum_{n \in \mathbb{N}_+} \frac{1}{n^2}

を計算する問題を数値計算でやってみるというものです。
よく知られているようにこの問題の答えは$\pi^2/6$です。

解法1: 素直にとく

数値計算でこの問題を解く場合、当然無限大までの和を考えることができませんから、どこかで打ち切ったもの、つまりある程度大きな数$N$を定義して、

S_n := \sum_{n = 1}^{N} \frac{1}{n^2}

という数を計算することになります。例えば、C++ではこのように計算できます。

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

const double target = M_PI * M_PI / 6;

void calc1(long n_max)
{
  double sum = 0;
  cout << "start calc1" << endl;
  for (long i = 1; i <= n_max; ++i){
    sum += 1.0 / (i * i);
  }
  cout << "result: " << fixed << setprecision(15) << sum << endl;
  cout << "diff: " << scientific << target - sum << endl;

}

int main () 
{
  const long n_max = 1000000000; // 10^9
  cout << "target: " << setprecision(15) << target << endl;
  calc1(n_max);
  return 0;
}

これを計算すると、

target: 1.64493406684823
start calc1
result: 1.644934057834575
diff: 9.013651380840315e-09

このように計算できるはずです。

方法２: ちょっと工夫する

少し工夫して、逆順に計算してみましょう。

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

const double target = M_PI * M_PI / 6;

void calc1(long n_max)
{
  double sum = 0;
  cout << "start calc1" << endl;
  for (long i = 1; i <= n_max; ++i){
    sum += 1.0 / (i * i);
  }
  cout << "result: " << fixed << setprecision(15) << sum << endl;
  cout << "diff: " << scientific << target - sum << endl;

}

void calc2(long n_max)
{
  double sum = 0;
  cout << "start calc2" << endl;
  for (long i = n_max; i >= 1; --i){
    sum += 1.0 / (i * i);
  }
  cout << "result: " << fixed << setprecision(15) << sum << endl;
  cout << "diff: " << scientific << target - sum << endl;
}

int main () 
{
  const long n_max = 1000000000;
  cout << "target: " << setprecision(15) << target << endl;
  calc1(n_max);
  calc2(n_max);
  return 0;
}

これを実行すると

target: 1.64493406684823
start calc1
result: 1.644934057834575
diff: 9.013651380840315e-09
start calc2
result: 1.644934065848226
diff: 1.000000082740371e-09

となります。
１桁程度精度が改善しています。

まとめ

これが起こる理由は、大きさの極端に異なる数を足し合わせることによる情報落ちによるものです。

大きい方から足し合わせた方法１の場合、足し合わせる量の大きさの差が徐々に大きくなっていくことになり、情報落ちが起きやすくなっていました。

一方、方法２では、小さいものから徐々に足したため、足す数同士の大きさは比較的同じぐらいになり、情報落ちが緩和されました。

情報落ちは、常に大きさが近いものを足し合わせることである程度回避できるため、例えば、

1. 足し合わせるリストをソートし
2. もっとも小さい2つを足し
3. 足し合わせたものを大きさ順で挿入し
4. 2.と3.を繰り返す

という計算を行うともう少し精度が改善します。

計算順序も数値計算には必要なこともあります。¹

1. とはいえ、数学的にはどちらも同じものです。順序のない世界に順序を導入しないといけないということに数値解析をやり始めた時は混乱しました。 ↩