二维阵列的 FFTW.jl：扩散仅发生在 1D

如何解决二维阵列的 FFTW.jl：扩散仅发生在 1D

从我读到的内容来看，当 A 是二维数组时，使用 FFTW.jl / AbstractFFTs.jl 的 fft(A) 应该在 2D 中执行 fft，而不是按列执行。知道为什么当（我认为）我将缩放的二阶空间导数添加到 u(t,x) 时，我只看到按列扩散，就像使用显式求解器计算时间一样？

谢谢！我对此很陌生。

code and heatmap screenshot

using Random
using FFTW
using Plots

gr()

N = (100,100)

# initialize with gaussian noise
u = randn(Float16,(N[1],N[2])).*0.4.+0.4
# include square of high concentration to observe diffusion clearly
u[40:50,40:50] .= 3

N = size(x)
L = 100
k1 = fftfreq(51)
k2 = fftfreq(51)
lap_mat = -(k1.^2 + k2.^2)

function lap_fft(x)
    lapF = rfft(x)
    lap = irfft(lap_mat.*lapF,100)
    return lap
end

# ode stepper or Implicit-Explicit solver
for i in 1:100000
    u+=lap_fft(u)*0.0001
end
# plot state
heatmap(u)

解决方法

仅仅因为您正在执行真正的 FFT，并不意味着您可以对结果进行真正的逆FFT。 rfft 来自 R -> C。但是你可以做的是：

function lap_fft(x)
    lapF = complex(zeros(100,100));        # only upper half filled
    lapF[1:51,1:100] = rfft(x) .* lap_mat; # R -> C
    return abs.(ifft(lapF));               # C -> R
end

实FFT到复频域（由于数据冗余只填充了上半部分），在频域乘以你的滤波器，逆FFT到复数图像域并获得幅度abs.()，实部real.()等
但老实说，为什么要使用真正的 fft 呢？

using Random
using FFTW
using Plots

gr()

N = (100,100)

# initialize with gaussian noise
u = randn(Float16,(N[1],N[2])).*0.4.+0.4;
# include square of high concentration to observe diffusion clearly
u[40:50,40:50] .= 3;

N = size(u);
L = 100;
k1 = fftfreq(100);
k2 = fftfreq(100);
tmp = -(k1.^2 + k2.^2);
lap_mat = sqrt.(tmp.*reshape(tmp,1,100));


function lap_fft(x)
    return abs.(ifftshift(ifft(fftshift(ifftshift(fft(fftshift(x))).*lap_mat))));
end

# ode stepper or Implicit-Explicit solver
for i in 1:100000
    u+=lap_fft(u)*0.001;
end
# plot state
heatmap(u)