2012-10-30 08:34:12
一元三次方程的求解比一元二次方程困难, 求根公式看起来也很复杂, 中文维基百科词条 上原先的三角函数解公式有误, 我已经修改过了. 但是我觉得还有必要将这些公式简化一下, 给出具体的算法, 以便编程使用.
首先需要明确, 一元三次方程 $f(x)=ax^3+bx^2+cx+d=0,(a\neq0)$ 至少有一个实根. 这是由于 $f(-\infty)=-\infty, f(\infty)=\infty$, 而 $f(x)$ 在整个区间上连续, 所以由微积分的介值定理, 必存在一点 $x_0$ 满足 $f(x_0)=0$ . 直观上也很容易理解, 曲线 $f(x)$ 取值范围为负无穷到正无穷, 那么它与x轴至少有一个交点, 这个交点就是方程的实根. 这个结论可以推广到所有奇数次的多项式方程.
对一元三次方程 $f(x)=ax^3+bx^2+cx+d=0, (a\neq0)$, 首先化为 $f(x)=x^3+{3b \over 3a} x^2+{6c \over 6a} x+{2d \over 2a}=x^3+3b’x^2+6c’x+2d’=0$ 的形式, 为简便我们仍记作 $f(x)=x^3+3bx^2+6cx+2d=0$ , 令
\[\begin{align} \alpha &=-b^3+3bc-d \\ \beta &=b^2-2c \\ \Delta &=\alpha^2-\beta^3 \\ R_1 &=\sqrt[3]{\alpha+\sqrt{\Delta}} \\ R_2 &=\sqrt[3]{\alpha-\sqrt{\Delta}}={\beta \over R_1} \end{align}\]则方程的三个根可写作：
\[\begin{align} x_1 &=-b+R_1+R_2 \\ x_2 &=-b-{(R_1+R_2) \over 2} +{\sqrt{3}(R_1-R_2) \over 2} i \\ x_3 &=-b-{(R_1+R_2) \over 2} -{\sqrt{3}(R_1-R_2) \over 2} i \end{align}\]其中 $\Delta$ 即为三次方程根的判别式, 由其容易看出：
- 当 $\Delta>0$ 时, 方程有一个实根 $x_1$ 和两个共轭复根 $x_2$, $x_3$
- 当 $\Delta=0$ 时, 有 $R_1=R_2=\sqrt[3]{\alpha}=R$ ,
- 若 $R=0$ , 方程只有一个三重复根 $x_1=x_2=x_3=-b$
- 若 $R\neq0$ , 方程有一个实根 $x_1=-b+2R$ 和一个二重复根 $x_2=x_3=-b-R$
- 当 $\Delta \lt 0$ 时, 方程有三个不等实根, 令 $\theta=\arccos( {\alpha \over \beta^{3/2}} )$ , 则
$x_1=-b+2\sqrt{\beta} \cos( {\theta \over 3} )
x_2=-b+2\sqrt{\beta} \cos( { {\theta+2\pi} \over 3} )
x_3=-b+2\sqrt{\beta} \cos( { {\theta-2\pi} \over 3} )$
很显然, $\Delta \lt 0$ 的情况最复杂, 也很让人迷惑, 因为实根必须借助于复数才能求得, 这也是最初要引入复数的原因. 由于我没有学过复变函数, 下面的说法可能有误, 说出来只是为了帮助大家理解.
当 $\Delta \lt 0$ 时, 记 $z=\alpha+\sqrt{|\Delta|}i, \bar z=\alpha-\sqrt{|\Delta|}i, R_1=\sqrt[3]{z}, R_2=\sqrt[3]{\bar z}$, 若我们认为 $\sqrt[3]{\bar z}=\overline {\sqrt[3]{z} }$ , 则 $R_2=\bar {R_1}$, 记 $R_1=R, R_2=\bar R, \omega= { {-1+\sqrt{3} i} \over 2}, \bar \omega={ {-1-\sqrt{3} i} \over 2}$ , 则三个根可写作:
\[\begin{align} x_1 &=-b+R_1+R_2=-b+R+\bar R \\ x_2 &=-b-{(R_1+R_2) \over 2} +{\sqrt{3}(R_1-R_2) \over 2} i =-b+\omega R+\bar \omega \bar R \\ x_3 &=-b-{(R_1+R_2) \over 2} -{\sqrt{3}(R_1-R_2) \over 2} i =-b+\bar \omega R+ \omega \bar R \end{align}\]很显然
\[\begin{align} x_1 &=-b+2\Re(R) \\ \Re(R) &=\Re \sqrt[3]{z} = \sqrt[3]{|z|} \cos({ {\theta+2k\pi} \over 3}), k=0, 1, -1 \\ |z| &=\sqrt{\alpha^2+|\Delta|}=\sqrt{\alpha^2-\alpha^2+\beta^3} = \beta^{3/2} \\ \theta &= \arccos({\alpha \over |z|}) = \arccos({\alpha \over \beta^{3/2}}) \end{align}\]由于复变函数的多值性, 我们一下就得到了所有的三个根, 正如上面所给出的. 另外由于 $\omega$ 对应的正是旋转 $2\pi \over 3$ 的复数, 所以 $x_2, x_3$ 并不能给出新的根.
代码
# Language: fortran
Subroutine getCubicRoot(P, X)
real*8, parameter:: TwoPi = 8.D0*atan(1.D0)
real*8 P(*), a, b, c, d, Alph, Beta, Delt, R1, R2, tht, X(3)
X = 0.D0
a = P(1)
b = P(2)/(3.D0*a)
c = P(3)/(6.D0*a)
d = P(4)/(2.D0*a)
Alph = -b*b*b + 3.D0*b*c - d
Beta = b*b - 2.D0*c
Delt = Alph*Alph-Beta*Beta*Beta
if(Delt>0.D0) then
tht = Alph+sqrt(Delt); R1 = sign(abs(tht)**(1.D0/3.D0), tht)
tht = Alph-sqrt(Delt); R2 = sign(abs(tht)**(1.D0/3.D0), tht)
X(1) = -b+R1+R2
else if(Delt==0.D0) then
R1 = sign(abs(Alph)**(1.D0/3.D0), Alph)
if(R1==0.D0) then
X(1) = -b
else
X(1) = -b+2.D0*R1
X(2) = -b-R1
end if
else if(Delt<0.D0) then
tht = acos(Alph/(sqrt(Beta)*Beta))
X(1) = -b+2.D0*sqrt(Beta)*cos(tht/3.D0)
X(2) = -b+2.D0*sqrt(Beta)*cos((tht+TwoPi)/3.D0)
X(3) = -b+2.D0*sqrt(Beta)*cos((tht-TwoPi)/3.D0)
end if
End Subroutine getCubicRoot
