Rはなぜここで合計を間違って実行するのですか？

Question

（答えはRが小さい小数を格納するためのバイナリ数値をどのように使用するかを実際にある場合は、まだこれは、この上の知識を持つ人からなぜ起こるかについての詳細を聞くのが大好きだ）

私は具体的には、合計を計算するためにRを使用していますこの1：

ここに私のRコードだ：

r = 21 #Number of times the loop executes 
n = 52 
sum = 0 #Running total 
for(k in 0:r){ 
    sum = sum + (1-(choose(2^(k), n)*(factorial(n)/(2^(n*k))))) 
    print(sum) 
} 

あなたは出力を見れば、あなたは気づくでしょう：

[1] 1 
[1] 2 
... 
[1] 11.71419 
[1] 11.71923 
[1] 11.72176 
[1] 12.72176 
[1] 13.72176 

なぜ19回目の反復後に1ずつ増加し始めますか？

他の自由に利用できる計算エンジンがこのタスクに適していますか？

Answer 1

すべての計算で浮動小数点数を使用していますが、これは一般に階乗や冪乗にはあまり適していません。

> factorial(52)/(2^(52*19)) 
[1] 3.083278e-230 
> factorial(52)/(2^(52*20)) 
[1] 0 

をし、巴里-GPと比較：

、試してみてください

? \p 256 
realprecision = 269 significant digits (256 digits displayed) 
? precision((52!)/2^(52*19) + . , 24) 
%1 = 3.08327794368108826958435659488643289724 E-230 
? precision((52!)/2^(52*20) + . , 24) 
%2 = 6.8462523287873017654100595727727116496 E-246 

Answer 2

はあなたが抱えているオーバーフロー/アンダーフロー問題を回避する方法を探しているなら、あなたは（中間計算のための合理的な大きさを保つために）ログを使用して、最後に累乗することができます

options(digits=20) 

for(k in 0:r){ 
    sum = sum + (1 - exp(lchoose(2^k, n) + log(factorial(n)) - k*n*log(2))) 
    print(paste0(k,": ",sum)) 
} 

[1] "0: 1" 
[1] "1: 2" 
[1] "2: 3" 
... 
[1] "19: 11.7217600143979" 
[1] "20: 11.7230238079842" 
[1] "21: 11.7236558993777" 

これが正しいことを確認するために、私はMathematicaで（ログをとらずに）元々の集計を実行し、小数点以下12桁まで同じ結果を出しました。

Rで問題を解決できますが、コンピュータ代数システム（記号計算と正確な計算を行うことができます）を使いたい場合は、Sageは無料でオープンソースです。