盛金公式解一元三次方程_用盛金公式求解一元三次方程

解一元三次方程一般用盛金公式求解，算法高效且求出来的解精确。

百度百科关于盛金公式有如下解释：

盛金公式

Shengjin's Formulas

一元三次方程aX^3+bX^2+cX+d=0，(a，b，c，d∈R，且a≠0)。

重根判别式：A=b^2－3ac；B=bc－9ad；C=c^2－3bd，

总判别式：Δ=B^2－4AC。

当A=B=0时，盛金公式①：

X⑴=X⑵=X⑶=－b/(3a)=－c/b=－3d/c。

当Δ=B^2－4AC>0时，盛金公式②：

X⑴=(－b－Y⑴^(1/3)－Y⑵^(1/3))/(3a)；

X(2,3)=(－2b+Y⑴^(1/3)+Y⑵^(1/3))/(6a)±i3^(1/2)(Y⑴^(1/3)－Y⑵^(1/3))/(6a)，

其中Y(1, 2)=Ab+3a(－B±(B^2－4AC)^(1/2))/2，i^2=－1。

当Δ=B^2－4AC=0时，盛金公式③：

X⑴=－b/a+K；X⑵=X⑶=－K/2，

其中K=B/A，(A≠0)。

当Δ=B^2－4AC<0时，盛金公式④：

X⑴=(－b－2A^(1/2)cos(θ/3))/(3a)；

X(2,3)=(－b+A^(1/2)(cos(θ/3)±3^(1/2)sin(θ/3)))/(3a)，

其中θ=arccosT，T=(2Ab－3aB)/(2A^(3/2))，(A>0，－1

盛金判别法

Shengjin's Distinguishing Means

①：当A=B=0时，方程有一个三重实根；

②：当Δ=B^2－4AC>0时，方程有一个实根和一对共轭复根；

③：当Δ=B^2－4AC=0时，方程有三个实根，其中有一个两重根；

④：当Δ=B^2－4AC<0时，方程有三个不相等的实根。

个人实现的代码实现如下：

#include

using namespace std;

#define ZERO                                             1e-6

#define ISZERO(Value)                                    ((-1.0f * ZERO < (Value)) && ((Value)  < 1.0f * ZERO))

#define LESSTHANZERO(Value)                              ((Value) < ZERO)

#define GREATERTHANZERO(Value)                           (ZERO < (Value))

#define ISBETWEEN(Value, Value1, Value2)

(GREATERTHAN(Value, MIN(Value1, Value2)) && LESSTHAN(Value,

MAX(Value1, Value2)))

typedef double RealNum;

typedef struct tagComplex

{

RealNum  real;

RealNum  image;

}Complex;

//rA * (X ^3) + rB * (X ^ 2) + rC * (X) + rD = 0

int Shengjin(RealNum rA, RealNum rB, RealNum rC, RealNum rD, RealNum* rX1, RealNum* rX2, RealNum* rX3)

{

RealNum fDA = rB * rB - 3.0f * rA * rC;

RealNum fDB = rB * rC - 9.0f * rA * rD;

RealNum fDC = rC * rC - 3.0f * rB * rD;

RealNum rTempA = ABSZERO;

RealNum rTempB = ABSZERO;

RealNum fDelta = fDB * fDB - 4.0f * fDA * fDC;

if (ISZERO(fDA) && ISZERO(fDB))

{

if (!ISZERO(rA))

{

*rX1 = -1.0f * rB / (3.0f * rA);

}

else if (!ISZERO(rB))

{

*rX1 = -1.0f * rC / rB;

}

else if (!ISZERO(rC))

{

*rX1 = -1.0f * rD / rC;

}

else// if (!ISZERO(rD))

{

return 0;

}

*rX2 = *rX1;

*rX3 = *rX1;

return 1;

}

if (GREATERTHANZERO(fDelta))

{

RealNum rY1 = fDA * rB + 3.0f * rA * (0.5f * (-1.0f * fDB + sqrtf(fDelta) ));

RealNum rY2 = fDA * rB + 3.0f * rA * (0.5f * (-1.0f * fDB - sqrtf(fDelta) ));

if (GREATERTHANZERO(rY1))

{

rTempA = powf(rY1, 1.0f / 3.0f);

}

else

{

rTempA = -1.0f * powf(-1.0f * rY1, 1.0f / 3.0f);

}

if (GREATERTHANZERO(rY2))

{

rTempB = powf(rY2, 1.0f / 3.0f);

}

else

{

rTempB = -1.0f * powf(-1.0f * rY2, 1.0f / 3.0f);

}

*rX1 = -1.0f * (rB + rTempA + rTempB) * (1.0f / (3.0f * rA));

Complex cX2;

cX2.real  = (1.0f / (6.0f * rA)) * (-2.0f * rB + rTempA + rTempB);

cX2.image = (1.0f / (6.0f * rA)) * sqrtf(3.0f) * (rTempA -rTempB);

Complex cX3;

cX3.real  = cX2.real;

cX3.image = -1.0f * cX2.image;

return 1;

}

else if (ISZERO(fDelta))

{

if (ISZERO(fDA))

{

return 0;

}

RealNum rK = fDB / fDA;

*rX1 = -1.0f * rB / rA + rK;

*rX2 = -0.5f * rK;

*rX3 = *rX2;

return 2;

}

else

{

if (LESSTHANZERO(fDA))

{

return 0;

}

if (GREATERTHANZERO(fDA))

{

rTempA = powf(fDA, 3.0f);

}

else

{

rTempA = -1.0f * powf(-1.0f * fDA, 3.0f);

}

RealNum rT = (2.0f * fDA * rB - 3.0f * rA * fDB) / (2.0f * sqrtf(rTempA));

if (!ISBETWEEN(rT, -1.0f, 1.0f))

{

return 0;

}

RealNum rTheta = acosf(rT);

*rX1 = -1.0f * (rB + 2.0f * sqrtf(fDA) * cos(rTheta / 3.0f)) * (1.0f / (3.0f * rA));

*rX2 = (-1.0f * rB + sqrtf(fDA) * (cos(rTheta / 3.0f) + sqrtf(3.0f) * sin(rTheta / 3.0f))) * (1.0f / (3.0f * rA));

*rX3 = (-1.0f * rB + sqrtf(fDA) * (cos(rTheta / 3.0f) - sqrtf(3.0f) * sin(rTheta / 3.0f))) * (1.0f / (3.0f * rA));

return 3;

}

return 0;

}