## Integrating functions in python

| categories: | tags: | View Comments

Problem statement

find the integral of a function f(x) from a to b i.e.

$$\int_a^b f(x) dx$$

In python we use numerical quadrature to achieve this with the scipy.integrate.quad command.

as a specific example, lets integrate

$$y=x^2$$

from x=0 to x=1. You should be able to work out that the answer is 1/3.

def integrand(x):
return x**2

ans, err = quad(integrand, 0, 1)
print ans
0.333333333333

## 1 double integrals

Integrate $$f(x,y)=y sin(x)+x cos(y)$$ over

$$\pi <= x <= 2\pi$$

$$0 <= y <= \pi$$

i.e.

$$\int_{x=\pi}^{2\pi}\int_{y=0}^{\pi}y sin(x)+x cos(y)dydx$$

The syntax in dblquad is a bit more complicated than in Matlab. We have to provide callable functions for the range of the y-variable. Here they are constants, so we create lambda functions that return the constants. Also, note that the order of arguments in the integrand is different than in Matlab.

import numpy as np

def integrand(y, x):
'y must be the first argument, and x the second.'
return y * np.sin(x) + x * np.cos(y)

ans, err = dblquad(integrand, np.pi, 2*np.pi,
lambda x: 0,
lambda x: np.pi)
print ans
-9.86960440109

we use the tplquad command to integrate $$f(x,y,z)=y sin(x)+z cos(x)$$ over the region

$$0 <= x <= \pi$$

$$0 <= y <= 1$$

$$-1 <= z <= 1$$

import numpy as np

def integrand(z, y, x):
return y * np.sin(x) + z * np.cos(x)

0, np.pi,  # x limits
lambda x: 0,
lambda x: 1, # y limits
lambda x,y: -1,
lambda x,y: 1) # z limits

print ans
2.0

## 2 Summary

scipy.integrate offers the same basic functionality as Matlab does. The syntax differs significantly for these simple examples, but the use of functions for the limits enables freedom to integrate over non-constant limits.