|
|
> No.1831[元記事へ]
エス・テーさんへのお返事です。
> Ferrariの公式で4次方程式を解くプログラム
!フェラーリ（Ferrari）の方法
OPTION ARITHMETIC COMPLEX !複素数を扱う
!a4*x^4+a3*x^3+a2*x^2+a1*x+a0=0、a4≠0の解
LET a4=1 !1±√-1,-1±√-3
LET a3=0
LET a2=2
LET a1=-4
LET a0=8
PRINT a4;a3;a2;a1;a0 !debug
!x=y-a3/(4*a4)を代入して、3次の項を消去すると、
!y^4+py^2+qy+r=0
!ただし、
! p=-3*a3^2/(8*a4^2)+a2/a4
! q=a3^3/(8*a4^3)-a2*a3/(2*a4^2)+a1/a4
! r=-3*a3^4/(256*a4^4)+a2*a3^2/(16*a4^3)-a1*a3/(4*a4^2)+a0/a4
LET p=-3*a3^2/(8*a4^2)+a2/a4
LET q=a3^3/(8*a4^3)-a2*a3/(2*a4^2)+a1/a4
LET r=-3*a3^4/(256*a4^4)+a2*a3^2/(16*a4^3)-a1*a3/(4*a4^2)+a0/a4
!!!PRINT p;q;r !debug
!q=0のとき、
! 複2次式 (y^2 - (-p-√{p^2-4r})/2)(y^2 - (-p+√{p^2-4r})/2)=0 と因数分解される。
!q≠0のとき、
! y^4=-p*y^2-q*y-rとして、両辺にy^2*z+z^2/4を加えて、変形すると、
! (y^2+z/2)^2=(z-p){y-q/(2*(z-p))}^2 + {1/(4*(z-p))}(z^3-p*z^2-4*r*z+4*p*r-q^2)
! ここで、z^3-p*z^2-4*r*z+4*p*r-q^2=0となるzを1つ求めると、
! 同様に、複2次式
! (y^2+z/2 - √(z-p){y-q/(2*(z-p))})(y^2+z/2 + √(z-p){y-q/(2*(z-p))})=0
! と因数分解できる。
IF q=0 THEN
LET t=p^2-4*r !判別式
CALL Solve2Equ(0,-(-p-SQR(t))/2 ,y1,y2)
CALL Solve2Equ(0,-(-p+SQR(t))/2 ,y3,y4)
ELSE
CALL Solve3Equ(-p,-4*r,4*p*r-q^2, z1,z2,z3)
LET t=SQR(z1-p)
CALL Solve2Equ(-t, q/(2*t)+z1/2, y1,y2)
CALL Solve2Equ( t,-q/(2*t)+z1/2, y3,y4)
END IF
LET x1=y1-a3/(4*a4) !x=y-a3/(4*a4)
LET x2=y2-a3/(4*a4)
LET x3=y3-a3/(4*a4)
LET x4=y4-a3/(4*a4)
!解を表示する
PRINT x1
PRINT x2
PRINT x3
PRINT x4
END
EXTERNAL SUB Solve2Equ(P,Q, x1,x2) !代数方程式 x^2+P*x+Q=0 の解
OPTION ARITHMETIC COMPLEX !複素数を扱う
LET D=P^2-4*Q !判別式
LET x1=(-P+SQR(D))/2 !解の公式より
LET x2=(-P-SQR(D))/2
END SUB
EXTERNAL SUB Solve3Equ(P,Q,R, x1,x2,x3)!代数方程式 x^3+P*x^2+Q*x+R=0 の解
OPTION ARITHMETIC COMPLEX !複素数を扱う
LET a=-P^2/3+Q !カルダノ（Cardano）の方法より
LET b=2*P^3/27-P*Q/3+R
LET t=b^2/4+a^3/27 !(b/2)^2+(a/3)^3
LET z1=-b/2+SQR(t)
LET z2=-b/2-SQR(t)
IF t>=0 THEN
!!!PRINT "1実根と2虚根" !debug
LET u=SGN(z1)*ABS(z1)^(1/3) !u=(-b/2+SQR(t))^(1/3)∈R
LET v=SGN(z2)*ABS(z2)^(1/3) !v=(-b/2+SQR(t))^(1/3)∈R
ELSE !不還元の場合
!!!PRINT "3実根" !debug
LET u=EXP(LOG(z1)/3) !u=(-b/2+SQR(t))^(1/3)∈C
LET v=EXP(LOG(z2)/3) !v=(-b/2-SQR(t))^(1/3)∈C
END IF
!!!PRINT u*v; -a/3 !debug
LET w=COMPLEX(-1/2,SQR(3)/2) !ω=(-1+(√3)i)/2は1の原始3乗根の1つ
LET y1=u+v !y^3+a*y+b=0
LET y2=w*u+w^2*v
LET y3=w^2*u+w*v
LET x1=y1-P/3
LET x2=y2-P/3
LET x3=y3-P/3
END SUB
|
|