計算誤差について

質問

こんばんは。去年までフォートランをやっていたため現在Ｃ言語に苦戦中です。以下のプログラムを書いたのですが計算誤差が出てしまいます。修正をお願いいたします。
#include <iostream.h>
#include <math.h>
int main (void)
{
double t,x,pi=4.0*atan(1.0);
for(t=0.0; t<=20.0; t+=0.1){
x=double(sin(pi/3*t));
cout<<t<<" "<<x<<endl;
}
return 0;
}
この文はｔの値を0.1ずつを変化させて正弦波を描こうと思って書いたプログラムです。t=3.0のときsin(180°)となりｘ＝0.0になるのはずなんですが、そうはならず10の-15乗くらいの誤差が出てしまいます。どうしたら誤差をなくせれるでしょうか？
また、フォートランではｔ＝0.1d0というような表現をして計算誤差をなくしていたのですがC(++)にはこういった方法はないのでしょうか？
ご教授よろしくお願いいたします。

回答No.5

wolv · ベストアンサー

piが１５桁めまでの精度までしかないのが原因で、
sin(pi) や sin(pi/6) (=1/2) が
仮に１５桁よりちょっと下の桁で誤差eがあった場合、

sin(pi)= 0 + e
sin(pi/6)=0.5 + e

になります。これが精度がdoubleの変数に入り、その精度の位置で四捨五入された場合、sin(pi)は誤差eがそのまま残りますが、sin(pi/6)は四捨五入されると0.5ぴったりになる、ということじゃないでしょうか。ちなみに、０．５は二進数でもキリのいい数字ですよ。

実際に、eがどの程度の値になるかは、
e = sin(x+δ）-sin(x)
で見積もれます。
e/δ＝（ sin(x+δ）-sin(x) )/ δ
は大体sinの微分になるので、
e/δ ＝ cos(x) 程度、つまり、
e = δ cos(x)程度でしょう。
x=piの時、cos(x)=-1なので、e=-δ程度、
x=pi/6の時、cos(x)= √3 / 2 = 1.732/2 = 0.866 なので、e= 0.87 δ程度、
よって 0.5＋e = 0.5 + 0.87 δ = 0.5 * (1+1.74δ）となり、この計算結果からだけでは(1+1.74δ）がdoubleの精度で四捨五入されて1になるかどうかは微妙であることしかしめせませんでしたが、結果がそうなのだから、そうなのでしょう。

お礼

前半は良く分かりました。でも後半は僕にはちょっと難しいようです^^;
ありがとうございました。

回答No.4

wolv

計算誤差ではなく、まるめ誤差ですよ。
sin(pi)の値を表示してみましたか？

doubleの精度で、10の-15乗程度まで０ならば、それは０ということです。たぶん。

piの値がdoubleの精度しかなければ、sin piの値もdoubleの精度までしかでません。
引数の単位がラジアンではないサイン関数を作ってください。

おまけ。
t=20のときにループの内容は実行されていますか？私の実行系では、0.1ずつ足すと、t=20.00000000000000001のような値になるので、t=約19.9の回までしか実行されません。
for(i=1;i<200;i++){
　t=i/10.0;
　/*処理*/
}
のようにしましょう。

お礼

遅くなってすみません。ご回答ありがとうございます。
御礼のところで申し訳ないのですが、なぜsin(pi/3.0*t)=0.0となるはずのところだけ誤差が出るのでしょうか？piがdoubleの制度ならば他のところも誤差が出ているはずだと思うのですが。例えばt=0.5のときはちゃんとx=0.5が表示されています。よろしければご教授ください。

回答No.3

nagare

10の-15とはしびあですね

pi/3*tで既に誤差がでていませんか？
一度定数で求めて確認してはいかがでしょうか？
また、pi/3を定数(詳しい値)にしてtを掛けたらいかがでしょうか？
もうひとつ、pi/3*tを整数で計算し、少数に直した方が誤差は減ると思います（割り算で誤差がでますので）
例
pi/3*t/100
(pai、tは100倍)
(精度上げるならば、倍率をあげる）

うまく動けばいいですね

お礼

遅くなってすみません。良い方法をありがとうございます。
ご回答ありがとうございました。

回答No.2

wogota

C言語(C++)では、浮動小数点はIEEE準拠がほとんどですので、2ビット表現による
丸め誤差が発生します。

どうしても回避しようとすると、2進化10進のライブラリを探してくるか、
作るしかないかと思います。BCDのライブラリを私は知りませんが、もしかしたら
あるかもしれませんので。

お礼

遅くなってすみません。たかが誤差にライブラリを持ってくるのは面倒ですね。
ご回答ありがとうございました。

回答No.1

sha-girl

コンピューター上では2進数であるため
誤差がでるのは仕方ありません。

処理系にもよると思いますが
例えば下記のプログラムではおそらくans=1にはならないでしょう。
int i;
double ans=0;
for (i=0;i<10000;i++){
ans += 1/10000;
}

さらに無理数が含まれればなおの事です。
誤差をなくすためには四捨五入してはどうでしょうか？

お礼

遅くなってすみません。Cには誤差をなくす方法はないんですね～
ご回答ありがとうございました。