Mafeni Alpha
Mafeni Alpha

Reputation: 308

Double Trapezoidal Integral in numpy

I have a two-dimensional function $f(x,y)=\exp(y-x)$. I would like to compute the double integral $\int_{0}^{10}\int_{0}^{10}f(x,y) dx dy$ using NumPy trapz. After some reading, they say I should just repeat the trapz twice but it's not working. I have tried the following

import numpy as np

def distFunc(x,y):
    f = np.exp(-x+y)
    return f

# Values in x to evaluate the integral.
x = np.linspace(.1, 10, 100)
y = np.linspace(.1, 10, 100)

list1=distFunc(x,y)
int_exp2d = np.trapz(np.trapz(list1, y, axis=0), x, axis=0)

The code always gives the error

IndexError: list assignment index out of range

I don't know how to fix this so that the code can work. I thought the inner trapz was to integrate along y first then we end by the second along x. Thank you.

Upvotes: 3

Views: 5956

Answers (1)

Chris Mueller
Chris Mueller

Reputation: 6680

You need to convert x and y to 2D arrays which can be done conveniently in numpy with np.meshgrid. This way, when you call distfunc it will return a 2D array which can be integrated along one axis first and then the other. As your code stands right now, you are passing a 1D list to the first integral (which is fine) and then the second integral receives a scalar value.

import numpy as np

def distFunc(x,y):
    f = np.exp(-x+y)
    return f

# Values in x to evaluate the integral.
x = np.linspace(.1, 10, 100)
y = np.linspace(.1, 10, 100)
X, Y = np.meshgrid(x, y)

list1=distFunc(X, Y)
int_exp2d = np.trapz(np.trapz(list1, y, axis=0), x, axis=0)

Upvotes: 5

Related Questions