memory-management - 迭代数组时内存使用量激增

Question

我正在使用 Julia 研究泊松方程的 Jacobi 求解器。求解器被迭代调用直到 err足够小(~1e-8)，我的 nx = ny = 80 的函数需要大约 25,000 次循环测试用例。分析显示大部分时间都花在了内部循环中(如预期的那样)，但内存分配似乎正在消失——@time 宏分配了 38 GB 以达到收敛，这似乎太多了，因为我不这样做不认为我正在为每个循环创建新数组。

function jacobi(P::Array{Float64,2}, maxiter::Int64)
    P_old = copy(P)
    for j = 2:ny-1
        # Main body loop
        for i = 2:nx-1 
             @inbounds P[i,j] = ((P_old[i+1,j] + P_old[i-1,j])*dx2 
                               + (P_old[i,j+1] + P_old[i,j-1])*dy2)/denom-Rmod[i,j]
        end
    end
    err = vecnorm(P::Array{Float64,2}-P_old::Array{Float64,2})/sqrt(nx+ny)
    return (P, err)
end

我已经为 1000 次循环计时了函数，从设置初始条件的函数包装器 ( methodwrap ) 调用:

function methodwrap(solver, maxiter::Int64) # (solver fn name, max # of iterations)
    P = copy(P0)
    iter = 1
    err  = 1.0
    maxerr = 1e-8
    prog = Progress(maxiter,.2, "Solving using $solver method", 10) # Show progress bar
    while (err > maxerr) && (iter < maxiter)
        P, err = solver(P, maxiter)
        next!(prog) # Iterates progress bar counter
        iter += 1
    end
    println()
    return (P, iter, err)
end

与我的愿望相反，看起来内存分配与循环次数成比例，所以我做错了。看起来好像每个 Jacobi pass 分配了大约 1.4 mb:

julia> @time methodwrap(jacobi,1000)
Solving using jacobi method 98%|##########|  ETA: 0:00:00
elapsed time: 4.001988593 seconds (1386549012 bytes allocated, 26.45% gc time)

我尝试将内部循环数组减少为向量子数组并使用@simd:

function jacobi2(P::Array{Float64,2}, maxiter::Int64)
    P_old = copy(P)::Array{Float64,2}
    for j = 2:ny-1
        # Main body loop
        Pojm = sub(P_old,:,j-1)
        Poj = sub(P_old,:,j)
        Pojp = sub(P_old,:,j+1)
        Pj = sub(P,:,j)
        Rmodj = sub(Rmod,:,j)
        @simd for i = 2:nx-1 
             @inbounds Pj[i] = ((Poj[i+1] + Poj[i-1])*dx2
                              + (Pojp[i] + Pojm[i])*dy2)/denom-Rmodj[i]
        end
    end
    err = vecnorm(P::Array{Float64,2}-P_old::Array{Float64,2})/sqrt(nx+ny)
    return (P, err)
end

但是，这似乎只会增加内存分配并降低速度，并且我收到了 @simd 警告:

julia> @time methodwrap(jacobi2,1000);
Warning: could not attach metadata for @simd loop.
Solving using jacobi2 method: 100%|##########|  ETA: 0:00:00
elapsed time: 4.947097666 seconds (1455818184 bytes allocated, 29.85% gc time)

这是我在 Julia 的第一个项目，所以我可能犯了一个非常明显的错误，但我还没有找到解决方案。我已经将全局变量定义为常量。我已经多次阅读了性能提示，我已经整理了文件，我已经使用 TypeCheck 来确保我的类型是一致的，并且一切在我看来都非常洁净。我究竟做错了什么？我已经发布了我的 full code如果您也想检查一下，请在 Gist 上查看。

Answer 1

事实证明，这个问题很微妙。我做了 3 处更改(见下文)。我确实使用了@IainDunning 建议的 --track-allocation=user ，它指向了有问题的行。这两个问题都来自使用全局变量。

这些变化之后

julia> @time methodwrap(jacobi,1000)
elapsed time: 0.481986712 seconds (116650236 bytes allocated)

更改 1 将 const 添加到 nx 和 ny

除了这 2 个变量之外，你到处都有 const，但是当离开非 const 和 global 时，会导致循环迭代器 i 不必要地分配。

nx=80 # Number of mesh points in the x-direction
ny=80 # Number of mesh points in the y-direction

改为

const nx=80 # Number of mesh points in the x-direction
const ny=80 # Number of mesh points in the y-direction

变化 2:避免 Array{Any,2} 类型的 Rmod

const Rmod = dx2*dy2*R/(2*(dx2+dy2))

改为

const Rmod = convert(Array{Float64,2},dx2*dy2*R/(2*(dx2+dy2)))