Julia实现2D信号高斯卷积性能低于MATLAB？技术问题问询

• 完全浪费了循环矩阵的结构优势，计算量爆炸
  你当前的实现本质是把3801×3801的循环矩阵当作普通稠密矩阵做乘法，单次行向量乘矩阵的复杂度是O(nt²)，全量计算的复杂度是O(nr×nt²)≈671×1444万≈969亿次浮点运算，远高于MATLAB版本实际的有效计算量（MATLAB虽然写了全长度kernal的点乘求和，但高斯核大部分区域数值趋近于0，实际等效计算量远小于这个量级）。
• 存储顺序错误导致缓存命中率极低
  Julia和MATLAB都是列优先存储，你把时间序列放在行维度，每次取sC[i,:]都会产生非连续内存访问，缓存命中率极低，还会产生大量临时切片拷贝，你给出的@time结果里193.7亿次分配、288GiB内存占用就是这个原因导致的。
• mapslices本身性能极差
  Julia的mapslices是为了易用性设计的高层函数，内部有大量动态派发和临时对象生成，性能远低于手写循环或者向量化运算。
• 矩阵乘法方向错误
  行向量×矩阵的计算效率远低于矩阵×列向量，列优先存储下后者的内存访问是连续的，性能差距可达数倍到数十倍。

方案1：直接调用成熟的一维卷积实现（最推荐） #

不要自己造轮子实现卷积，直接用DSP.jl的卷积函数，或者ImageFiltering.jl的imfilter函数，天然支持自定义填充方式，且内部已经做了极致的性能优化：

using DSP, ImageFiltering

# 构造高斯核，只保留3σ有效长度
sigma = 63.7
kern = Kernel.gaussian((sigma,), (round(Int, 3*sigma)*2+1,))
# 自定义填充逻辑：左补0，右补边缘值，符合你的需求
padded_sC = padarray(sC, Pad(:replicate, (0, 900), (0,0)))
padded_sC = padarray(padded_sC, Pad(:constant, (0,0), (0,900)))
conv = imfilter(padded_sC, kern, NoPad())

这个实现的时间复杂度只有O(nr×nt log nt)，预计耗时会降到1秒以内。

方案2：优化循环矩阵实现 #

如果你一定要用循环矩阵，首先确保你用的Circulant类型来自ToeplitzMatrices.jl，这个包已经重载了高效的矩阵乘法，利用FFT实现O(nt log nt)的乘法复杂度，同时调整运算顺序避免逐行计算：

using ToeplitzMatrices, LinearAlgebra

K = Circulant(G.Fx)
# 转置信号为(nt, nr)，一次完成全量乘法，不需要逐行/逐列遍历
conv = transpose(K * transpose(sC))

这个实现也可以把耗时降到几秒级别，且不会有巨量内存分配。

方案3：手写循环实现（和MATLAB逻辑对齐） #

Julia的原生循环性能远高于MATLAB，你可以直接对齐MATLAB的逻辑写循环，加上@inbounds和@views宏关闭边界检查、用视图代替拷贝，性能至少是MATLAB版本的3倍以上：

const dt = 0.5
const b = 1/(2*sigma^2)
const preb = 1/(sqrt(2π)*sigma)
const half_kern = round(Int, 3*sigma)

WWc = similar(WWc0)
@inbounds @views for it = 1:nt
    # 只计算高斯核的有效范围，不需要全长度计算
    start_t = max(1, it - half_kern)
    end_t = min(nt, it + half_kern)
    kern = preb .* exp.(-b .* (tc[start_t:end_t] .- tc[it]).^2)
    N = sum(kern) * dt
    kern ./= N
    for iq = 1:Nq
        WWc[iq,it] = dt * sum(WWc0[iq, start_t:end_t] .* kern)
    end
end

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