Метод Ньютона с заданными цифрами точности - PullRequest
Новая
       21

Метод Ньютона с заданными цифрами точности

3 голосов
/ 12 февраля 2011

Я пытаюсь написать функцию на Java, которая вычисляет n-й корень числа.Я использую метод Ньютона для этого.Тем не менее, пользователь должен иметь возможность указать, сколько цифр точности они хотят.Это та часть, с которой у меня возникают проблемы, поскольку мой ответ часто не совсем верен.Соответствующий код здесь: http://pastebin.com/d3rdpLW8. Как я могу исправить этот код, чтобы он всегда давал ответ с точностью не менее p цифр?(не делая больше работы, чем необходимо) 

import java.util.Random;

public final class Compute {

    private Compute() {
    }

    public static void main(String[] args) {
        Random rand = new Random(1230);
        for (int i = 0; i < 500000; i++) {
            double k = rand.nextDouble()/100;
            int n = (int)(rand.nextDouble() * 20) + 1;
            int p = (int)(rand.nextDouble() * 10) + 1;
            double math = n == 0 ? 1d : Math.pow(k, 1d / n);
            double compute = Compute.root(n, k, p);
            if(!String.format("%."+p+"f", math).equals(String.format("%."+p+"f", compute))) {
                System.out.println(String.format("%."+p+"f", math));
                System.out.println(String.format("%."+p+"f", compute));
                System.out.println(math + " " + compute + " " + p);
            }
        }
    }

    /**
     * Returns the n-th root of a positive double k, accurate to p decimal
     * digits.
     * 
     * @param n
     *            the degree of the root.
     * @param k
     *            the number to be rooted.
     * @param p
     *            the decimal digit precision.
     * @return the n-th root of k
     */
    public static double root(int n, double k, int p) {     
        double epsilon = pow(0.1, p+2);
        double approx = estimate_root(n, k);
        double approx_prev;

        do {
            approx_prev = approx;
            // f(x) / f'(x) = (x^n - k) / (n * x^(n-1)) = (x - k/x^(n-1)) / n
            approx -= (approx - k / pow(approx, n-1)) / n;
        } while (abs(approx - approx_prev) > epsilon);
        return approx;
    }

    private static double pow(double x, int y) {
        if (y == 0)
            return 1d;
        if (y == 1)
            return x;
        double k = pow(x * x, y >> 1);
        return (y & 1) == 0 ? k : k * x;
    }

    private static double abs(double x) {
        return Double.longBitsToDouble((Double.doubleToLongBits(x) << 1) >>> 1);
    }

    private static double estimate_root(int n, double k) {
        // Extract the exponent from k.
        long exp = (Double.doubleToLongBits(k) & 0x7ff0000000000000L);
        // Format the exponent properly.
        int D = (int) ((exp >> 52) - 1023);
        // Calculate and return 2^(D/n).
        return Double.longBitsToDouble((D / n + 1023L) << 52);
    }
}

Ответы [ 3 ]

3 голосов
/ 12 февраля 2011

Просто выполняйте итерацию до тех пор, пока обновление не станет меньше, чем, скажем, 0,0001, если вам нужна точность в 4 десятичных знака.

То есть установите для эпсилона значение Math.pow(10, -n), если вы хотите n цифры точности.

2 голосов
/ 12 февраля 2011

Давайте вспомним, что говорит анализ ошибок метода Ньютона. По сути, это дает нам ошибку для n-й итерации как функцию от ошибки n-1-й итерации.

Итак, как мы можем определить, что ошибка меньше k? Мы не можем, если не знаем ошибку при е (0). И если бы мы знали ошибку при e (0), мы бы просто использовали ее, чтобы найти правильный ответ.

То, что вы можете сделать, это сказать «e (0) <= m». Затем вы можете найти n такое, что e (n) <= k для желаемого k. Однако для этого необходимо знать максимальное значение f '' в вашем радиусе, что является (в общем) такой же сложной проблемой, как и поиск перехвата x. </p>

Что вы проверяете, так это если ошибка изменяется менее чем на k, что является вполне приемлемым способом сделать это. Но это не проверка, если ошибка меньше k. Как отметил Аксель и другие, существует множество других алгоритмов аппроксимации корней, некоторые из которых приведут к более простому анализу ошибок, и если вы действительно этого хотите, вам следует использовать один из них.

0 голосов
/ 16 февраля 2016

У вас есть ошибка в вашем коде.Последняя строка вашего pow() метода должна читатьсяreturn (y & 1) == 1 ? k : k * x;а не
return (y & 1) == 0 ? k : k * x;

...