Composite Simpson's Rule

def composite_simpsons(f, a, b, n):
    h = (b-a)/n
    s1 = 0
    s2 = 0
    for i in xrange(1, n/2):
        s1 += f(h*2*n)
    for i in xrange(1, n/2 + 1):
        s2 += f(h*(2*n - 1))
    return h/3 * (f(a) + 2*s1 + 4*s2 + f(b))