IDL에서 적분을 해봅시다 (Part 2)

IDL에서의 수치적분 방법에 대해서 지난번에 글을 올렸던게 벌써 한 달이 넘었군요. Part 2를 이제서야 올립니다. 지난 Part 1의 말미에서 언급했던 것처럼, 함수형태를 알 수 없는 그냥 데이터 포인트들의 형태로만 갖고 있는 경우에 대한 수치적분값을 어떻게 구하느냐에 관한 내용이 되겠습니다. 사실 내용 자체는 매우 간단합니다. 바로 INT_TABULATED라는 내장함수를 사용하면 됩니다. 예를 들어, 다음과 같은 형태의 X, Y 자료값 배열을 갖고 있다고 합시다. 그리고 이 자료값 점들을 플롯해보면 다음 그림과 같습니다.

x = [0, 1, 2, 4, 6, 9, 10]

y = [3.2, 4.3, 7.1, 11.6, 12.4, 8.9, 8.2]

win = WINDOW(DIMENSIONS=[500, 500])

pl = PLOT(x, y, SYMBOL='circle', YRANGE=[0, 15], /SYM_FILLED, SYM_COLOR='crimson', LINESTYLE=2, CLIP=0, /CURRENT)

즉 이러한 자료값 점들로 이루어지는 함수를 추정하고 그 적분값을 구해야 합니다. 특히 그 "함수"를 추정하는 방법이 관건인데, INT_TABULATED 함수에서는 내부적으로 5-point Newton-Cotes Formula라고 하는 기법을 사용합니다. 이 방법은 간단하게만 언급한다면, 연속되는 5개의 점들로 이루어지는 구간을 4차식으로 근사하여 매 구간마다의 적분값을 구하고 이를 총합하는 수치적분 기법입니다. 사실 위 그림에서는 자료값 점들을 직선 형태의 점선으로 이어놓았지만, 실제 수치적분 과정에서 가정되는 내삽 함수의 형태는 좀 더 울퉁불퉁한 곡선이 될겁니다. 어쨌든 이와 같은 방법으로 수치적분값을 얻는 과정은 아래와 같이 달랑 한 줄이면 됩니다.

S_5nc = INT_TABULATED(x, y)

너무 간단해서 좀 허무한 감이 없잖아 있는데, 어쨌든 다음과 같이 적분값을 출력해보면 약 93.67이라는 값을 얻게 됩니다.

PRINT, '5-point N-C Result :', S_5nc

5-point N-C Result :      93.6688

만약에 각 구간을 위의 그림처럼 무조건 직선으로 이어서 각각의 삼각형 또는 사다리꼴 면적을 계산하여 전부 합쳐본다면, 제가 해본 바로는 그 값은 92.65로 산출이 됩니다. 앞서 INT_TABULATED를 사용하여 얻은 값과 크게 다르지는 않습니다. 물론 이 함수를 사용하더라도, 원래 자료값 구간이 더 촘촘하다면 산출되는 적분값도 꽤 정확해질 것입니다. 그래서, 어떤 자료값 점들의 형태를 우리가 아는 함수의 형태로 맞추기가 어려운 상태에서 적분값을 얻어야 할 경우라면, 오늘 소개한 INT_TABULATED라는 내장함수를 사용하는 것이 좀 더 유용한 방법일 것으로 생각됩니다.