Numerical Integration




#Trapezoidal

def f(x):
    return 1/(1+x)
a,b,n=input('Input lower limit(a), upper limit(b) & no. of steps(n): ')
h=(b-a)/float(n)
x=[a+h*i for i in range(n+1)]
y=[f(x[i]) for i in range(n+1)]
trapint=(h/2.0)*(y[0]+ y[-1]+ 2*sum(y[1:-1:1]))
print 'Integral using trapezoidal method= ',trapint

Output:
Input lower limit(a), upper limit(b) & no. of steps(n): 0,1,4
Integral using trapezoidal method=  0.697023809524



#Simpsons:
def f(x):
    return 1/(1+x)
a,b,n=input('Input lower limit(a), upper limit(b) & no. of steps(even no.)(n): ')
h=(b-a)/float(n)
x=[a+h*i for i in range(n+1)]
y=[f(x[i]) for i in range(n+1)]
simpint=(h/3.0)*(y[0]+y[-1]+4*sum(y[1:-1:2])+2*sum(y[2:-1:2]))
print 'Integral using Simpson\'s 1/3 rule= ', simpint

Output:
Input lower limit(a), upper limit(b) & no. of steps(even no.)(n): 0,1,4
Integral using Simpson's 1/3 rule=  0.693253968254


Comments