資訊專欄INFORMATION COLUMN
摘要：首發于接上一篇讓代碼飛起來高性能學習筆記一，繼續整理高性能學習筆記。和都只能表示特定的整數范圍，超過范圍會。通用代碼一般會用，這就有可能導致性能問題。
首發于 https://magicly.me/hpc-julia-2/
接上一篇：讓代碼飛起來——高性能 Julia 學習筆記（一）， 繼續整理高性能 Julia 學習筆記。
數字Julia 中 Number 的 size 就跟 C 語言里面一樣， 直接依賴于底層的 CPU/OS， 32 位 OS 上 integer 默認是 32 位， 64 位 OS 上 integer 默認是 64 位。
可以用bitstring查看 number 的二進制表示：
julia > bitstring(3) "0000000000000000000000000000000000000000000000000000000000000011" julia > bitstring(-3) "1111111111111111111111111111111111111111111111111111111111111101"
有時候數字可能會被 boxed，Julia 的 compiler 會自動 boxing/unboxing。
Int64 和 Int32 都只能表示特定的整數范圍， 超過范圍會 overflow。
julia> typemax(Int64) 9223372036854775807 julia> bitstring(typemax(Int64)) "0111111111111111111111111111111111111111111111111111111111111111" julia> typemin(Int64) -9223372036854775808 julia> bitstring(typemin(Int64)) "1000000000000000000000000000000000000000000000000000000000000000" julia> typemax(Int64) + 1 -9223372036854775808 julia> typemin(Int64) - 1 9223372036854775807
如果要表示任意精度的話， 可以用BitInt。
julia> big(typemax(Int64)) + 1 9223372036854775808
float point 遵循IEEE 754標準。
julia> bitstring(1.5) "0011111111111000000000000000000000000000000000000000000000000000" julia> bitstring(-1.5) "1011111111111000000000000000000000000000000000000000000000000000"
無符號整數 Unsigned integer 可以用 UInt64 和 UInt32 表示， 不同類型可以轉， 但是如果超出可表示范圍會報錯。
julia> UInt64(UInt32(1)) 0x0000000000000001 julia> UInt32(UInt64(1)) 0x00000001 julia> UInt32(typemax(UInt64)) ERROR: InexactError: trunc(UInt32, 18446744073709551615) Stacktrace: [1] throw_inexacterror(::Symbol, ::Any, ::UInt64) at ./boot.jl:567 [2] checked_trunc_uint at ./boot.jl:597 [inlined] [3] toUInt32 at ./boot.jl:686 [inlined] [4] UInt32(::UInt64) at ./boot.jl:721 [5] top-level scope at none:0
有時候不需要考慮是否溢出， 可以直接用%符號取低位， 速度還更快：
julia> (typemax(UInt64) - 1) % UInt32 0xfffffffe
精度和效率平衡
有時候為了更高的效率， 可以使用@fastmath宏， 它會放寬一些限制， 包括檢查 NaN 或 Infinity 等， 類似于 GCC 中的-ffast-math編譯選項。
julia> function sum_diff(x)
n = length(x); d = 1/(n-1)
s = zero(eltype(x))
s = s + (x[2] - x[1]) / d
for i = 2:length(x)-1
s = s + (x[i+1] - x[i+1]) / (2*d)
end
s = s + (x[n] - x[n-1])/d
end
sum_diff (generic function with 1 method)
julia> function sum_diff_fast(x)
n=length(x); d = 1/(n-1)
s = zero(eltype(x))
@fastmath s = s + (x[2] - x[1]) / d
@fastmath for i = 2:n-1
s = s + (x[i+1] - x[i+1]) / (2*d)
end
@fastmath s = s + (x[n] - x[n-1])/d
end
sum_diff_fast (generic function with 1 method)
julia> t=rand(2000);
julia> sum_diff(t)
1124.071808538789
julia> sum_diff_fast(t)
1124.0718085387887
julia> using BenchmarkTools
julia> @benchmark sum_diff(t)
BenchmarkTools.Trial:
memory estimate: 16 bytes
allocs estimate: 1
--------------
minimum time: 4.447 μs (0.00% GC)
median time: 4.504 μs (0.00% GC)
mean time: 4.823 μs (0.00% GC)
maximum time: 17.318 μs (0.00% GC)
--------------
samples: 10000
evals/sample: 7
julia> @benchmark sum_diff_fast(t)
BenchmarkTools.Trial:
memory estimate: 16 bytes
allocs estimate: 1
--------------
minimum time: 574.951 ns (0.00% GC)
median time: 580.831 ns (0.00% GC)
mean time: 621.044 ns (1.04% GC)
maximum time: 65.017 μs (99.06% GC)
--------------
samples: 10000
evals/sample: 183
差異還是蠻大的， 差不多 8 倍差距！ 我們可以用macroexpand看看宏展開的代碼：
julia> ex = :(@fastmath for i in 2:n-1
s = s + (x[i+1] - x[i+1]) / (2*d)
end)
:(#= REPL[57]:1 =# @fastmath for i = 2:n - 1
#= REPL[57]:2 =#
s = s + (x[i + 1] - x[i + 1]) / (2d)
end)
julia> macroexpand(Main, ex)
:(for i = 2:(Base.FastMath).sub_fast(n, 1)
#= REPL[57]:2 =#
s = (Base.FastMath).add_fast(s, (Base.FastMath).div_fast((Base.FastMath).sub_fast(x[(Base.FastMath).add_fast(i, 1)], x[(Base.FastMath).add_fast(i, 1)]), (Base.FastMath).mul_fast(2, d)))
end)
可以看到， 主要是把普通的加減乘除換成了Base.FastMath中的方法。
本章最后介紹了K-B-N求和減少誤差，以及 Subnormal numbers， 感覺都是科學計算里面才會用到的， 暫時沒管。
數組科學計算以及人工智能里面有大量向量、矩陣運算， 在 Julia 里都直接對應到 Array。 Vector 和 Matrix 其實就是 Array 的特例：
julia> Vector
Array{T,1} where T
julia> Matrix
Array{T,2} where T
即 Vector 是一維數組， Matrix 是二維數組。從這里也可以看出， Julia 中類型的參數類型可以是具體的 value， 比如這里的 1 和 2。
Julia 中 Array 數據是連續存放的， 如下圖：
跟存放指針相比好處是少了一次內存訪問， 并且可以更好地利用 CPU 的 pipeline 和 cache，以及 SIMD。
Julia 中多維數組是按列優先存儲的(類似 Matlab 和 Fortran)，這點跟 C 語言中不一樣：
按照數據在內存中的布局來讀取數據， 能顯著提高效率。
julia> function col_iter(x)
s=zero(eltype(x))
for i = 1:size(x, 2)
for j = 1:size(x, 1)
s = s + x[j, i] ^ 2
x[j, i] = s
end
end
end
col_iter (generic function with 1 method)
julia> function row_iter(x)
s=zero(eltype(x))
for i = 1:size(x, 1)
for j = 1:size(x, 2)
s = s + x[i, j] ^ 2
x[i, j] = s
end
end
end
row_iter (generic function with 1 method)
julia> a = rand(1000, 1000);
julia> @benchmark row_iter(a)
BenchmarkTools.Trial:
memory estimate: 0 bytes
allocs estimate: 0
--------------
minimum time: 2.217 ms (0.00% GC)
median time: 2.426 ms (0.00% GC)
mean time: 2.524 ms (0.00% GC)
maximum time: 11.723 ms (0.00% GC)
--------------
samples: 1974
evals/sample: 1
julia> @benchmark col_iter(a)
BenchmarkTools.Trial:
memory estimate: 0 bytes
allocs estimate: 0
--------------
minimum time: 815.984 μs (0.00% GC)
median time: 910.121 μs (0.00% GC)
mean time: 917.850 μs (0.00% GC)
maximum time: 1.292 ms (0.00% GC)
--------------
samples: 5416
evals/sample: 1
Julia runtime 會對 array 訪問做 bound 判斷， 看是否超出邊界。 超出邊界的訪問會導致很多 bugs，甚至是安全問題。 另一方面， bound check 會帶來額外的開銷，如果你能很確定不會越界， 那可以用@inbound 宏告訴 Julia 不要做 bound check， 效率會提高不少。
julia> function prefix_bounds(a, b)
for i = 2:size(a, 1)
a[i] = b[i-1] + b[i]
end
end
prefix_bounds (generic function with 1 method)
julia> function prefix_inbounds(a, b)
@inbounds for i = 2:size(a, 1)
a[i] = b[i-1] + b[i]
end
end
prefix_inbounds (generic function with 1 method)
julia> a = rand(1000);
julia> b = rand(1000);
julia> @benchmark prefix_bounds(a, b)
BenchmarkTools.Trial:
memory estimate: 0 bytes
allocs estimate: 0
--------------
minimum time: 742.360 ns (0.00% GC)
median time: 763.264 ns (0.00% GC)
mean time: 807.497 ns (0.00% GC)
maximum time: 1.968 μs (0.00% GC)
--------------
samples: 10000
evals/sample: 125
julia> @benchmark prefix_inbounds(a, b)
BenchmarkTools.Trial:
memory estimate: 0 bytes
allocs estimate: 0
--------------
minimum time: 179.294 ns (0.00% GC)
median time: 181.826 ns (0.00% GC)
mean time: 185.124 ns (0.00% GC)
maximum time: 635.909 ns (0.00% GC)
--------------
samples: 10000
evals/sample: 690
慎用@inbound?。?！ 最好是限制 loop 直接依賴于 array 長度， 比如：for i in 1:length(array)這種形式。
可以在啟動的時候加上--check-bounds=yes/no來全部開啟或者關閉（優先級高于@inbound 宏）bound check！ 再說一次， 慎用！
Julia 內置了很多操作 Array 的函數， 一般都提供兩個版本， 注意可變和不可變版本差異很大?。?！
julia> a = rand(1000); julia> @benchmark sort(a) BenchmarkTools.Trial: memory estimate: 7.94 KiB allocs estimate: 1 -------------- minimum time: 16.050 μs (0.00% GC) median time: 17.493 μs (0.00% GC) mean time: 18.933 μs (0.00% GC) maximum time: 726.416 μs (0.00% GC) -------------- samples: 10000 evals/sample: 1 julia> @benchmark sort!(a) BenchmarkTools.Trial: memory estimate: 0 bytes allocs estimate: 0 -------------- minimum time: 4.868 μs (0.00% GC) median time: 4.997 μs (0.00% GC) mean time: 5.282 μs (0.00% GC) maximum time: 13.772 μs (0.00% GC) -------------- samples: 10000 evals/sample: 7
我們可以看到， 不可變版本sort無論時間還是內存占用和分配上都比可變版本sort!高。 這很容易理解， 不可變版本需要 clone 一份數據出來，而不是直接修改參數。 根據需要選擇最適合的版本。
類似的，合理利用預分配，可以有效降低時間和內存占用。
julia> function xpow(x)
return [x x^2 x^3 x^4]
end
xpow (generic function with 1 method)
julia> function xpow_loop(n) s= 0
for i = 1:n
s = s + xpow(i)[2]
end
s
end
xpow_loop (generic function with 1 method)
julia> @benchmark xpow_loop(1_000_000)
BenchmarkTools.Trial:
memory estimate: 167.84 MiB
allocs estimate: 4999441
--------------
minimum time: 68.342 ms (4.11% GC)
median time: 70.378 ms (5.21% GC)
mean time: 71.581 ms (6.04% GC)
maximum time: 120.430 ms (39.92% GC)
--------------
samples: 70
evals/sample: 1
julia> function xpow!(result::Array{Int, 1}, x)
@assert length(result) == 4
result[1] = x
result[2] = x^2
result[3] = x^3
result[4] = x^4
end
xpow! (generic function with 1 method)
julia> function xpow_loop_noalloc(n)
r = [0, 0, 0, 0]
s= 0
for i = 1:n
xpow!(r, i)
s = s + r[2]
end
s
end
xpow_loop_noalloc (generic function with 1 method)
julia> @benchmark xpow_loop_noalloc(1_000_000)
BenchmarkTools.Trial:
memory estimate: 112 bytes
allocs estimate: 1
--------------
minimum time: 7.314 ms (0.00% GC)
median time: 7.486 ms (0.00% GC)
mean time: 7.599 ms (0.00% GC)
maximum time: 9.580 ms (0.00% GC)
--------------
samples: 658
evals/sample: 1
Julia 中做 Array slicing 很容易，類似 python 的語法：
julia> a = collect(1:100);
julia> a[1:10]
10-element Array{Int64,1}:
1
2
3
4
5
6
7
8
9
10
語法容易使用很容易造成濫用， 導致性能問題， 因為：array slicing 會 copy 一個副本！ 我們來計算一個矩陣的每列的和， 簡單實現如下：
julia> function sum_vector(x::Array{Float64, 1})
s = 0.0
for i = 1:length(x)
s = s + x[i]
end
return s
end
sum_vector (generic function with 1 method)
julia> function sum_cols_matrix(x::Array{Float64, 2})
num_cols = size(x, 2)
s = zeros(num_cols)
for i = 1:num_cols
s[i] = sum_vector(x[:, i])
end
return s
end
sum_cols_matrix (generic function with 1 method)
julia> t = rand(1000, 1000);
julia> @benchmark sum_cols_matrix(t)
BenchmarkTools.Trial:
memory estimate: 7.76 MiB
allocs estimate: 1001
--------------
minimum time: 1.703 ms (0.00% GC)
median time: 2.600 ms (0.00% GC)
mean time: 2.902 ms (19.10% GC)
maximum time: 40.978 ms (94.27% GC)
--------------
samples: 1719
evals/sample: 1
x[:, j]是取第 j 列的所有元素。 算法很簡單， 定義一個函數sum_vector計算向量的和， 然后在sum_cols_matrix中取每一列傳過去。
由于x[:, j]這樣 slicing 會 copy 元素， 所以內存占用和分配都比較大。 Julia 提供了view函數，可以復用父數組的元素而不是 copy， 具體用法可以參考https://docs.julialang.org/en... 。
由于view返回的是SubArray類型， 所以我們先修改一下sum_vector的參數類型為AbstractArray：
julia> function sum_vector2(x::AbstractArray)
s = 0.0
for i = 1:length(x)
s = s + x[i]
end
return s
end
sum_vector2 (generic function with 1 method)
julia> function sum_cols_matrix_views(x::Array{Float64, 2})
num_cols = size(x, 2)
s = zeros(num_cols)
for i = 1:num_cols
s[i] = sum_vector2(view(x, :, i))
end
return s
end
sum_cols_matrix_views (generic function with 1 method)
julia>
julia> @benchmark sum_cols_matrix_views(t)
BenchmarkTools.Trial:
memory estimate: 7.94 KiB
allocs estimate: 1
--------------
minimum time: 812.209 μs (0.00% GC)
median time: 883.884 μs (0.00% GC)
mean time: 897.888 μs (0.00% GC)
maximum time: 6.552 ms (0.00% GC)
--------------
samples: 5474
evals/sample: 1
可以看到，提升還是蠻大的。
SIMD全稱Single Instruction Multiple Data， 是現代 CPU 的特性， 可以一條指令操作多條數據。 來加兩個向量試試：
julia> x = zeros(1_000_000); y = rand(1_000_000); z = rand(1_000_000);
julia> function sum_vectors!(x, y, z)
n = length(x)
for i = 1:n
x[i] = y[i] + z[i]
end
end
sum_vectors! (generic function with 1 method)
julia> function sum_vectors_inbounds!(x, y, z)
n = length(x)
@inbounds for i = 1:n
x[i] = y[i] + z[i]
end
end
sum_vectors_inbounds! (generic function with 1 method)
julia> function sum_vectors_inbounds_simd!(x, y, z)
n = length(x)
@inbounds @simd for i = 1:n
x[i] = y[i] + z[i]
end
end
sum_vectors_inbounds_simd! (generic function with 1 method)
julia> @benchmark sum_vectors!(x, y, z)
BenchmarkTools.Trial:
memory estimate: 0 bytes
allocs estimate: 0
--------------
minimum time: 1.117 ms (0.00% GC)
median time: 1.156 ms (0.00% GC)
mean time: 1.174 ms (0.00% GC)
maximum time: 1.977 ms (0.00% GC)
--------------
samples: 4234
evals/sample: 1
julia> @benchmark sum_vectors_inbounds!(x, y, z)
BenchmarkTools.Trial:
memory estimate: 0 bytes
allocs estimate: 0
--------------
minimum time: 1.118 ms (0.00% GC)
median time: 1.148 ms (0.00% GC)
mean time: 1.158 ms (0.00% GC)
maximum time: 1.960 ms (0.00% GC)
--------------
samples: 4294
evals/sample: 1
julia> @benchmark sum_vectors_inbounds_simd!(x, y, z)
BenchmarkTools.Trial:
memory estimate: 0 bytes
allocs estimate: 0
--------------
minimum time: 1.118 ms (0.00% GC)
median time: 1.145 ms (0.00% GC)
mean time: 1.155 ms (0.00% GC)
maximum time: 2.080 ms (0.00% GC)
--------------
samples: 4305
evals/sample: 1
這個測試結果， 無論用@inbounds還是@simd，都沒有提速（難道 Julia 現在默認會使用 SIMD？）， 測試了幾次都不行。 另外， 我在測試https://github.com/JuliaCompu... 的時候， 也發現有沒有 simd 差別不大， 如有知道原因的童鞋歡迎留言告知， 非常謝謝。
另外說Yeppp這個包也能提高速度， 還沒有測試。
如果我們設計的函數給別人用， 那么 API 會設計的比較通用（比如參數設計成 AbstractArray）， 可能會接受各種類型的 Array，這時候需要小心如何遍歷數組。
有兩種 indexing 方式。 一種是 linear indexing， 比如是一個三維 array， 每維度 10 個元素， 則可以用 x[1], x[2], ..., x[1000]連續地訪問數組。 第二種叫 cartesian indexing， 訪問方式為 x[i, j, k]。 某些數組不是連續的（比如 view 生成的 SubArray），這時候用 cartesian indexing 訪問會比用 linear indexing 訪問快， 因為 linear indexing 需要除法轉化成每一維的下標。 通用代碼一般會用 linear indexing， 這就有可能導致性能問題。
julia> function mysum_linear(a::AbstractArray)
s=zero(eltype(a))
for i = 1:length(a)
s=s+a[i]
end
return s
end
mysum_linear (generic function with 1 method)
julia> mysum_linear(1:1000000)
500000500000
julia> mysum_linear(reshape(1:1000000, 100, 100, 100))
500000500000
julia> mysum_linear(reshape(1:1000000, 1000, 1000))
500000500000
julia> mysum_linear(view(reshape(1:1000000, 1000, 1000), 1:500, 1:500) )
62437625000
julia> @benchmark mysum_linear(1:1000000)
BenchmarkTools.Trial:
memory estimate: 0 bytes
allocs estimate: 0
--------------
minimum time: 1.380 ns (0.00% GC)
median time: 1.426 ns (0.00% GC)
mean time: 1.537 ns (0.00% GC)
maximum time: 13.475 ns (0.00% GC)
--------------
samples: 10000
evals/sample: 1000
julia> @benchmark mysum_linear(reshape(1:1000000, 1000, 1000))
BenchmarkTools.Trial:
memory estimate: 0 bytes
allocs estimate: 0
--------------
minimum time: 1.379 ns (0.00% GC)
median time: 1.392 ns (0.00% GC)
mean time: 1.482 ns (0.00% GC)
maximum time: 23.089 ns (0.00% GC)
--------------
samples: 10000
evals/sample: 1000
julia> @benchmark mysum_linear(view(reshape(1:1000000, 1000, 1000), 1:500, 1:500))
BenchmarkTools.Trial:
memory estimate: 0 bytes
allocs estimate: 0
--------------
minimum time: 1.016 ms (0.00% GC)
median time: 1.047 ms (0.00% GC)
mean time: 1.071 ms (0.00% GC)
maximum time: 2.211 ms (0.00% GC)
--------------
samples: 4659
evals/sample: 1
可以看到最后一次測試， 元素總數更少了， 但是時間更長， 原因就是數組不是連續的， 但是又用了 linear indexing。 最簡單的解決方法是直接迭代元素， 而不是迭代下標。
julia> function mysum_in(a::AbstractArray)
s = zero(eltype(a))
for i in a
s=s+ i
end
end
mysum_in (generic function with 1 method)
julia> @benchmark mysum_in(view(reshape(1:1000000, 1000, 1000), 1:500, 1:500))
BenchmarkTools.Trial:
memory estimate: 0 bytes
allocs estimate: 0
--------------
minimum time: 224.538 μs (0.00% GC)
median time: 238.493 μs (0.00% GC)
mean time: 246.847 μs (0.00% GC)
maximum time: 413.477 μs (0.00% GC)
--------------
samples: 10000
evals/sample: 1
可以看到， 跟 linear indexing 相比， 效率是 4 倍多。 但是如果我們需要下標怎么辦呢？可以用eachindex()方法，每一種 array 都會定義一個優化過的eachindex()方法：
julia> function mysum_eachindex(a::AbstractArray)
s = zero(eltype(a))
for i in eachindex(a)
s = s + a[i]
end
end
mysum_eachindex (generic function with 1 method)
julia> @benchmark mysum_eachindex(view(reshape(1:1000000, 1000, 1000), 1:500, 1:500))
BenchmarkTools.Trial:
memory estimate: 0 bytes
allocs estimate: 0
--------------
minimum time: 243.295 μs (0.00% GC)
median time: 273.168 μs (0.00% GC)
mean time: 285.002 μs (0.00% GC)
maximum time: 4.191 ms (0.00% GC)
--------------
samples: 10000
evals/sample: 1
通過這兩篇文章介紹， 我們基本上掌握了 Julia 高性能的方法。 如果還不夠， 那就只能求助于并行和分布式了， 等著我們下一篇介紹吧。
歡迎加入知識星球一起分享討論有趣的技術話題。
文章版權歸作者所有，未經允許請勿轉載,若此文章存在違規行為，您可以聯系管理員刪除。
轉載請注明本文地址：http://specialneedsforspecialkids.com/yun/19846.html
摘要：前面兩篇讓代碼飛起來高性能學習筆記一讓代碼飛起來高性能學習筆記二，介紹了如何寫出高性能的代碼，這篇結合我最近的項目，簡單測試對比一下各種語言用算法計算的效率。下一篇，我們就來看一下中如何利用并行進一步提高效率。 前面兩篇讓代碼飛起來——高性能 Julia 學習筆記（一） 讓代碼飛起來——高性能 Julia 學習筆記（二）， 介紹了如何寫出高性能的 Julia 代碼， 這篇結合我最近的項...
摘要：精準截圖不再需要細調截圖框這一細節功能真心值無數個贊相比大多數截屏軟件只能檢測整個應用窗口邊界，對界面元素的判定讓你操作時可以更加精準快捷，下面的動圖就可以讓你直觀地感受到這個功能的強大之處。接著打開截屏的現象，將里面的顯示邊框寬度調整為。 ...
摘要：介紹性能號稱可以趕得上，我很好奇的執行速度，因為我一直用的是，所以就想把和做一下簡單的比較?？偨Y從上面的比較來看，確實比快很多，不過這里只做了簡單的比較，并沒有做嚴謹的測試，僅供參考。 1、介紹 Julia性能號稱可以趕得上c/c++，我很好奇Julia的執行速度，因為我一直用的是Java，所以就想把Julia和Java做一下簡單的比較。這次比較一下Julia和Java做一億次加法運算...
摘要：本文匯總了中數款儀表盤工具，簡單比較其使用場景學習難度，挑選趁手的來玩即可。進一步學習誰是中最強開發工具誰是中最強開發工具 儀表盤 (Dashboard)，可簡單的理解為一個交互式網頁，在其中，用戶可以不懂代碼，拖拖拽拽即可與數據交互、做數據探索建模分析、展示自己關注的結果。 ...
閱讀 3574·2019-08-30 15:55
閱讀 1372·2019-08-29 16:20
閱讀 3656·2019-08-29 12:42
閱讀 2660·2019-08-26 10:35
閱讀 1010·2019-08-26 10:23
閱讀 3405·2019-08-23 18:32
閱讀 897·2019-08-23 18:32
閱讀 2891·2019-08-23 14:55
