Python实现Simpson 1/3二重积分计算结果错误排查

问题背景 #

使用Simpson 1/3规则实现二重积分计算，待积分函数为：

已知该积分理论正确结果为0.025，但编写的代码运行输出仅为0.008869，与理论值偏差较大。
实现参考的Simpson 1/3二重积分计算公式为：

原始实现代码 #

a= 4 
b= 4.4
c= 2
d= 2.6
N1= 2
N2=2
h=(b-a)/N1 
k=(d-c)/N2

def f(x,y):
    return (1/(x*y))

for i in range(N1+1):
     r = a+i*h
   for j in range(N2+1):
       t = c + j*k
       first=f(a,c)+f(a,d)
       if i%2==0:
          Xf1 =first+2*f(a,t)
       else:
          Xf1 =first+0
       if i%2!=0:
           Xf2 = Xf1+ 4*f(a,t)
       else:
          Xf2 = Xf1+ 0
    
       last=f(b,c)+f(b,d)
       if i%2==0:
          Xl1 =last +2*f(b,t)
       else:
         Xl1 =last+0
       if i%2!=0:
          Xl2 = Xl1+ 4*f(b,t)
       else:
          Xl2 = Xl1+ 0
    
       if i%2==0:
           even = 2* f(r,t)
       else:
          odd = 4 *f(r,t)
    
    

  print("Xf=",Xf2 )
  print("Xl=",Xl2 )
  print("even=",even )
  print("odd",odd)
  ans = Xl2+Xf2+even+odd
  ans = ((h*k)/9)*ans
  print(ans)

代码存在的核心逻辑错误 #

• 缩进语法错误：Python依赖缩进判断代码块归属，当前代码内层j循环、循环后的打印与计算语句缩进层级不匹配，会直接触发语法错误，或导致执行顺序完全偏离预期。
• 累加变量未初始化且循环内反复覆盖：Xf2、Xl2、even、odd四个存储加权和的变量没有在循环外初始化为0，每次内层j循环都会直接覆写上一轮的计算值，最终仅保留最后一次网格点的计算结果，完全没有实现全网格点的加权求和逻辑。
• 角点重复累加：四个角点f(a,c)、f(a,d)、f(b,c)、f(b,d)按照公式权重为1，仅需计算1次，但代码在每一次j循环中都重复将这几个值计入和中，导致边界值计算完全错误。
• 权重判断逻辑缺失：Simpson 1/3二重积分的点权重需要同时根据x方向索引i、y方向索引j的奇偶性共同判定，当前代码仅通过i的奇偶性赋值，完全没有判断j的奇偶性，不符合二维Simpson公式的权重规则。
• 边界与内部点计算逻辑错误：x=a、x=b两条y方向边界，以及中间x位置的所有行，都需要沿y方向逐点加权累加，当前代码没有实现沿j方向的逐点累加逻辑，仅在单个点上做赋值，完全偏离公式的求和要求。

修正提示：实现时需先将总加权和初始化为0，优先计算四个角点的贡献，再分别计算四条边界上非角点的加权贡献，最后计算区域内部所有网格点的加权贡献，所有点的权重由所在位置的i、j奇偶性共同决定，最终乘以系数h*k/9即可得到正确计算结果。

内容的提问来源于stack exchange，提问作者mugdho nahian