整合一个二维向量场阵列（反转 np.gradient）

Question

我有以下问题：我想整合一个二维数组，所以基本上是反转一个梯度算子。

假设我有一个非常简单的数组，如下所示：

shape = (60, 60)
sampling = 1
k_mesh = np.meshgrid(np.fft.fftfreq(shape[0], sampling), np.fft.fftfreq(shape[1], sampling))

然后我将我的向量场构造为一个复值数组（x 向量 = 实部，y 向量 = 虚部）：

k = k_mesh[0] + 1j * k_mesh[1]

所以例如真实的部分看起来像这样

现在我取梯度：

k_grad = np.gradient(k, sampling)

然后我使用傅立叶变换来反转它，使用以下函数：

def freq_array(shape, sampling):

    f_freq_1d_y = np.fft.fftfreq(shape[0], sampling[0])
    f_freq_1d_x = np.fft.fftfreq(shape[1], sampling[1])
    f_freq_mesh = np.meshgrid(f_freq_1d_x, f_freq_1d_y)
    f_freq = np.hypot(f_freq_mesh[0], f_freq_mesh[1])

    return f_freq


def int_2d_fourier(arr, sampling):
    freqs = freq_array(arr.shape, sampling)

    k_sq = np.where(freqs != 0, freqs**2, 0.0001)
    k = np.meshgrid(np.fft.fftfreq(arr.shape[0], sampling), np.fft.fftfreq(arr.shape[1], sampling))

    v_int_x = np.real(np.fft.ifft2((np.fft.fft2(arr[1]) * k[0]) / (2*np.pi * 1j * k_sq)))
    v_int_y = np.real(np.fft.ifft2((np.fft.fft2(arr[0]) * k[0]) / (2*np.pi * 1j * k_sq)))

    v_int_fs = v_int_x + v_int_y
    return v_int_fs


k_int = int_2d_fourier(k, sampling)

不幸的是，结果是不其中位置非常精确的k具有急剧变化，如可在下面的曲线图，这里显示的一水平线轮廓可以看到k和k_int。

任何想法如何提高准确性？有没有办法让它完全一样？

Answer 1

我实际上找到了解决方案。积分本身会产生非常准确的结果。然而，numpy 中的梯度函数计算二阶精确中心差，这意味着梯度本身已经是一个近似值。

当您将上述问题替换为解析公式（例如二维高斯）时，可以通过解析方式计算导数。当对这个分析导出的函数进行积分时，误差约为 10^-10（取决于高斯的宽度，这可能导致混叠效应）。

长话短说：上面提出的集成功能按预期工作！