Langue de Julia. Comment battre les opérations vectorisées?

Question

On dit que les boucles for-Julia sont aussi rapides que les opérations vectorisées et même plus rapides (si elles sont utilisées correctement). J'ai deux morceaux de code. L'idée est de trouver une statistique d'échantillon pour une séquence 0-1 donnée, qui est x (dans ces deux exemples j'essaie de trouver une somme, mais il y a des exemples plus compliqués, j'essaie juste de comprendre une signification générale des pièges de performance dans mon code). Le premier ressemble à:

S = 2 * sum(x) - n 
s_obs_new = abs(S)/sqrt(2 * n) 
pval = erfc(s_obs_new) 

et le second est quelque chose « naïf » et classique:

S = 0 
for i in eachindex(x) 
    S += x[i] 
end 
S = 2 * S - n 
s_obs_new = abs(S)/sqrt(2 * n) 
pval = erfc(s_obs_new) 

En utilisant @benchmark j'ai trouvé que la durée du premier exemple est d'environ 11,8 ms, et pour la seconde est de 38 ms. Cet exemple est très important pour moi, car il y a beaucoup d'autres endroits, où la vectorisation n'est pas possible, donc je veux faire des calculs de manière "dévoralisée" aussi vite qu'en vectorisé.

Y a-t-il des idées pour lesquelles le code devectorisé est probablement 4x fois plus lent que vectorisé? La stabilité de type est OK, il n'y a pas d'allocation grande mémoire inutiles et etc.

Le code pour la première fonction est:

function frequency_monobit_test1(x :: Array{Int8, 1}, n = 0) 
# count 1 and 0 in sequence, we want the number 
# of 1's and 0's to be approximately the same 
# reccomendation n >= 100 
# decision Rule(at 1% level): if pval < 0.01 -> non-random 
if (n == 0) 
    n = length(x) 
end 
S = 2 * sum(x) - n 
s_obs_new = abs(S)/sqrt(2 * n) 
pval = erfc(s_obs_new) 
return pval 

Le second est:

function frequency_monobit_test2(x :: Array{Int8, 1}, n = 0) 
# count 1 and 0 in sequence, we want the number 
# of 1's and 0's to be approximately the same 
# reccomendation n >= 100 
# decision Rule(at 1% level): if pval < 0.01 -> non-random 
if (n == 0) 
    n = length(x) 
end 
S = 0 
@inbounds for i in eachindex(x) 
    S += x[i] 
end 
S = 2 * S - n 
s_obs_new = abs(S)/sqrt(2 * n) 
pval = erfc(s_obs_new) 
return pval 

Answer 1

Ceci est un curieux Cas. Il semble y avoir un problème de performances lors de l'accumulation de Int8 dans une variable Int64.

Essayons ces fonctions:

using SpecialFunctions, BenchmarkTools 

function frequency_monobit_test1(x, n=length(x)) 
    S = sum(x) 
    return erfc(abs(2S - n)/sqrt(2n)) 
end 

function frequency_monobit_test3(typ::Type{<:Integer}, x, n=length(x)) 
    S = zero(typ) 
    for i in eachindex(x) 
     @inbounds S += x[i] 
    end 
    return erfc(abs(2S - n)/sqrt(2n)) 
end 

Initialiser certains vecteurs

N = 2^25; 
x64 = rand(0:1, N); 
x8 = rand(Int8[0, 1], N); 
xB = rand(Bool, N); 
xb = bitrand(N); 

Analyse comparative:

Pour Int64 entrée:

julia> @btime frequency_monobit_test1($x64) 
    17.540 ms (0 allocations: 0 bytes) 
0.10302739916042186 

julia> @btime frequency_monobit_test3(Int64, $x64) 
    17.796 ms (0 allocations: 0 bytes) 
0.10302739916042186 

julia> @btime frequency_monobit_test3(Int32, $x64) 
    892.715 ms (67106751 allocations: 1023.97 MiB) 
0.10302739916042186 

Nous voyons que sum et une boucle explicite sont également rapides, et que l'initialisation avec un Int32 est une mauvaise idée.

Pour Int32 entrée:

julia> @btime frequency_monobit_test1($x32) 
    9.137 ms (0 allocations: 0 bytes) 
0.2386386867682374 

julia> @btime frequency_monobit_test3(Int64, $x32) 
    8.839 ms (0 allocations: 0 bytes) 
0.2386386867682374 

julia> @btime frequency_monobit_test3(Int32, $x32) 
    7.274 ms (0 allocations: 0 bytes) 
0.2386386867682374 

sum et la boucle sont similaires à la vitesse. Accumuler dans Int32 permet d'économiser un peu de temps.

entrée Int8:

julia> @btime frequency_monobit_test1($x8) 
    5.681 ms (0 allocations: 0 bytes) 
0.16482999123032094 

julia> @btime frequency_monobit_test3(Int64, $x8) 
    19.517 ms (0 allocations: 0 bytes) 
0.16482999123032094 

julia> @btime frequency_monobit_test3(Int32, $x8) 
    4.815 ms (0 allocations: 0 bytes) 
0.16482999123032094 

boucle explicite si un peu plus rapide lors de l'accumulation dans un Int32, mais vache sacrée! que se passe-t-il avec Int64? C'est super lent!

Qu'en est-il de Bool?

julia> @btime frequency_monobit_test1($xB) 
    9.627 ms (0 allocations: 0 bytes) 
0.7728544347518309 

julia> @btime frequency_monobit_test3(Int64, $xB) 
    9.629 ms (0 allocations: 0 bytes) 
0.7728544347518309 

julia> @btime frequency_monobit_test3(Int32, $xB) 
    4.815 ms (0 allocations: 0 bytes) 
0.7728544347518309 

boucle et sum ont même vitesse mais, dans l'accumulation Int32 permet d'économiser la moitié du temps.

Maintenant, nous allons essayer une BitArray:

julia> @btime frequency_monobit_test1($xb) 
    259.044 μs (0 allocations: 0 bytes) 
0.7002576522570715 

julia> @btime frequency_monobit_test3(Int64, $xb) 
    19.423 ms (0 allocations: 0 bytes) 
0.7002576522570715 

julia> @btime frequency_monobit_test3(Int32, $xb) 
    19.430 ms (0 allocations: 0 bytes) 
0.7002576522570715 

Ainsi, sum sur un BitArray est fou rapide, parce que vous pouvez effectuer des additions CHUNKED, mais extraire des éléments individuels dans une boucle a une surcharge.

Conclusions:

• Vous pouvez obtenir de meilleures performances similaires ou avec une boucle qu'avec sum, sauf pour BitArray s ce qui est un cas très particulier.
• Si vous connaissez la longueur de votre tableau, et sachez que Int32 est assez pour contenir votre somme, alors c'est un gain de temps.
• Il se passe quelque chose de très étrange lors de l'accumulation de Int8 dans un Int64. Je ne sais pas pourquoi la performance est si mauvaise.
• Si vous êtes uniquement intéressé par les zéros et les uns, utilisez un tableau de Bool s, pas de Int8 s, et peut-être accumuler dans un Int32.
• BitArray s peut être très rapide dans certains cas.
• sum sur un vecteur de Int8 s est étrangement rapide, deux fois plus vite que sum(::Vector{Bool}).