Как вычислить n-й корень из очень большого целого числа

Вы можете заставить его работать немного быстрее, избегая циклов while в пользу установки низкого значения на 10 ** (len (str (x)) / n) и высокого на низкое значение * 10. Вероятно, лучше заменить len (str (x) )) с битовой длиной и использованием сдвига битов. Основываясь на моих тестах, я оцениваю ускорение на 5% от первого и ускорение на 25% со второго. Если целые числа достаточно большие, это может иметь значение (и ускорения могут отличаться). Не доверяйте моему коду без тщательного тестирования. Я провел некоторое базовое тестирование, но, возможно, пропустил крайний случай. Кроме того, эти ускорения зависят от выбранного количества.

Если фактические данные, которые вы используете, намного больше, чем вы разместили здесь, это изменение может стоить.

from timeit import Timer

def find_invpow(x,n):
    """Finds the integer component of the n'th root of x,
    an integer such that y ** n <= x < (y + 1) ** n.
    """
    high = 1
    while high ** n < x:
        high *= 2
    low = high/2
    while low < high:
        mid = (low + high) // 2
        if low < mid and mid**n < x:
            low = mid
        elif high > mid and mid**n > x:
            high = mid
        else:
            return mid
    return mid + 1

def find_invpowAlt(x,n):
    """Finds the integer component of the n'th root of x,
    an integer such that y ** n <= x < (y + 1) ** n.
    """
    low = 10 ** (len(str(x)) / n)
    high = low * 10

    while low < high:
        mid = (low + high) // 2
        if low < mid and mid**n < x:
            low = mid
        elif high > mid and mid**n > x:
            high = mid
        else:
            return mid
    return mid + 1

x = 237734537465873465
n = 5
tests = 10000

print "Norm", Timer('find_invpow(x,n)', 'from __main__ import find_invpow, x,n').timeit(number=tests)
print "Alt", Timer('find_invpowAlt(x,n)', 'from __main__ import find_invpowAlt, x,n').timeit(number=tests)

Норма 0.626754999161

Alt 0.566340923309

Попробуйте преобразовать экспоненту в плавающее число, так как в Python поведение по умолчанию - целочисленное деление

п ** (1 / поплавок (3))

В более старых версиях Python 1/3равен 0. В Python 3.0 1/3равен 0.33333333333 (и 1//3равен 0).

Итак, либо измените свой код для использования, 1/3.0либо переключитесь на Python 3.0.

Gmpy - это C-кодированный модуль расширения Python, который оборачивает библиотеку GMP, чтобы обеспечить коду Python быструю арифметику с множественными точностями (целочисленные, рациональные и с плавающей точкой), генерацию случайных чисел, расширенные теоретические функции и многое другое.

Включает rootфункцию:

x.root (n): возвращает 2-элементный кортеж (y, m), такой что y является (возможно, усеченным) n-м корнем x; m, обычный Python int, равен 1, если корень точный (x == y ** n), иначе 0. n должен быть обычным Python int,> = 0.

Например, 20-й корень:

>>> import gmpy
>>> i0=11968003966030964356885611480383408833172346450467339251 
>>> m0=gmpy.mpz(i0)
>>> m0
mpz(11968003966030964356885611480383408833172346450467339251L)
>>> m0.root(20)
(mpz(567), 0)

Если это действительно большое число. Вы можете использовать бинарный поиск.

def find_invpow(x,n):
    """Finds the integer component of the n'th root of x,
    an integer such that y ** n <= x < (y + 1) ** n.
    """
    high = 1
    while high ** n <= x:
        high *= 2
    low = high/2
    while low < high:
        mid = (low + high) // 2
        if low < mid and mid**n < x:
            low = mid
        elif high > mid and mid**n > x:
            high = mid
        else:
            return mid
    return mid + 1

Например:

>>> x = 237734537465873465
>>> n = 5
>>> y = find_invpow(x,n)
>>> y
2986
>>> y**n <= x <= (y+1)**n
True
>>>
>>> x = 119680039660309643568856114803834088331723464504673392511960931441>
>>> n = 45
>>> y = find_invpow(x,n)
>>> y
227661383982863143360L
>>> y**n <= x < (y+1)**n
True
>>> find_invpow(y**n,n) == y
True
>>>

Что ж, если вас не особо беспокоит точность, вы можете преобразовать ее в строчку, отрубить несколько цифр, использовать функцию экспоненты, а затем умножить результат на корень того, сколько вы отрубили.

Например, 32123 примерно равно 32 * 1000, кубический корень примерно эквивалентен кубическому корню из 32 * кубического корня из 1000. Последний может быть рассчитан путем деления числа 0 на 3.

Это исключает необходимость использования модулей расширения.

О, для чисел , что большой, вы бы использовать модуль десятичного.

нс: ваш номер в виде строки

ns = "11968003966030964356885611480383408833172346450467339251196093144141045683463085291115677488411620264826942334897996389485046262847265769280883237649461122479734279424416861834396522819159219215308460065265520143082728303864638821979329804885526557893649662037092457130509980883789368448042961108430809620626059287437887495827369474189818588006905358793385574832590121472680866521970802708379837148646191567765584039175249171110593159305029014037881475265618958103073425958633163441030267478942720703134493880117805010891574606323700178176718412858948243785754898788359757528163558061136758276299059029113119763557411729353915848889261125855717014320045292143759177464380434854573300054940683350937992500211758727939459249163046465047204851616590276724564411037216844005877918224201569391107769029955591465502737961776799311859881060956465198859727495735498887960494256488224613682478900505821893815926193600121890632"
from decimal import Decimal
d = Decimal(ns)
one_third = Decimal("0.3333333333333333")
print d ** one_third

и ответ: 2.287391878618402702753613056E + 305

Т.З. указал, что это не точно ... и он прав. Вот мой тест.

from decimal import Decimal

def nth_root(num_decimal, n_integer):
    exponent = Decimal("1.0") / Decimal(n_integer)
    return num_decimal ** exponent

def test():
    ns = "11968003966030964356885611480383408833172346450467339251196093144141045683463085291115677488411620264826942334897996389485046262847265769280883237649461122479734279424416861834396522819159219215308460065265520143082728303864638821979329804885526557893649662037092457130509980883789368448042961108430809620626059287437887495827369474189818588006905358793385574832590121472680866521970802708379837148646191567765584039175249171110593159305029014037881475265618958103073425958633163441030267478942720703134493880117805010891574606323700178176718412858948243785754898788359757528163558061136758276299059029113119763557411729353915848889261125855717014320045292143759177464380434854573300054940683350937992500211758727939459249163046465047204851616590276724564411037216844005877918224201569391107769029955591465502737961776799311859881060956465198859727495735498887960494256488224613682478900505821893815926193600121890632"
    nd = Decimal(ns)
    cube_root = nth_root(nd, 3)
    print (cube_root ** Decimal("3.0")) - nd

if __name__ == "__main__":
    test()

Это примерно на 10 ** 891

Возможно, для вашего любопытства:

http://en.wikipedia.org/wiki/Hensel_Lifting

Это может быть техника, которую Maple будет использовать, чтобы найти n-й корень больших чисел.

Представьте тот факт x^n - 11968003.... = 0 mod p, и идти оттуда ...

Если вы ищете что-то стандартное, пишите быстро с высокой точностью. Я хотел бы использовать десятичную и настроить точность (getcontext (). Prec) по крайней мере до длины х.

Код (Python 3.0)

from decimal import *

x =   '11968003966030964356885611480383408833172346450467339251\
196093144141045683463085291115677488411620264826942334897996389\
485046262847265769280883237649461122479734279424416861834396522\
819159219215308460065265520143082728303864638821979329804885526\
557893649662037092457130509980883789368448042961108430809620626\
059287437887495827369474189818588006905358793385574832590121472\
680866521970802708379837148646191567765584039175249171110593159\
305029014037881475265618958103073425958633163441030267478942720\
703134493880117805010891574606323700178176718412858948243785754\
898788359757528163558061136758276299059029113119763557411729353\
915848889261125855717014320045292143759177464380434854573300054\
940683350937992500211758727939459249163046465047204851616590276\
724564411037216844005877918224201569391107769029955591465502737\
961776799311859881060956465198859727495735498887960494256488224\
613682478900505821893815926193600121890632'

minprec = 27
if len(x) > minprec: getcontext().prec = len(x)
else:                getcontext().prec = minprec

x = Decimal(x)
power = Decimal(1)/Decimal(3)

answer = x**power
ranswer = answer.quantize(Decimal('1.'), rounding=ROUND_UP)

diff = x - ranswer**Decimal(3)
if diff == Decimal(0):
    print("x is the cubic number of", ranswer)
else:
    print("x has a cubic root of ", answer)

Ответ

х представляет собой кубическое количество 22873918786185635329056863961725521583023133411 451452349318109627653540670761962215971994403670045614485973722724603798 107719978813658857014190047742680490088532895666963698551709978502745901 704433723567548799463129652706705873694274209728785041817619032774248488 2965377218610139128882473918261696612098418

Я придумал свой собственный ответ, который берет идею @Mahmoud Kassem, упрощает код и делает его более пригодным для повторного использования:

def cube_root(x):
    return decimal.Decimal(x) ** (decimal.Decimal(1) / decimal.Decimal(3))

Я протестировал его в Python 3.5.1 и Python 2.7.8, и, похоже, он работал нормально.

Результат будет иметь столько цифр, сколько указано в десятичном контексте, в котором выполняется функция, которая по умолчанию составляет 28 десятичных знаков. Согласно документации для powerфункции в decimalмодуле, « результат четко определен, но только« почти всегда правильно округлен ». ». Если вам нужен более точный результат, это можно сделать следующим образом:

with decimal.localcontext() as context:
    context.prec = 50
    print(cube_root(42))