2015年12月28日月曜日

第20回シェル芸勉強会でびっくりしたことまとめ

第20回シェル芸勉強会に参加しました。幾つか初めて知って(過去に見たことがあるかもしれないけど)びっくりしたことがあったので備忘録としてまとめてみた。基本的に午前中で燃焼したので午後のことは書いていない。manを読むべきということを再確認した。そもそもGNUの指針だとmanにいっぱい書いてあるのは良くないはずじゃなかったっけ、とか言う話は知らない。(本当はinfoに書くべき(?))

AWKのPattern Ranges

とりあえずman抜粋

A pattern range consists of two expressions separated by a comma; in this case, the action shall be performed for all records between a match of the first expression and the following match of the second expression, inclusive. At this point, the pattern range can be repeated starting at input records subsequent to the end of the matched range.
見たことあるかもしれないけど、知らなかった。というかパターンのとこにカンマをカンマ演算子のつもりで書いてうまくいかなかったことがあったかも。 つまり

awk 'NR==2,NR==4'

sed -ne '2,4 p'
と同じ意味になるらしい。ついでに言うと

awk 'NR==2,0'

sed -ne '2,$ p'
と同じ意味で使うらしい。

nawkとそれ以外の$0代入の挙動が微妙に違う

これはやばいやつ。MacOSXだとawkがnawkらしいので結構重要かも。

具体的には(i > 1)で


awk '$0=$i'
の挙動が怪しい。gawkとmawkでは全体の真偽値が$iの真偽値と一致すると思って良いがnawkではそうでもないらしい。manを読んでもよくわからなかったので適当にコードを実行して調べた。本を読むと良いのかも。少なくとも調べた限りでは$iの中身がstringの場合には思った通りの挙動をするが、numericの場合には常に偽を返すっぽく見える。

もちろん


nawk '$0=$i'

nawk '($0=$i)||($0)'
の結果は変わらないのだが、(numericだと偽)

nawk '($0=$2)||(($1=$1)&&$0)'
ということにして一回フィールドを書き換えると$0の真偽値が想定しているものに変わる。因みに$0の値はprintして見る限りでは特に変なことになっていない。この挙動が意図したものなのかバグなのかはよくわからないが、少なくともMacOSXでドヤ顔して$0への代入をすると冷や汗をかくことになりそう。

seqと{..}

seq 100 1はだめだけどecho {100..1}が動くの知らなかった。{1..10..2}ができることも知らなかった。

因みに{foo,bar}で色々やる話は知っているけど技術的に使えないやつなのでとりあえず置いておく。${foo:ここ色々}系列も見たことあるし使ったこともある(manを見ないと半分くらい使えない)から今回は置いておく。

午後も自分じゃない端末ってどうやって判別するのかとかよくわかっていない点はあった(今もちゃんとやるにはどうすればいいかよくわかっていない)が、サーバエンジニアじゃないしな、って言って忘却の彼方に投げやった。

2015年12月24日木曜日

どうしてGNU Grepは高速なのか

この記事は アルゴリズム Advent Calendar 2015 24日目の記事です。 あまりアルゴリズムらしい話題でもないですがGNU grepがどうして高速なのかという話について大して詳しくもないですが最近知ったので書きます。

決して最近の投稿ではないですがwhy GNU grep is fastという投稿がFreeBSDのメーリングリストにありました。これはGNU grepの最初の作者のMike Haertelによる投稿で、当時幾つか解説記事も書かれました。FreeBSDのgrepよりGNUのgrepの方がかなり速い(今の速度は知らないですが当時はだいぶ違ったみたいです)みたいで、どうしてGNU grepが速いかということを書いた投稿です。

とりあえず前提知識について確認しておきます。grepは文字列中で正規表現を探索するプログラムで、伝統的にはThompsonのアルゴリズムでオートマトンを求めてマッチします。fgrepは文字列中で文字列を探索するプログラムで、伝統的にはAho-Corasick法で求めます。これらのアルゴリズムについてはここでは特に説明しません。

さて、本題に戻ります。この投稿によるとGNU grepではこんなことが速さの秘訣みたいです。

  • すべてのバイトを読むことなく処理している
  • 各バイトについて実行される命令数を少なくしている

そりゃ速くなりそうということが書かれていますがどうやってこれを実行しているのかが問題です。アルゴリズム的でない点としては例えばstdio.hにある標準関数ではなくmmapなどのシステムコールを直接叩いて入力についてはバッファリングしていないようです。また、例えばa$のような正規表現はa\nに書き換えることで読み込みを行単位でなくて済むようにしています。他にもループアンローリングとかまあやりそうな最適化はやっているらしいです。

さて、これはアルゴリズムアドベントカレンダーの記事なのでアルゴリズム的な話に触れます。主にすべてのバイトを読むことなく処理しているという点に効きます。原文ではこんな感じの内容です。

GNU grep uses heuristics to look for a fixed string that any string matching the regex *must* contain, and uses that fixed string as the bases for its initial Boyer-Moore search. For example if your regex is /foo.*bar/, the initial Boyer-Moore search is (probably) searching for foo.

GNU grepではBoyer-Mooreアルゴリズムを使っています。Boyer-Mooreアルゴリズムとは、検索対象の文字列(長さn)の中からキーの文字列(長さm)を探すアルゴリズムです。ナイーブに全探索で求めるとO(mn)になります。Boyer-Mooreアルゴリズムではキーの文字列に対して事前にあるアルゴリズムでテーブルを作り、そのテーブルにしたがって適宜文字を読み飛ばすことで高速に探索を行うことができます。文字の種類が多いというそこそこ現実的な条件で大体O(n)とかO(n/m)位の性能が出るようです。Boyer-Mooreアルゴリズムの詳細についてもここでは扱わないので調べてください。

ところでBoyer-Mooreアルゴリズムはgrepとは違って、文字列の中から正規表現ではなく文字列にマッチさせるアルゴリズムです。なのでナイーブに考えるとgrepというよりfgrepに使えそうなアルゴリズムです。GNU grepでは正規表現の中で含まれなければならない断片の固定文字列にマッチさせて、そこから正規表現マッチを行うというヒューリスティックを使っているらしいです。どうも単純にこのヒューリスティックを使っているだけでもなく、どうやら本当は色々な実装を組み合わせているらしいです。また、もっと複雑な正規表現の場合はそもそも難しいし滅多に使われないから多少遅くてもいいよねというヒューリスティック(?)の元で実装をしているらしいです(超訳)。

本当にそんな感じになっているのか気になるのでベンチマークを取ってみます。

$ time grep test /usr/share/dict/words > /dev/null

real    0m0.009s
user    0m0.007s
sys     0m0.000s

まず固定文字です。まあ速いです。

$ time grep tes.*t /usr/share/dict/words > /dev/null

real    0m0.010s
user    0m0.007s
sys     0m0.000s

多分BMでtesを探してからtを探します。固定文字とあまり時間が変わりません。固定文字部分がもう少し長いともっと速くなります。

$ time grep [a-z] /usr/share/dict/words  > /dev/null

real 0m0.074s
user 0m0.073s
sys 0m0.000s

$ time grep [A-Z] /usr/share/dict/words  > /dev/null

real 0m0.142s
user 0m0.137s
sys 0m0.003s

$ time grep [A-Z].*[A-Z] /usr/share/dict/words  > /dev/null

real    0m0.168s
user    0m0.163s
sys     0m0.003s

$ time grep [a-z].*[A-Z] /usr/share/dict/words  > /dev/null

real 0m0.337s
user 0m0.337s
sys 0m0.000s

3つ目の例では[A-Z]がなかなか見つからないので時間がかかっています。ひとつ目の[A-Z]もあまり見つからないので[A-Z]のみの場合とそこまで時間が変わらないです。4つ目の例になると[a-z]はほぼ常に見つかるのに[A-Z]がなかなか見つからないので時間がかかっています。[a-z]にかかる時間や[A-Z]にかかる時間と比べてかなり時間がかかっています。

2015年12月17日木曜日

シェル芸初心者向け勉強会に後日自宅にてサテライト参加(参加とは言わない)【後半】

さて、ここからシェル芸勉強会の過去問です。シェル芸勉強会の方は「コワイ」のイメージが先行して敷居が高いのでawk縛りとかawk使わない縛りとかいうことは言わずに解きます。

シェル芸勉強会過去問第2回 問題1: 文字化けしたファイルの削除

問題

次のように、空ファイル abc, DEFG と、文字化けした空ファイルを作ってください。 【今回は問題データのQ1ディレクトリ以下に作成済み】

$ touch abc
$ touch DEFG
$ echo ほげ | nkf -s | xargs touch
$ ls

%82ق%B0        DEFG        abc

解法

rm $(ls | grep -a -v [A-z])

解説

今回は文字化けどころではなく[A-z]を含まないファイルを全部消しました。

シェル芸勉強会過去問第2回 問題2: 計算

問題

以下のファイル中の数字を全部足し算してください。

$ cat Q2/num

1
2 3 4
5 6
7 8 9 10

解法

cat /tmp/num | tr "\n " + | sed 's/+$/\n/' | bc

解説

前述の通りこういう問題は式をいじって最後に評価するタイプの解法が好きなのでこう解きました。

シェル芸勉強会過去問第2回 問題3: 条件でデータを取り出し

問題

下のhogeファイルから、aとbについてそれぞれ一番大きな数を求めましょう。

$ cat Q3/hoge

a 12
a 13
b 13
a 432
b 111
b 43

解法

awk '{if(s[$1] < $2){s[$1]=$2}}END{for(i in s){print i,s[i]}}'

解説

なんかここまで普通に書くと負けた気がするので怪しい解ですが

cat /tmp/hoge | sed 's/a/000/;s/b/999/' | tr -d " " | sort -rn | sed 's/000/a /;s/999/b /' | uniq -w 1
というものも作ってみました。

シェル芸勉強会過去問第2回 問題4: 計算

問題

以下のファイル中の数字を、キー(a,b)ごとに全部足し算してください。

$ cat Q4/num2

a 1
b 2 3 4
a 5 6
b 7 8 9 10

解法

cat /tmp/num2  | awk '{for(i=2;i<=NF;i++){s[$1]+=$i}}END{for (i in s){print i,s[i]}}'
cat /tmp/num2  | sort | sed 's/\([0-9]\) /\1+/g;s/ /+=/' | cat - <(echo a;echo b) | bc

解説

bcを使ったほうがたのしいよ!

シェル芸勉強会過去問第2回 問題5: 日付と曜日

問題

以下のようにファイルを作って、何曜日が何日あるか集計してください。 【今回は問題データのQ5ディレクトリ以下に作成済み】

$ seq 1990 2012 | awk '{print $1 "0101"}' > osyoga2
$ head -n 5 osyoga2

19900101
19910101
19920101
19930101
19940101

解法

cat /tmp/osyoga2 | xargs -I{} date --date={} +%u | sort | uniq -c

解説

date使って曜日を計算させました。さっきから集計のときにsort | uniq -cばっかり使ってますね。

シェル芸勉強会過去問第2回 問題6: ダミーデータの作成

問題

seq 1 100 の出力をランダムに並び替えてください。

解法

seq 1 100 | shuf

解説

shufという便利な(?)コマンドがcoreutilsにあります。awkで解くのは嫌です。

シェル芸勉強会過去問第2回 問題7: 検索

問題

大文字小文字を区別しない場合、以下の辞書ファイルから、重複する単語を検索してください。 asciiコード以外の字がありますが、 とりあえず気にしないでください。

$ head Q7/words

A
A's
AA's
AB's
ABM's
AC's
ACTH's
AI's
AIDS's
AM's

解法

cat words | tr '[:upper:]' '[:lower:]' | sort | uniq -d

解説

とりあえず全部小文字にして重複文字を出しましたがこれでいいんでしょうか?

シェル芸勉強会過去問第2回 問題8: ファイルの比較

問題

file2から、file1にない数字を抽出してください。

$ cat Q8/file1

1
3
4
6
9

$ cat Q8/file2

2
3
4
5
9

解法

grep -v -f file1 file2

解説

grepの逆マッチ、ファイル入力版です。

シェル芸勉強会過去問第2回 問題9: 形式変換

問題

gameから、resultように整形してください。

$ cat Q9/game

1 表 2
1 裏 2
2 表 0
2 裏 1
3 表 3
3 裏 0
4 表 3
4 裏 0
5 表 0
5 裏 0

$ cat Q9/results

   1 2 3 4 5
表 2 0 3 3 0
裏 2 1 0 0 0

解法

cat game | awk '$2=="表"{omote[$1]=$3}$2=="裏"{ura[$1]=$3}END{printf "  ";for (i=1;i<=NR/2;i++){printf " %d",i} printf "\n表";for (i=1;i<=NR/2;i++){printf " %d",omote[i]} printf "\n裏";for (i=1;i<=NR/2;i++){printf " %d",ura[i]} printf "\n"}'
cat game | awk '{score[$2=="裏"][$1]=$3}END{printf "  ";for (i=1;i<=NR/2;i++){printf " %d",i} for(j=0;j<2;j++){printf "\n%s",j?"裏":"表";for (i=1;i<=NR/2;i++){printf " %d",score[j][i]}} printf "\n"}'

解説

awkの配列に入れて最後にforで取り出しました。あまりの泥臭さに少しなんとかしましたがもっと良い解法があるはず。

シェル芸勉強会過去問第2回 問題10: ファイルの結合

問題

file1、file2 から、下の出力を得てください。

$ cat Q10/file1

001 うそ
002 笑止
003 矢追純

$ cat Q10/file2

001 800
002 10000000
003 1

欲しい出力

001 うそ 800
002 笑止 10000000
003 矢追純 1

解法

cut -d ' ' -f 2 file2 | paste -d ' ' file1 -

解説

横結合はpasteでできますが、流石に題意を取り違えているような気がする。

2015年12月12日土曜日

シェル芸初心者向け勉強会に後日自宅にてサテライト参加(参加とは言わない)【前半】

シェル芸初心者向け勉強会が福岡であったようなので問題を解いてみました。ウォーミングアップ問題はawkを封じて解いた解法とawkで解いた解法を書きました。因みにtukubaiは全く使えない(僕の技術的に)ので使っていません。

ウォーミングアップ問題1

問題

1から10まで合計を計算してください

解法

seq -s + 10 | bc
seq 10 | awk '{t+=$1}END{print t}'

解説

awk封じでもかなり色々な解法があると思いますが個人的にはseqで式を生成してbcに食わせるのが好きです。

awkの方はかなりシンプルな感じで書きました。

ウォーミングアップ問題2

問題

下記のファイルを2列目の数字の小さい順に並べ替えてください。

$ cat list

b 11
d 5
a 3
e 4
c 2

解法

cat list | sort -k2,2 -n
cat list | awk '{s[$2]=$0}END{for (i in s){print s[i]}}'

解説

sortのオプションで二列目についてソートするようにしていすると簡単です。

awkでソートといえば連想配列ですが、今回は最初にふたつ目のキーの範囲がわかっていないという体裁で解きました。gawkだとこれで整列されるみたいですがmawkだと別の順になったのであんまり正当派解法ではないです。

ウォーミングアップ問題3

問題

文字列の並びを逆にしてください。

$ echo 'a b c'

解法

echo 'a b c' | rev
echo 'a b c' | awk '{for (i = NF;i>1;i--){printf "%s%s",$i,FS}print $1}'

解説

revというまさにこの用途のコマンドがあるので使いました。ちなみに行レベルでひっくり返したいときはtacを使います。

awkの方は普通に書きました。

ウォーミングアップ問題4

問題

下記の日付一覧から、月毎の日の数を数えてください

$ cat date1

20150101
20150121
20150201
20150202
20150203
20150310

解法

cat date1 | cut -c1-6 | uniq -c
cat date1 | awk '{s[substr($0,1,6)]++}END{for (i in s){print i,s[i];}}'

解説

cutで月のデータだけを切り出してuniq -cで行数を数えました。ソート済みでない場合はuniqの前にsort -nとかをしないといけないです。

awkの方も大体同じ感じで解きました。uniq相当のものは配列で作ります。

ウォーミングアップ問題5

問題

日付と商品が売れた数があります。月毎の商品の売れた数を数えてください。

$ cat date2

20150101 1
20150121 2
20150201 5
20150202 2
20150203 3
20150310 2

解法

cat date2 | cut -c 1-6,9-10 | xargs -I{} bash -c "yes \$(echo {}| cut -d ' ' -f 1) | head -n \$(echo {} | cut -d ' ' -f 2)" | uniq -c
cat date2 | awk '{s[substr($1,1,6)]+=$2}END{for (i in s){ print i,s[i];}}'

解説

awkを使わない方はかなり気合を入れて書きました。回数を数えるのにuniqを使おうとすると月ごとの行を商品の個数回繰り返す必要があります。それをyesとheadで各月について繰り返して、xargsで回しました。awkもtukubaiも使わなくてももう少し良い解法があるような気がします。

awkの方は普通に足し算をしました。

ウォーミングアップ問題6

問題

下記のファイルの積集合と和集合を出力してください。(註: ここで言う積集合とは直積(product)集合ではなく結び(intersection)の方です。)

$ cat data1

a
b
c
d

$ cat data2

b
d
e

解法

cat data1 | grep -f data2
cat data1 data2 | sort | uniq
cat data1 |awk '++s[$1];END{while((getline < "data2")>0){if(!s[$1])print;}}'
cat data1 | awk '{s[$1]++}END{while((getline < "data2")>0){if(s[$1])print $1;}}'

解説

awkを使わない方は大したことないです。grep -fを知らないと辛いかもしれません。

awkの方はちょっと辛いです。awkで複数ファイルを扱うときにはgetlineを使います。あとは大体前に書いたuniq相当の連想配列の使い方と同じ感じです。

ここから本題とは関係ない補足です。このgetlineの使い方はそんなに問題ないですが、BEGIN/END以外の場所で特にファイル指定をしないでgetlineを使ったりすると何を言っているのかわかりづらい挙動をするので、回避できるならやめたほうが良いとされています。

途中で力尽きたので第二回シェル芸勉強会の問題は後半ということにします。

2015年12月10日木曜日

そこまで遅くないShellスクリプトの書き方

この記事は Shell Script Advent Calendar 2015 10日目の記事です。
9日目の記事はryoana14さんの麗しきawkの世界でした。

Shellスクリプトがいつまで経ってもまともに書けないMasWagです。書けないなりにも人の書いた(昔の自分が書いたものも多く含む)スクリプトを見てこれは遅いなと思うことはたまにあります。書き方のコツというか考え方が幾つかあると思うのでまとめてみようと思います。基本的な話なので多分Shellスクリプトをあんまり普段書かない人向けだと思います。ShellはBashを前提として書きます。zshだともっと色々できるのかもしれないです。細かい説明は(そんなに細かくなくても?)省いているので適宜manやinfoを参照すると良いでしょう。

forkを減らす

Shellでコマンドを使うということは多くの場合プロセスをforkしていることになります。Cygwinでプロセスを立てるコストが非常に高いということが問題になったりしますが、*nixでもプロセスを立てるコストは多くの場合無視できません。そうとなるとshellの組み込みの機能でできることはやりたくなります。 Bashの組込み機能で使用頻度が高そうなものを挙げていきます。

まず数式処理"(())"があります。これは大体exprの代替だと思って間違いないです。案外使える機能に((t++))でインクリメントできたり((t+=2))のようなものもあります。

下の結果を見ると解ると思いますが、やはりforkの回数を減らすと目に見えて速度が速くなります。一回だけ実行するタイプではbcもかなり健闘しています。

次にtest/[]の代替の[[]]があります。過去の記事にもあるように単にtestの代替として使うだけでもちょっと速くなりますが、実は[[]]は正規表現マッチができます。grepを使って正規表現マッチしている分岐があったら書き直した方がいいです。またこれは正規表現の話ですが、複数の正規表現のorで分岐するならひとつの正規表現にまとめた方がいいです。

下の結果のようにまあ[[]]を使うと速いです。testとの比較だとそこまで時間が変わっていないですがgrepとの比較だとかなり時間に差が出てしまいます。

 最後に変数展開時の文字列操作があります。変数に対してsedとかを使いたくなったらこっちを使えないか考えるといいと思います。幾つか種類があるので詳細はmanやinfoを参照してください

pipeで繋げられるものは積極的に繋ぐ

Shellスクリプトでコマンドをpipeでつなげると、マルチプロセスで前の処理がすべて終わるのを待たずに次の処理を行うことになるので大抵速いです。特にこのご時世マルチコアなので本当に並列化ができて速いです。良い時代になりました。

専用コマンドがあったらそれを使う

専用コマンドがあったらそれを使うというのはシンプルですが有効な方法です。Shellスクリプトで使うコマンドは基本的にCで実装されているものが多く、まあ速いです。機能が少ない(用途に特化している)物のほうが大抵速いので専用のコマンドがあればそっちを使ったほうがいいです。例えば高機能なcutとしてawkを使うことはあると思いますが、cutで足りる場合にはcutの方が速かったりします。awkが悪者みたいですがshellでloopの処理があったときにawkに置き換えて高速化できることは良くあります。

while readはawkで書きなおしたほうが大体速い(こんなこと書くと怒られうる)

あんまり滅多なこと書くとShell界隈のいろんな人に怒られますが、while readはawkで書きなおしたほうが簡潔にかけて速いことが多いと思っているので僕はawkを使います。僕は基本的にawkが好きな人です。awkはShellスクリプトではない派の人もいるので気をつけましょう。while内部の変数は扱いに注意が必要なので僕はawkを推します。


$ time seq 100000 | awk '$0=(t=$1+t)' | tail -n 1
5000050000

real 0m0.102s
user 0m0.103s
sys 0m0.007s

$ time seq 100000 | while read a ;do echo $((t+=$a));done | tail -n 1
5000050000

real 0m1.712s
user 0m1.490s
sys 0m0.557s

11日目はzayarwinttunさんです。

2015年11月12日木曜日

openboxでディスプレイの輝度を変える話

ノートパソコンに付いている輝度調整ボタンをopenboxで使うという話です。openboxではショートカットキーはrc.xmlに追記する形で設定をします。適当な輝度調整用のプログラムの呼び出しをショートカットキーに割り当てるという形になります。今回使う輝度調整用のプログラムはxbacklightです。

設定例はhttps://wiki.archlinuxjp.org/index.php/Openbox にありますが、何故か手元の環境ではxbacklightで+/-を使うと現在の値からの割合での調整になっています。使っているノートパソコンの仕様なのかxbacklightの仕様なのかよくわかりませんが、これでは一回0%まで輝度を下げてしまうと+を使っても0%から増えません。ということでこれをそのまま輝度ボタンに割り振るとまずいです。そこで次のような感じで書きました。

見ての通り非常に短いのでrc.xmlに直接書こうかと思いましたがワンライナーはrc.xmlに書けないようなので別にスクリプトを作って呼び出す形にしました。

2015年9月11日金曜日

chiselで色々なデータを扱ってみた話

chiselというscala baseのalt HDLがあります。それなりに速く動く回路を書けるらしいのでそれなりにいいみたいです。テストもscalaで書けて、C++のテストコードを生成したり、verilogの回路をchiselはscalaのSDLですが、scalaとは別にchiselのレイヤーで型推論器を持っていたりなかなかしっかりしているみたいです。scalaのインスタンスがchiselの型に対応して、chiselのデータもscalaのインスタンスの内部に格納されているっぽいですがその辺は正直良く理解していないので迂闊なことは書かないでおきます。

ラッチ

chiselにはRegというclassというか型があります。Regはハードウェアで言うところのラッチに相当するものを生成するらしく、クロックのライジングエッヂのタイミングで入力データが書き込まれるっぽいです。テストを含めて書いてみたものはこれです。

メモリ

chiselにはメモリを表わすMemという型があります。シーケンシャルかどうかはパラメタで決められてまあメモリっていう感じです。少なくともテストの上ではラッチと同様にライジングエッヂのタイミングで動作するみたいです。テストを含めて書いてみたものはこれです。

Vec型

ところでchiselにはVec型という配列のようなデータ構造があります。ハードウェア的にはROMなんでしょうか。良くわかりません。実はVecのラッチのReg (Vec)型を作ることもできます。これがなんとMemと大体同じような動作をします。コードはこれです。ちなみにverilogを生成すると結構違う感じになります。verilogは良くわからないのでどっちが速いかは何とも言えないですが、Memの方が綺麗なコードが吐かれました。

二次元配列のようなもの

じゃあMemとVecを組み合わせたりすると二次元配列 (rowとcolを指定するメモリ)の様なものを扱えるかという気がしますが、なぜかこれは上手くいきませんでした。 (これがコード)じゃあReg (Vec (Vec))は動くのかというとこっちは動きました(これがコード)。世の中よくわからないものですが、コードの感じから言ってMemを使って線形なメモリとして扱うのが良さそうです。

2015年8月18日火曜日

Runge-Kutta法のまとめ

差分法 (オイラー法) とからめてRunge-Kutta法についてまとめてみた。近似次数とか刻み幅制御とかは省いている。

差分法

一階微分方程式の近似解を求めることを考える。

前進差分法

hが十分小さい場合
$$ f'(x) \approx \frac{f (x+h) - f(x)}{h} $$
と近似することができる。これを用いると

$$ f (x+h) \approx f (x) + h f'(x) $$
と近似することができる。例えば$f' (x) = f (x),f (0) = 1$の場合

$$ f (x+h) \approx (1 + h)f (x) $$

となるので$f (nh) \approx (1+h)^n$となる。厳密解の$f (x) = e^x$と比較してplotするとこんな感じになる。

後退差分法

$f'$は
$$ f'(x) \approx \frac{f (x) - f(x-h)}{h} $$
と近似することもできる。これを用いると

$$ f (x+h) \approx f (x) + h f'(x+h) $$

と近似することができる。例えば$f' (x) = f (x),f (0) = 1$の場合

$$ f (x+h) \approx \frac{1}{1-h}f (x) $$

となる。この場合は容易に解けるが、右辺に$(x+h)$に依存した項が出てくるので一般には連立方程式を解く必要がある。このような解法を陰的解法と呼ぶ。 逆に前進差分法の様に順々に代入していけば求まる解法を陽的解法と呼ぶ。

Runge-Kutta法

差分法では傾きを一次で近似した解法であった。 (テイラー展開と比較すると示せる。)更に良い精度で傾きを近似する方法として (一般の) Runge-Kutta法がある。
Runge-Kutta法は、$f (x,y) = y'(x)$としたときにパラメタとして$s$次行列$A$,$s$次ベクトル$\mathbf{b},\mathbf{c}$を与えて

$$\begin{align*} y (x_0 + (n + 1)h) &= y (x_0 + nh) + h \sum_{i=1}^{s} b_i k_i\\ k_i &= f (x_n + c_ih,y_n + h\sum_{j=1}^{s}a_{ij}k_j) \end{align*}$$

として近似する方法である。通例$c_i = \sum_{j=1}^{s}a_ij,\sum_{i=1}^{s}b_i = 1$となるようにする。Butcher配列

$$\begin{equation*} \begin{array}{c | c} \mathbf{c}& A\\ \hline & \mathbf{b}\\ \end{array} \end{equation*}$$
を用いてRunge-Kutta法を表わすことができる。

例えば前進差分法は

$$\begin{equation*} \begin{array}{c | c} 0& 0\\ \hline & 1\\ \end{array} \end{equation*}$$

であり、後退差分法は

$$\begin{equation*} \begin{array}{c | c} 1& 1\\ \hline & 1\\ \end{array} \end{equation*}$$

である。

RK4 (古典的Runge-Kutta法)

Runge-Kutta法のなかでも4次のRunge-Kutta法

$$\begin{equation*} \begin{array}{c | c c c c} 0& 0 &0 &0 &0\\ \frac{1}{2} & \frac{1}{2}& 0& 0& 0\\ \frac{1}{2} & 0&\frac{1}{2}& 0& 0\\ 1& 0& 0& 1& 0\\\hline & \frac{1}{6}& \frac{1}{3}& \frac{1}{3}& \frac{1}{6} \end{array} \end{equation*}$$

精度が良く一般に良く用いられる。

RK4のように$A$が狭義下三角行列になる場合陽的に解くことができるので陽的Runge-Kutta法と呼ばれる。 そうでない場合を陰的Runge-Kutta法と呼ぶ。陰的Runge-Kutta法の方が陽的Runge-Kutta法よりも計算は煩雑だが少ない段数でより近似精度が良いという事実が知られている。

2015年8月1日土曜日

HaskellでGraphvizのデータをパースする

Graphvizとはオープンソースのグラフ描画ソフトでDOT言語というグラフを表現するための言語を読み込んでグラフを描画してpngやepsなどの形式で出力することができます。このあたりの詳細はWikipediaを参照したりするとわかると思います。

Haskellでグラフの処理をするときの入出力形式としてDOT言語が使えると便利です。Happyとかを使って自分でパーサを書いても問題ないですが、ライブラリにパーサがあるならそれを使ったほうが楽ですね。Haskell用のGraphvizのライブラリgraphvizがあります。グラフの描画などGrapvhvizを用いた諸々の処理を行えるらしいですが、今回はこれをパーサとして使います。使うのはparseDotGraphで、Textを入力として受け取ってParseDotReprを出力します。ParseDotReprはグラフの表現の型クラスで、例えばDotGraphでラベルとしてStringを用いる場合はparseDotGraphの出力の型がDotGraph Stringに推論されるように型を記述するとよいです。graphNodes (parseDotGraph $ pack "digraph {a -> b}" ::DotGraph String)とか書くとこのグラフのノードのリストが得られます。

2015年6月22日月曜日

PAPIでのFLOPSの計算の方法

PAPI(Performance API)というものがあります。キャッシュヒット率、ミス率、FLOPS、命令の呼び出し回数などを計測できるAPIで、CPUの機能を使います。命令の呼び出し回数を計測できるので、今回は浮動小数点演算の呼び出し回数を計測してFLOPSを計測するのに使おうと思います

PAPIには高レベルAPIと低レベルAPIがあります。高レベルAPIにはまさにFLOPSを計算するためのPAPI_flops関数があります。これを計測したい処理の前後に挟むだけでFLOPSを取得できます。ただし、複数回計測する場合は値がリセットされないので入力時に-1を入れておく必要があります。

普段使う分には高レベルAPIを使えばいいですが高レベルAPIは並列化に対応していないのでopenmpなどで並列化したときは低レベルAPIを使う必要があります。低レベルAPIを使う場合は時間と実行命令数からFLOPSを計算する過程を自分で実装する必要があります。

参考

2015年6月21日日曜日

シュトラッセンのアルゴリズムの計算量

シュトラッセンのアルゴリズムという行列積を求めるアルゴリズムがあります。行列積を求めるときに部分行列の積を求めることによって演算回数が少なくて済みます。通常の行列積の演算は$O(n^3)$ですが、シュトラッセンのアルゴリズムを使うと$O(n^{\log_2 7})$で行えます。巷には日本語だとシュトラッセンのアルゴリズムの説明はちらほらありますが計算量について言及している説明はほとんどないのでしてみようと思います。

サイズNの正方行列の掛け算を考える。因みにシュトラッセンのアルゴリズムを使うときは正方行列でなかったら端に0だけの行or列を追加して無理やり正方行列にする必要がある。Nが2のべき乗でない場合は漸化式が途中で計算できなくなるので途中で諦めるか最初から諦める必要があると思う。以下Nは2のべき乗ということにする。サイズNの正方行列のシュトラッセンのアルゴリズムによる計算時間を$T(N)$と書く。シュトラッセンのアルゴリズムでは半分のサイズの行列積を7回求める。行列和の計算量は$O(N^2)$なので $$ T(N) = 7 T(N/2) + O(N^2) $$ となる。マスターの定理を使うと$O(N^{\log_2 7})$となる(らしい)。

マスターの定理はよくわからないのでマスターの定理を使わないで$T(N) = O(N^{\log_2 7})$を示す。 $$2^n := N$$ として $$U(n) := T(N)$$ と書くと $$U(n) = 7 U(n-1) + O(2^{2n})$$ となる。漸化式を解く感じで整理すると \begin{align*} &\frac{U(n)}{2^{2n}} + \frac{4}{3} c \leq \frac{7}{4} \left(\frac{U(n-1)}{2^{2n-2}} + \frac{4}{3}c\right)\\ &\leq \left(\frac{7}{4}\right) ^ n \left( U (0) + \frac{4}{3}c\right)\\ &\leq \left(\frac{7}{4}\right) ^ n c' (c'はなんかしらの定数) \end{align*} となるので $$ \frac{T(N)}{N^2} \leq c'\left(\frac{7}{4}\right) ^ {\log_2 N} $$ \begin{align*} T(N) &\leq c' N^2 \left(\frac{7}{4}\right) ^ {\log_2 N} \\ &= c'7 ^ {\log_2 N}\\ &= c'N ^ {\log_2 7} \end{align*} となり $$ T(N) = O(N^{\log_2 7}) $$ が示せた。

計算量解析をしてみたけれど数値計算的にはシュトラッセンのアルゴリズムはほとんど使われない。というのも

  • オーダは小さいけど定数が大きいから普通はシュトラッセンのアルゴリズムの方が時間がかかる
  • 疎行列だともっといい方法がある
  • 部分行列の計算のところでメモリを一杯食う
  • 足し引きを一杯行うから浮動小数点演算の誤差が一杯乗る
だからだそうだ。

参考

2015年5月14日木曜日

Evernoteとmarkdownで幸せになった話

最近 (と言っても結構前からほそぼそと使っていましたが) ノートや簡単なレポートやこのブログを書く環境にevernoteを組み込んでみたので書いてみる。

今までは主にmarkdownで書いてpandocをかましてlatex経由でpdf化していました。mediawikiを出力したかったりhtmlが欲しかったりslideにしたかったりということもありましたが、全部pandocで変換していました。pandocは素晴らしい。基本的に入力はemacsで行っています。この方法は特に悪くないんですが、しいて言うなら全部localのパソコンで完結しているので電車の中でスマホで内容を確認したり推敲したりということが出来ません。以前はDropboxを使っていましたが、容量がいっぱいになってしまって最近使ってませんでした。

そこでevernoteです!!evernoteはスマホとパソコンとの間での連携が非常に簡単です。というかスマホのクライアントが優秀です。スマホからテキストをいじる用途ではDropboxよりずっと使いやすいです。以前からevernoteは使っていましたが、evernoteはemacsで編集できません。 (余談ですが、以前はevernote-modeというものがあって割と便利でしたが今はもうAPIが変わったとかで動きません。)emacsで編集できないとなるとストレスが貯まります。このままでは死んでしまいます。

そこでgeeknoteというものがあります。geeknoteはevernoteをcliで操作できるすぐれものです。好きなエディタでevernoteを編集することができます。勿論editorにはemacsclientを登録します。そして、なんとgeeknoteではevernoteをmarkdownにして扱います。つまり、

スマホ -- evernoteアプリ --- evernote --- emacs --- markdown --- latex,html,...

という風にevernoteで何でも書けます。素晴らしいこれで何でもできますね。ちなみにこのブログを書く際にはgeeknote経由でmarkdownで出力したあとstackedit経由でBloggerに投稿しています。

2015年5月13日水曜日

FlycheckでVHDLのシンタックスチェック

FlycheckはEmacsでコードを書いている途中にリアルタイムでシンタックスチェックやコーディングスタイルのチェックを行ってくれる便利な拡張機能です。Flymakeと似ていますが、デフォルトで多くのシンタックスチェックが書かれています。僕は以前はFlymakeを使っていましたが、大体の言語でデフォルトの設定割とつかえるということと、自分でruleを書くのもFlymakeほどではないがそこまで難しくないのでFlycheckに乗り換えました。

僕は学科の課題などでVHDLを書くことがありますが、VHDLは比較的マイナーな言語だからか、無料でない処理系が多いからかFlycheckのルールがありませんでした。だから作りました。GHDLのシンタックスチェック機能を使いました。

最後のadd-hookのところでvhdl-modeでデフォルトでvhdl-ghdlのcheckerが使われるようになります。

2015年3月7日土曜日

Stumpwmを使おうと思ったけどやめた話

先日 http://sourceforge.jp/projects/tilingwm/wiki/FrontPage やらawesomeを使っている友人やらに感化されてタイル型ディスプレイマネージャーを使ってみようとおもった。その中でもStumpwmが(Emacsみたいに)拡張性が高そうで良さそうに感じたのでしばらく軽く使ってみた。そこで起こったいくつかの問題とか思ったことを書いてみる。

日本語入力がうまく行かない -> cjkみたいなimを使わないと生きていけない人の情報が少ない

ibus+mozcをここ何年か使っている。別にibusにこだわりがあるわけではないがとりあえずibusを使ってみようと思ったらどうにも入力ができていない。stumpwmについてはwikiとかに色々な情報があるにはあるが如何せんimを使う文化圏の人の情報が多くなくてibusとかについては調べてもだいぶ古い情報しか見つからなかった。やはり利用者がそれなりにいるor日本語ユーザの書いたドキュメントがある環境でないと辛いなと感じたときでした。(そういえば僕はこんな感じの理由で永らくVine Linuxを使ってましたね。今はArch Linuxを使っていますが。)

Commom Lisp環境についての知識がないとそんなに幸せになれない

僕はLisperではないです。Lisperではないですが、学科の課題でScheme処理系をSchemeで書いたりEmacs Lispの設定を書いたりするくらいにはLispは書けます。Common Lispも昔若干書いたことがあるので書けなくはないです。だからStumpwmを使おうかと思いました。しかしQuicklispとかasdfとかそういったCommon Lispを取り巻く周辺事情については全く知りません。Stumpwmにmoduleを追加したりしようと思ったらどうしてもその辺の知識が必要になるようです。Lisperになりたいなら勿論こういったことを勉強すればいいと思いますが、僕はLisperではなくて結局Unix文化圏の人でLispも書けなくはないくらいの人です。Window Managerを使うためにCommon Lispの世界に足を突っ込むのは割に合うのかと考えてしまった瞬間に僕は colon quitしてgnomeを開きました。gnomeはJavaScriptを知らなくても動きます。

つまり何を思ったか

僕はawkとかCとかbashとかが好きなUnix文化圏の人であってLisp文化圏の人ではない(と思っている)ので慣れない世界のものを中途半端に取り入れるのは辛いなあと。(そのくせEmacsは使ってるんだよなあ…)面白そうではあるのでしっかり時間が取れるときにちゃんと勉強したいなあと。

2015年2月20日金曜日

.emacs.dを書き直す会 第一弾

ここのところ.emacs.dの秘伝のタレ化が進んで非常に汚く、重くなってきたので.emacs.dを書き直した。第一弾として特定のモードに依らない設定を書いた。

package.el(melpa)

Emacs 24からデフォルトで入るようになったpackage管理。非常に便利なので引き続き使う。公式レポジトリだけだと足りないのでmelpaも追加する。marmaladeは以前は入れていたがひとまず入れないでみる。

(require 'package)
(add-to-list 'package-archives '("melpa" . "http://melpa.milkbox.net/packages/") t)
(package-initialize) 

init-loader

設定ファイルを分散して書けるようなもの。.emacs.d/inits以下に設定ファイルを書くと、辞書順に呼ばれる。

(require 'init-loader)
(init-loader-load)

twittering-mode

emacsからtweetするために必要。普段も使っていた。いつでも呟ける様にキーバインドを設定しておく。

(require 'twittering-mode)
(global-set-key "\C-xt" 'twittering-update-status-interactive)

appearance

特定のモードではなく背景色とかそういうもの。以前はXのクリップボードとemacsのコピーが同期されるようなものも書いていたが気がついたら勝手に同期されるようになっていた。背景色はthemeとかを使おうかと思ったが良さそうな色がなかったので”秘伝のタレ”を取り出した。

;;背景色の設定
( when window-system 
  (custom-set-faces
   '(default ((t
           (:background "#000040" :foreground "#e0e0e0")
           )))
   '(cursor ((((class color)
           (background dark))
          (:background "#00AA00"))
         (((class color)
           (background light))
          (:background "#999999"))
         (t ())
         ))))

;; フレーム透過設定
(add-to-list 'default-frame-alist '(alpha . (0.75 0.75)))
;; tool-barを消す
(tool-bar-mode 0)

;; C-h でbackspaceにする
(global-set-key "\C-h" 'delete-backward-char)

;; 括弧を強調表示
(show-paren-mode t)

helm

色々と便利な子。こいつがいないともうファイルも開けない。設定がよくわからなかったのでここも”秘伝のタレ”。

(when (require 'helm-config nil t)
  (helm-mode 1)

  (define-key global-map (kbd "M-x")     'helm-M-x)
  (define-key global-map (kbd "C-x C-f") 'helm-find-files)
  (define-key global-map (kbd "C-x C-r") 'helm-recentf)
  (define-key global-map (kbd "M-y")     'helm-show-kill-ring)
  (define-key global-map (kbd "C-c i")   'helm-imenu)
  (define-key global-map (kbd "C-x b")   'helm-buffers-list)

  (define-key helm-map (kbd "C-h") 'delete-backward-char)
  (define-key helm-find-files-map (kbd "C-h") 'delete-backward-char)
  (define-key helm-find-files-map (kbd "TAB") 'helm-execute-persistent-action)
  (define-key helm-find-files-map (kbd "C-z") 'helm-select-action)
  (define-key helm-read-file-map (kbd "TAB") 'helm-execute-persistent-action)

  ;; Disable helm in some functions
  (add-to-list 'helm-completing-read-handlers-alist '(find-alternate-file . nil))

  ;; Emulate `kill-line' in helm minibuffer
  (setq helm-delete-minibuffer-contents-from-point t)
  (defadvice helm-delete-minibuffer-contents (before helm-emulate-kill-line activate)
    "Emulate `kill-line' in helm minibuffer"
    (kill-new (buffer-substring (point) (field-end))))

  (defadvice helm-ff-kill-or-find-buffer-fname (around execute-only-if-exist activate)
    "Execute command only if CANDIDATE exists"
    (when (file-exists-p candidate)
      ad-do-it))

  (defadvice helm-ff-transform-fname-for-completion (around my-transform activate)
    "Transform the pattern to reflect my intention"
    (let* ((pattern (ad-get-arg 0))
           (input-pattern (file-name-nondirectory pattern))
           (dirname (file-name-directory pattern)))
      (setq input-pattern (replace-regexp-in-string "\\." "\\\\." input-pattern))
      (setq ad-return-value
            (concat dirname
                    (if (string-match "^\\^" input-pattern)
                        ;; '^' is a pattern for basename
                        ;; and not required because the directory name is prepended
                        (substring input-pattern 1)
                      (concat ".*" input-pattern)))))))

日本語入力

日本語入力は普段\C-\ではなく[変換/無変換]を使っているのでその辺の設定をまずする。(2015/02/25修正)Emacs 24.3以降でinactive-input-methodがobsoleteになったので修正した。

;; IMEの設定
(define-key global-map [zenkaku-hankaku] 'toggle-input-method)
;; 変換キーでIME ON
(define-key global-map [henkan]
  (lambda ()
    (interactive)
    (if current-input-method (deactivate-input-method))
    (toggle-input-method)))

;; 無変換キーでIME OFF
(define-key global-map [muhenkan]
  (lambda ()
    (interactive)
    (deactivate-input-method)))

次にmozcの設定。

;; mozc
(require 'mozc)
(setq default-input-method "japanese-mozc")

auto-complete

コードを書くときとかにいい感じで補完してくれる子。emacs lispを書くにもこれがないとつらい。Cとか特定の言語に依存する設定は後で。

(require 'auto-complete)
(require 'auto-complete-config)
(ac-config-default)

magit

emacsからgitを操作するもの。そう言えば結構前は標準で入るVCからgitを使っていたが、VCはgitに特化しておらず色々とやりたいことが出来なかったのでmagitにした。git rebase -iに当たるものをmagit-statusじゃない場所からやるとか使用頻度の少ないキーバインドは消した。

(require 'magit)
(global-set-key "\C-xgs" 'magit-status)
(global-set-key "\C-xgb" 'magit-branch-manager)

iedit

共通部分を全部まとめて編集するもの。変数名とかを変えるときに便利。”\C-;”と書きたいがそれだと動かないのでこんな感じで書いた。

;; iedit
(global-set-key [?\C-;] 'iedit-mode)

Written with StackEdit.

2015年1月1日木曜日

AwkLike classを作った

Google Apps Scriptでawkみたいな行志向でパターンマッチして動かしたくなったのでAwkLikeを作った。 こっから先は後日追記する。