Computations with decimals (arbitrary precision)

New to Julia
1

Hello everybody !
I’m busy now with studying the speed of convergence of sequences of real numbers. I begin to study classical algorithms like Archimedes polygonal approximation of pi, and newton computations of a square root. I already illustrated with some simple python programs and I want to translate into Julia . Here’s an example :

from decimal import Decimal, getcontext 
  
def u(n): 
    if n == 0: 
        return Decimal(1) 
    return (u(n - 1) + (2 - u(n - 1) * u(n - 1)) / (2 * u(n - 1)))  #appel récursif 
 
 
def main(): 
    """affichage des 10 premiers termes avec 100 décimales""" 
    getcontext().prec = 101 
    for i in range(0, 8): 
        print("i=", i, ":", u(i), end='\n') 
 
if __name__ == '__main__': 
    main()

The output is this :

i= 0 : 1
i= 1 : 1.5
i= 2 : 1.4166666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666667
i= 3 : 1.4142156862745098039215686274509803921568627450980392156862745098039215686274509803921568627450980392
i= 4 : 1.4142135623746899106262955788901349101165596221157440445849050192000543718353892683589900431576443402
i= 5 : 1.4142135623730950488016896235025302436149819257761974284982894986231958242289236217849418367358303566
i= 6 : 1.4142135623730950488016887242096980785696718753772340015610131331132652556303399785317871612507104752
i= 7 : 1.4142135623730950488016887242096980785696718753769480731766797379907324784621070388503875343276416016

    √2=1.4142135623730950488016887242096980785696718753769480731766797379907324784621070388503875343276415727

I use colors to show at every step how many exact decimals we get. Theoretical discussion follows leading to acceleration processes.

The base for a Julia equivalent would be :

using Decimals
#setprecision(1600)-evaluation of a size in bits to contain such data
u(n)= n==0 ? decimal(1) : (u(n - 1) + (2 - u(n - 1) * u(n - 1)) / (2 * u(n - 1)))
for i in 0:7
    print("i=",i,":",u(i)); println("")
end

So my problem is:
What is the julia equivalent for python
getcontext().prec = 101
???
I found some discussions here and there about this but the concern was mainly to produce big integers and I am familiar with this technique even with julia.
First I couldn’t appreciate the effect of ‘setprecision’ whatever the value given I always have the same result , without any error.

i=0:1
i=1:1.5
i=2:1.41666666666666666667
i=3:1.41421568627450980392
i=4:1.41421356237468991063
i=5:1.4142135623730950488
i=6:1.4142135623730950488
i=7:1.4142135623730950488

which is not enough for this special case and others to come later.
Is it just a matter of parametrization of println using a default rounding , or something more serious.
Note: specifying initial term as decimal(1) instead of simple 1 (Int64) gives 3 additional places of decimals no more…

2

I would use the built-in BigFloat type, not the Decimals package — you just need arbitrary-precision floating-point for this, not specifically decimal floating-point. The setprecision function by default applies to the BigFloat type.

For example:

function u(n)
    val = BigFloat(1)
    for i = 1:n
       val = val + (2 - val^2) / 2val
    end
    return val
end

and you can then do e.g.:

julia> setprecision(1024) # about 300 decimal digits
1024

julia> u.(0:7) # results for u(0) .. u(7)
8-element Vector{BigFloat}:
 1.0
 1.5
 1.416666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666663
 1.414215686274509803921568627450980392156862745098039215686274509803921568627450980392156862745098039215686274509803921568627450980392156862745098039215686274509803921568627450980392156862745098039215686274509803921568627450980392156862745098039215686274509803921568627450980392156862745098039215686274509803921
 1.414213562374689910626295578890134910116559622115744044584905019200054371835389268358990043157644340231759948346756380195058959459000237876779828049070581438814693988513949774017059163353382947633126040710911747714683793794814286199748530261324633839671050395894926428110238896251741597852312502123899819893296
 1.41421356237309504880168962350253024361498192577619742849828949862319582422892362178494183673583035655831431067505948202657730831732002621577830222200111296051917731217809008141213275773401309092053579529061462415958396989506279589018810128467136090375399417852437387204925332371777716921286946425937818658863
 1.414213562373095048801688724209698078569671875377234001561013133113265255630339978531787161250710475216048375111261837640718851266168116167162101238843783439868142965883500768770461491093225893183284963099150622339263287815909281327778653897986858072121764773149281450019241702294721718190060345552963844612524
 1.414213562373095048801688724209698078569671875376948073176679737990732478462107038850387534327641601639785778384557829824991246370587677393279969296443529709905820330039117363037412424448038077011337853455199672296933820918084119855070984754438400660813977905907724148608431057295333625478998730127918890051192

See also this notebook for more Julia code and explanations about Newton’s method for square roots.

3 Likes
3

Thank you very much stevengj for your quick and accurate answer.
Exactly what I needed when precision set to 334.
And of course to compare with ‘real’ value of the limit

println(sqrt(BigFloat(2)))

with the same settings.
Thanks again for your help!