레이블이 수학인 게시물을 표시합니다. 모든 게시물 표시
레이블이 수학인 게시물을 표시합니다. 모든 게시물 표시

2024년 1월 27일 토요일

람베르트 W 함수(Lambert W Function)

[경고] 아래 글을 읽지 않고 "람베르트 W 함수"을 보면 바보로 느껴질 수 있습니다.


[그림 1] 람베르트 W 함수의 변화 특성(출처: wikipedia.org)

우리가 쓰는 대부분의 함수는 이미 오래전에 정립되어 대부분의 성질이 잘 규명되어 있다. 반면 람베르트 W 함수(Lambert W function)는 람베르트Johann Heinrich Lambert(1728–1777)가 1758년람베르트 30세, 조선 영조 시절에 시작했지만 1996년김영삼 정부 시절에 와서야 수학 함수로 명확히 인정받았다[1]. 람베르트 W 함수 $W(x)$는 지수 함수 $e^y$를 이용해서 정의한다.

                          (1)

여기서 $y$ = $W(x)$이다. 람베르트 함수에 알파벳 W를 쓴 이유는 기호 계산(symbolic computation) 프로그램인 메이플(Maple)에서 W를 쓰기 때문이다.[메이플에서 W를 쓴 이유는 불분명하다. 아마도 W는 잘 쓰지 않는 알파벳인 이유가 클 것이다. 혹은 이 함수에 기여한 라이트(Edward Maitland Wright) 교수[2]의 첫자를 따서 W라는 설도 있다.]

[그림 2] 가지 자름(branch cut)으로 표현한 람베르트 W 함수의 다가성(출처: wikipedia.org)

람베르트 W 함수는 주어진 $x$에 대해 답이 여러 개인 다가 함수(multi-valued function)이다. 이 개념을 이해하기 위해 식 (1)에 자연 로그 함수 $\log(x)$를 취한다.

                  (2a)

                  (2b)

여기서 $x$ = $x_0 e^{2 n \pi i}$, $x_0$의 편각(偏角, argument)은 $-\pi < \operatorname{arg}(x_0) \le \pi$, $n$ = $0, \pm 1, \pm 2, \cdots$, $n$은 가지를 구별하는 가지 번호(branch number)이다. 변수 $x$가 같더라도 편각은 $2 \pi$의 정수배만큼 달라질 수 있어서 $y$는 하나가 아니고 무한 개의 답이 나온다. 변수 $x$처럼 $y$도 $\log y$ = $\log y_0 + 2m\pi i$로 만들면, 식 (2b)에 따라 $n$ = $m$이다. 여기서 $\log y_0 + y_0 - \log x_0$ = $0$이다. 그래서 람베르트 W 함수의 정확한 표기는 [그림 1]에 나온 $W_n(x)$이다. 여기서 $n$은 [그림 2]와 같은 가지 자름(branch cut)을 가리키는 정수이다. 또한 식 (2)처럼 람베르트 W 함수는 $y$ 기준으로 $\log x$와 관계되고 곱 연산도 있어서 곱 로그(product logarithm)로 부를 수 있다.

[그림 3] 람베르트 W 함수의 등각 사상(출처: wikipedia.org)

등각 사상(conformal mapping) 관점에서 [그림 3]에 보인 람베르트 W 함수와 가지 자름의 특성을 고찰한다. 실수 영역에서 정의한 식 (1)을 복소수 영역으로 확장한다.

                  (3a)

여기서 $w$ = $u+ iv$, $z$ = $x +iy$이다. 람베르트 W 함수에는 식 (2)처럼 로그 특성이 있어서 $z$의 가지 자름은 음의 실수축으로 선택한다. 등각 사상에 따라 $z$에 대한 음의 실수축[$x < 0$, $y$ = $0$]은 $w$ 영역에서 다음과 같이 사상된다.

                  (3b)

식 (3b)에서 $v \ne k \pi$로 놓고 등식과 부등식을 푼다.

                  (3c)

여기서 $y$ = $0$, $k$ = $0, \pm 1, \pm 2, \cdots$이다. 그러면 $z$평면에서 음의 실수축[$x < 0$, $y$ = $0$]에 정의한 가지 자름이 $w$평면으로 넘어가 형성한 곡선을 [그림 3]처럼 그릴 수 있다. 이 곡선은 $v$에 대해 우함수(even function)라서 $v$축에 대칭이다. 먼저 표본화 함수(sampling function) $\operatorname{Sa}(v)$가 0보다 큰 범위인 $2n \pi < v < (2n+1) \pi$를 시작으로 $u$의 범위를 결정한다. 여기서 가지 번호 $n$은 $n \ge 0$으로 제한한다.
  • $2n \pi < v < (2n+0.5) \pi$: $\cos v > 0$이라서 $u$는 0보다 작음
  • $(2n+0.5) \pi \le v < (2n+1) \pi$: $\cos v \le 0$으로 인해 $u$는 0과 같거나 큼
이 결과를 이해하면서 [그림 3]을 보면, 가지 번호 $n$에 대해 $v$ = $2 n \pi$ 및 $u$는 음의 무한대에서 출발해 $v$ = $(2n+0.5) \pi$에서 $u$ = $0$ 되며, $v$ = $(2n+1) \pi$ 및 $u$가 무한대로 점근하는 방식으로 곡선이 끝난다. 가지 번호 $n$ = $0$ 혹은 $v$ = $0$인 경우는 가지 자름을 특별하게 선정한다. 왜냐하면 $v \mathbin{/} \sin v$ = $1 \mathbin{/} \operatorname{Sa}(v)$는 $v$ = $0$에서 잘 정의되기 때문이다. 식 (3b)에 $v$ = $0$을 넣으면, $u < 0$인 모든 값이 가능하다. 다만 $u$ = $-1$ 혹은 $x$ = $-e^{-1}$에서 기울기가 무한대로 가서 함수의 다가성이 생긴다. 그래서 [그림 3]과 같이 $v$ = $0$ 및 $u$ = $-1$에서 가지 자름을 다시 만들게 되어, $(u, v)$ = $(-1, 0)$은 가지점(branch point)이 된다. 이를 이해하기 위해 $W(x)$의 미분을 구한다.

                  (4)

식 (4)에서 분모가 0이 되는 경우는 $w$ = $-1$이 유일하다. 함수 $W(x)$의 다가성을 해석적으로 만들려고 검정색 곡선은 $n$ = $0$에 넣고, 파란색 곡선은 $n$을 하나 더 낮추어 $n$ = $-1$로 배정한다. [그림 1]의 제시처럼 $x$ = $-e^{-1}$을 기준으로 해 $W_0 (x)$와 $W_{-1}(x)$를 각각 정의한다. 이로 인해 $n < 0$인 경우에 대한 가지 자름을 만드는 기준은 $n \ge 0$ 조건과 약간 달라진다.
  • $(2n+1.5) \pi < v < (2n+2) \pi$: $\cos v > 0$이라서 $u$는 0보다 작음
  • $(2n+1) \pi < v \le (2n+1.5) \pi$: $\cos v \le 0$으로 인해 $u$는 0과 같거나 큼
가지 번호가 $n$ = $0$인 $W_0(x)$는 실수 범위의 정의역이 넓어서, $n$ = $0$은 주요 가지(principal branch)가 된다. 이상에 나온 논의를 바탕으로 [그림 3]에 나온 굵은 선 모양의 가지 자름은 $v$를 변화시키면서 $(u, v)$ = $(-v \cos v / \sin v, v)$인 궤적으로 그린다.
람베르트 W 함수는 고차 방정식의 해법을 찾는 과정에서 우연히 발견되었다[1]. 람베르트는 1758년에 $x^m - x + q$ = $0$의 답을 멱급수의 반전(inversion of power series)으로 찾았다. 여기서 $m$은 자연수, $q$는 상수이다. 요즘은 라그랑주 반전 정리(Lagrange inversion theorem)로 편하게 이 방정식을 풀 수 있다. 그 후 1783년오일러 76세, 조선 정조 시절에 오일러Leonhard Euler(1707–1783)는 변수 치환을 통해 더 쉽게 답을 얻는 방법을 제시했다. 오일러의 생각을 따라가려고 원래 고차 방정식에서 $m$ = $\alpha / \beta$, $q$ = $(\alpha - \beta) c$로 놓고 $x$ = $u^{-\beta}$로 치환한다.

                  (5a)

식 (5a)에 대해 $\beta \to \alpha$로 가는 극한을 적용해서 정리한다.

                  (5b)

람베르트 W 함수 형태로 만들기 위해 $\alpha$ = $1$로 두고 간략화한다.

                  (5c)

만약 $\alpha \ne 1$이 아니면 식 (5b)의 양변에 $\alpha$를 곱해서 $\alpha \log u$ = $\log u^\alpha$ = $\alpha c u^\alpha$로 만든다. 그러면 $t$ = $u^\alpha$인 치환을 통해 다시 식 (5b)와 같은 형태가 될 수 있다. 이 모두를 종합한 결과식을 보면 신기하게도 고차 방정식 $x^m - x + q$ = $0$은 람베르트 W 함수 $W(x)$를 내재하고 있다.
람베르트 W 함수가 도입됨으로 인해 지수와 일차 함수, 로그와 일차 함수, 밑수와 지수가 함께 변하는 방정식 등을 손쉽게 해결할 수 있다. 예를 들어, $e^{-ax}$ = $bx + c$를 식 (5c)와 같은 방식으로 풀어본다.

                  (6a)

지수 관계인 $2^4$ = $4^2$처럼 밑수와 지수가 서로 바뀌는 방정식의 결과도 람베르트 W 함수로 공식화된다.

                  (6b)

마찬가지로 밑수와 지수가 같은 지수 함수의 역함수도 람베르트 W 함수에 속한다.

   
               (6c)

람베르트 W 함수의 급수해는 라그랑주 반전 정리(Lagrange inversion theorem)로 손쉽게 획득한다. 식 (1)에 따라 $x$ = $f(y)$ = $ye^y$로 놓고 역함수를 써서 $y$ = $g(x)$ = $f^{-1}(x)$를 계산한다.

                          (7a)

                  (7b)

식 (7b)에 나온 무한 급수의 수렴 구간은 비율 판정(ratio test)으로 결정한다.

                  (7c)

따라서 $|x| < e^{-1}$라면 식 (7b)는 절대 수렴한다.

[참고문헌]
[1] R. M. Corless, G. H. Gonnet, D. E. G. Hare, D. J. Jeffrey, and D. E. Knuth, "On the Lambert W function," Adv. Comput. Math., vol. 5, pp. 329–359, Dec. 1996.
[2] E. M. Wright, "Solution of the equation $ze^z$ = $a$", Proceedings of the Royal Society of Edinburgh Section A: Mathematics, vol. 65, no. 2, pp. 193–203, 1959.

2024년 1월 14일 일요일

에어리 미분 방정식(Airy Differential Equation)

[경고] 아래 글을 읽지 않고 "에어리 미분 방정식"을 보면 바보로 느껴질 수 있습니다.


[그림 1] 에어리 함수의 변화 모습(출처: wikipedia.org)

빛 산란(light scattering)이나 렌즈 초점(lens focus) 변화를 분석할 때 쓰이는 에어리 함수(Airy function)에어리 미분 방정식(Airy differential equation)을 만족한다[1].

                          (1)

식 (1)의 해 중 하나인 제1종 에어리 함수(Airy function of the first kind) $\operatorname{Ai}(x)$는 무한 적분을 이용해서 정의한다.

                          (2)

식 (2)를 식 (1)에 직접 대입해서 $\operatorname{Ai}(x)$는 식 (1)의 타당한 해임을 보일 수 있다.

                          (3)

식 (3)에서 유도한 마지막 식이 0이 되는 이유는 복소 해석학 혹은 함수론(complex analysis or complex function theory) 때문이다.

[그림 2] 제1종 에어리 함수를 위한 닫힌 경로(출처: wikipedia.org)

[그림 2]처럼 실수축에 위치한 적분 경로를 복소 반평면으로 확대한다. 그러면 코쉬의 적분 정리(Cauchy's integral theorem)에 의해 적분값은 0이 된다.

                          (4)

여기서 $R$이 커짐에 따라 경로 $c_1$상의 적분값은 조르당의 보조 정리(Jordan's Lemma)에 의해 0이 된다.[∵ $z$의 허수부가 0보다 매우 커져서 $e^{iz}$는 지수 함수적으로 감쇠한다.]
식 (2)를 참고해서 제2종 에어리 함수(Airy function of the second kind) $\operatorname{Bi}(x)$도 비슷하지만 식 (2)와 독립되게 정의한다.

                          (5)

여기서 사인 함수는 식 (2)의 코사인 함수와 독립이어서 도입되며, 지수 함수는 $\operatorname{Bi}(x)$가 에어리 미분 방정식의 해가 되도록 돕는다. 식 (3)과 동일하게 $\operatorname{Bi}(x)$를 넣고 미분 방정식을 풀어쓴다.

                          (6)

식 (6)의 최종 결과도 복소 함수론으로 증명해야 한다. [그림 2]와 다르게 실수축과 허수축이 모두 포함되도록 적분 경로를 설정한다.
 
[그림 3] 제2종 에어리 함수에 쓰는 닫힌 경로

경로 $c_1, c_2, c_3$이 닫히도록 복소수 $z$를 정의해서 코쉬의 적분 정리에 넣는다.

                          (7)

적분 구간이 $[0, -R]$로 가는 경우도 식 (7)과 비슷하게 구해서 적분값을 $i$로 계산한다. 이 두 결과를 식 (6)에 넣으면 최종값은 0이 되어 유도가 완성된다.
다만 $t$가 커질 때 피적분 함수의 위상은 $t^3$ 크기로 빠르게 변화해서 적분이 존재하는지 확인해야 한다. 식 (2)에 부분 적분을 적용해서 제1종 에어리 함수가 수렴하는 특성을 증명할 수 있다.

                          (8)

여기서 $a$는 분모를 0이 되지 않기 위해 선택한 0보다 큰 적당한 양수이다. 아니면 $u_n$ = $(2n+1) \pi/2$ = $t_n^3/3 + xt_n$으로 두고, 적분 구간을 $\pi$로 잘라서 교대 급수(alternating series)를 만든다.

                          (9)

여기서 $u$ = $t^3 + xt$, $n$ = $1,2,\cdots$이다. 항 $a_n$으로 나타낸 적분은 $n$이 커질수록 구간이 계속 짧아져서 적분값의 크기는 줄어든다. 그래서 $a_n$은 각 항의 크기가 단조 감소하며 $0$으로 수렴하기 때문에, 라이프니츠 기준(Leibniz criterion)에 의해 에어리 함수는 수렴한다. 제2종 에어리 함수에 대해서도 동일한 논리를 적용해서 수렴성을 이끌어낼 수 있다.

[그림 4] 함수값 $\operatorname{Ai}(0)$ 계산에 사용되는 적분 경로

에어리 함수는 위상이 3차 함수로 변해서 함수값 계산이 쉽지 않다. 다행히 $x$ = $0$인 경우는 수월하게 답이 나온다. 식 (2)에 $x$ = $0$를 대입해서 지수 함수 형태로 만든다.

                          (10a)

식 (10a)에서 $e^{iu}$를 포함한 적분은 닫힌 경로를 [그림 4]처럼 선택해서 결과를 얻는다.

                          (10b)

여기서 $\Gamma(x)$는 감마 함수(gamma function), $c_2$상의 적분은 피적분 함수가 지수 함수적으로 감쇠해서 0, $c_4$를 가진 적분은 반지름이 너무 작아서 0이 된다. 피적분 함수가 $e^{-iu}$인 경우는 [그림 4]를 쓸 수 없고 허수부가 0보다 작은 닫힌 경로[그림 4에 나온 경로를 $x$축에 대해 대칭한 경로]를 선택한다.

                          (10c)

식 (10b)와 (10c)를 식 (10)에 넣어서 $\operatorname{Ai}(0)$를 결정한다.

                          (11a)

여기서 $i$ = $e^{i \pi/2}$이다. 식 (11a)에 다시 오일러의 반사 공식(Euler's reflection formula)으로 만든 $\Gamma(1/3) \sin(\pi/3) / \pi$ = $1 / \Gamma(2/3)$을 적용해서 최종 결과를 얻는다.

                          (11b)

식 (2)를 미분해서 에어리 함수의 도함수도 구한다.

                          (12)

함수값 $\operatorname{Ai}(0)$처럼, 식 (12)의 첫째식에 $x$ = $0$을 대입해서 $\operatorname{Ai}'(0)$을 계산한다.

                          (13a)

식 (13a)에 식 (10b)와 (10c)를 대입해서 정리한다.

                          (13b)

                          (13c)

식 (11b)와 (13c)의 유도 과정을 참고해서 $\operatorname{Bi}(0)$과 $\operatorname{Bi}'(0)$을 유도한다.

                          (14a)

                          (14b)

                          (15a)

                          (15b)

지금과 같이 식 (2)와 (5)를 적분해서 모든 $x$에 대한 에어리 함수값을 모두 구할 수 있지만, 계속 이런 방식으로 적분할 수는 없다. 그래서 에어리 함수의 근본인 에어리 미분 방정식으로 돌아가서, 에어리 함수와 베셀 함수(Bessel function) 사이의 관계식을 도출한다[1]. 이를 위해 $u$ = $-x$, $y(x)$ = $\sqrt{u} \phi(u)$로 변수 치환한다.

                          (16a)

                          (16b)

다시 $v$ = $(2/3) u^{3/2}$로 치환해서 식 (16b)를 다시 기술한다.

                          (16c)

식 (16c)의 마지막식은 베셀의 미분 방정식(Bessel's differential equation)이므로, 서로 독립적인 두 해는 $\phi$ = $J_{\pm 1/3}(v)$, $y$ = $\sqrt{-x}J_{\pm 1/3}[2/3\cdot(-x)^{3/2}]$이다. 따라서 입력이 음수인 에어리 함수를 베셀 함수의 선형 결합으로 표현한다.

                          (17a)

여기서 $c_1, c_2, c_3, c_4$는 결정해야 할 상수이다. 식 (17a)의 첫째식에 $\operatorname{Ai}(0)$, $\operatorname{Ai}'(0)$을 대입해서 $c_1, c_2$를 정한다.

                          (17b)

             (17c)

식 (17b)와 (17c)를 각각 풀어서 $c_1$ = $c_2$ = $1/3$을 얻어서 식 (17a)의 첫째식에 대입한다. 마찬가지 방식으로 $c_3, c_4$를 계산한다.

                  (17d)

식 (17d)에 따라 $c_3$ = $1/\sqrt{3}$, $c_4$ = $-1/\sqrt{3}$이다. 상수 $c_1, c_2, c_3, c_4$를 식 (17a)에 넣어서 공식을 완성한다.

                          (18)

여기서 $x \ge 0$이다. 에어리 함수의 입력이 0보다 크면, 식 (18)에 $x$ 대신 $-x$를 넣는다. 다만 베셀 함수 입력을 다룰 때는 [그림 4]에 있는 가지 자름(branch cut)처럼 음의 실수축에 주의해야 한다. 차수가 1/3인 제1종 베셀 함수는 다음과 같은 과정을 거쳐 제1종 변형 베셀 함수(modified Bessel function)가 된다.

             (19)

여기서 베셀 함수에 해석적 연속(analytic continuation) $J_\nu(e^{i \pi} z)$ = $e^{i \nu \pi} J_\nu( z)$을 적용한다. 다음 단계로 식 (19)를 식 (18)에 넣어서 간략화한다.

                          (20)

여기서 $x \ge 0$이다.

[그림 5] 물잔이 렌즈 역할해서 생긴 소작(燒灼, caustic) 현상(출처: wikipedia.org)

베셀 함수랑 비슷해 보이지만, 다루기가 매우 까다로운 에어리 함수는 도대체 어디에 사용될까? 에어리 함수는 파동(wave), 더 정확히는 광학(optics)에 주로 사용한다[2], [3]. 파동의 중요한 특징은 위상(phase)이라서, 파동이 특정 위치 $\bar r$ = $(x, y, z)$에서 가지는 위상 $\phi(r)$을 균일 평면파(uniform plane wave)처럼 정의한다.

                          (21)

여기서 $k_0$은 진공 중의 파수, $\phi_0$은 기준 위상이다. 식 (21)에 나온 제곱근 함수를 그대로 사용할 수 있으면 좋겠지만, 제곱근 함수는 적분하기 쉽지 않다. 그래서 식 (21)은 주로 테일러 급수(Taylor series)로 전개한 멱급수(power series)로 어림해서 사용된다. 예를 들어, 식 (21)을 $x$에 대해 3차 항까지 전개해서 위상을 $\phi(r)$ $\approx$ $ax^3 + bx^2 + cx+ d$로 근사한다. 이 3차 방정식(cubic equation)에 변수 치환을 해서 위축된 3차 방정식(depressed cubic equation)을 만든다[3].

                          (22)

여기서 $p,q$는 계수 $a,b,c,d$로 만드는 상수, $x$ = $u - b \mathbin{/}(3a)$이다. 만약 테일러 급수를 2차 항까지만 쓰면, 이 경우는 프레넬 근사(Fresnel approximation)가 되고, 프레넬 적분(Fresnel integral)이 관련된다. 위상 성분인 식 (22)가 파동 성질에 끼치는 기여는 연속 파수 $\zeta$를 쓰는 푸리에 변환 형태로 표현된다.

                          (23)

여기서 $p,q,s$는 $\zeta$에 대해 적당한 상수, $x$ = $p/(3s)^{1/3}$, 적분 핵심(integral kernel)인 $F(\zeta)$는 $e^{i \phi}$보다 매우 느리게 변한다고 가정한다.
식 (23)은 빛 산란을 분석할 때 에어리 함수가 등장하는 수학적 논지를 보여주지만, 물리적 이해는 또 다른 차원의 문제이다. 에어리 함수의 특성을 알기 위해 식 (1)의 $x$를 상수 $x_0$으로 가정한다. 그러면 식 (1)은 전형적인 상수 계수 선형 상미분 방정식이 되어서 해가 매우 쉽게 구해진다. 만약 $x_ 0 > 0$이면, 해는 [그림 1]처럼 지수적으로 감쇠하거나 발산한다. 하지만 $x_0 < 0$이면, 해가 삼각 함수의 선형 결합으로 바뀌어 $x_0$에 따라 ($+$)와 ($-$)를 진동한다. 이러한 에어리 함수의 수렴과 진동 특성을 보여주는 예는 [그림 5]에 보여준 소작(燒灼, caustic) 현상이다. 소작은 광선이 동위상으로 모여서 기하 광학(geometrical optics)이 발산하는 영역이다. 에어리 함수 관점에서는 $x$ = $0$을 만족하는 선이 바로 소작선(caustic line)이다. 소작선에 근접해 수렴하는 광선은 $x > 0$, 멀어져 발산하는 경우는 $x < 0$이 되어야 한다. 이 현상은 제1종 에어리 함수의 변화 특성인 [그림 1]이 잘 보여주고 있다. 광선의 수렴과 발산을 이해하는 출발점은 급속 하강 방법(method of steepest descent)에 나오는 안장점(saddle point) 유무이다. 안장점은 접선 기울기가 0이면서 극값을 가지지 않는 점이다. 식 (23)에서 $x$ = $0$일 때는 위상이 $t^3/3$으로 변하며, $t$ = $0$에서 기울기는 0이지만 이 점 근방에서 위상값은 계속 커진다. 그래서 $t$ = $0$은 안장점이 되고, 위상이 빠르게 변하는 적분은 잘 수렴한다. 더 구체적으로 분석하려면 위상 항을 미분한 $\delta \phi$를 적용한다.

                          (24)

만약 $x > 0$이면, 안장점은 아니지만 $t$ = $0$ 근처에서 기울기가 $x$만큼 더 커진 ($+$)라서 위상값은 안정적으로 더 빨리 커진다. 그래서 여전히 적분은 되지만, [그림 1]처럼 적분값이 감쇠하는 특성으로 나타난다. 하지만 $x < 0$에서는 기울기가 ($-$)도 될 수 있어서 위상값의 증감이 생기고, 적분값은 진동하는 방식으로 도출된다. 따라서 식 (24)는 위상 관점에서 본 페르마의 원리(Fermat's principle)이다[2]. 통상적인 페르마의 원리는 빛이 이동하는 시간이 최소가 되는 경로가 실체라는 뜻이다. 위상 기준으로는 위상 차이 $\delta \phi$의 크기가 최소 혹은 0이 되는 경로가 실재라고 생각하면 쉽다.

[참고문헌]
[1] V. Lakshminarayanan and L. S. Varadharajan, Chapter 4. Airy FunctionsSpecial Functions for Optical Science and Engineering, SPIE Press, 2015. (방문일 2024-01-12)
[2] R. D. Blandford and K. S. Thorne, "7. Geometric Optics", Applications of Classical Physics, 2012. (방문일 2024-01-12)
[3] N. C. Albertsen, P. Balling and N. E. Jensen, "Caustics and caustic corrections to the field diffracted by a curved edge," IEEE Trans. Antennas Propag., vol. 25, no. 3, pp. 297–303, May 1977.