Вычисление дробных показателей в C

Вроде твоя задача - посчитать

⎛ xN ⎞(aN / aD)
⎜⎼⎼⎼⎼⎟           where xN,xD,aN,aD ∈ ℤ,  xD,aD ≠ 0
⎝ xD ⎠

с использованием только умножения, деления, сложения и вычитания с методом Ньютона в качестве предлагаемого метода для воплощать в жизнь.

Уравнение, которое мы пытаемся решить (для y):

             (aN / aD)
y = (xN / xD)            where y ∈ ℝ

Метод Ньютона находит корень функции. Если мы хотим использовать его для решения вышеуказанной проблемы, мы вычитаем правую часть из левой, чтобы получить функцию, ноль которой дает нам желаемое y:

                  (aN/aD)
f(y) = y - (xN/xD)        = 0

Не очень помогло. Я думаю, это все, что вам нужно? Дело здесь в том, чтобы пока не формировать эту функцию, потому что у нас нет способа вычислить рациональную степень рационального числа!

Во-первых, давайте решим, что оба значения aD и xD положительны. Мы можем сделать это, просто отрицая и aN, и aD, если aD было отрицательным (так что знак aN / < em> aD не изменяется), и отрицание значений xN и xD, если xD было отрицательным. Помните, что по определению ни xD, ни aD не равны нулю. Затем мы можем просто возвести обе стороны в ад степени:

 aD            aN     aN     aN
y   = (xN / xD)   = xN   / xD

Мы даже можем исключить деление, умножив обе части на последний член:

 aD     aN     aN
y   × xD   = xN

Теперь это выглядит многообещающе! Функция, которую мы получаем из этого,

        aD   aN     aN
f(y) = y   xD   - xN

Метод Ньютона также требует производной, которая, очевидно,

f(y)            aD   aN
⎼⎼⎼⎼ = df(y) = y   xD   y / aD
 dy

Сам метод Ньютона основан на повторении

             f(y)
y    = y  - ⎼⎼⎼⎼⎼⎼
 i+1    i    df(y)

Если вы разберетесь с математикой, вы обнаружите, что итерация просто

                                 aD
                y[i]      y[i] xN
y[i+1] = y[i] - ⎼⎼⎼⎼ + ⎼⎼⎼⎼⎼⎼⎼⎼⎼⎼⎼⎼⎼⎼
                 aD           aD   aN
                       aD y[i]   xD

Вам не нужно хранить все значения y в памяти; достаточно запомнить последний и прекратить повторение, когда их разница достаточно мала.

У вас все еще есть возведение в степень выше, но теперь это только целочисленное возведение в степень, т.е.

  aD
xN   = xN × xN × .. × xN
       ╰───────┬───────╯
              aD times

что вы можете сделать очень просто, например, просто умножив аргумент сам на себя желаемое количество раз, например в C,

double ipow(const double base, const int exponent)
{
    double result = 1.0;
    int    i;
    for (i = 0; i < exponent; i++)
        result *= base;
    return result;
}

Существуют более эффективные методы для выполнения целочисленного возведения в степень, но указанная выше функция должна быть вполне приемлемой для этого.

Последняя проблема - выбрать начальное y, чтобы получить сходимость. Вы не можете использовать 0, потому что (степень) y используется в качестве знаменателя при делении; вы получите деление на нулевую ошибку. Лично я проверял, должен ли результат быть положительным или отрицательным, а также быть меньше или больше единицы по величине; два общих правила, чтобы выбрать безопасное начальное y.

Вопросы?

Вы можете использовать обобщенную биномиальную теорему. Заменить y=1 и x=a-1. Вы можете усечь бесконечный ряд после достаточного количества членов, исходя из желаемой точности. Чтобы иметь возможность связать количество терминов с точностью, вам необходимо убедиться, что x^r термины уменьшаются по абсолютной величине. Итак, в зависимости от значений a и n, вы должны применить формулу для вычисления одного из a^n и a^(-n) и использовать это для получения желаемого результата.

Решение возведения целого числа в степень:

int poweri (int x, unsigned int y)
{
    int temp;
    if (y == 0)
        return 1;

    temp = poweri (x, y / 2);
    if ((y % 2) == 0)
        return temp * temp;
    else
        return x * temp * temp;
}

Однако квадратный корень не дает такого чистого закрытого решения. На странице wikipedia-square root и Алгоритмы квадратного корня Wolfram Mathworks Оба предоставляют несколько методов, которые будут удовлетворить ваши потребности, вам просто нужно выбрать тот, который соответствует вашим целям.

С небольшими изменениями эта процедура из Википедии (измененная для возврата квадратного корня и повышения точности) возвращает удивительно точный квадратный корень. Да, будут вопли об использовании объединения, и это действительно только тогда, когда целочисленное хранилище и хранилище с плавающей запятой эквивалентны, но если вы взламываете свой собственный квадратный корень, это относительно эффективно:

float sqrt_f (float x)
{
        float xhalf = 0.5f*x;
        union
        {
            float x;
            int i;
        } u;
        u.x = x;
        u.i = 0x5f3759df - (u.i >> 1);
        /* The next line can be repeated any number of times to increase accuracy */
        // u.x = u.x * (1.5f - xhalf * u.x * u.x);
        int i = 10;
        while (i--)
            u.x *= 1.5f - xhalf * u.x * u.x;

        return 1.0f / u.x;
}