パフォーマンス Tips

以降の節では、Julia コードを可能な限り高速にするために役立つテクニックを簡単に紹介します。

グローバル変数を避ける

グローバル変数の値は好きなタイミングで変更できます。これは、グローバル変数の型が任意の時点で変わり得ることを意味します。このためコンパイラはグローバル変数が含まれるコードを上手く最適化できません。可能な場合は必ず、変数をローカルにするか、関数の引数として渡してください。

パフォーマンスが重要なコード、あるいはベンチマークされるコードは必ず関数の内部にあるべきです。

グローバル変数は定数であることも多いはずです。その場合は定数として宣言すると性能が大きく向上します:

const DEFAULT_VAL = 0

定数でないグローバル変数であっても、それを使うコードに型注釈を付ければ最適化が可能です:

global x = rand(1000)

function loop_over_global()
    s = 0.0
    for i in x::Vector{Float64}
        s += i
    end
    return s
end

この例では x を関数に引数として渡すスタイルの方が優れています。x を引数にすればコードの再利用性が高まり、入力と出力が明確になります。

情報

REPL に入力される全てのコードはグローバルスコープで評価されるので、トップレベルで定義・代入された変数はグローバル変数となります。モジュールのトップレベルスコープで宣言される変数もグローバル変数です。

また REPL への入力

julia> x = 1.0

は、次の等価です:

julia> global x = 1.0

このため、上述の性能劣化の問題は REPL でも生じます。

@time で性能を計測し、メモリのアロケートに注目する

性能を計測する便利なツールが @time マクロです。例として、上述のグローバル変数を使ったコードから型注釈を取り除いたものを考えます:

julia> x = rand(1000);

julia> function sum_global()
           s = 0.0
           for i in x
               s += i
           end
           return s
       end;

julia> @time sum_global()
  0.017705 seconds (15.28 k allocations: 694.484 KiB)
496.84883432553846

julia> @time sum_global()
  0.000140 seconds (3.49 k allocations: 70.313 KiB)
496.84883432553846

sum_global 関数は最初の呼び出し @time sum_global() でコンパイルされます。加えて、このセッションで @time を使ったことがなければ時間計測に必要な関数もコンパイルされるので、初回の実行の計測結果を信用してはいけません。

信用できる二回目の計測において、大量のメモリがアロケートされていることに注目してください。ここでは 64 ビット浮動小数点数のベクトルの和を計算しているだけなので、@time が計測するヒープメモリのアロケートは必要ないはずです。

予期せぬメモリのアロケートはほぼ間違いなくコードに問題があるサインです。たいていは不安定な (特定されていない) 型を使ったり、一時的な小さい配列を大量に作成することによって引き起こされます。このような場合、コードが表す処理がアロケートを含んでいるのに加えて、関数に対して生成されるコードが最適からほど遠いものになる可能性が非常に高くなります。こういったサインを決して見逃さず、以下のアドバイスに従ってください。

x を関数の引数として渡すようにすればメモリのアロケートが起こらなくなり、実行速度が格段に向上します:

julia> x = rand(1000);

julia> function sum_arg(x)
           s = 0.0
           for i in x
               s += i
           end
           return s
       end;

julia> @time sum_arg(x)
  0.007701 seconds (821 allocations: 43.059 KiB)
496.84883432553846

julia> @time sum_arg(x)
  0.000006 seconds (5 allocations: 176 bytes)
496.84883432553846

5 回のアロケートが計測されているのは、@time マクロをグローバルスコープで実行しているためです。計測を関数として行えば、アロケートが起こっていないことを確認できます:

julia> time_sum(x) = @time sum_arg(x);

julia> time_sum(x)
  0.000001 seconds
496.84883432553846

関数が処理の途中でメモリをアロケートしなければならない場合は性能の測定をこの例ほど単純には行えません。後述するツールを使って問題を診断するか、メモリのアロケートをアルゴリズム的な部分と切り離したバージョンの関数を作成してください (参照: 出力を事前にアロケートする)。

情報

さらに本格的なベンチマークを行うパッケージに BenchmarkTools.jl があります。このパッケージは関数を複数回評価してノイズを低減するといった処理を行います。

ツールを利用する

コードが持つ問題を診断し、コードの改善と性能向上を可能にするツールが Julia とそのパッケージエコシステムによって提供されます。

コンテナの型パラメータに抽象型を使わない

パラメータを持つ型 (例えば配列) を使うときは、パラメータに抽象型を使うのは避けた方が賢明です。

次のコードを考えます:

julia> a = Real[]
Real[]

julia> push!(a, 1); push!(a, 2.0); push!(a, π)
3-element Array{Real,1}:
 1
 2.0
 π = 3.1415926535897...

抽象型 Real の配列である a は任意の Real 型の値を保持できなければなりません。Real オブジェクトのサイズと構造に制限は無いので、a は個別にアロケートされた Real オブジェクトへのポインタの配列として表される必要があります。これに対して、同一の型 (例えば Float64) の数値だけを a に入れられるようにすれば、オブジェクトを効率良く格納できます:

julia> a = Float64[]
Float64[]

julia> push!(a, 1); push!(a, 2.0); push!(a,  π)
3-element Array{Float64,1}:
 1.0
 2.0
 3.141592653589793

aFloat64 の配列と宣言すると、a の要素に代入される数値は Float64 に変換されます。a は 64 ビット浮動小数点数が隙間無く並んだブロックに要素を格納でき、配列操作の効率が向上します。

関連する議論はパラメトリック型の節にもあります。

型宣言を (必要なら) 付ける

型宣言を省略できる多くの言語において、型宣言の追加はコードを高速化するための第一歩です。しかし Julia ではそうではありません。Julia では、基本的にコンパイラは全ての関数引数・ローカル変数・式の型を知っています。しかし、型宣言でコードが高速化する状況も一部存在します。

抽象型を持つフィールドを使わない

Julia ではフィールドの型を指定しない型宣言が可能です:

julia> struct MyAmbiguousType
           a
       end

こうすると a は任意の型になれます。この機能が役立つ場面も多くありますが、これには一つ欠点があります: MyAmbiguousType 型のオブジェクトが絡む処理では、コンパイラは高性能なコードを生成できません。なぜなら、コンパイラはコードを組み立てるときにオブジェクトの (値ではなく) 型を使うためです。残念なことに、MyAmbiguousType 型のオブジェクトから得られる情報はほとんどありません:

julia> b = MyAmbiguousType("Hello")
MyAmbiguousType("Hello")

julia> c = MyAmbiguousType(17)
MyAmbiguousType(17)

julia> typeof(b)
MyAmbiguousType

julia> typeof(c)
MyAmbiguousType

bc は同じ型を持ちますが、メモリ上での内部表現は大きく異なります。a に格納したのが数値であったとしても、UInt8Float64 の表現は異なるので、CPU は異なる種類の命令を使って別々に処理を行わなければなりません。行うべき処理を判断するための情報が型から取得できないので、その判断は実行時に回されます。これが性能を低下させます。

a の型を宣言すれば性能が向上しますが、ここでは a が複数の型になり得る状況をさらに考えます。このような状況では型パラメータを使うのが自然な解決法です。例えば次のようにします:

julia> mutable struct MyType{T<:AbstractFloat}
           a::T
       end

これは、次の構造体より優れています:

julia> mutable struct MyStillAmbiguousType
           a::AbstractFloat
       end

なぜなら、最初のバージョンでは a を包むオブジェクトの型が a の型を特定するからです。このことが分かる例を示します:

julia> m = MyType(3.2)
MyType{Float64}(3.2)

julia> t = MyStillAmbiguousType(3.2)
MyStillAmbiguousType(3.2)

julia> typeof(m)
MyType{Float64}

julia> typeof(t)
MyStillAmbiguousType

m の型を見れば m.a の型が分かりますが、t の型を見ても t.a の型は分かりません。実を言えば、t.a の型は変更できます:

julia> typeof(t.a)
Float64

julia> t.a = 4.5f0
4.5f0

julia> typeof(t.a)
Float32

これに対して、m.a の型は m が構築された後に変更できません:

julia> m.a = 4.5f0
4.5f0

julia> typeof(m.a)
Float64

m.a の型が m の型から分かるという事実 (および m.a の型が関数の中で変わらないという事実) があるので、コンパイラは m を使う処理に対して高度に最適化されたコードを生成できます。一方で t のようなオブジェクトに対してこれは不可能です。

もちろん、以上の議論は具象型を使って m を構築したときにだけ成立します。抽象型を使ったときは仮定が成り立たず、例えば m.a の型は変更できます:

julia> m = MyType{AbstractFloat}(3.2)
MyType{AbstractFloat}(3.2)

julia> typeof(m.a)
Float64

julia> m.a = 4.5f0
4.5f0

julia> typeof(m.a)
Float32

実際の処理において、MyType{AbstractFloat} のオブジェクトは MyStillAmbiguousType のオブジェクトと全く同一に振る舞います。

抽象型を使ったときに次の簡単なコードに対して生成される膨大なコードを確認しておくと理解に役立つでしょう:

func(m::MyType) = m.a+1

生成されるコードは次の関数で表示できます:

code_llvm(func, Tuple{MyType{Float64}})
code_llvm(func, Tuple{MyType{AbstractFloat}})

返り値は非常に長いのでここに載せませんが、自分で確認してみるべきです。型が完全に特定される一つ目の場合には実行時に型を解決するコードが必要ないので、コンパイラは短くて高速なコードを生成できます。

フィールドに抽象コンテナを使わない

同じことはコンテナ型に対しても言えます:

julia> struct MySimpleContainer{A<:AbstractVector}
           a::A
       end

julia> struct MyAmbiguousContainer{T}
           a::AbstractVector{T}
       end

構築されるオブジェクトの型を確認する例を示します:

julia> c = MySimpleContainer(1:3);

julia> typeof(c)
MySimpleContainer{UnitRange{Int64}}

julia> c = MySimpleContainer([1:3;]);

julia> typeof(c)
MySimpleContainer{Array{Int64,1}}

julia> b = MyAmbiguousContainer(1:3);

julia> typeof(b)
MyAmbiguousContainer{Int64}

julia> b = MyAmbiguousContainer([1:3;]);

julia> typeof(b)
MyAmbiguousContainer{Int64}

MySimpleContainer では型とそのパラメータが内部のオブジェクトを完全に特定するので、コンパイラは最適化された関数を生成できます。多くのケースで、こうしておけば十分なはずです。

特定される型を使えばコンパイラは問題なく仕事を行えますが、a要素型に応じて異なる処理を行いたい場合もあるでしょう。次の例における foo のように、要素型によって変わる処理を一つの関数を使ってラップするのが通常は一番の方法です:

julia> function sumfoo(c::MySimpleContainer)
           s = 0
           for x in c.a
               s += foo(x)
           end
           s
       end
sumfoo (generic function with 1 method)

julia> foo(x::Integer) = x
foo (generic function with 1 method)

julia> foo(x::AbstractFloat) = round(x)
foo (generic function with 2 methods)

こうすればコードは綺麗に保たれ、さらにコンパイラは全ての場合で最適化されたコードを生成できます。

しかし、ときには MySimpleContainer 型のフィールド aAbstractVector 型やその要素型に応じて異なるバージョンの関数を宣言しなければならない場合もあります。そのときは次のようにできます:

julia> function myfunc(c::MySimpleContainer{<:AbstractArray{<:Integer}})
           return c.a[1]+1
       end
myfunc (generic function with 1 method)

julia> function myfunc(c::MySimpleContainer{<:AbstractArray{<:AbstractFloat}})
           return c.a[1]+2
       end
myfunc (generic function with 2 methods)

julia> function myfunc(c::MySimpleContainer{Vector{T}}) where T <: Integer
           return c.a[1]+3
       end
myfunc (generic function with 3 methods)

実行例を示します:

julia> myfunc(MySimpleContainer(1:3))
2

julia> myfunc(MySimpleContainer(1.0:3))
3.0

julia> myfunc(MySimpleContainer([1:3;]))
4

型の無いコンテナから取得した値には型注釈を付ける

任意の型を格納できるデータ構造 (Array{Any} 型の配列など) が上手く利用できる状況もよくあります。しかし、そういったデータ構造を使うときに特定の要素の型が (たまたま) 分かっているなら、その情報をコンパイラに与えることでコードの最適化が可能です:

function foo(a::Array{Any,1})
    x = a[1]::Int32
    b = x+1
    ...
end

このコードは a の第一要素が Int32 で決まっていることを想定しています。こういった型注釈を付けると型が異なる場合に実行時エラーとなるので、型に関するバグを速く見つけられるという利点もあります。

a[1] の型は正確に分からないが、どんな型であれ Int32 に変換してから使う」ことが分かっているなら、x = convert(Int32, a[1])::Int32 という宣言が利用できます。convert 関数を使うことで Int32 に変換可能な任意のオブジェクト (例えば UInt8) をa[1] に使えるようになります。つまり、型の要件を緩めることでコードの汎用性が高められます。この文脈では convert にも型注釈を付けないと型が安定にならないことに注意してください。この理由は、コンパイラが関数の返り値の型を推論するには関数の全ての引数の型が分かっていなければならないためです (convert が特別扱いされることはありません)。

型が実行時に構築されるときは、型注釈を付けても性能は向上しません (逆に低下します)。型注釈があってもコンパイラは以降のコードを特殊化できず、型の確認に時間をとられるだけだからです。例えば

function nr(a, prec)
    ctype = prec == 32 ? Float32 : Float64
    b = Complex{ctype}(a)
    c = (b + 1.0f0)::Complex{ctype}
    abs(c)
end

において、c に付いている型注釈 ::Complex{ctype} は性能に悪影響を及ぼします。実行時に構築される型が絡むコードを最適化するには、後述する関数バリアと呼ばれるテクニックを使ってください。これは構築される型の値を引数に受け取るカーネル関数を別に用意して、どの型が構築されてもコンパイラによって適切に特殊化されたカーネル関数が実行されるようにするというテクニックです。

例えば上の関数では、b が構築されたらすぐにカーネル関数 k に渡すようにすれば関数バリアを張れます。関数 k の引数 b を型パラメータ T を使って Complex{T} 型と宣言すれば、k に次の形をした宣言を入れても性能が悪化しません (ただし向上もしません):

c = (b + 1.0f0)::Complex{T}

コンパイラは k をコンパイルするときに c の型を決定できるためです。

Julia が特殊化を行わない状況に気を付ける

Julia はヒューリスティックとして、関数の引数が Type, Function, Vararg のときは型パラメータに関するメソッドの自動的な特殊化を行いません。また Julia はメソッド内で引数が使われるときは特殊化しますが、引数が他の関数に渡されているだけなら特殊化を行いません。通常こうしても実行時の性能は低下せず、コンパイラの性能が向上します。この仕様が実行時の性能に悪影響を及ぼしていることが判明したら、メソッドの宣言に型パラメータを追加すれば特殊化できます。いくつか例を示します:

次の関数は特殊化されません:

function f_type(t)  # あるいは t::Type
    x = ones(t, 10)
    return sum(map(sin, x))
end

次の関数は特殊化されます:

function g_type(t::Type{T}) where T
    x = ones(T, 10)
    return sum(map(sin, x))
end

次の関数は特殊化されません:

f_func(f, num) = ntuple(f, div(num, 2))
g_func(g::Function, num) = ntuple(g, div(num, 2))

次の関数は特殊化されます:

h_func(h::H, num) where {H} = ntuple(h, div(num, 2))

次の関数は特殊化されません:

f_vararg(x::Int...) = tuple(x...)

次の関数は特殊化されます:

g_vararg(x::Vararg{Int, N}) where {N} = tuple(x...)

特殊化を起こすには型パラメータを一つ導入すれば十分であり、そのとき他の型が制約されていなくても構いません。例えば次の関数は特殊化されます。引数が全て同じ型でないときに使えます:

h_vararg(x::Vararg{Any, N}) where {N} = tuple(x...)

@code_typed 関数とその仲間を使えば、Julia が通常なら特殊化しないメソッド呼び出しに対しても特殊化したコードを確認できます。引数を変えたときに特殊化されたコードが生成されるかどうかを確認するには、メソッドの内部情報の確認が必要です。具体的には (@which f(...)).specializations が考えている引数に対する特殊化を含むことを確認してください。

関数を複数の定義に分割する

小さなメソッドの集合体として関数を書けば、コンパイラは最も適したコードを直接呼び出すことができ、インライン化も可能になります。

次に示すのは望ましくない “混ざった関数” の例です。この関数は複数の定義を使って書くべきです:

using LinearAlgebra

function mynorm(A)
    if isa(A, Vector)
        return sqrt(real(dot(A,A)))
    elseif isa(A, Matrix)
        return maximum(svdvals(A))
    else
        error("mynorm: invalid argument")
    end
end

次のように書けば明快かつ効率的になります:

norm(x::Vector) = sqrt(real(dot(x, x)))
norm(A::Matrix) = maximum(svdvals(A))

ただし、コンパイラは mynorm にあるような死んでいる枝 (到達不可能なブロック) を最適化して取り除くことに非常に長けています。

「型安定」 な関数を書く

可能なら、関数は常に同じ型の値を返すようにするべきです。例えば次の関数定義を考えます:

pos(x) = x < 0 ? 0 : x

何の変哲もないコードに見えるかもしれませんが、x は任意の型であり得るのに対して 0 が (Int 型の) 整数であることが問題になります。つまり、x の値によってはこの関数は二つ以上の型の値を返します。Julia でこの振る舞いは許されており、必要な場合もあるでしょう。ただ今の例では簡単に修正できます:

pos(x) = x < 0 ? zero(x) : x

同様の oneunit 関数もあります。また、より一般的な oftype(x, y) 関数を使えば、x の型に変換した y を取得できます。

変数の型を変更しない

型の安定性の問題は関数の中で何度も使われる変数でも発生します:

function foo()
    x = 1
    for i = 1:10
        x /= rand()
    end
    return x
end

最初ローカル変数 x は整数であり、反復を一度実行すると (/ 演算子の返り値の) 浮動小数点数になります。このためコンパイラはループの本体を上手く最適化できません。考えられる修正方法はいくつかあります:

カーネル関数を分離する (関数バリア)

多くの関数は「最初にセットアップ処理を行い、それから中心的な処理を何度も反復して行う」というパターンを持ちます。可能なら中心的な処理を個別の関数とするのが良いアイデアです。例えば、次に示すのはランダムに選択した型の配列を返す (わざとらしい) 関数です:

julia> function strange_twos(n)
           a = Vector{rand(Bool) ? Int64 : Float64}(undef, n)
           for i = 1:n
               a[i] = 2
           end
           return a
       end;

julia> strange_twos(3)
3-element Array{Float64,1}:
 2.0
 2.0
 2.0

この関数は次のように書くべきです:

julia> function fill_twos!(a)
           for i = eachindex(a)
               a[i] = 2
           end
       end;

julia> function strange_twos(n)
           a = Vector{rand(Bool) ? Int64 : Float64}(undef, n)
           fill_twos!(a)
           return a
       end;

julia> strange_twos(3)
3-element Array{Float64,1}:
 2.0
 2.0
 2.0

Julia のコンパイラは関数を境界として引数の型に関してコードを特殊化するので、最初の実装をコンパイルするときにループ中の a の型を知ることはできません (a の型はランダムに選ばれるためです)。二番目の実装ではループ中の処理を fill_twos! 関数として a の型を知った状態でコンパイルできるので、こちらの実装の方が一般に高速です。

二番目の形式はスタイル的にも優れており、再利用しやすいコードが生まれやすくなります。

このパターンは Julia Base でも使われています。例えば abstractarray.jlvcat, hcat, fill! がそうです (fill_tows! を自分で書く代わりに fill! を使うこともできたでしょう)。

strange_twos のような関数は不確かな型を持つデータを扱っているときによく書かれます。例えば入力ファイルから読み込んだデータに整数・浮動小数点数・文字列といったデータが混ざっている場合です。

値を型パラメータに持つ型の扱いに気を付ける

各次元のサイズが 3 の N 次元配列を作りたいとします。例えば、次のようにすればそういった配列を作成できます:

julia> A = fill(5.0, (3, 3))
3×3 Array{Float64,2}:
 5.0  5.0  5.0
 5.0  5.0  5.0
 5.0  5.0  5.0

コンパイラは AArray{Float64,2} 型だと判断できるので、このアプローチに問題はありません。コンパイラが A の型を判断できる理由は、fill する値の型 5.0::Float64 と次元 (3, 3)::NTuple{2,Int} を知っているからです。そのため、この後に同じ関数で A を使うと、コンパイラは非常に最適化されたコードを生成できます。

では、引数によって定まる次元 3×3×... を持つ配列を作成したいとしましょう。次の関数を書きたくなるかもしれません:

julia> function array3(fillval, N)
           fill(fillval, ntuple(d->3, N))
       end
array3 (generic function with 1 method)

julia> array3(5.0, 2)
3×3 Array{Float64,2}:
 5.0  5.0  5.0
 5.0  5.0  5.0
 5.0  5.0  5.0

これは正しい値を計算しますが、(@code_warntype array3(5.0, 2) を実行すれば分かるように) 出力の型を推論できないことが問題です: 引数 NInt 型のであり、コンパイル時の型推論では N の値は与えられません (そもそも実行時にならないと判明しません)。これは array3 関数の出力を使うコードが保守的ならざるを得ず、A に対するアクセスのたびに型を確認しなければならないことを意味します。そういったコードは非常に低速になります。

もちろん、こういった問題を解決する非常に優れた手段の一つは関数バリアのテクニックを使うことです。しかし、型の不安定性を完全に取り除きたい場合もあるはずです。そのような場合には、Val{T}() を使って次元をパラメータとして渡すのが一つのアプローチです (参照: "値型"):

julia> function array3(fillval, ::Val{N}) where N
           fill(fillval, ntuple(d->3, Val(N)))
       end
array3 (generic function with 1 method)

julia> array3(5.0, Val(2))
3×3 Array{Float64,2}:
 5.0  5.0  5.0
 5.0  5.0  5.0
 5.0  5.0  5.0

Julia には第二引数に Val{::Int} を受け取る特殊なバージョンの ntuple が存在します。N を型パラメータとして渡すことで、その “値” をコンパイラに伝えることができます。同様に、この array3 ではコンパイラが返り値の型を予測できます。

ただし、このテクニックを正しく使うには驚くほど細かな制御が必要になります。例えば、array3 を次のような関数から呼び出すのは正しくありません:

function call_array3(fillval, n)
    A = array3(fillval, Val(n))
end

こうすると、最初の問題がまた表れます: コンパイラは n の値を推測できず、Val(n)も分からないのです。Val を使おうとして使い方を間違えると、多くの場合で性能は簡単に悪化します (このようなコードを使うべき唯一の状況は、Val と関数バリアを効率的に組み合わせるときです)。

Val の正しい使い方の例を示します:

function filter3(A::AbstractArray{T,N}) where {T,N}
    kernel = array3(1, Val(N))
    filter(A, kernel)
end

この例では N がメソッドのパラメータとして渡されるので、N の値はコンパイラに伝わります。Val{T} が本質的な意味を持つのは、TVal(3) のようなハードコードされた値/リテラルのとき、あるいは型のドメインで既に指定されているときだけです。

多重ディスパッチの危険性を認識する

初心者が多重ディスパッチの使い方を学んで理解すると、(気持ちは分かるのですが) あらゆることに多重ディスパッチを使いまくる傾向が見られます。例えば、型パラメータに情報を格納するといった使い方です:

struct Car{Make, Model}
    year::Int
    ...その他のフィールド...
end

こうした上で Car{:Honda,:Accord}(year, args...) などとすればオブジェクトに対するディスパッチを行えます。

次のいずれかが成り立つなら、このやり方に意味があると言えます:

二つ目の条件が成り立つなら、同種の Car を収めた配列を処理する関数は生産的に特殊化できます: 配列の各要素が同じ具象型なので Julia は各要素の型を事前に知ることができ、関数をコンパイルするときに呼び出すべきメソッドを “検索” でき (実行時の型の確認は不要となり)、配列全体を処理する効率的なコードを生成できます。

もし二つの条件が成り立たないなら、このやり方で得られるものは無い可能性が高くなります: さらに悪いことに、型パラメータを追加することで起こる「組み合わせ論的な型の爆発」は生産性を害するでしょう。items[i+1]item[i] と異なる型を持っていると、Julia は実行時に型を取得し、適切なメソッドをメソッドテーブルから検索し、(型の交差を計算することで) 適合するのがどれかを判断し、JIT コンパイルされているかを調べ、JIT コンパイルされていなければコンパイルし、そして実際の呼び出しを行わなくてはなりません。つまり、これは型システムと JIT コンパイル機構に向かってコードに対する switch 文あるいは辞書の探索を行えと言っているのに等しいのです。

型のディスパッチ・辞書の探索・switch 文を比べた実行時ベンチマークがメーリングリストにあります。

おそらく実行時間への影響よりもさらに深刻なのが、コンパイル時間への影響です。Julia は異なる Car{Make, Model} のそれぞれに対して特殊化された関数をコンパイルするので、型の個数が数百あるいは数千に達するとき、オブジェクトをパラメータに受け取る全ての関数 (オブジェクト固有の get_year のような関数から、Julia Base に含まれる push! のような汎用な関数まで) が数百回あるいは数千回コンパイルされることになります。それぞれがコンパイル済みコードのキャッシュのサイズを増加させ、内部のメソッドリストを長くし、他にも様々な悪影響を及ぼします。「型パラメータに値を使って高速化」を突き詰めすぎると、途方もない量のリソースがすぐに必要になってしまいます。

配列にはメモリに沿った列優先の順序でアクセスする

Julia の多次元配列は各要素を列優先 (column-major) の順序で格納します。これはメモリ上に配列の各列が順に並ぶことを意味します。次に示すように、[:] 構文 (あるいは vec 関数) を使えばこの事実を確認できます。配列の要素が [1 2 3 4] ではなく [1 3 2 4] と並んでいることに注目してください:

julia> x = [1 2; 3 4]
2×2 Array{Int64,2}:
 1  2
 3  4

julia> x[:]
4-element Array{Int64,1}:
 1
 3
 2
 4

配列の順序に関するこの慣習は Fortran, MATLAB, R など多くの言語で使われています。列優先でない順序として使われるのは行優先 (column-major) の順序であり、これは C や Python (numpy) で採用されています。

配列に対して反復を行うとき配列の格納順序を意識しないと性能が大幅に低下することがあります。最も重要なこととして、一般に列優先の順序を使う配列では最初の添え字が一番速く変わることを覚えておいてください。これは Julia 配列のループ処理では一番内側のループ変数をスライスの最初に使うと高速になることを意味します。さらに、配列に対する : を使った添え字アクセスは特定の次元に含まれる全ての要素に対してアクセスする暗黙なループなのだと意識してください。例えば行を取り出すよりも列を取り出す方が高速になる可能性があります。

次の人工的な例を考えます。 Vector を受け取って、それを行または列に沿ってコピーした正方 Matrix を作る関数を書きたいとして、行と列のどちらにコピーするかは (他のコードは簡単に対応できるなどの理由で) 重要でないとします。この処理は少なくとも四つの方法で書けるはずです (実際のコードでは repeat を推奨します):

function copy_cols(x::Vector{T}) where T
    inds = axes(x, 1)
    out = similar(Array{T}, inds, inds)
    for i = inds
        out[:, i] = x
    end
    return out
end

function copy_rows(x::Vector{T}) where T
    inds = axes(x, 1)
    out = similar(Array{T}, inds, inds)
    for i = inds
        out[i, :] = x
    end
    return out
end

function copy_col_row(x::Vector{T}) where T
    inds = axes(x, 1)
    out = similar(Array{T}, inds, inds)
    for col = inds, row = inds
        out[row, col] = x[row]
    end
    return out
end

function copy_row_col(x::Vector{T}) where T
    inds = axes(x, 1)
    out = similar(Array{T}, inds, inds)
    for row = inds, col = inds
        out[row, col] = x[col]
    end
    return out
end

同じ 10000×1 の乱数ベクトルを入力として、この四つの関数の実行時間を計測します:

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_colscopy_rows よりずっと速いことが分かります。これは copy_colsMatrix のメモリが持つ列優先のレイアウトに沿って一列ごとにデータをコピーしているためであり、不思議なことではありません。また copy_col_rowcopy_row_col よりずっと速いのは、スライスの最初の添え字を一番内側のループ変数とするべきであるという一般的なアドバイスに従っているためです。

出力配列を事前にアロケートする

関数が Array や複雑な型を返す場合は、メモリをアロケートする必要があります。残念なことに、メモリのアロケートとその逆操作であるガベージコレクションは大きなボトルネックとなることが少なくありません。

出力を事前にアロケートすれば、関数内でアロケートを回避できる可能性があります。簡単な例として、次のコードを考えます:

julia> function xinc(x)
           return [x, x+1, x+2]
       end;

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

これは次のように書き換えられます:

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

julia> function loopinc_prealloc()
           ret = Vector{Int}(undef, 3)
           y = 0
           for i = 1:10^7
               xinc!(ret, i)
               y += ret[2]
           end
           return y
       end;

実行時間は次の通りです:

julia> @time loopinc()
  0.529894 seconds (40.00 M allocations: 1.490 GiB, 12.14% gc time)
50000015000000

julia> @time loopinc_prealloc()
  0.030850 seconds (6 allocations: 288 bytes)
50000015000000

事前のアロケートには他の利点もあります。一つはアルゴリズムが “出力” する型を呼び出し側が制御できることです。例えば上の xinc! では、必要なら Array ではなく SubArray を渡すことができます。

この考え方を過激に適用してあらゆるものを事前にアロケートしようとすると、コードは汚くなります。そのため性能の計測といくらかの割り切りが必要になるでしょう。ただしベクトル化された (要素ごとの) 関数に対しては、便利な構文 x .= f.(y) を使うことで融合されたループによるインプレースな処理を一時的な配列を使わずに呼び出せます (参照: 関数をベクトル化するドット構文)。

もっとドット: ベクトル化演算を融合する

Julia では特殊なドット構文を使って任意のスカラー関数を「ベクトル化された」関数として呼び出すことができ、同様に任意の演算子にドットを付ければ「ベクトル化された」演算子として呼び出せます。そのときネストされた「ドット呼び出し」には融合 (fusing) という特別な処理が行われます: 入れ子になったドット呼び出しは構文レベルで一つのループへと変換され、一時的な配列をアロケートしなくなります。.= などのベクトル化された代入演算子を使えば、計算結果は事前アロケートされた配列へインプレースに格納されます (上述)。

線形代数の文脈では、これは vector + vectorvector * scalar は定義されるものの、vector .+ vectorvector .* scalar を使う方が望ましいことを意味します。なぜなら最終的なループが他の計算と融合される可能性があるからです。例えば次の二つの関数を考えます:

julia> f(x) = 3x.^2 + 4x + 7x.^3;

julia> fdot(x) = @. 3x^2 + 4x + 7x^3; # 3 .* x.^2 .+ 4 .* x .+ 7 .* x.^3 と等価

ffdot が計算する値は同じです (fdot@. マクロを使って定義されています)。しかし配列に適用すると、fdot の方がはるかに高速になります:

julia> x = rand(10^6);

julia> @time f(x);
  0.019049 seconds (16 allocations: 45.777 MiB, 18.59% gc time)

julia> @time fdot(x);
  0.002790 seconds (6 allocations: 7.630 MiB)

julia> @time f.(x);
  0.002626 seconds (8 allocations: 7.630 MiB)

この測定結果からは、fdot(x)f(x) より十倍高速で、六分の一のメモリしかアロケートしないことが分かります。これは f(x) に含まれる *+ が一時的な配列を新しくアロケートし、異なるループとして処理を実行するためです。もちろん、この例では f.(x) とすれば fdot(x) と同程度の速度が得られますが、多くの状況においてドットを式の中に入れる方がベクトル化のために個別の関数を定義するよりも優れています。

スライスにビューを使えないか考える

Julia では、array[1:5, :] のような配列の「スライス」式がデータのコピーを作成します (代入の左辺にある array[1:5, :] = ... のようなスライスは例外で、このときは切り出された領域へインプレースな代入が起こります)。スライスに対してたくさんの処理を行うなら、データをコピーによって処理対象のデータ領域がメモリ上で小さな連続領域となるので、元の配列へ添え字アクセスを行うより効率が向上します。一方で、スライスに対して少量の処理を行うだけなら、アロケートとコピーのコストが割に合わない可能性があります。

このような場合のもう一つの選択肢が配列の「ビュー (view)」です。ビューは SubArray 型の配列として振る舞うオブジェクトですが、コピーを作成することなくアクセスのたびに元の配列のデータをその場で参照します (例えばビューへ書き込めば元の配列のデータも変更されます)。ビューは個別のスライスに対して view を呼ぶか、スライスを使う式あるいはブロックの前に @views マクロを付けることで利用できます。例を示します:

julia> fcopy(x) = sum(x[2:end-1]);

julia> @views fview(x) = sum(x[2:end-1]);

julia> x = rand(10^6);

julia> @time fcopy(x);
  0.003051 seconds (7 allocations: 7.630 MB)

julia> @time fview(x);
  0.001020 seconds (6 allocations: 224 bytes)

fview を付けたバージョンでは三倍の高速化に加えて消費メモリの削減が達成されています。

「データのコピーは悪い」と決めつけない

配列はメモリ上の連続領域に格納されるので、規則的なアクセスに対しては CPU によるベクトル化とキャッシュによるメモリアクセスの隠蔽が可能です。配列へのアクセスを列優先にすべき (上述) なのも同じ理由です。逆に規則的でないアクセスパターンや連続でないビューを使うと、逐次的でないメモリアクセスにより配列に対する計算が劇的に遅くなることがあります。

規則的にアクセスできないデータを処理の前に規則的にアクセスできる連続配列にコピーしておくことで、性能が大きく向上する場合もあります。例えば次の例では、行列とベクトルから 800,000 個の要素をランダムに切り出して乗算を行っています。最初にビューを通常の配列へコピーすると、コピー処理のコストがあるにもかかわらず乗算の速度は向上します:

julia> using Random

julia> x = randn(1_000_000);

julia> inds = shuffle(1:1_000_000)[1:800000];

julia> A = randn(50, 1_000_000);

julia> xtmp = zeros(800_000);

julia> Atmp = zeros(50, 800_000);

julia> @time sum(view(A, :, inds) * view(x, inds))
  0.412156 seconds (14 allocations: 960 bytes)
-4256.759568345458

julia> @time begin
           copyto!(xtmp, view(x, inds))
           copyto!(Atmp, view(A, :, inds))
           sum(Atmp * xtmp)
       end
  0.285923 seconds (14 allocations: 960 bytes)
-4256.759568345134

コピーが可能なだけのメモリが存在するなら、行列-ベクトル乗算を連続配列に対して行えることによる性能向上はビューを配列にコピーするコストを大きく上回ります。

IO で文字列補間を使わない

ファイルなどの IO デバイスへデータを書き込むときに文字列補間を使うと途中で不必要な文字列が作成されるので、オーバーヘッドになる可能性があります。次のようには書かないでください:

println(file, "$a $b")

次のように書いてください:

println(file, a, " ", b)

最初のバージョンは文字列を作成してからファイルへ書き込みますが、二番目のバージョンは値をファイルに直接書き込みます。また文字列補間を使うとコードが読みにくくなる場合もあります。例えば

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

と、次の呼び出しを比較してください:

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

並列処理におけるネットワーク IO を最適化する

関数のリモートコールを行う次の二つのコードでは、前者の方が高速です:

using Distributed

responses = Vector{Any}(undef, nworkers())
@sync begin
    for (idx, pid) in enumerate(workers())
        @async responses[idx] = remotecall_fetch(foo, pid, args...)
    end
end
using Distributed

refs = Vector{Any}(undef, nworkers())
for (idx, pid) in enumerate(workers())
    refs[idx] = @spawnat pid foo(args...)
end
responses = [fetch(r) for r in refs]

前者は全てのワーカーと一回のネットワークラウンドトリップを行いますが、後者では二回となります ──@spawnatfetch で一回ずつです (この後に wait があればさらにそこでもラウンドトリップが起こります)。fetchwait を逐次的に実行するときも同様に性能が低下します。

非推奨の警告を修正する

非推奨の関数を実行すると、警告を一度だけ出力する処理でメモリの読み込みが行われます。この読み込みが性能を大きく低下させる可能性があるので、非推奨の関数は警告に従って修正してください。

細かい修正を行う

性能が重要な内側のループで性能を向上させるかもしれない細かいポイントをいくつか示します:

性能のための注釈

プログラムの特性を約束することで高度な最適化が可能になる場合があります。

1:n という添え字による AbstractArray へのアクセスは配列が非慣習的な添え字を持つとき安全でないので、境界検査を無効化しているとセグメンテーションフォルトが発生する可能性があります。代わりに LinearIndices(x)eachindex(x) を使ってください (参照: 独自の添え字を持った配列)。

情報

@simd は一番内側の for ループの直前に付ける必要がありますが、@inbounds@fastmath は単一の式あるいはコードブロック全体に適用できます。例えば @inbounds begin あるいは @inbounds for ... とします。

@inbounds@simd のマークアップを使う例を示します (オプティマイザが賢くなりすぎてベンチマークの結果が分かりにくくならないように @noinline を使っています):

@noinline function inner(x, y)
    s = zero(eltype(x))
    for i=eachindex(x)
        @inbounds s += x[i]*y[i]
    end
    return s
end

@noinline function innersimd(x, y)
    s = zero(eltype(x))
    @simd for i = eachindex(x)
        @inbounds s += x[i] * y[i]
    end
    return 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/sec        = ", 2n*reps / time*1E-9)
    time = @elapsed for j in 1:reps
        s += innersimd(x, y)
    end
    println("GFlop/sec (SIMD) = ", 2n*reps / time*1E-9)
end

timeit(1000, 1000)

2.4 GHz の Intel Core i5 プロセッサでは次の結果となりました:

GFlop/sec        = 1.9467069505224963
GFlop/sec (SIMD) = 17.578554163920018

(GFlop/sec が測定された性能であり、大きい値が優れています)

上で示した三つのマークアップを全て使う例を示します。このプログラムはまず一次元配列の有限差分を計算し、その L2 ノルムを評価します:

function init!(u::Vector)
    n = length(u)
    dx = 1.0 / (n-1)

    # u を Vector と宣言したので、1 ベースの添え字を仮定できる
    @fastmath @inbounds @simd for i in 1:n
        u[i] = sin(2pi*dx*i)
    end
end

function deriv!(u::Vector, du)
    n = length(u)
    dx = 1.0 / (n-1)
    @fastmath @inbounds du[1] = (u[2] - u[1]) / dx
    @fastmath @inbounds @simd for i in 2:n-1
        du[i] = (u[i+1] - u[i-1]) / (2*dx)
    end
    @fastmath @inbounds du[n] = (u[n] - u[n-1]) / dx
end

function mynorm(u::Vector)
    n = length(u)
    T = eltype(u)
    s = zero(T)
    @fastmath @inbounds @simd for i in 1:n
        s += u[i]^2
    end
    @fastmath @inbounds return sqrt(s)
end

function main()
    n = 2000
    u = Vector{Float64}(undef, n)
    init!(u)
    du = similar(u)

    deriv!(u, du)
    nu = mynorm(du)

    @time for i in 1:10^6
        deriv!(u, du)
        nu = mynorm(du)
    end

    println(nu)
end

main()

2.7 GHz の Intel Core i7 プロセッサでは次の結果となりました:

$ julia wave.jl;
  1.207814709 seconds
4.443986180758249

$ julia --math-mode=ieee wave.jl;
  4.487083643 seconds
4.443986180758249

ここで使われているコマンドライン引数 --math-mode=ieee によって @fastmath マクロが無効化されるので、結果の比較が可能です。

この例では @fastmath による高速化が 3.7 倍であることが分かります。この値は異常に高く、通常はもっと小さくなります。この例で高速化が大きくなるのは、ベンチマークが処理するデータが小さくプロセッサの L1 キャッシュに全て載るので、メモリアクセスのレイテンシが問題にならず、計算時間が CPU 使用量だけに支配されるためです。しかし現実の問題ではこの条件は成り立ちません。またこの例では計算結果の数値が変化していませんが、一般には結果が少しだけ違う値になることがあり得ます。場合によっては、特に数値的に不安定なアルゴリズムを使ったときは、結果が大きく異なる値になることもあります。

注釈 @fastmath を付けると浮動小数点数演算を含む式が組み変わります。例えば評価の順序が変わったり、変数が特殊な値 (inf, nan) にならないことが仮定されたりします。この例 (をこのマシンで実行した場合) における大きな違いは、deriv! 関数に含まれる 1/(2*dx) がループの外に引っ張り出され、ループの外で計算された idx = 1/(2*dx) を使ってループの中の ... / (2*dx)... * idx と書き換えたような振る舞いになることです。これにより計算はずっと高速になります。

もちろん、コンパイラが実際に適用する最適化、そして最終的な高速化の度合いはハードウェアに大きく依存します。生成されたコードの変化は Julia の code_native 関数を使うと調べられます。

なお、@fastmath は計算で NaN が発生しないことも仮定します。この仮定が不思議な挙動につながることがあります:

julia> f(x) = isnan(x);

julia> f(NaN)
true

julia> f_fast(x) = @fastmath isnan(x);

julia> f_fast(NaN)
false

非正規数をゼロとして扱う

非正規化数 (subnormal number, かつては denormal number とも) が便利な状況もありますが、一部のハードウェアでは非正規化数があると性能が低下します。set_zero_subnormals(true) を呼び出すと浮動小数点数演算の入力および出力に含まれる非正規化数をゼロとして扱うことを許可でき、一部のハードウェアで性能が向上します。set_zero_subnormals(false) を呼び出せば非正規化数に関して厳密な IEEE の振る舞いを強制できます。

一部のハードウェアで起こる非正規化数による性能の低下が分かる例を示します:

function timestep(b::Vector{T}, a::Vector{T}, Δt::T) where T
    @assert length(a)==length(b)
    n = length(b)
    b[1] = 1                            # 境界条件
    for i=2:n-1
        b[i] = a[i] + (a[i-1] - T(2)*a[i] + a[i+1]) * Δt
    end
    b[n] = 0                            # 境界条件
end

function heatflow(a::Vector{T}, nstep::Integer) where T
    b = similar(a)
    for t=1:div(nstep,2)                # nstep は偶数と仮定する。
        timestep(b,a,T(0.1))
        timestep(a,b,T(0.1))
    end
end

heatflow(zeros(Float32,10),2)           # コンパイルを行わせる。
for trial=1:6
    a = zeros(Float32,1000)
    set_zero_subnormals(iseven(trial))  # 奇数回目の試行では厳密な IEEE 算術を使う。
    @time heatflow(a,1000)
end

次のような出力が得られます:

  0.002202 seconds (1 allocation: 4.063 KiB)
  0.001502 seconds (1 allocation: 4.063 KiB)
  0.002139 seconds (1 allocation: 4.063 KiB)
  0.001454 seconds (1 allocation: 4.063 KiB)
  0.002115 seconds (1 allocation: 4.063 KiB)
  0.001455 seconds (1 allocation: 4.063 KiB)

偶数回目の試行が高速であることは明らかです。

この例で a に含まれる値は指数的に減少する曲線を描き、時間の経過と共に少しずつ傾きが緩やかになるので、非正規化数が多く生成されます。

非正規化数をゼロと扱う設定は慎重に利用するべきです。「x-y = 0 ならば x == y」といった恒等命題がこの設定を使うと成り立たなくなります:

julia> x = 3f-38; y = 2f-38;

julia> set_zero_subnormals(true); (x - y, x == y)
(0.0f0, false)

julia> set_zero_subnormals(false); (x - y, x == y)
(1.0000001f-38, false)

一部の応用では、非正規化数をゼロと扱う代わりに少しだけノイズを入れるという対処もできます。例えば a をゼロで初期化するのではなく、次の値で初期化します:

a = rand(Float32,1000) * 1.f-9

@code_warntype を利用する

@code_warntype マクロ (および同じ処理を行う code_warntype 関数) は型に関する問題を診断するのに役立つことがあります。例を示します:

julia> @noinline pos(x) = x < 0 ? 0 : x;

julia> function f(x)
           y = pos(x)
           return sin(y*x + 1)
       end;

julia> @code_warntype f(3.2)
Variables
  #self#::Core.Compiler.Const(f, false)
  x::Float64
  y::UNION{FLOAT64, INT64}

Body::Float64
1       (y = Main.pos(x))
   %2 = (y * x)::Float64
   %3 = (%2 + 1)::Float64
   %4 = Main.sin(%3)::Float64
└──      return %4

@code_warntype マクロの出力の読解には少し慣れが必要です。これは親戚の @code_lowered, @code_typed, @code_llvm, @code_native といったマクロでも同様です。出力されるコードでは最終的なコンパイル結果の機械語を生成する処理がかなり進んでおり、大部分の式には ::T で示される型が付いています (TFloat64 といった値です)。@code_warntype の最も重要な特徴が、具象型でない型が赤く表示されることです: Markdown で書かれるこのドキュメントでは色を付けられないので、上の実行例では赤い文字を大文字で示しています。

Body::Float64 の行は推論された関数の返り値の型を示し、以降の行は f の本体を Julia の SSA IR 形式で表しています。一番左に付いている数字はラベルであり、goto によるジャンプ先を表します。関数の本体を見ると、最初に pos が呼ばれていることが分かります。pos の返り値 y の推論された型を Variables の欄で確認すると、UNION{FLOAT64, INT64} という Union 型であることが分かります。大文字で表示されるこの型は具象型ではありません。これは pos の返り値の正確な型が入力の型からは決定できないことを意味します。しかし、yFloat64Int64 のどちらであっても、y*x の結果は Float64 で変わりません。つまり f(x::Float64) の途中の処理では型が不安定になるにもかかわらず、f(x::Float64) は入力の型に対して型安定であるというのが最終的な結論です。

@code_warntype で得られる情報をどう使うかはあなた次第です。pos を型安定な関数として書き直すべきなのは間違いありません: そうすれば f に含まれる変数は全て具象型となり、性能は最適となるはずです。ただ、この種の一時的な型の不安定性があまり問題にならない場合もあります: 例えば pos が単独で使われないなら、f の出力の型が (入力 がFloat64 型のとき) 安定であるという事実が型の不安定性を以降のコードに伝播させない防波堤となります。これは型の不安定性の解決が難しいもしくは不可能な場合に特に重要です。そのような場合には、上述のアドバイス (型注釈を入れる/関数を分割する) が型の不安定性による “損害” を抑え込むための最良のツールです。

型が不安定な関数は Julia Base にさえ含まれることも知っておいてください。例えば findfirst 関数はキーが見つかれば配列に対する添え字を返し、見つからなければ nothing を返すので、この関数の型は明らかに不安定です。重要な型の不安定性を目立たせるために、Missing または Nothing を含む Union は赤ではなく黄色でハイライトされます。

葉型でない型を持つ式の発見・対処に慣れてもらうために、そういった式の例を示します:

変数のキャプチャを工夫する

関数の中で関数を定義する次のコードを考えます:

function abmult(r::Int)
    if r < 0
        r = -r
    end
    f = x -> x * r
    return f
end

abmult 関数は「r の絶対値と引数の乗算結果を返す関数」f を返します。f に代入されている内部関数は「クロージャ」と呼ばれます。クロージャは do ブロックやジェネレータ式でも利用されます。

このスタイルのコードは言語として解決すべき性能の問題を提示します。パーサーはこのコードをトップレベルの命令に変換するときクロージャを個別のブロックとして切り出すので、このコードは大きく形を変えます。クロージャと外側のスコープで共有される r のような「キャプチャされた」変数も切り出され、ヒープにアロケートされた「ボックス」に移されます。クロージャが持つ r と外側の r が同一であることを言語が定めているので、このボックスにはクロージャとその外側の両方からアクセスが可能です。そのためクロージャの定義の後に外側のスコープ (あるいは他のクロージャ) で r が変更されると、クロージャが参照する値も変更されます。

前の段落の議論で「パーサー」という言葉を使いました。パースは abmult を含むモジュールを最初に読み込んだ段階で行われるコンパイル処理であり、後になって関数を最初に呼び出したときに行われる処理ではありません。そのためパーサーは Int が固定された型であることを知らず、また r = -rIntInt に変える処理であることも知りません。型推論の魔法はコンパイルのもっと後の段階で起きます。

つまり、パーサーは r が固定された型 Int を持つことも、クロージャが作られた後に r が値を変えない (よってボックスが必要にならない) ことも理解できません。そのためパーサーはオブジェクトを Any のような抽象型かのように保持するボックスを使ったコードを生成し、実行時には r を使うたびにディスパッチが必要になります。これは @code_warntypeabmult が返す関数に適用すれば確認できます。ボックスと実行時の型ディスパッチの両方が性能を低下させます。

もしコードの性能が重要な部分がキャプチャされた変数を使っているなら、次のアドバイスを使えば性能をいくらか取り戻せるでしょう。まず、キャプチャされた変数の型が変わらないと分かっているなら、型注釈を使ってその型を明示的に宣言してください (宣言するのは右辺の式の型ではなく左辺の変数の型です):

function abmult2(r0::Int)
    r::Int = r0
    if r < 0
        r = -r
    end
    f = x -> x * r
    return f
end

パーサーがボックスに含まれるオブジェクトに具象型を関連付けられるようになるので、この型注釈を使うとキャプチャによる性能の低下を緩和できます。さらに、もし (作成後に再代入されないために) キャプチャされた変数をボックスを入れなくて構わないなら、let ブロックを使ってそのことを提示できます:

function abmult3(r::Int)
    if r < 0
        r = -r
    end
    f = let r = r
            x -> x * r
    end
    return f
end

このコードで let ブロックはスコープがクロージャと等しい新たな変数 r を作成します。この二つ目のテクニックを使うとキャプチャされた変数が存在する場合でも言語が提供する完全な性能を達成できます。

キャプチャされた変数に関する機能はコンパイラの中で高速に進化している領域であり、将来のリリースではプログラマからの注釈が無くとも性能を保てるようになる予定です。それまでの間は、ユーザーによって作られた FastClosures.jl などのパッケージが利用できます。このパッケージは abmult3 のような let の利用を自動化できます。

シングルトンを使った等価性の判定

ある値が特定のシングルトンかどうかを判定するときは、== ではなく === を使うと性能が向上します。同様に != ではなく !== を使うべきです。この種の判定は頻繁に起こります: 例えば反復プロトコルを実装するときは、iterate の返り値が nothing かどうかを判定するはずです。


日本語 Julia 書籍 (Amazon アソシエイト)
1 から始める Julia プログラミング
Julia プログラミングクックブック―言語仕様からデータ分析、機械学習、数値計算まで
スタンフォード ベクトル・行列からはじめる最適化数学