ベアストゥ法
ベアストゥ法（Bairstow’s method）は、高次多項式方程式 $P_n(x) = 0$ の根（解）を求める数値計算法の１つです。
最大の特徴は、実数係数の多項式から実数根だけでなく、複素根も含めて効率よく計算できる点にあります。実数係数の多項式において複素根は必ず共役のペア（$a \pm ib$）で出現することを利用し、２次式の多項式で割る除算を繰り返し行います。
基本的な仕組み
多項式 $P_n(x)$ を、試行的な２次式 $x^2-px-q$ で割ることを考えます。
$$P_n(x) = (x^2-px-q) Q_{n-2}(x) + R(x)$$
ここで、余り $R(x)$ は高々１次式となるため、$R(x) = fx+g$ の形に書くことができます。もし、$x^2-px-q$ が $P_n(x)$ の完全な約数であれば、余り $R(x)$ の係数はゼロ（ $f=0$ かつ $g=0$ ）となります。
ベアストゥ法は、この条件を満たすパラメータ $(p, q)$ を求める数値計算アルゴリズムです。多項式から２次の多項式で因数分解し、これを繰り返すことで、次数を順次下げることで全ての根を求めていきます。
メリットとデメリット
- メリット
- 全ての計算を実数範囲で完結できる。複素数演算を直接使わずに、複素数根を求めることができます。
- 収束速度が速い。ニュートン法をベースにしているため、初期値が適切であれば２次収束します。
- デメリット
- 初期値依存性。初期値の取り方によっては収束しなかったり、振動・発散することがあります。
- 重根に弱い。重根や重根に近い根が存在する場合、収束が遅くなったり精度が落ちたりします。
アルゴリズム
ベアストゥ法のアルゴリズム（手順）は以下になります。
多項式の合成除算
設定した $(p, q)$ を使って多項式の係数から数列 $b_k$ を順次計算し、余りの項 $f, g$ を求めます。これは、次のような一般的な $n$ 次の多項式について、
$$a_nx^n+a_{n-1}x^{n-1}+a_{n-2}x^{n-2}+\cdots+a_1x+a_0=0 -①$$
以下のように２次式で除算します。
$$(b_{n-2}x^{n-2}+\cdots+b_1x+b_0)(x^2+px+q)+fx+g=0 -②$$
②を展開すると、各係数は①と等しくなるため、
$a_n=b_{n-2}$
$a_{n-1}=b_{n-3}+pb_{n-2}$
$a_{n-2}=b_{n-4}+pb_{n-3}+qb_{n-2}$
・・・・
$a_2=b_0+pb_1+qb_2$
$a_1=pb_0+qb_1+f$
$a_0=qb_0+g$
$a_0\sim a_n$ は既知であるため、これらより $b_0\sim b_{n-2}$ を求めると、$f$ と $g$ は $p$、$q$ の関数として表すことができます。
$f=f(p,q)$
$g=g(p,q)$
n=4 の場合の計算例
$n=4$ の場合の②は、
$$(b_2x^2+b_1x+b_0)(x^2+px+q)+fx+g=0$$
これと①と同じ次数同士の係数を比べると、
$a_4=b_2 -(1)$
$a_3=pb_2+b_1 -(2)$
$a_2=qb_2+pb_1+b_0 -(3)$
$a_1=pb_0+qb_1+f -(4)$
$a_0=qb_0+g -(5)$
(1)～(3) から $b_0\sim b_2$ を求めると、
$b_2=a_4$
$b_1=a_3-a_4p$
$b_0=a_2-a_4q-a_3p+a_4p^2$
これより (4)と(5) から $f$ と $g$ を求めると、
$f=a_1-a_2p-a_3q+2a_4pq+a_3p^2+a_4p^3$
$g=a_0-a_2q+a_4q^2+a_3pq-a_4p^2q$
初期値の設定
２次多項式の係数 $(p, q)$ の初期値を適当に設定します。例えば、 $p=0, q=0$ です。
ニュートン法による近似更新
剰余 $fx+g$ を０にするような $p$、$q$ を求める問題は、$f$ と $g$ の根 $\alpha$、$\beta$ を求める問題に帰着します。
$f(\alpha,\beta)=0$
$g(\alpha,\beta)=0$
これらの根をニュートン法を使って解きます。（$p_1$、$q_1$）の周りでテイラー展開すると、
$$f(p_2,q_2)=f(p_1,q_1)+(p_2-p_1)\frac{\partial f}{\partial p}+(q_2-q_1)\frac{\partial f}{\partial q}+\cdots$$$$g(p_2,q_2)=g(p_1,q_1)+(p_2-p_1)\frac{\partial g}{\partial p}+(q_2-q_1)\frac{\partial g}{\partial q}+\cdots$$
右辺の２次の項を無視し、左辺を０、および、
$$f_p\equiv\frac{\partial f}{\partial p} , f_q\equiv\frac{\partial f}{\partial q}$$$$g_p\equiv\frac{\partial g}{\partial p} , g_q\equiv\frac{\partial g}{\partial q}$$
と置いて解くと、
$$p_2=p_1-\frac{g_qf(p_1,q_1)-f_qg(p_1,q_1)}{f_pg_q-f’\dot{g}} -③$$$$q_2=q_1+\frac{g_pf(p_1,q_1)-f_pg(p_1,q_1)}{f_pg_q-f_qg_p} -④$$
以上により、新しい（$p_2$、$q_2$）が求められます。
係数の収束判定
③と④を繰り返すことにより（$p_2$、$q_2$）$\cdots$（$p_n$、$q_n$）を計算し、変化量が許容誤差範囲 $\epsilon$ 未満になれば、（$p_n$、$q_n$）を根と考えます。
$|p_n-p_{n-1}|\lt\epsilon_p$ かつ、$|q_n-q_{n-1}|\lt\epsilon_q$
ならば、
$fx+g\simeq0$
従って、元の多項式①は以下のように書き換えることができます。
$$(b_{n-2}x^{n-2}+\cdots+b_1x+b_0)(x^2+px+q)=0$$
根の算出と次数ダウン
剰余 $fx+g$ が０になるような $p$、$q$ を見つけることができれば、①の根を求める問題は、次の２つの多項式の根を求める問題に帰着します。
$$b_{n-2}x^{n-2}+\cdots+b_1x+b_0=0 -⑤$$$$x^2+px+q=0 -⑥$$
⑥は解の公式を使って根を求めることができます。
$$ \alpha,\beta = \frac{-p\pm\sqrt{p^2-4q}}{2} $$
一方、⑤が高次の多項式の場合は、この手順を繰り返し、次数を下げていきます。
