Jamie M
Jamie M

Reputation: 127

Sympy precision using rationals

In sympy using large rationals, ie


xa2=2869041017039531/549755813888

xtotatives=5221

print(Rational(2869041017039531/549755813888)*5221)

For a larger set of similar fractions I am getting unexpected output when I check the number of distinct fractions. Is there a way to increase the precision in sympy for rationals using simple code like this?

Edit1:

import math
import numpy as np
from collections import Counter
from sympy import Symbol, Rational, fraction
#from decimal import *

#getcontext().prec = 10000 #digits of precision for decimal


A060753=[1, 2, 3, 15, 35, 77, 1001]
A038110=[1, 1, 1, 4, 8, 16, 192]

totatives=[]

#primesList=[2,3,5,7]
#primesList=[2,3,5,7,11]
primesList=[2,3,5,7,11,13]

primeProduct=math.prod(primesList)
nToUse=len(primesList)

valuesToCheck=range(1,primeProduct)
for n in valuesToCheck:
    if math.gcd(n, primeProduct)==1 and n<primeProduct:
        totatives.append(n)

print(len(totatives))

correctOutput=[]
incorrectOutput=[]
a2=[]

k=0
totativeCount=len(totatives)
while k<totativeCount:
    #if k%round(totativeCount/50)==1:
        #print(f"loop {k+1} of {totativeCount}")
    
    a2.append(A060753[nToUse]/A038110[nToUse]*(k+1))

    correctOutput.append(A060753[nToUse]*(k+1)-A038110[nToUse]*totatives[k])
        
    incorrectOutput.append((Rational(a2[k])-Rational(totatives[k]))*A038110[nToUse])

    k+=1


print(f"count of distinct correct: {len(np.unique(correctOutput))}")
print(f"sum(correctOutput): {sum(correctOutput)}")

print(f"count of distinct incorrect: {len(np.unique(incorrectOutput))}")
print(f"sum(incorrectOutput): {float(sum(incorrectOutput))}")

Output:

count of distinct correct: 2422
sum(correctOutput): 2882880
count of distinct incorrect: 3408
sum(incorrectOutput): 2882880.000000864

Cheers, Jamie

Upvotes: 1

Views: 244

Answers (1)

asmeurer
asmeurer

Reputation: 91620

SymPy doesn't use a global precision like decimal. Every SymPy Float object stores its own precision (defaulting to 15 digits).

If you want to use rational numbers, your best bet is to just avoid floats entirely. This can be done by ensuring your integers are SymPy integers. This can be done by calling sympify() on values that start out as Python integers, like

A060753 = sympify([1, 2, 3, 15, 35, 77, 1001])
A038110 = sympify([1, 1, 1, 4, 8, 16, 192])

And also using sympy.prod and sympy.gcd instead of math.prod and math.gcd. This avoids the gotcha where dividing two Python int objects produces a float, which loses information about the exact rational number it represents.

If you do this, the results will be rational numbers. You can then convert these to floats with as many digits as you want with number.evalf()

Upvotes: 1

Related Questions