Генерация цифр квадратного корня из 2

Question

Я хочу сгенерировать цифры квадратного корня из двух до 3 миллиона цифр.

Мне известно о Ньютоне-Рафсоне , но я не очень разбираюсь, как реализовать его в C или C ++ из-за отсутствия поддержки biginteger. Может ли кто-нибудь указать мне правильное направление?

Кроме того, если кто-нибудь знает, как это сделать на python (я новичок), я также был бы признателен за это.

Answer 1

Вы можете попробовать использовать сопоставление:

a/b -> (a+2b)/(a+b), начиная с a= 1, b= 1.Это сходится к sqrt (2) (фактически дает его непрерывные дробные представления).

Теперь ключевой момент: это можно представить в виде умножения матриц (аналогично Фибоначчи)

Еслиa_n и b_n - это n-е числа в шагах, тогда

[1 2] [a_n b_n] ^T = [a_ (n + 1) b_ (n + 1)] ^T [1 1]

, который теперь дает нам

[1 2] ⁿ [a_1 b_1] ^T = [a_ (n + 1)B_ (п + 1)] ^T [1 1]

Таким образом, если матрица 2x2 равна A, нам нужно вычислить A ⁿ , что можно сделать путем многократного возведения в квадрат и использовать только целочисленную арифметику (поэтому вам не нужнобеспокоиться о проблемах точности).

Также обратите внимание, что a / b, который вы получите, всегда будет в уменьшенной форме (как gcd (a, b) = gcd (a + 2b, a + b)), так что есливы думаете об использовании класса дроби для представления промежуточных результатов, не надо!

Поскольку n-й знаменатель подобен (1 + sqrt (2)) ^ n, для получения 3 миллионов цифр вам, скорее всего, понадобитсявычислять до 3671656 ^th term.

Примечание. Даже если вы ищете ~ 3,6-миллионный термин, повторное возведение в квадрат позволит вам вычислить n-й член в O (Log n)умножения и сложения.

Кроме того, это легко сделать параллельным, в отличие от итеративных, таких как Ньютон-Рафсон и т. д.

Answer 2

РЕДАКТИРОВАТЬ : мне нравится эта версия лучше, чем предыдущая. Это общее решение, которое принимает как целые, так и десятичные дроби; при n = 2 и точности = 100000 это занимает около двух минут. Спасибо Полу Макгуайру за его предложения и другие предложения!

def sqrt_list(n, precision):
    ndigits = []        # break n into list of digits
    n_int = int(n)
    n_fraction = n - n_int

    while n_int:                            # generate list of digits of integral part
        ndigits.append(n_int % 10)
        n_int /= 10
    if len(ndigits) % 2: ndigits.append(0)  # ndigits will be processed in groups of 2

    decimal_point_index = len(ndigits) / 2  # remember decimal point position
    while n_fraction:                       # insert digits from fractional part
        n_fraction *= 10
        ndigits.insert(0, int(n_fraction))
        n_fraction -= int(n_fraction)
    if len(ndigits) % 2: ndigits.insert(0, 0)  # ndigits will be processed in groups of 2

    rootlist = []
    root = carry = 0                        # the algorithm
    while root == 0 or (len(rootlist) < precision and (ndigits or carry != 0)):
        carry = carry * 100
        if ndigits: carry += ndigits.pop() * 10 + ndigits.pop()
        x = 9
        while (20 * root + x) * x > carry:
                x -= 1
        carry -= (20 * root + x) * x
        root = root * 10 + x
        rootlist.append(x)
    return rootlist, decimal_point_index

Answer 3

Наилучшим способом, вероятно, является использование продолжения расширения дроби [1; 2, 2, ...] квадратный корень из двух.

def root_two_cf_expansion():
    yield 1
    while True:
        yield 2

def z(a,b,c,d, contfrac):
    for x in contfrac:
        while a > 0 and b > 0 and c > 0 and d > 0:
            t = a // c
            t2 = b // d
            if not t == t2:
                break
            yield t
            a = (10 * (a - c*t))
            b = (10 * (b - d*t))
            # continue with same fraction, don't pull new x
        a, b = x*a+b, a
        c, d = x*c+d, c
    for digit in rdigits(a, c):
        yield digit

def rdigits(p, q):
    while p > 0:
        if p > q:
           d = p // q
           p = p - q * d
        else:
           d = (10 * p) // q
           p = 10 * p - q * d
        yield d

def decimal(contfrac):
    return z(1,0,0,1,contfrac)

decimal((root_two_cf_expansion()) возвращает итератор всех десятичных цифр. t1 и t2 в алгоритме - это минимальное и максимальное значения следующей цифры. Когда они равны, мы выводим эту цифру.

Обратите внимание, что это не обрабатывает некоторые исключительные случаи, такие как отрицательные числа в непрерывной дроби.

(Этот код является адаптацией кода на Haskell для обработки непрерывных дробных фрагментов.)

Answer 4

на работу?Используйте библиотеку!

Для развлечения?Хорошо для вас:)

Напишите программу, которая будет имитировать то, что вы будете делать с карандашом и бумагой.Начните с 1 цифры, затем 2 цифры, затем 3, ..., ...

Не беспокойтесь о Ньютоне или о ком-либо еще.Просто сделай это по-своему.

Answer 5

Что касается произвольных больших чисел, вы можете взглянуть на Арифметическая библиотека множественной точности GNU (для C / C ++).

Answer 6

Ну, вот код, который я написал. Он сгенерировал миллион цифр после десятичного числа для квадратного корня из 2 примерно за 60800 секунд для меня, но мой ноутбук спал, когда запускал программу, это должно быть быстрее. Вы можете попытаться сгенерировать 3 миллиона цифр, но это может занять пару дней.

def sqrt(number,digits_after_decimal=20):
import time
start=time.time()
original_number=number
number=str(number)
list=[]
for a in range(len(number)):
    if number[a]=='.':
        decimal_point_locaiton=a
        break
    if a==len(number)-1:
        number+='.'
        decimal_point_locaiton=a+1
if decimal_point_locaiton/2!=round(decimal_point_locaiton/2):
    number='0'+number
    decimal_point_locaiton+=1
if len(number)/2!=round(len(number)/2):
    number+='0'
number=number[:decimal_point_locaiton]+number[decimal_point_locaiton+1:]
decimal_point_ans=int((decimal_point_locaiton-2)/2)+1
for a in range(0,len(number),2):
    if number[a]!='0':
        list.append(eval(number[a:a+2]))
    else:
        try:
            list.append(eval(number[a+1]))
        except IndexError:
            pass
p=0
c=list[0]
x=0
ans=''
for a in range(len(list)):
    while c>=(20*p+x)*(x):
        x+=1
    y=(20*p+x-1)*(x-1)
    p=p*10+x-1
    ans+=str(x-1)
    c-=y
    try:
        c=c*100+list[a+1]
    except IndexError:
        c=c*100
while c!=0:
    x=0
    while c>=(20*p+x)*(x):
        x+=1
    y=(20*p+x-1)*(x-1)
    p=p*10+x-1
    ans+=str(x-1)
    c-=y
    c=c*100
    if len(ans)-decimal_point_ans>=digits_after_decimal:
            break
ans=ans[:decimal_point_ans]+'.'+ans[decimal_point_ans:]
total=time.time()-start
return ans,total

Answer 7

Вот краткая версия для вычисления квадратного корня из целого числа a до цифр точности.Он работает, находя целочисленный квадратный корень из a после умножения на 10, повышенного до 2 x цифр .

def sqroot(a, digits):
    a = a * (10**(2*digits))
    x_prev = 0
    x_next = 1 * (10**digits)
    while x_prev != x_next:
        x_prev = x_next
        x_next = (x_prev + (a // x_prev)) >> 1
    return x_next

Всего несколько предостережений.

Вам нужно будет преобразовать результат в строку и добавить десятичную точку в правильном месте (если вы хотите, чтобы десятичная точка была напечатана).

Преобразование очень большого целого числа в строку не являетсяочень быстро.

Деление очень больших целых чисел тоже не очень быстро (в Python).

В зависимости от производительности вашей системы, для вычисления квадратного корня может потребоваться час или большеот 2 до 3 миллионов знаков после запятой.

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

Answer 8

Вот более эффективная целочисленная функция квадратного корня (в Python 3.x), которая должна завершаться во всех случаях.Он начинается с числа, расположенного намного ближе к квадратному корню, поэтому требуется меньше шагов.Обратите внимание, что для int.bit_length требуется Python 3.1+.Ошибка проверки пропущена для краткости.

def isqrt(n):
    x = (n >> n.bit_length() // 2) + 1
    result = (x + n // x) // 2
    while abs(result - x) > 1:
        x = result
        result = (x + n // x) // 2
    while result * result > n:
        result -= 1
    return result

Answer 9

Python уже поддерживает большие целые числа из коробки, и если это единственное, что сдерживает вас в C / C ++, вы всегда можете написать быстрый контейнерный класс самостоятельно.

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