忍者ブログ
日々の研究生活のメモ書きなど
なんか自分の書いたコードが思っているのと違う挙動をしていてちょっと困った・・・・

fprintf(stdout, "hoge ");
みたいな感じで改行を挟まずにループを何個も回して、最後に
fprintf(stdout, "\n");
と改行を入れるコード

実際に動かしてみると、待ってもhogeが表示されない・・・
しかし計算はどんどん進んでいっていることがstderrの情報からはわかる
stdoutは内部に出力バッファを持っていて、printfやfprintfでは改行文字を出力するか出力バッファが一杯になるまで溜め込みます。
改行文字がなくても強制的に吐き出させるのが fflush です。
■ 参考 : fflush(stdout)の意味

とのことらしい
なるほどね・・・
fprintfは、\nが来るまで実際に画面表示がされないというのは今回初めて知った・・・・

今回のコードでstderrは毎回\nが入っていたので、その都度表示されていたのも納得がいった
なおかつ、stderrとstdoutはバッファが違うから、今回のようなことになったんだろうなぁ〜














PR
printf("%%\n");
でOK

\%でできるかな〜と思ったけどダメだった








判別する方法は、
#include <math.h>

if ( isnan(x) ) {
 // これはNaNの場合
} else {
 // NaNでない場合
}

if ( isinf(x) ) {
 // これはInfの場合
} else {
 // Infでない場合
}
じゃあ判別できているのかのテストの方法は、

■ 参考 : infやnanになる条件

log(0)とか、5/0.0とかやってみたけどすべてInfになってしまった・・・

■ 参考 : 8.3 math.h

このままではNaNのテストができないので、
math.hにNANってマクロがあるらしいので、それに対して調べてみた

あとでもうちょっと調べてみたら、sqrt(-1.0)って作り方もあるらしい




以前軽く調べた知識によると、IEEE_754では
NaNもInfも指数部のビット全て1(つまり255)
NaNは仮数部が0以外の任意の値
Infは仮数部が0
らしい

そんなことを知らなくても、上記のisnanとかで調べられるからまぁいっか・・・・










前回書いたMPIの記事はあまりにめちゃくちゃだった・・・
今回もめちゃくちゃです

一応、既存のコードの並列化自体はできて、並列処理はできました
(そもそも既存のコードが遅いので、、まだまだ改善の余地がありますが)

もっと早くする方法について色々と調べて行くと「master/slave方式」というのが見つかった




通常並列化するときは

例えば何かある値をn個計算するとして逐次コードだと、
for(i=0; i<n; i++){
// 何か計算
}
みたいにしてやる


これをncore個のcoreで並列化するときは、サイクリック配置
for(i=0; i<n; i=i+(ncore)){
// 何か計算
// MPI_Gatherで集める または printfとか
}
または、ブロック配置では
for(i=(myrank*ncore); i<(myrank+1)*ncore; i++){
// myrankはそのプロセスのrank
// 何か計算
// MPI_Gatherで集める または printfとか
}
みたいな感じか?(上ので合ってるかはわからんけど)
他にもこれらを組み合わせたブロックサイクリック配置というのがあるらしい

どうやって逐次コードを並列化するかはその計算処理自体に依存するのでケースバイケース





一方でmaster/slaveでは、myrank==0のプロセスは他の監視するためだけに使う
myrank!=0のプロセスは何かの処理や計算をする

言うのは簡単だけど、これを実装するのがなかなか大変・・・
サンプルコードも色々と探したけど、そのまま活用できそうなのが見つからない
もうちょっと書くと、計算に必要なパラメーターをmasterがsend(MPI_Send)して、slaveがreceive(MPI_Recv)して計算を走らせる
計算が終わったらslaveがmasterに値をsendして、masterが値をreceiveする
そのslaveは次の計算ができるので、また走らせて・・・というのを管理していく
masterはいつでも受信できるように
MPI_Recv(&result, 1, MPI_DOUBLE, MPI_ANY_SOURCE, MPI_ANY_TAG, MPI_COMM_WORLD, &status);
のようにANYを使う

たぶんこういう理解で合ってるはず・・・

大まかなmaster/slave方式は以下のコードを見てもらえるとわかるはず
myidという変数で場合分けしてる
■ 参考 : code_mpi_exp/test_mpi5.c

他の例もある(githubでも割と探した)
■ 参考 : Parallel programming: an MPI example

■ 参考 : pi-pp/serie_4/master_slave_template.c

■ 参考 : pi-pp/serie_4/mpi-master-slave-template.c

一つ上のとすごく似ているが微妙に違う

■ 参考 : ValiTDS/valiTDS/gcd.c

ポルトガル語で説明が書いてあるけど、要するに3つの整数の最小公倍数を探す計算らしい


■ 参考 : Master/slave programs in MPI

動くコードがすべて書いてくれている
ありがたい


■ 参考 : /* 1からこの数までの和を求める */

どういう経路で見つけたページか忘れたけどすべてのコードが書いてある
コピペで完全に動いてビビった





他のサンプルコードは、「MPI 並列化 "サンプルコード" C言語 -pdf」とかで検索しても全然見つからない
その検索で見つけた興味深い記事をメモしておく

■ 参考 : Message Passing Interface (MPI) 統合スレ

ネット上にマジでまとまった情報がないから2chに行ってみた
2004年からある老舗だが、あまり活発ではないらしい
とりあえず一度全部さっと目を通しておくとMPIのソースコードには慣れられる


■ 参考 : link集/MPI

2009年で更新が止まっている


■ 参考 : MPI Programming samples

リンク集の中にあった東工大のサンプルコードページ
見つけたサイトの中ではかなり有意義なページ
実際にコピペでは動かないものがあるので注意



■ 参考 : MPI Documents

日本語のpdfもあるが、1997年のものだった
少しだけsample codeがある

■ 参考 : Web pages for MPI and MPE

上記のページのdocumentは少し検索しずらかったり、読みにくい
web page版があった






色々と調べて感じたのが、MPIによる並列化ってまだ今も使われている技術なんだろうか・・・
並列化って今はもっと他のライブラリとかを使ってたりするのだろう
ちょっと不安だけど、まぁもうちょっとやりこんでみる・・・・




(2016/10/10 追記)

今更だけど、すごくよくまとまったページを見つけた

■ 参考 : PCクラスタ超入門 講習会テキスト


この中で特に「MPIによる並列プログラミングの基礎」というのが役に立ちそう
他のドキュメントはまだ見てないけど、高速な並列化をするために重要なことが書かれていそう・・・ まだ見てないけど




(2017/06/18 追記)

たまたま↑のページを見に行ったら、ページがなくなっていた
オーマイガッ

















■ 参考 : MPI_Wtime
double MPI_Wtime( void )











計算に11時間かかるッスって言ったら並列化すればいいじゃん
なんならコードの並列化をやってあげる
という優しい言葉を受けて初めての並列化に手を出してしまった・・・・
わからないことだらけなのでメモしておく



まずはMPIのインストール
省略




■ 参考 : MPI Hello World

上記のサイトを参考にしてテストコードを書いてみる
途中まで読んで英語だったからめんどくなった・・・


■ 参考 : MPI「超」入門(C言語編) - 東京大学
↑難しいすぎる、猿の私には上っ面しかわからんかった

■ 参考 : MPIを用いた並列プログラミングの概要


こっちに日本語の文章があったのでこっちを見て行く

大事なのは「各プロセスは同じことをやるがデータがそれぞれ異なる、大規模なデータを分割し、各部分について各プロセスが計算をする」ということだと思う
並列のためにでかいデータを分割する、というのは初めてなので工夫しないといけないなぁ〜
スライドの中盤にも、「領域分割」としてそのことが書いてあった

・PE : Processing Elementの略、ここではプロセスという意味で使っているらしい

・MPIのプロセス番号は0から始まる、なので8プロセス走らせる場合は0~7が割り当てられる

・MPIのコードのコンパイラはgccではなくて、mpiccというのを使う

・MPIに関連した変数名はMPIで始まっている
 なので、ユーザー独自の変数としてMPIなんとかは使わない方がいい
#include <mpi.h>
まぁ、こういうヘッダーファイルがいるんだろう
インストールがうまくいってPATHがきちんと通っていれば問題ないはず
MPI_Init(&argc,&argv);

MPI_Finalize();
すべての並列化処理は MPI_Init から MPI_finalize の間で書かれなければならない
MPI_Comm_size(MPI_COMM_WORLD,&numprocs);
必須ではないらしい
これで総プロセス数の取得をしているらしい
numprocsがプロセス数
MPI_Comm_rank(MPI_COMM_WORLD,&myid);
起動した各プロセスにはそれぞれ「ランク」と呼ばれる値が割り当てられ、それを取得する関数
ランク値は他のプロセスを指定して通信するときとかに使うらしい


便利そうなのは
MPI_Barrier
すべてのプロセスの同期を取ってくれる関数
すべてのプロセスがこのサブルーチンを通らない限り、次のステップには進まない

と思ったけど「主にデバッグに使う、オーバーヘッドが大きいので実用計算には使わない方が無難」と書いてあるわ・・・・
sprintf(FileName, "%d.txt", MyRank);
のように各プロセスの通し番号を使うこともできる


MPIはCPU core毎にMPIプロセスを実行する
Rank数はコアの数に等しい

hostfileのフォーマットは
node01 slots=4
node02 slots=4
node03 slots=4
node04 slots=4
node05 slots=4
みたいな感じらしい
slotsには各ノードのCPU数を指定する


プロセスごとに通信(送信/受信)とかもできるらしいが、そんな難しいのはまだいらないと思うので・・・
今やりたいのは、単純に1つのコードを複数のプロセスに分担させて並列で走らせたいだけ





(2016/10/1 追記)

■ 参考 : 技術コラム

この記事の第50回〜58回あたりはMPIの並列化について説明されていて、わかりやすかった















ほんとたまに使うだけなので、いつも関数名を忘れるのでメモ
(そうだこのブログはポケモンGOブログではなくて備忘録だった!)
#include <complex.h>
double complex hoge;

double real_part = creal(hoge); //実数部のみ
double imaginary_part = cimag(hoge); //虚数部のみ

double absolute = cabs(hoge); //複素数の絶対値
double phase = carg(hoge); //複素平面上の偏角

double complex conjugate = conj(hoge); //複素数の共役複素数

double complex conjugate = cpow(hoge); //複素数のベキ乗(使うのか?)









3ヶ月前から動かないコードが手元にあった
元々は3年以上前に自分が書いたコードなのだが、なぜか一括で
% make
とするとコンパイルがそこで止まる・・・
特に手を加えたわけでもないのになーぜー

Dropboxで過去のファイルに戻ろうとするも、古すぎて履歴が残っていない(もしかしたら、ファイルをコピーして作ったからとか移動させたからかもしれない)

とりあえず急務なものでもなかったので、放置してた
(実はそれがけっこう便利なプログラムで、周波数と周波数スペクトルの値を渡すと横軸質量、縦軸inspiral rangeのグラフを書いてくれる)

コンパイルが通らない理由がやっとわかって、コードの中からカンマがなくなっていた
「いや、なんでだよ」と言わざるを得ない・・・




コードに一つ改良を加えた
今までpdfしかグラフ出力できなかったが、与えた出力ファイル名から拡張子を判断してグラフを書かせる
今までこういうグラフ処理までさせるときは実行部分をC言語で書いて、グラフ描画をシェルスクリプトと分担していたが、今回のコードはかなり前のなので全部C言語でやっていた

そのためにはまず与えた文字列から拡張子を判断する必要がある
strrchrを使う
strchrも似たような文字列検索だけど、こっちは後ろから検索するので拡張子を見つけられる
そもそも.で検索すれば一発なはずだから、どっちを使っても同じ結果になるはずだけど後ろの情報が欲しいので実行時間を減らすためにもsttchrで
(.が2個以上あったら困るんじゃない?と思うかもしれないけど、たぶんMacでは.が2個あるファイルが作れない)

■ 使い方
cFindが検索する対象の文字列
cKeywordが検索する文字列(int型)
#include
char *strrchr(const char *cFind, int cKeyword);
例えば
cFind="hoge.png";
cKeyword=".";
だとすると、返り値は .png になる
拡張子にするためには 1つポインタを進めればよい(っていう表現で合ってる・・・?)


■ 例文
初めて使って動いたのでたぶんこんな感じなんだろう、知らん
char *cTerm=""; //初期化大事
int cKeyword = '.'; //検索する文字列
cTerm = strrchr( output, cKeyword); // cTermになんか入る

if( cTerm == NULL ) {
// こっちはkeywordが見つからなかったとき
//なんか例外処理
}else{
// keywordが見つかったとき
cTerm++; // これで.の後ろの文字列からなので拡張子になるはず
//なんか処理
}






非常に解りやすい記事があったのでそのメモ

■ 参考 : 2 関数のデータの渡し方




関数に変数を渡しても、その計算はmain関数のメモリとは別の場所で計算される
なので関数に何か変数を渡しても、その値が書き換えられることはない

では逆に書き換えて返してほしい場合はどうするか?
1. 関数の返り値に指定する
2. 関数に変数のアドレスを渡す

1は
int hoge(int a, int b);
みたいな感じ
けど、これじゃあ1つしか変数が返せない・・・


2は
void hoge(int a, int b, int *c, int *d);
みたいにして、aとbはただ関数に渡すだけで、cとdに欲しい値が入って返ってくる
この関数を呼ぶときは
hoge(a, b, &c, &d);
みたいにして、変数のアドレスで渡さないといけない
また、関数hogeの中ではcとdに何かを代入するときは
*c = a+b;
*d = a-b;
みたいにしないとだめ


(もしかしたら間違った理解をしているかもしれないので注意・・・・)









# 今更C言語のことを書くとは・・・

scanf系の関数でなぜオーバーフローが起こる可能性があるのかは、ここに書いてあります
■ 参考 : [迷信] scanf ではバッファオーバーランを防げない


また、その解決方法はこちら・・・
■ 参考 : scanf 系関数での文字列長指定方法

マクロ関数を使って、うまいこと回避しているようです



この記事を読んで初めてマクロ関数というのに触れました
一応知識では知っていましたが
#define SCNs(n) SCNs_(n)
#define SCNs_(n) #n "s"
の#nが何をしているのかわからず調べました
nという渡された変数に#を付けると、その値をそのままそこに展開できるらしい

これを使うと、
#define N 10
char s[N-1];
scanf("%"SCNs(N), s);
で、SCANs(N)の部分が"10s"に置き換わって、うまいこと文字配列の上限値を指定できると・・・

これ、char s[N+1]; の間違い?
よくわからん・・・・



これで
scanf without field width limits can crash with huge input data on some versions of libc.
というwarningがなくなればいいのが・・・・






最近、関数ポインタについてあれこれ理解する必要に迫られたので調べた

関数ポインタを利用して呼び出す関数を動的に変更する

このページの説明が割としっくり来た
関数Aに、さらに関数Bを渡すときにポインタでそれを渡してしまうことで
関数Aの中で使われる関数Bを簡単に置き換えることができる

例えば積分する関数などで(今思い返せばnumerical receipeの中で積分する関数もこれを用いていたような・・・)
積分するための関数をAとする
被積分関数をBとする
関数Aの引数にBを渡す

このときBをxとかcos(x)とか変えれば、積分値もそれにともなって変わるというのが関数ポインタ


という理解であってるのかなぁ?














やりたいことはタイトルの通り
forkというのを使ってみる
詳細な説明は参考にした記事をどうぞ

■ 参考
http://linuxjm.sourceforge.jp/html/LDP_man-pages/man2/wait.2.html



とりあえずサンプルコードを書いてみた
本当に必要なヘッダーがどれかわからん・・・
もしかしたら不要なものを呼んでるかもしれない

このforkを使えば色々と計算させるプログラムで役に立つかも・・・?









これまで4年ほど計算機に触れてきたけど、こんなに便利なものを見過ごしてきたなんて・・・・

とりあえず設定ファイルにPKG_CONFIG_PATHを設定する
自分の場合は
setenv PKG_CONFIG_PATH /usr/local/lib/pkgconfig:/opt/local/lib/pkgconfig/
の2つくらい

これで
$ pkg-config --cflags fftw3
-I/usr/local/include
$ pkg-config --libs fftw3
-L/usr/local/lib -lfftw3 -lm
あとはこれをmakefileの中でうまいこと使えばおーけー。





某ライブラリを使ってプログラムを書いていた
無事コンパイルまでは到達したが(ここまで4時間・・・)
実行してみると、以下のようなエラーが出てきた・・・
一瞬とても嫌な汗をかいた
error while loading shared libraries: libhoge.so.1: cannot open shared object file: No such file or directory
ググってみると、
「どうやらダイナミックリンクを使っているプログラムを実行はしたものの、そのダイナミックリンクを見つけられませんでした」
というエラーらしい

解決方法は、
・ダイナミックリンクをデフォルトで探しに行く場所(/usr/libとか?)にシンボリックリンクを貼る
・ダイナミックリンカーが参照する環境変数 LD_LIBRARY_PATH に追記する
らしい。

とりあえず LD_LIBRARY_PATH にライブラリの場所を追記したら、無事動きました・・・・。
setenv LD_LIBRARY_PATH /home/hoge/foo/lib

ほんと新しい計算機で、プログラムを動かすだけなのに一苦労・・・・

■ 参考
LD_LIBRARY_PATH とは








今後C++を使う機会が多そうなので、早めに習得をしていく
とりあえずC言語との違うところをまず抑えて、次によく使うライブラリーを使ってみようと思う。






C++で画面出力するときは次のようにする
int a = 5;
cout << "a = " << a << "\n"
これで次のように表示されるはず
a = 5
これがもし科学の計算ならいろいろと注文したいことがある
1, 桁数をもっと多く表示してほしい(例えば15桁全部)
2, 12345e+5とか1.2345とか表示形式を指定したい
(3, 桁数を揃えたい)




とりあえず1と2について
double e = 1e-15;
cout << fixed << setprecision(16) << 1e-15 << "\n";
cout << scientific << setprecision(16) << 1e-15 << "\n";
これで、
0.0000000000000010
1.0000000000000001e-15
こんな感じに表示されるはず


ここでfixとかscientficは操作子と呼ばれるらしいです
他にもあるので、今後使う可能性があるものを簡単にまとめておきます

■基数
dec 整数の入出力を10進法で行う
hex 整数の入出力を16進法で行う
oct 整数の入出力を8進法で行う
setbase(n) 整数の入出力をn進法で行う


■ 浮動小数点
fixed 浮動小数点の数を23.45という風に表示する
scientific 浮動小数点の数を23.45e+6という風に表示する
setprecision(n) n桁の精度で表示する
showpoint 小数点を付けて表示する
noshowpoint 小数点を付けずに表示する
ただし setw, setprecision, setw他を使う場合は #include とプログラムの上の方に書いておく必要がある









現在新しい言語を模索中
Haskellを触ってみたが、いろいろと挫折

とりあえず昔買ったけどずっと放置していたc++の本を読んで、
ソースをそのまま書き写してみたりして勉強中

いろいろと便利すぎてやべぇw


今回は変数の上限値を表示させてみる
よくある例のヤツ

char   -128~127
signed char   -128~127
unsighed char   0~255

int   -2147483648~2147483647
unsigned int   0~4294967295

long   -9223372036854775808~9223372036854775807
unsighed long   0~18446744073709551615

double   2.22507e-308~1.79769e+308
unsigned double    0~1.79769e+308



include すべきヘッダーはclimitsとcfloat
climitsにはint型、cfloatにはfloat型の変数の上限値が書かれている

■ 参考リンク
(limits.h)、  (float.h)

最近、gnuplotでアニメを作って遊んでいる
探したら、ネットにいろいろなアニメを作るソースが落ちている

例えばものを投げたときの放物線のアニメとか
最初に解析式を求めるわけじゃなくて、微分方程式をその点について解いて点をどんどんうっていくような感じ

その中で「ソフトウェアタイマ」と呼ばれるものが使われていた

つまり
int i, k;
for(i=0; i<1e+6; i++) k++;
みたいなの
単純な演算を繰り返して、時間を無駄に使う。要するにここで足踏み

これは実行するマシンによって、時間が変わるしよろしくないらしい

そこで
unsigned int usec;
usleep( usec );
もしくは
unsigned int sec;
sleep( sec );
ってのを使うらしい

usecはミリ秒、secは秒だけ足踏み可能。

ブックマークを整理してたら見つけた

■リンク
計算物理のためのC/C++言語入門



ソートは色々と種類があって、それぞれ並べ替える方式の違いから実行ステップ数が変わってくる
その結果、かかる時間がデータの種類によってマチマチ
本当はそこんところをしっかり考えないといけない・・・・が、

なんでもいいからソートしたいという時もあるものです


そういうときは、

#include <gsl/gsl_sort.h>
#include <gsl/gsl_sort_vector.h>
をインクルードして、

gsl_sort(hoge, 1, 5);
とすればおーけー

hogeは配列の名前
1hoge[i]のiをいくつずつずらしていくか
5は配列の大きさ

とりあえずこれでソートはしてくれるようです
タイトルのわかりにくさ・・・・



要するに、2桁の数字を出力するときに、3桁で余ってるところを0埋めしたいと
i=10;
printf("%03d¥n", i);

で出力したら、
010
的な。


%03dの0がポイントです


■参考サイト
リーディングゼロの指定

危うく逆にするところだったのでメモ

atan2(y, x)で
949b4810
■参考
atan2






(2017/4/8 追記)

atan2の返り値は [-pi : pi]
acosの返り値は[0:pi]
asinの返り値は[-pi/2:pi/2]









とみながだいすけさんが日本語訳してくれているこちらのドキュメント によると、

17章に積分ルーチンが載っています。
(執筆時はversion 1.15が最新バー)

今までは、台形積分やらシンプソン則の積分しかC言語では書いて使ったことがなかったので
どうやって使うのかもさっぱり・・・・

とりあえずサンプルコードが一つあったので、そのまま書いて実行してみます。



が、エラー
gcc -c sample_integral.c
sample_integral.c: In function ‘main’:
sample_integral.c:37: warning: format ‘%ld’ expects type ‘long int’, but argument 2 has type ‘int’
gcc -o sample_integral sample_integral.o -lgsl -gslcblas
コンパイルは通って、実行ファイルも出来てますが・・・・

原因は構造体wのメンバsizeを呼び出して出力するときの型の問題。
sizeはsize_t型ですが、
同期曰く、size_t型は一般にはintと同じ型だけど、ここではlong int型として扱ってしまっているらしく

printf("%d¥n", w->size);

printf("%ld¥n", w->size);

と変更すればエラーは消えました。



しばらくはこれを使ってみようかと・・・
できれば中でどういうことをしているのかを知りたいですね。
double xx,yy,zz;
xx = 10.0/11.0;
yy = 10.0/11;
zz = 10/11;
printf("%lf %lf %lf\n", xx, yy, zz);



計算結果は
0.909091 0.909091 0.000000


doubleがどっかに入ってたら、結果もdoubleになるらしい

けど、怖いので、
yy = 10.0/(double)11;
とかにしたいですね
cal_eventrate.c: In function ‘cal_hogehoge’:
cal_hoge.c:12: error: nested functions are disabled, use -fnested-functions to re-enable
cal_hoge.c:92: error: expected declaration or statement at end of input


「ネストしちゃってるよ(ある関数の中でその関数を呼んでる、みたいな感じか?)」
「括弧が不足もしくは多い」

ってエラーで1時間ほど詰まった。

原因はヘッダーファイル

なんとなく苦労メモ
FFTWの公式サイトから落としてきたファイルを解凍して、中に入る

いつものcd fftw-3.3.1/
./configure
make


これでbenchが使えるようになる。



使い方と実行結果$bench 16                  #16はファイルサイズ
Problem: 16, setup: 18.00 us, time: 81.74 ns, ``mflops'': 3914.8


あとはファイルのサイズを変えてbenchに食わせるようなシェルスクリプトをちょいちょいと書いて実行。

Mac mini とMac Book airではこんな感じの結果になりました。

プロフィール
HN:coffee
職業:物理屋(自称)
趣味:映画鑑賞
自己紹介:
#include <stdio.h>
#include "MyProfile.h"

#define TWITTER coffee_pote

#define WISH_LIST
amazonのほしい物リスト
#ifdef RICH_FLAG
// ↑いつも支援いただきありがとうございます m(_ _)m
#endif


int main(void){

printf("\n");
printf("D論・・・? あぁそんな子もいましたね(執筆中)\n");
printf("\n");
printf("猿でもわかるgnuplot を執筆中(こっちの執筆は半年以上何も進んでいない・・・・)\n");
/* 最終更新 2017/07/19 */
return 0;

}
カウンター
ブログ内検索
ツイートするボタン
Flickr

Template "simple02" by Emile*Emilie
忍者ブログ [PR]