CORDIC(Volder's Algorithm)
자주 있는 일은 아니지만, 살면서 적어도 한 번쯤은 삼각함수를 사용하는 코드를 작성해야 할 일이 있을 것입니다. C++의 math.h 헤더나 Python의 math 모듈같이, 대부분의 언어에서 삼각함수를 비롯한 여러 기능을 지원하는 수학 라이브러리가 기본으로 제공됩니다.
>>> from math import *
>>> sin(1) # usage of sin function in Python
0.8414709848078965
그런데, 주어진 각도가 $ \pi/6 $, $\pi/4$, $\pi/3$과 같이 딱 떨어지는 특수각이 아닐 때에도 컴퓨터는 어떻게 삼각함수의 값을 구할 수 있는 걸까요? 위 예시에서 $\sin 1$의 값은 도대체 어떤 과정을 통해 얻어진 걸까요?
근사값 구하기
사실 특수각이 아닌 경우 삼각함수의 값을 정확히 구하기란 불가능합니다. 대신, 여러가지 근사 기법을 이용해 실제 값에 매우 가까운 근사값을 계산해내는 것입니다. 대표적인 예시로, 수학에 관심이 있는 분이라면 아래와 같은 테일러 급수를 이용해서 삼각함수의 값을 꽤 정확하게 근사해낼 수 있다는 것을 알고 계실 것입니다.
$ \sin x = x - \frac{x^3}{3!} + \frac{x^5}{5!} - \frac{x^7}{7!} + \cdots $
$ \cos x = 1 - \frac{x^2}{2!} + \frac{x^4}{4!} - \frac{x^6}{6!} + \cdots $
자세한 설명은 생략하겠지만, 위 급수는 빠르게 수렴하며 실제로 사용하기에도 무리가 없습니다. 하지만 꼭 이것만이 유일한 방법인 것은 아닙니다.
이번 글에서는 삼각함수의 값을 구할 때 쓰이는 또 다른 전통적인 근사 알고리즘인 CORDIC에 대해 알아보도록 하겠습니다.
CORDIC이란 COordinate Rotation DIgital Computer의 약자로, 말 그대로 좌표의 회전을 이용해 값을 계산해내는 알고리즘입니다. 좌표의 회전이라는 특성 상, CORDIC에서는 2차원에서의 회전 변환이 중요한 역할을 합니다.
회전 변환 행렬
어떤 2차원 벡터 $(x,y)$를 원점 기준 반시계 방향으로 $\theta$만큼 회전시키는 행렬 $R_{\theta}$는 아래와 같이 표현됩니다. (이 글에서 각도 $\theta$는 항상 라디안 단위입니다)
$ R_{\theta} = \begin{bmatrix} {\cos \theta} & {-\sin \theta} \ {\sin \theta} & {\cos \theta} \end{bmatrix} $
예를 들어, 벡터 $(2,4)$를 $\pi/4$만큼 회전시킨 결과는 아래처럼 회전 변환 행렬을 왼쪽에 곱해서 구할 수 있습니다.
\[R_{\frac{\pi}{4}} \begin{bmatrix} 2 \\ 4 \end{bmatrix} = \begin{bmatrix} \cos {\frac{\pi}{4}} & -\sin {\frac{\pi}{4}} \\ \sin {\frac{\pi}{4}} & \cos{\frac{\pi}{4}} \end{bmatrix} \begin{bmatrix} 2 \\ 4 \end{bmatrix} = \begin{bmatrix} \frac{\sqrt{2}}{2} & -\frac{\sqrt{2}}{2} \\ \frac{\sqrt{2}}{2} & \frac{\sqrt{2}}{2} \end{bmatrix} \begin{bmatrix} 2 \\ 4 \end{bmatrix} = \begin{bmatrix} -\sqrt{2} \\ 3\sqrt{2} \end{bmatrix}\]또한, 어떤 점을 $\theta_1$만큼 회전시킨 후 다시 $\theta_2$만큼 회전시킨다면, 결과적으로 $\theta_1 + \theta_2$만큼 회전시킨 것이 되므로 회전 변환 행렬끼리의 곱으로 또 다른 회전 변환 행렬을 나타낼 수 있습니다.
$ R_{\theta_1} \cdot R_{\theta_2} = R_{\theta_1 + \theta_2} $
회전 변환 행렬에 대한 더 자세한 정보가 궁금하시다면 위키백과를 참조하시길 바랍니다.
아이디어
먼저, 편의를 위해 구하는 각도가 $[-\pi/2, \pi/2]$의 범위에 있다고 가정합시다. 삼각함수의 주기성에 의해, 모든 각도에 대한 삼각함수의 값은 이 범위 안에 있는 삼각함수의 값으로 변환될 수 있으므로 일반성을 잃지 않습니다.
하나의 예시로, $ \cos 1 $, $ \sin 1 $의 값을 근사하고 싶다고 해봅시다. 만약 벡터 $(1,0)$를 여러 번 회전시켜 각도가 $1$인 벡터, 즉 $( \cos 1, \sin 1 )$에 아주 비슷하게 만들 수 있다면, 해당하는 회전 결과를 회전 변환 행렬로 직접 계산함으로써 $\cos 1$과 $\sin 1$의 근사값을 얻을 수 있을 것입니다.
이제, 다음과 같은 아이디어를 생각할 수 있습니다.
위와 같이 초기 벡터 $v_0 = (1,0)$을 잡습니다. 이 때, 위에서 $target$으로 표시된 벡터 $( \cos 1, \sin 1 )$에 최대한 가깝게 회전시키는 것이 목표입니다.
$v_0$의 각도가 목표 벡터보다 작으므로, $\frac{\pi}{4}$만큼 회전시킵시다.
그러면 새로운 벡터 $v_1 = R_{\frac{\pi}{4}} v_0 $가 됩니다.
이제 $v_1$의 각도를 보면, 여전히 목표 벡터보다 작으므로 $\frac{\pi}{8}$만큼 회전시킵시다.
그러면 새로운 벡터 $v_2 = R_{\frac{\pi}{8}} v_1$이 됩니다.
이제 $v_2$의 각도가 목표 각도보다 더 커졌으므로, 이번에는 $-\frac{\pi}{16}$만큼 회전시킵시다.
그러면 새로운 벡터 $v_3 = R_{-\frac{\pi}{16}} v_2$이 됩니다.
이 과정을 계속 반복하다 보면, $v_n$의 각도가 원하는 각도 $\theta$에 매우 가까워질 것임을 직관적으로 알 수 있습니다. 따라서,
\[\begin{bmatrix} \cos 1 \\ \sin 1 \end{bmatrix} \approx R_{\frac{\pi}{4}} R_{\frac{\pi}{8}} R_{-\frac{\pi}{16}} R_{\frac{\pi}{32}} \cdots \begin{bmatrix} 1 \\ 0 \end{bmatrix}\]처럼 근사할 수 있으며, 단계를 더 진행하면 진행할수록 더 정확한 값을 얻을 수 있습니다.
이 아이디어를 사용하면, $R_{\frac{\pi}{4}}$, $R_{\frac{\pi}{8}}$, $R_{\frac{\pi}{16}}$, $\cdots$, $R_{-\frac{\pi}{4}}$, $R_{-\frac{\pi}{8}}$, $R_{-\frac{\pi}{16}}$, $\cdots$ 에 해당하는 행렬만 알고 있다면 어떤 각도 $\theta$에 대해서도 효과적으로 $\cos \theta$와 $\sin \theta$를 근사할 수 있습니다.
CORDIC의 구현
실제 CORDIC은 위와 비슷한 아이디어를 사용하지만, 계산의 편의와 효율을 향상하기 위해 더하거나 빼는 각도가 살짝 다릅니다. 그 수학적 배경은 다음과 같습니다.
삼각함수에서 $\theta$가 $[-\pi/2, \pi/2]$에 있을 때, 아래와 같은 항등식이 성립합니다.
$ \cos \theta = \frac{1}{\sec \theta} = \frac{1}{\sqrt{1 + \tan^2 \theta}} $
$ \sin \theta = \tan \theta \cos \theta = \frac{\tan \theta}{\sqrt{1 + \tan^2 \theta}} $
이를 이용하면, 회전 변환 행렬 $R_{\theta}$를 다음과 같이 나타낼 수 있습니다.
\[R_{\theta}= \begin{bmatrix} {\cos \theta} & {-\sin \theta} \\ {\sin \theta} & {\cos \theta} \end{bmatrix} = \frac{1}{\sqrt{1 + \tan^2 \theta}} \begin{bmatrix} 1 & -\tan \theta \\ \tan \theta & 1 \end{bmatrix}\]이 때, $i$번째 단계에서 더할 각도 $\gamma_i$를 $\tan \gamma_i = 2^{-i}$가 되도록 설정하면, $i$번째 단계에서 곱해야 할 회전 변환 행렬 $R_i$는
$$ R_i = R_{\pm \gamma_i} = \frac{1}{\sqrt{1 + \tan^2 \gamma_i}} \begin{bmatrix} 1 & \mp \tan \gamma_i \ \pm \tan \gamma_i & 1 \end{bmatrix} = \frac{1}{\sqrt{1 + 2^{-2i}}} \begin{bmatrix} 1 & \mp 2^{-i} \ \pm 2^{-i} & 1 \end{bmatrix} $
처럼 표현됩니다. 이 행렬들을 서로 곱할 때, 실제로 행렬 곱셈에 관여하는 부분은 $2^{-i}$로 컴퓨터가 오차 없이 정확히 표현할 수 있는 수인 동시에 계산도 특별한 과정 없이 쉽게 가능하므로, 위에서 제시한 방법보다 훨씬 효율적입니다.
이 경우, 알고리즘은 $i$번째 단계에서 $\frac{\pi}{2^{i+2}}$ 대신 $\gamma_i$ 크기만큼 회전시키는 것으로 바뀌게 되며, 이 때 각 $\gamma_i = \arctan(2^{-i})$의 값은 미리 구해서 저장해 놓은 값을 사용합니다.
이 때, $\gamma_i$의 무한합이 $1.74329\cdots$로 수렴하기 때문에 $\gamma_i$의 합으로 표현 가능한 각도의 범위는 $(-1.74329\cdots, 1.74329\cdots)$이고, 이는 $[\pi/2, \pi/2]$를 포함하기 때문에 $\gamma_i$를 사용해도 원하는 범위에서는 문제가 없음을 알 수 있습니다.
회전 변환을 어느 방향으로 적용할지는 마찬가지로 $v_n$의 각도가 목표 각도 $\theta$ 보다 큰지 작은지에 따라 결정되며, $i$번째 단계에서 $\gamma_i$를 어느 방향으로 회전할 지를 $\sigma_i$를 이용하여 나타냅시다. (반시계 방향일 경우 $\sigma_i = 1$, 시계 방향일 경우 $\sigma_i = -1$) 또한, 앞에 붙은 계수를 $K_i$로 나타내면 $i$번째에 적용할 회전 변환 행렬을 아래처럼 간단하게 나타낼 수 있습니다.
$ R_i = \frac{1}{\sqrt{1 + 2^{-2i}}} \begin{bmatrix} 1 & \mp 2^{-i} \ \pm 2^{-i} & 1 \end{bmatrix} = K_i \begin{bmatrix} 1 & - \sigma_i 2^{-i} \ \sigma_i 2^{-i} & 1 \end{bmatrix}$
결과적으로, $(\cos \theta , \sin \theta)$는 다음과 같은 식으로 근사 가능합니다.
\[\begin{bmatrix} \cos \theta \\ \sin \theta \end{bmatrix} \approx R_n \cdots R_1 R_0 \begin{bmatrix} 1 \\ 0 \end{bmatrix}\] \[=K_n \begin{bmatrix} 1 & - \sigma_i 2^{-n} \\ \sigma_n 2^{-n} & 1 \end{bmatrix} \cdots K_1 \begin{bmatrix} 1 & - \sigma_1 /2 \\ \sigma_1 /2 & 1 \end{bmatrix} K_0 \begin{bmatrix} 1 & - \sigma_0 \\ \sigma_0 & 1 \end{bmatrix} \begin{bmatrix} 1 \\ 0 \end{bmatrix}\] \[=K_n \cdots K_1 K_0 \begin{bmatrix} 1 & - \sigma_i 2^{-n} \\ \sigma_n 2^{-n} & 1 \end{bmatrix} \cdots \begin{bmatrix} 1 & - \sigma_1 /2 \\ \sigma_1 /2 & 1 \end{bmatrix} \begin{bmatrix} 1 & - \sigma_0 \\ \sigma_0 & 1 \end{bmatrix} \begin{bmatrix} 1 \\ 0 \end{bmatrix}\]이 때, 각각의 n에 대해
\[K_n \cdots K_1 K_0 = \prod_{i=0}^{n} K_i = \prod_{i=0}^{n} \frac{1}{\sqrt{1+2^{-2i}}}\]값 역시 미리 구해서 저장해 놓으면 더 효율적으로 근사값을 구할 수 있습니다.
아래는 C++로 구현한 CORDIC 알고리즘 코드의 예시입니다.
void cordic(double target_angle, int n){ // works for -PI/2 <= target_angle <= PI/2, n < 28
double gamma[28]={
0.78539816339745, 0.46364760900081, 0.24497866312686, 0.12435499454676,
0.06241880999596, 0.03123983343027, 0.01562372862048, 0.00781234106010,
0.00390623013197, 0.00195312251648, 0.00097656218956, 0.00048828121119,
0.00024414062015, 0.00012207031189, 0.00006103515617, 0.00003051757812,
0.00001525878906, 0.00000762939453, 0.00000381469727, 0.00000190734863,
0.00000095367432, 0.00000047683716, 0.00000023841858, 0.00000011920929,
0.00000005960464, 0.00000002980232, 0.00000001490116, 0.00000000745058
}; // pre-processed gamma_i
double K[28]={
0.70710678118655, 0.63245553203368, 0.61357199107790, 0.60883391251775,
0.60764825625617, 0.60735177014130, 0.60727764409353, 0.60725911229889,
0.60725447933256, 0.60725332108988, 0.60725303152913, 0.60725295913894,
0.60725294104140, 0.60725293651701, 0.60725293538591, 0.60725293510314,
0.60725293503245, 0.60725293501477, 0.60725293501035, 0.60725293500925,
0.60725293500897, 0.60725293500890, 0.60725293500889, 0.60725293500888,
0.60725293500888, 0.60725293500888, 0.60725293500888, 0.60725293500888
}; // pre-proecessed K[i] = K_0 * K_1 * ... * K_i
double current_angle = 0.0;
double x = 1.0;
double y = 0.0;
double pow2 = 1.0; // 1, 1/2, 1/4, ...
for(int i=0;i<=n;i++){
int sigma; // sign of rotation
if(current_angle < target_angle) sigma = 1;
else sigma = -1;
current_angle += sigma * gamma[i];
double nx = x - sigma * pow2 * y;
double ny = sigma * pow2 * x + y; // matrix multiplication
x = nx;
y = ny;
pow2 /= 2;
}
x *= K[n];
y *= K[n]; // scalar multiplication for K_0 * K_1 * ... * K_n
// (x,y) results in (cos(target_angle), sin(target_angle))
printf("cos(%lf) = %lf\n", target_angle, x);
printf("sin(%lf) = %lf\n", target_angle, y);
}
CORDIC의 활용
삼각함수의 값을 근사하는 다른 알고리즘들과 비교하여, CORDIC은 비교적 간단한 연산들로 구현됩니다. 특히, 고정 소수점을 사용한다면 곱셉을 전혀 사용하지 않고 덧셈, 뺄셈, 비트 시프트만으로 CORDIC 구현이 가능합니다. 다른 연산에 비해 더 많은 자원을 요구하는 곱셉을 사용하지 않고도 삼각함수의 값을 구할 수 있다는 것은 매우 큰 장점입니다. 비슷한 이유로, CORDIC은 소프트웨어로 코딩되는 경우보다 전기 회로 자체에서 CORDIC을 수행하도록 하드웨어 코딩이 되는 경우가 많다고 합니다. 예를 들어, 공학용 계산기에는 이 같은 회로가 탑재되어 있어 효율적으로 삼각함수의 값을 구할 수 있게 됩니다.
또한, 이번 글에서는 CORDIC을 이용해 $\sin \theta$와 $\cos \theta$, 그리고 이 둘을 나눠서 얻을 수 있는 $\tan \theta$까지 기본적인 삼각함수만을 구하는 법을 알아보았습니다. 하지만, CORDIC을 응용하면 $\arccos$, $\arcsin$, $\arctan$과 같은 역삼각함수를 비롯하여 다양한 초월함수의 값을 구할 수 있다고 합니다. 이에 대해서는 나중에 기회가 되면 다뤄보도록 하겠습니다. 마지막으로, 긴 글 읽어주셔서 감사합니다.