pythondecimalprecision

How to make python use more than 15 decimals in calculations?


In this question I made a claim that this construction of pentagons will always converge.

I drew the first 50 iteration using Geogebra here. However, zooming in on the 50th iteration resulted in significant lag and imprecise rendering, hindering further exploration. Additionally, manually coding the constructions in GeoGebra was time-consuming (over 6 hours for 50 iterations).

I needed to use another tool than GeoGebra to calculate the pentagons so I wrote this Python code:

import math 
# enter your coordinates here
A=[0,0]
B=[0,0]
C=[0,0]
D=[0,0]
E=[0,0]
print("Number of pentagons")
n=int(input() )
print("Enter A_1")
A[0]=float (input() )
A[1]=float (input() )
print("Enter B_1")
B[0]=float (input() )
B[1]=float (input() )
print("Enter C_1")
C[0]=float (input() )
C[1]=float (input() )
print("Enter D_1")
D[0]=float (input() )
D[1]=float (input() )
print("Enter E_1")
E[0]=float (input() )
E[1]=float (input() )
X=[0,0]
Y=[0,0]
Z=[0,0]
W=[0,0]
T=[0,0]


#defining some useful functions
def v(a,b):
    result= [b[0]-a[0], b[1]-a[1]]
    return result 
def crs(a, b):
    result = [-0.5*((a[0]*b[1])-(a[1]*b[0]))]
    return result
def dis(a,b):
    result =[math.sqrt((a[0]-b[0])**2+(a[1]-b[1])**2)]
    return result 




for i in range (2,n+1):
    AB= dis(A,B)[0]
    BC= dis(B,C)[0]
    CD= dis(C,D)[0]
    DE= dis(D,E)[0]
    EA= dis(E,A)[0]



    #Here to calculate A_n
    x=   CD*crs(v(D,C),v(D,E))[0]
    y=  -CD*crs(v(E,D),v(E,B))[0]+  2*DE*crs(v(C,B),v(C,D))[0]
    z=  -CD*crs(v(C,B),v(C,E))[0]+  2*BC*crs(v(D,C),v(D,E))[0]
    t=   CD*crs(v(C,B),v(C,D))[0]
    X[0]=(x *B[0]+y*C[0]+z*D[0]+t*E[0])/(x+y+z+t)
    X[1]=(x *B[1]+y*C[1]+z*D[1]+t*E[1])/(x+y+z+t)



    #Here to calculate B_n
    x=   DE*crs(v(E,D),v(E,A))[0]
    y=  -DE*crs(v(A,E),v(A,C))[0]+  2*EA*crs(v(D,C),v(D,E))[0]
    z=  -DE*crs(v(D,C),v(D,A))[0]+  2*CD*crs(v(E,D),v(E,A))[0]
    t=   DE*crs(v(D,C),v(D,E))[0]
    Y[0]=(x *C[0]+y*D[0]+z*E[0]+t*A[0])/(x+y+z+t)
    Y[1]=(x *C[1]+y*D[1]+z*E[1]+t*A[1])/(x+y+z+t)



    #Here to calculate C_n
    x=   EA*crs(v(A,E),v(A,B))[0]
    y=  -EA*crs(v(B,A),v(B,D))[0]+  2*AB*crs(v(E,D),v(E,A))[0]
    z=  -EA*crs(v(E,D),v(E,B))[0]+  2*DE*crs(v(A,E),v(A,B))[0]
    t=   EA*crs(v(E,D),v(E,A))[0]
    Z[0]=(x *D[0]+y*E[0]+z*A[0]+t*B[0])/(x+y+z+t)
    Z[1]=(x *D[1]+y*E[1]+z*A[1]+t*B[1])/(x+y+z+t)



    #Here to calculate D_n
    x=   AB*crs(v(B,A),v(B,C))[0]
    y=  -AB*crs(v(C,B),v(C,E))[0]+  2*BC*crs(v(A,E),v(A,B))[0]
    z=  -AB*crs(v(A,E),v(A,C))[0]+  2*EA*crs(v(B,A),v(B,C))[0]
    t=   AB*crs(v(A,E),v(A,B))[0]
    W[0]=(x *E[0]+y*A[0]+z*B[0]+t*C[0])/(x+y+z+t)
    W[1]=(x *E[1]+y*A[1]+z*B[1]+t*C[1])/(x+y+z+t)



    #Here to calculate E_n
    x=   BC*crs(v(C,B),v(C,D))[0]
    y=  -BC*crs(v(D,C),v(D,A))[0]+  2*CD*crs(v(B,A),v(B,C))[0]
    z=  -BC*crs(v(B,A),v(B,D))[0]+  2*AB*crs(v(C,B),v(C,D))[0]
    t=   BC*crs(v(B,A),v(B,C))[0]
    T[0]=(x *A[0]+y*B[0]+z*C[0]+t*D[0])/(x+y+z+t)
    T[1]=(x *A[1]+y*B[1]+z*C[1]+t*D[1])/(x+y+z+t)
    A=X
    B=Y
    C=Z
    D=W
    E=T

    
    print(f"A_{i}= {A}")
    print(f"B_{i}= {B}")
    print(f"C_{i}= {C}")
    print(f"D_{i}= {D}")
    print(f"E_{i}= {E}")  

The problem is that after about 20 or 30 iterations, Python gave me this error:

ZeroDivisionError: float division by zero.

I figured that this has to do with the fact that after 20-30 iterations the first 15 digits will be the same, so I need Python to use more than 50 decimals in all of these calculations.

I tried to google "how to make Python use more decimals" I tried some solutions but none of them worked for me (I'm a beginner in coding).


Solution

  • HThe builtin decimal package provides support for this precisely. A simple example adapted from the docs:

    >>> from decimal import Decimal
    >>> import decimal
    >>> decimal.getcontext().prec
    28
    >>> Decimal(1)/7
    Decimal('0.1428571428571428571428571429')
    >>> decimal.getcontext().prec = 6
    >>> Decimal(1)/7
    Decimal('0.142857')
    

    You will need to adapt your library to box numbers with Decimal though. Operations must also be amended to use the library:

    math.sqrt(x)
    

    becomes

    x.sqrt()
    

    assuming x was buzzed in a Decimal.

    If you are using Numpy which is always a good idea for vectorized implementations you can't choose arbitrary precision but there are data types that are more accurate:

    https://numpy.org/doc/stable/user/basics.types.html#relationship-between-numpy-data-types-and-c-data-data-types

    for example

    https://numpy.org/doc/stable/reference/arrays.scalars.html#numpy.longdouble (np.float128)

    While this is not as good as arbitrary precision, this uses platform supported precision so is very fast in computations.