IDL/Math

INT_2D를 이용한 2-변수(Bivariate) 함수의 적분

이상우_idl 2020. 4. 2. 15:45
728x90
반응형

IDL에서 적분 계산을 하는 방법에 관해서는 제가 예전에 관련 게시물을 올렸던 적이 있습니다. 여기서 언급했던 내용은 1-변수 함수, 즉 F(x)의 형태로 존재하는 수식에 대한 적분값을 계산하는데 있어서 IDL의 관련 내장 함수인 QROMB 또는 QSIMP를 사용하면 된다는 내용이었습니다. 그런데 오늘 언급하고자 하는 것은 2-변수 함수, 즉 F(x, y)의 형태를 갖는 수식에 대한 2차원 적분값을 IDL에서 계산하는 방법에 관한 것입니다. 이러한 수식을 도식적으로 나타낸다면 XYZ 공간상에서서 존재하는 서피스(Surface)과 같은 형태를 갖습니다. 그리고 적분의 구간은 X축 방향 및 Y축 방향에 대하여 각각 설정됩니다. 예를 들면 다음 그림과 같은 상황일 것입니다.



물론 IDL에서는 이러한 경우를 위한 내장함수가 제공되고 있습니다. 바로 INT_2D라는 함수인데요. 그래서 이 INT_2D 함수를 이용하여 적분값을 계산하는 예제를 한번 살펴보겠습니다. 사실 INT_2D 함수가 사용 방법이 약간 독특해서 도움말의 설명만으로 이해가 쉽지 않은 측면도 있긴 한데, 제가 아는 범위 내에서 최대한 쉽게 소개해보도록 노력하겠습니다. 우선 예제로 사용한 F(x, y) 형태의 수식을 하나 가정해야 합니다. 일단 다음과 같이 매우 간단한 형태의 수식을 가정해보기로 하겠습니다.


F(x, y) = 2*x


그리고 이 수식에 대한 2차원 적분에 있어서 적분 구간은 X 방향으로는 0~2 그리고 Y 방향으로는 0~4로 가정합니다. 이 상황을 도식적으로 나타내본다면 다음 그림과 같습니다.



사실 이와 같은 경우에는 2차원 적분값은 결국 X, Y, Z 방향의 길이가 각각 2, 4, 4인 육면체의 부피의 절반에 해당됩니다. 그 값은 16이므로 결국 2차원 적분값 역시 16이 나올 것이라는 예상을 충분히 해볼 수 있을 것입니다. 그러면 이러한 결과를 얻기 위한 IDL 코딩이 필요한데요. 이 코딩에서 가장 먼저 필요한 것은 적분의 대상이 될 위의 수식을 함수형 부프로그램으로 나타내는 것입니다. 그 내용은 예를 들면 다음과 같습니다.


FUNCTION Fxy, x, y

  RETURN, 2*x

END


그리고 또 필요한 것이 있는데, 바로 Y축 방향의 적분 구간을 나타내는 함수형 부프로그램을 별도로 만들어줘야 한다는 것입니다. 앞서 언급했듯이 Y 방향의 적분 구간은 0~4이기 때문에 다음과 같은 방식으로 코딩을 합니다.


FUNCTION PQ_Limits, x

  RETURN, [0.0, 4.0]

END


그런데 Y축 방향의 적분 구간을 이렇게 굳이 함수형 부프로그램의 형태로 정의해야 한다는 것이 좀 특이하게 보일 수 있습니다. 그것은 2-변수 적분의 특성 때문인데요. 2차원적인 적분이라는 그 특성상 적분의 '구간'은 기하학적으로는 XY 평면 상에서 존재하는 마치 사각형(Rectangle)의 면과 같은 형태를 띄게 됩니다. 그래서 적분 '구간'이라기보다는 적분 '면'이라고 부르는 것이 더 나을지도 모르겠습니다. 지금 우리가 다루고 있는 예제에서 가정한 사각형 형태의 적분 '면'을 도식적으로 나타낸다면 다음 그림과 같을 것입니다.



하지만 2차원 적분에서 어떤 경우에는 적분 '면'이 사각형이 아닌 경우도 있습니다. 예를 들면 Y 방향 구간이 X 값에 따라 변동되는 경우도 있다는 것입니다. 이 경우에는 적분 '면'의 형태는 사각형외의 다양한 형태들이 될 수도 있습니다. 이렇게 적분 '면'의 형태에 대한 일반적인 경우들까지 다룰 수 있으려면 Y 방향 구간을 X에 대한 함수로 나타내야 할 필요가 있기 때문에 굳이 위와 같은 문법을 따른다 정도로 이해하시면 될 것 같습니다. 어쨌든 INT_2D를 사용하려면 위와 같은 방식의 준비가 먼저 필요합니다. 이제는 주프로그램에 해당되는 내용에 대한 코딩을 해줘야 하는데, 먼저 X축 방향의 적분 구간을 다음과 같이 정의합니다.


AB_Limits = [0.0, 2.0]


여기까지 하면 INT_2D를 사용할 준비는 다 끝난 것입니다. 여기서 잠시 INT_2D 함수의 사용 문법을 IDL 도움말에서 참조해보면 다음과같습니다.


Result = INT_2DFxyAB_LimitsPQ_LimitsPts [, /DOUBLE] [, /ORDER] )


여기서 첫번째 인수로는 적분의 대상이 되는 2-변수 수식에 해당되는 함수형 부프로그램의 이름이 문자값으로 제공되어야 합니다. 그리고 두번째 인수로는 X 방향의 적분 구간을 나타내는 두 값을 묶어서 배열의 형태로 제공하면 됩니다. 세번째 인수로는 Y 방향 적분 구간을 나타내는 함수형 부프로그램의 이름이 문자값으로 제공되어야 합니다. 그리고 네번째 인수에 대해서는 도움말에서 다음과 같이 설명되어 있습니다.


The number of transformation points used in the computation. Possible values are: 6, 10, 20, 48, or 96.


저도 정확한 의미는 완전히 파악을 못하였지만 아마도 숫자가 높을수록 더 정밀한 계산 결과를 산출한다는 의미로 봐도 무방할 것 같습니다. 따라서 이왕이면 가장 높은 값인 96을 사용하는게 나을 것 같습니다. 어차피 계산 시간이 눈에 띄게 더 걸리는 것도 아닙니다. 따라서 앞서 코딩했던 내용을 바탕으로 하면 다음과 같이 INT_2D 함수를 사용하면 됩니다.


result = INT_2D('Fxy', AB_Limits, 'PQ_Limits', 96)

PRINT, 'Result :', result


그러면 지금까지의 코딩 내용을 모두 모아서 정리를 해보겠습니다. 사실 오늘의 예제에서는 함수형 부프로그램들이 사용되어야 하기 때문에 위의 과정들을 하나의 프로그램에 담는데 있어서 그 순서에 유의해야 합니다. 만약 이 모든 과정을 test_int_2d라는 이름의 프로그램에 담는다고 한다면, 이 프로그램에는 위의 전 과정이 다음과 같은 형태로 순서대로 작성되어야 합니다.


FUNCTION Fxy, x, y

  RETURN2*x

END


FUNCTION PQ_Limits, x

  RETURN, [0.04.0]

END


PRO test_ind_2d


AB_Limits = [0.02.0]


result = INT_2D('Fxy', AB_Limits, 'PQ_Limits'96)

PRINT'Result :', result


END


이것은 IDL 프로그램이 단일 형태가 아니라 주프로그램/부프로그램이 혼재된 형태일 경우, 부프로그램들이 앞쪽에 배치되고 주프로그램이 맨 뒤에 배치되어야 하기 때문입니다. 그래서 이 내용을 test_2d_int.pro라는 프로그램 파일로 저장하여 컴파일 및 실행을 하면 온전하게 결과를 얻을 수 있습니다. 그러면 다음과 같이 16.0이라는 결과값이 출력되는 것을 확인할 수 있습니다.


Result :      16.0000


위에서 봤듯이 코딩 문법에 있어서의 약간의 독특함만 감안한다면 사용법 자체는 그리 복잡하진 않습니다. 그래서 다양한 경우들에 대해서 비교적 쉽게 시도를 해볼 수 있습니다. 예를 들어서 적분 구간은 그대로 두고 수식만 다음과 같이 약간 더 복잡한 형태로 바꿔봅시다.


FUNCTION Fxy, x, y

  RETURN, SQRT(x)*EXP(-(y-2)^2)

END


이렇게 하여 얻은 결과값은 아마 3.3265 정도가 될 것입니다. 이 수식의 모습을 도식적으로 보면 다음 그림과 같습니다.



비슷한 요령으로 여러가지 경우들에 대하여 적용할 수 있으므로 여러분들도 다양하게 시도해보시기 바랍니다. 마지막으로 INT_2D에 관하여 한가지만 더 언급해야 할 것 같은데, 바로 ORDER 키워드에 관한 것입니다. 이 키워드의 값은 0 또는 1이 될 수 있고 디폴트 값은 0입니다. 0은 dydx 순서의 적분에 해당되고 1은 dxdy 순서의 적분에 해당됩니다. 위의 예제들에서는 이 키워드를 따로 명시하지는 않았기 때문에 디폴트 값인 0으로 간주되었습니다. 즉 dydx 순서의 적분이 수행된 셈입니다. 그런데 만약 위의 예제들에서 ORDER 키워드의 값을 1로 따로 설정한다면 결과값은 어떻게 될까요? 즉 INT_2D가 사용된 부분을 다음과 같이 수정할 경우입니다.


result = INT_2D('Fxy', AB_Limits, 'PQ_Limits'96, ORDER=1)


실제로 이렇게 하여 다시 실행해보면 결과값이 바뀌는 것을 바로 확인할 수 있습니다. 첫번째 수식에 대한 결과의 경우 원래는 16.0이 나왔는데 ORDER 키워드를 1로 설정하면 32.0이 나옵니다. 두번째 수식에 대한 결과의 경우 원래는 3.3265 정도가 나왔는데 ORDER 키워드를 1로 설정하면 4.7044라는 전혀 다른 값이 나옵니다. 그런데 이와 같이 ORDER 키워드의 값을 1로 설정하여 얻은 값들은 모두 틀린 결과입니다.


사실은 이와 같이 ORDER 키워드의 값을 1로 설정할 경우에는 각별히 주의해야 할 것이 있습니다. 바로 적분 구간에 대한 설정에 있어서도 X, Y의 순서를 서로 바꿔주어야 한다는 것입니다. 즉 이 경우에는 X 구간을 0~4로 설정하고 Y 구간을 0~2로 설정해야 합니다. 왜냐하면 ORDER=1 설정에 의하여 적분의 순서가 dydx에서 dxdy로 바뀌게 되는데, 위의 코드에서 AB_LIMITS 및 PQ_LIMITS에 대한 설정도 이에 맞춰서 맞바꿔주어야 하기 때문입니다. 이러한 부분은 자칫 불필요한 혼란을 일으킬 가능성이 큽니다. 따라서 제 개인적인 의견으로는 ORDER 키워드에 대한 별도의 설정은 가급적 피하고 그냥 dydx의 순서로 작업이 되도록 코딩하는 것이 혼란도 줄이고 명료하게 작업하는 방법이 될 것으로 보입니다.


그러면 INT_2D 함수를 이용하여 2-변수 수식에 대한 적분 결과를 얻는 방법에 관한 소개는 이 정도로 마무리하겠습니다. 관련하여 추가적으로 언급할만한 이슈가 나중에라도 또 나온다면 다시 다뤄보도록 하겠습니다.


* 참고로 맨 처음에 제시되었던 그림에 해당되는 수식은 다음과 같습니다. 이 경우에 대해서도 계산 결과를 한번 얻어보시기 바랍니다. 제가 얻은 적분값은 35.69입니다.


FUNCTION Fxy, x, y

  RETURN, y*COS(x^2)+4

END

반응형