p2_erf_erfc

誤差関数

モジュールのインポートや初期設定など

import matplotlib.pyplot as plt  
import numpy as np
import scipy as sp
from scipy.special import erf, erfc # 誤差関数，余誤差関数
from scipy.integrate import quad
np.set_printoptions(precision=3)

import scienceplots
#Warning : As of version 2.0.0, you need to add import scienceplots before setting the style (plt.style.use('science')).
plt.style.use(['science', 'notebook'])

# rcPramsを使ってTeXのフォントをデフォルトにする
#plt.rcParams['font.family'] = 'Times New Roman' # font familyの設定
plt.rcParams['font.family'] = 'serif' # font familyの設定
plt.rcParams['mathtext.fontset'] = 'cm' # math fontの設定

誤差関数および余誤差関数

　誤差関数 $\mathrm{erf} (x)$ ，余誤差関数 $\mathrm{erfc} (x)$ は， \begin{gather} \mathrm{erf} (x) = \frac{2}{\sqrt{\pi}} \int_0^x e^{-t^2} dt, \ \ \ \ \ \mathrm{erfc} (x) = \frac{2}{\sqrt{\pi}} \int_x^{\infty} e^{-t^2} dt \end{gather}

uu = np.linspace(0.0,4.0,9) 
print(f"{'u':>3s} {'erf(u)':>17s} {'erfc(u)':>17s} {'erf(u)+erfc(u)':>17s}")
for j,u in enumerate(uu):
    f_erf = erf(u) # 誤差関数
    f_erfc = erfc(u) # 余誤差関数
    print(f"{u:>.1f} {f_erf:>.15f} {f_erfc:>.15f} {f_erf+f_erfc:>.15f}")

  u            erf(u)           erfc(u)    erf(u)+erfc(u)
0.0 0.000000000000000 1.000000000000000 1.000000000000000
0.5 0.520499877813047 0.479500122186953 1.000000000000000
1.0 0.842700792949715 0.157299207050285 1.000000000000000
1.5 0.966105146475311 0.033894853524689 1.000000000000000
2.0 0.995322265018953 0.004677734981047 1.000000000000000
2.5 0.999593047982555 0.000406952017445 1.000000000000000
3.0 0.999977909503001 0.000022090496999 1.000000000000000
3.5 0.999999256901628 0.000000743098372 1.000000000000000
4.0 0.999999984582742 0.000000015417258 1.000000000000000

数値積分によって求めると，

def fe(x):
    return np.exp(-x**2)

uu = np.linspace(0.0,4.0,9) 
print(f"{'u':>3s} {'quad:erf(u)':>17s} {'quad:erfc(u)':>17s} {'erf(u)+erfc(u)':>17s}")
for j,u in enumerate(uu):
    q_erf = quad(fe,0,u)[0]*2/np.sqrt(np.pi) # 誤差関数
    q_erfc = quad(fe,u,np.inf)[0]*2/np.sqrt(np.pi) # 余誤差関数
    print(f"{u:>.1f} {q_erf:>.15f} {q_erfc:>.15f} {q_erf+q_erfc:>.15f}")

  u       quad:erf(u)      quad:erfc(u)    erf(u)+erfc(u)
0.0 0.000000000000000 1.000000000000000 1.000000000000000
0.5 0.520499877813047 0.479500122186954 1.000000000000000
1.0 0.842700792949715 0.157299207050285 1.000000000000000
1.5 0.966105146475311 0.033894853524689 1.000000000000000
2.0 0.995322265018953 0.004677734981047 1.000000000000000
2.5 0.999593047982555 0.000406952017438 0.999999999999993
3.0 0.999977909503002 0.000022090496999 1.000000000000000
3.5 0.999999256901628 0.000000743098362 0.999999999999990
4.0 0.999999984582742 0.000000015417258 1.000000000000000

積分範囲の上限が無限大の場合， \begin{gather} \int_0^{\infty} e^{-t^2} dt = \frac{\sqrt{\pi}}{2} \end{gather} あるいは， \begin{gather} \frac{2}{\sqrt{\pi}} \int_0^{\infty} e^{-t^2} dt = 1 \end{gather}

quad(fe, -np.inf, np.inf)[0]*2/np.sqrt(np.pi)/2

1.0

被積分関数

def fe(x):
    return np.exp(-x**2)
x = np.linspace(-4, 4, 201)
y = fe(x)
fig = plt.figure() # グラフ領域
plt.grid()
plt.plot(x, y)
fig.savefig('erf_integrand.pdf')
plt.show()