数值积分|自适应辛普森积分公式

import math
###自适应辛普森公式求积分
### y = 1/( 1+x^2 )


def Func(x):
    return 1/( 1+pow(x,2) )
    

def AdaptiveSimpsonCtrl(Func, a, b, eps = 1e-6):
    kmax = 9000   #最大迭代步数
    h = b-a       # 积分区间
    n = 1
    t0 = 0.5*h*(Func(a) + Func(b))  
    
    for k in range(1,kmax+1):
        sumf = 0
        for i in range(1,n+1): sumf += Func(a+(i-0.5)*h)
        
        t = 0.5*(t0 + h*sumf)
        
        s = (4*t - t0)/3
        
        if (k > 1):
            if (math.fabs(s-s0) <= eps * math.fabs(s)): break
            if (math.fabs(s) <= eps and math.fabs(s) <= math.fabs(s-s0)): break
            
        h *= 0.5
        n *= 2
        s0 = s
        t0 = t
        
    if (k >= kmax): print("已经达到最大迭代步数!")
    
    return s
    

S = AdaptiveSimpsonCtrl(Func, 0, 1, eps = 1e-6)
print(S)