Xuan Shi
Xuan Shi

Reputation: 11

Why sympy gives different/wrong answer when integrate fraction of power?

I am running a definite integral but the sympy gives me a different answer. I don't know what went wrong.

The equation integrated by Wolfram, the result is as expected 1.097

import sympy as sp
y=sp.Symbol('y')
fy=(y-2.4)**2*sp.sqrt(y)*0.1875
print(sp.integrate(fy,(y,0,4)))

This is giving me 0.485255751108899

Upvotes: 1

Views: 212

Answers (1)

Oscar Benjamin
Oscar Benjamin

Reputation: 14530

This looks like a bug to do with floats that can be problematic in integration:

In [15]: from sympy import *

In [16]: y = symbols('y')

In [17]: fy=(y-2.4)**2*sp.sqrt(y)*0.1875

In [18]: fy
Out[18]: 
                                 2
1.08⋅√y⋅(0.416666666666667⋅y - 1) 

In [19]: integ = Integral(fy, (y, 0, 4))

In [20]: integ
Out[20]: 
4                                      
⌠                                      
⎮                                  2   
⎮ 1.08⋅√y⋅(0.416666666666667⋅y - 1)  dy
⌡                                      
0   

In [21]: integ.doit() # symbolic integration with floats
Out[21]: 0.485255751108899

In [22]: integ.evalf() # numeric integration
Out[22]: 1.09714285714286

In [23]: integ_rational = Integral(nsimplify(fy), (y, 0, 4))

In [24]: integ_rational # same integral with Rational instead of Float
Out[24]: 
4                    
⌠                    
⎮                2   
⎮       ⎛5⋅y    ⎞    
⎮ 27⋅√y⋅⎜─── - 1⎟    
⎮       ⎝ 12    ⎠    
⎮ ──────────────── dy
⎮        25          
⌡                    
0                    

In [25]: integ_rational.doit() # symbolic integration with rationals
Out[25]: 
192
───
175

In [26]: _.n() # numerically evaluate the last result
Out[26]: 1.09714285714286

Firstly don't use floats with sympy unless you have a good reason e.g. use Rational('0.1') or S(1)/10 etc rather than 0.1. The nsimplify function can be used to convert float to Rational (it tries to guess what rational you really meant).

Secondly this is a bug so it should be reported to github rather than SO: https://github.com/sympy/sympy/issues

Upvotes: 1

Related Questions