scipy

Numerical Integration

quad()

In [1]:
import numpy as np
from scipy.integrate import quad
In [2]:
I=quad(lambda x: np.exp(-x**2),1,2)
I
Out[2]:
(0.13525725794999466, 1.5016572202374808e-15)
In [3]:
I[0]
Out[3]:
0.13525725794999466
In [4]:
I=quad(lambda x: np.exp(-x),0,np.inf)
I[0]
Out[4]:
1.0000000000000002

dblquad()

In [5]:
from scipy.integrate import dblquad
I=dblquad(lambda x,y:x*y,0,0.5,lambda y:0,lambda y:1-2*y)
In [6]:
I[0]-1.0/96
Out[6]:
1.734723475976807e-18
In [7]:
I=dblquad(lambda x,y: np.exp(-x**2-y**2),-np.inf,np.inf,lambda y:-np.inf,lambda y:np.inf)
In [8]:
I[0]
Out[8]:
3.141592653589777

newton()

In [9]:
from scipy.optimize import newton
newton(lambda x: x**3-1,0.5)
Out[9]:
0.9999999999999973

Polynomials

poly1d()

In [10]:
p=np.poly1d([3,4,5])
In [11]:
print p
   2
3 x + 4 x + 5
In [12]:
print p*p
   4      3      2
9 x + 24 x + 46 x + 40 x + 25
In [13]:
print p+p
   2
6 x + 8 x + 10
In [14]:
print p.deriv()
 
6 x + 4
In [15]:
print p.integ()
   3     2
1 x + 2 x + 5 x
In [16]:
print p.integ(k=6)
   3     2
1 x + 2 x + 5 x + 6
In [17]:
print p
   2
3 x + 4 x + 5
In [18]:
p(2)
Out[18]:
25
In [19]:
p([3,4,5])
Out[19]:
array([ 44,  69, 100])

np.roots()

In [20]:
np.roots(p)
Out[20]:
array([-0.66666667+1.1055416j, -0.66666667-1.1055416j])