Julia の基礎

julia> collect(1:3)
3-element Vector{Int64}:
 1
 2
 3

スコープ

for ループの中で定義された変数は、ループ変数を含めて for ループの外では使えない。

julia> for i = 1:5
           println(i)
       end
1
2
3
4
5

julia> i
ERROR: UndefVarError: `i` not defined

while ループ

while ループは、条件が成り立つあいだ繰り返す。

julia> i = 1
0

julia> while i <= 5
           println(i)
           i += 1
       end
1
2
3
4
5

continue と break

continue でループを次に進めることができる。

julia> for i = 1:10
           if i % 2 == 0
               continue
           end
           println(i)
       end
1
3
5
7
9

break でループを抜けることができる。

julia> for i = 1:10
           if i == 5
               break
           end
           println(i)
       end
1
2
3
4

goto

一応、goto 文を実現するマクロがある。たとえば、入れ子のループを一気に抜けたい場合はは、次のようにする。

begin
    for i = 1:5
        for j = 1:5
            if i == 3
                @goto loop_end
            end
            println(i, " ", j)
        end
    end
    @label loop_end    
end

実行結果

$ julia goto.jl
1 1
1 2
1 3
1 4
1 5
2 1
2 2
2 3
2 4
2 5

ただし、この場合、次のようにループを構成すれば、ふつうに break で抜けられる。

for i = 1:5, j = 1:5
  ...

例外処理

関数に変な値を渡すとエラーになる。

julia> sqrt("a")
ERROR: MethodError: no method matching sqrt(::String)

これをユーザー側で受け止めることもできる。

julia> try
           sqrt("a")
       catch e
           println("Error")
       end
Error

try...catch には else 節をつけることができる。これは、例外が起こらなかった場合に実行される。

julia> x = 0
0

julia> try
           x = sqrt(2)
       catch e
           println("Error")
       else
           println(x)
       end
1.4142135623730951

上の例では最初に変数を定義しているが、try...catch の節はそれぞれローカルなスコープをもつためである (else の中から try の中の変数が見えない)。

finally 節もある。この節は、例外が起ころうが起こるまいが実行される。

例外を自分で投げる場合は、throw() を用いる。

julia> throw(ArgumentError("Error"))
ERROR: ArgumentError: Error

自分で例外を定義する場合は、次のようにする。

julia> struct MyException <: Exception end

プログラムの停止

プログラムを止めるには、exit() を用いる。

julia> exit(2)
$ echo $?
2

引数に終了コードを指定できる。引数を省略すれば 0 を指定したことになる。

関数

関数は、引数を受け取って値を返す処理のまとまりである。次のように定義する。

julia> function f(x, y)
           x + y
       end

julia> f(1, 2)
3

最後の式の値が返り値になる。この場合、次のようにも表現できる。

julia> f(x, y) = x + y

return も使える。

julia> function f(x, y)
           return x + y
       end

複数の値を返したければ、タプルを使う。

julia> function f(x, y)
           x + y, x - y
       end

julia> f(1, 2)
(3, -1)

値を何も返したくない場合は、nothing を返せばよい。

julia> function f(x)
           x[1] = 0
           nothing
       end

関数を変数に入れることもできる。

julia> g = f

julia> g(1, 2)
3

Julia の関数は、同じ名前でも引数の数や型が異なれば別の関数とみなされる。

引数

関数内で引数の値を変更したとき、単純な型であれば呼び出し元の変数の値は変わらない。

julia> function f(x)
           x = 0
       end

julia> x = 1
1

julia> f(x)
0

julia> x
1

ただし、引数のオブジェクトの値 (内部状態) を変更した場合、呼び出し元のオブジェクトの値も変わってしまう。

julia> function f(x)
           x[1] = 0
       end

julia> x = [1, 2, 3]
3-element Vector{Int64}:
 1
 2
 3

julia> f(x)
0

julia> x
3-element Vector{Int64}:
 0
 2
 3

Julia の変数は値 (オブジェクト) に結び付けられた名前なので、変数が「持っている」のはオブジェクトの参照であり、引数で渡されるのもそれである。

キーワード引数

; 以降の引数は、キーワード引数になる。

julia> function f(x ; y)
           x - y
       end

julia> f(3, y=2)
1
# f(3, 2) だと怒られる

; 以降の引数以外はキーワード指定できない。; の位置の違いだけでも別の関数と区別されるので注意。

オプション引数

引数に値を与えておけば、その引数は指定してもしなくてもよい。

julia> function f(x=1, y=2)
           x + y
       end

julia> f()
3

オプション引数は、引数リストの最後に位置する必要がある (; の直前でもよい)。

可変長引数

可変長引数にするには、"x..." の形にする。

julia> function f(x...)
           sum(x)
       end

julia> f(1, 2, 3)
6

ここで、x の型はタプルになっている。

"..." は引数リストの最後に位置する必要がある (; の直前でもよい)。

この "..." は、ベクトルの要素を引数に展開するのにも使える。

julia> f([1, 2, 3]...)
6

型指定

引数および返り値には型を指定できる。

julia> function f(x::Vector{Int64})::Int64
           x[1] = 0
       end

型に合わないものを渡すと怒られる。

julia> x = 1
1

julia> f(x)
ERROR: MethodError: no method matching setindex!(::Int64, ::Int64, ::Int64)

再帰呼び出し

関数の再帰呼び出しはふつうにできる。

function fibonacci(n::Int64)::Int64    
    if n == 0    
        return 0    
    elseif n == 1    
        return 1    
    else    
        return fibonacci(n - 2) + fibonacci(n - 1)    
    end    
end 

無名関数

無名関数も定義できる。

julia> f = (x, y) -> x + y

julia> f(1, 2)
3

=> ではなくて -> であることに注意 (=> はペア)。

無名関数は、map() のようにベクトルの要素に関数を適用するなどのときに使う。

julia> map(x -> 2x, [1, 2, 3])
3-element Vector{Int64}:
 2
 4
 6

do 構文

map() のように関数が渡せる関数には、do 構文というものが使える。

julia> map([1, 2, 3]) do x
           2x
       end
3-element Vector{Int64}:
 2
 4
 6

複合関数とパイプ

∘ (\circ) で複合関数を作ることができる。

julia> (sqrt ∘ +)(1, 2)
1.7320508075688772

これは以下と同じである。

julia> sqrt(1 + 2)
1.7320508075688772

関数の実行結果を別の関数に流し込む「パイプ」 |> もある。

julia> 1 + 2 |> sqrt
1.7320508075688772

.|> でベクトルのそれぞれの要素にパイプを適用できる。

julia> [1, 2, 3] .|> [x -> 2x, sqrt, sin]
3-element Vector{Real}:
 2
 1.4142135623730951
 0.1411200080598672

ドット構文

関数の最後にドット . をつけると、ベクトルやタプルの要素にそれぞれ関数を適用するようになる。

julia> function f(x, y)
           x + y
       end

julia> f([1, 2], [3, 4])
2-element Vector{Int64}:
 4
 6

julia> f.([1, 2], [3, 4])
2-element Vector{Int64}:
 4
 6

julia> f((1, 2), (3, 4))
ERROR: MethodError: no method matching +(::Tuple{Int64, Int64}, ::Tuple{Int64, Int64})

julia> f.((1, 2), (3, 4))
(4, 6)

ベクトルはベクトル同士の足し算に対応しているため、ドットを付けても付けなくても変わらない。タプルはタプル同士の足し算に対応していないため、ドットを付けないとエラーになる。

@. でも同じような効果がある。

julia> @. f((1, 2), (3,4))
(4, 6)

クロージャ

次のように、呼ぶたびに 1 増えた数値を返す関数を定義したい。

julia> f()
1

julia> f()
2

julia> f()
3

そのためには、ローカルの環境にオブジェクトを閉じ込める形で関数を定義すればよい。この仕組みをクロージャという。目的の関数は、次のように実現できる。

julia> f = let x = 0; () -> x += 1 end

julia> f()
1

julia> f()
2

julia> f()
3

f には関数が入っていて、その関数は let ブロック内の x を参照し続けるため、こうなる。

上の例だと、x を初期化できない。代わりに、次のような実装も考えられる。

julia> setup = (x) -> () -> x += 1

julia> f = setup(0)

julia> f()
1

julia> f()
2

julia> f()
3

上の考え方を使うと、最初の例は、次のようにも書ける。

julia> f = ((x) -> () -> x += 1)(0)

メソッド

C++ などのオブジェクト指向プログラミング対応言語では、データと手続きをひとまとまりとして、オブジェクトの設計図であるクラスを構成することができる。そのときの手続きをメンバー関数とかメソッドなどという。Julia では、これに似たことを "ディスパッチ" (あるいは "多重ディスパッチ") という仕組みで実現する。あるひとつの抽象的機能 (たとえば「足し算」など) を実行する同じ名前の関数群を用意して、その中から引数のデータ型や数に応じた適切な関数を選び出す。この関数群をメソッドと呼ぶ。

たとえば、+ 演算子は、これも実は関数だが、メソッドとしての関数のリストは methods() で得られる。

julia> methods(+)
# 207 methods for generic function "+" from Base:
   [1] +(x::T, y::T) where T<:Union{Int128, Int16, Int32, Int64, Int8, UInt128, UInt16, UInt32, UInt64, UInt8}
     @ int.jl:87
   [2] +(x::T, y::T) where T<:Union{Float16, Float32, Float64}
     @ float.jl:408
   [3] +(c::Union{UInt16, UInt32, UInt64, UInt8}, x::BigInt)
     @ Base.GMP gmp.jl:539

       ...

 [207] +(a, b, c, xs...)
     @ operators.jl:578

数値でもベクトルでも同じように + を使えるのは、この仕組みによる。また、タプルで + が使えないのも、この中にタプル用の関数が定義されていないからである。

多重ディスパッチとは

多重ディスパッチとは？ 多重でないディスパッチがあるのか？ そもそもディスパッチとは？

ディスパッチとは、オブジェクト指向プログラミングのための機構で、同じ名前の関数群から、オブジェクトのクラスに応じた関数を「動的に」決定して呼び出す仕組みのことを言うようである。「動的に」がポイントで、C++ の関数の多重定義のことをディスパッチとは言わない (どの関数を呼び出すかはコンパイル時に決定されているため)。C++ の例でいうと、ディスパッチと呼べるのは仮想関数 (多相性) の仕組みである。

C++ の仮想関数の仕組みをディスパッチとみなしたとして、これは「多重」(multiple) とは呼ばれない。というのも、「どのクラスを選ぶか」というところは動的に決まるが、そのあとは選ばれたクラスの関数から呼び出す関数が静的に決まるためである。これは、実質的には関数の第 1 引数の型だけで関数を選んでいることになるため (Python のメソッドの self のイメージ)、「単一」(single) と呼ばれる。

一方、多重ディスパッチは、CLOS (Common Lisp Object System) の仕組みのことをそう呼ぶようであるが、呼び出される関数が複数の引数の型によって決まる。

パラメータメソッド

データ型でパラメータ型 (parametric type) というものがあったが、関数でもデータ型をパラメータにすることができる。

julia> function f(x::T, y::T) where T
           x + y
       end

julia> f(1, 2)
3

上記の例は、x と y は任意の型でよいが、同じ型でなければならない、ということを意味している。型が違っていると怒られる。

julia> f(1, 2.9)
ERROR: MethodError: no method matching f(::Int64, ::Float64)

数値型なら OK、という指定をしたい場合は、次のようにする。

julia> function f(x::T, y::T) where T<:Number
           x + y
       end

一見、f(x::Number, y::Number) とどう違うのか、と思うかもしれないが、前者は相変わらず x と y の型が異なると怒られるが、後者は型が異なっても問題ない。

引数の型は揃っていてほしいが、型が揃っていなくても適切な変換で型を揃えることができる場合は、promote() で対応できる。上の例では、別途次のような関数を定義しておけばよい。

julia> function f(x::Number, y::Number)
           f(promote(x, y)...)
       end

julia> f(1, 2.9)
3.9

明示的に変換を記述する方法もある。

julia> function f(x::Int64, y::Float64)
           f(convert(Float64, x), y)
       end
f (generic function with 2 methods)

julia> f(1, 2.9)
3.9

関数オブジェクト

関数のように振る舞うオブジェクト (function-like object) を定義できる。

julia> mutable struct Point
           x
           y
       end

julia> function (p::Point)(x, y)
           p.x = x
           p.y = y
           p
       end

julia> p = Point(0, 0)
Point(0, 0)

julia> p(1, 2)
Point(1, 2)

コンストラクタ

構造体 (複合型) のオブジェクトは、その型に引数を与えて生成する。

julia> mutable struct Point
           x
           y
       end

julia> p = Point(1, 2)
Point(1, 2)

これは次のように書いても同じである。

julia> mutable struct Point
           x
           y
           Point(x, y) = new(x, y)
       end

これを用いて、コピーコンストラクタを作るには、次のようになる。

julia> mutable struct Point
           x
           y
           Point(x, y) = new(x, y)
           Point(p) = new(p.x, p.y)
       end

julia> p = Point(1, 2)
Point(1, 2)

julia> p2 = Point(p)
Point(1, 2)

julia> p2.x = 2
2

julia> p
Point(1, 2)

julia> p2
Point(2, 2)

パラメータ型を用いる場合は、new() にも型指定がいる。

julia> mutable struct Pont{T}
           x::T
           y::T
           Point(x::T, y::T) where T = new{T}(x, y)
       end

コンストラクタは外で定義してもよい (new() は使えない)。

julia> mutable struct Point                                                      
           x                                                                     
           y                                                                     
       end

julia> Point(p) = Point(p.x, p.y)

型変換

型変換は convert() で行うが、自前の型変換を定義したい場合は、convert() を用意すればよい。たとえば、次のような Point 型があり、タプルからのインスタンスの生成にも対応しているものとする。

julia> mutable struct Point{T<:Number}
           x::T
           y::T
           Point(x::T, y::T) where T = new{T}(x, y)
           Point(t::Tuple{T, T}) where T = new{T}(t...)
       end

julia> Point(1, 2)
Point{Int64}(1, 2)

julia> Point((1, 2))
Point{Int64}(1, 2)

このとき、タプルから Point 型への変換は、次のように書ける。

julia> convert(::Type{Point{T}}, x::Tuple{T, T}) where T = Point(x)

julia> convert(Point{Int64}, (1, 2))
Point{Int64}(1, 2)

プロモーション

promote() は引数の型を揃えてタプルで返す。これを自前で用意するには、promote() を定義すればよい (以下の例はファイルの形で実行)。

mutable struct Point{T<:Number}
   x::T
   y::T
   Point(x::T, y::T) where T = new{T}(x, y)
   Point(t::Tuple{T, T}) where T = new{T}(t...)
end

convert(::Type{Point{T}}, x::Tuple{T, T}) where T = Point(x)
promote(x::Point{T}, y::Tuple{T, T}) where T = (x, convert(Point{T}, y))

println(promote(Point(1, 2), (1, 2)))

代わりに promote_rule() を定義する方法もある。

Base.convert(::Type{Point{T}}, x::Tuple{T, T}) where T = Point(x)
Base.promote_rule(::Type{Point{T}}, ::Type{Tuple{T, T}}) where T = Point{T}

コードの分割

include() で、別ファイルのコードを取り込むことができる。

julia> include("geometry.jl")

コードを機能ごとに分けて、メインのコードで取り込むといったことができる。

モジュール

モジュールは、プログラムを分割管理するのに使う。モジュールは名前空間を区切ることができる。

julia> module Geometry
       struct Point
           x
           y
       end
       end

julia> p = Geometry.Point(1, 2)
Main.Geometry.Point(1, 2)

モジュールを使う方法として、次のような方法がある。

julia> module Geometry
       export Point
       struct Point
           x
           y
       end
       end

julia> using .Geometry

julia> p = Point(1, 2)
Point(1, 2)

using でモジュール内にアクセスできるようにする。ここでは、export されている Point を直接 (モジュール名なしで) 扱えるようになる。

使う項目を直接指定する方法もある。

julia> using .Geometry: Point

この場合、export 指定にかかわらず、Point を直接扱える。

項目の名前を付け直すこともできる。

julia> using .Geometry: Point as pt

サブモジュールも扱える。

julia> module A
       module B
       struct BS
           x
       end
       end
       using .B: BS
       struct AS
           AS(x) = BS(x)
       end
       end
Main.A

julia> using .A: AS

julia> as = AS(1)
Main.A.B.BS(1)

using に似たものとして import がある。using と異なるのは、import は export を無視することである。using だと、使わないものも含めて export されているものが使える状態になってしまう。それが気持ち悪い場合は、import で運用したほうがよいのだろう。

パッケージ

Julia では、様々なライブラリがパッケージの形で配布されている。Julia のコマンドラインで "]" を押すとパッケージモードに入るが、これが Julia のパッケージマネージャに当たる。

(@v1.9) pkg>

パッケージの情報は、status コマンドで見ることができる (最初はグローバルの環境が作成される)。

(@v1.9) pkg> status
Status `~/.julia/environments/v1.9/Project.toml` (empty project)

パッケージは、グローバルの環境に追加してもよいし、特定の環境に追加してもよい。この環境をプロジェクトと呼ぶ。作業環境として test プロジェクトを作ってみる (";" でシェルモードに入る。行頭バックスペースで元に戻る)。空のディレクトリを作って移動。

shell> mkdir test

shell> cd test

プロジェクトを有効化。

(@v1.9) pkg> activate .
  Activating new project at `~/work/julia/test`

(test) pkg>

パッケージを追加。

(test) pkg> add CSV

パッケージは using や import で使うことができる。

julia> import CSV

パッケージモードの status でパッケージが追加されていることがわかる。

(test) pkg> status
Status `~/work/julia/test/Project.toml`
  [336ed68f] CSV v0.10.11

パッケージを消す場合は、rm を使う。

(test) pkg> rm CSV

グローバル環境に戻りたければ、単に activate を打てばよい。

(test) pkg> activate
  Activating project at `~/.julia/environments/v1.9`

プロジェクトをちゃんと作る

generate コマンドでちゃんとしたプロジェクトを作ることができる。

(@v1.9) pkg> generate test
  Generating  project test:
    test/Project.toml
    test/src/test.jl

(@v1.9) pkg> activate test
  Activating project at `~/work/julia/test`

shell> cd test

ソースコードができている。

shell> vi src/test.jl

test.jl

module test
    
greet() = print("Hello World!")
    
end # module test

これはすでにパッケージの形をしており、import で使ってみることができる。

julia> import test

julia> test.greet()
Hello World!

サブプロジェクトのパッケージのインポート

サブプロジェクトなどのパッケージのインポートを行うには、LOAD_PATH にパスを追加する方法がある。

julia> push!(LOAD_PATH, "./test")

julia> import test

あるいは、パッケージのディレクトリに移動してインポートする方法もある (ドットがいる)。

julia> cd("test")

julia> import .test

ドキュメント

モジュールや構造体、関数の上に文字列を置いておけば、それがドキュメントになる。

geometry.jl

"""
  Geometry module
"""
module Geometry

"""
  Point
"""
struct Point
    x
    y
end

end

ドキュメントは "?" で確認できる ("?" と入れるとプロンプトが "help?>" に変わる)。

julia> include("geometry.jl")
Main.Geometry

help?> Geometry
search: Geometry

  Geometry module

help?> Geometry.Point
  Point

マクロ

マクロは、言語を拡張するための機能である。たとえば、以下のように使う。

julia> macro sayhello()
           return :( println("Hello") )
       end

julia> @sayhello
Hello

引数を取る場合は、次のようになる。

julia> macro saymsg(msg)
           return :( println($msg) )
       end

julia> @saymsg("Hello")
Hello

あるいは

julia> @saymsg "Hello"
Hello

これだけだと関数と変わらないが、マクロでは、次のようなことができる。

julia> macro saymsg(msg)
           return :( println($(string(msg))) )
       end

julia> @saymsg 1 + 2
1 + 2

関数の場合、引数に 1 + 2 の結果が渡ってしまうため、上のようなことができない。これを用いて、コードのチェックを行う assert マクロが定義されている。

julia> @assert 1 == 0
ERROR: AssertionError: 1 == 0

論理式をチェックして、結果が false だったら例外を投げてその論理式を表示する。実装は、たとえば次のようになる。

julia> macro assert(ex)
           return :( $ex ? nothing : throw(AssertionError($(string(ex)))) )
       end

マクロが返している :(...) はシンボルである。Meta.parse() で、コードの文字列からシンボルを得ることができる。

julia> ex = Meta.parse("""println("Hello")""")
:(println("Hello"))

これは、eval() で実行できる。

julia> eval(ex)
Hello

引数込みの関数呼び出しを 2 回行う、ということができる。

julia> x = """println("Hello")"""
"println(\"Hello\")"

julia> ex = Meta.parse("$x; $x")
:($(Expr(:toplevel, :(println("Hello")), :(println("Hello")))))

julia> eval(ex)
Hello
Hello

これをマクロ化すると、次のようになる。

julia> macro double(x)
           return Meta.parse("$x; $x")
       end

julia> @double println("Hello")
Hello
Hello

マクロが生成するシンボルを確認するには、macroexpand マクロを用いる。

julia> @macroexpand @double println("Hello")
:($(Expr(:toplevel, :(Main.println("Hello")), :(Main.println("Hello")))))

コマンドライン引数

コマンドラインからスクリプトと引数を指定した場合、スクリプト名は PROGRAM_FILE、引数は ARGS に入る。

args.jl

@show PROGRAM_FILE
@show ARGS

$ julia args.jl 1 2 3
PROGRAM_FILE = "args.jl"
ARGS = ["1", "2", "3"]

文字列

文字列の変換

文字列は文字を変更できない。文字列を編集したければ変換するしかない。そのためには、一旦ベクトルに変換したほうがよいかもしれない。

julia> collect("abc")
3-element Vector{Char}:
 'a': ASCII/Unicode U+0061 (category Ll: Letter, lowercase)
 'b': ASCII/Unicode U+0062 (category Ll: Letter, lowercase)
 'c': ASCII/Unicode U+0063 (category Ll: Letter, lowercase)

ただ、どうせループを回すなら、ベクトルに変換すせずに直接回すことができる。

julia> for c in "abc"
           println(c)
       end
a
b
c

たとえば、小文字を大文字に変換する関数を雑に作るとすると、次のようになるだろう。

julia> function f(x)
           r = []
           for c in x
               if 'a' <= c <= 'z'
                   push!(r, c + ('A' - 'a'))
               else
                   push!(r , c)
               end
           end
           string(r...)
       end

julia> f("abc")
"ABC"

julia> f("aB6")
"AB6"

関数の最後で、ベクトルを文字列に変換している。

文字列を繰り返す

文字列を繰り返すには、repeat() を用いる。

julia> repeat("=", 20)
"===================="

連結と分割

ベクトルを文字列として連結したり、文字列をベクトルに分割したりするには、それぞれ join()、split() を用いる。

julia> join(["sabatora", "chatora", "mike"], ",")
"sabatora,chatora,mike"

julia> split("sabatora,chatora,mike", ",")
3-element Vector{SubString{String}}:
 "sabatora"
 "chatora"
 "mike"

文字を探す

文字列から文字を探すには、findfirst()、findnext()、findprev() を使う。

julia> findfirst('n', "penguin")
3

julia> findnext('n', "penguin", 4)
7

julia> findprev('n', "penguin", 6)
3

文字列がある文字列を含んでいるかどうかは、occursin() で調べられる。

julia> occursin("gui", "penguin")
true

julia> occursin("cui", "penguin")
false

文字数を数えるだけなら、count() が使える。

julia> count('n', "penguin")
2

正規表現

r"..." の形で正規表現を表すことができる。

julia> typeof(r"\s*$")
Regex

match() で条件に合致しているものを取り出したりできる。

julia> m = match(r"\s(p.*)[\s\.]", "This is a penguin.")
RegexMatch(" penguin.", 1="penguin")

julia> m[1]
"penguin"

括弧でくくったところ (この場合は "(p.*)") を取り出せる。

取り出す部分に名前をつけることもできる。

julia> m = match(r"(?\d+):(?\d+)", "12:45")
RegexMatch("12:45", hour="12", minute="45")

julia> m[:hour]
"12"

replace() で文字の置き換えができる。

julia> replace("This is a penguin.    ", r"\s*$" => "")
"This is a penguin."

ベクトル

大きさだけ確保する

とりあえずベクトルの大きさだけ確保するには、次のようにする。

julia> a = Vector{Float64}(undef, 3)
3-element Vector{Float64}:
 0.0
 0.0
 0.0

ここで undef は未初期化を意味する。

zeros() でもよい。

julia> zeros(3)
3-element Vector{Float64}:
 0.0
 0.0
 0.0

fill() でも同様のことができる。

julia> a = fill(0.0, 3)
3-element Vector{Float64}:
 0.0
 0.0
 0.0

ある値で埋める

全要素を特定の値にしたい場合などは、fill!() が使える。

julia> a = [1, 2, 3]
3-element Vector{Int64}:
 1
 2
 3

julia> fill!(a, 0)
3-element Vector{Int64}:
 0
 0
 0

以下のようにしてもよい。

julia> a .= 0
3-element Vector{Int64}:
 0
 0
 0

結合する

ベクトルを結合するには、次のようにする。

julia> a = [1, 2, 3]
3-element Vector{Int64}:
 1
 2
 3

julia> b = [4, 5, 6]
3-element Vector{Int64}:
 4
 5
 6

julia> [a; b]
6-element Vector{Int64}:
 1
 2
 3
 4
 5
 6

ベクトルへの変換

範囲型などからベクトルを作るには、collect() を使う。

julia> a = collect(1:5)
5-element Vector{Int64}:
  1
  2
  3
  4
  5

最初と最後を取り出す

ベクトルから最初と最後の要素を取り出すには、それぞれ first()、last() がある。

julia> first(collect(1:100), 3)
3-element Vector{Int64}:
 1
 2
 3

julia> last(collect(1:100), 3)
3-element Vector{Int64}:
  98
  99
 100

要素が含まれるかどうか

要素が含まれるかどうかは、in でわかる。

julia> "a" in ["a", 2, "?"]
true

julia> 2 in ["a", 2, "?"]
true

julia> 1 in ["a", 2, "?"]
false

∈ (\in) も使える。

julia> "a" ∈  ["a", 2, "?"]
true

含まれないかどうかを調べる場合は、! で逆にするか、∉ (\notin) を用いる。

julia> !(1 in ["a", 2, "?"])
true

julia> 1 ∉ ["a", 2, "?"]
true

要素の並び替え

要素を並び替えるには sort!() を用いる。

julia> a = [1, 2, 3, 5, 6, 4]
6-element Vector{Int64}:
 1
 2
 3
 5
 6
 4

julia> sort!(a)
6-element Vector{Int64}:
 1
 2
 3
 4
 5
 6

sort() (! がないバージョン) もある。

引数 rev で逆順にもできる。

julia> sort(a, rev=true)
6-element Vector{Int64}:
 6
 5
 4
 3
 2
 1

単純に逆にするだけなら reverse() でよい。

julia> a = [1, 2, 3, 5, 6, 4]
6-element Vector{Int64}:
 1
 2
 3
 5
 6
 4

julia> reverse(a)
6-element Vector{Int64}:
 4
 6
 5
 3
 2
 1

重複要素の削除

重複要素を削除するには unique() を使う。

julia> a = [1, 2, 2, 5, 4, 6, 3]
7-element Vector{Int64}:
 1
 2
 2
 5
 4
 6
 3

julia> unique(a)
6-element Vector{Int64}:
 1
 2
 5
 4
 6
 3

添字がほしい場合は、次のような方法がある。

julia> unique(i -> a[i], eachindex(a))
6-element Vector{Int64}:
 1
 2
 4
 5
 6
 7

最大値と最小値

要素の最大値、最小値を得るには、それぞれ maximum()、minimum() を用いる。

julia> a = rand(5)
5-element Vector{Float64}:
 0.03623237283710057
 0.39652922861963347
 0.18058660806946192
 0.035906856264285025
 0.9215430627026505

julia> maximum(a)
0.9215430627026505

julia> minimum(a)
0.035906856264285025

最大値と最小値を両方ほしい場合は、extrema() を使う方法もある。

julia> extrema(a)
(0.035906856264285025, 0.9215430627026505)

最大値、最小値のそれぞれの添字がほしい場合は、argmax()、argmin() がある。

julia> argmax(a)
5

julia> argmin(a)
4

最大値、最小値とそれぞれの添字の両方ほしい場合は、findmax()、findmin() を用いる。

julia> findmax(a)
(0.9215430627026505, 5)

julia> findmin(a)
(0.035906856264285025, 4)

各要素に関数を適用した結果の最大値、最小値を得ることもできる。

julia> a = ["annko", "omochi", "mizore"]
3-element Vector{String}:
 "annko"
 "omochi"
 "mizore"

julia> maximum(length, a)
6

julia> minimum(length, a)
5

extrema() は次のような使い方もできる。

julia> extrema(sin, 0:0.001:pi)
(0.0, 0.9999999792586128)

絶対値

abs() で絶対値を得られるが、ベクトル用の関数はない。ただし、ドットを使えばよいだけである。

julia> a = [1, -2, 3, -4, 5]
5-element Vector{Int64}:
  1
 -2
  3
 -4
  5

julia> abs.(a)
5-element Vector{Int64}:
 1
 2
 3
 4
 5

合計

sum() で合計を取ることができる。

julia> a = [1, 2, 3, 4, 5]
5-element Vector{Int64}:
 1
 2
 3
 4
 5

julia> sum(a)
15

足していく途中がほしい場合は cumsum() を使う。

julia> cumsum(a)
5-element Vector{Int64}:
  1
  3
  6
 10
 15

積

prod() で各要素の積を取ることができる。

julia> a = [1, 2, 3, 4, 5]
5-element Vector{Int64}:
 1
 2
 3
 4
 5

julia> prod(a)
120

積の途中がほしい場合は cumprod() を用いる。

julia> cumprod(a)
5-element Vector{Int64}:
   1
   2
   6
  24
 120

演算の積算

reduce() で、ある演算子を要素に適用して積算した結果が得られる。

julia> a = [1, 2, 3, 4, 5]
5-element Vector{Int64}:
 1
 2
 3
 4
 5

julia> reduce(+, a)
15

julia> reduce(*, a)
120

初期値を指定できる。

julia> reduce(*, a, init=-1)
-120

演算の順番を指定する方法もある。

julia> foldl(=>, a)
(((1 => 2) => 3) => 4) => 5

julia> foldr(=>, a)
1 => (2 => (3 => (4 => 5)))

積算の途中がほしい場合は、accumulate() を用いる。

julia> accumulate(+, [1, 2, 3, 4, 5])
5-element Vector{Int64}:
  1
  3
  6
 10
 15

階差

diff() で階差を作る。

julia> diff(cumprod([1, 2, 3, 4, 5]))
4-element Vector{Int64}:
  1
  4
 18
 96

ベクトル同士の演算

ベクトルの要素同士を足し合わせるには、+ でよい。

julia> a = [1, 2, 3]
3-element Vector{Int64}:
 1
 2
 3

julia> b = [4, 5, 6]
3-element Vector{Int64}:
 4
 5
 6

julia> a + b
3-element Vector{Int64}:
 5
 7
 9

- も同様である。* は使えない。ベクトル同士の積がほしい場合は、.* を使う。

julia> a .* b
3-element Vector{Int64}:
  4
 10
 18

スカラー相手では、* は機能するが、+、- については機能しない。代わりに .+、.- を使う。

各要素への関数の適用

別のベクトルを作る

ベクトルの各要素へ関数を適用して別のベクトルを作るには、map() を用いる。

julia> a = collect(1:5)
5-element Vector{Int64}:
  1
  2
  3
  4
  5

julia> map(x -> 2x + 1, a)
5-element Vector{Int64}:
  3
  5
  7
  9
 11

あるいは次のようにも書ける。

julia> map(a) do x
           2x + 1
       end
5-element Vector{Int64}:
  3
  5
  7
  9
 11

map() を使わない方法もある。

julia> (x -> 2x + 1).(a)
5-element Vector{Int64}:
  3
  5
  7
  9
 11

map() で cumsum()、cumprod()、diff() を表現できる。

julia> map((s -> v -> s += v)(0), [1, 2, 3, 4, 5])
5-element Vector{Int64}:
 1
 3
 6
 10
 15

julia> map((s -> v -> s *= v)(1), [1, 2, 3, 4, 5])
5-element Vector{Int64}:
  1
  2
  6
 24
 120

julia> map(((a, b) -> v -> (b = a; a = v; a - b))(0, 0), cumprod([1, 2, 3, 4, 5]))[2:end]
4-element Vector{Int64}:
 1
 4
 18
 96

マップ後の演算の積算

map() 後に reduce() したい場合があるとする。

julia> reduce(+, map(x -> 2x + 1, 1:5))
35

これを行うための mapreduce() がある。

julia> mapreduce(x -> 2x + 1, +, 1:5)
35

同様のものとして mapfoldl()、mapfoldr() もある。

julia> mapfoldl(x -> 2x + 1, =>, 1:5)
(((3 => 5) => 7) => 9) => 11

julia> mapfoldr(x -> 2x + 1, =>, 1:5)
3 => (5 => (7 => (9 => 11)))

単純な適用

単純に各要素に関数を適用したい場合は、for を使えばよい。

julia> b = []
Any[]

julia> for x in 1:5
           push!(b, 2x + 1)
       end

julia> b
5-element Vector{Any}:
  3
  5
  7
  9
 11

これと同じことを行うのに、foreach() がある。

julia> foreach(x -> push!(b, 2x + 1), 1:5)

julia> b
5-element Vector{Any}:
  3
  5
  7
  9
 11

map() との違いは、関数で値を返す必要がないところである。

上の例の場合は、内包表記でよい。

julia> b = [2x + 1 for x in 1:5]
5-element Vector{Int64}:
  3
  5
  7
  9
 11

抽出

filter() で条件に合った要素だけを取り出すことができる。

julia> a = collect(1:10)
10-element Vector{Int64}:
  1
  2
  3
  4
  5
  6
  7
  8
  9
 10

julia> filter(x -> x % 2 == 0, a)
5-element Vector{Int64}:
  2
  4
  6
  8
 10

map() を使って Bool 型の添字リストを作る方法もある。

julia> a[map(x -> x % 2 == 0, a)]
5-element Vector{Int64}:
  2
  4
  6
  8
 10

もっと直接的な書き方もできる。

julia> a[a .% 2 .== 0]
5-element Vector{Int64}:
  2
  4
  6
  8
 10

条件に合う場所の添字がほしい場合は、findall() を用いる。

julia> a = collect(2:10)
9-element Vector{Int64}:
  2
  3
  4
  5
  6
  7
  8
  9
 10

julia> findall(x -> x % 2 == 0, a)
5-element Vector{Int64}:
 1
 3
 5
 7
 9

julia> findall(a .% 2 .== 0)
5-element Vector{Int64}:
 1
 3
 5
 7
 9

条件判定

どれかの要素で条件が成り立つかを調べるには any() が使える。また、すべての要素で条件が成り立つかどうかは all() でわかる。

julia> a = [1, 3, 5, 8, 9]
5-element Vector{Int64}:
 1
 3
 5
 8
 9

julia> any(x -> x % 2 == 0, a)
true

julia> all(x -> x % 2 == 1, a)
false

julia> a = [1, 3, 5, 7, 9]
5-element Vector{Int64}:
 1
 3
 5
 7
 9

julia> all(x -> x % 2 == 1, a)
true

条件に合うものが何個あるかを数えるには、count() を使う。

julia> a = [1, 3, 5, 8, 9, 10]
6-element Vector{Int64}:
  1
  3
  5
  8
  9
 10

julia> count(x -> x % 2 == 0, a)
2

要素の置き換え

replace() で条件に合った要素を別の値で置き換えることができる。

julia> a = replace(x -> x % 2 == 0 ? missing : x, Any[1, 2, 3, 4, 5])
5-element Vector{Any}:
 1
  missing
 3
  missing
 5

また、特定の値を別の特定の値に置き換えることもできる。

julia> replace(a, missing => 0)
5-element Vector{Any}:
 1
 0
 3
 0
 5

replace() は map() と似ているが、若干結果が異なる。

julia> a = map(x -> x % 2 == 0 ? missing : x, Any[1, 2, 3, 4, 5])
5-element Vector{Union{Missing, Int64}}:
 1
  missing
 3
  missing
 5

この場合、もっと直接的に以下のようにもできる。

julia> a = Vector{Any}([1, 2, 3, 4, 5]);

julia> a[a .% 2 .== 0] .= missing;

julia> a
5-element Vector{Any}:
 1
  missing
 3
  missing
 5

julia> b = copy(a);

julia> b[ismissing.(b)] .= 0;

julia> b
5-element Vector{Any}:
 1
 0
 3
 0
 5

行列

大きさだけ確保する

とりあえず行列の大きさだけ確保するには、次のようにする。

julia> a = Array{Float64}(undef, 2, 3)
2×3 Matrix{Float64}:
 0.0  0.0  0.0
 0.0  0.0  0.0

zeros() でもよい。

julia> a = zeros(2, 3)
2×3 Matrix{Float64}:
 0.0  0.0  0.0
 0.0  0.0  0.0

fill() でもよい。

julia> a = fill(0., 2, 3)
2×3 Matrix{Float64}:
 0.0  0.0  0.0
 0.0  0.0  0.0

多次元配列がほしければ、さらに引数を指定すればよい。

julia> a = Array{Float64}(undef, 2, 3, 2)
2×3×2 Array{Float64, 3}:
[:, :, 1] =
 0.0  0.0  0.0
 0.0  0.0  0.0

[:, :, 2] =
 0.0  0.0  0.0
 0.0  0.0  0.0

ある値で埋める

全要素を特定の値にしたい場合などは、fill!() が使える。

julia> a = [1 2; 3 4]
2×2 Matrix{Int64}:
 1  2
 3  4

julia> fill!(a, 0)
2×2 Matrix{Int64}:
 0  0
 0  0

以下のようにしてもよい。

julia> a .= 0
2×2 Matrix{Int64}:
 0  0
 0  0

形状を変える

reshape() で形状を変えられる。

julia> [1:3 4:6]
3×2 Matrix{Int64}:
 1  4
 2  5
 3  6

julia> a = reshape([1:3 4:6], 1, 6)
1×6 Matrix{Int64}:
 1  2  3  4  5  6

julia> [1:6;]
6-element Vector{Int64}:
 1
 2
 3
 4
 5
 6

julia> a = reshape([1:6;], 3, 2)
3×2 Matrix{Int64}:
 1  4
 2  5
 3  6

[1:6;] の ; は、ベクトルにするために入れている。

ベクトル化するだけであれば、vec() が使える。

julia> a = [1 2; 3 4]
2×2 Matrix{Int64}:
 1  2
 3  4

julia> vec(a)
4-element Vector{Int64}:
 1
 3
 2
 4

結合する

行列を結合するには、以下のようにすればよい。

julia> a = [1 2; 3 4]
2×2 Matrix{Int64}:
 1  2
 3  4

julia> b = [5 6; 7 8]
2×2 Matrix{Int64}:
 5  6
 7  8

julia> [a; b]
4×2 Matrix{Int64}:
 1  2
 3  4
 5  6
 7  8

julia> [a b]
2×4 Matrix{Int64}:
 1  2  5  6
 3  4  7  8

行列の部分への参照を得る

@view で行列の部分への参照を得られる。

julia> a = [1 2; 3 4]
2×2 Matrix{Int64}:
 1  2
 3  4

julia> b = @view a[:, 1]
2-element view(::Matrix{Int64}, :, 1) with eltype Int64:
 1
 3

これで b は a の 1 列目を指している。b の値を変更すると a の値も変わる。

julia> b .= 0
2-element view(::Matrix{Int64}, :, 1) with eltype Int64:
 0
 0

julia> a
2×2 Matrix{Int64}:
 0  2
 0  4

ループ

行ごとのループ

julia> for x in eachrow(a)
           @show x
       end
x = [1, 2]
x = [3, 4]

列ごとのループ

julia> for x in eachcol(a)
           @show x
       end
x = [1, 3]
x = [2, 4]

要素ごとのループ

julia> for x in a
           @show x
       end
x = 1
x = 3
x = 2
x = 4

julia> for i in eachindex(a)
           @show a[i]
       end
a[i] = 1
a[i] = 3
a[i] = 2
a[i] = 4

要素ごとのループの順番からすると、列方向が自然な動きのようである。添字 (i, j) でアクセスするなら、次のようになる。

julia> for j in 1:size(a, 2), i in 1:size(a, 1)
           @show a[i, j]
       end
a[i, j] = 1
a[i, j] = 3
a[i, j] = 2
a[i, j] = 4

転置

行列の転置は引用符 (') で得られる。

julia> a = [1 2; 3 4]
2×2 Matrix{Int64}:
 1  2
 3  4

julia> a'
2×2 adjoint(::Matrix{Int64}) with eltype Int64:
 1  3
 2  4

これで得られるものは、正確には随伴行列である。

julia> a = [1 + 2im; 3 + 4im]
2-element Vector{Complex{Int64}}:
 1 + 2im
 3 + 4im

julia> a'
1×2 adjoint(::Vector{Complex{Int64}}) with eltype Complex{Int64}:
 1-2im  3-4im

adjoint() でも同じ結果が得られる。

julia> adjoint(a)
1×2 adjoint(::Vector{Complex{Int64}}) with eltype Complex{Int64}:
 1-2im  3-4im

本当の意味での転置は transpose() で得られる。

julia> transpose(a)
1×2 transpose(::Vector{Complex{Int64}}) with eltype Complex{Int64}:
 1+2im  3+4im

複素共役行列

conj() で複素共役行列が得られる。

julia> a = [1 + 2im; 3 + 4im]
2-element Vector{Complex{Int64}}:
 1 + 2im
 3 + 4im

julia> conj(a)
2-element Vector{Complex{Int64}}:
 1 - 2im
 3 - 4im

行列同士の演算

行列同士の和は + で計算できる。

julia> a = [1 2; 3 4]
2×2 Matrix{Int64}:
 1  2
 3  4

julia> b = [5 6; 7 8]
2×2 Matrix{Int64}:
 5  6
 7  8

julia> a + b
2×2 Matrix{Int64}:
  6   8
 10  12

アスタリスク (*) だと行列の積になる。

julia> a*b
2×2 Matrix{Int64}:
 19  22
 43  50

要素同士の掛け算をしたい場合は、.* を用いる。

julia> a.*b
2×2 Matrix{Int64}:
  5  12
 21  32

累乗も行列の累乗になる。

julia> a = [1 2 3; 4 5 6; 7 9 8]
3×3 Matrix{Int64}:
 1  2  3
 4  5  6
 7  9  8

julia> a^2
3×3 Matrix{Int64}:
 30   39   39
 66   87   90
 99  131  139

julia> a^-2
3×3 Matrix{Float64}:
  3.74074   -3.85185    1.44444
 -3.25926    3.81481   -1.55556
  0.407407  -0.851852   0.444444

線形代数ライブラリ

標準ライブラリのひとつとして、線形代数ライブラリ LinearAlgebra がある。次のように使う。

using LinearAlgebra

単位行列

julia> a = [1 2 3; 4 5 6; 7 8 9]
3×3 Matrix{Int64}:
 1  2  3
 4  5  6
 7  9  8

julia> a*I
3×3 Matrix{Int64}:
 1  2  3
 4  5  6
 7  9  8

julia> I(3)
3×3 Diagonal{Bool, Vector{Bool}}:
 1  ⋅  ⋅
 ⋅  1  ⋅
 ⋅  ⋅  1

ノルム

julia> a = [1, 2, 3]
3-element Vector{Int64}:
 1
 2
 3

julia> norm(a)
3.7416573867739413

正規化

julia> normalize(a)
3-element Vector{Float64}:
 0.2672612419124244
 0.5345224838248488
 0.8017837257372732

内積

julia> a = [1, 2, 3]
3-element Vector{Int64}:
 1
 2
 3

julia> b = [4, 5, 6]
3-element Vector{Int64}:
 4
 5
 6

julia> dot(a, b)
32

外積

julia> cross(a, b)
3-element Vector{Int64}:
 -3
  6
 -3

トレース

julia> a = [1 2 3; 4 5 6; 7 9 8]
3×3 Matrix{Int64}:
 1  2  3
 4  5  6
 7  9  8

julia> tr(a)
14

対角成分

julia> diag(a)
3-element Vector{Int64}:
 1
 5
 8

行列式

julia> det(a)
8.999999999999996

逆行列

julia> inv(a)
3×3 Matrix{Float64}:
 -1.55556    1.22222   -0.333333
  1.11111   -1.44444    0.666667
  0.111111   0.555556  -0.333333

固有値

julia> a = [1 2; 3 4]
2×2 Matrix{Int64}:
 1  2
 3  4

julia> f = eigen(a)
Eigen{Float64, Float64, Matrix{Float64}, Vector{Float64}}
values:
2-element Vector{Float64}:
 -0.3722813232690143
  5.372281323269014
vectors:
2×2 Matrix{Float64}:
 -0.824565  -0.415974
  0.565767  -0.909377

julia> f.values
2-element Vector{Float64}:
 -0.3722813232690143
  5.372281323269014

julia> f.vectors
2×2 Matrix{Float64}:
 -0.824565  -0.415974
  0.565767  -0.909377

julia> eigvals(a)
2-element Vector{Float64}:
 -0.3722813232690143
  5.372281323269014

julia> eigvecs(a)
2×2 Matrix{Float64}:
 -0.824565  -0.415974
  0.565767  -0.909377

分解

julia> a = [1 2 3; 4 5 6; 7 9 8]
3×3 Matrix{Int64}:
 1  2  3
 4  5  6
 7  9  8

julia> factorize(a)
LU{Float64, Matrix{Float64}, Vector{Int64}}
L factor:
3×3 Matrix{Float64}:
 1.0        0.0  0.0
 0.142857   1.0  0.0
 0.571429  -0.2  1.0
U factor:
3×3 Matrix{Float64}:
 7.0  9.0       8.0
 0.0  0.714286  1.85714
 0.0  0.0       1.8

factorize() でどういう手法が用いられるかは、与えた行列の種類による。

LU 分解

julia> f = lu(a)
LU{Float64, Matrix{Float64}, Vector{Int64}}
L factor:
3×3 Matrix{Float64}:
 1.0        0.0  0.0
 0.142857   1.0  0.0
 0.571429  -0.2  1.0
U factor:
3×3 Matrix{Float64}:
 7.0  9.0       8.0
 0.0  0.714286  1.85714
 0.0  0.0       1.8

julia> f.L
3×3 Matrix{Float64}:
 1.0        0.0  0.0
 0.142857   1.0  0.0
 0.571429  -0.2  1.0

julia> f.U
3×3 Matrix{Float64}:
 7.0  9.0       8.0
 0.0  0.714286  1.85714
 0.0  0.0       1.8

Cholesky 分解

julia> a = [4 12 -16; 12 37 -43; -16 -43 98]
3×3 Matrix{Int64}:
   4   12  -16
  12   37  -43
 -16  -43   98

julia> f = cholesky(a)
Cholesky{Float64, Matrix{Float64}}
U factor:
3×3 UpperTriangular{Float64, Matrix{Float64}}:
 2.0  6.0  -8.0
  ⋅   1.0   5.0
  ⋅    ⋅    3.0

QR 分解

julia> a = [12 -51 4; 6 167 -68; -4 24 -41]
3×3 Matrix{Int64}:
 12  -51    4
  6  167  -68
 -4   24  -41

julia> f = qr(a)
LinearAlgebra.QRCompactWY{Float64, Matrix{Float64}, Matrix{Float64}}
Q factor:
3×3 LinearAlgebra.QRCompactWYQ{Float64, Matrix{Float64}, Matrix{Float64}}:
 -0.857143   0.394286   0.331429
 -0.428571  -0.902857  -0.0342857
  0.285714  -0.171429   0.942857
R factor:
3×3 Matrix{Float64}:
 -14.0   -21.0   14.0
   0.0  -175.0   70.0
   0.0     0.0  -35.0

特異値分解

julia> a = [1 0 0 0 2; 0 0 3 0 0; 0 0 0 0 0; 0 2 0 0 0]
4×5 Matrix{Int64}:
 1  0  0  0  2
 0  0  3  0  0
 0  0  0  0  0
 0  2  0  0  0

julia> svd(a)
SVD{Float64, Float64, Matrix{Float64}, Vector{Float64}}
U factor:
4×4 Matrix{Float64}:
 0.0  1.0   0.0  0.0
 1.0  0.0   0.0  0.0
 0.0  0.0   0.0  1.0
 0.0  0.0  -1.0  0.0
singular values:
4-element Vector{Float64}:
 3.0
 2.23606797749979
 2.0
 0.0
Vt factor:
4×5 Matrix{Float64}:
 -0.0        0.0  1.0  -0.0  0.0
  0.447214   0.0  0.0   0.0  0.894427
  0.0       -1.0  0.0   0.0  0.0
  0.0        0.0  0.0   1.0  0.0

連立方程式を解く

julia> a = [1 2 3; 4 5 6; 7 9 8]
3×3 Matrix{Int64}:
 1  2  3
 4  5  6
 7  9  8

julia> b = [1; 2; 3]
3-element Vector{Int64}:
 1
 2
 3

julia> x = a\b
3-element Vector{Float64}:
 -0.11111111111111116
  0.22222222222222232
  0.22222222222222215

julia> a*x
3-element Vector{Float64}:
 1.0
 2.0
 3.0

時間の計測

関数の実行時間を計るには、@time を用いる。

julia> @time fibonacci(30)
  0.005647 seconds
832040

プロファイリングするには、Profile を用いる。

julia> using Profile

julia> @profile fibonacci(30)
832040

julia> Profile.print(format=:flat)
 Count  Overhead File                    Line Function
 =====  ======== ====                    ==== ========
    35         0 @Base/boot.jl            370 eval
    35         0 @Base/client.jl          421 (::Base.var"#1017#1019"{Bool, B...
    35         0 @Base/client.jl          522 _start()
    35         0 @Base/client.jl          322 exec_options(opts::Base.JLOptions)
    35         0 @Base/client.jl          405 run_main_repl(interactive::Bool...
    35         0 @Base/essentials.jl      816 #invokelatest#2
    35         0 @Base/essentials.jl      813 invokelatest
     1         1 @Base/int.jl              86 -
    35         0 ....9/REPL/src/REPL.jl   379 run_repl(repl::REPL.AbstractREP...
    35         0 ....9/REPL/src/REPL.jl   234 start_repl_backend(backend::REP...
    35         0 ....9/REPL/src/REPL.jl   153 eval_user_input(ast::Any, backe...
    35         0 ....9/REPL/src/REPL.jl   249 repl_backend_loop(backend::REPL...
    35         0 ....9/REPL/src/REPL.jl   365 run_repl(repl::REPL.AbstractREP...
    35         0 ....9/REPL/src/REPL.jl   231 kwcall(::Any, ::typeof(REPL.sta...
     6         6 .../julia/fibonacci.jl     ? fibonacci(n::Int64)
    17        17 .../julia/fibonacci.jl     1 fibonacci(n::Int64)
     3         3 .../julia/fibonacci.jl     2 fibonacci(n::Int64)
    34         7 .../julia/fibonacci.jl     7 fibonacci(n::Int64)
Total snapshots: 36. Utilization: 100% across all threads and tasks. Use the `groupby` kwarg to break down by thread and/or task.

並列計算

マルチスレッド

Julia にスレッド数を指定すれば、マルチスレッドにできる。

$ julia -t 4

...

julia> Threads.nthreads()
4

スレッド並列を行うには、for ループに Threads.@threads をつけるか、関数を Threads.@spawn で呼ぶ方法がある。たとえば、配列の和を取ることを考える。

threads.jl

const n = 10000000
const m = Threads.nthreads()
const dn = round(Int64, n/m)
const a = collect(1:n)

function f()
    s = 0
    for i = 1:m
        j1 = dn*(i - 1) + 1
        j2 = dn*i
        for j = j1:j2
            s += a[j]
        end
    end
    s
end

function ft()
    s = Threads.Atomic{Int64}(0)
    Threads.@threads for i = 1:m
        j1 = dn*(i - 1) + 1
        j2 = dn*i
        for j = j1:j2
            Threads.atomic_add!(s, a[j])
        end
    end
    s.value
end

function ft2()
    s = zeros(Int64, m)
    Threads.@threads for i = 1:m
        j1 = dn*(i - 1) + 1
        j2 = dn*i
        for j = j1:j2
            s[i] += a[j]
        end
    end
    sum(s)
end

function ft3()
    tasks = map(1:m) do i
        Threads.@spawn (i -> begin
            s = 0
            j1 = dn*(i - 1) + 1
            j2 = dn*i
            for j = j1:j2
                s += a[j]
            end
            s
        end)(i)
    end
    sum(fetch.(tasks))
end

f() は並列計算しないもの、ft()、ft2() は Threads.@threads の for ループで並列計算を行うもの、ft3() は Threads.@spawn で並列で関数を呼ぶものである。ft() では Threads.Atomic などを用いているが、並行して同じ変数に足し込むので、ふつうの変数だと競合が起こるためである。ft2() では配列を用いて競合を避けている。

julia> include("threads.jl")

julia> for i = 1:3; @time f(); end
  0.004412 seconds
  0.004410 seconds
  0.004435 seconds

julia> for i = 1:3; @time ft(); end
  0.174059 seconds (11.32 k allocations: 830.396 KiB, 41.33% compilation time)
  0.183646 seconds (26 allocations: 2.328 KiB)
  0.157423 seconds (28 allocations: 2.391 KiB)

julia> for i = 1:3; @time ft2(); end
  0.025975 seconds (11.56 k allocations: 818.498 KiB, 316.88% compilation time)
  0.005410 seconds (26 allocations: 2.406 KiB)
  0.005359 seconds (27 allocations: 2.438 KiB)

julia> for i = 1:3; @time ft3(); end
  0.045985 seconds (55.50 k allocations: 3.910 MiB, 161.96% compilation time)
  0.003005 seconds (38 allocations: 2.578 KiB)
  0.002870 seconds (38 allocations: 2.578 KiB)

最初はコンパイル時間が入るため、3 回計算時間を計測している。ft() は遅くなっているが、Atomic のせいだろうか。ft2() は少しましになっているが、まだ少し遅い。ft3() では少し速くなっている。

分散計算

Julia にプロセス数を指定すれば、分散計算できる。

$ julia -p 4

...

julia> nworkers()
4

プロセス数を指定すると、ライブラリの Distributed が有効になる。nworkers() でワーカーの数がわかる。似たもので nprocs() があるが、こちらは指定のプロセス数より 1 多いので、おそらく管理プロセスみたいなのを含む。

配列の和を取ることを考える。

procs.jl

const n = 10000000
const m = nworkers()
const dn = round(Int64, n/m)
const a = collect(1:n)

function f()
    s = 0
    for i = 1:m
        j1 = dn*(i - 1) + 1
        j2 = dn*i
        for j = j1:j2
            s += a[j]
        end
    end
    s
end

function fp()
    @distributed (+) for i = 1:m
        s = 0
        j1 = dn*(i - 1) + 1
        j2 = dn*i
        for j = j1:j2
            s += a[j]
        end
        s
    end
end

function fp2()
    tasks = map(1:m) do i
        @spawnat :any (i -> begin
            s = 0
            j1 = dn*(i - 1) + 1
            j2 = dn*i
            for j = j1:j2
                s += a[j]
            end
            s
        end)(i)
    end
    sum(fetch.(tasks))
end

println("loaded")

f() は並列計算しないもの、fp() は @ditributed (+) for ループで足し合わせるもの、fp2() は @spawnat で並列に関数を呼び出すものである。

この場合、include() 自体も特殊である。

julia> @everywhere include("procs.jl")
loaded
      From worker 3:	loaded
      From worker 5:	loaded
      From worker 4:	loaded
      From worker 2:	loaded

@everywhere ですべてのワーカーで include() する。

julia> for i = 1:3; @time f(); end
  0.004352 seconds
  0.004292 seconds
  0.004445 seconds

julia> for i = 1:3; @time fp(); end
  0.134823 seconds (12.01 k allocations: 824.546 KiB, 11.72% compilation time)
  0.007191 seconds (292 allocations: 13.438 KiB)
  0.007004 seconds (291 allocations: 13.141 KiB)

julia> for i = 1:3; @time fp2(); end
  0.018105 seconds (1.28 k allocations: 75.156 KiB)
  0.003484 seconds (484 allocations: 20.359 KiB)
  0.003252 seconds (485 allocations: 20.656 KiB)

fp() は遅くなっている。fp2() は少し速くなっている。

ファイル

ファイルを開いたり閉じたりするには open()、close() を使う。

julia> f = open("test.txt", "w")
IOStream()

julia> write(f, "test")
4

julia> close(f)

julia> f = open("test.txt", "r")
IOStream()

julia> readline(f)
"test"

julia> close(f)

あるいは、次のような書き方もできる。

julia> open("test.txt", "r") do f
           readline(f)
       end
"test"

全体を読み込むには readlines() が使える。

julia> open("test.txt", "w") do f
           write(f, "aaa\n")
           write(f, "bbb\n")
           write(f, "ccc\n")
       end
4

julia> lines = open("test.txt", "r") do f
           readlines(f)
       end
3-element Vector{String}:
 "aaa"
 "bbb"
 "ccc"

1行ずつ読み込んでファイルの最後まで処理するなら、次のような方法もある。

julia> open("test.txt", "r") do f
           for line in eachline(f)
               println(line)
           end
       end
aaa
bbb
ccc

行列を読み書きすることもできる。

julia> a = [1 2; 3 4]
2×2 Matrix{Int64}:
 1  2
 3  4

julia> open("test.txt", "w") do f
           write(f, a)
       end
32

julia> a = zeros(Int64, 2, 2)
2×2 Matrix{Int64}:
 0  0
 0  0

julia> open("test.txt", "r") do f
           read!(f, a)
       end
2×2 Matrix{Int64}:
 1  2
 3  4

ただし、バイナリデータになる。

書式

Printf パッケージで @printf などが使える。C 言語の printf() のような書式が使えるようになる。

julia> using Printf

julia> a = 1.23456
1.23456

julia> @printf "a = %.3f" a
a = 1.235

@sprintf もある。

julia> s = @sprintf "a = %.3f" a
"a = 1.235"