在矩阵的列上有记忆效率的分拣精子

Question

我有一个大型矩阵,我想应用于sortperm该矩阵的每一列.天真的事情是

order = sortperm(X[:,j])

制作副本.这似乎是一种耻辱,所以我想我会尝试SubArray:

order = sortperm(sub(X,1:n,j))

但那甚至更慢.为了笑,我试过了

order = sortperm(1:n,by=i->X[i,j])

但当然那太可怕了.最快的方法是什么？

这是一些基准代码:

getperm1(X,n,j) = sortperm(X[:,j])
getperm2(X,n,j) = sortperm(sub(X,1:n,j))
getperm3(X,n) = mapslices(sortperm, X, 1)
n = 1000000
X = rand(n, 10)
for f in [getperm1, getperm2]
    println(f)
    for it in 1:5
        gc()
        @time f(X,n,5)
    end
end
for f in [getperm3]
    println(f)
    for it in 1:5
        gc()
        @time getperm3(X,n)
    end
end

结果:

getperm1
elapsed time: 0.258576164 seconds (23247944 bytes allocated)
elapsed time: 0.141448346 seconds (16000208 bytes allocated)
elapsed time: 0.137306078 seconds (16000208 bytes allocated)
elapsed time: 0.137385171 seconds (16000208 bytes allocated)
elapsed time: 0.139137529 seconds (16000208 bytes allocated)
getperm2
elapsed time: 0.433251141 seconds (11832620 bytes allocated)
elapsed time: 0.33970986 seconds (8000624 bytes allocated)
elapsed time: 0.339840795 seconds (8000624 bytes allocated)
elapsed time: 0.342436716 seconds (8000624 bytes allocated)
elapsed time: 0.342867431 seconds (8000624 bytes allocated)
getperm3
elapsed time: 1.766020534 seconds (257397404 bytes allocated, 1.55% gc time)
elapsed time: 1.43763525 seconds (240007488 bytes allocated, 1.85% gc time)
elapsed time: 1.41373546 seconds (240007488 bytes allocated, 1.82% gc time)
elapsed time: 1.42215519 seconds (240007488 bytes allocated, 1.83% gc time)
elapsed time: 1.419174037 seconds (240007488 bytes allocated, 1.83% gc time)

凡mapslices版本是10倍的getperm1版本,如你所期望.

值得指出的是,至少在我的机器上,copy + sortperm选项并不比只有相同长度的矢量上的sortperm慢得多,但是不需要内存分配,因此避免使用它会很好.

Answer 1

在一些非常特殊的情况下（例如连续查看Array），您可以使用指针魔法来击败 SubArray 的性能：

\n\n

function colview(X::Matrix,j::Int)\n    n = size(X,1)\n    offset = 1+n*(j-1) # The linear start position\n    checkbounds(X, offset+n-1)\n    pointer_to_array(pointer(X, offset), (n,))\nend\n\ngetperm4(X,n,j) = sortperm(colview(X,j))\n

\n\n

该函数colview将返回一个成熟的，Array与原始X. 请注意，这是一个糟糕的想法，因为返回的数组引用了 Julia 仅通过 跟踪的数据X。这意味着如果X在列“视图”之前超出范围，数据访问将因段错误而崩溃。

\n\n

结果：

\n\n

getperm1\nelapsed time: 0.317923176 seconds (15 MB allocated)\nelapsed time: 0.252215996 seconds (15 MB allocated)\nelapsed time: 0.215124686 seconds (15 MB allocated)\nelapsed time: 0.210062109 seconds (15 MB allocated)\nelapsed time: 0.213339974 seconds (15 MB allocated)\ngetperm2\nelapsed time: 0.509172302 seconds (7 MB allocated)\nelapsed time: 0.509961218 seconds (7 MB allocated)\nelapsed time: 0.506399583 seconds (7 MB allocated)\nelapsed time: 0.512562736 seconds (7 MB allocated)\nelapsed time: 0.506199265 seconds (7 MB allocated)\ngetperm4\nelapsed time: 0.225968056 seconds (7 MB allocated)\nelapsed time: 0.220587707 seconds (7 MB allocated)\nelapsed time: 0.219854355 seconds (7 MB allocated)\nelapsed time: 0.226289377 seconds (7 MB allocated)\nelapsed time: 0.220391515 seconds (7 MB allocated)\n

\n\n

我没有研究为什么 SubArray 的性能更差，但它可能只是由于每次内存访问时额外的指针取消引用。值得注意的是，分配实际上花费的时间很少 - getperm1 的时间变化更大，但它仍然偶尔会击败 getperm4！我认为这是由于Array共享数据的内部实现中的一些额外的指针数学造成的。还有一些疯狂的缓存行为\xe2\x80\xa6 getperm1 在重复运行时变得明显更快。

\n