Scipy nquad gives a wrong output

user_4859331 0
Viewed 38

I'm trying to write a code to solve

integral from 1 to 2 integral from 0 to x, xy^2 dydx = 31/15

I verified the manual solution with Wolfram|Alpha as well. But the code gives the output as 3.1.

from scipy.integrate import nquad

def f(x, y):
    return x * y**2
def limits_y(x):
    return [0, x]
def limits_x():
    return [1, 2]

ans,err = nquad(f, [limits_y, limits_x])
print(round(ans,3))

Any idea why? Am I missing something?

scipy python
1 Answers

nquad expects the first argument of the integrand function to be the variable corresponding to the innermost integral. If you change your definition of f to

def f(y, x):
    return x * y**2

then nquad returns the expected result:

In [73]: def f(y, x):
    ...:     return x*y**2 
    ...:

In [74]: ans, err = nquad(f, [limits_y, limits_x])

In [75]: ans
Out[75]: 2.066666666666667
user_12173580
Related