梯形法数值积分高斯函数模平方时出现y轴偏移的原因及解决方法咨询

你好呀，我看了你的问题和代码，很快就发现了问题核心——你把解析解和数值积分的积分范围搞反了！咱们一步步拆解清楚：

问题根源：积分下限不一致 #

1. 你的u_integral函数用梯形法计算的是从time_arr的起始点（也就是-20）到当前时间t的数值积分，数学上是：
   $$\int_{-20}^t |u(\tau)|^2 d\tau$$
2. 而你的gaussian_squared_integral函数实现的是从负无穷到t的解析积分，也就是：
   $$\int_{-\infty}^t |u(\tau)|^2 d\tau$$

这两个积分的下限差了一个区间，自然会产生固定偏移——这个偏移量就是$\int_{-\infty}^{-20} |u(\tau)|^2 d\tau$。因为高斯函数在t=-20处的值几乎为0，这个积分值非常接近0.5（毕竟整个高斯模平方从-∞到+∞的积分是1，左半部分从-∞到0的积分就是0.5），这正好对应你看到的y轴偏移。

解决方法：统一积分范围 #

有两种简单的修正思路，任选其一即可：

思路一：修正解析解的积分下限

把解析解改成和数值积分一样，从time_arr的起始点开始计算，这样两者的范围完全匹配：

from scipy.integrate import trapezoid
from scipy.interpolate import CubicSpline
import numpy as np
from scipy.special import erf
import matplotlib.pyplot as plt

def u_time(t_0):
    prefactor = 1 / (np.pi**0.25)
    def func(time):
        return prefactor * np.exp(-(time - t_0) ** 2 / 2)
    return func

def u_integral(t_0, times):
    u_values = u_time(t_0)(times)
    u_abs_2 = np.abs(u_values) ** 2
    # 改成float类型即可，不需要复数
    u_abs_2_int = np.zeros_like(times, dtype=np.float64)
    for k in range(1, len(times)):
        u_abs_2_int[k] = trapezoid(u_abs_2[:k], times[:k])
    return CubicSpline(times, u_abs_2_int)

# 修正解析解：传入数值积分的起始点，统一积分范围
def gaussian_squared_integral(t_0, start_time):
    return lambda time: 0.5 * (erf(time - t_0) - erf(start_time - t_0))
    
t_0 = 0
time_arr = np.linspace(-20, 100, 1000)
numeric = u_integral(t_0, time_arr)(time_arr)
# 传入数值积分的起始点-20
analytic = gaussian_squared_integral(t_0, time_arr[0])(time_arr)

plt.plot(time_arr, numeric, label="Numeric")
plt.plot(time_arr, analytic, label="Analytic")
plt.legend()
plt.show()

思路二：调整数值积分的起始点

如果想让数值积分尽量接近解析解的负无穷下限，可以把time_arr的起始点设得足够小（比如-100），此时$\int_{-\infty}^{-100} |u(\tau)|^2 d\tau$几乎为0，数值结果和解析解的偏移就可以忽略不计。

对你后续问题的建议 #

针对你提到的频率域反射脉冲积分的情况，只要记住一个原则：确保你对比的两个积分（数值/解析，或者不同函数的数值积分）的上下限完全一致，偏移问题就不会再出现啦。

备注：内容来源于stack exchange，提问作者Frostman