FC2ブログ

# 誤差関数 Error Function
# 引数 値 ($X)
# 戻り値 誤差関数 ($ErrorFunction)
sub ERRORFUNCTION{
my ($X) = @_;
my $ErrorFunction = 0;
my $RombergsMethod = 0;
my @Romberg = ();
my $Sum = 0;
my $a = 0;
my $b = $X;
my $h = 0;
my $N = 10;
my $Power4 = 0;
my $Epsilon = 1.0e-10;
my $Pi = atan2(1, 1) * 4;

if($X == 0){
# 誤差関数 Error Function
$ErrorFunction = 0;

return $ErrorFunction;
}

# ロンバーグ積分 Romberg's Method
# (0, 0)
$Romberg[0] = (($b - $a) / 2) * (&INTEGRAND($a) + &INTEGRAND($b));

for(my $i = 1; $i < $N; $i++){
$Power4 = 1;

for(my $j = $i; $j >= 0; $j--){
if($i == $j){
$Sum = 0;
$h = ($b - $a) / (2 ** $j);

for(my $k = 1; $k <= (2 ** ($i - 1)); $k++){
$Sum += &INTEGRAND($a + (((2 * $k) - 1) * $h));
}

# ($i, 0)
$Romberg[$j] = ($Romberg[$j - 1] / 2) + ($Sum * $h);
}else {
# 累乗
$Power4 = 4 * $Power4;

# リチャードソンの補外 Richardson's Extrapolation
$Romberg[$j] = (($Power4 * $Romberg[$j + 1]) - ($Romberg[$j])) / ($Power4 - 1);

# 収束判定
if(abs($Romberg[$j] - $Romberg[$j + 1]) < $Epsilon){
$Romberg[0] = $Romberg[$j];
$i = $N;
last;
}
}
}
}

# ロンバーグ積分 Romberg's Method
$RombergsMethod = $Romberg[0];

# 誤差関数 Error Function
$ErrorFunction = (2 / sqrt($Pi)) * $RombergsMethod;

return $ErrorFunction;
}

# 被積分関数 Integrand
# 引数 変数 被積分関数 ($X, \@Integrand)
# 戻り値 被積分関数 ($Integrand)
sub INTEGRAND{
my ($X) = @_;
my $Integrand = exp(-($X * $X));

return $Integrand;
}


参考URL
誤差関数 - Wikipedia
誤差関数 - 高精度計算サイト
オンライン コンパイラ/インタプリタ
テクニカル分析
プロフィール

Author:雨宮
Firefoxを使用しているので気づかなかったけど、IE6でソースコードを上手くコピーできない

5/3
携帯用ならIE6でもソースコードをコピーできる
携帯用

検索フォーム


あわせて読みたいブログパーツ
一寸先は闇 RSS