数值积分 python代码实现
老规矩，数学原理什么的就不写了。
直接贴代码和实例演示，以下代码基于python和numpy。
在这里，我将用代码实现复化梯形算法、复化 Simpson 算法、Romberg 积分算法和三点 Gauss-Legendre求积算法。
往期博客：
线性方程组的迭代法 python代码实现
函数插值法之牛顿插值法 python代码实现
数值积分复化梯形算法定义函数参数说明实例运行复化 Simpson 算法定义函数参数说明实例运行Romberg 积分算法定义函数参数说明实例运行三点 Gauss-Legendre求积算法定义函数参数说明实例运行总结 复化梯形算法复化梯形公式如下图所示：
首先
import numpy as np
定义函数
以下便是我定义的函数：
def tx_fh(x0,f,n):
a=x0[0]
b=x0[1]
x=np.linspace(a,b,n+1)
y=2*np.sum(f(x))-f(a)-f(b)
tn=((b-a)/n)*y/2
return tn
参数说明
“x0”指的是自变量的定义域。
“f”指的是函数（以后再实现不是函数，而是确定值的）。
“n”指的是将定义域分为n等份。
实例运行x=np.array([0,1])
f=lambda x:4/(1+x**2)
tx_fh(x,f,2**3)
得出以下结果：
3.1389884944910893
复化 Simpson 算法复化 Simpson 公式如下图所示：
以下便是我定义的函数：
def simpson_fh(x0,f,n):
a=x0[0]
b=x0[1]
x=np.linspace(a,b,2*n+1)
x1=x[1::2].copy()
x2=x[0::2].copy()
y=4*np.sum(f(x1))+2*np.sum(f(x2))-f(a)-f(b)
sn=((b-a)/n)*y/6
return sn
参数说明
“x0”指的是自变量的定义域。
“f”指的是函数（以后再实现不是函数，而是确定值的）。
“n”指的是将定义域分为2n等份。
实例运行x=np.array([0,1])
f=lambda x:4/(1+x**2)
simpson_fh(x,f,2**2)
得出以下结果：
3.141592502458707
综上可知，在相同的等份数下，相比于复化梯形，复化 Simpson更为精确。
Romberg 积分算法Romberg 积分算法的计算过程如下图所示：
以下便是我定义的函数：
def romberg(x0,f,n):
k=0
xlb=np.zeros((4,n+3))
for i in range(n+3):
xlb[0,i]=tx_fh(x0,f,2**i)
for i in range(n+2):
xlb[1,i+1]=4*xlb[0,i+1]/3-xlb[0,i]/3
for i in range(n+1):
xlb[2,i+2]=16*xlb[1,i+2]/15-xlb[1,i+1]/15
while k<n:
xlb[3,k+3]=64*xlb[2,k+3]/63-xlb[2,k+2]/63
k+=1
k=np.arange(n+3)
xl=np.vstack([k,xlb])
print(xl.T)
return xlb[3][3:]
参数说明
“x0”指的是自变量的定义域。
“f”指的是函数（以后再实现不是函数，而是确定值的）。
“n”指的是得到n个Romberg 积分值，也就是得到R1~Rn。
实例运行x=np.array([0.0000001,1])
f=lambda x:np.sin(x)/x
romberg(x,f,1)
得出以下结果（下列数组对照上一张图，可知每一个值所对应的数学意义）：
[[0. 0.9207354 0. 0. 0. ]
[1. 0.93979319 0.94614578 0. 0. ]
[2. 0.94451342 0.94608683 0.9460829 0. ]
[3. 0.94569076 0.94608321 0.94608297 0.94608297]]
array([0.94608297])
上面这个“array([0.94608297])”就是我们要的结果。
三点 Gauss-Legendre求积算法三点 Gauss-Legendre求积公式如下图所示：
以下便是我定义的函数：
def gauss_l3(x0,f):
a=x0[0]
b=x0[1]
t1=-(b-a)*0.6**0.5/2+(b+a)/2
t2=(b+a)/2
t3=(b-a)*0.6**0.5/2+(b+a)/2
y=((5*f(t1)+8*f(t2)+5*f(t3))/9)*(b-a)/2
return y
参数说明
“x0”指的是自变量的定义域。
“f”指的是函数（以后再实现不是函数，而是确定值的）。
实例运行x=np.array([0.0000001,1])
f=lambda x:np.sin(x)/x
gauss_l3(x,f)
得出以下结果：
0.9460830340784266
可以看出，相比于Romberg 积分算法，三点 Gauss-Legendre求积算法不用迭代，所以能够更快的得到答案。
总结在相同的等份数下，相比于复化梯形，复化 Simpson更为精确。
相比于Romberg 积分算法，三点 Gauss-Legendre求积算法不用迭代，所以能够更快的得到答案。
作者：一只不爱晒太阳的猫
