この計算は、浮動小数点数だと誤差がたまってくるが、幸い2進数では有限小数になるので、ある程度の桁までは正確に計算できる。
Arrayで、2048進数で考えるとプログラミングがやさしくなる。
例えば、3 + 11/2048 + 13/2048^2 → [ 3, 11, 13 ]とする。
10万項まで計算したところ、
[ 2045, 2047, 2047, 2047, 2047, 2046, 439, 386, 1912, 1451, ... ]
すなわち、2046に4.96e-17程度足りない数値となった。
p(100000,10)は2.4e-25程度だったようだ。
var n = 10;
var N = 100000;
var base = 1 << (n + 1);
var p = [ ];var E = [ 0 ]; // 期待値
p[n] = [ 0, 2 ];
update_E(E, n, p[n]);
for(var m = n + 1; m <= n * 2; m++) {
p[m] = [ 0, 1 ];
update_E(E, m, p[m]);
}var sum = [ 0, 1 ];
for(var m = n * 2 + 1; m <= N; m++) {
update_sum(sum, p[m-n-1]);
p[m-n-1] = null; // 捨て
p[m] = sum.slice(); // コピー
update_E(E, m, p[m]);
}function update_sum(a, b) {
for(var i = b.length - 1; i > 0; i--) {
var tmp = a[i+1];
if(tmp == undefined)
tmp = 0;
tmp -= b[i];
if(tmp < 0) {
a[i]--;
a[i+1] = tmp + base;
}
else {
a[i+1] = tmp;
}
}
}function update_E(E, m, p) {
for(var i = p.length - 1; i > 0; i--) {
var tmp = E[i];
if(tmp == undefined)
tmp = 0;
tmp += m * p[i];
var a = Math.floor(tmp / base);
E[i-1] += a;
E[i] = tmp - a * base;
}
}