Как вычислить цифры иррационального числа одну за другой?

Question

Я хочу прочитать di git по di git десятичные дроби квадрата 5 в C. Квадрат root из 5 равен 2,23606797749979 ..., так что это будет ожидаемый результат:

2
3
6
0
6
7
9
7
7
...

Я нашел следующий код :

#include<stdio.h>

void main()
{
    int number;

    float temp, sqrt;

    printf("Provide the number: \n");

    scanf("%d", &number);

    // store the half of the given number e.g from 256 => 128
    sqrt = number / 2;
    temp = 0;

    // Iterate until sqrt is different of temp, that is updated on the loop
    while(sqrt != temp){
        // initially 0, is updated with the initial value of 128
        // (on second iteration = 65)
        // and so on
        temp = sqrt;

        // Then, replace values (256 / 128 + 128 ) / 2 = 65
        // (on second iteration 34.46923076923077)
        // and so on
        sqrt = ( number/temp + temp) / 2;
    }

    printf("The square root of '%d' is '%f'", number, sqrt);
}

Но этот подход сохраняет результат в переменной с плавающей запятой, и я не хочу зависеть от ограничений типов с плавающей запятой, например, я хотел бы извлечь, например, 10 000 цифр. Я также пытался использовать встроенную функцию sqrt () и приводить ее к строковому номеру, используя этот метод , но я столкнулся с той же проблемой.

Answer 1

То, о чем вы спрашивали, это очень сложная проблема, и возможно ли вообще делать это «один за другим» (т.е. без требования к рабочему пространству, которое зависит от того, насколько далеко вы хотите go) зависит как от конкретного иррационального числа , так и от базового значения , в котором вы хотите его представить. Например, в 1995 году, когда была обнаружена формула для числа пи, позволяющая вычислить n-й двоичный код di git в O (1) пробел , это было действительно большое дело. Это было не то, чего ожидали люди.

Если вы готовы принять O (n) пробел, то некоторые случаи, подобные упомянутому вами, довольно просты. Например, если у вас есть первые n цифр квадрата root числа в виде десятичной строки, вы можете просто попробовать добавить каждую ди git 0 к 9, а затем возвести в квадрат строку с длинным умножением (так же, как вы узнали в начальной школе), и выбирая последний, который не промахивается. Конечно, это очень медленно, но это просто. Самый простой способ сделать это намного быстрее (но все же асимптотически так же плохо) - использовать математическую библиотеку произвольной точности вместо строк. Достижение значительно лучших результатов требует более продвинутых подходов и в целом может оказаться невозможным.

Answer 2

Как уже отмечалось, вам нужно изменить алгоритм на di git -by-di git (на странице Wikipedia есть несколько примеров методов вычисления квадратных корней ) и использовать библиотеку c произвольной точности для выполнения вычислений (например, GMP ).

В следующем фрагменте я реализовал вышеупомянутый алгоритм, используя GMP (но не квадратная root функция, которую предоставляет библиотека). Вместо вычисления одного десятичного числа di git за один раз, эта реализация использует большее основание, наибольшее кратное 10, которое помещается в unsigned long, так что оно может производить 9 или 18 десятичных знаков на каждой итерации.

Он также использует адаптированный метод Ньютона, чтобы найти фактическое "di git".

#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <gmp.h>

unsigned long max_ul(unsigned long a, unsigned long b)
{
    return a < b ? b : a;   
}

int main(int argc, char *argv[])
{
    // The GMP functions accept 'unsigned long int' values as parameters.
    // The algorithm implemented here can work with bases other than 10,
    // so that it can evaluate more than one decimal digit at a time.
    const unsigned long base = sizeof(unsigned long) > 4
                             ? 1000000000000000000
                             : 1000000000;
    const unsigned long decimals_per_digit = sizeof(unsigned long) > 4 ? 18 : 9;

    // Extract the number to be square rooted and the desired number of decimal
    // digits from the command line arguments. Fallback to 0 in case of errors.
    const unsigned long number = argc > 1 ? atoi(argv[1]) : 0;
    const unsigned long n_digits = argc > 2 ? atoi(argv[2]) : 0;

    // All the variables used by GMP need to be properly initialized before use.
    // 'c' is basically the remainder, initially set to the original number
    mpz_t c;
    mpz_init_set_ui(c, number);

    // At every iteration, the algorithm "move to the left" by two "digits"
    // the reminder, so it multplies it by base^2.
    mpz_t base_squared;
    mpz_init_set_ui(base_squared, base);
    mpz_mul(base_squared, base_squared, base_squared);

    // 'p' stores the digits of the root found so far. The others are helper variables
    mpz_t p;
    mpz_init_set_ui(p, 0UL);    
    mpz_t y;
    mpz_init(y);
    mpz_t yy;
    mpz_init(yy);
    mpz_t dy;
    mpz_init(dy);
    mpz_t dx;
    mpz_init(dx);
    mpz_t pp;    
    mpz_init(pp);

    // Timing, for testing porpuses
    clock_t start = clock(), diff;

    unsigned long x_max = number;
    // Each "digit" correspond to some decimal digits
    for (unsigned long i = 0,
         last = (n_digits + decimals_per_digit) / decimals_per_digit + 1UL;
         i < last; ++i)
    {
        // Find the greatest x such that:  x * (2 * base * p + x) <= c
        // where x is in [0, base), using a specialized Newton method

        // pp = 2 * base * p
        mpz_mul_ui(pp, p, 2UL * base);

        unsigned long x = x_max;
        for (;;)
        {            
            // y = x * (pp + x)
            mpz_add_ui(yy, pp, x);
            mpz_mul_ui(y, yy, x);

            // dy = y - c
            mpz_sub(dy, y, c);

            // If y <= c we have found the correct x
            if ( mpz_sgn(dy) <= 0 )
                break;

            // Newton's step:  dx = dy/y'  where  y' = 2 * x + pp            
            mpz_add_ui(yy, yy, x);
            mpz_tdiv_q(dx, dy, yy);

            // Update x even if dx == 0 (last iteration)
            x -= max_ul(mpz_get_si(dx), 1);
        }        
        x_max = base - 1;

        // The actual format of the printed "digits" is up to you       
        if (i % 4 == 0)
        {
            if (i == 0)
                printf("%lu.", x);
            putchar('\n');
        }
        else
            printf("%018lu", x);

        // p = base * p + x
        mpz_mul_ui(p, p, base);
        mpz_add_ui(p, p, x);

        // c = (c - y) * base^2
        mpz_sub(c, c, y);
        mpz_mul(c, c, base_squared);
    }

    diff = clock() - start;
    long int msec = diff * 1000L / CLOCKS_PER_SEC;
    printf("\n\nTime taken: %ld.%03ld s\n", msec / 1000, msec % 1000);

    // Final cleanup
    mpz_clear(c);
    mpz_clear(base_squared);
    mpz_clear(p);
    mpz_clear(pp);
    mpz_clear(dx);
    mpz_clear(y);
    mpz_clear(dy);
    mpz_clear(yy);
}

Здесь вы можете увидеть выведенные цифры .

Answer 3

Ваш заголовок гласит:

Как вычислить цифры иррационального числа одну за другой?

Иррациональные числа не ограничены большинством квадратных корней. Они также включают числа вида log(x), exp(z), sin(y), et c. ( трансцендентные числа ). Однако есть некоторые важные факторы, которые определяют, насколько быстро вы можете вычислить цифры данного иррационального числа одну за другой (то есть слева направо).

• Не все иррациональные числа вычислимы; то есть никто не нашел способа приблизить их к любой желаемой длине (будь то выражение в виде замкнутой формы, ряд или иным образом).
• Существует множество способов выражения чисел, например, их двоичные или десятичные разложения, как непрерывные дроби, как ряды, et c. И существуют разные алгоритмы для вычисления цифр данного числа в зависимости от представления.
• Некоторые формулы вычисляют цифры данного числа в конкретной базе (например, в базе 2), а не в произвольной базе.

Например, помимо первой формулы для извлечения цифр числа пи без вычисления предыдущих цифр, существуют другие формулы этого типа (известные как формулы типа BBP ), которые извлекают цифры некоторых иррациональных чисел. Однако эти формулы работают только для конкретной базы, не все формулы типа BBP имеют формальное доказательство, и, что более важно, не все иррациональные числа имеют формулу типа BBP (по сути, только некоторые логарифмические и арктановые константы, а не числа форма exp(x) или sqrt(x)).

С другой стороны, если вы можете express иррациональное число в виде непрерывной дроби (которое есть у всех действительных чисел), вы может извлекать свои цифры слева направо и в любой требуемой базе, используя алгоритм c. Более того, этот алгоритм работает для любой действительной числовой константы, включая квадратные корни, экспоненты (e и exp(x)), логарифмы и т. Д. c., Если вы знаете, как express это как продолжение доля. Для реализации см. " Цифры пи и Python генераторы ". См. Также Код для генерации одного Di git за один раз .