数值积分方法问题 python

numerical integration methods question python

我正在尝试整合以下公式:

下面是我尝试使用 f(a) = sin(a) 的自制方案执行此集成。

def func(x):
    return math.sin(x)

def integration (f, n, r, a, dtheta ):

    summation = 0
    theta = 0
    while theta <= 2*np.pi:

        f_arg = a + r*np.exp(1j*theta)
        second = np.exp(-1j*theta*n)

        summation += f(f_arg) * second * dtheta
        theta += dtheta

    return math.factorial(n)*summation / (2*np.pi*r**n)

integration(func, n=1, r=1, a=0, dtheta=2*np.pi/10000)

f(a) = sin(a) 的一阶导数 (n=1) 是 f'(a) = cos(a)。当在 a = 0 处计算时,这应该给出 cos(0) = 1,但是,它没有。我哪里错了?

看来你的问题是 math.sin 函数,它不支持复杂参数:

i = np.exp(.5j * np.pi)
math.sin(i), np.sin(i)
(6.123233995736766e-17, (9.44864380126377e-17+1.1752011936438014j))

它还会发出警告(但不是错误...):

ComplexWarning: Casting complex values to real discards the imaginary part

改用 np.sin 可以解决问题。

一般来说,如果更多地使用 numpy,实现可能更容易表达(也更容易调试),例如

def integration(func, a, n, r, n_steps):
    z = r * np.exp(2j * np.pi * np.arange(0, 1, 1. / n_steps))
    return math.factorial(n) * np.mean(func(a + z) / z**n)

np.allclose(1., integration(np.sin, a=0., n=1, r=1., n_steps=100))

True