数理解析研究所講究録 1295 巻 2002 年 209-212209

近似根の係数の計算について
北本卓也 1
山口大学教育学部

Abstract
2 変数多項式 $g(x, y)$ が与えられたとき、 その$x$に関するを根 $h(y)$ とする。 このとき、 ある条件下での解析関数であり、 のべき級数としてあらわすことが可能である。記号的ニュートン法を使えばこのべき級数の係数を効率よく計算できるが、記号的ニュートン法では、高次の項の係数の計算に低$h(y)$は$y$$y$次の項の係数が使われるため誤差の伝播が起こる可能性がある。 そこで、本稿では高次の項係数を低次の項の係数を用いずに直接計算する新しい算法を提案する。 この算法は効率の点で記号的ニュートン法に劣るので、計算誤差の確認といった目的に使用することが考えられる。

1

序論
著者は 2 変数多項式 $f(x, y)=0$ の展開を$k$次の近似根 (以下、 $\phi^{(k)}(y)$がべき級数展開出来るとき、 次までのべき級数と記す) として定義し、 その計算法を提案して来た。 これまでの算法は$x$に関する根$k$$\phi(y)$基本的に記号的ニュートン法と呼ばれるものの拡張であり、 $\phi^{(0)}(y)arrow\phi^{(p)}(y)arrow\phi^{(p^{2})}(y)arrow\cdotsarrow\phi^{(p^{k})}(y)$のように、近似根を少しずつ次数を上げながら計算を行う。 この算法は計算効率が良い反而、高次のべきを持つ項の係数の計算にそれより低次のべきを持つ項の係数が使われており、次数力吠きくなるに従い少しずつ計算誤差が高次の項に蓄積してしまうという性質を持っている。

そこで、本稿では高次の項を低次の項を用いずに直接計算する新しい算法を提案する。 この算法は効率の点で記号的ニュートン法に劣るので、 計算誤差の確認といった目的に使うことが考えられる。

2

記号的ニュートン法とその問題点$g(x,y)$ を $y$ 1 こ関してモニックな二変数多項式とする.

を $g(x, y)$ の$y$

このとき、 次式を満たす $m$ 次の多項式 $f^{(m)}(x)$} こ関する $m$ 次の近似根と呼ぶ.
$g(x, f^{(m)}(x))\equiv 0(\mathrm{m}\mathrm{o}\mathrm{d} x^{m+1})$(1)の次数 $m$ 次までべき級数展開したものに等しく、 次の記号的ニュートン法を用いて、 容易に計算することが可能である (詳細は $[5][6]$ を参照)。

近似根は正確な根を$x$

記号的ニュートン法(1) $g(0, y)=0$ の根の 1 つを、 $f^{(0)}(x)$ と置く.

.yamaguchi-u.ac$*\mathrm{k}\mathrm{i}\mathrm{t}\mathrm{a}\mathrm{m}\mathrm{o}\mathrm{t}\mathrm{o}@\mathrm{p}\mathrm{o}.\mathrm{c}\mathrm{c}$$.\mathrm{j}\mathrm{p}$XXXI-I

210(2) 次式を用いて、 $f^{(m)}(x)$ ( $f^{(k)}(x)$ は 次の近似根を示す) を計算する. ただし、 右辺の除算はべき級数の除算を用い、 計算結果をべき級数とする.
$k$(2)$f^{(2^{l})}(x) \equiv f^{(2^{1-1})}(x)-\frac{g(x,f^{(2^{\mathrm{I}-1})}(x))}{\Delta,\partial y\partial(x,f^{()}2^{l-1}(x))}(\mathrm{m}\mathrm{o}\mathrm{d} x^{2^{l}+1})$

上の記号的ニュートン法は効率的に近似根を計算するが、高次の項の係数を計算するために低次の項の係数を用いており、誤差の伝播が起こりうる。 そこで、本稿では高次の項の係数を直接計算する算法を与える。本稿の計算法は計算効率の点では記号的ニュートン法に劣るが、記号的ニュートン法で求めた近似根の精度の検証などに有用である。

3

基本的な考え方2 変数多項式 $h(x,y)$ の$y$} こ関する $m$ 次の近似根 $\phi^{(m)}(y)$ が(3)$\phi^{(m)}(y)=\alpha_{0}+\alpha_{1}y+\alpha_{2}y^{2}+\cdots+\alpha_{m}y^{m}$で与えられるとする。 Cauchy の積分公式より $\phi^{(m)}(y)$ の$k$

次の項の係数$\alpha_{k}$は(4)$\alpha_{k}=\frac{1}{2\pi}\int\frac{\phi(y)}{y^{k+1}}dy$で与えられる。上式で$y=\rho e^{i\theta}$

とおいて整理すると(5)$\alpha_{k}=\frac{1}{2\pi\rho^{k}}\int_{0}^{2\pi}\phi(\rho e^{i\theta})e^{-ik\theta}dy$を得る。実部と虚部を分けるために$\alpha_{k}=\frac{1}{2\pi\rho^{k}}\int_{0}^{2\pi}[f(\theta)$$f(\theta)={\rm Re}(\phi(\rho e^{i\theta})),$$g(\theta)={\rm Im}(\phi(\rho e^{i\theta}))$とおくと上式はcos(k の十 $g(\theta)\sin(k\theta)+i\{g(\theta)\cos(k\theta)-f(\theta)$ sin(k の月(6)$d\theta$となる。

(6) 式より、を数値積分法で計算すればよいのであるが、実際に計算してみると 力吠きくなったときには、 Gauss 積分などでは計算誤差が多く、実用的でないことがわかった。 これはcos(k のやsin(k のの関数が激しく振動し、数値積分が難しい関数になっているからである。そこで、
cos(k のやsin(k のなども積分できる数値積分法が必要になった。

4$k$$\alpha_{k}$$f(\theta)$$g(\theta)$$f(\theta)$$g(\theta)$

Filon の方法

sin(k のなども積分できる数値積分法として、 Filon の方法と呼ばれる方法に着目しcos(7c のやた。 Filon の方法では、を次のようにして計算する (も同様に計算す( $N$ は分割数) とおき、ることができるが、 ここでは省略する)。 まず、$C_{2n-1}$ を次のよう$f(\theta)$$g(\theta)$$\int_{a}^{b}f(\theta)\cos(k\theta)d\theta$$\int_{a}^{b}f(\theta)\sin(k\theta)d\theta$$h= \frac{b-a}{2N}$$C_{2n},$に求める。
$C_{2n}$$=$$\frac{1}{2}f(a)\cos(ka)+f(a+2h)\cos k(a+2h)+f(a+4h)\cos k(a+4h)+\cdots+$$\frac{1}{2}f(b)\cos(kb)$$C_{2n-1}$$=$$f(a+h)\cos k(a+h)+f(a+3h)\cos k(a+3h)+\cdots+f(b-h)\cos k(b-h)$
XXXI-2

211

表 1: 数値実験の結果

次に$\alpha,$$\beta,\gamma$を下のように定めると$=$$kh$$\alpha$$=$$(\theta^{2}+\theta\sin\theta cos\theta-2\sin 2\theta)/\theta^{3}$$\beta$$=$2 $[\theta(1+\cos^{2}\theta)-2\sin\theta\cos\theta]/\theta^{3}$$=$$4(\sin\theta-\theta\cos\theta)/\theta^{3}$$\theta$$\gamma$$\int_{a}^{b}f(\theta)\cos(k\theta)d\theta$の近似値はそれぞれ下のように与えられる。

$\int_{a}^{b}f(\theta)\cos(k\theta)d\theta\cong h$

5{$\alpha[f(b)\sin$kb-f(a) $\sin ka]+\beta C_{2n}+\gamma C_{2n-1}$ }(7)

アルゴリズム以上をまとめると、 (3) の$\alpha_{k}$$\alpha_{k}$を求めるアルゴリズムとして、 次のものを得る。の計算法(1),${\rm Res}(g,\partialx$$x)\theta\Delta=0$の根 (特異点) を求め、$\rho$を決定する ($\rho e^{:\theta}$が特異点を囲まないように、$\rho$を十分小さく取る)。

(2)$f(\theta)={\rm Re}(\phi(\rho e^{i\theta})),$ $g(\theta)={\rm Im}(\phi(\rho e^{i\theta}))$とおき、$\int_{0}^{2\pi}f(\theta)\sin k\theta d\theta,$$\int_{0}^{2\pi}g(\theta)\sin k\theta d\theta$を Filon の方法で計算する。

(3)6$\alpha_{k}$を (7) より計算する。

数値実験$h(x,y)=x^{2}+xy+y^{2}-\cdot 3y-2$ の近似根の $k$ 次の係数$\alpha_{k}$を前節のアノレゴリズムで求め、 それに含まれる相対誤差を調べた。 その結果を表 1 に示す。

7

結論

Cauchy の積分公式を用いて、 べき級数根の高次の係数を直接求める計算式を導 ‘た。 この計算式で ま、振動の激しい関数を数値積分する必要がある。 そこで、 Filon の方法を用いてこの数値積分を行 ‘、数値実験を行った。 その結果、 100 次の項の係数までは十分な精度で求められたが、 200 次の項の係数には実用的でないほどの誤差が含まれていた。今後は、 より良い数値積分法を探すことと誤差の解析を課題としたい。
$\mathrm{A}$$\mathrm{A}$XXXI-3212
参考文献[1] T.Kitamoto “Approximate Eigenvalues, Eigenvectors and Inverse of aMatrix with Polynomial Entries,” Jpn. Indus. Appl. Math., Vol. 11, No. 1, pp. 73-85, 1994.
[2] T.Kitamoto “Hensel Construction with an Arbitrary Degree of Convergence,”
Math., Vol. 13, No. 2, pp. 203-215, 1996.

$\mathrm{J}\mathrm{p}\mathrm{n}.$Indus. Appl.

[3] T.Kitamoto “記号的ニュートン法の高次収束への拡張について” 電子情報通信学会論文誌掲載予定.
[4] D.E.Knuth “準数値算法/算術演算” サイエンス社, 東京, 1985.

. Kung and1978.

[5]$\mathrm{H}.\mathrm{T}$[6]$\mathrm{J}.\mathrm{D}$$\mathrm{J}.\mathrm{F}$. Traub “All algebraic functions can be computed fast,” J. ACM, 25,$\mathrm{P}\mathrm{P}$. 245-260,. Lipson “Newton’s method: A great algebraic algorithm,” Proc. of ACM symposium on symbolicand algebraic computations, pp. 260-270, 1976.

[7] M.Shaw and J.F.Traub “On the Number of Multiplication for the Evaluation of aPolynomial and. 21, pp. 161-167, 1974.
Some of Its Derivatives,” JACM,$\mathrm{V}\mathrm{o}\mathrm{l}$XXXI-4