davo36
davo36

Reputation: 714

Arbitrary Precision Arithmetic in Julia

This has kinda been asked, but not in this way. I have a little Python program which finds continued fractions for square roots of n (1 <= n <= 10000).

I have been trying to do this in Julia and I can't see how to. Mainly because it deals with irrational numbers (sqrt(x) is irrational if x is not a perfect square, e.g. sqrt(2) = 1.414213...). So I don't think I can use the rational class.

It says here https://docs.julialang.org/en/latest/manual/integers-and-floating-point-numbers/#Arbitrary-Precision-Arithmetic-1 that Julia can do arbitrary precision arithmetic with BigFloats. But they don't seem to be accurate enough.

I have also tried to use PyCall and the Decimals package in Python (from Julia), but get weird errors (I can post them if they would be helpful).

Here is my Python program which works. And my question is how to do this in Julia please?

def continuedFractionSquareRoots():
''' 
  For each number up to 100, get the length of the continued fraction 
  of the square root for it.
'''

decimal.getcontext().prec = 210 # Need lots of precision for this problem.

continuedFractionLengths = []
for i in range(1, 101):

    # For perfect squares, the period is 0
    irrationalNumber = decimal.Decimal(i).sqrt()
    if irrationalNumber == int(irrationalNumber):
        continue        

    continuedFractionLength = 0
    while True:

        intPart = irrationalNumber.to_integral_exact(rounding=decimal.ROUND_FLOOR)
        if continuedFractionLength == 0:
            firstContinuedFractionTimes2 = int(intPart*2)

        continuedFractionLength += 1
        if intPart == firstContinuedFractionTimes2:
            # Have reached the 'period' end for this fraction
            break

        fractionalPart = irrationalNumber - intPart
        irrationalNumber = 1 / fractionalPart

continuedFractionLengths.append(continuedFractionLength)
return continuedFractionLengths

So as you can see, I need a way of computing an exact square root and also a way to get the integer part of the number. And that's all really, apart from lots and lots of precision!

Guys, I didn't post my Julia code because I didn't want to have a small manuscript for an answer! But here it is, it works. As I said in the comments below, I used the setprecision function to set the precision to a high value and it works. I got the value 711 empirically.

setprecision(711)

function continuedFractionSquareRoots()
#=
  For each number up to 100, get the length of the continued fraction 
  of the square root for it.
=#

continuedFractionLengths = Int[]
for i=1:100

    # For perfect squares, the period is 0
    irrationalNumber = BigFloat(sqrt(BigFloat(i)))
    if irrationalNumber == floor(irrationalNumber)
        continue
    end

    continuedFractionLength = 0
    while true

        intPart = floor(irrationalNumber)
        if continuedFractionLength == 0
            firstContinuedFractionTimes2 = intPart*2
        end

        continuedFractionLength += 1
        if intPart == firstContinuedFractionTimes2
            # Have reached the 'period' end for this fraction
            break
        end

        fractionalPart = irrationalNumber - intPart
        irrationalNumber = BigFloat(1) / fractionalPart

    end

    push!(continuedFractionLengths, continuedFractionLength)

end


return continuedFractionLengths
end

So anyway, user2357112 solved it, thanks very much.

Upvotes: 4

Views: 3639

Answers (2)

Tasos Papastylianou
Tasos Papastylianou

Reputation: 22225

user2357112's answer is at the heart of the problem and the correct answer to this question.

However, for the sake of completion, the 'literal' question of "how would I get this python script to run in julia" is an interesting question in itself, as it is not as straightforward as it might first seem due to this issue, so I'll show how to do this as well.

Assuming you have your python script called "testo.py" in your current directory (and with the correct import decimal statements etc), then here's how to import it as a python module and run the relevant function:

using PyCall
unshift!(PyVector(pyimport("sys")["path"]), ""); # as per https://github.com/JuliaPy/PyCall.jl#troubleshooting
testo = pyimport("testo");
testo[:oddPeriodSquareRoots]()  # will output '1322'

Upvotes: 1

user2357112
user2357112

Reputation: 280733

Just like with Python's decimal.Decimal, you can configure the precision of Julia's BigFloat:

setprecision(however_many_bits)

Note that this is in bits, unlike decimal.Decimal, because BigFloat doesn't use decimal.

Upvotes: 4

Related Questions