Calculate numerical integration of normal distribution

ガウス積分の数値計算.md

概要

以下の積分値を数値計算により導出する

手法

x = (1 - t) / tという変数変換を行うことで、下式を得る。

ここで、

とおくと、区分求積法により、Iの値は以下のように近似できる。

ここで、

にであることより（証明略）、

これを求める言語のプログラムは以下のようになる。

#include <stdio.h>
#include <math.h>
#include <stdlib.h>

// 区分求積法の分割数
#define N 1000

// 変数変換を行った被積分関数
double g(double x) {
    double x_squared = x * x;
    double one_minus_x = 1 - x;
    return exp(-(one_minus_x * one_minus_x) / (2.0 * x_squared)) / x_squared;
}

int main() {
    double sum = 0.5;
    for (int k = 1; k < N; ++k) {
        sum += g((double)k / N);
    }

    printf("%lf\n", sum / N);
    return 0;
}

実行結果

$ ./a.out
1.252312