Babylonian Method, Heron's Method
바빌로니아의 제곱근 알고리즘은 $\sqrt S$의 값을 초기값으로 설정된 임의의 양수 $x_0$의 값으로 다음과 같은 계산을 반복해 나가는 알고리즘으로 알려져 있다.
$$x_{n+1} = \frac 1 2 \left( x_n +\frac S {x_n} \right)$$
계산해 보면 알겠지만, 이 방법으로 $\sqrt S$의 값을 상당히 빠르게 얻어낼 수 있다. 하지만, 많은 수학사 책을 통해 살펴보면 바빌로니아 시대에는 제곱근 계산에 대한 기록은 있으나, 위와 같은 계산법은 아직 발견된 바가 없다. 바빌로니아의 유물로 알려진 것 중 제곱근에 관련해서 가장 유명한 것은 다음 점토판이다.
이 점토판의 중앙에는 $\sqrt 2$의 근삿값이 60분법으로 기록되어 있으며, 정밀도가 꽤 높다. 하지만, 바빌로니아인이 제곱근 계산에 사용한 방법을 현대 수학 기호를 사용하여 표현하면 다음과 같다:
$$\sqrt{a^2 +b^2} \approx a + \frac {b^2}{2a}$$
이 글 시작 부분에 언급한 알고리즘은 현재까지 확인된 바로는 AD 60년 경 그리스 수학자이자 엔지니어였던 헤론(Heron of Alexandria)이 그의 책 Metrica에 소개한 것이 처음이다.[^1] 헤론의 Metrica는 넓이와 부피를 구하고 이 양들을 다양하게 사용하는 방법에 대해 총 3권에 걸쳐 다룬 책으로, 그 책 1권에는 세 변이 주어진 삼각형의 넓이를 구하는, 이른바 헤론의 공식이 적혀 있다.
만약 $A$가 세 변이 $a$, $b$, $c$인 삼각형의 넓이라 하고 $s= \frac 1 2 (a+b+c)$라 하면 다음 식이 성립한다: $$ A^2 = s(s-a)(s-b)(s-c)$$
이 공식은 아르키메데스가 이미 알고 있었던 것을 헤론이 기록해 놓은 것이다.[^2] 헤론이 이 공식을 소개하는 부분에는 $A^2$에서 $A$를 얻어내는 과정이 설명되어 있는데, 이 과정이 여기서 말하려는 바빌로니아 방법, 혹은 헤론의 방법이다. 또한, 헤론은 Metrica 3권에서 세제곱근을 구하는 독특한 알고리즘도 제시했는데, 그 내용은 여기서 논하지 않는다.
헤론 이후에 제곱근을 구하는 여러 가지 방법이 개발되었는데, 뉴턴도 미분을 활용해 계산하는 방법을 개발하게 된다. 그런데, 뉴턴이 개발한 방법을 제곱근에 적용하면 헤론의 방법과 동일한 식이 얻어진다. 헤론은 무슨 생각으로 그런 알고리즘을 알아냈을까?
헤론의 아이디어가 기록된 것은 없지만, 헤론이 만든 식을 잘 뜯어보면 기하학적 의미가 숨어있음을 알 수 있다. 기하학적으로 알고리즘을 해석하면 다음과 같다.
- 넓이가 $S$인 직사각형을 임의로 설정하고 한 변의 길이를 $x_1$이라 하자. 그러면 다른 변의 길이는 $\dfrac {S} {x_1}$이 된다.
- 이제 두 변의 산술평균을 구해 새로운 사각형의 한 변의 길이 $x_2$로 설정한다. 그러고 넓이를 유지하기 위해 남은 한 변의 길이를 $\dfrac S {x_2}$로 잡는다.
- 2번 과정을 정사각형에 충분히 가까워질 때까지 반복한다.
위 그림은 헤론의 방법을 통해 $\sqrt 2$의 값을 초기값 $x_1 =4$를 가지고 계산하는 과정을 그림으로 나타낸 것이다. 사각형의 모양이 금방 정사각형에 가까워지고 있음을 알 수 있다.
산술평균, 기하평균, 조화평균
다시 헤론의 방법을 나타내는 식을 보자:
$$x_{n+1} = \frac 1 2 \left( x_n + \frac S {x_n} \right) $$
$x_n >0$이라 할 때, $x_{n+1}$은 $x_n$과 $\dfrac S {x_n}$의 산술평균이고, $\sqrt S$는 기하평균이다. 또한, $x_{n+1} \times \dfrac S {x_{n+1}} = S$이므로 $\dfrac S {x_{n+1}}$은 $x_n$과 $\dfrac S {x_n}$의 조화평균이다. (산술평균과 조화평균의 곱은 기하평균의 제곱이다.) 이 해석을 통해 $x_n < \sqrt S$이면(반대의 경우도 비숫함)
$$ x_n < x_{n+1} < \sqrt S < \dfrac S {x_{n+1}} < \frac S {x_n} $$
임과
$$ \left| x_{n+1} - x_n \right| \leq \left(\frac 1 2 \right)^{n-1} \left| x_1 - \frac S {x_1} \right|$$
을 직관적으로 이해할 수 있고, 수렴성도 간단히 설명 가능하다.
세제곱근으로 확장하면...
이 아이디어를 $S$의 세제곱근으로 확장해보자. 초기 설정을 어떻게 하는가에 따라 두 가지 방법을 생각할 수 있다.
만약 초깃값으로 $x_1$, $x_2$ 두 수를 설정한다면 다음과 같은 알고리즘을 얻을 수 있다.
- 세 변의 길이가 $x_1$, $x_2$, $\dfrac S {x_1 x_2}$로 부피가 $S$인 직육면체 세 변의 산술평균을 새로운 변의 길이 $x_3$으로 한다.
- 이제, 부피를 유지하기 위해 세 변의 길이가 $x_2$, $x_3$, $\dfrac S {x_2 x_3}$인 직육면체를 생각하고 이 세 변의 산술평균을 $x_4$라 한다.
- 2번과 같은 작업을 반복한다.
이 알고리즘에서 나오는 식은 다음과 같다.
$$ x_{n+2} = \frac 1 3 \left( x_{n+1} + x_n + \frac S {x_{n+1} x_n} \right)$$
한편, 초깃값으로 $x_1$ 하나만을 이용한다면 다음과 알고리즘을 얻을 수 있다.
- 세 변의 길이가 $x_1$, $x_1$, $\dfrac S {(x_1)^2}$로 부피가 $S$인 직육면체 세 변의 산술평균을 새로운 변의 길이 $x_2$라 한다.
- 이제, 부피를 유지하기 위해 세 변의 길이가 $x_2$, $x_2$, $\dfrac S {(x_2)^2}$인 직육면체를 생각하고 이 세 변의 산술평균을 $x_3$이라 한다.
- 2번과 같은 작업을 반복한다.
이 알고리즘에서 나오는 식은 다음과 같다.
$$x_{n+1} = \frac 1 3 \left( 2x_n + \frac S {(x_n)^2} \right) $$
어느 방법이 더 빠르게 수렴할까? 다음 표는 $\sqrt[3]{2}$를 구한 과정을 표로 나타낸 것이다. 새로운 값의 계산은 n이 2일 때부터 시작된다. 이 경우를 보면 후자가 더 빠름을 알 수 있다.
뉴턴의 방법으로 식을 만들어보면 두 번째 방법의 식이 얻어진다.
이제, $n$-제곱근
차원을 더 높여서 생각하면 $n$-제곱근은 다음과 같은 식을 통해 얻어내는 것이 합리적이다.
$$x_{n+1} = \frac 1 n \left( (n-1)x_n + \frac S {(x_n)^{n-1}} \right)$$
[^1]: Heath, Sir Thomas L. (1921). A History of Greek Mathematics, Vol. 2. Oxford: Clarendon Press. pp. 323–324.
[^2]: G Deslauriers and S Dubuc, Le calcul de la racine cubique selon Héron, Elem. Math. 51 (1) (1996), 28-34.