import numpy as np
from scipy import integrate
import matplotlib.pyplot as plt


def computePi(x):
    return 4/(1+x**2)

res = integrate.quad(computePi, 0.0, 1.0)
print("(res, error)=", res)
res = integrate.romberg(computePi, 0, 1)
print("(res, error)=", res)
res = integrate.quad(lambda x: x**5*np.exp(-x)*np.sin(x), 0, np.inf)
print("(res, error)=", res)

x = np.linspace(0, 2, num=2**4+1)
y = x**4

y_int = integrate.cumtrapz(y, x, initial=0)

plt.plot(x, y_int, 'ro', x, y[0] + 0.2 * x**5, 'b-')
plt.show()
