C87 お疲れ様でした

コミケに参加した皆様、お疲れ様でした。 3日目に私のサークルに来てくださった方々はありがとうございます。 想定よりちょっと多く捌けました。

次回以降

次回は外部ライブラリの利用の話を軽くした後、 マクロ、もしくは並列化のどちらかを解説するつもりです。 毎回一冊にマージしちゃうとどんどん分厚くなって、 しかもそれを毎回買わせるというのも流石に気がひけるので、(一段落つくまで)マージはしません。 今回の本が改めて1巻、次回が2巻、その次が3巻、...、と続いていきます。

今回の本についてはまだ在庫があるので、次回以降はこれを刷り増します。 本そのものの訂正や、あともちろんJulia のアップデートへの追従も必要ですが、 これらは別途小冊子/web 形式で頒布することで対応しようかと思います。

以下今回の話や反省点など。

今覚えている、現地で受けた質問とその時の解答、そしてそれに対する今思った補足

  • 特徴を一言で言うと?
    • 「ぼくのかんがえた最強の言語」ですかねぇ……
      • 無茶ぶりに近いけれど、答えとしてはあながち間違ってはいないかなと。
      • 元の意味を考えるとネガティブ寄りだけれど、まぁそんな最強言語実在するとは思ってないですしー(ぇー
  • 関数型なの?
    • 少なくとも関数は第一級の値だし、無名関数リテラルを使って動的に作れるけれど、それだけでいわゆる関数型って言っていいんですかね?
      • まぁ関数型プログラミング自体は(きっと)できる。高階関数とか部分適用とかで無名関数を使いまくった結果のプログラムが速いかどうかは別として。
      • とりあえず、関数ごとにJIT コンパイルが効くので*1、いろいろ関数を作ってそれを組み合わせろ、という事にはなっている
  • どんなものが得意?
    • 割となんでもこなせるから逆にこれってものはないけれど、スクリプト言語なので、ゲームとかはインタプリタ込で配布しないといけないから面倒臭いんじゃない?あと起動時間が他のスクリプト言語の数倍から十数倍はあるので、シェルのユーティリティとかフィルタとかみたいな実行時間がそもそも短いのには向かない
      • 機械学習とかは、ユーザが多いから目立っているだけで、特にとりたてて得意だというわけではないんじゃないかなと。苦手なものを挙げるしかない
      • 仕様が固まっていないので、広く長く使うプログラムにはまだまだ辛いかもー

(過去と比べて相対的に)よさげだったところ

表紙絵の導入 (Julia-tan)

やっぱりワンポイントでも表紙絵がつくと違うなぁと。 絵を用いた表紙デザインがまっとうにできているかどうかはわかりませんが……

1からまとめたこと

去年の夏に出した1巻、完全な入門部分は無事完売して、 一方で冬に出した2巻であるJulia の型の話はあんまりでした。 当時はコピー製本だったので、1巻をまた一から作るのが超絶面倒臭かった*2 ので2巻だけ持っていったのです。 Julia なんて正直いって(少なくとも去年の冬の段階では)文法どころか名前すらも知らない人が多くて、そんな状況で入門すっ飛ばした2巻だけ売られても買わねーよ、っていうのが最大の原因だったかなぁと。

そんなわけで今回は1巻も刷り直そうとして、せっかくだからと両方一冊に合わせて印刷・製本を外注したのでした。 1つにまとめ直して最新版に追従しただけで印刷所に伝えたページ数にピッタリ到達しちゃったよ……

POP の改良

正直言って人見知りよりのコミュ障気味なので、こちらから声かけるのが辛いのです。 今回はPOP にQ&A *3を書いた事によって、 最低限の情報を言うまでもなく伝えることで面倒臭さ(ぇー)を軽減できたかなぁと。 あと自分以外に店番頼んだ時にも役立つし。 ただし、これだけだと相互コミュニケーションにならないので、誤解を生む可能性があるのが難点かな。

反省点

「(科学)技術計算」

いつも特に何も考えずにふつーに使っちゃってるんですが、「科学技術計算」って今では何を指しているんですかね。 方程式解いたり数値シミュレーションしたり*4でいいんでしょうか。 少なくとも今年の日本では、データ・統計処理とか機械学習の分野でJulia が注目されているので、そのあたりを押した方がキャッチーかなぁ、と思いつつも、自分の専門分野でもないのでつらいところ。

参考資料について

参考資料2冊(Julia を扱った市販本)も並べておいておいたけれど、 "Seven More ..." の方は特に話題にもならず、また「データサイエンティスト ... 」の方は市販本であることに気づかなかった人もそれなりにいたようです。 説明など何の工夫もせずにおいただけなのが悪いというのもあるけれど、 あまり余計なことはしないほうがいいかな……。

*1:これ自体は関数型プログラミングどうこうとは関係ないとは思うけれど

*2:あと単純に2巻を刷りすぎた

*3:Julia ってなに / どんな言語なの

*4:この2つだいたい同じだけれど

C87 情報

JuliaLang Advent Calendar 2014 は楽しかったですね(挨拶

明日からビッグサイトコミックマーケットが開催されますが、1年ぶりにJulia 本を出します。

3日目(火)西 く-14b EREX工房

です。去年の2冊を1冊にしてv0.4 向けに加筆修正した感じです。 B5、60ページで500円です。 表紙のJulia-tan が目印です。 今回ついに印刷・製本を外注して、現物は当日現地に届くので、どうなっているのかちょっとドキドキですね。 加筆修正とはいうものの、ひとつにまとめて整理しただけで予定したページ数に到達してしまったので本質的な加筆はありません……いいかげんマクロどうにかしたいですね……

ついでに、ここ1,2ヶ月の間に出た*1 Julia 関連の書籍2冊(データサイエンティスト養成読本 R活用編, Seven More Languages in Seven Weeks のうちJulia 部分 )も持っていくので、中身をちょっと見てみたいという人もどうぞ。私がスペースに居るとき*2、かつ周りにあまり迷惑のかからないぐらい空いている時、だけですが。 ちなみにこの二冊のJulia 部分、結構長くて*3、立ち読みだけでは済まないので、興味持ったら買うといいですよ。電子版もあるし。

当日はよろしくお願いいたします。

*1:むしろ今まで出てなかったけれど

*2:11時までと、13-15 時ぐらいには多分いるんじゃないかと思いますが、ふらふらと買い物に出ているかも。

*3:後者はもちろん、前者も思ったより長かった

ifelse() と関数の分離による高速化 -- Base.randn() を題材にして --

JuliaLang Advent Calendar 2014 の 18 日目の記事です。遅れてすみません。 今回も速度向上の記事です。

概要

Base.randn() に最近なされた速度向上のための更新(#9126, #9132)を通じて、 ifelse() 関数や関数分離の有用性を示します。

サンプルコード

https://gist.github.com/yomichi/62d5f121ab11831b0759

Base.randn() について

正規分布はホワイトノイズ の分布として誤差解析によく使われ、 中心極限定理によって平均値の分布にもなり、 さらには拡散方程式の解として直接現れるなど、 自然科学や工学の分野では一様分布と同じぐらい重要な分布となっています。 また、手で積分の計算ができるほどに非常に性質が良いため、 様々なランダム性を含む理論モデルで、そのランダム性の分布として使われています。 モンテカルロ法などの数値シミュレーションにもよく使われるので、 正規分布に従うガウス乱数を高速に生成することはとても大事になります。

Base.randn() は 平均 μ = 0, 分散 σ^2 = 1正規分布 N(0,1) に従うガウス乱数を生成する関数です。 任意の平均、分散のガウス乱数を使いたい場合、μ + σ * randn() で変換できます。

ガウス乱数を生成するアルゴリズムとして、Box-Muller 法が(簡単なので)有名ですが、 Julia ではより洗練されて高速な Ziggurat 法 が用いられています *1四辻「計算機シミュレーションのための確率分布乱数生成法」によると*2、 Ziggurat 法は Box-Muller のおよそ5-6 倍速いようです。 アルゴリズムの詳細は本記事の範囲ではないので、 興味のある方はWikipedia や原著論文、上記教科書を参照してください。 もちろん、Julia の実装を読み解くというのもありだと思います:D

最近、#9126#9132 によって、 Base.randn() が更新されました。 これは恐ろしく単純な変更なのですが、 環境にもよりますが、2倍近くまで加速しています。 ここで使われた考えは決して Base.randn() に特有なものではなく、 普段のパフォーマンスチューニングにも役立つものとなっているので*3、本記事ではそれを紹介します。

初期バージョン

Ziggurat 法には必要な定数や関数がいくつかあるのですが、 あまりにも多いので Base.Random モジュールの中から拝借してきます *4

## 必要なパラメータや関数を持ってくる
import Base.Random: GLOBAL_RNG, rand_ui52
import Base.Random: ki, wi, fi, ke, we, fe  # それぞれ256 個のFloat64 の配列
import Base.Random: ziggurat_nor_r, ziggurat_nor_inv_r, ziggurat_exp_r

さて、加速する前の最初のバージョンは次のとおりです。

## 初期バージョン
@inline function randn_first(rng::AbstractRNG=GLOBAL_RNG)
    @inbounds begin
        r = rand_ui52(rng)
        rabs = int64(r>>1)
        idx = rabs & 0xFF

        # Point 1
        # 分岐命令
        x = (r % Bool ? -rabs : rabs)*wi[idx+1]

        rabs < ki[idx+1] && return x # 99.3% の確率で、ここでreturn

        # Point 2
        # ここに来る確率は 0.7 %
        @inbounds if idx == 0
          while true
            xx = -ziggurat_nor_inv_r*log(rand(rng))
            yy = -log(rand(rng))
            yy+yy > xx*xx && return (rabs >> 8) % Bool ? -ziggurat_nor_r-xx : ziggurat_nor_r+xx
          end
        elseif (fi[idx] - fi[idx+1])*rand(rng) + fi[idx+1] < exp(-0.5*x*x)
          return x # return from the triangular area
        else
          return randn_first(rng)
        end
    end
end

アルゴリズムの詳細は重要ではないので述べません。 この randn_first()N = 10000*10000 回呼び出すのにかかる時間を計測すると、

julia> @time for i in 1:N; randn_first(); end
elapsed time: 1.835435398 seconds (0 bytes allocated)

となりました*5。 加速にあたって改造する場所は2箇所、Point 1 とPoint 2 です。

Point 1 : ifelse(c :: Bool, x, y)

まずPoint 1 ですが、ここでは三項条件演算子による分岐が生じます *6。 ここを次のように ifelse(c :: Bool ,x, y) で 書き換えます。

## Point 1:
## 分岐命令(if ... end / ?: ) をifelse() 関数に置き換え

@inline function randn_ifelse(rng::AbstractRNG=GLOBAL_RNG)
    @inbounds begin
        r = rand_ui52(rng)
        rabs = int64(r>>1)
        idx = rabs & 0xFF

        # Point 1
        # ifelse() 関数(分岐命令なしで値が選択される!)
        x = ifelse(r % Bool, -rabs, rabs)*wi[idx+1]

        rabs < ki[idx+1] && return x # 99.3% の確率で、ここでreturn

        @inbounds if idx == 0
          while true
            xx = -ziggurat_nor_inv_r*log(rand(rng))
            yy = -log(rand(rng))
            yy+yy > xx*xx && return (rabs >> 8) % Bool ? -ziggurat_nor_r-xx : ziggurat_nor_r+xx
          end
        elseif (fi[idx] - fi[idx+1])*rand(rng) + fi[idx+1] < exp(-0.5*x*x)
          return x # return from the triangular area
        else
          return randn_ifelse(rng)
        end
    end
end

これについて計算時間を測定すると、

julia> @time for i in 1:N; randn_ifelse(); end
elapsed time: 1.1361549 seconds (0 bytes allocated)

と、たったの1行を書き換えただけで、およそ6 割程度も速くなりました。

今回使った ifelse(c::Bool, x, y) 関数 は、基本的には if ... else ... end?: を再現した関数で、 c == true ならば x を、 c == false ならば y を返すのですが、

  1. 関数なので xy はあらかじめ評価される
  2. LLVM IR における select instructionコンパイルされる

という2点が他の構文と大きく異なります。 特に2つ目が非常に重要で、これにより ifelse 関数は(LLVM IR のレベルで)分岐命令を生成しなくなり、 特に最内ループで使うことで性能が向上する 可能性 があります。 これらの特徴は、生成されたLLVM IR を実際に眺めれば確認できて、 まず条件構文は

julia> condition_operator(N) = iseven(N) ? div(N,2) : 3N+1
condition_operator (generic function with 1 method)

julia> code_llvm(condition_operator, (Int,) )

define i64 @julia_condition_operator_64603(i64) {
top:
  %1 = and i64 %0, 1, !dbg !577
  %2 = icmp eq i64 %1, 0, !dbg !577
  br i1 %2, label %pass2, label %L, !dbg !577

pass2:                                            ; preds = %top
  %3 = sdiv i64 %0, 2, !dbg !577
  ret i64 %3, !dbg !577

L:                                                ; preds = %top
  %4 = mul i64 %0, 3, !dbg !577
  %5 = add i64 %4, 1, !dbg !577
  ret i64 %5, !dbg !577
}

となります。br が分岐命令で、変数 %2 の真偽で %pass2:%L: へとジャンプします。 返り値はジャンプしてから計算していることがわかりますね。 一方のifelse() 関数は

julia> condition_ifelse(N) = ifelse(iseven(N), div(N,2), 3N+1)
condition_ifelse (generic function with 1 method)

julia> code_llvm(condition_ifelse, (Int,) )

define i64 @julia_condition_ifelse_64604(i64) {
pass2:
  %1 = mul i64 %0, 3, !dbg !580
  %2 = add i64 %1, 1, !dbg !580
  %3 = and i64 %0, 1, !dbg !580
  %4 = icmp ne i64 %3, 0, !dbg !580
  %5 = sdiv i64 %0, 2, !dbg !580
  %6 = select i1 %4, i64 %2, i64 %5, !dbg !580
  ret i64 %6, !dbg !580
}

となります。こちらでは先に条件の真偽値 %4 と2通りの返り値%2%5 とを計算しておいて、 select にこれらを渡して結果を返しています *7

ifelse を使うと速くなるケースが割とあるのですが、 xy の両方の式をかならず評価しないといけないため、 稀にしか選ばれない方の式が重い場合、 例えばPoint 2 で示すような rand() < 0.993 ? x : heavy_function(x) などでは、逆に遅くなります。 また、再帰関数で使うとスタックオーバーフローで落ちますし、 さらには前段の条件によってはどちらかでエラーを吐くような場合、 例えば isdefined(:x) ? x : 0 なんかはifelse の呼び出しにすらたどりつかずにエラーで落ちます。 何事も適材適所があるのです。。。

Point 2 : 関数の分離

次にPoint 2 ですが、この直前にある

rabs < ki[idx+1] && return x # 99.3% の確率で、ここでreturn

のおかげで、この部分にはおよそ0.7 % という小さな確率でしか到達しません *8。 しかしながら、この稀にしか到達しない分岐には、まだまだ長い処理が待っています *9。 こういう時には、長い部分を別の関数に分けて、関数呼び出しに置き換えることで、 性能が向上する 場合が あります。

## Point 2:
## およそ0.7% でしか到達しないのにやたら長い部分を
## 別関数に分離

@inline function randn_separate(rng::AbstractRNG=GLOBAL_RNG)
    @inbounds begin
        r = rand_ui52(rng)
        rabs = int64(r>>1)
        idx = rabs & 0xFF
        x = (r % Bool ? -rabs : rabs)*wi[idx+1]
        rabs < ki[idx+1] && return x # 99.3% の確率で、ここでreturn

        # Point 2
        # ここに来る確率は0.7 %
        # 関数呼び出しに置き換えて関数の長さを減らす
        return separate_unlikely(rng, idx, rabs, x)
    end
end

# Point 2
# たまにしか使われないくせに長い部分
function separate_unlikely(rng, idx, rabs, x)
    @inbounds if idx == 0
        while true
            xx = -ziggurat_nor_inv_r*log(rand(rng))
            yy = -log(rand(rng))
            yy+yy > xx*xx && return (rabs >> 8) % Bool ? -ziggurat_nor_r-xx : ziggurat_nor_r+xx
        end
    elseif (fi[idx] - fi[idx+1])*rand(rng) + fi[idx+1] < exp(-0.5*x*x)
        return x # return from the triangular area
    else
        return randn_separate(rng)
    end
end

こいつの実行時間を計測すると

julia> @time for i in 1:N; randn_separate(); end
elapsed time: 1.629525055 seconds (0 bytes allocated)

となり、最初のバージョンから比べて1割程度速くなります *10

Point 1+2 : 現在のrandn()

上記2つの加速を合わせたのが次のコードで、実際にBase.randn() として使われているものになります。

## 20141218 現在のrandn()
## Point 1+2
@inline function randn_final(rng::AbstractRNG=GLOBAL_RNG)
    @inbounds begin
        r = rand_ui52(rng)
        rabs = int64(r>>1)
        idx = rabs & 0xFF

        # Point 1
        x = ifelse(r % Bool, -rabs, rabs)*wi[idx+1]

        rabs < ki[idx+1] && return x # 99.3% の確率で、ここでreturn

        # Point 2
        return randn_unlikely(rng, idx, rabs, x)
    end
end

function randn_unlikely(rng, idx, rabs, x)
    @inbounds if idx == 0
        while true
            xx = -ziggurat_nor_inv_r*log(rand(rng))
            yy = -log(rand(rng))
            yy+yy > xx*xx && return (rabs >> 8) % Bool ? -ziggurat_nor_r-xx : ziggurat_nor_r+xx
        end
    elseif (fi[idx] - fi[idx+1])*rand(rng) + fi[idx+1] < exp(-0.5*x*x)
        return x # return from the triangular area
    else
        return randn_final(rng)
    end
end

実行時間は

julia> @time for i in 1:N; randn_final(); end
elapsed time: 1.014320666 seconds (0 bytes allocated)

となり、最初と比べるとおよそ8割、2倍近くまで速度が向上します。

まとめ

ifelse() 関数による分岐の書き換えはうまくハマるとかなり速くなります。 もちろん大した効果の得られない時もありますし、場合によっては 遅くなることもありますが、 パフォーマンス測定をしながら程々に使っていくとよいでしょう。

関数の分解をすることで可読性や保守性、再利用性が高まりますが、 今回のように速度面でも有利になることがあるので、 がんがん切り分けていきましょう。 もちろん、あまりやり過ぎると可読性や保守性が落ちるし 関数呼び出しのオーバーヘッドも現れうる*11ので、 やはりパフォーマンス測定をしつつ、何事も程々に*12

*1:ちなみに他の言語やライブラリでは、GSL (C) やnumpy (Python) がBox-Muller, Boost (C++) がZiggurat でした。

*2:どこもかしこも在庫切れのようですが。。。

*3:それと単純に感動したので

*4:これらはexport されていないので、いつの日か名前が変わったり別のものになったりなくなったりするかもしれません。あまりマネしないようにしましょう:P

*5:計算時間は同時に走っている関係ないプロセスから かなり影響を受けるため、本当に正確に測定するならば、 できるだけ他のプロセスを殺してから測定したり、 その上で何度も測定して最小値を取ったりする必要があります。 ただし、今回はそこまでやっていないので、 速度向上の倍率などは参考程度に捉えてください。

*6: 大多数の人は見慣れないかと思いますが、 r % Bool isodd(r) と同じです。 生成されるLLVM IR も全く同じです。

*7:アセンブラは私がよくわかっていないので見ないことにします(ぇー

*8: ここでは && 演算子によるショートカットを利用しています。 生成されるLLVM IR は if 構文と同じで、 標準ライブラリにおいては一行でif 構文を書くのに好んで用いられるようです。

*9:つまり、Point 1 のようにifelse 関数を使うと確実に遅くなります。 再帰っぽいですが、再帰するかしないかは乱数で決まるので、一応帰ってきます。

*10:結構揺らぎますが。。。

*11:一応LLVM は自動インライン化もあるようですが

*12:大事なことなので2回言いました

NumericExtensions.jl / NumericFuns.jl

JuliaLang Advent Calendar の10日目の記事です。

Julia 言語を使う動機として、開発速度と実行速度を共に向上させたい、 というものがやはり大きくなるのではないでしょうか。 NumericFuns.jlNumericExtensions.jl は、 Julia の実行速度、特に配列を用いた数値計算の速度を更に向上させる事のできるパッケージです。

予備知識 -- 時間測定

実行速度を向上させる、というためには、実際に実行速度(実行時間)を計測する必要があります。 Julia では@time を使うことで簡単に実行時間を計測できます。

# @time は与えた式を評価して、それと同時に標準出力に計算時間と必要メモリを表示
julia> @time sleep(0.5)
elapsed time: 0.501986776 seconds (472 bytes allocated)

なお、Julia には他にもコードカバレッジ測定ツールやプロファイラなどが 標準で 用意されています。詳しくはJuliaTokyo #2 での 発表スライド をご参照ください*1

なお、Julia では関数単位のJIT コンパイルがあるために、それぞれの関数を予め(小さなデータでもいいので)実行しておくことでより高速化を図れます。 逆に言えばパフォーマンス測定のためには予め空回ししておく必要があります。 以下では省スペース化のため省略していますが、時間測定前に一度は同じ計算を行っています。

何のパッケージ?

主に配列演算を高速化するパッケージです。 例えば2つのベクトルの差の自乗を取りたいときは、一番単純には

julia> N = 1000000
1000000

julia> xs, ys = rand(N), rand(N);

julia> @time sum(abs2(xs-ys))
elapsed time: 0.010056548 seconds (16000312 bytes allocated)
166877.0598863598

となります。これはこれで何をやっているのかわかりやすいのですが、 xs-ysabs2 とで2回走査してそのたびに新しいベクトルを作り、 更に最後に sum でもう1回操作するという流れになっています。 NumericExtensions.jl を使うと

julia> using NumericExtensions

# mapreduce
# foldl 部分で走査が1回減るが xs-ys がまだ残っている
julia> @time foldl(Add(), 0.0, Abs2Fun(), xs - ys)
elapsed time: 0.003482581 seconds (8000264 bytes allocated)
166877.05988635868

# 2つの配列の差に関して走る
julia> @time foldl_fdiff(Add(), 0.0, Abs2Fun(), xs, ys)
elapsed time: 0.001559327 seconds (176 bytes allocated)
166877.05988635868

# 縮約操作を和に特殊化
julia> @time sumfdiff(Abs2Fun(), xs, ys)
elapsed time: 0.001133061 seconds (176 bytes allocated)
166877.05988635944

# まさにやりたかったことを行う、そのものずばりな関数
julia> @time sumsqdiff(xs, ys)
elapsed time: 0.001081251 seconds (176 bytes allocated)
166877.05988635944

と言った具合に、様々なレベルで特殊化を行った関数が用意されているので、 より高速に配列操作を行うことができます。

関数オブジェクトによる高速化

数直線上原点からスタートして、単位時間あたりに ある決められた確率 p で+1 方向に進み、 または確率 1-p で-1 方向に進む粒子の運動(ランダムウォーク)の数値シミュレーションを行います。 同時にたくさんの粒子についてシミュレーションを行うことで、 粒子の存在位置の確率分布を(経験的に)求めることができます。

まず、次の関数は確率 px+1 を、1-px-1 を返す関数です。

julia> update(x::Integer, p::Real) = ifelse(rand()<p, x+1, x-1)
update (generic function with 1 method)

これを使って M 個の粒子を N ステップだけランダムウォークさせる シミュレーションは次のようにかけます *2 *3

julia> N,M = 100_000, 1000;

julia> xs = zeros(Int, M);

julia> for i in 1:N; xs = map(x->update(x,0.5), xs); end

julia> xs
1000-element Array{Int64,1}:
  -14
 -778
 -218
  374
 -570
  224
 -356
 -296
 -598
 -250
   96
    ⋮
 -316
 -208
  114
  390
   10
   -2
  228
 -144
   26
  154
   32

# 平均
julia> mean(xs)
-0.274

# 分散(不偏分散)
# なお、キーワード引数 corrected に false を渡すことで
# 不偏分散ではなく標本分散を計算できる
julia> var(xs)
93933.3422662663

問題がかなり簡単だという事情もありますが、 単純でわかりやすい部品の組立てによって簡単に開発することができます*4。 さて、こんな感じにもっと大きなプログラムを書いて、 解きたい問題のサイズも大きくなっていくと、実行時間も伸びてくるので、 高速化をしたくなります。 その時には、どこがホットスポットなのか、つまりどこが一番時間がかかっているかを プロファイラなどを使って測定するのですが、 とりあえずこの更新部分に最も時間がかかっていたとします。 そこでこの部分の高速化を行ってみましょう。

最初期のバージョンでは

julia> @time for i in 1:N; xs = map(x->update(x, 0.5), xs); end
elapsed time: 5.398709802 seconds (936531816 bytes allocated, 10.51% gc time)

という実行時間になりました。 実行時間のおよそ1/10 がガベージコレクタにかかる時間なのですが、 map は非破壊な関数、つまり引数として与えた配列に触らずに新しい配列をつくる関数なので、 毎回メモリ確保をすることが問題なのではないかと考えるわけです。 そこで代わりに破壊的な関数、つまり引数を直接書き換える関数であるmap! を使ってみます。

julia> @time for i in 1:N; map!(x->update(x, 0.5), xs); end
elapsed time: 7.686995863 seconds (873712320 bytes allocated, 7.36% gc time)

使用メモリがちょっと減って、ガベコレ時間の割合も減っているのですが、 実行時間が5割増しぐらいになってしまっています。 使用メモリも大して減っていないし、ぶっちゃけバグなんじゃないかとも思いますが、 とにかく別の手段を考える必要があります。

そもそも最初の例が遅い*5のは、 map 関数に渡す関数が無名関数リテラルで表現されており、 map の呼び出しごとに新しい無名関数を作り出しているため、 JIT コンパイラが機能しないからです *6。 しかしそもそもなぜ無名関数を作っているかというと、 もともと2引数関数だったものを部分適用して1引数関数にしたいからです。 今回の例では、部分適用する引数である確率p は (大体のケースで)シミュレーション全体を通して固定なので、 予め部分適用した関数を定義しておけますが、 そういうことができない場合もあります。 こういう時には、あらかじめ部分適用した関数をラッピングした型(関数オブジェクト)を作ってしまい、 束縛する値を型のフィールドとして持っておくことで、 部分適用関数にJIT コンパイルをかけることができます。 関数オブジェクトやそれに関する高階関数を定義しているのがNumericExtensions.jlNumericFuns.jl です。 NumericFuns にはN 個の引数を取る関数オブジェクトを表す抽象型 Functor{N} があります*7N=1,2 についてはそれぞれ UnaryFunctor, BinaryFunctor という別名があります。 update を1引数関数オブジェクトにしてみましょう。 もともと引数を2つ持つので、残った1つ(p) に束縛される実引数をフィールドとして保存しておく必要があります。

julia> using NumericFuns, NumericExtensions

julia> type UpdateFun <: NumericFuns.UnaryFunctor
       p :: Float64
       end

関数適用はNumericFuns.evaluate で行うので、メソッドを定義しておきます。 また型推論のための関数である NumericFuns.result_type も必要です。 なおNumericExtensions は これらの関数を import および export しています。

julia> NumericExtensions.evaluate(uf::UpdateFun, x::Integer) = update(x, uf.p)
evaluate (generic function with 143 methods)

julia> NumericExtensions.result_type{T<:Integer}(::UpdateFun, ::Type{T}) = T
result_type (generic function with 100 methods)

julia> evaluate(UpdateFun(0.5), 0)
-1

ここまでくれば後はNumericExtensins.map に渡すだけです。

julia> @time for i in 1:N; xs = map(UpdateFun(0.5), xs); end
elapsed time: 1.695740634 seconds (822383688 bytes allocated, 29.14% gc time)

最初のケースと比べて3倍以上早くなりました。 毎回関数オブジェクトが作られますが、呼ばれる関数は既に名前のついた NumericExtensions.evaluate なので、JIT コンパイルが働いて性能が向上します。 重要なポイントとして、UpdateFunコンストラクタに渡す値が変わっても、 JIT コンパイル自体されていることに変わりはないという点です *8。 実際にやってみると、

julia> @time for i in 1:N; xs = map(UpdateFun(rand()), xs); end
elapsed time: 1.620748121 seconds (822383688 bytes allocated, 30.23% gc time)

あれ、rand() の時間を考えるとむしろ早くなっているような……

なお、UpdateFunimmutable にするとむしろ性能が悪化します。

# immutable UpdateFun2
julia> @time for i in 1:N; xs = map(UpdateFun2(0.5), xs); end
elapsed time: 12.022501026 seconds (4637697080 bytes allocated, 23.65% gc time)

要求メモリが格段に増えているのが原因でガベコレが発生しまくっているのですが、 これは素朴に考えると、immutable 型ではC/C++ でいう参照渡しから値渡しに切り替わって、 コピーのコストがかかるのではないかと思います*9

さて、NumericExtensions.jl 版と初期版とを比較すると、 必要メモリ量は殆ど変わらず、従ってガベコレの時間はほとんど変わらないので、相対的にガベコレが3倍ほど重くなっています。 そのため、破壊的な関数を使うことでメモリ使用量を減らしてみましょう。 map! ではなく map1! となっていることに注意してください。

julia> @time for i in 1:N; map1!(UpdateFun(0.5), xs); end
elapsed time: 1.592824116 seconds (7983688 bytes allocated)

使用メモリ量が大幅に減った結果、ガベコレが発生しなくなりました。 ただ、実行実時間はちょっとしか減っていないので、本質的には遅くなっているようです。 この辺りのトレードオフは、どのくらいのサイズの問題を解くかにもよるので、 実際に測定するまでは判断できないでしょう。

その他の拡張たち

NumericFuns.jl

標準ライブラリにない数学関数がいくつか追加されています。 何があるかはREADME を参照してください。

また、標準ライブラリにある数学関数も含めて、関数オブジェクト化されています。 型名は$(capitalize(funcname))Fun です。例えばsin => SinFun など。

NumericExtensions.jl

前述のmap など、主に配列操作の拡張がメインです。 標準ライブラリにあるものでも、関数オブジェクト NumericFuns.Functor{N} を受け付けるようになっていたり、 単純に実行速度が上がっていることがあります。 この場合、using NumericExtensions とするだけで高速化がなされることがあります *10。 非常に沢山あるのでドキュメント を見てください。

  • mapping
    • (いくつかの)配列から配列を作る操作 例えば `map(abs2, 1:3) ==> [1,4,9]
    • map(f, xs)
      • f として2引数、3引数関数を用いる map(f,xs,ys)map(f,xs,ys,zs) もあります
    • map!(f, dst, xs) ==> dst = map(f,xs)
    • map1!(f, xs) ==> x = map(f,xs)
    • mapdiff(f, xs, ys) ==> map(f, xs-ys)
    • add!(xs,ys) ==> xs += ys
      • 他にもsubtract!, multiply!divide! など
    • absdiff(xs,ys) ==> abs(xs-ys)
    • sqrdiff(xs,ys) ==> abs2(xs-ys)
    • fma(xs,ys,zs) ==> xs + ys .* zs
  • reduction
    • いわゆる畳み込み操作
      • 例えば sum(1:3) ==> 6
    • sum(xs) ==> reduce(+,xs)
    • sum(f, xs) ==> mapreduce(f,+,xs)
      • 複数引数パターンもあります
    • sumdiff(f, xs, ys) ==> mapreduce(f, +, xs-ys)
    • sumabs(xs) ==> sum(abs(xs))
    • sumsq(xs) ==> sum(abs2(xs))
  • scan
    • 畳み込みながら、途中経過を逐次残すことで配列を作る
      • 例えば cumsum(1:3) ==> [1,3,6]

Julia 0.4.0-dev がでたけれど

2日目の記事でも書きましたが、 Julia の次のバージョンである v0.4.0 では Base.callメソッドを書くことで、 関数呼び出し構文 a() を任意の型の値 a に適用できます。 普通の関数と同じ記法を使えるので、NumericFuns.Functor{N} よりも扱いやすいのですが、 残念なことに現状ではこうして作った関数オブジェクトは Base.map などに渡すことができません*11。 とは言え流石に標準ライブラリの高階関数は、将来的にはユーザ定義型を受け入れるようになるんじゃないかとは思います。 また、性能が上がっている拡張関数たちも、いつかは標準ライブラリに引き受けられるだろうので、 いつの日にかこれらの2パッケージは役目を終えるかもしれません。 しかし、どうやら今はまだその時ではないようなので、是非是非使ってみてください。

おまけ

twitter における、日本語話者な Julia ユーザのリストを作りました。 まだまだまったくもって不完全だと思うので、追加とか削除とかの希望を(twitter で直接)いただければと思います。

https://twitter.com/yomichi_137/lists/julialangjp

*1:ちなみにこの記事は当時割愛した「チューニングの実例」の補足となっています。

*2:ifelse(cond, T, F)cond が真ならT を、偽ならF を返す関数。 関数なのでTFともに必ず評価されるが、 その代わりif や三項条件演算子?: と比べて速い

*3:実際には毎ステップごとに平均や分散などの計算を行うことになるでしょうが、本記事の本質ではないのでおもいっきり簡便化しました。

*4:言うまでもなく、他の現代的な言語でもこのぐらい簡単だったりしますが

*5:正確に言えば、NumericExtensions.jl のような、より高速な比較対象があって初めて「遅い」といえるのですが

*6: 忘れがちですが、do 記法は関数の第一引数として無名関数リテラルを流し込むという仕様になっているため、便利だからと乱用すると確実に性能を落とします。何度も実行されるものはきちんと外に関数を書きましょう。

*7:C++ 的な意味でのFunctor

*8:例えば複雑な分岐がある関数において、めったに起きないコードパスを通るパラメータを渡してしまった場合などでは影響があったりするかもしれませんが、コンパイラの癖も含めて調べていないのでよくわかりません

*9:それにしても上がり過ぎな気もするけれど

*10:ただ、単純な高速化の場合、既に標準に取り入れられていることもあるので、その意味で過度な期待はしないこと

*11:引数の型がUnion(Function, DataType) なので。

Julia v0.4.0-dev

Julia Advent Calendar 2014 の2日目*1の記事です。 前回の記事X 分で学ぶJulia -- りんごがでている では、Julia の最新リリースバージョンであるv0.3 の文法や特徴について、 ざっくりと、しかしながら仕事ができるぐらいには十分に解説がなされました。

この記事では、Julia の次期リリースに向けた開発バージョンである v0.4.0-dev について、v0.3 と大きく異なる点を中心に解説していきます。

開発バージョンについて

Julia は言語そのものの開発が活発であり、まだまだ仕様や文法が固まっていません。 しかしながら、ユーザを確保するためには、ある程度の期間にわたって仕様が変化しない安定バージョンが必要になります。 そのため、Julia には安定バージョン(最新リリースバージョン)と開発バージョンの2つが同時に存在します。 現在の安定バージョンはv0.3.3 (以下v0.3)で、開発バージョンはv0.4.0-dev (以下v0.4)です。

開発バージョンにはバグフィクスから性能向上、時には文法の変更まで、 大小様々な変更が日々なされていきます。 常々言われるように、Julia は非常に開発が活発であるため、開発バージョンは本当に毎日(毎時)変わります。 そのため開発バージョンはNightly build とも呼ばれます。

適当な時期になると、開発バージョンはバグフィクス以外の更新をストップして、 RC(Release Candidate: リリース候補)版へと移行します。 旧安定版を使っていたユーザは、RC 版が出ている間に、RC 版を実際に使うことで 次世代バージョンへと移行することが望まれます。 RC 版を出してしばらく様子を見て、特に問題がないようだと RC 版が次の安定リリース版になります。 なお、Github に登録されているマイルストーンによれば、 v0.4.0-RC は2015 年の2月10日頃に出る予定のようです*2

この記事は何?

RC 版が出てから次世代バージョンに移行するには数週間の期間が空きますが、 しつこいようですがJulia は言語自体の開発が非常に活発なので、 バージョン毎の変更点が結構大きなものとなっています。 実際に今回の開発バージョンをみるに、v0.3 とv0.4 との間とでは よく使われる文法レベルでもギャップが開いています。

この記事は、現在の開発版v0.4.0-devにおける最新リリース版v0.3.3 からの目立った変更点を列挙・概説することで、 現在v0.3.x を使っているだろう大多数のユーザが、 将来的に大きな混乱なく移行できることを目指します*3

面倒くさい人は「後方互換性を破る大きな変更」という項目と、 後は新機能の最初にある「ユーザ定義ドキュメント機能」だけ読めばいいかもしれません。

警告

(追記:2014-12-16

「開発版は毎日ビルドしてナンボだろJK」とか思っていたんですが、 ふと冷静に考えてみたら全くそんなことはないことに気づきました。 なので、適当なバージョンで固定して開発版をガンガン使っていくといいと思います(ぇー

本当にNightly Build する人には↓)

開発版を日常的にビルドして使っていると、ある日いきなり 今まで動いていたプログラムが動かなくなるといったことが往々にして起こります。 さらに、動くには動くけれど、警告が大量に表示されて可読性が著しく低下する、 と言った事態は日常茶飯事です。 むしろ今現在そんな状況です。 (20141217 追記: deprecated function を使うことで表示される警告を、Julia のコマンドラインオプションで抑制できるようになりました。 julia --depwarn={yes|no} で切り替え可能です。) 開発版は、このような正直いって面倒くさい状況をすら楽しめる人や、 文法の変更にかこつけてパッケージにプルリクエストを投げるのが 好きな人にしかおすすめしません。 特に、会社や組織にJulia を広げよう、と考えている人は、 絶対に 開発版 Nightly Buildを薦めてはいけないでしょう。

後方互換性を破る大きな変更

v0.4 では、後方互換性を破る変更がv0.3 の時以上にあるのですが、 その中でもかなり影響の広そうな破壊的変更が(私見では)3つあります。

です。これらは他の破壊的な変更点と比べても使っている事が多いので、 そのままv0.4-dev に移行するとまず確実に警告が出ます *4。 詳細に行く前に、v0.3 を使いながらでもできる準備があるので先に触れておきます。

Compat.jl パッケージ

Compat.jl パッケージを using すると、 開発版で使われている名前で古い版の型や変数を扱うことができます。 さらに @compat マクロを用いることで、 古い版からでも開発版で導入される新しい文法を使えるようになります *5

# Compat パッケージのインストール
julia-0.3> Pkg.update()
julia-0.3> Pkg.add("Compat")

# v0.4 ではString がAbstractString にリネームされる
# v0.3 では当然使えない
julia-0.3> AbstractString
ERROR: AbstractString not defined

# v0.4 における辞書型リテラル
# v0.3 では当然使えない
julia-0.3> Dict("Kagerou" => 1, "Shiranui" => 2) 
ERROR: unsupported or misplaced expression =>

# 魔法のパッケージ
julia-0.3> using Compat

# AbstractString がString の別名として定義される
julia-0.3> AbstractString
String

# v0.3 でも @compat を使うとv0.4 の辞書リテラルを使える
julia-0.3> @compat Dict("Kagerou" => 1, "Shiranui" => 2)
Dict{ASCIIString,Int64} with 2 entries:
  "Kagerou"  => 1
  "Shiranui" => 2

一度書き直す手間はかかりますが、 一旦書きなおしてしまえばそのまま移行することができます。 年末進行忙しいかもですが暇を見つけて是非やっておきましょう。

辞書型リテラルの変更

例えばアスキー文字列を整数に結びつける辞書型 Dict{ASCIIString, Int} を作るリテラルは、v0.3 までは

julia-0.3> (ASCIIString => Int)["Kagerou" => 1, "Shiranui" => 2]
Dict{ASCIIString,Int64} with 2 entries:
  "Kagerou"  => 1
  "Shiranui" => 2

のように配列リテラルの特別版のようなものになっていました。 この表記は自然なように思えますが、 =>予約語になってしまうという問題を抱えています。

そこで、v0.4 では、 まず有向グラフの辺を表現する型である Pair{A,B} 型が定義され、 =>Pair{A,B} コンストラクタの別名として オブジェクトとして定義されました。 そして、Dict{K,V}コンストラクタPair{K,V} な値を 可変長引数として渡すことで、辞書を作る事になりました。 言葉にすると面倒くさい感じがしますが、 実際には (K=>V)[]Dict{K,V}() と変えるだけです。

julia-0.4> Dict{ASCIIString, Int}( "Kagerou" => 1, "Shiranui" => 2 )
Dict{ASCIIString,Int64} with 2 entries:
  "Kagerou"  => 1
  "Shiranui" => 2

型推論があるので、大抵のケースでは辞書型の型パラメータは省略可能です。

julia-0.4> Dict( "Kagerou" => 1, "Shiranui" => 2 )
Dict{ASCIIString,Int64} with 2 entries:
  "Kagerou"  => 1
  "Shiranui" => 2

リスト内包表記は

julia-0.4> Dict{Int,Int}( [ x => x*x for x in 1:10] )
Dict{Int64,Int64} with 10 entries:
  7  => 49
  4  => 16
  9  => 81
  10 => 100
  2  => 4
  3  => 9
  5  => 25
  8  => 64
  6  => 36
  1  => 1

のように行います*6

Any 型配列リテラルの削除

v0.3 では、 [] の代わりに {} を用いて配列リテラルを書くことで 型推論をさせずに Any 型の配列、つまりなんでも入れられる配列を作ることが 出来ました。 しかし、同じことは型注釈を用いて Any[] とできるので、 この構文は冗長なものとなっています。 v0.4 ではこの構文がなくなることとなりました。

文字列抽象型 String の名前の変更

文字列を表す具体型 ASCIIStringUTF8String などは 抽象型 String の subtype でした。 この抽象型の名前が AbstractString に変更されます。

この変更により、 AbstractString が抽象型であることがはっきりし、 また将来的に具体型の名前として String を使えるようになりました。

新機能

ひとつでも気になるものがあったら開発版を触ってみましょう!(悪魔のささやき)

ユーザ定義ドキュメント機能

Julia ではhelp 関数に名前を渡したり、 REPL で? キーを叩いてヘルプモードに移行して調べたい名前を入力することで オンラインヘルプを見ることができます。

# REPL で'?' キーを押すとヘルプモードに移行する
help?> println
Base.println(x)

   Print (using "print()") "x" followed by a newline.

標準ライブラリの関数*7には大抵ドキュメントが付いているので 非常に便利なのですが、 各種パッケージを含むユーザ定義型では使えないため、歯痒い思いをすることが多くありました。

v0.4 では、そんなオンラインドキュメントにユーザが自分で追記をする機能が 追加されました。

# 最後のセミコロンは評価した値の表示の抑制
julia-0.4> const planck = 6.62606957e-34;

julia> @doc "Planck constant in SI Units" planck;

help?> planck
  Planck constant in SI Units

julia-0.4> @doc """
       function `hello` says "Hello, @doc!"
       """ hello() = println("Hello, @doc!")
hello (generic function with 1 method)

help?> hello
  function hello says "Hello, @doc!"

planck のように定義済みのオブジェクトにドキュメントを書くこともできるし、 hello のようにドキュメントを書くのと同時に定義することもできます。 ドキュメント文字列にはデフォルトでMarkdown 文法が使えます*8。 今回再現していないので伝わらないのですが、 バッククォートで囲った部分は違う色で表示されます。 また、今回のようにダブルクォーテーションを""" ... """ と3重にしておくと、 文字列中に出てくるダブルクォーテーションマークをエスケープしなくて良くなります。

ドキュメント文字列の次の行に関数や型をおく場合、 -> が必要です。 変数(定数)や関数にとどまらず、型にもドキュメントを書くことができます。

julia> @doc "
       type `Point` represents a point in 2D space
       " ->
       type Point
       x :: Float64
       y :: Float64
       end

help?> Point
  type Point represents a point in 2D space

Julia では同じ名前の関数を、別の型をもつ引数の組み合わせで定義することで 多重定義でき、それぞれの定義を引数の型に対応した「メソッド」と呼びます。 これらのメソッドにそれぞれ別のドキュメントを書くことができます。 ヘルプモードにおいて引数もしくは引数の型も同時に与えることで、 対応するメソッドのドキュメントを表示できます。

julia-0.4> @doc """
       `hello(x::Any)` greets `x`
       """ hello(x::Any) = println("Hello, $(x)!")
hello (generic function with 2 methods)

help?> hello
  hello() says "Hello, @doc!"

  hello(x::Any) greets the argument x

help?> hello()
  hello() says "Hello, @doc!"

help?> hello(Any)
  hello(x::Any) greets x

help?> hello(1)
  hello(x::Any) greets x

ちなみに、ドキュメント文字列を書いている場所には、 文字列以外の任意のオブジェクトが置けます。 また、@doc マクロはドキュメントを定義するだけでなく、 定義したドキュメントを取り出すこともできます。

julia> @doc 42 answer() = 0
answer (generic function with 1 method)

julia> @doc answer
1-element Array{Int64,1}:
 42

help?> answer
1-element Array{Int64,1}:
 42

お行儀は良くないけれど、何か面白いことができるかもしれませんね。

なお、書法などはまだまだ変更する可能性があるようです。 詳しくはREPL で実際に@doc のヘルプを見てください。 ちなみに @doc のヘルプは実際に @doc を用いて書かれているので、 base/docs.jl を読むと更に学べるでしょう。

関数呼び出し演算子 () をユーザ定義型で使えるようになった

a(x...) とかくと自動的にcall(a,x...) となるようになりました。 つまり、ユーザ定義型 A について、Base.call(a::A, x) を 定義することで、ユーザ定義型についても 関数呼び出し構文が使えるようになりました。

今まで関数の引数を部分適用したり高階関数を書くときは、 ラムダ式を使って無名関数を用いていました。 概念的にはこれでも全く問題ないのですが、 無名関数にはJIT コンパイルが働かないために、 性能の面からは使いづらいものとなっていました。 この新機能を用いると、例えば関数の部分適用は、 元の関数と束縛する仮引数の番号と束縛する実引数との3つ組を保存して、 call 関数のメソッドとして部分適用して得られた関数の適用を行うものを 実装することで実現できます。

キーワード引数が使いやすくなった

関数にキーワード引数を渡す方法は従来では次のとおりでした。

julia> foo(; x=0, y=0, z=0) = 100x + 10y + z

# 普通の呼び方
julia> foo(; x=1, y=2, z=3)
123

# 辞書のアンパック代入
julia> params = Dict( :x => 1, :y => 2, :z => 3) ;

julia> foo(; params...)
123

キーワード引数の名前を動的に変えたい場合(変数で渡したい場合)は、 次のように変数名シンボルと値のタプルを配列に入れてアンパック代入する必要がありました。

julia> a, b, c = :x, :y, :z ;

# (Symbol, Any) 配列のアンパック代入 ...
julia> foo(; [(a,1), (b,2), (c,3)]...)
123

v0.4 からは次のように簡便にかけます。

# 新記法
julia-0.4> foo(; a => 1, (b, 2), [c, 3] )
12

最初の a の例は、実際にはPair{Symbol,Int} を作っていることになります。 他にも普通使うだろう記法としてタプルと配列を例示しましたが、 実際にはイテレーティブなオブジェクト、 つまり for x in xsxs の位置にかけるオブジェクトならなんでもかけます。 この時、最初の2つの要素が名前と値として使われます。

@inline マクロの導入

まず、関数にコンパイラ向けのメタデータを付けられるようになりました。 詳細はよくわかっていませんが、何もしない場合は特に何も付けないようなので、 従来の使い方でのオーバーヘッドはないようです。たぶん。

さて、@inline マクロを関数定義の頭につけることで、 その関数にインライン化用のメタデータをつけられます。 インライン関数 inner を中で呼ぶ関数 outerJIT コンパイルされるとき、 inner の関数本体が outer の該当箇所にベタ書きされます。 これにより、速度が向上する場合があります。

注意としては、外側の関数がコンパイルされて初めて効果が出るということです。 そのため、コンパイルされていない最初の outer の呼び出しや、 トップレベルからの inner の呼び出しではインライン化されずに 普通に inner が関数として呼び出されます。 この時、(多分メタデータの解釈にオーバーヘッドがかかるために) 普通の関数を呼び出すよりも性能が低下するようです。 使い方の例と、実行速度の測定結果は以下のgist を参照ください。 使い方と結果

型パラメータに任意のisbits オブジェクトを取れるようになった

isbits な型とは、変更不可能で、 他の非isbits なオブジェクトへの参照をメンバに持たない、 真に変更不可能なオブジェクトです。 例えば UInt8Float64 等の数値型や、 これらだけを含む immutable 型である Complex{Float64} などが該当します *9

v0.3 までは型パラメータとして通常の型の他に、Int 型や Bool 型の値、 つまり整数 1 や真偽値 true など、を指定出来ました。 例えば配列型 Array{T,N}T は要素の型を、N は配列の添字の数*10を表します。

v0.4 では型パラメータとして他のisbits 型の値、例えば浮動小数点数 1.38e-23 や 分数 1 // 2 などが指定できるようになります。 これにより型の表現力が強化されます。 例えば、科学技術計算における次元解析を行う SIUnits.jl というパッケージが 将来的に拡張されるはずですし、 SI 単位系以外への変換も行えるようになるかもしれません。

staged function の導入

文法にはまだまだ考慮の余地があるようですが、 staged function を導入する予定のようです。 staged function については私は専門外でよくわからないのですが、 ひじょーーーーーーに大雑把にいえば、 複数引数(パラメータ)を取る関数に関して、 広範囲、例えばプログラム全体、にわたって変わらないパラメータを固定して、 残りの引数に関する関数を静的に作る技術、のようです。 staged function をうまいこと導入することによって メタプログラミングDSL などがより効率的に書けるようになるということです。

正直いって良くわかっていないので、是非ともプログラミング言語の専門家の方(もしくはそれに準じる方)に、 Julia を用いて具体的解説していただきたいところです。

後方互換性を破る小さな変更

最初に述べた3つに比べれば影響は小さいと思います。

符号なし整数型のスペル変更

符号なし整数型は今まで UintUint64 のように書かれましたが、 v0.4 ではUIntUInt64 といった具合に、大文字の I を使うようになりました。

なお、AppleSwift 言語における符号なし整数型を見て このスペルの合理性 (U(nsigned) I(nt))に気づいたそうですw

NothingNone の名前変更

Nothing

Julia の関数では、明示的に return を書かずとも最後に評価した式の値が返り値となります。 基本的に関数というのは、引数を元にしてなにかを返すものですが、 時には副作用だけを期待して何も返さない関数、というのを書きたくなることがあります。 呼び出し元に、でてきた結果をそのまま捨ててもらっても良いのですが、 返り値には何の意味もないのだということを明示したほうがユーザには優しいでしょう *11 *12。 そんな時に使えるのが nothing です。その型は Nothing です。

None

抽象型を使うことでいくつかの具体型をひとまとめにすることができますが、 既存の具体型に外から抽象型を与えることはできません。 その時に使うのが Union です。

julia> Int64 <: Integer
true

julia> Int64 <: Union(Int64, ASCIIString)
true

julia> Union( Union(Int64, UInt64), ASCIIString)
Union(ASCIIString, Int64, UInt64)

Union はこのように型の「足し算」となっているのですが、 その零元が None です。

Nullable

さて、「無」に関連づいた型が既に2つあるわけですが、 v0.4 では新たに無効な値になりうる値の型として、 Nullable{T} が導入されました。

julia-0.4> invaid = Nullable{Int}()
Nullable{Int64}()

julia-0.4> isnull(invalid)
true

julia-0.4> get(invalid)
ERROR: NullException()
 in get at nullable.jl:26

julia-0.4> get(invalid, 137)
137

julia-0.4> valid = Nullable{Int}(42)
Nullable(42)

julia-0.4> isnull(valid)
false

julia-0.4> get(valid)
42

julia-0.4> get(valid, 137)
42

「無」多すぎ問題

こうして似たような名前が増えてきたので、 NothingVoid に、NoneUnion() に変わりました *13

余談ですがもう1つ重要な「無」があって、それは () です。 主に0-引数関数 の引数を表すのに使います。これ自体の型も () です。

CharInteger ではなくなった

今までは 文字型 CharInteger のsubtype だったので、 加減や比較などのやりたい演算の枠を超えてなんでも *14 できる状況でした。

v0.4 では Char が直接 Any のsubtype になり、 文字コード操作のための比較や加減など、限られた演算のみがサポートされます。

他の変更と比べると、普通のバグフィックスですね。

convert で数値型を変更するときにオーバーフローチェックをするようになった

数値型をより狭い数値型に変換するとき、元の値によってはあふれることがあります。 今までは無視していたのですが、v0.4 では例外を投げるようになりました。

その他 deprecated な関数達

base/deprecated.jl に大量にならんだ関数のうち、 # 0.3 deprecations# 0.4 deprecations の間にあるものはv0.4.0rc ですべて消されます。 これらの関数を呼ぶと警告が出るので、さっさと直しましょう:D

機能改善

線形代数や乱数を中心として、様々な機能改善*15 がなされています。 たくさんあるので全部は挙げきれないで、個人的に気になったもののみ。

ガウス乱数生成関数 randn の高速化

速度が5割増しから2倍近くぐらいまで向上しました。 乱数生成はモンテカルラーにとっては最適化の律速になりうるので、 非常に喜ばしい改善です。

結構単純な変更で劇的な速度向上をしているので、 後日別の記事にて改めて触れたいと思います。 (20141220 : 触れました。http://yomichi.hateblo.jp/entry/2014/12/20/020748

Unicode7 のサポート

文字種がふえるよ!やったn(ry

最後に

どれか1つでも「これは!」という項目があれば、 実際に最新開発版を手元において実行してみることをおすすめします。 最新安定バージョンとの共存は、 きっと次の記事にてchezou さんが解説してくださるはずです *16。 開発版は本当に毎日、それもドキュメントなどの更新も含めて何回も更新されるので、 何が変わったのかわかりづらいかと思いますが、 基本的にはNEWS.md だけを追っていれば問題無いと思います。

*1:まだ27:40 ですよ!!!

*2:v0.3.0-RC は思いっきりオーバーランしましたが……

*3:私自身の備忘録という面もあります。

*4:v0.4.0rc までは警告を出しつつよしなにしてくれるが、rc 版では容赦なく落としてくると思われる。

*5:新文法を旧文法に変換することで実現している

*6:今までのやり方でも警告なしで動くが、仕様かどうかは不明

*7:型名 = コンストラクタも含む

*8:今後増やすのかどうかは不明

*9:isbits 関数のヘルプより。

*10:テンソルの階数

*11:実装を変更することで返り値が変わったりするので

*12:また、ある関数を同じ型で呼んだ時には、常に同じ型の値が返るようにしたほうがパフォーマンスが向上する

*13:None はもともと Union() の別名だったが、それをやめた形になる。

*14: 指数関数に入れたり、文字型の複素数とか有理数を作ったりもできる。 ちなみに 'a' // 'b' とかやるとスタックオーバーフローで落ちるw

*15:主に速度向上

*16:念のためさらっとコマンドだけ述べておくと、 git checkout v0.3.3 した後に make cleanall && make install して、作られたディレクトリ直下の bin/julia にパスを通せばよい。

スクリプト言語の起動時間

以前から、Julia は起動が超絶遅いのがネック、とか言っていたんですが、Julia0.3 になって起動が速くなったようです。

今まではおよそ2 秒かかっていて、明らかに待たされていたんだけれど、今は0.2 秒となり、一桁改善されました。

他のスクリプト言語はさらに一桁速く、流石にこれらと比べちゃうとワンテンポ間が開いてしまうように感じるけれど、Julia だけを単体で使う分には対して気にならないレベルになりました。

何が変わったのかを見ると楽しいかもしれないけれど、ちょっと今はそんな時間がないのでスルー。

initial time of several script languages

Julia 本Vol.1 の公開

2013年夏コミ(C84) で出したJulia 本のWeb 公開です*1

https://dl.dropboxusercontent.com/u/61544927/JuliaBook-Vol1-140105.pdf

C84 とくらべて一部加筆修正があります。具体的には

  • インストール方法をもう少し細かく書いた
    • コンパイラの指定などビルドオプション
    • 過去バージョンのインストール方法や複数バージョンの共存方法
  • 1ウォーカーの問題において、update! 関数 が書けないことについて訂正
  • Julia のバージョンアップ(0.2->0.3)に伴う微修正
    • Julia の起動速度について
    • methods 関数について

といったところです。

 

*1:冬コミにだしたVol.2 では「1/3 頃に公開するよー」とか書きましたが、遅れてしまいました。ごめんなさい。