您可以使用
Descartes’ rule of signs和
Sturm’s theorem给出根数的上限.可以想象这可以告诉你预期有多少根.笛卡尔的规则很快,因为它只涉及
比较符号变化.
但是,您可以在当前代码中添加一些分支,以避免分支,从而获得复杂的结果.
public static double[] solveRealQuarticRoots(double a, double b, double c, double d, double e) {
double s1 = 2 * c * c * c - 9 * b * c * d + 27 * (a * d * d + b * b * e) - 72 * a * c * e,
q1 = c * c - 3 * b * d + 12 * a * e;
double discrim1 = -4 * q1 * q1 * q1 + s1 * s1;
if(discrim1 >0) {
double s2 = s1 + Math.sqrt(discrim1);
q2 = Math.cbrt(s2 / 2),
s3 = q1 / (3 * a * q2) + q2 / (3 * a),
discrim2 = (b * b) / (4 * a * a) - (2 * c) / (3 * a) + s3;
if(discrim2>0) {
double s4 = Math.sqrt(discrim2);
double s5 = (b * b) / (2 * a * a) - (4 * c) / (3 * a) - s3;
double s6 = (-(b * b * b) / (a * a * a) + (4 * b * c) / (a * a) - (8 * d) / a) / (4 * s4);
double discrim3 = (s5 - s6),
discrim4 = (s5 + s6);
// actual root values, may not be set
double r1, r2, r3, r4;
if(discrim3 > 0) {
double sqrt1 = Math.sqrt(s5-s6);
r1 = -b / (4 * a) - s4/2 + sqrt1 / 2;
r2 = -b / (4 * a) - s4/2 - sqrt1 / 2;
} else if(discrib3 == 0) {
// repeated root case
r1 = -b / (4 * a) - s4/2;
}
if(discrim4 > 0) {
double sqrt2 = Math.sqrt(s5+s6);
r3 = -b / (4 * a) + s4/2 + sqrt2 / 2;
r4 = -b / (4 * a) + s4/2 - sqrt2 / 2;
} else if(discrim4 ==0) {
r3 = -b / (4 * a) + s4/2;
}
if(discrim3 > 0 && discrim4 > 0)
return {r1,r2,r3,r4};
else if( discrim3 > 0 && discrim4 == 0 )
return {r1,r2,r3};
else if( discrim3 > 0 && discrim4 < 0 )
return {r1,r2};
else if( discrim3 == 0 && discrim4 > 0 )
return {r1,r3,r4};
else if( discrim3 == 0 && discrim4 == 0 )
return {r1,r3};
else if( discrim3 == 0 && discrim4 < 0 )
return {r1};
else if( discrim3 < 0 && discrim4 > 0 )
return {r3,r4};
else if( discrim3 < 0 && discrim4 == 0 )
return {r3};
else if( discrim3 < 0 && discrim4 < 0 )
return new double[0];
}
}
return new double[0];
}
再看下数学交流帖的答案,Don Herbison-Evans,Michel Daoud Yacoub等人提到了“寻找四方的真实根源”. Gustavo Fraidenraich.你可以找到here.在那篇论文中,他考虑了数值问题.
不要打折数值方法.牛顿非常快,可以快速收敛.
扫一扫
9983
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
