Performance#
Julia avec la compilation Just-In-Time est un langage naturellement performant. Il n’est pas allergique aux boucles comme le sont les langages Python et R. Les opérations vectorisées fonctionnent également très bien à condition d’être attentifs aux allocations mémoire et aux vues explicites.
Allocations#
using Random, LinearAlgebra, BenchmarkTools
function test(A, B, C)
C = C - A * B
return C
end
A = rand(1024, 256); B = rand(256, 1024); C = rand(1024, 1024)
@btime test($A, $B, $C);
2.903 ms (4 allocations: 16.00 MiB)
Dans l’appel de la macro @benchmark on interpole les arguments avec le signe $ pour être sur que les fonctions
rand aient déjà été evaluées avant l’appel de la fonction test. La matrice C est modifiée dans la fonction suivante donc par convention on ajoute un ! au nom de la fonction. Par convention également, l’argument modifié se placera en premier. Comme dans la fonction push! par exemple.
function test!(C, A, B)
C .= C .- A * B
return C
end
@btime test!( $C, $A, $B);
2.119 ms (2 allocations: 8.00 MiB)
En effectuant une opération “en place”, on supprime une allocation mais celle pour effectuer l’opération A * B est toujours nécessaire. On peut supprimer cette allocation en utilisant la bibliothèque BLAS, cependant le code perd en lisibilité ce qu’il a gagné en performance.
function test_opt!(C, A, B)
BLAS.gemm!('N','N', -1., A, B, 1., C)
return C
end
@btime test_opt!($C, $A, $B);
954.786 μs (0 allocations: 0 bytes)
function test_opt_mul!(C, A, B)
mul!(C, A, B, -1, 1) # mul!(C, A, B, α, β) -> ABα + Cβ
return C
end
@btime test_opt_mul!($C, $A, $B);
955.358 μs (0 allocations: 0 bytes)
Alignement de la mémoire#
Les opérations le long des premiers indices seront plus rapides.
using FFTW
T = randn(1024, 1024)
@btime fft($T, 1);
10.439 ms (7 allocations: 32.00 MiB)
@btime fft($T, 2);
28.566 ms (7 allocations: 32.00 MiB)
Voici un autre exemple ou on calcule la dérivée de la quantité \(f\) suivant la coordonnée \(y\) en passant par l’espace de Fourier
using FFTW
xmin, xmax, nx = 0, 4π, 1024
ymin, ymax, ny = 0, 4π, 1024
x = LinRange(xmin, xmax, nx+1)[1:end-1]
y = LinRange(ymin, ymax, ny+1)[1:end-1]
ky = 2π ./ (ymax-ymin) .* [0:ny÷2-1;ny÷2-ny:-1]
exky = exp.( 1im .* ky' .* x)
function df_dy( f, exky )
ifft(exky .* fft(f, 2), 2)
end
f = sin.(x) .* cos.(y') # f is a 2d array created by broadcasting
@btime df_dy($f, $exky);
59.388 ms (14 allocations: 64.00 MiB)
En utilisant les “plans” de FFTW qui permettent de pré-allouer la mémoire nécessaire et le calcul “en place”. On peut améliorer les performances. On réutilise le même tableau pour la valeur de \(f\) et sa transformée de Fourier. On prend soin également de respecter l’alignement de la mémoire en transposant le tableau contenant \(f\) pour calculer la FFT. On utilise plus de mémoire, on fait plus de calcul en ajoutant les transpositions mais finalement le calcul va 3 fois plus vite car on évite les allocations et on limite les accès mémoire.
f = zeros(ComplexF64, (nx,ny))
fᵗ = zeros(ComplexF64, reverse(size(f)))
f̂ᵗ = zeros(ComplexF64, reverse(size(f)))
f .= sin.(x) .* cos.(y')
fft_plan = plan_fft(fᵗ, 1, flags=FFTW.PATIENT)
function df_dy!( f, fᵗ, f̂ᵗ, exky )
transpose!(fᵗ,f)
mul!(f̂ᵗ, fft_plan, fᵗ)
f̂ᵗ .= f̂ᵗ .* exky
ldiv!(fᵗ, fft_plan, f̂ᵗ)
transpose!(f, fᵗ)
end
@btime df_dy!($f, $fᵗ, $f̂ᵗ, $exky );
21.431 ms (2 allocations: 112 bytes)
Vues explicites#
@btime sum(T[:,1]) # Somme de la première colonne
1.879 μs (3 allocations: 8.16 KiB)
19.351104431816317
@btime sum(view(T,:,1))
242.526 ns (3 allocations: 80 bytes)
19.351104431816317
Eviter les calculs dans l’environnement global.#
v = rand(1000)
function somme()
acc = 0
for i in eachindex(v)
acc += v[i]
end
acc
end
@btime somme()
77.261 μs (3978 allocations: 77.77 KiB)
514.0748698860377
@code_lowered somme()
CodeInfo(
1 ─ acc = 0
│ %2 = Main.eachindex(Main.v)
│ @_2 = Base.iterate(%2)
│ %4 = @_2 === nothing
│ %5 = Base.not_int(%4)
└── goto #4 if not %5
2 ┄ %7 = @_2
│ i = Core.getfield(%7, 1)
│ %9 = Core.getfield(%7, 2)
│ %10 = acc
│ %11 = Base.getindex(Main.v, i)
│ acc = %10 + %11
│ @_2 = Base.iterate(%2, %9)
│ %14 = @_2 === nothing
│ %15 = Base.not_int(%14)
└── goto #4 if not %15
3 ─ goto #2
4 ┄ return acc
)
Il faut écrire des fonctions avec les variables utilisées en argument
function somme( x )
acc = 0
for i in eachindex(x)
acc += x[i]
end
acc
end
@btime somme( v )
3.461 μs (1 allocation: 16 bytes)
514.0748698860377
@code_lowered somme(v)
CodeInfo(
1 ─ acc = 0
│ %2 = Main.eachindex(x)
│ @_3 = Base.iterate(%2)
│ %4 = @_3 === nothing
│ %5 = Base.not_int(%4)
└── goto #4 if not %5
2 ┄ %7 = @_3
│ i = Core.getfield(%7, 1)
│ %9 = Core.getfield(%7, 2)
│ %10 = acc
│ %11 = Base.getindex(x, i)
│ acc = %10 + %11
│ @_3 = Base.iterate(%2, %9)
│ %14 = @_3 === nothing
│ %15 = Base.not_int(%14)
└── goto #4 if not %15
3 ─ goto #2
4 ┄ return acc
)
Pour comprendre pourquoi l’utilisation de variable global influence les performances, prenons un exemple simple d’une fonction additionnant deux nombres:
variable = 10
function addition_variable_globale(x)
x + variable
end
@btime addition_variable_globale(10)
16.569 ns (0 allocations: 0 bytes)
20
Comparons la performance avec cette fonction qui retourne la somme de ses deux arguments
function addition_deux_arguments(x, y)
x + y
end
@btime addition_deux_arguments(10, 10)
1.286 ns (0 allocations: 0 bytes)
20
On remarque que la deuxième fonction est 300 fois plus rapide que la première. Pour comprendre pourquoi elle est plus rapide, on peut regarder le code généré avant la compilation. On s’appercoit que le code est relativement simple avec une utilisation unique de l’instruction add.
@code_llvm addition_deux_arguments(10, 10)
; @ In[16]:1 within `addition_deux_arguments`
define i64 @julia_addition_deux_arguments_1863(i64 signext %0, i64 signext %1) #0 {
top:
; @ In[16]:2 within `addition_deux_arguments`
; ┌ @ int.jl:87 within `+`
%2 = add i64 %1, %0
; └
ret i64 %2
}
Si on regarde le code généré utilisant la variable globale, on comprend rapidement pourquoi c’est plus long. Pourquoi le code est-il si compliqué ? Ici le langage ne connait pas le type de variable, il doit donc prendre en compte le fait que ce type puisse être modifié à tout moment. Comme tous les cas sont envisagés, cela provoque un surcoût important.
@code_llvm addition_variable_globale(10)
; @ In[15]:3 within `addition_variable_globale`
define nonnull {}* @julia_addition_variable_globale_1886(i64 signext %0) #0 {
top:
%1 = alloca [2 x {}*], align 8
%gcframe2 = alloca [4 x {}*], align 16
%gcframe2.sub = getelementptr inbounds [4 x {}*], [4 x {}*]* %gcframe2, i64 0, i64 0
%.sub = getelementptr inbounds [2 x {}*], [2 x {}*]* %1, i64 0, i64 0
%2 = bitcast [4 x {}*]* %gcframe2 to i8*
call void @llvm.memset.p0i8.i64(i8* align 16 %2, i8 0, i64 32, i1 true)
%thread_ptr = call i8* asm "movq %fs:0, $0", "=r"() #6
%tls_ppgcstack = getelementptr i8, i8* %thread_ptr, i64 -8
%3 = bitcast i8* %tls_ppgcstack to {}****
%tls_pgcstack = load {}***, {}**** %3, align 8
; @ In[15]:4 within `addition_variable_globale`
%4 = bitcast [4 x {}*]* %gcframe2 to i64*
store i64 8, i64* %4, align 16
%5 = getelementptr inbounds [4 x {}*], [4 x {}*]* %gcframe2, i64 0, i64 1
%6 = bitcast {}** %5 to {}***
%7 = load {}**, {}*** %tls_pgcstack, align 8
store {}** %7, {}*** %6, align 8
%8 = bitcast {}*** %tls_pgcstack to {}***
store {}** %gcframe2.sub, {}*** %8, align 8
%variable = load atomic {}*, {}** inttoptr (i64 139910041751984 to {}**) unordered, align 16
%9 = getelementptr inbounds [4 x {}*], [4 x {}*]* %gcframe2, i64 0, i64 2
store {}* %variable, {}** %9, align 16
%10 = call nonnull {}* @ijl_box_int64(i64 signext %0)
%11 = getelementptr inbounds [4 x {}*], [4 x {}*]* %gcframe2, i64 0, i64 3
store {}* %10, {}** %11, align 8
store {}* %10, {}** %.sub, align 8
%12 = getelementptr inbounds [2 x {}*], [2 x {}*]* %1, i64 0, i64 1
store {}* %variable, {}** %12, align 8
%13 = call nonnull {}* @ijl_apply_generic({}* inttoptr (i64 139909873737088 to {}*), {}** nonnull %.sub, i32 2)
%14 = load {}*, {}** %5, align 8
%15 = bitcast {}*** %tls_pgcstack to {}**
store {}* %14, {}** %15, align 8
ret {}* %13
}
Il est donc possible d’améliorer la performance en fixant la valeur de la variable globale avec l’instruction const.
const constante = 10
function addition_variable_constante(x)
x + constante
end
@btime addition_variable_constante(10)
1.271 ns (0 allocations: 0 bytes)
20
On peut également fixer le type de cette variable. C’est mieux mais cela reste éloigné, en terme de performance, du résultat précedent.
function addition_variable_typee(x)
x + variable::Int
end
@btime addition_variable_typee(10)
2.277 ns (0 allocations: 0 bytes)
20
Pour régler notre problème de performance avec une variable globale, il faut la passer en argument dans la fonction.
function addition_variable_globale_en_argument(x, v)
x + v
end
addition_variable_globale_en_argument (generic function with 1 method)
@btime addition_variable_globale_en_argument(10, $variable)
2.781 ns (0 allocations: 0 bytes)
20
Instabilité de type#
Une fonction est de type stable lorsque vous pouvez déduire ce que doit être la sortie de la fonction. L’exemple ci-dessous rendra les choses plus claires. En règle générale, les fonctions de type stable sont plus rapides.
function carre_plus_un(v::T) where T <:Number
g = v * v
return g+1
end
carre_plus_un (generic function with 1 method)
v = rand()
0.2692504830309297
@code_warntype carre_plus_un(v)
MethodInstance for carre_plus_un(::Float64)
from carre_plus_un(v::T) where T<:Number @ Main In[23]:1
Static Parameters
T = Float64
Arguments
#self#::Core.Const(carre_plus_un)
v::Float64
Locals
g::Float64
Body::Float64
1 ─ (g = v * v)
│ %2 = (g + 1)::Float64
└── return %2
w = 5
5
@code_warntype carre_plus_un(w)
MethodInstance for carre_plus_un(::Int64)
from carre_plus_un(v::T) where T<:Number @ Main In[23]:1
Static Parameters
T = Int64
Arguments
#self#::Core.Const(carre_plus_un)
v::Int64
Locals
g::Int64
Body::Int64
1 ─ (g = v * v)
│ %2 = (g + 1)::Int64
└── return %2
Sur les deux exemples précedents on peut déduire le type de sortie de la fonction.
function carre_plus_un(v::T) where T <:Number
g = v*v # Type(T * T) ==> T
return g+1 # Type(T + Int)) ==> "max" (T,Int)
end
Le type de la valeur de retour peut être différent: Float64 ou Int64. Mais la fonction est toujours stable.
Créons maintenant un nouveau type:
mutable struct Cube
length
width
height
end
volume(c::Cube) = c.length*c.width*c.height
volume (generic function with 1 method)
mutable struct Cube_typed
length::Float64
width::Float64
height::Float64
end
volume(c::Cube_typed) = c.length*c.width*c.height
volume (generic function with 2 methods)
mutable struct Cube_parametric_typed{T <: Real}
length :: T
width :: T
height :: T
end
volume(c::Cube_parametric_typed) = c.length*c.width*c.height
volume (generic function with 3 methods)
c1 = Cube(1.1,1.2,1.3)
c2 = Cube_typed(1.1,1.2,1.3)
c3 = Cube_parametric_typed(1.1,1.2,1.3)
@show volume(c1) == volume(c2) == volume(c3)
volume(c1) == volume(c2) == volume(c3) = true
true
using BenchmarkTools
@btime volume(c1) # not typed
@btime volume(c2) # typed float
@btime volume(c3) # typed parametric
23.127 ns (1 allocation: 16 bytes)
6.572 ns (1 allocation: 16 bytes)
17.361 ns (1 allocation: 16 bytes)
1.7160000000000002
c4 = Cube_parametric_typed{Float64}(1.1,1.2,1.3)
@btime volume(c4)
17.300 ns (1 allocation: 16 bytes)
1.7160000000000002
The second and the third function calls are faster! Let’s call @code_warntype and check type stability
@code_warntype volume(c1)
MethodInstance for volume(::Cube)
from volume(c::Cube) @ Main In[29]:1
Arguments
#self#::Core.Const(volume)
c::Cube
Body::Any
1 ─ %1 = Base.getproperty(c, :length)::Any
│ %2 = Base.getproperty(c, :width)::Any
│ %3 = Base.getproperty(c, :height)::Any
│ %4 = (%1 * %2 * %3)::Any
└── return %4
@code_warntype volume(c2)
MethodInstance for volume(::Cube_typed)
from volume(c::Cube_typed) @ Main In[30]:6
Arguments
#self#::Core.Const(volume)
c::Cube_typed
Body::Float64
1 ─ %1 = Base.getproperty(c, :length)::Float64
│ %2 = Base.getproperty(c, :width)::Float64
│ %3 = Base.getproperty(c, :height)::Float64
│ %4 = (%1 * %2 * %3)::Float64
└── return %4
@code_warntype volume(c3)
MethodInstance for volume(::Cube_parametric_typed{Float64})
from volume(c::Cube_parametric_typed) @ Main In[31]:6
Arguments
#self#::Core.Const(volume)
c::Cube_parametric_typed{Float64}
Body::Float64
1 ─ %1 = Base.getproperty(c, :length)::Float64
│ %2 = Base.getproperty(c, :width)::Float64
│ %3 = Base.getproperty(c, :height)::Float64
│ %4 = (%1 * %2 * %3)::Float64
└── return %4
Conclusion: Les types en Julia sont importants donc si vous les connaissez, ajoutez-les, cela peut améliorer les performances.
function zero_or_val(x::Real)
if x >= 0
return x
else
return 0
end
end
@code_warntype zero_or_val(0.2)
MethodInstance for zero_or_val(::Float64)
from zero_or_val(x::Real) @ Main In[38]:1
Arguments
#self#::Core.Const(zero_or_val)
x::Float64
Body::Union{Float64, Int64}
1 ─ %1 = (x >= 0)::Bool
└── goto #3 if not %1
2 ─ return x
3 ─ return 0
function zero_or_val_stable(x::Real)
T = promote_type(typeof(x),Int)
if x >= 0
y = x
else
y = zero(T)
end
return y
end
@code_warntype zero_or_val_stable(0.2)
MethodInstance for zero_or_val_stable(::Float64)
from zero_or_val_stable(x::Real) @ Main In[39]:1
Arguments
#self#::Core.Const(zero_or_val_stable)
x::Float64
Locals
y::Float64
T::Type{Float64}
Body::Float64
1 ─ Core.NewvarNode(:(y))
│ %2 = Main.typeof(x)::Core.Const(Float64)
│ (T = Main.promote_type(%2, Main.Int))
│ %4 = (x >= 0)::Bool
└── goto #3 if not %4
2 ─ (y = x)
└── goto #4
3 ─ (y = Main.zero(T::Core.Const(Float64)))
4 ┄ return y
Conclusion: promote_type peut vous permettre de supprimer une instabilité de type en utilisant la réprésentation la plus haute dans l’abre des types.
Je vous propose le jeu suivant: Soit un vecteur de nombres. Calculons la somme comme suit.
Pour chaque nombre du vecteur, on lance une pièce de monnaie (rand()), si c’est face (>=0.5), vous ajoutez 1. Sinon, vous ajoutez le nombre lui-même.
function flipcoin_then_add(v::Vector{T}) where T <: Real
s = 0
for vi in v
r = rand()
if r >=0.5
s += 1
else
s += vi
end
end
end
function flipcoin_then_add_typed(v::Vector{T}) where T <: Real
s = zero(T)
for vi in v
r = rand()
if r >=0.5
s += one(T)
else
s += vi
end
end
end
myvec = rand(1000)
@show flipcoin_then_add(myvec) == flipcoin_then_add_typed(myvec)
flipcoin_then_add(myvec) == flipcoin_then_add_typed(myvec) = true
true
@btime flipcoin_then_add(rand(1000))
@btime flipcoin_then_add_typed(rand(1000))
7.114 μs (1 allocation: 7.94 KiB)
1.247 μs (1 allocation: 7.94 KiB)
Conclusion: Think about the variables you are declaring. Do you know their types? If so, specify the type somehow.
@code_XXX#
Nous avons vu durant ce chapitre que regarder le code généré peut nous aider à améliorer les performances. Voici toutes les macros à votre disposition:
# @code_llvm
# @code_lowered
# @code_native
# @code_typed
# @code_warntype
function flipcoin(randval::Float64)
if randval<0.5
return "H"
else
return "T"
end
end
flipcoin (generic function with 1 method)
@code_lowered flipcoin(rand()) # syntax tree
CodeInfo(
1 ─ %1 = randval < 0.5
└── goto #3 if not %1
2 ─ return "H"
3 ─ return "T"
)
@code_warntype flipcoin(rand()) # try @code_typed
MethodInstance for flipcoin(::Float64)
from flipcoin(randval::Float64) @ Main In[42]:7
Arguments
#self#::Core.Const(flipcoin)
randval::Float64
Body::String
1 ─ %1 = (randval < 0.5)::Bool
└── goto #3 if not %1
2 ─ return "H"
3 ─ return "T"
@code_llvm flipcoin(rand()) # this and code_warntype are probably the most relevant
; @ In[42]:7 within `flipcoin`
define nonnull {}* @julia_flipcoin_2470(double %0) #0 {
top:
; @ In[42]:8 within `flipcoin`
; ┌ @ float.jl:536 within `<`
%1 = fcmp uge double %0, 5.000000e-01
%. = select i1 %1, {}* inttoptr (i64 139910097979208 to {}*), {}* inttoptr (i64 139910097979152 to {}*)
; └
; @ In[42] within `flipcoin`
ret {}* %.
}
@code_native flipcoin(rand())
.text
.file "flipcoin"
.section .rodata.cst8,"aM",@progbits,8
.p2align 3 # -- Begin function julia_flipcoin_2486
.LCPI0_0:
.quad 0x3fe0000000000000 # double 0.5
.text
.globl julia_flipcoin_2486
.p2align 4, 0x90
.type julia_flipcoin_2486,@function
julia_flipcoin_2486: # @julia_flipcoin_2486
; ┌ @ In[42]:7 within `flipcoin`
# %bb.0: # %top
push rbp
mov rbp, rsp
movabs rax, offset .LCPI0_0
vmovsd xmm1, qword ptr [rax] # xmm1 = mem[0],zero
; │ @ In[42]:8 within `flipcoin`
; │┌ @ float.jl:536 within `<`
vucomisd xmm1, xmm0
movabs rcx, 139910097979208
movabs rax, 139910097979152
cmovbe rax, rcx
; │└
; │ @ In[42] within `flipcoin`
pop rbp
ret
.Lfunc_end0:
.size julia_flipcoin_2486, .Lfunc_end0-julia_flipcoin_2486
; └
# -- End function
.section ".note.GNU-stack","",@progbits