Decimal k-th root

Ideas
1

Looking at this stackoverflow thread, it seems like Context.power() does not work for very large base values. Particularly when the exponent is non-integral. I understand that this rounding is expected, but for very very large values of X, the answer returned by power(X, Decimal(1.0)/(Decimal(3.0))) is very far off of the actual value. Is there a way that this can be rectified? Maybe with a k-th root algorithm that uses binary search/newton’s method to calculate the integer and decimal components of the root?

2

The decimal module can do this but you need to use a high precision for the 1/3 exponent:

In [11]: n = 1196800396603096435688561148038340883317234645046733925119609314414104568346308529111567
    ...: 74884116202648269423348979963894850462628472657692808832376494611224797342794244168618343965
    ...: 22819159219215308460065265520143082728303864638821979329804885526557893649662037092457130509
    ...: 98088378936844804296110843080962062605928743788749582736947418981858800690535879338557483259
    ...: 01214726808665219708027083798371486461915677655840391752491711105931593050290140378814752656
    ...: 18958103073425958633163441030267478942720703134493880117805010891574606323700178176718412858
    ...: 94824378575489878835975752816355806113675827629905902911311976355741172935391584888926112585
    ...: 57170143200452921437591774643804348545733000549406833509379925002117587279394592491630464650
    ...: 47204851616590276724564411037216844005877918224201569391107769029955591465502737961776799311
    ...: 85988106095646519885972749573549888796049425648822461368247890050582189381592619360012189063
    ...: 2

In [12]: import decimal

In [13]: ctx = decimal.getcontext()

In [14]: ctx2 = ctx.copy()

In [15]: ctx.prec = 400

In [16]: ctx2.prec = 10000

In [17]: third = ctx2.divide(1, 3)

In [18]: r = ctx.power(n, third)

In [19]: r
Out[19]: Decimal('228739187861856353290568639617255215830231334114514523493181096276535406707619622159719944036700456144859737227246037981077199788136588570141900477426804900885328956669636985517099785027459017044337235675487994631296527067058736942742097287850418176190327742484882965377218610139128882473918261696612098418.0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000')

For comparison some exact resultd using SymPy:

In [20]: sympy.integer_nthroot(n, 3)
Out[20]: 
(228739187861856353290568639617255215830231334114514523493181096276535406707619622159719944036700456144859737227246037981077199788136588570141900477426804900885328956669636985517099785027459017044337235675487994631296527067058736942742097287850418176190327742484882965377218610139128882473918261696612098418,
 True)

In [22]: print(sympy.Integer(n) ** sympy.Rational(1, 3))
228739187861856353290568639617255215830231334114514523493181096276535406707619622159719944036700456144859737227246037981077199788136588570141900477426804900885328956669636985517099785027459017044337235675487994631296527067058736942742097287850418176190327742484882965377218610139128882473918261696612098418

Obviously an nth root function could take care of ensuring high enough precision and could be more efficient/accurate than using general pow (SymPy’s evalf already takes care of ensuring high enough precision if approximate roots are needed). Other libraries like SymPy, mpmath, gmpy2 etc can handle this better than decimal and provide many different alternatives depending on whether exact or approximate integer or non-integer roots are wanted. This sort of thing likely goes beyond the original scope of the decimal module since it does not really have anything to do with decimal vs binary floating point.

1 Like