Adaptive Simpson's Method - Sample Implementations - Python

Python

Here is an implementation of adaptive Simpson's method in Python. Note that this is explanatory code, without regard for efficiency. Every call to recursive_asr entails six function evaluations. For actual use, one will want to modify it so that the minimum of two function evaluations are performed.

def simpsons_rule(f,a,b): c = (a+b) / 2.0 h3 = abs(b-a) / 6.0 return h3*(f(a) + 4.0*f(c) + f(b)) def recursive_asr(f,a,b,eps,whole): "Recursive implementation of adaptive Simpson's rule." c = (a+b) / 2.0 left = simpsons_rule(f,a,c) right = simpsons_rule(f,c,b) if abs(left + right - whole) <= 15*eps: return left + right + (left + right - whole)/15.0 return recursive_asr(f,a,c,eps/2.0,left) + recursive_asr(f,c,b,eps/2.0,right) def adaptive_simpsons_rule(f,a,b,eps): "Calculate integral of f from a to b with max error of eps." return recursive_asr(f,a,b,eps,simpsons_rule(f,a,b)) from math import sin print adaptive_simpsons_rule(sin,0,1,.000000001)