Search code examples
pythonsympyarbitrary-precisionmpmath

Python Sympy and mpmath give different results

I did an evaluation of mpmath and sympy and found inconsistencies, which I cannot resolve. I believe mpmath gives correct results. Please, can somebody help in resolving this?

import sympy as sp
import mpmath as mp

# Set the desired precision to 120 decimal places using mpmath
mp.mp.dps = 30               

print('mp: ', mp.sin('0.4'))

# Define the variable and calculate the sine function
x = sp.symbols('x')

x = '0.4'
print('sp: ', (sp.sin(x).evalf(30)))
print('N : ', sp.N(sp.sin(x), 30))

x = sp.symbols('x')
sin_func = sp.sin('0.4')

# Convert the sine function to mpmath format
sin_func_mpmath = sp.lambdify(x, sin_func, modules='mpmath')

# Evaluate the sine function at a specific value of x
result = sin_func_mpmath(mp.mpf('0.4'))
print('mo: ', result)

with output:

mp: 0.389418342308650491666311756796

sp: 0.389418342308650522465285348517

N : 0.389418342308650522465285348517

mo: 0.389418342308650522465285348517

Solution

  • In SymPy when you pass a float to a function it evaluates automatically using the float precision. Passing a string is no different because if the string has a decimal point then it is converted to a float:

    In [13]: sympy.sin(0.4)
Out[13]: 0.389418342308651

In [14]: sympy.sin('0.4')
Out[14]: 0.389418342308651

    Since sin('0.4') has already been evaluated with 15 decimal digit precision, calling evalf(30) cannot increase the accuracy of the expected result. You get more digits but they are only more digits representing the exact value of the 15 digit approximation. With enough digits you will see the exact value of that approximation:

    In [15]: sympy.sin('0.4').evalf(1000)
Out[15]: 
0.38941834230865052246528534851677250117063522338867187500000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
000000000000000000000000000000000000000000000000000000000000000000000000000

    With SymPy if you want exact values or you want to use more than the default 15 decimal digits (53 binary bits) of precision then you should create exact rational numbers rather than floats:

    In [16]: sympy.sympify('0.4')
Out[16]: 0.400000000000000

In [17]: sympy.Rational('0.4')
Out[17]: 2/5

In [18]: sympy.sin(sympy.Rational('0.4'))
Out[18]: sin(2/5)

In [19]: sympy.sin(sympy.Rational('0.4')).evalf(30)
Out[19]: 0.389418342308650491666311756796

    An exact expression like sin(2/5) above is needed for evalf to be able to compute an answer that is accurate to the specified number of digits. This last result agrees with mpmath for 30 digits and is in fact computed by mpmath which SymPy uses for evalf internally.