Python实现Simpson积分算法（小白也能学会的数值积分方法）

Simpson公式的数学原理

对于区间 \([a, b]\) 上的函数 \(f(x)\)，Simpson公式为：

\[ \int_a^b f(x) \, dx \approx \frac{b - a}{6} \left[ f(a) + 4f\left(\frac{a+b}{2}\right) + f(b) \right] \]

为了提高精度，通常将区间 \([a, b]\) 分成偶数个子区间（比如 \(n\) 段，\(n\) 为偶数），然后在每两个相邻子区间上应用上述公式，这就是复合Simpson公式。

用Python实现Simpson积分

下面我们用Python编写一个通用的Simpson积分函数。这个实现适合初学者理解，也便于在实际项目中使用。

import numpy as npdef simpson_integration(func, a, b, n):    """    使用复合Simpson法则计算定积分        参数:        func: 被积函数（Python函数）        a: 积分下限        b: 积分上限        n: 子区间数量（必须为偶数）            返回:        积分近似值    """    if n % 2 != 0:        raise ValueError("n 必须是偶数！")        h = (b - a) / n  # 步长    x = np.linspace(a, b, n + 1)  # 生成n+1个点    y = func(x)        # Simpson公式：首尾项 + 4*奇数项 + 2*偶数项    result = y[0] + y[-1] + 4 * np.sum(y[1:-1:2]) + 2 * np.sum(y[2:-2:2])    result *= h / 3        return result# 示例：计算 sin(x) 从 0 到 π 的积分if __name__ == "__main__":    import math        def f(x):        return np.sin(x)        a, b = 0, math.pi    n = 100  # 偶数        integral = simpson_integration(f, a, b, n)    print(f"Simpson积分结果: {integral:.6f}")    print(f"理论值 (2): {2:.6f}")    print(f"误差: {abs(integral - 2):.8f}")

代码详解

• 输入检查：确保子区间数 n 是偶数，因为Simpson法要求每两个小区间组成一组。
• 步长计算：h = (b - a) / n 是每个小区间的宽度。
• 节点生成：使用 np.linspace 生成等间距的采样点。
• 权重分配：Simpson法对端点、奇数索引点和偶数索引点分别赋予1、4、2的权重。
• 最终求和：按公式加权求和后乘以 h/3 得到积分近似值。

为什么选择Simpson积分法？

相比梯形法（误差阶为 \(O(h^2)\)），Simpson法的误差阶为 \(O(h^4)\)，这意味着当步长减半时，误差大约减少16倍！因此，在Python科学计算中，Simpson法是平衡精度与效率的优秀选择。

进阶建议

在实际项目中，你可以使用SciPy库中的 scipy.integrate.simpson 函数，它已经高度优化。但理解手动实现有助于你掌握数值分析教程中的核心思想。

现在你已经掌握了Simpson积分法的基本原理和Python实现！快去试试计算其他函数的积分吧。如果你觉得这篇文章对你有帮助，欢迎分享给更多学习Python数值积分的朋友！