호너 방법

1. 개요

호너 방법 Horner's method
하는 일	다항식 $p(x)=\sum a_k x^k$ 의 값 계산, 일차식 $x-r$ 로 나누기
연산량	곱셈 $n$ 번 · 덧셈 $n$ 번 (일반 다항식에 최적)
오차	후방 안정 — 계수 상대 섭동 $\gamma_{2n}\approx 2nu$ 이내
이름	호너 (1819), 중국 진구소 (1247), 뉴턴·라그랑주도 사용
변형	조립제법, 2차 호너(도함수), 에스트린 방식, 보정 호너

괄호를 다르게 치는 것만으로 곱셈이 반으로 줄고 오차 분석이 깔끔해진다. 수치해석에서 이보다 가성비 좋은 아이디어는 드물다.

호너 방법(Horner’s method, Horner’s scheme)은 $n$ 차 다항식

$$p(x) = a_n x^n + a_{n-1}x^{n-1} + \cdots + a_1 x + a_0$$

을 안쪽부터 괄호를 쳐서

$$p(x) = \Big(\cdots\big((a_n x + a_{n-1})x + a_{n-2}\big)x + \cdots + a_1\Big)x + a_0$$

로 계산하는 방법이다. 코드로는 b = a[n]; for k = n-1 … 0: b = b*x + a[k] 한 줄이다. 곱셈 $n$ 번, 덧셈 $n$ 번으로 끝나고, 거듭제곱 $x^k$ 를 따로 만들지 않으니 메모리도 스칼라 하나다.

이름은 1819년 왕립학회에 근사 근 찾기 논문을 낸 윌리엄 조지 호너에서 왔지만, 같은 계산 배치는 13세기 중국 진구소(秦九韶)의 『수서구장』에 이미 있고 뉴턴도 썼다. 호너의 원래 기여는 값 계산이 아니라 근을 한 자리씩 확정하며 다항식을 평행이동해 나가는 근 찾기 절차였다¹. 오늘날 이 이름은 거의 전적으로 위의 중첩 곱셈 평가와, 같은 점화식으로 하는 일차식 나눗셈(조립제법)을 가리킨다.

윌킨슨 다항식·파울하버 공식·반데르몽드 행렬이 호너를 도구로 썼다. 이 문서는 그 도구 자체 — 왜 최적인지, 오차가 어떻게 나는지, 그 오차를 어떻게 줄이는지, 그리고 병렬 하드웨어에서 왜 다른 괄호가 이기는지를 다룬다.

2. 조립제법 — 같은 점화식의 두 얼굴

호너 점화식의 중간값 $b_k$ 를 버리지 않고 모으면 나눗셈이 된다. $b_n = a_n$, $b_k = b_{k+1}x_0 + a_k$ ($k = n-1,\dots,0$) 로 두면

$$p(x) = (x - x_0)\,q(x) + b_0,\qquad q(x) = b_n x^{n-1} + b_{n-1}x^{n-2} + \cdots + b_1$$

이다. 즉 $b_0 = p(x_0)$ 는 나머지 정리의 나머지이고, 나머지 $b_k$ 들은 몫 $q$ 의 계수다. 중학교에서 배우는 조립제법이 바로 이것이다.

여기서 두 가지가 공짜로 나온다.

• 도함수. $p'(x_0) = q(x_0)$ 이므로 $q$ 에 호너를 한 번 더 돌리면 $p'(x_0)$ 가 나온다. 두 점화식을 한 루프에 섞으면 $p$ 와 $p'$ 를 곱셈 $2n$ 번에 동시에 얻는다. 뉴턴-랩슨법으로 다항식 근을 찾을 때의 표준 내부 루프다. 반복하면 $p^{(j)}(x_0)/j!$, 곧 $x_0$ 에서의 테일러 계수 전부가 $O(n^2)$ 에 나온다 — 이게 호너가 1819년에 실제로 한 「평행이동」이다.
• 수축(deflation). 근 $r$ 을 하나 찾으면 $q(x) = p(x)/(x-r)$ 로 차수를 낮추고 나머지 근을 $q$ 에서 찾는다.

수축에는 순서 문제가 있다. 찾은 근 $\tilde r$ 은 근사값이라 몫 $q$ 의 계수에 오차가 섞이는데, 위의 점화식(높은 차수 계수부터 내려가는 전진 수축)은 절댓값이 작은 근부터 걷어낼 때 오차 증폭이 작다는 것이 윌킨슨의 고전적 결론이다. 실제로 근이 $1,2,\dots,10$ 인 10차 다항식에서 근 하나를 상대오차 $10^{-8}$ 로 틀리게 잡아 수축하면, 남은 9개 근의 최대 상대오차는 근 1을 걷어낸 경우 $2.5\times10^{-7}$, 근 10을 걷어낸 경우 $2.5\times10^{-6}$ 로 열 배 차이가 났다. 어느 순서든 실무의 처방은 같다 — 수축된 다항식에서 얻은 근을 원래 다항식 $p$ 에 대고 뉴턴 몇 번으로 다시 다듬는다(root polishing). 수축 오차가 다음 근으로 번지지 않게 끊는 장치다.

3. 왜 최적인가

$x^k$ 를 따로 계산해 더하는 순진한 방법은 거듭제곱을 누적해도 곱셈이 $2n$ 번 든다. 호너는 $n$ 번이다. 그럼 더 줄일 수는 없을까.

오스트롭스키(1954)는 일반 계수 다항식의 값을 구하는 데 덧셈 $n$ 번이 필요함을, 판(1966)은 곱셈·나눗셈 $n$ 번이 필요함을 증명했다. 즉 계수가 미리 정해지지 않은 다항식에 대해 호너는 연산 횟수 하한을 정확히 달성한다. 단서는 「계수를 미리 모른다」는 것이다. 같은 다항식을 수없이 평가한다면 계수를 미리 가공(preconditioning)해 곱셈을 약 $n/2$ 번까지 줄이는 방법이 있다 — 판(1966)의 결과가 이것까지 다룬다. 다만 가공된 계수는 수치적으로 나쁜 경우가 많아, 부동소수점에서는 거의 쓰지 않는다.²

4. 오차 — 후방 안정, 그러나 조건수는 그대로

호너는 후진 오차 해석의 교과서 예다. 각 단계의 곱셈·덧셈 반올림을 계수 쪽으로 밀어 넣으면, 계산값 $\hat p(x)$ 는 계수를 $|\Delta a_k| \le \gamma_{2n}|a_k|$ ($\gamma_m = mu/(1-mu)$, $u$ 는 단위 반올림)만큼 흔든 다항식의 정확한 값이다. 그래서 순방향 오차는

$$\frac{|\hat p(x) - p(x)|}{|p(x)|} \;\le\; \gamma_{2n}\,\frac{\sum_k |a_k||x|^k}{\left|\sum_k a_k x^k\right|} \;=\; \gamma_{2n}\;\mathrm{cond}(p, x)$$

로 묶인다. 분수 부분이 다항식 평가의 조건수다. 알고리즘은 잘못이 없고, 계수끼리 크게 상쇄하는 점에서는 문제 자체가 나쁘다는 뜻이다.

전형적인 시험대는 $(x-1)^n$ 을 펼친 계수 $\binom{n}{k}(-1)^{n-k}$ 를 근 근처에서 평가하는 것이다. 배정밀도에서 실제로 재면:

$n$	$x$	정확한 값	$\mathrm{cond}(p,x)$	호너 상대오차	보정 호너 상대오차
5	1.1	$1.0\times10^{-5}$	$4.1\times10^{6}$	$5.1\times10^{-11}$	0
10	1.1	$1.0\times10^{-10}$	$1.7\times10^{13}$	$2.4\times10^{-4}$	0
10	1.01	$1.0\times10^{-20}$	$1.1\times10^{23}$	$8.9\times10^{4}$	$5.7\times10^{-11}$
15	1.1	$1.0\times10^{-15}$	$6.8\times10^{19}$	$9.5\times10^{2}$	$1.2\times10^{-13}$
15	1.01	$1.0\times10^{-30}$	$3.5\times10^{34}$	$8.0\times10^{15}$	$8.5$

($x$ 는 해당 배정밀도 수 자체를 정확한 입력으로 보고 50자리 계산과 비교했다.) 호너의 오차는 거의 정확히 $\mathrm{cond}\times u$ 를 따라간다. $n=10$, $x=1.01$ 에서는 답의 부호조차 틀린다($-8.9\times10^{-16}$). 범인은 호너가 아니라 $(x-1)^n$ 을 단항식 계수로 적은 것이다. 이 다항식은 $(x-1)^n$ 그대로 계산하면 상대오차 $nu$ 수준으로 끝난다. 윌킨슨 다항식이 같은 교훈을 근 쪽에서 보여 준다.

같은 이유로 구간 해석에서도 호너 형태가 대개 펼친 형태보다 좁은 구간을 준다. $x$ 가 식에 덜 여러 번 등장하기 때문이다.

5. 보정 호너 — 정밀도를 두 배로

조건수가 큰 점에서 꼭 정확한 값이 필요하면 혼합 정밀도로 올리는 대신 오차 없는 변환(error-free transformation)을 쓸 수 있다. 부동소수점 곱 $a\cdot b$ 와 합 $a+b$ 는 각각 결과 $s$ 와 정확히 표현 가능한 오차 $e$ 로 쪼갤 수 있다(TwoProduct는 FMA 한 번 또는 데커 분할, TwoSum은 덧셈 6번). 그라야·랑글루아·루베(2005, 정리 2009)의 보정 호너(compensated Horner)는 호너의 각 단계에서 이 오차들을 모아 두 번째 다항식으로 같이 굴린 뒤 마지막에 더한다. 결과는

$$\frac{|\hat p_{\mathrm{comp}}(x) - p(x)|}{|p(x)|} \;\le\; u + \gamma_{2n}^2\,\mathrm{cond}(p,x)$$

로, 두 배 정밀도로 호너를 돌리고 배정밀도로 반올림한 것과 같은 정확도다. 위 표의 오른쪽 열이 그것이다. $\mathrm{cond}\lesssim 1/u \approx 10^{16}$ 이면 마지막 비트까지 맞고, $\mathrm{cond}\approx 10^{23}$ 에서도 열 자리가 산다. $\mathrm{cond}\approx10^{34}$ 처럼 $1/u^2$ 를 넘으면 보정 호너도 무너진다. 비용은 호너의 대략 서너 배 연산인데, 연산이 서로 독립적이라 파이프라인이 잘 채워져 실측 시간은 그보다 덜 늘어난다는 것이 저자들의 보고다.

6. 에스트린 방식 — 병렬에서는 괄호를 다시 친다

호너의 약점은 완전히 직렬이라는 것이다. $b_k$ 는 $b_{k+1}$ 이 나와야 계산할 수 있으니, 곱셈-덧셈 지연이 4사이클인 FMA 유닛에서 20차 다항식은 최소 80사이클을 기다린다. 유닛이 몇 개 놀고 있어도 소용없다.

에스트린(1960)의 방식은 계수를 짝지어 트리로 줄인다.

$$p(x) = (a_0 + a_1x) + (a_2 + a_3x)x^2 + \big((a_4 + a_5x) + (a_6 + a_7x)x^2\big)x^4 + \cdots$$

첫 층에서 일차식들을 동시에 만들고, 다음 층에서 $x^2$ 로, 그다음 $x^4$ 로 묶는다. 깊이가 $O(\log_2 n)$ 이 되는 대신 연산 수가 조금 늘고($x^{2^j}$ 를 따로 만든다) 레지스터를 더 쓴다. SIMD·GPU 의 초월함수 구현(지수·로그의 최소최대 근사 다항식)은 대개 호너와 에스트린을 섞은 형태로 짠다. 낮은 차수 몇 단계는 에스트린으로 병렬화하고, 나머지는 호너로 잇는 식이다. 오차 상한은 에스트린도 호너와 같은 꼴($\gamma_{c\log n}$ 정도의 계수 섭동)이라 정확도 손해는 거의 없다.

7. 단항식 밖의 호너

호너는 「$x$ 를 곱하고 계수를 더한다」는 구조라 단항식 기저에 묶여 있다. 그런데 앞의 표가 보여 주듯 단항식 계수 자체가 나쁜 표현인 경우가 많다. 다른 기저에도 같은 계열의 점화식이 있다.

• 뉴턴 형식. 보간 절점 $x_0,\dots,x_{n-1}$ 에 대한 분할차분 계수 $c_k$ 로 $p(x) = c_0 + (x-x_0)\big(c_1 + (x-x_1)(c_2 + \cdots)\big)$ 를 같은 중첩 곱셈으로 계산한다. 곱하는 인자만 $x - x_k$ 로 바뀐 호너다.
• 직교다항식. 체비쇼프 다항식 급수 $\sum c_k T_k(x)$ 는 세 항 점화식을 거꾸로 돌리는 클렌쇼 점화식으로 계산한다. 구조가 호너의 직교다항식판이고, 체비쇼프 계수는 대개 빠르게 감쇠해 조건수 문제가 훨씬 덜하다(클렌쇼-커티스 구적 참고). 일반 직교다항식의 세 항 점화식도 같은 방식이 통한다.
• 바리센트릭 형식. 다항식을 계수가 아니라 절점에서의 값으로 들고 다니면 바리센트릭 보간 공식이 $O(n)$ 에 평가하고, 체비쇼프 절점에서는 수치적으로 안정하다.

결론적으로 호너는 「단항식 계수가 주어졌을 때의 최선」이고, 그 계수가 주어진 것 자체가 문제라면 기저를 바꾸는 쪽이 근본 처방이다.

8. 실무 메모

• 다항식 값은 언제나 호너로 평가한다. a[k]*pow(x,k) 를 더하는 코드는 느리고 덜 정확하다³.
• 근 근처에서 부호 판정이 중요하면(근 찾기 종료 조건, 기하 술어) 위 오차 상한으로 동적 오차 한계를 같이 계산해 「부호를 믿을 수 있는가」를 판정한다. 합 $\sum |a_k||x|^k$ 도 호너로 한 번 더 돌리면 된다.
• 반복 근 찾기에서 수축을 쓰면 반드시 원래 다항식에 대고 다듬는다. 범용 근 계산은 오늘날 동반 행렬의 고유값이나 젠킨스-트라우브 알고리즘이 표준이고, 호너는 그 안의 부품으로 쓰인다.
• 계수가 정수이고 값을 정확히 원하면 정수 산술 호너가 가장 빠르다. 모듈러 산술의 다항식 해시(라빈-카프 등)가 정확히 이 루프다.

9. 관련 문서

• 윌킨슨 다항식 · 동반 행렬 · 뮬러 방법 · 뉴턴-랩슨법 · 젠킨스-트라우브 알고리즘
• 후진 오차 해석 · 조건수 · 부동소수점 연산 · 혼합 정밀도 · 구간 해석
• 체비쇼프 다항식 · 바리센트릭 보간 · 보간과 근사 · 반데르몽드 행렬 · 파울하버 공식
• 뉴턴 보간 · 클렌쇼 점화식 · 오차 없는 변환

10. Footnotes

1. 호너 본인은 이 방법의 우선권을 두고 당대에 이미 시비에 휘말렸다. 이탈리아의 파올로 루피니가 1804년에 비슷한 절차를 냈기 때문이다. 그래서 대륙 쪽 교과서에는 「루피니-호너」라고 적힌다. 그 루피니가 5차 방정식 근의 공식이 없다는 걸 (불완전하게) 처음 증명한 사람이라는 게 함정 카드. ↩

2. 「이 알고리즘이 연산 횟수 하한을 정확히 달성한다」는 증명이 붙은 실용 알고리즘은 생각보다 드물다. 행렬 곱셈은 아직도 지수 $\omega$ 를 두고 소수점 셋째 자리 싸움 중인데, 다항식 평가는 1966년에 이미 게임이 끝났다. ↩

3. pow(x, k) 는 일반 실수 지수용 함수라 내부에서 로그·지수를 돌리는 구현도 있다. 정수 차수 다항식에 쓰면 곱셈 몇 번이면 될 걸 초월함수 수십 번으로 갈아넣는 셈. 컴파일러가 알아서 바꿔 주길 기대하는 건 신에게 맡기는 수렴이다. ↩