必要変数など列挙

SDPA-C のうち、いくつか必要な変数などを列挙しておいた。
ここからは、sdpa_struct.cpp, sdpa_algebra.cpp, sdpa_parts.cpp の順番に進めていく形になると考えている。


これでようやくプログラミングに入れそうである。
今まで長かったが、頭の中にだいたいのデザインができたので、あとは組んでいくだけである。

今日の作業内容：SDPA-C 3h
今日のランチ：つかさ ぶりの照り焼き
明日の予測作業時間： 4h

データ構造を検討

SDPA-C の実装だが、入力行列のデータ構造を大きく変えた方が良さそうである。
SDPA の場合には、入力行列を列単位で計算することがないので特に問題ないが、SDPA-C の場合はボトルネックである Schur 補完行列の計算をするときに入力行列を列単位でアクセスする必要がある。

ここで問題なのが、メモリを多く取って高速化するか、あるいは計算時間は多少遅くなってもメモリを節約するか、である。
今回は、メモリ消費量を多くしても高速化する方向で行こうかと考えてもいる。

今日の作業内容：データ構造検討 2h
今日のランチ：らく　鶏の照り焼き定食
明日の予測作業時間：4h

アルゴリズム再検討

CHOLMOD のほうもだいぶ理解できてきたので、SDPA-C のアルゴリズムの再検討を始めた。
できるだけ無駄をなくすように計算すべきなので、どこが無駄になるかが検討の中心である。

たとえば、inv(X)は sparse で X は dense になるので、X を全く持たないで計算できれば、そのほうが効率的であることも分かるが、Xを全く持たない場合には、計算途中で inv(X) の逆計算を繰り返す必要があり、結局のところ X を clique で保持した方が速い。
こういった感じの部分がまだまだたくさんあるようなので、そういったところを再検討している。

おそらく、入力行列はこれまでの下三角だけとは異なり、上下の両方を持つことになると考えている。

ただ、ステップ長の計算のところがまだ不確定要素なので、このあたりはまだまだ考える必要がある。


今日の作業内容：SDPA-C 再検討 3h
今日のランチ：サイゼリヤ きのこのスパゲッティセット
明日の予測作業時間：2h

clique は取り出せるようになった

CHOLMOD からの情報で clique は取り出せるようになった。
問題は、ordering のほうで、行列の各要素について何番目に移動されているかを保持しておくと計算効率が向上するが、その分メモリを大きく使ってしまう。
このあたりのバランスが難しそうだ。

なにはともあれ、これで必要な計算を列挙して、データ構造を再設定することになる。

今日の作業内容：CHOLMOD 検討
今日のランチ：やぶ A定食(さわらの西京焼き）
明日の予測作業時間：5h

CHOLMOD から clique 情報を取り出せるようになった

ついに CHOLMOD から clique 情報を計算する方法が分かった。

そのまえに、Matlab のいくつかのテクニックについて。
etreeplot で elimination tree の表示ができるが、このままだとノード番号が表示されず、若干使いづらい。
これは treelayout で情報を得ると表示できる。

p = amd(A);
[parent,post] = etree(A(p,p));
figure(1);
etreeplot(A(p(post),p(post)));
hold on;

[x,y] = treelayout(parent,post);
for i=1:n
    text(x(post(i)),y(post(i)),sprintf('%d',i));
end
hold off;
とすれば、elimination tree のノード番号を表示できる。


また、Matlab では、system 関数でコマンドを実行できるが、この2つめで出力を取り込むことができる。
たとえば、a.txt の内容を取り込むのであれば、fopen などをしなくても
[status,result] = system('cat a.txt');
とすればOKである。
a.txt に
A = [1,2; 3,4];
と書いておいてから

[status,result] = system('cat a.txt');

eval(result);

とすると、もちろん、 A に行列が格納される。

このあたりは、他の C 言語のプログラムや awk などで処理した結果を取り込むのに便利である。

今日の本題の clique 情報だが、基本的に supernode ひとつにつき clique をひとつ生成すればよく、この場合は supernode に含まれる頂点と、それらに隣接する頂点で clique を構成することができる。

ただし、maximal clique をまじめに取った場合に比較して、より fill-in が多くなってしまう。

このあたりは BLAS で高速化されるはずだが、どの程度の高速化になるかはやってみないと分からないかもしれない。

今日の作業内容： CHOLMOD 確認 3h

今日のランチ：たちばな カツオ丼

明日の予測作業時間：4h

supernode と maximal cliques

今回は、
http://www.cs.odu.edu/~pothen/Papers/ds.pdf
にあった、「Elimination Structures in Scientific  Computing」という論文をチェックして supernode と maximal cliques の関係を確認しようとしたが、これだとあまり有益な情報を得ることができなかった。
むしろ metis 関係をみたほうがいいのかもしれない。

あと、乱数で行列を生成して、それで cholmod_demo から探っていくのもありかもしれない。

今日の作業内容：CHOLMOD 確認 3h
今日のランチ：角笛 ビーフステーキ丼
明日の予測作業時間：5h

CHOLMOD の permutation

CHOLMOD の結果から clique の情報を取り出すのが次の方向だが、これが意外と難しい。
今日までに分かった点では、cholmod_factor の中に入っている Perm の順番と実際にコレスキー分解に用いている順番が異なる、という点である。

Perm で求まるのは（AMD で計算していると仮定したとき）、Matlab 表記であれば
Perm = amd(A);
である。
これに対して、コレスキー分解に用いられているのは、
Perm = amd(A);
[parent, post] = etree(A(Perm,Perm));
Perm2 = Perm(post);
である。
つまり、elimination tree の post-ordering で並び変えた後にコレスキー分解を行っている。

この順番に並び変えた後に、単純に極大クリークを求めてみたが、この情報とCHOLMOD から得られる supernode の情報はかなり異なる。
post-ordering から supernode を生成する過程と、その過程と極大クリークの関係を理解する必要がありそうだ。


今日の作業内容：CHOLMOD 確認 3h
今日のランチ：らく 焼き魚定食
明日の予測作業時間：2h

supernode の情報取り出し

CHOLMOD の supernode の情報を直接取り出す方法がマニュアルに書かれていないため、ソースを参照しながら解読を進めている。
cholmod_common* cm; の場合は、cm->print を 4 以上にすると cholmod_print_factor で詳細な情報が得られるので、これで確認をしながら理解を進めている。

supernode の方法については、
http://www.cise.ufl.edu/tr/DOC/REP-2006-109.pdf
にざっくりと書かれている。

cholmod_print_factor で得られた情報と、Matlab の etree で表示される etree の情報は対応させることができるので、このあたりからどのように clique を取りだすかに検討が移る。

今日の作業内容：CHOLMOD 確認 3h
今日のランチ：味庵　鶏のから揚げ
明日の予測作業時間：2h

CHOLMODで supernode に切り替える

CHOLMOD は実際にコレスキー分解をする前にパラメータを変更すると、強制的にsupernode に切り替えることができる。
これは、cholmod_common にある supernode 変数を２以上にすると強制できる。

今日は、これがうまく行かず、２時間程度手こずった。
結局、cholmod_start で初期化した後に supernode を指定したにもかかわらず、もう一度 cholmod で初期化しなおしてしまっていた。
バグも分かると簡単だが、分かるまでは結構時間がかかる。

今日の作業内容：CHOLMOD 3h
今日のランチ：シッダルータ　ベジタブルカレー
明日の予測作業時間：5h

mingw64 と mingw32 だと設定が違うのか？

昨日の続きで SDPA-M の 64bit Windows Matlab 用のものをコンパイルしているが、以前の 32bit 用と同じように 64bit でコンパイルしたら、libgcc の DLL が必要だと言われてしまった。

ひょっとしたら、mingw32 と mingw64 で異なるのかもしれない。

（他にも gcc のバージョンアップに伴う変化の可能性もあり。）

いずれにしても、mingw64 で DLL なしで実行できる exe を作るためには、 gcc, g++ などでのリンクの時点で -static をつければよい。

これで、

$ x86_64-w64-mingw32-objdump -p sdpa.exe | grep DLL

などとすると、libgcc などの DLL がリストアップされなくなる。

ちなみに、このコマンドは exe や dll で必要な DLL をリストアップするためのものであり、KERNEL などWindows 標準の DLL なども表示される。

また、Matlab の mex ファイル拡張子である mexw32, mexw64 がついているファイルは、実質的には DLL なので、 gcc に -shared -static で DLL として作成して拡張子を変更すれば Debian のクロスコンパイラでも作れるのである（だから SDPA-M がコンパイルできている）。

mexw32, mexw64 と異なった拡張子になっているのは、たとえば、sample.dll だと 64bit,32bit の両方を同じファイルに格納できないが sample.mexw32, sample.mexw64 だと両方を同時に格納できるためと考えられる。

このあたりは、linux 上の mex ファイルの拡張子が so でないのも同じ理由と推測できる。

今日コンパイルした SDPA-M は別途動作確認をしてもらう予定である。

今日の作業内容：SDPA-M 3h + CHOLMOD 2h

今日のランチ：四川 鶏ネギ炒め

明日の予測作業時間：4h

SDPAM 64bit-Windows版に挑戦中

前にやろうとしてそのままになっていた 64bit-Windows 用の SDPA-M のコンパイルに再挑戦している。
よく分からないのは、 pthread-win32, OpenBLAS が 64bit でも対応できるのか、という点と、Debian 上でのクロスコンパイルで Matlab の mex ファイルを生成できるか、というあたりである。

OpenBLAS については、コマンドを間違えたりなどして、何回かコンパイルに失敗しているので、今日の夜をコンパイルの時間に割り当てる予定である。

今日の作業内容：SDPA-M 2h
今日のランチ：角笛 かじきの生姜焼き
明日の予測作業時間：5h

CHOLMOD デモプログラムのコンパイル

CHOLMOD のマニュアルは、目を通し終わった。
今後の実装でキーとなる項目としては、

16. Check Modules
17. Cholesky Modules (特に etree, rowcolcounts)
18. Supernodal Modules
のあたりである。
ただ、Cholesky 分解をする 17 のあたりの関数にプロトタイプが書かれていないので、これはヘッダーファイルで確認が必要そうだ。

ソースファイルの中にあるデモファイルを Debian でコンパイルしたが、必要なパッケージが入っている場合には、
$ cd CHOLMOD/Demo
$ gcc -c -I/usr/include/suitesparse/ cholmod_demo.c

$ gcc -o cholmod_demo cholmod_demo.o -lcholmod
だけでコンパイルできたので、簡単であった。(-lcholmod はかなりの関数を含んでいるのかもしれない。）
実行は
$ ./cholmod_demo Matrix/bcsstk01.tri
で実行できる。
次は、このデモの内容を確認する段階となりそうだ。

今日の作業内容：CHOLMOD チェック 3h
今日のランチ：味庵 バンバンジー
明日の予測作業時間：1h

CHOLMODのマニュアルを８３ページまで

今日は、CHOLMOD のマニュアルの８３ページまでをチェックした。
この中で重要そうなのは、cholmod_factor_structのsupernodalの情報といったあたりである。
この情報が完全に内部の情報なのか、あるいはユーザ側でアクセスできるのか、について、小規模プログラムを作って確認したほうがいいかと考えている。

今日の作業内容：マニュアル確認 3h
今日のランチ：らく　鶏の照り焼き定食
明日の予測作業時間：4h

CHOLMOD のマニュアル読み（３６ページまで）

再度、CHOLMODのマニュアル読みを進めている。
基本的に必要な計算がある程度分かっているので、そこを集中的に読むつもりである。
現在は３６ページまで来たが、ここまでは Matlab での使い方になっているが、このあたりから C 言語などの呼び出しとなっており、ここからが本番ともいえる。

今日の作業内容：資料確認 2h + CHOLMOD マニュアルチェック 2h
今日のランチ：つかさ　マグロのづけ丼
明日の予測作業時間：3h

資料作成

１５分程度の発表を急きょすることになったので、資料作成を行った。
１５分でも、分かりやすくするためにいろいろと工夫などを考えたので、結局５時間以上かかっての資料作成となった。
これで、自分で発表練習をしてみて、必要な時間などを考えるのがよさそうだ。

今日の作業内容：資料作成 5h
今日のランチ：サイゼリヤ ハンバーグとパンチェッタ
明日の予測作業時間：5h

SDPA-C に出てくる変数のまとめ

最近読んできていた SDPA-C だが、変数の扱い方が特殊なので、そのあたりをまとめておく。

X: clique で保持
dX~tilde : sparse で保持（ただし対称ではない）
dX: clique で保持（だから X と足せる）
inv(X)のcholesky: extended sparsity で保持: clique で作成した後、sparse に変換

Y: sparse で保持
dY: sparse で保持
Yのcholesky: extended sparsity

今日の作業内容：SDPA-C チェック 3h
今日のランチ：りきちゃん うどんセット
明日の予測作業時間：5h

I-Lin 打ち合わせの続き

昨日に続いて I-Lin との打ち合わせをした。
やはり、noise があるケースをどのように扱うか、という点がポイントになりそうである。
あと、タンパク質データについても、どの情報を活用すべきか、というのが重要でもある。


今日の作業内容：I-Lin 打ち合わせ 2h + SDPA-C チェック 2h
今日のランチ：食堂 カキフライ
明日の予測作業時間:3h

I-Lin 打ち合わせ

台湾から来ている I-Lin Wang と研究についての打ち合わせをした。
今までの研究内容の確認と、そこでどういったことができそうなのかをいろいろと話した。
ただ、今日はあまり時間が取れなかったので、明日も続きの打ち合わせをすることになりそうだ。

今日の作業内容：SDPA-C のソース読み 3h + 打ち合わせ 2h
今日のランチ：味庵 豚角煮とジャガイモの煮込み
明日の予測作業時間：4h

main loop の直前までは来た

SDPA-C のソースの理解は、とりあえず main loop の直前までは来た。
アルゴリズムの関係で aggregate と extended があるが、aggregate は足し算を高速化するために自前データで行っている現状を維持した方がよさそうである。
一方で、extended については、内部のデータを直接アクセスする機会が意外と少ないため、数値演算ライブラリになどに任せるのも一つの手かもしれない。
ただ、どのようなライブラリを利用すべきかは、今後検討する必要があるが。

今日の作業内容：SDPA-C チェック 2h
今日のランチ：らく　焼き魚定食
明日の予測作業時間：4h

SDPA-C の解析中

SDPA-C については、ソースを読むだけでは理解しづらいところばかりが残ってきたので、小さい入力ファイルを作って、それで計算過程を解析している。

ところどころで改良できそうなところもあるが、とりあえずは一通り理解するのを優先することにする。
いまは、2x2 のクリークが２個、という小規模なので、clique tree の効果がよく分かっていない。
このあたりは、一通り理解した後で、もう一回試すことになると考えている。

今日の作業内容：SDPA-C 2h
今日のランチ：つかさ サーモンハラス焼き
明日の予測作業時間：3h