java-使用法拉利方法查找四次方程的实根
作者:互联网
我目前正在尝试使用Ferrari’s method from Wikipedia求解四次方程.我只想获取真实的根,而舍弃虚的根.我的实现没有为真正的根返回良好的价值.我找不到公式中的错误.
我的三次方程式按预期工作,我的二次方程式也可以.现在,我只是错过了要完成的法拉利方法,但我无法使它起作用!
这是我的课:
public class QuarticFunction {
private final double a;
private final double b;
private final double c;
private final double d;
private final double e;
public QuarticFunction(double a, double b, double c, double d, double e) {
this.a = a;
this.b = b;
this.c = c;
this.d = d;
this.e = e;
}
public final double[] findRealRoots() {
if (a == 0) {
return new CubicEquation(b, c, d, e).findRealRoots();
}
if (isBiquadratic()) { //This part works as expected
return solveUsingBiquadraticMethod();
}
return solveUsingFerrariMethodWikipedia();
}
private double[] solveUsingFerrariMethodWikipedia() {
// http://en.wikipedia.org/wiki/Quartic_function#Ferrari.27s_solution
// ERROR: Wrong numbers when Two Real Roots + Two Complex Roots
// ERROR: Wrong numbers when Four Real Roots
QuarticFunction depressedQuartic = toDepressed();
if (depressedQuartic.isBiquadratic()) {
return depressedQuartic.solveUsingBiquadraticMethod();
}
double y = findFerraryY(depressedQuartic);
double originalRootConversionPart = -b / (4 * a);
double firstPart = Math.sqrt(depressedQuartic.c + 2 * y);
double positiveSecondPart = Math.sqrt(-(3 * depressedQuartic.c + 2 * y + 2 * depressedQuartic.d / Math.sqrt(depressedQuartic.c + 2 * y)));
double negativeSecondPart = Math.sqrt(-(3 * depressedQuartic.c + 2 * y - 2 * depressedQuartic.d / Math.sqrt(depressedQuartic.c + 2 * y)));
double x1 = originalRootConversionPart + (firstPart + positiveSecondPart) / 2;
double x2 = originalRootConversionPart + (-firstPart + negativeSecondPart) / 2;
double x3 = originalRootConversionPart + (firstPart - positiveSecondPart) / 2;
double x4 = originalRootConversionPart + (-firstPart - negativeSecondPart) / 2;
Set<Double> realRoots = findAllRealRoots(x1, x2, x3, x4);
return toDoubleArray(realRoots);
}
private double findFerraryY(QuarticFunction depressedQuartic) {
double a3 = 1;
double a2 = 5 / 2 * depressedQuartic.c;
double a1 = 2 * Math.pow(depressedQuartic.c, 2) - depressedQuartic.e;
double a0 = Math.pow(depressedQuartic.c, 3) / 2 - depressedQuartic.c * depressedQuartic.e / 2
- Math.pow(depressedQuartic.d, 2) / 8;
//CubicEquation works as expected! No need to worry! It returns either 1 or 3 roots.
CubicEquation cubicEquation = new CubicEquation(a3, a2, a1, a0);
double[] roots = cubicEquation.findRealRoots();
for (double y : roots) {
if (depressedQuartic.c + 2 * y != 0) {
return y;
}
}
throw new IllegalStateException("Ferrari method should have at least one y");
}
public final boolean isBiquadratic() {
return Double.compare(b, 0) == 0 && Double.compare(d, 0) == 0;
}
private double[] solveUsingBiquadraticMethod() {
//It works as expected!
QuadraticLine quadraticEquation = new QuadraticLine(a, c, e);
if (!quadraticEquation.hasRoots()) {
return new double[] {};
}
double[] quadraticRoots = quadraticEquation.findRoots();
Set<Double> roots = new HashSet<>();
for (double quadraticRoot : quadraticRoots) {
if (quadraticRoot > 0) {
roots.add(Math.sqrt(quadraticRoot));
roots.add(-Math.sqrt(quadraticRoot));
} else if (quadraticRoot == 0.00) {
roots.add(0.00);
}
}
return toDoubleArray(roots);
}
public QuarticFunction toDepressed() {
// http://en.wikipedia.org/wiki/Quartic_function#Converting_to_a_depressed_quartic
double p = (8 * a * c - 3 * Math.pow(b, 2)) / (8 * Math.pow(a, 2));
double q = (Math.pow(b, 3) - 4 * a * b * c + 8 * d * Math.pow(a, 2)) / (8 * Math.pow(a, 3));
double r = (-3 * Math.pow(b, 4) + 256 * e * Math.pow(a, 3) - 64 * d * b * Math.pow(a, 2) + 16 * c * a
* Math.pow(b, 2)) / (256 * Math.pow(a, 4));
return new QuarticFunction(1, 0, p, q, r);
}
private Set<Double> findAllRealRoots(double... roots) {
Set<Double> realRoots = new HashSet<>();
for (double root : roots) {
if (!Double.isNaN(root)) {
realRoots.add(root);
}
}
return realRoots;
}
private double[] toDoubleArray(Collection<Double> values) {
double[] doubleArray = new double[values.size()];
int i = 0;
for (double value : values) {
doubleArray[i] = value;
++i;
}
return doubleArray;
}
}
我尝试了一个更简单的实现found here,但是遇到了一个新问题.现在,一半的根源是好的,但是另一半是错误的.这是我的尝试:
private double[] solveUsingFerrariMethodTheorem() {
// https://proofwiki.org/wiki/Ferrari's_Method
CubicEquation cubicEquation = getCubicEquationForFerrariMethodTheorem();
double[] cubicRoots = cubicEquation.findRealRoots();
double y1 = findFirstNonZero(cubicRoots);
double inRootOfP = Math.pow(b, 2) / Math.pow(a, 2) - 4 * c / a + 4 * y1;
double p1 = b / a + Math.sqrt(inRootOfP);
double p2 = b / a - Math.sqrt(inRootOfP);
double inRootOfQ = Math.pow(y1, 2) - 4 * e / a;
double q1 = y1 - Math.sqrt(inRootOfQ);
double q2 = y1 + Math.sqrt(inRootOfQ);
double x1 = (-p1 + Math.sqrt(Math.pow(p1, 2) - 8 * q1)) / 4;
double x2 = (-p2 - Math.sqrt(Math.pow(p2, 2) - 8 * q1)) / 4;
double x3 = (-p1 + Math.sqrt(Math.pow(p1, 2) - 8 * q2)) / 4;
double x4 = (-p2 - Math.sqrt(Math.pow(p2, 2) - 8 * q2)) / 4;
Set<Double> realRoots = findAllRealRoots(x1, x2, x3, x4);
return toDoubleArray(realRoots);
}
private CubicEquation getCubicEquationForFerrariMethodTheorem() {
double cubicA = 1;
double cubicB = -c / a;
double cubicC = b * d / Math.pow(a, 2) - 4 * e / a;
double cubicD = 4 * c * e / Math.pow(a, 2) - Math.pow(b, 2) * e / Math.pow(a, 3) - Math.pow(d, 2)
/ Math.pow(a, 2);
return new CubicEquation(cubicA, cubicB, cubicC, cubicD);
}
private double findFirstNonZero(double[] values) {
for (double value : values) {
if (Double.compare(value, 0) != 0) {
return value;
}
}
throw new IllegalArgumentException(values + " does not contain any non-zero value");
}
我不知道我在想什么.我花了数小时进行调试以尝试查看错误,但是现在我完全迷失了(加上一些头痛).我究竟做错了什么?我不在乎要使用哪个公式,因为我只希望其中一个起作用.
解决方法:
我有两个错误.
>我所有的“整数常量”都应该加倍.
>在Wikipedia方法中,当沮丧的四次方是双二次方时,我们必须将根重新转换为原始四次方.我正在归还沮丧的根源.
这是我的最终实现
public class QuarticFunction {
private static final double NEAR_ZERO = 0.0000001;
private final double a;
private final double b;
private final double c;
private final double d;
private final double e;
public QuarticFunction(double a, double b, double c, double d, double e) {
this.a = a;
this.b = b;
this.c = c;
this.d = d;
this.e = e;
}
public final double[] findRealRoots() {
if (Math.abs(a) < NEAR_ZERO) {
return new CubicFunction(b, c, d, e).findRealRoots();
}
if (isBiquadratic()) {
return solveUsingBiquadraticMethod();
}
return solveUsingFerrariMethodWikipedia();
}
private double[] solveUsingFerrariMethodWikipedia() {
// http://en.wikipedia.org/wiki/Quartic_function#Ferrari.27s_solution
QuarticFunction depressedQuartic = toDepressed();
if (depressedQuartic.isBiquadratic()) {
double[] depressedRoots = depressedQuartic.solveUsingBiquadraticMethod();
return reconvertToOriginalRoots(depressedRoots);
}
double y = findFerraryY(depressedQuartic);
double originalRootConversionPart = -b / (4.0 * a);
double firstPart = Math.sqrt(depressedQuartic.c + 2.0 * y);
double positiveSecondPart = Math.sqrt(-(3.0 * depressedQuartic.c + 2.0 * y + 2.0 * depressedQuartic.d
/ Math.sqrt(depressedQuartic.c + 2.0 * y)));
double negativeSecondPart = Math.sqrt(-(3.0 * depressedQuartic.c + 2.0 * y - 2.0 * depressedQuartic.d
/ Math.sqrt(depressedQuartic.c + 2.0 * y)));
double x1 = originalRootConversionPart + (firstPart + positiveSecondPart) / 2.0;
double x2 = originalRootConversionPart + (-firstPart + negativeSecondPart) / 2.0;
double x3 = originalRootConversionPart + (firstPart - positiveSecondPart) / 2.0;
double x4 = originalRootConversionPart + (-firstPart - negativeSecondPart) / 2.0;
Set<Double> realRoots = findOnlyRealRoots(x1, x2, x3, x4);
return toDoubleArray(realRoots);
}
private double[] reconvertToOriginalRoots(double[] depressedRoots) {
double[] originalRoots = new double[depressedRoots.length];
for (int i = 0; i < depressedRoots.length; ++i) {
originalRoots[i] = depressedRoots[i] - b / (4.0 * a);
}
return originalRoots;
}
private double findFerraryY(QuarticFunction depressedQuartic) {
double a3 = 1.0;
double a2 = 5.0 / 2.0 * depressedQuartic.c;
double a1 = 2.0 * Math.pow(depressedQuartic.c, 2.0) - depressedQuartic.e;
double a0 = Math.pow(depressedQuartic.c, 3.0) / 2.0 - depressedQuartic.c * depressedQuartic.e / 2.0
- Math.pow(depressedQuartic.d, 2.0) / 8.0;
CubicFunction cubicFunction = new CubicFunction(a3, a2, a1, a0);
double[] roots = cubicFunction.findRealRoots();
for (double y : roots) {
if (depressedQuartic.c + 2.0 * y != 0.0) {
return y;
}
}
throw new IllegalStateException("Ferrari method should have at least one y");
}
private double[] solveUsingBiquadraticMethod() {
QuadraticFunction quadraticFunction = new QuadraticFunction(a, c, e);
if (!quadraticFunction.hasRoots()) {
return new double[] {};
}
double[] quadraticRoots = quadraticFunction.findRoots();
Set<Double> roots = new HashSet<>();
for (double quadraticRoot : quadraticRoots) {
if (quadraticRoot > 0.0) {
roots.add(Math.sqrt(quadraticRoot));
roots.add(-Math.sqrt(quadraticRoot));
} else if (quadraticRoot == 0.00) {
roots.add(0.00);
}
}
return toDoubleArray(roots);
}
public QuarticFunction toDepressed() {
// http://en.wikipedia.org/wiki/Quartic_function#Converting_to_a_depressed_quartic
double p = (8.0 * a * c - 3.0 * Math.pow(b, 2.0)) / (8.0 * Math.pow(a, 2.0));
double q = (Math.pow(b, 3.0) - 4.0 * a * b * c + 8.0 * d * Math.pow(a, 2.0)) / (8.0 * Math.pow(a, 3.0));
double r = (-3.0 * Math.pow(b, 4.0) + 256.0 * e * Math.pow(a, 3.0) - 64.0 * d * b * Math.pow(a, 2.0) + 16.0 * c
* a * Math.pow(b, 2.0))
/ (256.0 * Math.pow(a, 4.0));
return new QuarticFunction(1.0, 0.0, p, q, r);
}
private Set<Double> findOnlyRealRoots(double... roots) {
Set<Double> realRoots = new HashSet<>();
for (double root : roots) {
if (Double.isFinite(root)) {
realRoots.add(root);
}
}
return realRoots;
}
private double[] toDoubleArray(Collection<Double> values) {
double[] doubleArray = new double[values.size()];
int i = 0;
for (double value : values) {
doubleArray[i] = value;
++i;
}
return doubleArray;
}
}
标签:equation-solving,trigonometry,polynomial-math,java,math 来源: https://codeday.me/bug/20191028/1955957.html