2009年7月11日土曜日

粒子法のプログラム最終回(粒子の出力)


最終回となる今回は、計算した粒子の出力について説明します。

このプログラムでは、計算した粒子をpov-rayのフォーマットで出力しています。

第1回に載せたムービーは、このプログラムでタイムステップごとに出力した.povファイルをレンダリングし、さらにffmpegで動画にしたものです。

粒子を出力するコードは以下のようになっています。

void output_particles( Particles* ps )
{
static int num = 0;
std::string file_name = \
(boost::format( "result%08d.pov" ) % num++ ).str();
std::cout << "processing " << file_name << " ..." << std::endl;
std::ofstream f( file_name.c_str() );
if ( ! f )
{
std::cerr << "cannot open " << file_name << std::endl;
exit(1);
}

f << "#include ¥"colors.inc¥"¥n"
<< "camera {¥n"
<< "  location <10,>¥n"
<< "  look_at <10,>¥n"
<< "}¥n"
<< "light_source { <0,> color White }¥n";
FOR_EACH_PARTICLES( ps, p )
{
f << "sphere {¥n"
<< "  <" << r(0) << ", " << r(1) << ", " << r(2) << ">, 0.5¥n"
<< "  texture {¥n"
<< "    pigment { color Yellow }¥n"
<< "  }¥n"
<< "}¥n";
}
f << std::endl;
}
タイムステップごとに.povファイルを出力しています。

さらに、この投稿に載せたムービーでは、計算した粒子から求めた陰曲面をマーチングキューブでポリゴン化し、反射や屈折も考慮してレンダリングしています。

これでひととおりの説明が終わりました。このあとの発展としては、3次元にしたり、物体との相互作用を入れたり、GPUで処理したりといったことを考えています。

粒子法のプログラム
粒子法のプログラム第1回(概要)
粒子法のプログラム第2回(プログラムの大枠)
粒子法のプログラム第3回(データ構造)
粒子法のプログラム第4回(密度と圧力の計算)
粒子法のプログラム第5回(力の計算)
粒子法のプログラム第6回(境界条件と粒子位置の更新)
粒子法のプログラム最終回(粒子の出力)

その他の解説エントリ
SPHによる巻き波のシミュレーション1
SPHによる巻き波のシミュレーション2
SPHによる巻き波のシミュレーション3
このあとやりたいこと

固液連成シミュレーションに関するエントリ
粒子法による固液連成シミュレーション

流体シミュレーションに関するエントリ
【粒子法】粒子を流体としてレンダリング
3次元の粒子法シミュレーション
粒子法のシーンを2倍のサイズにしてみたが…
粒子法のシーンを2倍のサイズにしてみた
Haskell、OCamlでSPH法
カメラ位置を変えて流体をレンダリング
Bunny-shaped fluid simulation

剛体シミュレーションに関するエントリ
粒子ベース剛体シミュレーション(プレビュー)
粒子ベース多体衝突シミュレーション
引き続き、粒子ベース剛体シミュレーション
Falling Rigid Bunnies

動画
シミュレーションの結果をレンダリングして作った動画です。流体シミュレーションや剛体シミュレーションの動画を見ることができます。

動画の一覧


--

2009年6月20日土曜日

粒子法のプログラム第6回(境界条件と粒子位置の更新)


今回は、粒子の境界条件と粒子位置の更新について説明します。

境界条件と粒子位置の更新のコードは以下のようになります。

void calculate_position( Particles* ps )
{
vec accel(3), g(3), norm(3);
double speed, diff, adj;

g(0) = 0.0; g(1) = -9.8; g(2) = 0.0;
FOR_EACH_PARTICLES( ps, p )
{
accel = p->f * SPH_PMASS;

speed = norm2( accel );
if ( speed > SPH_LIMIT*SPH_LIMIT ) {
accel *= SPH_LIMIT / sqrt(speed);
}

// Z-axis walls
diff = 2.0 * RADIUS - ( p->r(2) - MIN(2) ) * SPH_SIMSCALE;
if ( diff > EPSILON )
{
norm(0) = 0.0; norm(1) = 0.0; norm(2) = 1.0;
adj = SPH_EXTSTIFF * diff - SPH_EXTDAMP * inner_prod( norm, p->v );
accel += adj * norm;
}
diff = 2.0 * RADIUS - ( MAX(2) - p->r(2) ) * SPH_SIMSCALE;
if ( diff > EPSILON )
{
norm(0) = 0.0; norm(1) = 0.0; norm(2) = -1.0;
adj = SPH_EXTSTIFF * diff - SPH_EXTDAMP * inner_prod( norm, p->v );
accel += adj * norm;
}

// X-axis walls
diff = 2.0 * RADIUS - ( p->r(0) - MIN(0) ) * SPH_SIMSCALE;
if ( diff > EPSILON )
{
norm(0) = 1.0; norm(1) = 0.0; norm(2) = 0.0;
adj = SPH_EXTSTIFF * diff - SPH_EXTDAMP * inner_prod( norm, p->v );
accel += adj * norm;
}
diff = 2.0 * RADIUS - ( MAX(0) - p->r(0) ) * SPH_SIMSCALE;
if ( diff > EPSILON )
{
norm(0) = -1.0; norm(1) = 0.0; norm(2) = 0.0;
adj = SPH_EXTSTIFF * diff - SPH_EXTDAMP * inner_prod( norm, p->v );
accel += adj * norm;
}

// Y-axis walls
diff = 2.0 * RADIUS - ( p->r(1) - MIN(1) ) * SPH_SIMSCALE;
if ( diff > EPSILON )
{
norm(0) = 0.0; norm(1) = 1.0; norm(2) = 0.0;
adj = SPH_EXTSTIFF * diff - SPH_EXTDAMP * inner_prod( norm, p->v );
accel += adj * norm;
}
diff = 2.0 * RADIUS - ( MAX(1) - p->r(1) ) * SPH_SIMSCALE;
if ( diff > EPSILON )
{
norm(0) = 0.0; norm(1) = -1.0; norm(2) = 0.0;
adj = SPH_EXTSTIFF * diff - SPH_EXTDAMP * inner_prod( norm, p->v );
accel += adj * norm;
}

accel += g;
p->v += accel * DT;
p->r += p->v * DT / SPH_SIMSCALE;

// stick on x-y plane
p->r(2) = 0.0;

}
}
それぞれの粒子について、粒子にかかる力から加速度を計算したあと、各壁面について境界条件を適用し、重力加速度を加えます。そうやって求まった加速度から、粒子の速度を求め、粒子の位置を更新します。

境界条件には、ペナルティ法を用いています。ペナルティ法とは、粒子が壁面へめりこんだ場合に、そのめり込みの程度によって壁面との衝突による反発力を与える方法です。

このようにして、境界条件を処理し、粒子の位置を更新しています。

次回は計算結果の出力です。最終回です。

粒子法のプログラム
粒子法のプログラム第1回(概要)
粒子法のプログラム第2回(プログラムの大枠)
粒子法のプログラム第3回(データ構造)
粒子法のプログラム第4回(密度と圧力の計算)
粒子法のプログラム第5回(力の計算)
粒子法のプログラム第6回(境界条件と粒子位置の更新)
粒子法のプログラム最終回(粒子の出力)

その他の解説エントリ
SPHによる巻き波のシミュレーション1
SPHによる巻き波のシミュレーション2
SPHによる巻き波のシミュレーション3
このあとやりたいこと

固液連成シミュレーションに関するエントリ
粒子法による固液連成シミュレーション

流体シミュレーションに関するエントリ
【粒子法】粒子を流体としてレンダリング
3次元の粒子法シミュレーション
粒子法のシーンを2倍のサイズにしてみたが…
粒子法のシーンを2倍のサイズにしてみた
Haskell、OCamlでSPH法
カメラ位置を変えて流体をレンダリング
Bunny-shaped fluid simulation

剛体シミュレーションに関するエントリ
粒子ベース剛体シミュレーション(プレビュー)
粒子ベース多体衝突シミュレーション
引き続き、粒子ベース剛体シミュレーション
Falling Rigid Bunnies

動画
シミュレーションの結果をレンダリングして作った動画です。流体シミュレーションや剛体シミュレーションの動画を見ることができます。

動画の一覧


--

2009年6月14日日曜日

Vivienne Westwoodの長財布


財布って、どんなのを使ってます?

私は父親からアメリカ土産にもらった革の三つ折り財布を10年以上使っているのですが、最近、長財布もいいなぁと思っていました。

それで今日、表参道に行って、Vivienne Westwoodの長財布を買ってきました。

買ったのは、写真にあるADVANのネイビー。表面に味のあるアドバン加工を施したレザーをベースに、全面に大胆なORB型の型押しを施してあります。ORBはVivienneのいろんなアイテムに施されていて、そこには「伝統をもって未来を創る」という彼女のコンセプトが込められています。

あと、外側のネイビーと内側のイエローのコンビネーションがとってもクール!この絶妙なコンビネーションが素晴らしい!!

ブティックで買い物をしたのが初めてだったので、ドキドキしたー

大事に長く使おう!

2009年5月30日土曜日

粒子法のプログラム第5回(力の計算)


今回は、粒子法のプログラムでどのように力を計算しているかを説明します。

力の計算は次のコードになります。


void calculate_force( Particles* ps )
{
double pterm, vterm, r, c;
vec dr(3), force(3), fcurr(3);

FOR_EACH_PARTICLES( ps, p )
{
force(0) = 0.0; force(1) = 0.0; force(2) = 0.0;
FOR_EACH_PARTICLES( ps, pj )
{
if ( p == pj ) continue;
dr = ( p->r - pj->r ) * SPH_SIMSCALE;
r = norm( dr );
if ( H > r )
{
c = H - r;
pterm = -0.5 * c * SpikyKern * ( p->p + pj->p ) / r;
vterm = LapKern * SPH_VISC;
fcurr = pterm * dr + vterm * ( pj->v - p->v );
fcurr *= c * p->rho * pj->rho;
force += fcurr;
}
}
p->f = force;
}
}
それぞれの粒子について、ナビエ=ストークス方程式の圧力項、粘性項にあたる力を計算しています。

影響半径内にある自分以外の粒子に対して、距離に応じて重みを足し合わせています。SpikyKernとLapKernがカーネルの係数です。

力の計算も、密度と圧力の計算と同様に、近傍粒子探索を工夫したり対称性を利用することで、高速化できます。

次回は、どのように粒子の位置を更新するかという説明になります。

粒子法のプログラム
粒子法のプログラム第1回(概要)
粒子法のプログラム第2回(プログラムの大枠)
粒子法のプログラム第3回(データ構造)
粒子法のプログラム第4回(密度と圧力の計算)
粒子法のプログラム第5回(力の計算)
粒子法のプログラム第6回(境界条件と粒子位置の更新)
粒子法のプログラム最終回(粒子の出力)

その他の解説エントリ
SPHによる巻き波のシミュレーション1
SPHによる巻き波のシミュレーション2
SPHによる巻き波のシミュレーション3
このあとやりたいこと

固液連成シミュレーションに関するエントリ
粒子法による固液連成シミュレーション

流体シミュレーションに関するエントリ
【粒子法】粒子を流体としてレンダリング
3次元の粒子法シミュレーション
粒子法のシーンを2倍のサイズにしてみたが…
粒子法のシーンを2倍のサイズにしてみた
Haskell、OCamlでSPH法
カメラ位置を変えて流体をレンダリング
Bunny-shaped fluid simulation

剛体シミュレーションに関するエントリ
粒子ベース剛体シミュレーション(プレビュー)
粒子ベース多体衝突シミュレーション
引き続き、粒子ベース剛体シミュレーション
Falling Rigid Bunnies

動画
シミュレーションの結果をレンダリングして作った動画です。流体シミュレーションや剛体シミュレーションの動画を見ることができます。

動画の一覧


--

2009年5月24日日曜日

目覚めを良くするための最良の方法


最近、枕を変えてみた。

これまではふつうのやわらかい枕。それを、四角くて固めのい草枕に変えた。

すると、目覚めがすごくすっきりするようになった。

枕でここまで変わるのかー

びっくり。これから蒸し暑くなるので、その点でもグッド!


2009年5月23日土曜日

粒子法のプログラム第4回(密度と圧力の計算)


今回は、粒子法のプログラムでどのように密度と圧力を計算しているかを説明します。

粒子の密度と圧力の計算は次のコードになります。
void calculate_density_and_pressure( Particles* ps )
{
double r2, c, sum, H2;
vec dr;

H2 = H*H;

FOR_EACH_PARTICLES( ps, p )
{
sum = 0.0;
FOR_EACH_PARTICLES( ps, pj )
{
if ( p == pj ) continue;
dr = (p->r - pj->r) * SPH_SIMSCALE;
r2 = norm2( dr );
if ( H2 > r2 )
{
c = H2 - r2;
sum += c * c * c;
}
}
p->rho = sum * SPH_PMASS * Poly6Kern;
p->p   = ( p->rho - SPH_RESTDENSITY ) * SPH_INTSTIFF;
p->rho = 1.0 / p->rho;  // take inverse for later calculation
}
}

それぞれの粒子について、影響半径内にある自分以外の粒子に対して、距離に応じて重みを足し合わせています。そのあと、足し合わせた重みに、粒子質量についての係数とカーネル係数をかけています。密度が計算された後、圧力を密度から計算します。

上のやり方だと、計算量のオーダーは、各粒子について他の粒子それぞれとの計算になるので、粒子の数をnとするとO(n^2)となります。しかし、影響半径外の粒子との計算は不要なので、他粒子すべてと計算する必要は本来ありません。なので、近傍粒子をあらかじめマップに入れておき、それを参照するようにすれば、計算量をO(n)にすることができます。

また、粒子pと粒子pjについて対称なので、その対称性を利用すれば、計算量をさらに半分にすることが可能です。

密度と圧力は、このように計算されています。いったん素直な実装になっていますが、上の2つの点で高速化することができます。

次回は、粒子にかかる力をどのように計算しているかという説明です。

粒子法のプログラム
粒子法のプログラム第1回(概要)
粒子法のプログラム第2回(プログラムの大枠)
粒子法のプログラム第3回(データ構造)
粒子法のプログラム第4回(密度と圧力の計算)
粒子法のプログラム第5回(力の計算)
粒子法のプログラム第6回(境界条件と粒子位置の更新)
粒子法のプログラム最終回(粒子の出力)

その他の解説エントリ
SPHによる巻き波のシミュレーション1
SPHによる巻き波のシミュレーション2
SPHによる巻き波のシミュレーション3
このあとやりたいこと

固液連成シミュレーションに関するエントリ
粒子法による固液連成シミュレーション

流体シミュレーションに関するエントリ
【粒子法】粒子を流体としてレンダリング
3次元の粒子法シミュレーション
粒子法のシーンを2倍のサイズにしてみたが…
粒子法のシーンを2倍のサイズにしてみた
Haskell、OCamlでSPH法
カメラ位置を変えて流体をレンダリング
Bunny-shaped fluid simulation

剛体シミュレーションに関するエントリ
粒子ベース剛体シミュレーション(プレビュー)
粒子ベース多体衝突シミュレーション
引き続き、粒子ベース剛体シミュレーション
Falling Rigid Bunnies

動画
シミュレーションの結果をレンダリングして作った動画です。流体シミュレーションや剛体シミュレーションの動画を見ることができます。

動画の一覧


--

2009年5月19日火曜日

【粒子法】粒子を流体としてレンダリング


粒子法では流体を粒子として表現します。

そのため、そのまま描画すると、この投稿のように球のあつまりとして描画することになります。

これを、流体として見えるようにレンダリングしてみました。

マーチングキューブで陰曲面をポリゴン化し、レイトレーシングでレンダリングしています。

フォトンマップも使っているのでコースティクスも表現されています。

屈折率を高めにしたのでキラキラ。




■固液連成シミュレーションに関するエントリ
粒子法による固液連成シミュレーション

■関連する記事
粒子法のプログラム第1回(概要)
粒子法のプログラム第2回(プログラムの大枠)
粒子法のプログラム第3回(データ構造)
粒子法のプログラム第4回(密度と圧力の計算)
粒子法のプログラム第5回(力の計算)
粒子法のプログラム第6回(境界条件と粒子位置の更新)
粒子法のプログラム最終回(粒子の出力)
【粒子法】粒子を流体としてレンダリング
3次元の粒子法シミュレーション
粒子法のシーンを2倍のサイズにしてみたが…
粒子法のシーンを2倍のサイズにしてみた

粒子法(SPH)のプログラム
粒子法(SPH)のプログラムを解説したシリーズです。ソースコードも公開しています。

粒子法(SPH)のプログラム一覧



動画
シミュレーションの結果をレンダリングして作った動画です。流体シミュレーションや剛体シミュレーションの動画を見ることができます。

動画の一覧


--

2009年4月19日日曜日

粒子法のプログラム第3回(データ構造)


今回は、粒子法のプログラムのデータ構造を説明します。

使用したデータ構造は、Particle構造体とParticlesベクタの2つです。

Particles構造体は以下です。
typedef boost::numeric::ublas::vector<double> vec;
typedef struct
{
vec    r;
vec    v;
double rho;
double p;
vec    f;
} Particle;
boostライブラリのublasを使い、3次元ベクトルを表現しています。
rは粒子座標、vは粒子速度、rhoは粒子地点での密度、pは粒子地点での圧力、fは粒子が受ける力です。

Particlesベクタは以下です。
typedef std::vector<particle> Particles;
STLのvectorコンテナでParticle構造体のベクタを表現しています。

使用しているデータ構造は、この2つのみです。このように、粒子法はデータ構造が非常に簡単なのが特徴です。

次回からは、粒子法のアルゴリズムの説明になります。


粒子法のプログラム
粒子法のプログラム第1回(概要)
粒子法のプログラム第2回(プログラムの大枠)
粒子法のプログラム第3回(データ構造)
粒子法のプログラム第4回(密度と圧力の計算)
粒子法のプログラム第5回(力の計算)
粒子法のプログラム第6回(境界条件と粒子位置の更新)
粒子法のプログラム最終回(粒子の出力)

その他の解説エントリ
SPHによる巻き波のシミュレーション1
SPHによる巻き波のシミュレーション2
SPHによる巻き波のシミュレーション3
このあとやりたいこと

固液連成シミュレーションに関するエントリ
粒子法による固液連成シミュレーション

流体シミュレーションに関するエントリ
【粒子法】粒子を流体としてレンダリング
3次元の粒子法シミュレーション
粒子法のシーンを2倍のサイズにしてみたが…
粒子法のシーンを2倍のサイズにしてみた
Haskell、OCamlでSPH法
カメラ位置を変えて流体をレンダリング
Bunny-shaped fluid simulation

剛体シミュレーションに関するエントリ
粒子ベース剛体シミュレーション(プレビュー)
粒子ベース多体衝突シミュレーション
引き続き、粒子ベース剛体シミュレーション
Falling Rigid Bunnies

動画
シミュレーションの結果をレンダリングして作った動画です。流体シミュレーションや剛体シミュレーションの動画を見ることができます。

動画の一覧


--

2009年3月29日日曜日

さくらの根元に・・・


今日は日本橋にいってきました。

通りの両側にさくらが並んでいて、3分咲きくらいになっていたのですが、その根元に、今の時期絶対ないはずのものが。

なんと、雪があったんです。

もちろん降った雪が積もっていたわけではなく、たぶん、桜の開花を遅らせるためなんだろうなと思います。

聞くこと


最近、「聞く」ことに注意している。

これまでも、もちろん人の話を聞いて相手が何を言おうとしているのか推測して返すようにはしていたのだけれど、どうも上手くいっていなかった。

なんで上手くいかないんだろうと思っていると、その原因は推測がずれていることだということに気づいた。

それ以来、人の話を聞くときには、自分の推測があっているか、「それってこういうこと?」と確認するようにしている。

そうすると、びっくりするくらい会話がスムーズになった!
すばらしい。

リボンが流行り


今日、新宿に行ったのでついでにルミネを見てきた。

ルミネでお店を見てみると、リボン付きの服をよくみかけた。

つい最近日経MJでリボンが流行っているよという読んだことがあったので、「確かに流行ってる!」と納得。

2009年3月28日土曜日

フィール・ザ・カライ・センセーション!


シンガポールのマクドナルドのテレビCM。
http://www.mcdonalds.com.sg/tvcs/wasabi_tvc.mov

台詞が日本語なんだけど、そのあとのメッセージがまるっきりカタカナ英語で、めっちゃウケる。ワサビ=日本=カタカナ英語。ちょっと!バカにされてるよ!

あと「カライ・センセーション」というのも謎だけど、確かに、ワサビの辛さって、hotでもspicyとも違って、"karai"んだなぁと。

2009年3月17日火曜日

チンチョウゲ


オフィスの前にチンチョウゲの花が咲いています。
チンチョウゲが子供のころから好きです。

子供のころ住んでいたマンションの下に植わっており、2月3月になるとあの特徴的な香りといっしょに咲いて、

けど、1年に1度しか咲かないので、
「ちんちょうげ」という名前がなかなか覚えられず、

何か不思議なにおい、として記憶に残っています。

2009年3月15日日曜日

はるみ


はるみというみかんが売ってたので、食べてみた。

薄皮をむかなくても食べられて、ふつうのみかんより大きいから食べごたえもあって、ちゃんと甘いので、すごくおいしい!

2009年3月9日月曜日

3月8日 東京都美術館〜上野公園〜大崎



昨日は昼から東京都美術館へ行って、上野公園を散歩したあと、大崎のスタバで新卒採用のブログを書いた。

東京都美術館では、「アート&クラフツ」展をみた。展示のテーマは、19席後半にイギリスで起こった「アーツ&クラフツ」というデザイン運動。ウィリアム・モリスが主導した「アーツ&クラフツ」というのは、手仕事の良さを見直し、自然や伝統から美を再発見し、シンプルなライフスタイルを提案する、という運動で、その背景には、産業革命以降の工業化によって粗悪な大量生産品が増えたことがあったそう。

東京都美術館から出たあと、上野公園を散歩した。あいにくの曇り空だったが、ふと見上げると、さくらが咲いていて、鳥が花の周りを飛び回っていた。そうか、3月か、と思った。空のくもり色と、逆光の枝の黒、鳥の影、そのなかに、咲き始めのさくらのほのかなピンク色が印象的で、ノートにそれを描いた。

そのあとは、大崎へいって、新卒採用のブログを書いた。就活中の学生に、うちの会社へ興味を持ってもらうためのブログ。ターゲットとシナリオを考えた結果、どんな人間がこの会社にいるのかを知ってもらうのがベストだと考え、徹底的に自分のこれまで歩んできた道を話すことにした。ちょっとした自叙伝?

充実した1日だった。感謝!

2009年3月6日金曜日

白いたいやき


白いたいやきが流行っているらしい。

白い理由は、タピオカでんぷんを使っているから。

大福みたいなもちもちな食感だって。

そういえば、白いイチゴというのも、こないだみたなー。

白いイチゴは、突然変異でたまたまできたのを、たくさん作れるようにしたんだって。

赤いイチゴとならべて、紅白で売ってたりするのがいい感じ。

2009年2月28日土曜日

粒子法のプログラム第2回(プログラムの大枠)


プログラムの大枠が上の図です。

まず初期化をします。

そのあと密度と圧力を計算し、力を計算し、粒子の位置を更新します。

粒子の位置を更新したら、そのステップでの結果をファイルに出力し、また密度と圧力の計算へ戻ります。

これを、指定ステップ数だけ繰り返します。


粒子法のプログラム
粒子法のプログラム第1回(概要)
粒子法のプログラム第2回(プログラムの大枠)
粒子法のプログラム第3回(データ構造)
粒子法のプログラム第4回(密度と圧力の計算)
粒子法のプログラム第5回(力の計算)
粒子法のプログラム第6回(境界条件と粒子位置の更新)
粒子法のプログラム最終回(粒子の出力)

その他の解説エントリ
SPHによる巻き波のシミュレーション1
SPHによる巻き波のシミュレーション2
SPHによる巻き波のシミュレーション3
このあとやりたいこと

固液連成シミュレーションに関するエントリ
粒子法による固液連成シミュレーション

流体シミュレーションに関するエントリ
【粒子法】粒子を流体としてレンダリング
3次元の粒子法シミュレーション
粒子法のシーンを2倍のサイズにしてみたが…
粒子法のシーンを2倍のサイズにしてみた
Haskell、OCamlでSPH法
カメラ位置を変えて流体をレンダリング
Bunny-shaped fluid simulation

剛体シミュレーションに関するエントリ
粒子ベース剛体シミュレーション(プレビュー)
粒子ベース多体衝突シミュレーション
引き続き、粒子ベース剛体シミュレーション
Falling Rigid Bunnies

動画
シミュレーションの結果をレンダリングして作った動画です。流体シミュレーションや剛体シミュレーションの動画を見ることができます。

動画の一覧


--

2009年2月23日月曜日

粒子法のプログラム第1回(概要)


粒子法のプログラムを実装しようとしたのですが、論文や本だと数式レベルでしか説明がなく、それをプログラムにして動くようにするまで苦労しました。

数式をプログラムに落とすこと自体は簡単なのですが、各定数の値が分からなかったりして、粒子が上手く動いてくれるまで大変でした。

なので、同じところでつまづく人も多いのではと思い、自分が書いた粒子法のプログラムを紹介しようと思います。下の動画を作ったときのプログラムです。SPH法です。



こちらがソースコードです。

※2009.11.4追記 こちらに整理したC++版のほか、Haskell版、OCaml版のコードがあります。

まずはシンプルにと考え、近傍粒子の探索は工夫せず、n(O^2)ですがわかりやすいアルゴリズムを採用しています。まず動くコードを書き、そこから高速化していくアプローチです。

また、計算量を減らした方が試しやすいのと、粒子法は3次元への拡張が簡単なことから、まずは2次元で実装しています。

プログラムの中身を別の記事で少しずつ説明していこうと思います。

粒子法のプログラム
粒子法のプログラム第1回(概要)
粒子法のプログラム第2回(プログラムの大枠)
粒子法のプログラム第3回(データ構造)
粒子法のプログラム第4回(密度と圧力の計算)
粒子法のプログラム第5回(力の計算)
粒子法のプログラム第6回(境界条件と粒子位置の更新)
粒子法のプログラム最終回(粒子の出力)

その他の解説エントリ
SPHによる巻き波のシミュレーション1
SPHによる巻き波のシミュレーション2
SPHによる巻き波のシミュレーション3
このあとやりたいこと

固液連成シミュレーションに関するエントリ
粒子法による固液連成シミュレーション

流体シミュレーションに関するエントリ
【粒子法】粒子を流体としてレンダリング
3次元の粒子法シミュレーション
粒子法のシーンを2倍のサイズにしてみたが…
粒子法のシーンを2倍のサイズにしてみた
Haskell、OCamlでSPH法
カメラ位置を変えて流体をレンダリング
Bunny-shaped fluid simulation

剛体シミュレーションに関するエントリ
粒子ベース剛体シミュレーション(プレビュー)
粒子ベース多体衝突シミュレーション
引き続き、粒子ベース剛体シミュレーション
Falling Rigid Bunnies

動画
シミュレーションの結果をレンダリングして作った動画です。流体シミュレーションや剛体シミュレーションの動画を見ることができます。

動画の一覧


--

2009年2月18日水曜日

AdWordsの検索クエリレポート

リスティングキーワードのプラニングを行うときには、部分一致で入札したキーワードが実際にどのようなキーワードで検索されているのかを調べると役に立ちます。

AdWordsで実際にどのようなキーワードで検索されているのかというレポートが出せるのかどうか調べてみると、出せることがわかりました。
レポートセンター >
新しいレポート >
レポートタイプを「検索クエリの掲載結果」
とすれば、検索クエリのレポートが出力されます。

■検索クエリレポートの活用方法
検索クエリレポートには、以下のような活用方法があります。
  • 入札キーワードが効果的にユーザーの検索クエリや潜在顧客をとらえているか確認可能
  • どの検索クエリが効果的なのかを把握し、入札キーワードとして登録されているかを確認可能
  • 入札キーワードに適したマッチタイプを選択するのに利用可能
  • 新規登録が必要なキーワードがないか、もしくは削除したり除外キーワード登録が必要なものがないか確認可能

■検索クエリレポートの注意点
検索クエリレポートをみる際には、以下のような点に注意しよう。

  • 広告のクリックにつながったクエリのみが表示される(インプレッションのみのクエリは表示されない)
  • 検索ネットワークのデータのみが含まれる(モバイル検索のデータは含まれない)
  • 確認可能なデータは2007年5月2日以降のもの
  • 検索クエリに対応する入札キーワードは表示されない
  • 検索回数の少ないものや、個人情報に関わる場合、[検索クエリ キーワードマッチ]には「X個の固有なクエリ」と表示され、Xは実際に含まれているユニーククエリ数に置き換えられる

2008年11月30日日曜日

ElasticFoxが使えない…

Amazon EC2がなぜか使えない。ElasticFox経由で使おうとしているんだけど、
EC2 responded with an error for DescribeInstances
AuthFailure: AWS was not able to validate the provided access credentials
と言われる。

Access KeyもSecret Access Keyもちゃんといれているつもりなんだが…

2008年11月29日土曜日

ログデータから、特定項目についてユニークなものを調べる

スペース区切りのログデータの3〜5項目めについて、
ユニークなものを調べる。
check unique rows on 3 - 5 columns of space-separated log data.
awk '{ print $3, $4, $5; }' file | sort -u
awk便利だ。
awk is very useful!

2008年6月18日水曜日

【Hadoop】JAVA_HOMEの設定


Hadoopを動かすには、conf/hadoop-env.shの中にJAVA_HOMEという値を設定する必要があります。JAVA_HOMEには、Javaをインストールしたディレクトリのパスを指定します。以下に、一般的なJavaホームディレクトリのパスを記載します。
※1.5.0の部分は、使用するJavaのバージョンによって変えてください。

Windowsでの例

C:\jdk1.5.0

MacOSXでの例

/System/Library/Frameworks/JavaVM.framework/Versions/1.5.0/Home/

Linuxでの例

/usr/java/j2sdk1.5.0


2008年6月16日月曜日

コロー展とスナップスケッチ


西洋美術館のコロー展に行ってきた。

美術館で絵を見るとき、いつもなら他の人と同じように1つ1つ絵を見ていくのだが、今日はちょっと違ったことにチャレンジしてみた。

絵を見ながらその場でメモ帳にエンピツでスナップスケッチを描いてみたのである。

実際に自分でスケッチしてみると、驚いた。

ただ見るだけでは気づかなかったものが次々と見えてきたり、その絵をあとからでもはっきりと頭に思い描くことが出来るのである。

模写して画家の描いたものをなぞることによって、はじめてその絵を見たことになるのかもしれない。

2008年6月5日木曜日

グーグルとアマゾンのクラウドの違い


@ITより

GAEはWebアプリケーション作成に特化したプラットフォームを提供
AWSはより汎用的なコンピューティングプラットフォームを提供

GAEで利用できる開発言語は今のところPythonだけ。
SDKには、Djangoが含まれている。
その他、グーグルが提供するAPI群も利用できる
- Datastore API
- Images API
- Mail API
- Memcache API
- URL Fetch API
- Users API
GAE上で稼働するWebサービスはグーグルのインフラ上で大規模利用
に耐えるほどスケールする、というのがグーグルの売り文句。

一方、Amazon EC2はよりPCに近い仮想環境を提供する。
開発言語どころか、どのOSを使うかも自由。

AWSにはAmazon EC2のほかにも、仮想ストレージ(Amazon S3)、
仮想データベース(SimpleDB、AmazonSQS)などがある。

AWSに出来てGWSに出来ないのは、次のような利用方法。
つまり、Hadoopを使った分散処理。
例として挙っているのは、ニューヨークタイムズの担当者が
過去の記事をデジタルアーカイブとして公開するために合計
4TBにもなるTIFF画像をPDF化するためにHadoopを使ったという例。
Amazon EC2に100インスタンスを生成。24時間以下で1100万枚の
画像から1.5TB分のPDFを作ることに成功したという。

2008年5月14日水曜日

iPod Touch 1.1.4でmouseDownイベントを取得するには


iPod Touch 1.1.4でmouseDownイベントを取得する。

・__GSEventの宣言が必要。

・__GSEvent.i4がx座標、__GSEvent.i5がy座標。
 いずれもfloat型。

・UIResponderのサブクラスに定義する。

・短辺がx軸、長辺がy軸。
typedef struct __GSEvent {
long i0;
float i1;
long eventType;
float i3;
float i4;
float i5;
} __GSEvent;

- (void)mouseDown:(struct __GSEvent *)fp8
{
NSLog(@"(%f,%f)", (*fp8).i4, (*fp8).i5);
}

2008年5月11日日曜日

ブラウザとサーバで双方向通信


古橋貞之の日記より
Ajaxはクライアントの都合でサーバーに通信を仕掛けるpull型の通信ができ、Cometはサーバーが好きなタイミングでクライアントへデータを送りつけるpush型の通信ができるわけですが、新たに双方向の通信ができる技術を開発しました。...
とのこと。

ブラウザとサーバで双方向通信ってどうやるんだ?
と思っていると、
...仕組みはどうなっているかと言うと、通信部分にFlashを使っています。Flashとサーバーとの通信には、ActionScriptのSocketを使っています。...
だって。
Flashってソケット通信できたのか!
なるほどなぁ

2008年5月7日水曜日

フォルダ内のファイルの拡張子を一括削除


OSXのAutomatorを使えば、
フォルダ内のすべてのファイルに一括で拡張子を付けることができる。
しかし、逆に一括で削除することが出来なくて困った。

そこで、シェルスクリプトで拡張子を一括削除する1行を書いてみた。
for i in `ls`; do mv $i ${i:0:${#i}-4};done
${variable:offset:length}で変数variableをスライスできる。

さらに、${#i}の部分でファイル名の文字数を取得し、
そこから"."+拡張子3文字分を引いている。

カレントディレクトリ内のすべてのファイルから、
3文字の拡張子を取り除く。

例:hoge.txt -> hoge

2008年5月4日日曜日

Cocoaで"~/"を展開するには


Cocoaで、"~/Documents"などを展開して

ホームディレクトリの絶対パスを得るには、
- [NSString stringByStandardizingPath]
を使う。たとえば、
NSString* d;
d = [@"~/Documents" stringByStandardizingPath];
NSLog(d);
のようにする。

2008年5月2日金曜日

【Haskell】Networkモジュール


Networkモジュールを使ってみた。
import Network
import System.IO

main = withSocketsDo (
do s <- listenOn $ PortNumber 8080
(h, n, p) <- accept s
print ( h, n, p )
getLineAndPrint s h
sClose s
)

getLineAndPrint s h = do l <- hGetLine h
print l
exitTest l
where exitTest "exit\r" = return ()
exitTest _ = getLineAndPrint s h

PortNumber 8080 をlistenOnして、acceptする。

hGetLineソケットから読み込み、printする。

もし入力されたものが"exit\r"なら終了。

そうでなければ、引き続き読み込む。

【Python】OSからのシグナルをハンドルする


signalモジュールを使うと、OSからのシグナルをハンドルできます。

たとえば、プログラムを停止するときの"Control + C"。

これもシグナルの1つです。

こんなコードを書くことで、"Control + C"をハンドルできます:
import signal

def handler( signum, frame ):
print 'Signal handler called with signal', signum

signal.signal( signal.SIGINT, handler )

while True:
pass
このコードを実行すると・・・"Control + C"を押しても停止しません!

決まった停止手順を踏ませたい場合など、

"Control + C"を無効にしたいときに使えます。

2008年4月25日金曜日

【SQL】FROM句でサブクエリを使う


WHERE句、SELECT句でのサブクエリは結果として値を利用する。
 
FROM句でサブクエリを利用すると、結果を1つのテーブルとして扱うことができる。FROM句でのサブクエリの利用は非常に応用性が高い。

FROM句にサブクエリを使うことで、集計関数を組み合わせた計算が出来る。たとえば、COUNTのAVGを計算することができる。

サブクエリには";"をつけない。

2008年4月23日水曜日

The Haskell School of Expression: Lerning Functional Programming Through Multimedia


こんな本を買ってみました。

Haskellの入門書なのですが、他の本とはひと味違います。

Haskellの本は、数学的な文脈で語られることが多いのですが、この本では、グラフィクス、アニメーション、コンピュータミュージックといったマルチメディアアプリケーションを題材として使っています。

レベルとしては、初級者と上級者の両方に向けて書かれています。
簡単なイントロダクションからはじまり、もっと進んだ話題へと進みます。Haskellがまったく初めての方にはややテンポが速いかもしれません。

けどさ・・・グラフィクスとかアニメーションって、純粋関数言語を使う意味があんまりないんじゃ?だって、副作用のかたまりやん;)

2008年4月20日日曜日

【Haskell】Haskellの低レイヤグラフィクスライブラリまとめ


Haskellの低レイヤグラフィクスライブラリをまとめました。
各プラットフォームのグラフィクスシステムへのバインディングです。
GLUT
OpenGL GLUT Libraryへのバインディング。
OpenGL GLUT Libraryは、プラットフォーム固有の処理を吸収してOpenGLを使いやすくするためのライブラリです。
GLFW
GLFWへのバインディングです。
GLFW(OpenGL Framework)は、GLUTと同様、プラットフォーム固有の処理を吸収してOpenGLアプリケーションを使いやすくするためのライブラリです。
GLUTよりも新しいです。
TclHaskell
TclHaskellは、Tck/tkへのバインディング。
Tck/Tkは、Unix、Windows、Macintosh上で動作するGUIスクリプティング環境です。
Win32
Win32ライブラリは、Win32APIのうちグラフィクスに関係する部分へのバインディングです。
Win32APIは、Windowsの標準APIです。
X11
X11ライブラリは、X11へのバインディングです。
X11は、Unix系OSのウィンドウシステムであるX Window Systemのバージョン11です。


2008年4月11日金曜日

【MySQL】カーソルの使い方


カーソルの宣言
DECLARE [カーソル名] CURSOR FOR [SELECT文];
カーソルを開く
OPEN [カーソル名]
1行ごとにデータを取り出す
FETCH NEXT FROM [カーソル名] INTO [変数リスト];
WHILE (@@FETCH_STATUS = 0)
BEGIN
FETCH NEXT FROM [カーソル名] INTO [変数リスト];
END;
カーソルを閉じ、リソースを解放する
CLOSE [カーソル名]
DEALLOCATE [カーソル名]

2008年4月6日日曜日

複数のファイルを1つ1つコマンドに読ませるには


複数のファイルを1つ1つコマンドに読ませる方法。

たとえば、POV-Rayのサンプルファイル(*.pov)を

一度にpovrayに処理させるには、
for i in $(ls *.pov); do povray $i; done
とする。

POV-Ray


POV-RayをPowerBookにいれた。

POV-RayはレイトレースでCGをレンダリングするためのプログラム。

OSX用のバイナリだとコマンドラインから実行できないみたいなので、ソースからコンパイル。

ここの下の方からソースコードをダウンロードして展開。

カレントディレクトリを移動して
$ cd povray-3.6.1
$ ./configure
$ make
$ sudo make install
これでok.
$ which povray
/usr/local/bin/povray
ちゃんと入っている。

2008年4月3日木曜日

MySQLでDBをコピーするには


MySQLのデータベースをコピーするには、mysqldumpを使います。

mysqldumpはデータベースをダンプするためのユーティリティですが、他のMySQLサーバにデータをコピーするのに使うことができます。

ローカルで動いているMySQLの中で"source_database"から"dest_database"へデータをコピーするには、以下のようにします。
$ mysqldump -u user_name -p -n source_database | ¥
> mysql -u user_name -p dest_database
-nオプションは、CREATE DATABASEを出力に含めないようにするためのオプションです。これを付けることで、dest_databaseへデータをコピーできるようになります。

また、コピーしたいデータベースが巨大なときもあるかと思います。そのようなときに、以下のオプションを使うと、さらに効率良くデータをコピーできます。

--no-autocommit
これをつけると、ダンプしたファイルの全体を囲うようにautocommit=0がつきます。つまり、後でリストア時にはファイル全体が1トランザクションになり、うまくいけば劇的に早くインポートできます。

--single-transaction
このオプションを付けてダンプすると--lock-tablesを付けずにシングルトランザクションでダンプするので、一貫性のあるデータをダンプできます。

--opt
これはデフォルトで有効になっているはずです。このオプションには--lock-tablesが含まれますがこれはダンプ時ではなくてインポート時にロックする指定です。


■関連エントリ

粒子法による固液連成シミュレーション


--

2008年3月31日月曜日

粒子法シミュレーション—物理ベースCG入門

東大の越塚さんによる、粒子法の入門書。

つい今月出版されたみたい。

Amazonで注文したが、まったく内容を見ずに買ったので少し不安。

粒子法というのは、流体をコンピュータでシミュレーションする方法の一つ。

流体をコンピュータでシミュレーションする方法には、大きく分けて2つある。

メッシュベースのアプローチとパーティクルベースのアプローチ。

粒子法というのは、その名の通りパーティクルベースのアプローチ。

パーティクルベースのアプローチの特徴は、波や雫といった液体表面の状態を表現しやすいという点。

この記事にソースコードを置いています。

■関連する記事
粒子法のプログラム第1回(概要)
粒子法のプログラム第2回(プログラムの大枠)
粒子法のプログラム第3回(データ構造)
粒子法のプログラム第4回(密度と圧力の計算)
粒子法のプログラム第5回(力の計算)
粒子法のプログラム第6回(境界条件と粒子位置の更新)
粒子法のプログラム最終回(粒子の出力)
【粒子法】粒子を流体としてレンダリング
3次元の粒子法シミュレーション
粒子法のシーンを2倍のサイズにしてみたが…
粒子法のシーンを2倍のサイズにしてみた



2008年3月26日水曜日

lifeclock


この時計は、人の一生を刻む時計です。

短針は、時の代わりに、その人が生まれてから何年経ったかを表します。

長針は、分の代わりに、直前の誕生日から何週間がたったかを表します。53週間で1週します。

長針の1週が1年にあたり、短針の1週がその人の一生にあたるのです。

(いまのところInternet Explorerでは動きません・・・)

使い方

あなたの生年月日を入れて以下のURLにアクセスすると、あなたのlifeclockが表示されます。

http://kamonama.sakura.ne.jp/lifeclock/lifeclock.html?yyyymmdd

yyyymmddを、あなたの生年月日で置き換えてください。

たとえば、1983年10月14日生まれの場合は、

http://kamonama.sakura.ne.jp/lifeclock/lifeclock.html?19831014

となります。

ブログパーツとして使うには

iframeを使えば、ブログパーツとして使うことができます。

たとえば、上のlifeclockは、以下のようなiframeで書かれています。
<iframe 
src="http://kamonama.sakura.ne.jp/lifeclock/lifeclock.html?19831014"
marginheight="0"
marginwidth="0"
width="172"
height="172"
frameborder="0">

2008年3月24日月曜日

【Python】sqlobject0.10


SQLObjectが0.10にバージョンアップした。
SQLObjectはPythonのO/Rマッパ。

以下、変更点の一部を紹介。
  • SQLBuilderのSelectがreversed, distinct, joinsをサポート

  • SelectResults.througnToが追加された。オブジェクトを経由せずに、SQLでリレーションをトラバースできる

  • ViewSQLObjectが追加された。

  • MySQLConectionがSSLに関係するパラメータを受け付けるようになった。ssl_key, ssl_cert, ssl_ca, ssl_capath。

  • etc.

Mac OSX トラブル対処法


Mac OSXがトラブったときの対処法。

電源を入れたあと、
  • shiftを押し続ける → セーフモード

  • Cを押し続ける → インストールディスクから起動

  • Dを押し続ける → ハードウェアテスト

  • マウスボタンを押し続ける → ディスクを取り出す

2008年3月19日水曜日

Erlang

Yaws
 HTTPサーバ

ElryWeb
 Webフレームワーク
 ORマッパつき

Mnesia
 データベース
 分散

2008年3月17日月曜日

【Python】メモリ使用量のいろんな調べ方


Pythonからプロセスのメモリ使用量を調べるにはいろいろ方法がある。
  • /proc(linux限定)
  • resourceモジュール
  • psコマンド
  • PSI(Python System Infomation)


2008年3月12日水曜日

【Python】Heapy: A Memory Profiler and Debugger for Python


http://guppy-pe.sourceforge.net/

Python用のメモリプロファイラ。
メモリに関する問題のデバッグ、最適化に使う。

・Data gathering
  オブジェクトヒープの情報を集める
・Data processing
  参照点からの差分を計算する
・Presentation
  プロファイリング結果を表示する

といった機能がある。

作者のmaster's thesisのpdfをダウンロードできる。
背景、設計、実装などが記述してある。

"Heapy: A Memory Profiler and Debugger for Python"
http://guppy-pe.sourceforge.net/heapy-thesis.pdf


drawing: hare

2008年2月28日木曜日

【Python】for文よりリスト内包記法のほうが効率的?


中で関数を呼ぶときは、for文を使うよりリスト内包記法を使う方が処理が速い。たとえば、以下の2つのスクリプトを試してみる。

■1.py
def a():
return "a"
s = ""
for _ in range( 0, 30000000 ):
s += a()
■2.py
def a():
return "a"
s = "".join( a() for _ in range( 0, 30000000 ) )

結果は以下。

■1.py
$ time python 1.py
real 1m45.563s
user 0m19.425s
sys 0m0.974s
■2.py
$ time python 2.py
real 1m4.801s
user 0m11.880s
sys 0m1.139s
後者の方が約60%の時間で処理できている。


ちなみに、以下のように関数呼び出しをなくしてみると・・・

■1.py
s = ""
for _ in range( 0, 30000000 ):
s += "a"
■2.py
s = "".join( "a" for _ in range( 0, 30000000 ) )

結果はこれ。

■1.py
$ time python 1.py
real 0m14.904s
user 0m13.254s
sys 0m0.972s
■2.py
$ time python 2.py
real 0m35.092s
user 0m6.368s
sys 0m1.125s
逆に、後者の方が遅いという結果に。

「文字列の連結は遅くて.joinメソッドは速い」という問題では
ないみたい。ちゃんと調べるか・・・?

splitコマンド


splitコマンドを使うと任意の行数やバイト数でファイルを分割できる。

使用例

split -l 1000 original.txt file-
original.txtを1000行ごとに分割する。
分割したものはfile-**という名前のファイルに保存される。

出力されるファイル

file-aa
file-ab
...

2008年2月27日水曜日

Yahoo! Open Search Platform



サードパーティが検索結果にコンテンツを追加できる。

たとえば、画像や構造化されたデータなど。

立ち上げのパートナーには、YelpやNewYorkTimesが含まれている。


manページの後ろの番号って何?


manページの後ろには番号がついている。
たとえば、grep(1)。
これはセクション番号を表している。
セクション番号とその意味は以下(ソース













1だれもが実行できるユーザコマンド
2システムコール、つまり、カーネルが提供する関数
3サブルーチン、つまり、ライブラリ関数
4デバイス、つまり、/dev ディレクトリのスペシャルファイル
5ファイルフォーマットの説明、例 /etc/passwd
6ゲーム(説明不要だろうネ)
7その他 例: マクロパッケージや取り決め的な文書
8システム管理者だけが実行できるシステム管理用のツール
9Linux 独自のカーネルルーチン用のドキュメンテーション
n新しいドキュメンテーション:よりふさわしい場所に移動されるだろう
o古いドキュメンテーション 猶予期間として保存されているもの
l独自のシステムについてのローカルなドキュメンテーション


・・・あれぇ?へんな空白ができてしまう。

2月の投稿は20本


2月は結局20本の投稿。

なかなか1日1本にならへんなぁ・・・

Mac OSXからLinuxにNFSマウントしたい


Mac OSXからLinuxにNFSマウントしたい。

ここによると、

ファイダーのサーバへつなげるでは駄目。

ターミナルから、
sudo mount_nfs -P 接続先:ディレクトリ マウント先

とする。
はずす時は、
sudo umount マウント先

とのこと。

やってみると、
$ sudo mount_nfs -P 192.168.1.10:/home/someone  /tmp
mount_nfs: can't access /home/someone: Permission denied
といわれてしまう。うむ・・・
どうすれば?

結局、sshfsを使いました。

2008年2月25日月曜日

パリに親子丼の店を開いてみては?


パリに親子丼の店を開いてみてはどうだろう?

親子丼を食べたあとにエスプレッソがでてくるような。

そう。店の前には、親子丼のスケッチを飾る。

だって、パリの人にとっては未知の食べ物だから。

日本にあるような見本は良くない気がする。

色鉛筆で色づけしたスケッチが良い。

写真でもいいかも?

うーん、やっぱりスケッチがいいでしょ。

白地に青の模様が入ったどんぶりに、

黄色いふわふわ卵が入っていて、

ワンポイントに緑のみつ葉がのっかってるの。

おぉ。いい感じだ。

2008年2月15日金曜日

Google docsからBloggerに投稿

おお!Google docsからBloggerに投稿できるではないか。

Load averageって?


Load averageとは、処理を待っているプロセスの平均数。

1以下なら、処理を待っているプロセスがないということ。

1を超えると、処理を待っているプロセスがあるということ。

muninのcpuspeedモジュール


muninには、cpuspeedというモジュールがある。

CPUの動作周波数をレポートする。

単位はMHz。たとえば、3.0kと表示されれば3GHz。

CPUの動作周波数は、/proc/cpuinfoから取得される。

2008年2月11日月曜日

消息筋ってダレ?


ニュースによく出てくる「消息筋」

いったい誰?と思って調べてみると、

「ある方面の事情に詳しい人々や関係機関」

のことだそうな。大辞泉より。

要は、報道の世界で情報源を匿名にするときの決まり文句ということか。

2008年2月10日日曜日

Gaucheのリスト操作ライブラリ SRFI-1


リストを操作するためのライブラリ集。

たとえば、1から100までの整数のリストを作ってみるとこんな感じ。
gosh> (iota 100 1)
(1 2 ... 100)

スペル修正プログラムはどう書くか


「スペル修正プログラムはどう書くか」Peter Norvig / 青木靖 訳
http://www.aoky.net/articles/peter_norvig/spell-correct.htm

Googleでたとえば"speling"を検索すると、
もしかして:"spelling"じゃないかといってくる。

このようなスペル修正プログラムを実際にかいてみせている。

IMF Data Mapper



各国のGDP成長率をグラフ&時系列で見ることができる。

IMF Data Mapper
http://www.imf.org/external/datamapper/index.php

iPhone/iPod touchアプリ 「Balls」


iPhone/iPod touch用に、アプリケーションを作ってみました。
上からボールがぽこぽこ落ちてきます。

ダウンロード

http://kamonama.sakura.ne.jp/download/Balls.zip

使い方

1.ダウンロードして解凍する
2.makeする
3.iPhone/iPod touchに転送する

パッチョの日記


http://home.tokyo-gas.co.jp/pa-cho/diary/index.html

火グマのパッチョがおもしろい。

2008年2月8日金曜日

hypertable


http://www.hypertable.org/index.html

Bigtableクローン。

Hadoopの上で動かす方法がドキュメントに書いてある。

2008年2月7日木曜日

ケニアの民族対立


大統領ムワイ・キバキ率いるのキクユ人と
最大野党党首ライラ・オディンガ率いるルオ人との対立。

きっかけは、昨年12月27日の大統領選。
ずさんな選挙で、国連視察団が勝者を決められないと判断するほど。

オディンガ側が大統領に求めているのは、
大統領選での敗北をみとめて辞任し、直ちに再選挙を行うこと。

一方、大統領側ではそれを拒否している。

最近はやりのアーティスト?


最近のはやり?
  • 9mm Parabellum Bullet
  • いきものがかり
  • Perfume
  • たむらぱん
  • RADWIMPS
9mm Parabellum Bullet、読みにくい!覚えられん!
Prefumeは何度か聞いたことがあった。女の子3人のグループ。
たむらぱんはmyspaceがオフィシャルサイトになってる。へぇ。

オバマとジョン・F・ケネディ


オバマとジョン・F・ケネディは似ているところがある。

オバマは黒人初の大統領を目指している一方、
ケネディはカトリックで初めての大統領だった。

似ているところもあるが、異なるところもある。

ジョン・F・ケネディの裏では、弟のロバート・ケネディが
政治に必要な裏の仕事を一手に引き受けていた。

オバマにはそういう存在がいまのところいないようだ。

2008年2月5日火曜日

スパーチューズデーって?

アメリカ大統領選で、多数の州が集中的に
予備選・党員集会を行う火曜日のこと。

その結果が選挙戦の行方を左右することから
「スーパーチューズデー」と呼ぶ。

なんで火曜日?

米国では多くの選挙が火曜日に行われる。
大統領選の本選挙も、連邦法で
「11月の第1月曜日の次の火曜日」と
決められている。(今年は11月4日)

もともと米国では休日である土、日曜は
「私的な時間」「安息日」と考えられ、
公的行事の選挙を行うことには抵抗があった。

また、国土が広大なため、有権者は
遠く離れた投票所に向けて月曜日に出発し、
火曜日に投票するのは通例になった、といわれる。

 

2008年2月3日日曜日

SICP第1章を読了

SICP第1章を読了。

プログラミング言語は以下を備えていることが大事:

  • プリミティブな式

  • 式を組み合わせられること

  • 式を抽象化できること

また、反復型のプロセスと再帰型のプロセスという2種類の計算プロセスについてと、時間と空間のオーダーという考え方。

最後の方でいろんな数学の手法を抽象化してみせた部分は感動ものだった。

節分


今日は節分。節分ってわける感がすごいな。節に分けるのか。

2008年1月31日木曜日

TracのBrowse Sourceをシンタックスハイライトするには

tracpygmentsをインストール

$ svn export http://trac-hacks.org/svn/tracpygmentsplugin/ \
branches/matt_good ./tracpygments
$ cd tracpygments
$ sudo python setup.py install

trac.iniに設定

[components]
tracpygments.* = enabled

Apacheを再起動

$ apache2ctl graceful

http://trac-hacks.org/

2008年1月30日水曜日

OSXのグローバルホットキーを使うには


OSXにはグローバルホットキーという機能があります。

グローバルホットキーに登録しておくと、
フォアグラウンドにでなくてもキー入力を受け取ることができます。

フレームワークの登録

Carbon.frameworkをリンクします。

コールバック関数を用意する

グローバルホットキーが押されたときに呼ばれるコールバック関数を定義します。
OSStatus MyHotKeyHandler( EventHandlerCallRef nextHandler,
EventRef theEvent,
void *userData )
{
// ホットキーが押されたときの処理
return noErr;
}

ホットキーを登録する

アプリケーションを初期化するときに以下を実行します。cmd+opt+spaceを登録してみます。
EventHotKeyRef gMyHotKeyRef;
EventHotKeyID gMyHotKeyID;
EventTypeSpec eventType;

eventType.eventClass = kEventClassKeyboard;
eventType.eventKind = kEventHotKeyPressed;
InstallApplicationEventHandler( &MyHotKeyHandler,
1,
&eventType,
NULL,
NULL );

gMyHotKeyID.signature = 'htk1';
gMyHotKeyID.id=1;
RegisterEventHotKey( 49, /* space key */
cmdKey + optionKey,
gMyHotKeyID,
GetApplicationEventTarget(),
0,
&gMyHotKeyRef );


2008年1月27日日曜日

Starbucks in Osaki Gate City is pretty good!


i love to read books at Starbucks.

today i went to the Starbucks in Osaki Gate City for the first time, to which i have really taken because of it's size and silence. additionally, it opens until 10:30 PM!

by the way, i heard that Starbucks of the United States has a tight game recently. McDonald begins to open cafe and take Starbucks' customers. McDonald?

in Japan, Starbucks is place to enjoy the atmosphere with time and space, not only some coffee. On the other hand, McDonald is place to eat hamburgers. it's impossible for McDonald to be the alternative to Starbucks in Japan!

大崎ゲートシティのスタバが、so good!


よくスタバで本を読むんですが、

今日行った大崎のゲートシティにあるスタバがいい感じ。

広くて静かで。

しかも、夜10:30まで開いてる。

次からここにしよ。


それはそうと、アメリカの本家スタバはここ1年大苦戦らしい。

米マクドがカフェを展開し始めて、

そっちにお客をとられているんやって。

えー、マクドが?


日本やと、スタバには単にコーヒーだけでなくて、

時間とか場所といった雰囲気を買いに行っている感じがする。

マクドはやっぱりハンバーガーを食べるところ。

スタバの代わりにマクドなんてありえん(笑)

2008年1月26日土曜日

Applescriptで現在の日付を取得するには


Applescriptで現在の日付を取得するには、以下のようにする
set d to current date
文字列として使いたい場合は以下のようにする
set d to ( current date ) as string

2008年1月20日日曜日

<iframe>タグ覚え書き

<iframe>タグを使う時に、余計な枠や余白をなくしたいことがよくあります。そんなときのための覚え書き。
<iframe src="http://kamonama.sakura.ne.jp/bar/bar.html" 
marginheight="0"
marginwidth="0"
width="155"
height="30"
frameborder="0">
</iframe>

marginheight・・・縦マージン
marginwidth・・・横マージン
width・・・横幅
height・・・縦幅
frameborder・・・枠のあり/なし

SICP (Structure and Interpretation of Computer Programs)を読み始めた


SICP (Structure and Interpretation of Computer Programs)を読み始めた。

書いてあることはプログラミングとしてはほとんど分かる。けど、情報科学という意味でちゃんとやったことがないので、最後まで読みたい。

章の前書きも含め、小節が全部で121個ある。1つ1つの小節はせいぜい5ページ程度。これならなんとかなりそう。

できるだけ、1日1小節は読み進めたい。そんなわけで、はてなグラフを貼付けました。

しし鍋を食べにいった

しし鍋を食べにいった。野獣鍋、¥6800なり。

しし鍋と、鹿刺しと、鹿肉の立田揚げ、しし肉の煮込み、熊汁を食べた。

猪と鹿と熊。

1回食べたらもういっか、という味。しょっぱい。

昔のごちそうやってんな。

2008年1月18日金曜日

CocoaからAppleScriptを送信するには

CocoaからAppleScriptを送信するには、NSAppleScriptを使う。
NSDictionary*  errorDict;
NSAppleScript* scriptObject;
NSString* source;

source = @" \
tell Application \"CotEditor\"\n \
make new document\n \
activate\n \
end tell";

scriptObject = [[[NSAppleScript alloc] initWithSource: source]
autorelease];
[scriptObject executeAndReturnError: &errorDict];

【Mochikit】itemgetter

itemgetter(name)
という関数がある。
obj[name]を返すような関数function(obj)を返す。
へぇ〜。

例:
a = b = {};
a.x = 1;
b.x = 3;
x_getter = itemgetter( "x" );
assert( x_getter( a ) == 1 );
assert( x_getter( b ) == 3 );

2008年1月13日日曜日

【Hadoop】1つのノードで疑似的に分散処理

Hadoopのドキュメント"Hadoop Quickstart"の日本語訳です。
【Hadoop】クイックスタートの続きです。

設定ファイルを用意する

すでに用意されているconf/hadoop-site.xmlを使います。
<configuration>
<property>
<name>fs.default.name</name>
<value>localhost:9000</value>
</property>
<property>
<name>mapred.job.tracker</name>
<value>localhost:9001</value>
</property>
<property>
<name>dfs.replication</name>
<value>1</value>
</property>
</configuration>

パスフレーズなしでsshできるようにする

パスフレーズなしでローカルホストにsshできるか確認します。
$ ssh localhost
もしパスフレーズなしでsshできない場合は、以下を実行します。
$ ssh-keygen -t dsa -P '' -f ~/.ssh/id_dsa
$ cat ~/.ssh/id_dsa.pub >> ~/.ssh/authorized_keys

実行する

単一プロセスと同じことを疑似分散でやってみます。

新しい分散ファイルシステムを用意し、Hadoopデーモンを開始します。
$ bin/hadoop namenode -format
$ bin/start-all.sh
入力ファイルを分散ファイルシステムにコピーします。
$ bin/hadoop dfs -put conf input
分散モードの場合も、単一プロセスの場合と同じようにgrepを実行することができます。
$ bin/hadoop jar hadoop-*-examples.jar grep input output ¥
> 'dfs[a-z.]+'

結果を表示する

出力結果を分散ファイルシステムからローカルのファイルシステムにコピーして表示してみます。
$ bin/hadoop dfs -get output output
$ cat output/*
分散ファイルシステム上の出力結果を直接見ることもできます。
$ bin/hadoop dfs -cat output/*

終了する

終わったら、デーモンを止めます。
$ bin/stop-all.sh

おまけ

Hadoopデーモンのログは${HADOOP_LOG_DIR}に書き込まれます。
デフォルトでは${HADOOP_HOME}/logsです。

【Hadoop】クイックスタート

Hadoopのドキュメント"Hadoop Quickstart"の日本語訳です。

ダウンロードする

ここからHadoopをダウンロードし、展開します。

準備をする

conf/hadoop-env.shを編集しJAVA_HOMEを指定します。

動作確認のために、Hadoopスクリプトの使用方法を表示してみます。
$ bin/hadoop
Hadoopスクリプトの使用方法が表示されれば、正しく動作しています。

単一プロセスで実行してみる

まず、単一プロセスでHadoopを実行してみます。デフォルトでは、Hadoopは単一プロセスで実行するようになっています。

入力されたファイルの中から、与えられた正規表現にマッチする行をすべて表示してみます。入力にはconfディレクトリのコピーを使い、出力先にはoutputディレクトリを指定してみます。
$ mkdir input
$ cp conf/*.xml input
$ bin/hadoop jar hadoop-*-examples.jar grep input output \
> 'dfs[a-z.]+'
$ cat output/*

2008年1月10日木曜日

2008年1月8日火曜日

【Python】いとこモジュールのインポート

モジュールが以下のような構成になっているときに、module.a.a1からmodule.b.b1をインポートしたいことがある。
module
- a
- a1
- a2
- b
- b1
- b2

■2.4まで
2.4まではこのケースを相対パスで書くことはできない。
なので、おとなしくmoduleからのパスを書くしかない:
# module.a.a1.py
from module.b.b1 import hoge
 
■2.5から
2.5から、相対パスでモジュール内をたどってインポートできるようになった。このように書くことができる:
# module.a.a1.py
from ..b.b1 import hoge

モジュール階層が深くなってくるとうれしいはず。

2008年1月5日土曜日

iPod touchでrogue、可不可?

iPod touchでrogueが動かないかと思い、試してみた。

コンパイル


  1. ここからダウンロード
    Rogue 3.6.3 Source Code (r1043) 

  2. Makefileを以下のように変更してコンパイル

    ■変更前
    CC=gcc
    ■変更後
    CC=/usr/local/bin/arm-apple-darwin-gcc

試してみた


コンパイルしたファイルをiPod touchに送って、実行。

・・・

「横幅がたりないー!」って言われた・・・
横70以上、縦22以上が必要みたい。
残念。

rogue


http://rogue.rogueforge.net/

プラハ国立美術館展「ルーベンスとブリューゲルの時代」

200余年の歴史を誇るプラハ国立美術館は、ヨーロッパ有数の美術館として西洋絵画の名品を所蔵しています。

本展では同美術館の所蔵品のうち、16世紀から17世紀にかけて描かれたフランドル絵画黄金期の秀作70点を展覧するものです。農民の日常生活を数多く絵にしたブリューゲル一族や、ヨーロッパ各国の王侯貴族から絵の注文が殺到した人気画家ルーベンスなど、宗教画、風景画、風俗画、静物画と、あらゆるジャンルの作品が出品されます。それらの絵には、比喩や教訓的な意味合いが含まれたものも多く、絵を読むという視点からも楽しめる作品群といえるでしょう。

ヨーロッパの古典絵画の名品と、描かれた絵に込められた謎に触れる絶好の機会です。

プラハ国立美術館展「ルーベンスとブリューゲルの時代」
そごう美術館
2008年1月2日(水)〜1月23日(水)

cocopad所感

cocopad

さくさく動く感じが良い。
ラインや方眼、ガラスの切り替えが気持ちいい
タイプできればいいのに。全部手書きはキツい。
手書きするとベクターの円や四角になるというのはどう?
自由度を下げて使いやすくする。
名前をつけて保存は面倒。勝手に保存しといて欲しい
規定のサイズにぱっとできるのも良い
http://projects.digitalwaters.net/index.php?q=cocopad
iPhone版ができたら良い?

2008年1月4日金曜日

ぬば玉

ヒオウギの果実。

赤い美しい花が漆黒の玉となる。

「ぬば玉」もしくは「うば玉」と呼ばれ、

古来より黒の枕語として使われてきた。

Hello jiggy!

I've made a HelloWorld with "jiggy", the javascript-based IDE for iPhone/iPod touch. I'll show the steps.

Login to jiggy

  1. launch jiggy in your iPhone/iPod touch.
  2. access to jiggy with the browser of your Mac/PC. the IP to access is on the screen.
  3. login to jiggy with the User name and the Password on the screen.

Create new application -- HelloJiggy

  1. click "Create new application" on the left navigation.
  2. input "Title" and "Identifier." Below is an example. ("Description", "Version", "Author", "URL", "Icon" are your options)
    Title:      HelloJiggy
    Identifier: jp.kamonama(.hellojiggy)
  3. click "Save" below the form.
  4. it's ok that your iPhone/iPod touch has been locked.

Hello jiggy!

  1. unlock your iPhone/iPod touch and launch jiggy.
  2. access to jiggy.
  3. click "HelloJiggy" come up on the left navigation.
  4. click "main.js" in the "HelloJiggy" folder.
  5. write the code below in main.js.
    Plugins.load( "UIKit" );

    var rect = UIHardware.fullScreenApplicationContentRect
    rect[0] = rect[1] = 0;

    var window = new UIWindow( rect );
    window.setHidden( false );
    window.orderFront();
    window.makeKey();
    window.backgroundColor = [ 1.0 , 1.0 , 1.0, 1.0 ];

    var hello = new UITextLabel( rect );
    hello.text = "Hello jiggy!";
    hello.backgroundColor = [ 0.6 , 0 , 0 , 1.0 ];
    hello.setFont( new Font( "Trebuchet MS" , 2 , 46 ) );
    hello.color = [ 1 , 1 , 1 , 1 ];
    hello.centersHorizontally = true;

    window.setContentView( hello );
  6. click the floppy-disk icon to save main.js.
  7. click the right-pointing triangle icon to run HelloJiggy. (or tap the HelloJiggy icon on the spring board of your iPhone/iPod touch.)
  8. Congratulation! HelloJiggy is running, now!

jiggyでHelloWorld

iPhone/iPod touchのjavascriptベースの開発環境「jiggy」で、HelloWorldを作ってみます。ただ、HelloWorldはあんまりなので、HelloJiggyに。

jiggyにログインする

  1. iPhone/iPod touchでjiggyを起動する
  2. 手元のMac/PCのブラウザからjiggyにアクセスする
    (jiggyに表示されているIPを使う)
  3. jiggyに表示されているUser name/Passwordでログインする

新しいアプリケーションを作る

  1. 左ナビの「Create new application」をクリックする
  2. TitleとIdentifierを入力する(Description, Version, Author, URL, Iconは任意)
    Title:      HelloJiggy
    Identifier: jp.kamonama(.hellojiggy)
  3. 画面下部の「Save」をクリックする
  4. iPhone/iPod touchがロック状態になればOK

Hello jiggy!

  1. ロックを解除し、jiggyを起動する
  2. jiggyにアクセスする
  3. 左ナビに「HelloJiggy」が追加されているので、それをクリックする
  4. HelloJiggyフォルダ中のmain.jsをクリックする
  5. 以下のコードをmain.jsに記述する
    Plugins.load( "UIKit" );

    var rect = UIHardware.fullScreenApplicationContentRect
    rect[0] = rect[1] = 0;

    var window = new UIWindow( rect );
    window.setHidden( false );
    window.orderFront();
    window.makeKey();
    window.backgroundColor = [ 1.0 , 1.0 , 1.0, 1.0 ];

    var hello = new UITextLabel( rect );
    hello.text = "Hello jiggy!";
    hello.backgroundColor = [ 0.6 , 0 , 0 , 1.0 ];
    hello.setFont( new Font( "Trebuchet MS" , 2 , 46 ) );
    hello.color = [ 1 , 1 , 1 , 1 ];
    hello.centersHorizontally = true;

    window.setContentView( hello );
  6. フロッピーディスクのアイコンをクリックし、main.jsを保存する
  7. 右向き三角のアイコンをクリックし、HelloJiggyが起動すればOK
    (スプリングボードのHelloJiggyアイコンをタップしても起動できます)

2008年1月3日木曜日

Bloggerで投稿のスタイルを変更するには

Bloggerでは投稿のスタイルを自由に変更することができます。

手順

  1. 管理画面を表示する

  2. 「テンプレート」タブをクリックする

  3. 「HTMLの編集」リンクをクリックする

  4. テキストボックスにテンプレートが入っている。
    <b:skin>〜</b:skin>の中にスタイルが定義されているので、必要に応じて編集する。


具体例


<pre>〜</pre>タグのボーダーを破線にし、背景を薄い黄色にしてみる。

■変更前
1 -module(echo).
2 -export([echo/0]).
3
4 echo() ->
5 receive
6 Msg -> io:format("~necho: ~p", [Msg]),
7 echo()
8 end.

■変更後
1 -module(echo).
2 -export([echo/0]).
3
4 echo() ->
5 receive
6 Msg -> io:format("~necho: ~p", [Msg]),
7 echo()
8 end.

Erlangで関数の実行にかかった時間を計るには

Erlangで関数の実行にかかった時間を計るには、timer:tcを使います。
timer:tc( パッケージ名, 関数名, 引数リスト ).
出力は、{ かかった時間[microsec]、関数の出力 }

実行例

1> fib(40).
165580141
2> timer:tc(fib, fib, [40]).
{8464307,165580141}


ノンデザイナーズ・デザインブック




ノンデザイナーズ・デザインブック

お薦め本。
センスのいいレイアウトのための"コツ"が紹介されています。

1.まず、悪い(やってしまいがちな)例を載せて、
2.次にそれを良くした例を載せる。

という形で書いてあるので、自分ですぐに試すことができます。

「で、具体的には、どうやればいいの?」がわかる、良書です。

太助で牛タン

明日(2008.1.4)は虎ノ門の太助で牛タンを食べよう。

Pieter Brueghel the Elder

Pieter Brueghel the Elder

16世紀のフランドルの画家

Pieter Paul Rubens

Pieter Paul Rubens

17世紀、バロック期のフランドルの画家及び外交官

2008年1月2日水曜日

The Erlang mode for Emacs

EmacsのErlangモードを使うには、.emacsに以下のコードを追記します。
(setq load-path (cons "/usr/local/otp/lib/tools-<ToolsVer>/emacs"
load-path))
(setq erlang-root-dir "/usr/local/otp")
(setq exec-path (cons "/usr/local/otp/bin" exec-path))
(require 'erlang-start)


  • /user/local/otpは、Erlangのインストールパスで置き換えてください。
  • 1行目の<ToolsVer>には、インストールされているErlangのバージョンが入ります。実際のディレクトリ構成に合わせてください。

詳しくは:
http://www.erlang.org/doc/apps/tools/erlang_mode_chapter.html

Erlangでメッセージパッシング

Erlangは言語レベルでメッセージパッシングの仕組みを持っている。
そこで、Erlangでメッセージパッシングのコードを書いてみた。

サンプルコード(echo.erl)

1 -module(echo).
2 -export([echo/0]).
3
4 echo() ->
5 receive
6 Msg -> io:format("~necho: ~p", [Msg]),
7 echo()
8 end.

1行目
 モジュールの宣言

2行目
 echo関数のエクスポート

4〜8行目が本体の関数です。メッセージを受け取り、標準出力に出力します。再帰することでなんどもメッセージを受け取るようにしています。具体的には以下。

4行目
 関数の開始。関数名は"echo"。

5行目
 receive式でメッセージを受け取る。

6行目
 Msgにメッセージが入る。io:formatで、受け取ったMsgを標準出力に出力する。

7行目
 echo関数を再帰的に呼ぶ。こうすることで、なんどもメッセージを受け取ることができる。

8行目
 receive式を閉じる。

実行例

$ erl
1> c(echo.erl).
{ok,echo}
2> Pid = spawn( echo, echo, [] ).
<0.40.0>
3> Pid ! "Hello!".
"Hello!"
echo: "Hello!"
1行目
 echo.erlをコンパイル。

2行目
 echo関数をプロセスとして呼び出し。
 プロセスIDが返ってくるのでPidに入れる。
 <0.40.0>は返ってきたプロセスID。

3行目
 !演算子を使って、プロセスPidに"Hello!"メッセージを送る。
 echo関数がメッセージを受け取り、「echo: "Hello!"」を出力する。
 なお、その前の「"Hello!"」は!演算の評価結果。