[Numerical Analysis] Numerical Integration

• 수치적인 방법을 이용하여 적분을 수행하는 과정에 대한 포스팅이다.
• 'Numerical Methods and Optimizaiton : An Introduction (Sergiy Butenko)' - chapter 6 를 참고하였다
• 본 포스팅은 polynomial interpolation 학습을 가정한다. 생소하다면 아래의 포스팅을 참고하자
• 2021.07.06 - [[컴퓨터] 전산생물학/Modeling & Simulation] - [Numerical Analysis] Polynomial Interpolation

[Numerical Analysis] Polynomial Interpolation

함수 $f(x)$가 복잡하여 특별한 형태의 form이 존재하지 않거나, analytical 하게 풀 수 없는 경우 수치적인 방법을 통해 적분을 수행한다. 

$$\int_a^b f(x)\ dx$$

 

이러한 경우 $f(x)$를 적분이 상대적으로 효율적인 interpolating polynomial ($p_n(x)$) 형태로 변형해 준 뒤, 이를 적분하는 방법을 사용한다. 이때 lagrange method를 사용한 polynomial 형태는 아래와 같다. 생소하다면 위의 Polynomial interpolation 포스팅을 참고하자.

$$p_n(x) = \sum_{k=0}^n f(x_k)l_k(x)$$

 

하지만 interpolating polynomial 형태는 이미 approximate 한 형태이기 때문에 error term $(e_n)$을 가지게 된다.

$$f(x) = p_n(x) + e_n(x)$$

 

이를 적분하면 아래와 같다

$$\int_a^b f(x)\ dx = \int_a^b p_n(x)\ dx + \int_a^b e_n(x)\ dx$$

 

위의 lagrange method를 이용한 식을 $p_n(x)$에 대입하여 정리하면 다음과 같은 식을 얻는다.

$$\int_a^b p_n(x)\ dx = \int_a^b \left( \sum_{i = 0}^{n} f(x_i)l_i(x) \right) \ dx = \sum_{i =0}^n f(x_i) \int_a^b l_i(x)\ dx$$

 

그렇다면 이 식의 의미는 무었일까? n=1 , n=2의 경우를 각각 생각하여 보자

 

1. Trapezoidal Rule

n=1 인 경우의 $x_0 = a, x_1 = b$일 때 $p_1(x)$을 구하면

$$p_1(x) = f(a)l_0(x)+f(b)l_1(x)$$

$$= f(a)\frac{x-b}{a-b} + f(b)\frac{x-a}{b-a}$$

이를 적분하면

$$\int_a^b f(x)\ dx \approx  f(a) \int_a^b \frac{x-b}{a-b}\ dx + f(b)\int_a^b \frac{x-a}{b-a}\ dx$$

$$= \frac{b-a}{2}(f(a)-f(b))$$

 

이를 시각적으로 표현하면 아래와 같다.

즉, (a, f(a)), (b, f(b)를 이은 평행사변형의 넓이를 구하는 것과 같음을 알 수 있다. 

 

2. Simpson's Rule

그렇다면 n=2 인 경우를 살펴보자.

$x_0 = a, x_1 = \frac{a+b}{2}, x_2 = b$  3 points 을 이용하여 $p_2(x)$구하면 아래와 같다

$$p_2(x) = f(a)\frac{(x-\frac{a+b}{2})(x-b)}{(x_0 -\frac{a+b}{2})(x_0-b)}+ f\left(\frac{a+b}{2}\right)\frac{(x-a)(x-b)}{(x_1-a)(x_1-b)}+ f(b)\frac{(x-a)(x-\frac{a+b}{2})}{(x_2-a)(x_2-\frac{a+b}{2})}$$

이를 적분하면

$$\int_a^b p_2(x)\ dx = \frac{b-a}{6}\left( f(a) + 4f\left(\frac{a+b}{2}\right) + f(b) \right)$$

이는 시각적으로 표현하면 아래와 같다.

Trapezoidal Rule 에 비하여 integration error가 상대적으로 적어졌음을 알 수 있다.

하지만 두 method 모두 상대적으로 큰 a 와 b 사이 간격을 가지는 경우 error가 매우 커지는 것을 확인해 볼 수 있다.

이를 계산해 볼 수 있는데 정리는 생략하기로 하고 최종 에러에 관한 식만 적어보도록한다.

 

trapezoida rule error = $e_T = -\frac{f''(e)}{12}(b-a)^3$

simpson's rule error = $e_S = -\frac{f^{4}(e)}{90}\left(\frac{b-a}{2}\right)^5$

 

그렇다면 a 와 b 가 큰 경우는 어떤 방법을 사용하는가? 

 

3. Composite Rule

위 의 두 방법을 통하여 적분해주려는 a b의 거리가 먼 경우 큰 에러를 가지는 경우를 알 수 있었다. 이러한 경우 [a,b] 를 N 개의 subinterval 로 나우어 준 뒤에 각각의 interval 에 대하여 trapezoidal rule과 simpson's rule을 적용하는 방법을 생각해 볼 수 있다. 

우선 N 개의 동등한 길이의 subinterval을 나누고 각각의 길이를 h라 하자

$$h = \frac{b-a}{N}$$

따라서 N-1 개의 새로은 $x_i$들이 생겨나게 된다.

$$a_0 = x_0 < x_1 < x_2 < \dots < x_n = b$$

 

이러한 subinterval들의 적분을 trapezoidal rule을 이용하여 나타내면 아래와 같다.

$$T_N(f) = \sum_{i =0}^{N-1} \frac{h}{2}(f(x_i)+f(x_{i+1})) = h\sum_{i = 1}^{N-1} f(x_i) + \frac{h}{2}(f(x_0)+f(x_N))$$

 

composite trapezoidal rule을 (N=4)인 경우를 시각적으로 표현하면 아래와 같다.

에러가 상대적으로 적어짐이 확연하게 보인다.

이는 error term이

trapezoida rule error = $e_T = -\frac{f''(e)}{12}(h)^3$로 (b-a)에 비해 h가 현저하게 작은 값을 지니기 때문에 전제 error가 낮아짐을 수학적으로 확인가능하다

 

$$\lim_{N \to \infty} e_T = 0$$

 

Composite Simpson's rule 도 같은 방법으로 생각해 볼 수 있다.

$$S_N(f) = \sum_{i = 0}^{N-1} \frac{h}{6} \left( f(x_i) + 4f\left( \frac{x_i + x_{i+1}}{2} \right) + f(x_{i+1}) \right)$$