Julia 代碼性能優(yōu)化

代碼性能優(yōu)化

以下幾節(jié)將描述一些提高 Julia 代碼運(yùn)行速度的技巧。

避免全局變量

全局變量的值、類型，都可能變化。這使得編譯器很難優(yōu)化使用全局變量的代碼。應(yīng)盡量使用局部變量，或者把變量當(dāng)做參數(shù)傳遞給函數(shù)。

對(duì)性能至關(guān)重要的代碼，應(yīng)放入函數(shù)中。

聲明全局變量為常量可以顯著提高性能：

const DEFAULT_VAL = 0

使用非常量的全局變量時(shí)，最好在使用時(shí)指明其類型，這樣也能幫助編譯器優(yōu)化：

global x
y = f(x::Int + 1)

寫函數(shù)是一種更好的風(fēng)格，這會(huì)產(chǎn)生更多可重復(fù)和清晰的代碼，也包括清晰的輸入和輸出。

使用 @time 來衡量性能并且留心內(nèi)存分配

衡量計(jì)算性能最有用的工具是 @time 宏。下面的例子展示了良好的使用方式 :

  julia> function f(n)
             s = 0
             for i = 1:n
                 s += i/2
             end
             s
          end
  f (generic function with 1 method)

  julia> @time f(1)
  elapsed time: 0.008217942 seconds (93784 bytes allocated)
  0.5

  julia> @time f(10^6)
  elapsed time: 0.063418472 seconds (32002136 bytes allocated)
  2.5000025e11

在第一次調(diào)用時(shí) (@time f(1)), f 會(huì)被編譯. (如果你在這次會(huì)話中還 沒有使用過 @time, 計(jì)時(shí)函數(shù)也會(huì)被編譯.) 這時(shí)的結(jié)果沒有那么重要. 在 第二次調(diào)用時(shí), 函數(shù)打印了執(zhí)行所耗費(fèi)的時(shí)間, 同時(shí)請(qǐng)注意, 在這次執(zhí)行過程中 分配了一大塊的內(nèi)存. 相對(duì)于函數(shù)形式的 tic 和 toc, 這是 @time 宏的一大優(yōu)勢(shì).

出乎意料的大塊內(nèi)存分配往往意味著程序的某個(gè)部分存在問題, 通常是關(guān)于類型 穩(wěn)定性. 因此, 除了關(guān)注內(nèi)存分配本身的問題, 很可能 Julia 為你的函數(shù)生成 的代碼存在很大的性能問題. 這時(shí)候要認(rèn)真對(duì)待這些問題并遵循下面的一些個(gè)建 議.

另外, 作為一個(gè)引子, 上面的問題可以優(yōu)化為無內(nèi)存分配 (除了向 REPL 返回結(jié) 果), 計(jì)算速度提升 30 倍 ::

  julia> @time f_improved(10^6)
  elapsed time: 0.00253829 seconds (112 bytes allocated)
  2.5000025e11

你可以從下面的章節(jié)學(xué)到如何識(shí)別 f 存在的問題并解決。

在有些情況下, 你的函數(shù)可能需要為本身的操作分配內(nèi)存, 這樣會(huì)使得問題變得 復(fù)雜. 在這種情況下, 可以考慮使用下面的 :ref:工具? 之一來甄別問題, 或者將函數(shù)拆分, 一部分處理內(nèi)存分配, 另一部分處理算法 (參見 :ref:?預(yù)分配內(nèi)存）。 

 工具 

 Julia 提供了一些工具包來鑒別性能問題所在 :

•  profiling)可以用來衡量代碼的性能, 同時(shí)鑒別出瓶頸所在。對(duì)于復(fù)雜的項(xiàng)目, 可以使用 ?ProfileView  ?擴(kuò)展包來直觀的展示分析結(jié)果. 
•  出乎意料的大塊內(nèi)存分配, -- ?@time?, ?@allocated?, 或者 -profiler - 意味著你的代碼可能存在問題. 如果你看不出內(nèi)存分配的問題, -那么類型系統(tǒng)可能存在問題. 也可以使用 ?--track-allocation=user? 來 -啟動(dòng) Julia, 然后查看? *.mem? 文件來找出內(nèi)存分配是在哪里出現(xiàn)的. 
• ?TypeCheck ?_ 擴(kuò)展包可以指出程序一 些問題. 

避免包含一些抽象類型參數(shù) 

 當(dāng)運(yùn)行參數(shù)化類型時(shí)候，比如 arrays，如果有可能最好去避免使用抽象類型參數(shù)。 思考下面的代碼： 

a = Real[] # typeof(a) = Array{Real,1} 
if (f = rand()) x = [1 2; 3 4] 
2x2 Array{Int64,2}:
 1 2
 3 4 
julia> x[:]
4-element Array{Int64,1}: 
1
3 
2 
4

 這種給數(shù)組排序的約定在許多語言中都是常見的，比如 Fortran ， Matlab ，和 R 語言(舉幾個(gè)例子來說)。以列為主序的另一選擇就是以行為主序，其它語言中的 C 語言和 Python 語言(?numpy?)就是選用了這種方式。記住數(shù)組的順序?qū)?shù)組的查找有著至關(guān)重要的影響。要記住的一個(gè)查找規(guī)則就是對(duì)于基于列為順序的數(shù)組，第一個(gè)指針是變化最快的。這基本上就意味著如果在一段代碼中，循環(huán)指針是第一個(gè)，那么查找速度會(huì)更快。 

 我們來看一下下面這個(gè)人為的例子。假設(shè)我們想要實(shí)現(xiàn)一個(gè)功能，接收一個(gè) ?Vector ?并且返回一個(gè)方形的 ?Matrix?，且行或列為輸入矢量的復(fù)制。我們假設(shè)是行還是列為數(shù)據(jù)的復(fù)制并不重要（或許剩下的代碼可以相應(yīng)地更容易的適應(yīng)）。我們可以想到有至少四種方法可以實(shí)現(xiàn)這一點(diǎn)（除了建議的回訪正建的 ?repmat ?功能）： 

function copy_cols{T}(x::Vector{T}) 
    n = size(x, 1) 
    out = Array(eltype(x), n, n) 
    for i=1:n 
        out[:, i] = x 
    end 
    out 
end 
function copy_rows{T}(x::Vector{T}) 
    n = size(x, 1) 
    out = Array(eltype(x), n, n) 
    for i=1:n 
        out[i, :] = x 
    end 
    out 
end 
function copy_col_row{T}(x::Vector{T}) n = size(x, 1)
    out = Array(T, n, n) 
    for col=1:n, row=1:n 
        out[row, col] = x[row] 
    end 
    out 
end 
function copy_row_col{T}(x::Vector{T}) n = size(x, 1) 
    out = Array(T, n, n) 
    for row=1:n, col=1:n 
        out[row, col] = x[col] 
    end 
    out 
end 

 現(xiàn)在我們使用同樣的輸入向量 ?1? 產(chǎn)生的隨機(jī)數(shù) 10000? 給每個(gè)功能計(jì)時(shí)：

julia> x = randn(10000); 

julia> fmt(f) = println(rpad(string(f)*": ", 14, ' '), @elapsed f(x)) 

julia> map(fmt, {copy_cols, copy_rows, copy_col_row, copy_row_col}); 
copy_cols:    0.331706323 
copy_rows:    1.799009911 
copy_col_row: 0.415630047 
copy_row_col: 1.721531501 

 注意到 ?copy_cols ?比?copy_rows?快很多。這是意料之中的，因?yàn)??copy_cols ?遵守 ?Matrix?`界面的基于列的存儲(chǔ)，并且一次就填滿一列。除此之外，?copy_col_row?比 ?copy_row_col ?快很多，因?yàn)樗衔覀兊牟檎乙?guī)則，即在一段代碼中第一個(gè)出現(xiàn)的元素應(yīng)該是與最內(nèi)部的循環(huán)相聯(lián)系的。

 輸出預(yù)先分配

 如果你的功能返回了一個(gè) Array 或其它復(fù)雜類型，它可能不得不分配內(nèi)存。不幸的是，時(shí)常分配和它的相反事件，垃圾區(qū)收集，是有實(shí)質(zhì)性瓶頸的。 

 有時(shí)候，你可以在訪問每個(gè)功能時(shí)通過預(yù)先分配輸出來避開分配內(nèi)存的需要。作為一個(gè)很小的例子，比較一下 

function xinc(x) 
    return [x, x+1, x+2] 
end 
function loopinc() 
    y = 0 
    for i = 1:10^7 
        ret = xinc(i) 
        y += ret[2] 
    end y 
end

 和 

function xinc!{T}(ret::AbstractVector{T}, x::T) 
    ret[1] = x 
    ret[2] = x+1 
    ret[3] = x+2 
    nothing 
end 

function loopinc_prealloc() 
    ret = Array(Int, 3) 
    y = 0 
    for i = 1:10^7 
        xinc!(ret, i) 
        y += ret[2] 
    end 
    y 
end 

計(jì)時(shí)結(jié)果：

 julia> @time loopinc() 
elapsed time: 1.955026528 seconds (1279975584 bytes allocated) 
50000015000000 

julia> @time loopinc_prealloc() 
elapsed time: 0.078639163 seconds (144 bytes allocated) 
50000015000000

 預(yù)先分配有其他好處，比如，允許訪問者通過算法控制“輸出”類型。在上面的例子中，我們可以按照自己希望的，通過一個(gè)?SubArray ?而不是 ?Array?。

 按著最極端的來想，預(yù)先分配可以讓你的代碼看起來丑點(diǎn)，所以需要一些表達(dá)方式和判斷。 

避免輸入/輸出時(shí)的串插入 

 把數(shù)據(jù)寫入文件（或者其他輸入/輸出設(shè)備）時(shí)，中間字符串的形成是額外的開銷。而不是：

println(file, "$a $b")

使用：

println(file, a, " ", b) 

第一種代碼形成了一個(gè)字符串，然后把它寫入了文件，而第二種代碼直接把值寫入了文件。同樣也注意到在某些情況下，字符串的插入很難讀出來?？紤]一下：

 println(file, "$(f(a))$(f(b))")

 對(duì)比：

println(file, f(a), f(b))

 處理有關(guān)舍棄的警告

 被舍棄的函數(shù)，會(huì)查表并顯示一次警告，而這會(huì)影響性能。建議按照警告的提示進(jìn)行對(duì)應(yīng)的修改。

 小技巧

 注意些有些小事項(xiàng)，能使內(nèi)部循環(huán)更緊致。

• 避免不必要的數(shù)組。例如，不要使用 ?sum([x,y,z])?，而應(yīng)使用 ?x+y+z?
•  對(duì)于較小的整數(shù)冪，使用 ?*? 更好。如 ?x*x*x? 比 ?x^3? 好 
• 針對(duì)復(fù)數(shù) ?z? ，使用 ?abs2(z) ?代替 ?abs(z)^2 ?。一般情況下，對(duì)于復(fù)數(shù)參數(shù)，盡量用 ?abs2? 代替 ?abs ?
• 對(duì)于整數(shù)除法，使用 ?div(x,y)? 而不是 ?trunc(x/y)?, 使用 ?fld(x,y)? 而不是 ?floor(x/y)?, 使用 ?cld(x,y)? 而不是 ?ceil(x/y)?.

 性能注釋 

有時(shí)你可以設(shè)定某些項(xiàng)目屬性來獲得更好的優(yōu)化。

•  在檢查公式時(shí)，使用 ?@inbounds?來消除數(shù)組界限。一定要在這之前完成。如果下標(biāo)越界了，你可能會(huì)遇到崩潰或不執(zhí)行的問題。 
• 在 ?for?循環(huán)之前寫上 ?@simd?，這個(gè)可以幫你檢驗(yàn)。這個(gè)特征是試驗(yàn)性的而且在之后的 Julia 版本中可能會(huì)改變會(huì)消失。

 這里有一個(gè)包含兩種形式審定的例子： 

function inner( x, y ) 
    s = zero(eltype(x)) 
    for i=1:length(x) 
        @inbounds s += x[i]*y[i] 
    end 
    s 
end 

function innersimd( x, y ) 
    s = zero(eltype(x)) 
    @simd for i=1:length(x) 
        @inbounds s += x[i]*y[i] 
    end 
    s 
end 

function timeit( n, reps ) 
    x = rand(Float32,n) 
    y = rand(Float32,n) 
    s = zero(Float64) 
    time = @elapsed for j in 1:reps 
        s+=inner(x,y) 
    end 
    println("GFlop = ",2.0*n*reps/time*1E-9) 
    time = @elapsed for j in 1:reps 
        s+=innersimd(x,y) 
    end 
    println("GFlop (SIMD) = ",2.0*n*reps/time*1E-9) 
end 

timeit(1000,1000) 

在配有 2.4GHz 的 Intel Core i5 處理器的電腦上，產(chǎn)生如下結(jié)果：

GFlop = 1.9467069505224963 
GFlop (SIMD) = 17.578554163920018 

?@simd for? 循環(huán)應(yīng)該是一維范圍的。縮減變數(shù)是用于累積變量的，比如例子中的 ?s?。通過使用 ?@simd?，你可以維護(hù)循環(huán)的幾種性能：

•  -有縮減變數(shù)的特殊考慮后，在任意的或重疊的順序中執(zhí)行迭代都是安全的。 
•  減少變量的浮點(diǎn)操作可以被重復(fù)執(zhí)行，但是可能會(huì)比沒有 ?@simd? 產(chǎn)生不同的結(jié)果。 
• 不會(huì)有一個(gè)迭代在等待另一個(gè)迭代，以實(shí)現(xiàn)前進(jìn)。 

使用 ?@simd? 僅僅是給了編譯器矢量化的通行證。它是不是真的會(huì)這樣做還取決于編譯器。要真正從當(dāng)前的實(shí)現(xiàn)中獲益，你的循環(huán)應(yīng)該有如下額外的性能：

•   循環(huán)必須是內(nèi)部循環(huán)。 
• 循環(huán)主題必須是無循環(huán)程序。這就是為什么當(dāng)前所有的數(shù)組訪問都需要 ?@inbounds?的原因了。 
• 訪問必須有一個(gè)跨越模式，而且不能“聚集”（隨機(jī)指針讀?。┗蛘摺胺稚ⅰ保S機(jī)指針寫入）。 
• 跨越應(yīng)該是單元跨越。 
• 在一些簡(jiǎn)單的例子中，例如一個(gè) 2-3 數(shù)組訪問的循環(huán)中，LLVM 自動(dòng)矢量化可能會(huì)自動(dòng)生效，導(dǎo)致無需? @simd? 的進(jìn)一步加速。