我正在尝试在 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),所以显然我的脚本有问题。但我看不出它。
你落后了两倍。确实,梯形法则正如数学课上所教的那样,会使用像这样的增量
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
本文内容由网友自发贡献,版权归原作者所有,本站不承担相应法律责任。如您发现有涉嫌抄袭侵权的内容,请联系:hwhale#tublm.com(使用前将#替换为@)