14. 科学技術計算用言語 Julia

juliia web

Julia とは何か

ここで紹介する Julia とは、MIT とか NYU あたりの人たちが開発中の,2012年に公開が始まった,科学技術計算の専用言語のこと. 本家 web にマニュアルやらチュートリアルやら動画やら、大変多くの様々な資料があるので、困ったらまずはここをみると良い.

そしてその言語の能力であるが、 開発者らが欲張りだと自認するだけあって,後発の強みも活かして,

  1. 科学技術計算に便利な文法になっている

  2. 科学技術計算に便利なコマンドやライブラリが最初からある

  3. *計算が速い
    注: 実行時間は(JIT コンパイルにかかる時間は除いて)ほんとに高速で,C や Fortran とほぼ同等だ. Julia Micro-Benchmark に他のコンピュータ言語等での計算速度を比較したベンチマークのグラフが載っているので見てみよう.

  4. 言語設計が近代的なので,プログラミングが随分楽

という、大変多くの利点がある言語になっている.これからが楽しみな言語だ.

今のところ?JIT compile, つまり「実行時にコンパイルして動かすスタイル」をとっているので,使い勝手的にはインタプリタと言える プログラム開発中や,いろいろ試しながら動かすにも,これは便利.

よって,「たいへん気楽に使える高速な科学技術計算言語」として使える. また,ライセンスはMIT ライセンスなので安心して無料で使えるし,ソースも全部公開されているので中身をいじりたければ自由にいじれる.


ぜひ、この言語を自分のものとして、便利に使えるようになってもらいたい.

知っておくと混乱しないで済む?基本知識

少し、大雑把なことを書いておこう.

  1. Julia は(上に書いたように) JIT コンパイラがその場でプログラムソースをコンパイルして動かす仕組み.
    つまり,プログラムを「本当に」実行するマシン上には Julia 環境(Julia カーネル)がインストールされている必要がある.

  2. プログラムやルーチンをそのカーネルセッションで「初めて動かすとき」は JIT コンパイルにかかる時間が余計にかかる.
    「ちょっと 1回だけ触ってみた」という使い方だとこのせいで「Julia って遅くね?」と誤解するかもしれないので注意しよう.

  3. Julia プログラムを動かす「カーネル」と,人間用インターフェイスは分離できるし,違うコンピュータに置くこともできる.
    よって「使い勝手の良いインターフェイス」を自分で選ぶことが推奨されるし,インターフェイスは手元の PC で動かしつつ,実際にプログラムが動く環境は計算用の大きなサーバに置くということもできる.

  4. 緩いデフォルト的な存在1はあるけれども,グラフを描くライブラリや文法は自由に選べる.
  1. インターフェイスとしては,以下の様なものが有名だ.好きなものが選べるぞ.


web browser で操作するインターフェイス. セル入力して送信した内容をそのたびに動作する形で,試行錯誤するのにぴったり.

あと,Julia だけではなくて,Python や R でも同じインターフェイス jupyter が使えるぞ2


インターフェイスとして Julia との関係は上図の通り. わかりやすいし,資料も多くておすすめ. これを使うと,Mathematica や Matlab のような見た目と使い勝手になる.
今回はこれを用いよう



リリースから1年ぐらいの,かなり新しいインターフェイス. jupyter 同様,web browser で操作するものだが,セッション中の全入力が同期(一箇所の変更がただちにセッション中全ての入力と結果に反映される = 計算し直しが起きる,ということ)して動作する. わかりやすさはピカイチかな.

規模が大きい処理には向かないだろうが,比較的小さな実験を繰り返すとか,インタラクティブな操作をスムーズに行いたい時などは特に威力を発揮するだろう.


Julia for VSCode
VS code エディタのプラグインと Julia のパッケージをうまく使った IDE 環境. かつての,Atom エディタ上で動く Julia 用の統合開発環境だった Juno の後継環境のようだ3

大人気の VS code で使えるし,文法に沿ったいろいろなアシスト機能があって大変便利のようだ. ただし,トラブルが起きた時は Julia package, VS code, このプラグインの 3者のどこに原因があるかを調べねばならないため解決はやや面倒かもしれないね.


その他にも以下のようなインターフェイスもある.

  • 文字端末ソフトで直接対話する REPL

    Windows だと Powershell で,Unix 系だと shell から Julia を立ち上げて,いわゆる「黒い画面」で操作するもので, マニュアルなんかでは read-eval-print loop の略で REPL と呼ばれている.
    悪くない.悪くないが,カーネルが落ちるとそれまでのやりとりも画面から消えてしまうので,長時間扱うならやっぱインターフェイスは別にしておくといいかも.
    ちなみに,REPL でもグラフや画像の表示は思ったよりきちんとできる.画像がぽこっと別ウィンドウや web browser 呼び出しで出てくる様子はちょっとかわいいので一度試してみるといいかも.

  • Emacs での Julia 環境

    julia-mode 用の elisp (syntax highlight4 も行ってくれるぞ)や Emacs で動く統合開発環境 ESS (Emacs Speaks Statistics) などが古参の定番かな.

    あと,Julia の文法対応をしてくれる LanguageServer と組み合わせて使う eglot なんかは個人的にオススメだ.
  • Julia-vim

    vim の julia-mode らしい. vi使いは調べておくと良さそうだね.

それぞれの学習環境での状況

今回は Julia およびその web インターフェイス jupyter が使えないといけない. 詳細は下記に記載したが,なるべく自分の PC/Mac を用いるのが良いだろう.

環境 状況
阪大情報教育システム
Windows10
Julia 1.3.1 がインストール済み.ただしfirewall 下にあるため, 様々なファイルを取得するための git コマンドや curl コマンドがそのままだとインターネットと通信できないなど,Julia を使うには現状では不向きな環境だ.
どうしてもという場合は Julia 利用のために: 0. 準備編. ファイヤーウォール中の場合の設定 と, Julia 利用のために: 2. jupyter インストール編 の作業をしてからなら使えるはず…だったのだが,現在はこれもブロックされてしまうようだ.
今回はこれを用いるのは諦めざるを得ないかな.
阪大情報教育システム
CentOS7
これも上と同様だ.
今回はこれを用いるのは諦めざるを得ないかな.
各自PC CentOS7
obtained at OSBoxes
Julia 利用のために: 1. Julia インストール編Julia 利用のために: 2. jupyter インストール編 の作業を行えば良い.
今回はこれを用いても良い.
各自PC Windows Julia 利用のために: 1. Julia インストール編Julia 利用のために: 2. jupyter インストール編 の作業を行えば良い.
今回はこれを用いても良い.
各自PC Mac OS Julia 利用のために: 1. Julia インストール編Julia 利用のために: 2. jupyter インストール編 の作業を行えば良い.
今回はこれを用いても良い.
上のどの環境も使えない人の場合 クラウド上で julia が使える環境を提供している会社があり,無料プランもある.
たとえば, CoCalc などが挙げられる.こうした環境を使ってみるのも良いだろう.

Julia and jupyter をインストールして使えるようにするには

まあ,下記を順番に読んでいけば良い.多くの人は準備編は関係ないだろう.

Julia 利用のために: 0. 準備編. ファイヤーウォール中の場合の設定

Julia のパッケージと呼ばれるライブラリ環境のインストールと更新は内部的に git, curl といった unix のネットワーク系コマンドを用いて行われる. そのため,firewall の中(たとえば阪大情報教育システムはそうである)では,適切な設定をしないとパッケージをうまく扱えない.

そこで,たとえば阪大情報教育システム(windows)上の Julia を用いるのであれば,(Julia インストール直後ぐらいに)次のような設定をしておく必要がある.

1.git コマンドのための設定: Windows の メニュー → "Command Prompt" として command prompt を起動し,次の3つのコマンドを実行する.すると,とある設定ファイルが作られて下記の内容が書き込まれ, 以降のgit コマンドは proxy サーバを介して外へつながるようになる.

git config --global http.proxy http://toyonakans.ecs.cmc.osaka-u.ac.jp:3128
git config --global https.proxy http://toyonakans.ecs.cmc.osaka-u.ac.jp:3128
git config --global url."https://".insteadOf git://

2.curl コマンドのための設定:

2-a.おそらく c:\Users(もしくはユーザー)\アカウント名\.julia というディレクトリがあるので,その中に config というディレクトリを作る.

2-b.上で作ったディレクトリに startup.jl というファイルを作り,以下の2行を書き込んでおく.

1
2
ENV["HTTP_PROXY"] = "http://toyonakans.ecs.cmc.osaka-u.ac.jp:3128"
ENV["HTTPS_PROXY"] = "http://toyonakans.ecs.cmc.osaka-u.ac.jp:3128"

Julia 利用のために: 1. Julia インストール編

CentOS 7 に Julia をインストールする場合

下記のようにすれば最新バージョンの julia がインストールできるはずだ.

1
2
3
4
sudo yum update
sudo reboot  ← 念の為のリブート
sudo yum-config-manager --add-repo https://copr.fedorainfracloud.org/coprs/nalimilan/julia/repo/epel-7/nalimilan-julia-epel-7.repo
sudo yum install julia

注: 上の 3行目でなぜかエラーになることがある.エラーになった場合は,リブートしてやりなおしたりなど何回か試してみると良い.

通常の Windows に Julia をインストールする場合

1.Windows のアカウント名に
      半角アルファベット, 半角数字, 半角ハイフン, 半角ピリオド, 半角アンダースコア
以外の文字を使っていないか確認する(これ以外の文字を使っていると上記の jupyter 環境が使えない).

もしもこれらの文字をアカウント名の中で使っている人は,アカウント名をこの条件に沿う様に変えるか,この条件に沿ったアカウントを別に用意してそちらにログインし直し,下記の作業へ進もう. これらの文字をアカウント名の中に使っていない人は,そのまま次の作業へ進もう.

2.ver.1.6.1(for 64bit windows) をダウンロードしてダブルクリックしよう. 素の Julia がインストールできる.

通常の Mac に Julia をインストールする場合

Julia 1.6.1 for Mac 64bit をダウンロードし,help に従ってインストールしよう.

Julia 利用のために: 2. jupyter インストール編

jupyter インストールにはいくつかの方法があるが,ここでは最も簡単で,かつ,OS に依存しない方法を紹介しよう. 次の手順で一回だけ作業を行えば良い.

1.まず,Julia を起動する.Julia がインストールされているならば,どの OS でもメニューのどこかに Julia があるはずなのでそれを使おう. すると,

といった感じの window が起動するだろう(これは julia の REPL インターフェイス). ただし,window の色は環境依存なのでそこは気にしないように.

2.この Julia の window 中で ] キーを押す. すると,"julia>" というプロンプトが "(v1.6) pkg>" という感じのプロンプトに変わるだろう. これは Julia のパッケージ管理モードに入ったことを意味する.

3.このパッケージ管理モードのまま,

1
add IJulia

と入力する(大文字小文字の区別があるので注意). これは jupyter を使うための "IJulia" というパッケージのインストール作業である. しばらくかかるので待とう.

4.上の作業が終わったら,backspace キーを押す.すると通常モードに戻り "julia>" というプロンプトに戻るはずだ. そうしてから,次の 2つのコマンドを Julia に実行させよう.

1
2
using IJulia
notebook()

すると上の2つ目のコマンド notebook() を実行した時に jupyter 関連のソフトウェア一式をインストールしてくれる. このときインストール方式等について尋ねてくるが,よくわからない場合は y キーを押しておけば良い.

しばらく待つと,作業が終わる. もしかしたら web browser が起動したりもするが,この段階ではいったんその web browser も閉じ,この julia 環境も終了させてしまおう. これで一応,Julia と jupyter が使えるようになったはずだ.

Julia & jupyter を起動/終了するには

起動の手順

Julia を jupyter インターフェイスで利用するには,基本的には

1.Julia を起動.
2.Julia の中で以下の2行を行う.

1
2
using IJulia
notebook()

だけで良い. こうするだけで web browser が呼び出され,そこに下のような jupyter インターフェイスが表れるはずだ(見えている細かい内容は環境によって異なる).

これは jupyter のホーム画面とでも言うべき画面で,様々なファイル操作等はここを経由して行う.

次に Julia カーネルと接続して Julia を実際に使う画面を起動する. それには,この画面で右上に小さく見えている   New ▼   というボタンを押して,

というメニューを出現させ,ここで "Julia 1.6.1" という部分を選べば良い. すると

といった画面が起動する. ここで入手力を行うことで,Julia を使う,Julia プログラムを動かす,といったことができる.

終了の手順

実際の使い方はこのあとすぐ説明するとして,次にここから終了する手順を示しておこう.

まず,Julia を使っている画面の左上を見て,

上図の赤い丸がついている部分のアイコンを押し,念の為に現状をセーブしておく.

次に,その上の "File" という部分を押すとメニューが出るので,その一番下の "Close and Halt" を選ぶ.

次に,jupyter のホーム画面に移り, 右上の "Quit" ボタンを押す.これで jupyter インターフェイスが終了する. あとは Julia 環境を終了させれば良い.

ともかく試してみる

文法云々の堅い話はおいておいて、まあ、まずは触ってみよう.

  実習
インストール不要, 無料で Julia をすぐ使える JuliaBox が公開されている
へ飛ぼう.書いてあることの半分ぐらいが今回の貴殿らの環境に合致するだろうから,そうした部分について、jupyter の操作の練習を兼ねてごく簡単に Julia を触ってこよう.
(注: 上の web で言及されている JuliaBox というサービスは現在は廃止されている)

文法等のマニュアル、解説

日本語で書かれた解説が載っている M. Hiroi's Home Page: Julia Language Programming や、英語でよければ 本家のマニュアル あたりをみると詳細な文法の解説があるので、困ったらこれらを読もう.

また、できれば本を一冊買っておくと良い. おすすめは 書籍: Julia 1.0 Programming (2nd ed.) で、英語だけれども大変読み易い.

もう少し試してみる

実際に自分でプログラムを組む前に、もう少し試してみたほうが良いだろう. というわけで、下記の実習をやっておこう.

  実習

  1. 超入門
  2. 魔法の加速法

自分でやってみる

さて、次に、ごく簡単な問題に実際に Julia で取り組んでみるのが良いだろう. というわけで、次の

Collatz の予想(角谷の問題)

正の整数 $n$ に対して、次のような写像 $f$ を考える.

\begin{equation} f(n) = \left\{ \begin{array}{lcl} 3n + 1 & \mbox{     } & : n \mbox{ が奇数の場合} \newline n/2 & & : n \mbox{ が偶数の場合} \end{array} \right. \end{equation} このとき、任意の正の整数に対して、$f$ を複数回適用するといずれ必ず結果が 1 になる.

という、現在なお未解決の数学の予想がある. この問題を対象にして少しいろいろ練習してみよう.

まず最初は例を示そう.

上の写像 $f$ を関数としてプログラムする.

やり方はいろいろだ.例を示しておこう.

素朴な実装例

1
2
3
4
5
6
7
function f(n)
  if isodd(n)       # isodd は、中身が奇数かどうかを判定する
    return 3n+1     # Julia は定数係数と変数の間の積の記号を省略可能.
  else
    return div(n,2) # n/2 の結果は浮動小数なので、整数を返す div を使う.
  end
end

三項演算子を使ってシンプルに書いた例

1
2
3
function f(n)
  isodd(n) ? 3n+1 : div(n,2)
end

関数宣言もシンプルに書いた例

1
f(n) = isodd(n) ? 3n+1 : div(n,2)

  実習
自力で、もしくは上のように $f$ を実装して、$f(3)$ が 10 になること、$f(10)$ が 5 になることを、以下のように確認しておこう.

1
2
f(3)  # 結果は 10 のはず.
f(10) # 結果は 5 のはず.


$f$ を何回適用すると結果が 1 になるかをカウントする関数 count を作る

例を一つだけ示しておこう.しかし, なるべく例を見ないで、マニュアルを読むなどして自力で組んでみよう.

1
2
3
4
5
6
7
8
9
function count(n)
    m = n
    number = 0
    while m > 1
        m = f(m)
        number += 1
    end
    return number
end

  実習
自力で、もしくは上のように count を実装して、count(3) が 7 になること、count(17) が 12 になることを、以下のように確認しておこう.

1
2
count(3)  # 結果は 7 のはず.
count(17) # 結果は 12 のはず.


奇数 n = 3,5,7,... に対して、count(n) を求めてグラフを描く

  実習
これも例を示しておくが、なるべく自力で頑張ってみよう.

1.まず、count(n) の結果の配列を作る

1
y = [count(n) for n in 3:2:2001] # 1000個ほどデータを用意する.

2.適当なグラフ描画パッケージ(デフォルトは Plots だ)を使ってグラフを描く.まあ,ここでは Plots を使ってみよう.

1
2
3
4
using Plots # セッション中で一回宣言すれば良い.

plot(y, st = :scatter, leg = false )   
# st = :scatter は点描画,leg = false は legend表示なし,の指定.

すると、
collatz points
という感じの,なにか規則性があるような無いようなグラフが見て取れる. この規則性,数式にできる?


n を大きくしていくと,あたりまえだけど「大雑把には」count(n) の最大値が大きくなる傾向があるね

$i = 3,5,7,...,n$ での count(i) の最大値を v(n) とすると(式だと次のようになるね)

\begin{equation} v(n) = \max_{i \leq n} \{ \mbox{count}(i) \} \end{equation}

この挙動はどういうものだろうか.

すぐに数式で求められそうにはないので,プログラムで求めてグラフを描いてみるのが良さそうだ.

例えば次のような感じにできるだろう.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
N = 200001

n_sq = [ 3 ]

max_c = count(3)
v_sq = [ max_c ]

for i in 5:2:N
    global max_c, n_sq, v_sq
    
    v = count(i)
    if v > max_c
        max_c = v
        push!(n_sq, i)
        push!(v_sq, v)
    end
    
end

とすると,n_sq に最大値 $v(n)$ が更新されたときの $n$ のリストが,そして v_sq にはその $v(n)$ のリストが計算されて記録される. だから,グラフを描いてみるとこんな感じ.

1
2
plot(n_sq, v_sq, marker = :auto, leg = false)
# marker = :auto はマークを自動的に決める指定.

d

う~む,どうも $\log$ 関数に似ているような気もするので,適当に並べて描いてみて比べてみよう.

1
2
3
plot(n_sq, v_sq, marker = :auto, leg = false)
plot!(n_sq, 30 .* log.(n_sq) )
# plot! はグラフを重ね描く.

d

う~むむむむ… 似ているような似ていないような… もっと N を大きくすると様子がわかるかな?

N を大きくして,x 軸の対数をとってグラフにしてみよう.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
N = 10000001

n_sq = [ 3 ]

max_c = count(3)
v_sq = [ max_c ]

for i in 5:2:N
    global max_c, n_sq, v_sq
    
    v = count(i)
    if v > max_c
        max_c = v
        push!(n_sq, i)
        push!(v_sq, v)
    end
    
end

plot(n_sq, v_sq, marker = :auto, leg = false, xaxis = :log)

d

う~むむむむむむむむ… 直線に近いがちょっとだけ強い…?

両対数だとどうだろうか.

1
plot(n_sq, v_sq, marker = :auto, leg = false, xaxis = :log, yaxis = :log)

d

おお! なにか掴んだ気がするぞ! $n > 10^2$ ぐらいでの $v(n)$ をおおよその数式にしてみよう!(学生のみんなへの課題だね)

レポート

以下の課題について能う限り賢明な調査と考察を行い,
2021-AppliedMath7-Report-14
という題名をつけて e-mail にて教官宛にレポートとして提出せよ. なお,レポートを e-mail の代わりに TeX で作成した書面にて提出してもよい.

課題

  1. 上の内容を、一回自力でやり通してみよう. そして,最後のグラフの後に書いたように, $n > 10^2$ ぐらいでの $v(n)$ をおおよその数式の形で推測して,あらためてデータのグラフとその数式のグラフを同時に描いて比較してみよう.

  2. 「正 $n$ 角形の周長 $L(n)$ を定規で測ったデータ」をプログラム計算で擬似的に生成して5 、$n \rightarrow \infty$ で確かに $L(n) \rightarrow 2\pi r$ (ただし$r$ は正 $n$ 角形の外接円の半径相当)となることを数値データで確かめよう.

  3. 問2 のデータに対して、 魔法の加速法 で解説されている Aitken 加速法を適用したら収束先が「早い段階で」「より精密にわかるか」どうかを試して、その状況を報告しよう.さらに、多重に適用してみよう.

  1. 今の Julia では Plots package がグラフ描画のデフォルト的な存在だね. ↩︎

  2. そもそも jupyter という名前が JUlia + PYThon + R なのだ. ↩︎

  3. Atom エディタの開発が 2019年に実質的に停止した様子.そのため,Atom に依存していた Juno も先がないということで移行したのではないかと推測される. ↩︎

  4. 開発初期の syntax hilighting についてのブログ なんかは参考になるぞ. ↩︎

  5. 三角関数などを使って良いということだ.はっきり書くと、$L(n) = 2rn\sin(\pi/n)$ として良い、ということになるね. ↩︎