lec11_integrate - errorsimp=(abs(Isimp-Itrue)/Itrue*100.)...

Info iconThis preview shows page 1. Sign up to view the full content.

View Full Document Right Arrow Icon
from pylab import * def trap(f,xmin,xmax,numPanels): Isum=0. h=(xmax-xmin)/float(numPanels) for i in range(numPanels): x=xmin+i*h Isum+=(f(x+h)+f(x))*h/2. return Isum def simp13(f,xmin,xmax,numPanels): Isum=0. h=(xmax-xmin)/float(numPanels) for i in range(0,numPanels,2): x=xmin+i*h Isum+=(f(x)+4*f(x+h)+f(x+2*h))*h/3 return Isum #def f(x): return exp(-x) #x**2-0.1*x**3 def f(x): return sin(x) numPanels=10 Itrap=trap(f,0,np.pi,numPanels) Isimp=simp13(f,0,np.pi,numPanels) #Itrue=-(exp(-1)-exp(0.)) Itrue=2.0 errortrap=(abs(Itrap-Itrue)/Itrue*100.)
Background image of page 1
This is the end of the preview. Sign up to access the rest of the document.

Unformatted text preview: errorsimp=(abs(Isimp-Itrue)/Itrue*100.) print "Integral by trapeziod with %i panels = %3.7f" % (numPanels,Itrap) print "Error for Trap is %3.4f percent" % errortrap print "Integral by Simpson 1/3 with %i panels = %3.7f" % (numPanels,Isimp) print "Error for Simp is %3.4f percent" % errorsimp print "True value = ", Itrue x=linspace(0,1,50) fill_between(x,exp(-x),y2=0,color='yellow') plot(x,exp(-x),'b-',lw=2) xlim(-0.2,1.2) xlabel('$x$',fontsize=18) ylabel('$f(x)$',fontsize=18)...
View Full Document

This note was uploaded on 10/05/2010 for the course PHYS phy503 taught by Professor Gladden during the Spring '09 term at Ole Miss.

Ask a homework question - tutors are online