代数方程式

組み立て除法
m次の多項式P(x)を(x+β)で割った時、
商はm-1次の多項式Q(x)
余りをrとおける。
PとQは次の通り
P(x) = xm+ p1xm-1+ p2xm-2+ …+ pm-1x+ pm
Q(x) = xm-1+ q1xm-2+ q2xm-3+ …+ qm-2x+ qm-1
p0=q0=1と考えれる。

余りrの計算方法は
p0〜m、m、βが分かってる状態で
kを1からm-1まで増やしながら
→qk=pk-βk-1
r=pm-βqm-1
このアルゴリズムを組立除法(synthetic division)という。

#include <stdio.h>
main(){
	int p[]={1,-6,11,6};	//Pの係数
	int q[]={1,0,0,0};	//Qの係数になる予定
	int m=3,beta=-4,r;	//定数
	int i;

	for(i=1;i<m;i++)q[i]=p[i]-beta*q[i-1];
	r=p[m]-beta*q[m-1];

	//出力
	for(i=0;i<=m;i++)printf("%dx^%d %s",p[i],m-i,i!=m?"+ ":" = 0");
	printf("を(x+%d)で割った余りは\nr=%d",beta,r);
}


三次方程式
カルダノの公式(Cardano's formula)を使う。
xの三次式
x3 + a･x2 + b･x + c = 0

これにy=x+a/3を適用すると
y3 + q･y + r = 0 ; 但し q=a2+b, r=(2/27)a3-(1/3)ab+c
と、２次の係数が0になる。

さらに、y=(u+v)を代入して展開すると、
u3 + 3uv(u+v) + v3 + q(u+v) + r = 0
⇔u3+v3+r + (3uv+q)(u+v) = 0

この式が成り立つ"十分"条件として
u3+v3+r = 0　…�@
(3uv+q) = 0　…�A
が成り立てばよい。
�Aからv=-q/(3u)
これを�@に代入してu3倍すれば
u6 + r･u3 - (1/27)q3=0
これはU=(u3)の二次方程式として解ける。
U = -r/2 ± √(r2/4-a3/27);
u=3√U

uが出れば�Aより、
v = -q/(3u);
∴y = u+v
　　= u-q/(3u)
∴x = y-a/3
　　= u-q/(3u)-a/3

但し、U=u3からuを求める際に、
z=3√u3
としたら、
ωは1の3乗根の一つ (√3･i-1)/2 を使って
　zω、zω2
も、3乗すればu3になる。
よって(u,v) = (z,-q/(3z)),(zω,-qω2/(3z)),(zω2,-qω/(3z) )という3通り存在し、
xは３つ見付かる。


ニュートン法、再び
５次以上の多項式の解は代数的に計算することは出来ないから近似するしかない。
整式f(x)についてはf'(x)が容易に求められるから、ニュートン法が利用できる。
f(x) = Σ[k=0,n]akxk
f'(x) = Σ[k=1,n]k･akxk-1
x0=(適当な近似値);
次を反復する。
xk+1 = xk - f(xk)/f'(xk)