Hello, i'm currently trying to write some code to find the integral of exp(-x*x).
I have been given code for one trapezium rule (trap0) which is imperfect, and a second set of code for a new trapezium rule (trap1), which uses the old code and some new conditions to 'perfect' it. Here is my attempt.
import numpy import matplotlib import pylab as P import math from math import exp def f(x): return exp(-x*x) def trap0 (f,a,b,n): """Basic trapezium rule. Integrate f(x) over the 30 4 Numerical Integration interval from a to b using n strips.""" h= float (b-a)/n s =0.5*( f(a)+f(b)) for i in range (1,n): s=s+f(a+i*h) return s*h def trap1 (f,a,b,delta , maxtraps=512): """ Improved trapezium rule. Integrate f(x) over interval from a to b, trying to get relative accuracy delta. Optional last argument is maximum allowed number of trapezia.""" n=8 inew= trap0(f,100,-100,n) iold=- inew # make sure we don’t terminate immediately while ( fabs(inew - iold)>delta * fabs( inew )): iold= inew n=2*n if n> maxtraps: print " Cannot reach requested accuracy with"
The integral rapidly becomes negligible if x is large so i understand you would simply use the limits of -5 to 5 to calculate this integral. However, i'm still unsure of how to do this.
I am using python 2.7
Thanks for any help