gnuplotで遊ぼう 〜 階乗, フィボナッチ数, Brainfuck処理系まで

※この記事は古い記事です。最新のgnuplotではもっと簡単にコードを書くことができます。動くことは動きますが、まわりくどいコードが載っていますので、あまり参考にしないで下さい。

以下は自分の過去記事です.

どうやらこういうの好きみたい.


今回は, gnuplotでプログラミングをします.
最近, gnuplotの利用率が半端ないほど増えていて, 結構自由自在に扱えるようになってきたのです.
大量のグラフをgnuplotで捌いている時の快感は良いものです.


でも, この記事ではgnuplotなのにグラフは一切描きません!!!
他のサイトで勉強してください.
Makefileは足し算や引き算すら苦労したのですが, gnuplotは元々計算してプロットするためのものなので, 計算は簡単にできるのです.
しかし, ループを作る方法が特徴的です.

Hello, world!

まずは定番のハロー・ワールドから.

# main.gnuplot
print "Hello, world!"
~$ gnuplot main.gnuplot
Hello, world!
~$

見たまんまです.
printコマンドを使ってください.


文字列の結合は, 「.」です.

# main.gnuplot
s = "Hello, "."world!"
print s
~$ gnuplot main.gnuplot
Hello, world!
~$

数字の比較演算子は==やら!=, >=など, 普通なのですが, 文字列の比較はeq, neを使うことに注意です.

# main.gnuplot
s = "Hello"
if (s eq "Hello") print s; else print "No!"
s = s.", world!"
if (s eq "Hello") print s; else print "No!"
~$ gnuplot main.gnuplot
Hello
No!
~$

ループ

gnuplotは, 驚くべきことに標準的なforループを備えておりません(でした この記事は古い記事です).
ではループをしたい時どうすれば良いのでしょうか?


これまた驚くことなのですが, gnuplotでループをするには, rereadというコマンドでスクリプトを再読み込みするのが慣習的な方法です(でした この記事は古い記事です).
スクリプトを再読み込み」ですよ.
なんと原始的な!!!

# main.gnuplot
if (exists('i')) i = i + 1; else i = 0;
if (i < 10) print i; reread;
~$ gnuplot main.gnuplot
0
1
2
3
4
5
6
7
8
9

iという変数がなければ, そのスクリプトの最初ということで, 初期化します. そうでなければ1足します.
簡単ですね.


gnuplotは一行ごとに実行します(歴史を感じさせます).
ですから, 通常の言語のように改行して書きたい時は, バックスラッシュで続けます.

# main.gnuplot

if (exists('i'))\
  i = i + 1;\
else\
  i = 0;

if (i < 10)\
  print i;\
  reread;

どちらが綺麗か分からないですね.

階乗

ループのやり方が分かれば, 階乗を計算することができます.

# fact.gnuplot

if (exists('i'))\
  ans = ans * i;\
else\
  i = 0;\
  n = 10;\
  ans = 1;

i = i + 1;

if (i <= n)\
  reread;

print ans;
~$ gnuplot fact.gnuplot
3628800

rereadは, それ以降のコードを全て破棄します.
つまり, コードの再帰呼び出しをしているのとは異なり, 単純にスクリプトを最初から読み直す, ということです.
ですから, 上のコードで, 3628800から小さい方へansが出力されることはありません.

# fact.gnuplot

if (exists('i'))\
  ans = ans * i;\
else\
  i = 0;\
  n = 10;\
  ans = 1;

i = i + 1;

if (i <= n)\
  print i - 1, "! = ", ans;\
  reread;

print i - 1, "! = ", ans;
~$ gnuplot fact.gnuplot
0! = 1
1! = 1
2! = 2
3! = 6
4! = 24
5! = 120
6! = 720
7! = 5040
8! = 40320
9! = 362880
10! = 3628800

このrereadの挙動によって, エラーの行番号がスクリプトの行番号と一致しなくなります.

# error.gnuplot

if (exists('i'))\
  i = i + 1;\
else\
  i = 0;

if (i < 20)\
  reread;

error???
~$ gnuplot error.gnuplot

error???
^
"error.gnuplot", line 191: invalid command

スクリプトの中では11行目なのに, エラーの表示は191行目とあります.
スクリプトの9行めで20回スクリプトをrereadした後に, おかしな行に来るので, 9 * 20 + 11 = 191という事ですね.


これくらい簡単なスクリプトなら行数を自分で計算できるのですが, rereadがスクリプトの中で何度も出てくるようになると, エラーの行数からのデバッグが不可能になります.

九九表

ループをするにはrereadでした.
では, 二重ループをするにはどうすればいいのでしょうか?
答えは, 内側のループをするスクリプトを別ファイルに作り, それをloadする, です.

# kuku.gnuplot

if (exists('i'))\
  i = i + 1;\
  j = 1;\
else\
  i = j = 1;\
  s = '';

if (i < 10)\
  s = '';\
  load "kukusub.gnuplot";\
  print s;\
  reread;
# kukusub.gnuplot
if (j < 10)\
  s = s . ' ' . (i * j);\
  j = j + 1;\
  reread;
~$ gnuplot kuku.gnuplot
1 2 3 4 5 6 7 8 9
2 4 6 8 10 12 14 16 18
3 6 9 12 15 18 21 24 27
4 8 12 16 20 24 28 32 36
5 10 15 20 25 30 35 40 45
6 12 18 24 30 36 42 48 54
7 14 21 28 35 42 49 56 63
8 16 24 32 40 48 56 64 72
9 18 27 36 45 54 63 72 81

要は, まずメインのスクリプトがあって, その中にループがあるたびに補助的なファイルを作れ, そういうことなのです.


しかし, どうしても別ファイルのするのには抵抗がある, そういう人もいるでしょう.
一つのスクリプトで二重ループを作りたい!
それも, そんなに難しくないです.

# kuku.gnuplot

if (exists('i'))\
  j = j + 1;\
else\
  i = j = 1;\
  s = '';

if (j < 10)\
  s = s . ' ' . (i * j);\
  reread;

if (i < 10)\
  print s;\
  i = i + 1;\
  j = 0;\
  s = '';\
  reread;
~$ gnuplot kuku.gnuplot
1 2 3 4 5 6 7 8 9
2 4 6 8 10 12 14 16 18
3 6 9 12 15 18 21 24 27
4 8 12 16 20 24 28 32 36
5 10 15 20 25 30 35 40 45
6 12 18 24 30 36 42 48 54
7 14 21 28 35 42 49 56 63
8 16 24 32 40 48 56 64 72
9 18 27 36 45 54 63 72 81

どの条件の時にrereadして変数がどう変化するかを, きちんと追えばいいんです.
早期returnならぬ, 早期rereadですか.
或いは, ループの深さを表す変数を作るか, ですね.

フィボナッチ数

ループが出来ればフィボナッチ数も簡単ですね.

# fib.gnuplot

if (exists('i'))\
  i = i + 1;\
else\
  i = 1;\
  ans = 0;\
  ansb = 1;\
  n = 10;

if (i <= n)\
  print 'fib(', i - 1, ') = ', ans;\
  c = ansb + ans;\
  ans = ansb;\
  ansb = c;\
  reread;

print 'fib(', i - 1, ') = ', ans;
~$ gnuplot fib.gnuplot
fib(0) = 0
fib(1) = 1
fib(2) = 1
fib(3) = 2
fib(4) = 3
fib(5) = 5
fib(6) = 8
fib(7) = 13
fib(8) = 21
fib(9) = 34
fib(10) = 55

Brainfuckインタープリタ

ここまでできるなら, Brainfuckインタープリタだって作れるよね.

「やりましょう」

やった結果がこちら.

# bf.gnuplot
# Brainfuck in gnuplot (supported by sed, wc)
if (exists('i'))\
  i = i + 1;\
else\
  i = 1;\
  code = "+++++++++[>++++++++>+++++++++++>+++++<<<-]>.>++.+++++++..+++.>-.------------.<++++++++.--------.+++.------.--------.>+.";\
  sindex(s,i) = i < 0 ? '' : i > 255 ? sindex(sslice(s,255),i-255) :\
    system('echo "'.s.'" | sed "s/.\{'.i.'\}\(.\).*/\1/"');\
  iindex(s,i) = i < 0 ? '' : i > 255 ? iindex(islice(s,255),i-255) :\
    system('echo "'.s.'" | sed "s/\([^,]*,\)\{'.i.'\}\([^,]*\).*/\2/"');\
  sslice(s,i) = i < 0 ? s : i > 255 ? sslice(sslice(s,255),i-255) : \
    system('echo "'.s.'" | sed "s/\(.\)\{0,'.i.'\}//"');\
  islice(s,i) = i < 0 ? s : i > 255 ? islice(islice(s,255),i-255) : \
    system('echo "'.s.'" | sed "s/\([^,]*,\)\{0,'.i.'\}//"');\
  take(s,i) = i < 0 ? '' : i > 255 ? (take(s,255).take(islice(s,255),i-255)) :\
    system('echo "'.s.'" | sed "s/\(\([^,]*,\)\{'.i.'\}\).*/\1/"');\
  slen(s) = int(system('echo "'.s.'" | sed "s/\(.\)/\1 /g" | wc -w'));\
  ilen(s) = int(system('echo "'.s.'" | sed "s/,*\([^,]*\),*/\1 /g" | wc -w'));\
  string(i) = sprintf('%d', i);\
  replace(s,i,r) = take(s,i).r.','.islice(s,i+1);\
  replicate(s,i) = i < 1 ? '' : i < 2 ? s : replicate(s,i/2).replicate(s,i/2).(i&1?s:'');\
  pointerlen = 10;\
  pointer = replicate('0,', pointerlen);\
  codelen = slen(code);\
  cindex = 0;\
  pindex = 0;\
  outputstr = '';\
  state = 0;\
  counter = 0

if (state == 1)\
  cindex = cindex + 1; codechr = sindex(code, cindex)
if (state == 1 && codechr eq '[')\
  counter = counter + 1; reread
if (state == 1 && codechr eq ']' && counter > 0)\
  counter = counter - 1; reread
if (state == 1 && codechr ne ']')\
  reread
if (state == 1)\
  state = 0; counter = 0;

if (state == 2)\
  cindex = cindex - 1; codechr = sindex(code, cindex)
if (state == 2 && codechr eq ']')\
  counter = counter + 1; reread
if (state == 2 && codechr eq '[' && counter > 0)\
  counter = counter - 1; reread
if (state == 2 && codechr ne '[')\
  reread
if (state == 2)\
  state = 0; counter = 0;

codechr = sindex(code, cindex)
pointerval = int(iindex(pointer, pindex))

if (codechr eq '>')\
  pindex = pindex + 1
if (codechr eq '<')\
  pindex = pindex - 1
if (codechr eq '+')\
  pointer = replace(pointer, pindex, string(pointerval + 1))
if (codechr eq '-')\
  pointer = replace(pointer, pindex, string(pointerval - 1))
if (codechr eq '[' && pointerval == 0)\
  state = 1; reread
if (codechr eq ']' && pointerval != 0)\
  state = 2; reread
if (codechr eq '.')\
  outputstr = outputstr.sprintf('%c', pointerval);

if (pindex >= pointerlen)\
  pointer = pointer . '0,';\
  pointerlen = pointerlen + 1

if (pindex < 0)\
  print 'out of memory! (negative pointer region access)';\
  exit 1

formatptr = replicate('  %3d ',pindex).' [%3d]'.replicate('  %3d ',pointerlen - pindex - 1)
eval("print sprintf('%3d  %s   ".formatptr."  %s',cindex,codechr,".pointer."outputstr)")

cindex = cindex + 1
if (cindex < codelen) reread

print outputstr

はい、完全にsedです、ありがとうございました。

※この記事は古い記事です。最新のgnuplotではもっと簡単にコードを書くことができます。あまり参考にしないで下さい。というかsedを呼ばなくても普通に文字列の操作できるようになっているそうです。誰か書いておいて下さい。

それはともかく, 出力は美しい.
gnuplotの処理速度が遅いせいか, 或いはsedを何度も読んでいるために遅いのか, ポインターが動く様子が目に見えるくらいで, 良い感じ(gnuplotが入っている方は実行してみて欲しい).

~$ gnuplot bf.gnuplot
  0  +    [  1]    0     0     0     0     0     0     0     0     0
  1  +    [  2]    0     0     0     0     0     0     0     0     0
  2  +    [  3]    0     0     0     0     0     0     0     0     0
  3  +    [  4]    0     0     0     0     0     0     0     0     0
  4  +    [  5]    0     0     0     0     0     0     0     0     0
  5  +    [  6]    0     0     0     0     0     0     0     0     0
  6  +    [  7]    0     0     0     0     0     0     0     0     0
  7  +    [  8]    0     0     0     0     0     0     0     0     0
  8  +    [  9]    0     0     0     0     0     0     0     0     0
  9  [    [  9]    0     0     0     0     0     0     0     0     0
 10  >       9  [  0]    0     0     0     0     0     0     0     0
 11  +       9  [  1]    0     0     0     0     0     0     0     0
 12  +       9  [  2]    0     0     0     0     0     0     0     0
 13  +       9  [  3]    0     0     0     0     0     0     0     0
(中略)
 97  +       0    72  [113]   32     0     0     0     0     0     0   Hello, wo
 98  +       0    72  [114]   32     0     0     0     0     0     0   Hello, wo
 99  .       0    72  [114]   32     0     0     0     0     0     0   Hello, wor
100  -       0    72  [113]   32     0     0     0     0     0     0   Hello, wor
101  -       0    72  [112]   32     0     0     0     0     0     0   Hello, wor
102  -       0    72  [111]   32     0     0     0     0     0     0   Hello, wor
103  -       0    72  [110]   32     0     0     0     0     0     0   Hello, wor
104  -       0    72  [109]   32     0     0     0     0     0     0   Hello, wor
105  -       0    72  [108]   32     0     0     0     0     0     0   Hello, wor
106  .       0    72  [108]   32     0     0     0     0     0     0   Hello, worl
107  -       0    72  [107]   32     0     0     0     0     0     0   Hello, worl
108  -       0    72  [106]   32     0     0     0     0     0     0   Hello, worl
109  -       0    72  [105]   32     0     0     0     0     0     0   Hello, worl
110  -       0    72  [104]   32     0     0     0     0     0     0   Hello, worl
111  -       0    72  [103]   32     0     0     0     0     0     0   Hello, worl
112  -       0    72  [102]   32     0     0     0     0     0     0   Hello, worl
113  -       0    72  [101]   32     0     0     0     0     0     0   Hello, worl
114  -       0    72  [100]   32     0     0     0     0     0     0   Hello, worl
115  .       0    72  [100]   32     0     0     0     0     0     0   Hello, world
116  >       0    72   100  [ 32]    0     0     0     0     0     0   Hello, world
117  +       0    72   100  [ 33]    0     0     0     0     0     0   Hello, world
118  .       0    72   100  [ 33]    0     0     0     0     0     0   Hello, world!
Hello, world!


sedをなぜ使ったかって?
それは, 文字列からインデックスアクセスする方法が分からなかったので, system関数でsedちゃんを呼びました.
ちなみにsedはTuring completeらしい*1ので, これだけではgnuplotがそうである証明にはなれませんね*2.
まぁ入力は別に平らな文字列じゃなくてもいいって説もあるけど...
取り敢えず面倒なので, gnuplotオンリーでのBrainfuckインタープリタは, 読者への課題とします.
evalがあれば配列要らないし, なんでもできる気がしますね.
tableを使ってその都度ファイルに書き出し&読み込みを行なっていけばできるのかな...???


ハロー・ワールドくらい簡単だからできたんじゃとか思われる方は, bfのコードを

  code = "+++++++++++>+>>>>++++++++++++++++++++++++++++++++++++++++++++>++++++++++++++++++++++++++++++++<<<<<<[>[>>>>>>+>+<<<<<<<-]>>>>>>>[<<<<<<<+>>>>>>>-]<[>++++++++++[-<-[>>+>+<<<-]>>>[<<<+>>>-]+<[>[-]<[-]]>[<<[>>>+<<<-]>>[-]]<<]>>>[>>+>+<<<-]>>>[<<<+>>>-]+<[>[-]<[-]]>[<<+>>[-]]<<<<<<<]>>>>>[++++++++++++++++++++++++++++++++++++++++++++++++.[-]]++++++++++<[->-<]>++++++++++++++++++++++++++++++++++++++++++++++++.[-]<<<<<<<<<<<<[>>>+>+<<<<-]>>>>[<<<<+>>>>-]<-[>>.>.<<<[-]]<<[>>+>+<<<-]>>>[<<<+>>>-]<<[<+>-]>[<+>-]<<<-]";\

にしてみてください.
http://esoteric.sange.fi/brainfuck/bf-source/prog/fibonacci.txt から取ってきたのですが, 確かにフィボナッチ数を計算できることが分かります (変数iを確かめると, 33万1117回スクリプトがrereadされます).
それにしても, こういうBrainfuckのコードすらすら書ける人, 頭イッてると思う (いい意味ですよ).


あと, gnuplotBrainfuckインタープリタを書いた自分もイッてると思う.
もしかしたら, 自分が初めてなんじゃないかなぁ...
Makefileの時然りだけど, 「なぜBrainfuckインタープリタを作るのですか?」と問われれば, 「そこに言語があるから」と答えますね.
Turing completeな言語ね.

終了しないgnuplotスクリプト

永遠に終了しないスクリプトを書くことはできるのでしょうか.

reread
~$ gnuplot main.gnuplot
^C
~$

..... 簡単でした.
終了しないスクリプトの中では, あらゆる言語の中でも短い部類なんじゃないかな.
C言語だと「main(){while(1);}」で17文字, Haskellでも「main=main」で9文字, gnuplotの「reread」の6文字には叶わないです.


無限ループがあれば, 出来るだけ早く止めるゲーム出来ますね.

if (exists('i')) i = i + 1; else i = 0
print i
reread
~$ gnuplot main.gnuplot
0
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
^C
~$

僕は15が最高でした.
Ctrlキーは最初から押しておいて, Enter->Cと連続して押すのがコツですね.
何の話や...

whileあるよ

はい(白目).
最近のgnuplotではwhile文があります.
これを使うと, 階乗だって,

# fact.gnuplot

i = 0
n = 10
ans = 1
while (i < n) {
  i = i + 1
  ans = ans * i
}
print ans

フィボナッチ数だって,

# fib.gnuplot

i = 1
ans = 0
ansb = 1
n = 10

while (i <= n) {
  print 'fib(', i - 1, ') = ', ans
  c = ansb + ans
  ans = ansb
  ansb = c
  i = i + 1
}

print 'fib(', i - 1, ') = ', ans

九九表だって

# kuku.gnuplot

i = j = 1
s = ''
while (i < 10) {
  j = 1
  s = ''
  while (j < 10) {
    s = s . ' ' . (i * j)
    j = j + 1
  }
  i = i + 1
  print s
}


さらに言えば, do for [=:] { ... } という構文もあります.
これを使えば, 階乗はこうなるし,

# fact.gnuplot
n = 10
ans = 1
do for [i=1:n] {
  ans = ans * i
}
print ans

九九表は

# kuku.gnuplot
s = ''
do for [i=1:9] {
  s = ''
  do for [j=1:9] {
    s = s . ' ' . (i * j)
  }
  print s
}

ああ... 簡単ですね... 極めて自然なコードに見えます... 深い悲しみ...
ただし, whileもdo forも実装されたのがホントつい最近なので(2012年8月のバージョン4.6から), まだまだ動かない環境があるので注意です.
個人的には, rereadでループのほうがトリッキーだし歴史が感じられてなんかいいなぁ...

まとめ

階乗, フィボナッチ数を計算できたどころか, Brainfuckの処理系も書けた(ただしチートあり)ので, gnuplotは立派なプログラミング言語です. 何でも出来ます. 多分複雑な計算の処理が速いはずだから, Brainfuckインタープリタとかあほな事やってないで, gnuplotgnuplotらしい使い方を心がけましょう. グラフを描いて幸せになりましょう.

*1:http://www.catonmat.net/blog/proof-that-sed-is-turing-complete/

*2:ただ, rereadとifがあるのでTuring complete. sedを使ったのはコードへのアクセスの問題.