编程语言
首页 > 编程语言> > java-使用法拉利方法查找四次方程的实根

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