python - Python中的梯形规则

Question

我正在尝试在 Python 2.7.2 中实现梯形规则。我写了以下函数：

def trapezoidal(f, a, b, n):
    h = float(b - a) / n
    s = 0.0
    s += h * f(a)
    for i in range(1, n):
        s += 2.0 * h * f(a + i*h)
    s += h * f(b)
    return s

但是， f(lambda x:x**2, 5, 10, 100) 返回 583.333（它应该返回 291.667），所以很明显我的脚本有问题。不过我看不出来。

Answer 1

你差了两倍。事实上，数学课上教授的梯形规则会使用一个增量，比如

s += h * (f(a + i*h) + f(a + (i-1)*h))/2.0

(f(a + i*h) + f(a + (i-1)*h))/2.0是对网格上两个相邻点的函数高度进行平均。

由于每两个相邻的梯形都有一条公共边，因此上面的公式需要根据需要对函数进行两次评估。

更有效的实现（更接近您发布的内容）将结合来自相邻迭代的常用术语for-loop：

f(a + i*h)/2.0 + f(a + i*h)/2.0 =  f(a + i*h) 

到达：

def trapezoidal(f, a, b, n):
    h = float(b - a) / n
    s = 0.0
    s += f(a)/2.0
    for i in range(1, n):
        s += f(a + i*h)
    s += f(b)/2.0
    return s * h

print( trapezoidal(lambda x:x**2, 5, 10, 100))

产生

291.66875

Answer 2