2020년 9월 20일 일요일

이산 푸리에 변환(DFT: Discrete Fourier Transform)

[경고] 아래 글을 읽지 않고 "이산 푸리에 변환"을 보면 바보로 느껴질 수 있습니다.


이산 푸리에 변환(discrete Fourier transform, DFT)푸리에 변환(Fourier transform)을 컴퓨터 상에서 구현할 때 사용하는 중요 개념이다. 이산(離散, discrete) 푸리에 변환은 명칭에도 있듯이 연속적이지 않고 이산적으로 변하는 함수에 대한 푸리에 변환이다. 컴퓨터는 원론적으로 연속 함수(continuous function)를 다룰 수 없기 때문에, 수치 해석에서는 컴퓨터에 저장하기 쉬운 이산화(離散, discretization)를 적극적으로 채택한다. 

[그림 1] 이산 푸리에 변환의 적용 예시(출처: wikipedia.org)

연속적으로 변하는 함수 $f(t)$를 이산적으로 바꾸려면, 연속 함수에 [그림 2]와 같은 표본화(sampling)를 적용한다. 표본화하는 주기를 $\Delta t$, 표본 개수를 $M$이라 하면, $f(t)$를 이산화한 전체 표본화 시간은 $T$ = $M \Delta t$가 된다. 하지만 표본화 시간인 $0 \le t < M\Delta t$에서는 함수 $f(t)$의 값을 알지만, 나머지 시간에서는 함수값을 알 수 없게 된다. 이때는 $f(t)$를 과감하게 주기 $T$를 가진 주기 함수(periodic function)라 가정한다. 그래서 표본화된 함수 $f_\text{sa}(t)$는 주기가 $T$이면서 $t$ = $m \Delta t$[$m$ = $0, \pm 1, \pm 2, \cdots$]에서만 값이 이산적으로 존재하는 주기 함수가 된다.

[그림 2] 연속 함수의 표본화 예시(출처: wikipedia.org)

DFT 관계식을 유도하기 위해, 먼저 복소 푸리에 급수(complex Fourier series)푸리에 변환(Fourier transform) 공식을 고려한다.

                  (1)

                  (2)

여기서 $\omega_0$ = $2 \pi/T$이다. 푸리에 변환은 $F(\omega)$ = $\lim_{T_\text{tot} \to \infty} \widetilde{F}_m T_\text{tot}$로 정의한다. 여기서 $T_\text{tot}$는 $f(t)$가 존재하는 전체 시간이다. 함수 $f(t)$가 전체 표본화 시간 $T$를 주기로 반복되고 이상적인 표본화는 디랙 델타 함수(Dirac delta function)라고 가정하면, 푸리에 변환보다는 식 (1)에 있는 푸리에 급수가 표본화된 함수 $f_\text{sa}(t)$의 해석에 더 적절하다. 그러면 $f_\text{sa}(t)$에 대해, 푸리에 변환에서 유추한 새로운 DFT의 푸리에 계수 $F_k$를 다음처럼 얻을 수 있다.

                  (3)

여기서 $\delta(t)$는 디랙 델타 함수, 표본화 계수 $f_m$ = $f(m \Delta t)$이다. 식 (3)에 $k$ = $k + M$을 대입하면, $F_k$ = $F_{k+M}$이 항상 성립한다. 이 결과는 $F_k$가 주기 $M$을 가지고 반복되는 계수임을 뜻한다. 푸리에 계수 $F_k$를 식 (1)에 대입해서 $F_k$의 무한 합으로 $f_\text{sa}(t)$를 다음과 같이 표현한다.

                  (4a: 푸리에 급수 이용한 정의)

                  (4b)

여기서 $T$ = $M \Delta t$, $f_\text{sa}(t)$의 정의역을 $0 \le t < T$로 제한, $\delta(2 \pi t / \Delta t)$ = $\Delta t / (2 \pi) \delta(t)$이다; 식 (4b)에서 첫째 줄의 무한한 $k$[= $0, \pm1, \pm2, \cdots$]를 나눗셈의 원리를 이용해 둘째 줄에서 유한한 $0 \le k + mM < M$[$m$ = $0, \pm1, \pm2, \cdots$]으로 바꾸며, 둘째 줄의 $m, k$는 각각 몫(quotient)과 나머지(remainder)가 된다. 식 (3)의 첫째식과 식 (4)의 마지막식을 비교해서 다음 결과를 유도한다.

                  (5)

따라서 이산 푸리에 변환쌍(DFT pair)을 다음처럼 표현할 수 있다.

                  (6)

이산 푸리에 역변환(inverse discrete Fourier transform, IDFT)에 $1/M$이 나오는 이유는 적분의 이산화로 쉽게 설명할 수 있다. 미분소 $d \omega$의 이산화는 $\Delta \omega$ = $\omega_0$ = $2 \pi / T$ = $2 \pi / (M \Delta t)$이다. 이 결과를 식 (2)에 대입하고 $\Delta t$ = $1$로 두면, 최종적으로 $1/M$만 남는다. 식 (6)을 한 차원 더 확장해서 2차원 이산 푸리에 변환쌍을 정의한다.

                  (7)

식 (6)의 첫째식을 둘째식에 대입하거나 둘째식을 첫째식에 대입하면, 다음 항등식이 유도된다.

                  (8a)

                  (8b)

여기서 $\delta_{ml}$은 크로네커 델타(Kronecker delta)이다.
식 (8)의 최종식은 등비 급수(geometric series)를 적용해서 간단히 얻을 수도 있다.

                  (9)

또한 식 (9)는 인수 분해 공식 $z_m^M - 1$ = $(z_m - 1)$$(1 + z_m + z_m^2 + \cdots + z_m^{M-1})$에 $z_m$ = $e^{i2 \pi m / M}$, $z_m^M$ = $1$을 대입해 0이 됨을 보일 수 있다. 여기서 $m \ne 0$이다. 서로 다른 모든 기저 항(basis term)을 더하면 0이 된다는 뜻인 식 (9)는 기저 벡터(basis vector)선형 독립(linear independence)을 직관적으로 보여준다. 즉, $[1~z_m~z_m^2~\cdots~z_m^{M-1}]$인 기저 벡터는 모든 성분이 1인 기저 벡터에 직교한다. 식 (9) 대신 식 (8)의 첫째식을 고려하면, $[1~z_m~z_m^2~\cdots~z_m^{M-1}]$과 $[1~z_{m'}~z_{m'}^2~\cdots~z_{m'}^{M-1}]$는 켤레 전치 행렬(conjugate transpose) 연산에서 직교하므로, 이 기저 벡터로 구성한 행렬은 복소 영역의 직교 행렬인 유니터리 행렬(unitary matrix)이 된다. 여기서 $z_m^k$ = $z_k^m$이다. 따라서 식 (6)을 행렬 곱으로 바꾼 결과에는 유니터리 행렬이 등장하기 때문에, 이 행렬의 역행렬은 켤레 전치 행렬로 쉽게 구한다. 다시 말해 식 (9)는 식 (6)이 유니터리 관계임을 말하고 있어서, DFT를 보는 또 다른 관점은 유니터리 연산이다.
푸리에 변환(Fourier transform)의 대표 연산인 길쌈(convolution)을 이산화해서 DFT용 길쌈을 정의할 수 있다.

                  (10)

여기서 $f_n, g_n$의 지표(index) $n$이 가진 주기는 $M$이다. 길쌈의 DFT는 푸리에 변환과 마찬가지로 각 DFT의 곱이 된다.

                  (11a)

                  (11b)

여기서 $p_k$ = $2 \pi k \mathbin{/} M$이다. 식 (8b)를 사용해서 함수 곱의 DFT는 각 DFT의 길쌈으로 표현됨도 유도한다.

                  (12a)

                  (12b)

식 (11b)와 (12b)를  종합해서, IDFT에 대한 길쌈 관계식도 추가로 얻는다.

                  (13)

이상에서 논의한 기본적인 이산 푸리에 변환을 바탕으로 더욱 다양한 이산 변환을 생성할 수 있다.


   1. 이산 코사인 변환(discrete cosine transform, DCT)   

[그림 1.1] 여러 가지 이산 코사인 변환(출처: wikipedia.org)

이산 코사인 변환(discrete cosine transform, DCT)은 DFT의 현실적인 변형이다. DFT는 푸리에 변환을 실제 수치 해석에 적용할 때 매우 유용하다. 하지만 DFT의 계수 $F_k$는 복소수라서 다루기 불편하다. 만약 DFT 대신 DCT를 쓰면, 변환 결과는 실수가 되기 때문에 아주 편해진다. 현실적인 신호 처리에 필수적인 DCT는 아메드Nasir Ahmed(1940–)가 1972년아메드 32세, 박정희 정부 시절에 제안했다. [그림 1.1]에 있는 여러 DCT 중에서 가장 쉬운 형태인 DCT-I을 보자. [그림 1.1]처럼 모든 DCT는 우함수(even function) 특성을 가지기 때문에, DFT 기준으로 보면 $f_m$ = $f_{-m}$이 항상 성립한다. 또한 표본화된 함수는 계속 반복되므로, DCT는 주기적인 특성이 나타난다. 예를 들어, DCT-I의 주기는 $2M-2$로 나온다.[그림 1.1에서 빨간색 표본이 $M$개, 이어지는 검정색 표본은 $M-2$개가 된다.] 하지만 DCT를 적용하는 함수가 반드시 우함수이어야 할까? 아니다. 함수 $f(t)$는 우함수일 필요는 없고, 어떠한 $f(t)$이든지 [그림 1.1]과 같이 대칭적으로 함수값을 배정하여 우함수를 강제적으로 만들 수 있다. 즉, [그림 1.1]에서 빨간색은 직접 측정해서 대상의 정보를 가진 실제 표본이고, 검정색은 각 DCT 방법에 맞도록 실제 표본으로 생성한 부차적인 계산 표본이다. [그림 1.1]의 DCT-I처럼 $0 \le m < M$에서 표본화한 함수값[그림 1.1에서 빨간색]을 $x_m$ = $f(m \Delta t)$라 한다. 그러면 표본화한 범위를 넘는 $M \le m < 2M-2$에 대해, 우함수 특성인 $f_m$을 다음처럼 표현할 수 있다.

                  (1.1)

여기서 $f_m$은 $m$ = $M-1$ 기준으로도 우함수이다. 의도적으로 주기 $2M-2$를 가진 우함수로 만든 DFT의 표본화 계수 $f_m$을 식 (6)에 대입해 본다.

                  (1.2)

여기서 $F_k$의 주기도 $2M-2$이다. DCT-I의 표본화 및 푸리에 계수를 각각 $x_m$ = $f_m$, $X_k = F_k /2$로 정의해 정리하면, DCT-I의 변환식을 얻을 수 있다.

                  (1.3a)

                  (1.3b)

                  (1.3c)

여기서 $\alpha_{m, M}$ = $1+\delta_{m0} + \delta_{m,M-1}$, $\Sigma'$은 으뜸 합산(prime summation) 기호[1], [2]로서 처음과 끝 항은 $1/2$만 더하라는 합 연산이다. DCT-I의 표본화 계수 $x_m$이 실수라서 $X_k$도 당연히 실수가 된다. 또한 식 (1.3)에 의해 $X_k$는 우함수이면서 $k$ = $M-1$ 기준으로도 우함수가 된다. 즉, $X_k$ = $X_{-k}$ = $X_{2M-2-k}$[$k$ = $0, 1, \cdots, M-1$]가 성립한다. 그러면 $k$의 범위를 $0 \le k < M$으로 제한해도 문제는 없다. 이 특성을 식 (6)에 대입해서 $x_m$을 위한 DCT-I의 역변환식도 얻는다.

                  (1.4a)

                  (1.4)

여기서 $0 \le m < M$이다. 이산 코사인 역변환에 등장한 $2/(M-1)$은 푸리에 코사인 역변환(inverse Fourier cosine transform, IFCT)으로 쉽게 설명한다. IFCT에서 $\pi/2 \cdot d \omega$는 $2/\pi \cdot 2 \pi /[(2M-2)\Delta t]$ = $2/(M-1)$으로 이산화되어서 식 (1.4b)와 같은 상수가 나온다. 식 (1.3)와 (1.4)에서 새로운 DCT-I 공식을 얻었지만, DCT의 표현식이 간결해 보이지는 않는다. 그래서 DFT의 표본화 계수 $f_m$을 다음처럼 다시 정의해서 DCT-I과는 다른 공식 전개를 해본다.

                  (1.5)

여기서 $f_m$의 주기는 $4M$이다. 식 (1.2)와 비슷하게 식 (1.5)를 식 (6)에 대입한다.

                  (1.6)

여기서 $F_k$의 주기도 $4M$이다. DFT의 푸리에 계수를 바꾸어 $X_k = F_k /2$로 정의하면, 식 (1.6)에 의해 [그림 1.1]에 있는 DCT-II의 푸리에 계수를 얻는다.

                  (1.7)

여기서 $0 \le k < M$이다. 식 (1.3)과 비교하면 식 (1.7)은 매우 깔끔한 DCT-II의 표현식이다. 그래서 보통 DCT라고 하면 DCT-II를 의미한다. DCT-II 푸리에 계수 $X_k$의 주기 특성을 알기 위해 $X_{2M \pm k}$를 식 (1.7)에 대입하면, [그림 1.1]에 있는 DCT-III와 비슷하게 $X_{2M \pm k}$ = $-X_k$를 얻을 수 있다. 푸리에 계수 $X_k$의 주기는 분명 $4M$이지만, $X_k$와 $X_{2M \pm k}$는 서로 기함수 관계이기도 하다. 그러면 식 (1.7)에서 $k$의 범위를 $0 \le k < M$로 제한해도 문제없다. 이를 이용해 DCT-II의 역변환식을 구해보자.

                  (1.8)

여기서 $X_{2M}$ = $-X_0$. 기함수 관계 $X_{2M \pm k}$ = $-X_k$를 다시 식 (1.8)에 적용해서 간결하게 정리한다.

                  (1.9)

여기서 $0 \le m < M$, $X_M$ = $0$이다. 노이만 수(Neumann number) $\varepsilon_k$를 이용해 식 (1.9)을 더 간단하게 표현한다.

                  (1.10)

여기서 $\varepsilon_m$ = $2 - \delta_{m0}$, $\delta_{ml}$은 크로네커 델타(Kronecker delta)이다. 식 (1.3)와 (1.10)에서 유추해서 [그림 1.1]에 있는 DCT-III의 변환식도 정의할 수 있다.

                  (1.11)

여기서 $0 \le k < M$. DCT-II와 DCT-III은 비율 상수가 약간 다르지만 서로 역변환 관계를 가진다. 식 (7)에 있는 2차원 DFT와 비슷하게 2차원 이산 코사인 변환쌍을 정의할 수 있다. 예를 들어, DCT-II의 변환쌍은 다음과 같다.

                  (1.12)

                  (1.13)

[그림 1.1]에 4가지 종류의 DCT를 소개했지만, DCT의 표본화 계수 $x_m$에 기반을 둔 DFT의 표본화 계수 $f_m$을 다양하게 정의해서 다채로운 우리만의 DCT를 정의할 수도 있다. DCT를 어떻게 정의하든지 우함수 혹은 기함수의 대칭을 사용하므로, 표본화된 함수의 주기는 원래 경우보다 배 이상 늘어난다. 주기 $T$가 늘어나면, $\omega_0$ = $2 \pi/T$의 관계에 의해 주파수 영역의 계수 간격이 줄어든다. 따라서 DFT보다는 DCT의 푸리에 계수가 저주파에 더 많이 모이므로, 고주파 성분을 자르는 손실 있는 영상 압축 등에 DCT를 많이 사용한다.
DCT-I으로 얻어지는 풍성한 결과 중의 하나는 다음과 같은 항등식이다. 식 (1.4b)에 식 (1.3b)를 넣고 간단히 정리하면 곧 바로 증명이 된다.

                  (1.14a)

                  (1.14b)

여기서 $m$과 $m'$가 변하는 범위는 $0, 1, \cdots, M-1$, $\Sigma'$은 으뜸 합산(prime summation)이다. DCT-I을 쓰는 대신 식 (9)에 실수부를 선택하고 첨자 $k$를 두 영역으로 정리해도 유도가 가능하다.


                  (1.15a)

                  (1.15b)

                  (1.16)

여기서 $k$의 영역은 $k_1$ = $0, 1, \cdots, M-2$와 $k_2$ = $M-1, M, \cdots, 2M-3$으로 분해된다. 식 (1.15a)에서 $m$이 짝수라면, $k$와 $k + M-1$ 항은 서로 같은 값이어서[$\because$ $\pi m / (M-1) \cdot (M-1)$ = $m \pi$] $k_1$과 $k_2$를 가진 급수로 나누어진다. 다만 $m$이 홀수인 경우는 약간 다르게 접근한다. 새로운 첨자 $k_1$과 $k_2$로 만든 급수 합은 부호가 반대라서[$\because$ 이 두 급수를 합한 값은 자동적으로 0], $k_1$ 전체[= $0, 1, \cdots, M-2$]와 $M-1$을 가진 급수를 다시 $k_1$과 $M - 1 - k_1$인 항으로 정리해서 식 (1.15b)를 증명한다. 이때 $m$ = $0$인 급수는 합이 $M$이 되어 식 (1.15b)와 다르므로, $k$ = $0$과 $M-1$의 값을 반으로 줄여서 식 (1.15b)가 되게 한다.


   2. 이산 사인 변환(discrete sine transform, DST)   

[그림 2.1] 여러 가지 이산 사인 변환(출처: wikipedia.org)

이산 코사인 변환과 동일한 방법으로 푸리에 사인 변환(Fourier sine transform)의 이산화인 이산 사인 변환(discrete sine transform DST)을 정의할 수 있다. 표본화 계수가 $m$ = $0$과 $M-1$을 기준으로 기함수라 가정해서 DFT를 적용한다.

                  (2.1)

                  (2.2)

식 (2.1)의 둘째식에 강제로 $m$ = $M-1$을 넣으면, $f_{M-1}$ = $-f_{M-1}$이 되어서 $f_{M-1}$ = $0$이 되어야 한다. 비슷하게 $m$ = $2M-2$인 경우에 $f_{2M-2}$ = $f_0$ = $-f_{0}$으로 인해 $f_0$ = $0$으로 나온다. 그래서 식 (1.3)처럼 $x_m$ = $f_m$, $X_k$ = $i F_k/2$로 놓고 식 (2.2)를 깔끔히 정리한다.

                  (2.3)

여기서 $x_0$ = $x_{M-1}$ = $0$이다. 식 (1.3)과 같은 형태로 식 (2.3)을 다시 기술한 공식은 DST-I이라 부른다.

                  (2.4)

여기서 $\Sigma'$은 식 (1.3c)에 도입한 으뜸 합산(prime summation)이다[1], [2]. 표본화 계수 $x_m$과 유사하게 식 (2.4)에 의해 $X_0$ = $X_{M-1}$ = $0$이다. 또한 $X_k$는 $k$ = $0$과 $M-1$에 대해 기함수를 이룬다. 즉, $x_m$처럼 $X_k$ = $-X_{2M-2-k}$가 성립한다. 식 (2.4)에 대응하는 이산 사인 역변환(inverse discrete sine transform, IDST)은 다음과 같다.

                  (2.5a)

                  (2.5b)

이산 사인 변환인 식 (2.4)에서 $x_0$과 $x_{M-1}$은 0이라서 의미가 없으므로, 이 두 값을 제외해서 전체 표본수가 $M$이 되도록 식 (2.4)와 (2.5)의 $M$을 $M+1$로 늘린다.

                  (2.6)

                  (2.7)

여기서 $m$ = $0, 1, \cdots, M-1$, $k$ = $0, 1, \cdots, M-1$이다. 때때로 식 (2.6)을 DST-I이라 명하기도 한다.
식 (2.4)를 식 (2.5b)에 대입해서 DCT와 비슷한 항등식을 만들 수 있다.

                  (2.8)

여기서 $m$과 $m'$의 범위는 $1, 2, \cdots, M-2$로 한정한다. 혹은 식 (2.8)의 좌변에 삼각 함수 항등식을 적용하고 식 (1.15b)를 넣어서 더 쉽게 증명할 수도 있다.

                  (2.9)

여기서 계산을 쉽게 하기 위해 급수 범위를 $0 \le k \le M-1$로 늘린다.


   3. 이산 혼합 푸리에 변환(discrete mixed Fourier transform, DMFT)   

이산 혼합 푸리에 변환(discrete mixed Fourier transform, DMFT)은 이산 사인 및 코사인 변환을 섞어서 정의한다[1], [2]. DMFT를 이산화하기 전의 적분 변환은 바로 혼합 푸리에 변환(mixed Fourier transform, MFT)이다.

                  (3.1a)

                  (3.1b)

여기서 $x_m$ = $m \Delta x$, $p_k$ = $k \pi / [(M-1)\Delta x]$ = $k \Delta p$, $\Delta p\Delta x$ = $\pi /(M-1)$, $\alpha_{m, M}$ = $1+\delta_{m0} + \delta_{m,M-1}$, $\Sigma'$은 으뜸 합산(prime summation)[1], [2], $\alpha$는 $x_0$ 위치의 경계 조건을 결정하는 복소수인 혼합 계수(mixed coefficient)이다. 조건 $\Re[\alpha] \le 0$에 대한 DMFT의 역변환은 다음처럼 공식화된다.

             (3.2a)

             (3.2b)

혼합 푸리에 역변환(inverse mixed Fourier transform, IMFT)의 특성으로 인해, $\Re[\alpha] > 0$인 역변환은 다소 복잡해진다.

                  (3.3)

만약 $\alpha$ = $0$이면, DMFT는 전형적인 이산 코사인 변환으로 간략화된다. 하지만 혼합 푸리에 변환을 만드는 적분을 단순히 이산화해서 표현하기 때문에, 아래 적분과 다르게 $e^{-\alpha x_m}$과 $\alpha \sin(p_k x_m) - p_k \cos (p_k x_m)$의 급수에는 직교성이 성립하지 않는다.

                  (3.4)

                  (3.5)

결국 식 (3.3)에 변형을 가해서 $\Re[\alpha] > 0$인 경우에도 이산화한 두 함수 사이에 직교성이 생기도록 $F_k$를 재정의한다.

                  (3.6)

여기서 $\operatorname{Sa}(x)$ = $\sin x \mathbin{/} x$, $F_0$ = $F_{M-1}$ = $0$, $M$이 충분히 클 때에 $\operatorname{Sa}(p_k \Delta x)$ $\approx$ $1$, $\operatorname{Sa}(x)$를 도입함으로써 $p_k$와 $\cos(p_k x_m)$의 곱이 $\sin (p_k x_m)$의 차분(差分, difference)으로 바뀐다. 식 (3.6)에 등장한 차분은 중심 $x_m$에 대해 대칭적으로 차이를 만들어서 대칭 차분(symmetric difference) $\Delta_s$가 된다. 여기에 곱의 차분 공식을 적용해서 차분 관계를 $f_m$으로 바꾼다.

             (3.7)

여기서 중점 근사에 따라 $f_0 \sin(p_k x_0)$ = $0$ $\approx$ $[f_1 \sin(p_k x_1) + f_{-1} \sin(p_k x_{-1})]/2$, $f_{M-1} \sin(p_k x_{M-1})$ = $0$ $\approx$ $[f_M \sin(p_k x_M) + f_{M-2} \sin(p_k x_{M-2})]/2$이다. 식 (3.7)은 근사식이지만 $M$이 충분히 크면 항등식에 가까워진다. 근사이지만 타당한 식 (3.7)의 최종 결과를 식 (3.6)에 대입해서 $f_m$에 대한 2계 차분 방정식(the second order difference equation)이 생기도록 한다.

                  (3.8)

만약 $M$을 충분히 키우면, 근사식 (3.8)을 등식처럼 생각할 수 있다. 또한 식 (3.8)은 정확히 혼합 푸리에 변환[= $\int_0^\infty [df(x)/dx + \alpha f(x)]\sin(px)\,dx$]의 이산화가 된다. 상이한 관점으로 식 (3.8)은 이산 사인 변환이라서, 식 (2.5b)에 따라 $g_m \mathbin{/} (2 \Delta x)$의 역변환이 바로 얻어진다.

                  (3.9)

여기서 $g_0$ = $g_{M-1}$ = $0$, $m$ = $1,2,\cdots,M-2$이다. 미분 방정식의 특수해와 일반해 개념을 써서, 식 (3.9)를 $f_m$에 대해 풀어서 나오는 해는 $f_m$ = $f_{m,p} + f_{m,g}$로 둔다. 여기서 $f_{m,p}$와 $f_{m,g}$는 각각 차분 방정식의 특수해와 일반해이다. 식 (3.9)를 만족시키는 특수해 $f_{m,p}$는 식 (3.6)을 참고해서 설정한다.

                  (3.10)

여기서 $A_k$는 결정해야 하는 상수이다. 식 (3.9)에 나온 $f_m$에 식 (3.10)을 대입한 후, 그 결과의 좌변과 우변을 비교해서 $A_k$를 구한다.

                  (3.11)

여기서 $m$ = $1,2,\cdots,M-2$이다. 따라서 식 (3.11)에서 구한 $A_k$를 이용해서 $f_{m,p}$를 완전히 공식화한다.

             (3.12)

식 (3.9)에 있는 2계 차분 방정식의 일반해 $f_{m,g}$는 특성 방정식을 풀어서 정한다.

                  (3.13)

여기서 양끝점인 $m$ = $0$과 $M-1$에서 $\Delta_s (r^m) \mathbin{/} (2\Delta x)$ = $-\alpha r^m$을 만족하며, $c_1, c_2$는 $m$과 무관한 경계 조건으로 계산하는 상수[간단하게 $f_0$과 $f_{M-1}$로 $c_1, c_2$를 결정 가능]이다. 식 (3.12)와 (3.13)을 합쳐서 이산 혼합 푸리에 역변환(inverse discrete mixed Fourier transform, IDMFT)을 확정한다.

                  (3.14)

여기서 일반해를 표현할 때 식 (3.13)의 $(-r_1)^{-m}$ 대신 $m$에 대해 대칭적인 $(-r_1)^{M-m-1}$을 쓴다.[이때 일반해가 달라 식 (3.13)과 (3.14)의 $c_2$는 같지 않다.] 연립 방정식 개수 관점에서 $c_1, c_2$의 필요성을 이해할 수 있다. 식 (3.9)로 기술되는 차분 방정식의 개수는 $M-2$개밖에 되지 않아, $M$개인 미지수 $f_m$을 명확히 결정하기에는 방정식 2개가 부족하다. 이 결함은 $f_m$ 중의 어느 2개를 알고 있다는 가정으로 해결하며, 이런 타개책을 수학에서는 일반해라 명한다. 
상수 $\alpha \Delta x$가 매우 작은 경우에 $r_1, r_2$를 지수 함수로 근사할 수 있다.

                  (3.15)

식 (3.15)에 따라 $r_1$을 감소하는 지수 함수처럼 간주할 수 있어서 식 (3.5)처럼 급수 합을 정의한다.

                  (3.16)

신기하게도 식 (3.5)와 다르게, 식 (3.16)은 합산이 항상 0이 나오는 직교성(orthogonality)이 성립한다.

[DMFT 기저 함수의 직교성]

                          (3.17)

여기서 $k$ = $0, 1, 2, \cdots, M-1$, $\Sigma'$은 식 (1.3c)에 나온 으뜸 합산(prime summation) 기호이다.

[증명]
기저 함수를 삼각 함수의 합차 공식으로 분해한다.

                  (3.18a)

식 (3.18a)를 합산하고 첨자 $m$을 맞추어 정리한다.

                  (3.18b)

여기서 $m$ = $1, 2, \cdots, M-2$이다. 처음 항인 $m$ = $0$ 경우는 식 (3.18b)가 성립하지 않지만, 으뜸 합산의 정의[처음 항은 반만 더함]로 인해 $m$ = $0$의 결과는 $m$ = $1$에 대한 식 (3.18b)의 둘째 항이 됨으로써 $m$ = $1$의 합에 포함된다.

                  (3.18c)

마찬가지로 식 (3.17)의 끝 항인 $m$ = $M-1$은 $m$ = $M-2$ 합산의 셋째 항으로 간주된다. 또한 $r_2$ = $(-r_1)^{-1}$은 식 (3.13)의 또 다른 근이므로, 식 (3.18b)에 의해 식 (3.17)의 둘째 줄도 유도된다.
______________________________

DMFT의 기저 함수 직교성을 식 (3.14)에 적용하여 식 (3.3)과 비슷한 방식으로 $c_1, c_2$를 도출할 수 있다[1].

                  (3.19a)

                  (3.19b)

                  (3.19c)

모든 $f_m$에 대한 합산이 번거로울 때는 $f_0, f_{M-1}$을 안다는 조건을 사용해서 $c_1, c_2$를 풀어도 식 (3.19c)와 동일한 결과를 얻을 수 있다. 어떤 경우든지 $F_k$만을 연산해서 모든 $f_m$을 일의적으로 결정하기는 불가능하다. 따라서 DMFT는 수학적으로 완전한 변환이 아니고 결점을 수정해서 써야 하는 불완전 변환이다.


[참고문헌]
[1] G. D. Dockery and J. R. Kuttler, "An improved impedance-boundary algorithm for Fourier split-step solutions of the parabolic wave equation," IEEE Trans. Antennas Propag., vol. 44, no. 12, pp. 1592–1599, Dec. 1996.
[2] J. R. Kuttler and R. Janaswamy, "Improved Fourier transform methods for solving the parabolic wave equation," Radio Sci., vol. 37, no. 2, pp. 5-1–5-11, Apr. 2002.

[다음 읽을거리]

2020년 9월 19일 토요일

리우빌의 정리(Liouville's Theorem)

[경고] 아래 글을 읽지 않고 "리우빌의 정리"를 보면 바보로 느껴질 수 있습니다.


정의된 복소 영역(complex domain)에서 항상 미분 가능한 복소 함수(complex function) $f(z)$는 정칙 함수(正則函數, holomorphic function)라 부른다. 정의역에서 다음과 같은 테일러 급수(Taylor series)가 존재해서 항상 수렴하면, $f(z)$는 해석 함수(analytic function)가 된다.

                  (1)

복소 함수 $f(z)$의 정칙성 혹은 미분 가능성은 코쉬–리만 방정식(Cauchy–Riemann equation)으로 간단하게 확인할 수 있다.

                  (2)

여기서 $f(z)$ = $u(x, y) + i v(x, y)$이다. 복소 함수의 테일러 급수(Taylor series of complex function)에 대한 성질을 이용하면, 특정 영역에서 코쉬–리만 방정식을 만족하는 정칙 함수는 당연히 해석 함수가 됨을 증명할 수 있다.
정칙 함수의 정의역을 특정 영역이 아닌 복소 평면 전체로 확장하면, 이 복소 함수는 전해석 함수(全解析 函數, entire function)가 된다. 전해석 함수는 모든 위치에서 미분 가능하지만, $z$가 무한대로 갈 때 $f(z)$가 발산할 수도 있다. 만약 전해석이면서 모든 위치에서 유계(有界, bounded)라고 복소 함수의 특성을 제한하면 어떤 결과가 얻어질까? 복잡하리라는 우리 예상과는 다르게 이 복소 함수는 간단한 상수 함수가 된다. 복소 해석 함수(complex analytic function)와 유계 특성을 연결한 결과는 리우빌의 정리(Liouville's theorem)라 한다. 복소 함수론을 여러 문제에 적용할 때, 리우빌의 정리는 매우 유용한 도구가 될 수 있다.

[리우빌의 정리]
모든 복소 영역에서 정칙(holomorphic)이며 유계(bounded)인 복소 함수는 항상 상수이다.

[증명]
정칙인 함수는 해석적이어서 복소 함수 $f(z)$를 $z$ = $0$ 근방에서 다음과 같은 테일러 급수로 표현할 수 있다.

                  (3)

식 (3)에 등장하는 계수 $a_m$은 코쉬의 미분 공식(Cauchy's differentiation formula)을 이용해서 결정한다.

                  (4)

식 (4)에 절대값을 적용해보면, $a_m$은 다음과 같은 부등식 관계를 만족해야 한다.

                  (5)

여기서 $|f(z)| \le M$, 적분 경로 $c$는 복소 평면에서 반지름이 $R$인 원이다. 전영역에서 정칙인 복소 함수 $f(z)$를 위한 적분 경로의 반지름 $R$은 계속 커질 수 있다. 그러면 $a_0$을 제외한 $a_m$은 모두 $0$으로 수렴하므로, $f(z)$는 상수 함수가 되어야 한다.
______________________________

[그림 1] 어둠을 비추는 등대(출처: wikipedia.org)

리우빌의 정리는 명제를 구성하는 단어가 단순하고 증명 자체도 매우 쉬워서 리우빌 정리의 유용성을 놓치는 경우가 많다. 하지만 복소 함수의 특성을 진지하게 고민할 때, 리우빌의 정리는 우리 마음에 높이 솟은 밝은 등대가 된다. 수학 정리는 우리 지성의 한계를 더 넓은 영역으로 인도해준다. 리우빌의 정리를 사용하는 대표적인 예는 대수학의 기본 정리(代數學 基本定理, fundamental theorem of algebra)에 대한 엄밀한 증명이 될 수 있다.

[대수학의 기본 정리]
제$n$차 방정식의 해는 복소 영역에서 $n$개 존재한다.

[증명]
제$n$차 방정식을 위한 제$n$차 다항식을 다음처럼 정의한다.

                  (6)

여기서 $a_n \ne 0$, $f(z)$ = $0$을 만족하는 $z$가 방정식의 해이다. 만약 $f(z)$가 모든 복소 영역에서 $0$이 아니라면, $f(z)$의 역수를 취한 $g(z)$[= $1/f(z)$]를 생각할 수 있다. 그러면 $g(z)$는 모든 영역에서 정칙이면서 유계이므로, 리우빌의 정리에 의해 $g(z)$는 상수 함수가 되어야 한다. 하지만 $a_n \ne 0$인 조건과 배치되는 결과라서 $f(z)$는 어떤 점에서 반드시 $0$이 되어야 한다. 어떤 점 $z_1$이 $f(z_1)$ = $0$을 만족하는 해라면, $f(z)$에서 $(z-z_1)$을 나눈 함수를 $f_1 (z)$ = $f(z)/(z-z_1)$라 둔다. 그러면 $f_1 (z)$는 제$n-1$차 다항식이 된다. 복소 함수 $f(z)$와 유사하게 $f_1 (z)$도 전영역에서 해를 가지지 않는다면, 리우빌의 정리에 의해 $f_1 (z)$는 상수가 되어야 한다. 이는 $f_1 (z)$가 제$n-1$차 다항식이라는 조건에 모순이다. 이러한 적용을 계속 반복하면, $f(z)$는 복소 영역에서 $n$개의 해를 반드시 가진다.
______________________________

방정식의 해 혹은 근은 복소 함수론에서 영점(零點, zero: 함수가 $0$으로 수렴하는 점)이라 부른다. 영점은 극점(極點, pole: 함수가 무한대로 발산하는 점)에 대비되는 개념이다. 복소 함수 $f(z)$의 역수 함수를 $g(z)$라 하면, $f(z)$의 영점은 $g(z)$의 극점 혹은 $f(z)$의 극점은 $g(z)$의 영점이다.

2020년 9월 18일 금요일

사원수와 회전(Quaternion and Rotation)

[경고] 아래 글을 읽지 않고 "사원수와 회전"을 보면 바보로 느껴질 수 있습니다.


사원수(四元數, quaternion)는 물리학에 벡터(vector)란 개념을 선물해준 고마운 존재이다. 사원수가 아름답기 때문에 그동안 수많은 찬사를 받았지만, 사원수는 치명적인 약점이 존재한다. 우리 직관을 너무 벗어난 사원수의 대수적 특성은 사실 멀리하고 싶은 그리움이다. 우리만 이런 딜레마를 느낄까? 당연히 아니다. 사원수가 널리 퍼진 19세기말부터 대수 기반의 사원수가 아닌 직관적인 벡터 개념을 만들기 위한 경쟁이 시작되었다. 사원수에 벡터가 이미 포함되어 있었지만, 완벽한 사원수 대수를 버리고 어딘지 부실하게 벡터만 강조한 벡터 해석학(vector analysis)이 1881년기브스 42세, 조선 고종 시절에 출현했다[1], [2]. 벡터 해석학은 미국 최초의 공학 박사이자 예일대학교(Yale University) 교수인 기브스Josiah Willard Gibbs(1839–1903)가 만들었다. 사원수라는 엄밀한 수학 체계를 어려워하는 예일대 학생들을 위해 벡터 개념을 간단히 사용할 수 있도록 기브스 교수는 좌표계 기반 벡터에 대한 자체 교재를 만들어서 학생들을 가르쳤다. 좌표계 기반 벡터 교재는 영국에 있는 헤비사이드Oliver Heaviside(1850–1925)에게 1888년 무렵에 전달되었다. 깐깐한 헤비사이드가 기브스의 벡터 개념을 칭찬했지만, 헤비사이드는 이미 1884년헤비사이드 34세, 조선 고종 시절에 기브스와는 독립적으로 사원수로 기술된 맥스웰 방정식(Maxwell's equations)을 자신만의 벡터 기반 맥스웰 방정식으로 간략화했다. 그뒤 기브스는 너무 바빠서 새로운 벡터 개념을 다듬을 시간이 없었지만, 기브스의 제자인 윌슨이 벡터 해석학[3]이란 멋진 책을 써서 1901년에 출판했다. 이 교재로 인해 사원수라는 개념은 물리학자의 손을 떠나 원래 있어야 할 수학자에게 돌아갔다. 요즘 물리학자는 기브스와 헤비사이드가 제안한 좌표계 기반 벡터를 사용해서 사물의 움직임을 계산한다.

[그림 1] 회전축 $\hat e$에 대한 3차원 공간의 회전(출처: wikipedia.org)

[그림 2] 공간 회전의 사원수 표현식을 위한 3차원 좌표계(출처: wikipedia.org)

사원수는 수학자가 발견한 교환 법칙이 성립하지 않는 최초의 대수 체계라서 존재 가치가 분명히 있다. 하지만 사원수와 경쟁하는 벡터 개념이 너무 직관적이라서 사원수는 다수의 사랑을 다시 받기는 어렵다. 그럼에도 불구하고 3차원 회전 연산(rotation operation)만 보면, 사원수의 회전 표현식이 벡터나 행렬 공식보다 확실히 아름답다.

[3차원 공간 회전을 위한 사원수 표현식]
단위 벡터 $\hat k$를 회전축으로 놓고, 3차원 벡터 $\bar v$를 $\theta$만큼 회전시킨 벡터 $\bar v_\text{rot}$는 다음처럼 표현된다.

                  (1)

여기서 단위 벡터(unit vector)는 스칼라가 $0$인 단위 사원수(unit quaternion)이다.

[증명]
사원수 표현식을 증명하기 위한 사원수의 벡터 항등식은 다음과 같다.

                              (2)

                              (3)

회전을 표현하기 위해 사용한 사원수 $\bf q$의 크기는 $\theta$에 관계없이 항상 $1$이다.

                              (4)

그래서 식 (1)처럼 ${\bf q}^{-1}$ = ${\bf q}^*$이 성립한다. 벡터 $\bar v_\text{rot}$의 크기도 $|{\bf q} \bar v{\bf q}^*|$ = $|{\bf q}| |\bar v| |{\bf q}^*|$ = $|\bar v|$처럼 보존된다. 사원수 $\bf q$를 식 (1)에 대입한 후, 식 (2)와 (3)을 이용해 정리한다.

                              (5)

여기서 삼각 함수의 합차 공식에 의해 $\cos^2 (\theta/2) - \sin^2 (\theta/2)$ = $\cos \theta$, $2 \sin^2 (\theta/2)$ = $1 - \cos \theta$이다. 식 (5)는 로드리그의 회전 공식(Rodrigues' rotation formula)과 동일하므로, 3차원 공간의 회전 표현식이 증명된다.
______________________________

공간 회전에 대한 사원수 표현식을 증명할 때, 로드리그의 회전 공식과 비교한 부분이 약간 어색해 보일 수 있다. 하지만 역사적으로 보면, 로드리그의 회전 공식이 나온 직후에 사원수가 제안되었으므로 우리의 접근 방식은 타당하다. 

[참고문헌]
[1] M. J. Crowe, "A history of vector analysis," University of Louisville, 2002. (방문일 2020-09-18)
[2] E. B. Wilson, "Reminiscences of Gibbs by a student and colleague," Bull. Amer. Math. Soc., vol. 37, no. 6, pp. 401–416, 1931.