144. FDPS on Crystal (2018/8/21)
FDPS に Crystal 言語へのインターフェースをでっち上げて、
プログラムを全部 Crystal で書く、というのが
一応 動いた ので、
ちょっとだけ。
Crystal は、非常に Rubyに似た(同じプログラムが多くの場合にどっちでも動くほど)
文法をもつ、静的に型付けされたコンパイル言語です。
静的に型付けされたコンパイル言語であることのメリットは、
まずなんといっても実行速度が速いというか普通であることで、
普通に書いたプログラムは C/C++ で書いて最適化オプションを
つけてコンパイルしたのと同じ程度の速度で動きます。Ruby とか
Python では100倍遅いわけで、だいぶ幸せになれます。
「同じ程度」ならメリットない気もしますが、メリットは Ruby なみに簡潔に
書けることです。例えば、 FDPS のサンプルにあるリープフロッグ積分
の関数は、Crystal だと
def kick(psys_num,dt)
ptcl = return_psys_cptr(psys_num)
FDPS.get_nptcl_loc(psys_num).times{|i|
q = ptcl+i
q.value.vel = Vec_Float64.new(q.value.vel) +
Vec_Float64.new(q.value.acc)* dt
}
end
となります。まあその普通なら
q.value.vel = Vec_Float64.new(q.value.vel) +
Vec_Float64.new(q.value.acc)* dt
は
q.vel += q.acc* dt
ですんで欲しいところですが、ここは q の実体が C 言語構造体でCrystal 本
来のクラスではないため、一旦Crystal 本来のクラスに変換してから演算する
ので少し煩雑な記法になってます。また、FDPS の側にあるデータのポインタ
や粒子数へのアクセスも関数なので、これも煩雑な表記になっています。全
て Crystal の中なら、例えば
def kick(sys,dt)
sys.each{|q| q.vel += q.acc* dt}
end
ですむところです。逆にいうと、複雑な処理をするなら一旦
純粋 Crystal なデータ構造にコピーすればより簡潔に書けることになります。
C++ (FDPSのサンプル)だと
template<class Tpsys>
void kick(Tpsys & system,
const PS::F64 dt) {
PS::S32 n = system.getNumberOfParticleLocal();
for(PS::S32 i = 0; i < n; i++) {
system[i].vel += system[i].acc * dt;
}
}
で、FDPS の粒子型が C++ ネイティブで実装されているので普通に
system[i].vel += system[i].acc * dt;
と書けるのはいいのですが、それ以外はなかなか煩雑な記述になっている
ことがわかります。特に、関数宣言部の
template<class Tpsys>
void kick(Tpsys & system,
const PS::F64 dt) {
と
def kick(psys_num,dt)
の違いはなかなか大きいです。C++ は謎言語なので記述量が多いのでは、
という向きのために Fortran と C も用意してみました。まず
Fortran です。これは現行の FDPS の Fortran API のサンプルです。
subroutine kick(fdps_ctrl,psys_num,dt)
use fdps_vector
use fdps_module
use user_defined_types
implicit none
type(fdps_controller), intent(IN) :: fdps_ctrl
integer, intent(IN) :: psys_num
double precision, intent(IN) :: dt
integer :: i,nptcl_loc
type(full_particle), dimension(:), pointer :: ptcl
nptcl_loc = fdps_ctrl%get_nptcl_loc(psys_num)
call fdps_ctrl%get_psys_fptr(psys_num,ptcl)
do i=1,nptcl_loc
ptcl(i)%vel = ptcl(i)%vel + ptcl(i)%acc * dt
end do
nullify(ptcl)
end subroutine kick
えーと、その、 Fortran という名前がユーザーが使い始めるにあたって重要なのは
理解していますし、そのために Fortran API を開発したのですが、まあ、その、
ちょっと行数と文字数が多いような気もします。
今の Fortran にはクラスがあってクラスの演算、例えば
ベクトルの加算とか乗算も定義できるので、
ptcl(i)%vel = ptcl(i)%vel + ptcl(i)%acc * dt
と書けるのは悪くはないですが、現行の Fortran の大きな問題は、C++ でい
うところのテンプレートクラスが(まだ?)ないことで、例えば 3 次元ベクト
ルクラスは倍精度実数、単精度実数といった型毎に定義するソースコードが必
要になります。
Crystal ではテンプレートクラスと普通のクラスの区別はなくて、コンパイル
時の型推論でなんとかなればいいので、プログラムを書いたり読んだりする上
では邪魔な型宣言がほとんどないにもかかわらず、コンパイル時に全ての型
チェックはされていて、割合読みやすくて煩雑な繰り返しもないコードになる
傾向があります。
C ではどうかというと、例えば
void kick(int psys_num, double dt)
{
PFull_particle ptcl;
fdps_get_psys_cptr(psys_num,&ptcl);
int n = fdps_get_nptcl_loc(psys_num);
for (int i=0;i < n; i++){
PFull_particle pi = ptcl+i;
Cvec_Float64 *pv, *pa;
pv = &(pi->vel);
pa = &(pi->acc);
pv->x += pa->x * dt;
pv->y += pa->y * dt;
pv->z += pa->z * dt;
}
}
という感じです。これは、FDPS の Fortran インターフェースが実際には C
言語レベルでの API を定義していることを利用して、直接 その API を呼ぶ
ことで C 言語でプログラムを書くことを可能にしたものです。こんなのに使
い道はあるかというと、一つにはさらに他の言語への API を簡単につくれる
ようになるし、また C++ 難しい、という人が FDPS を使う Fortran 以外の方
法をもうひとつ用意できたことにもなります。C ではクラスやメンバー関数の
概念はないので、3次元ベクトルを表す演算子を定義とかできなくて、 C++ で
は
system[i].vel += system[i].acc * dt;
ですむところを手で展開して
pv->x += pa->x * dt;
pv->y += pa->y * dt;
pv->z += pa->z * dt;
としています。もちろんこれは気のきいたマクロを定義するとかで簡潔に書けない
わけではないですが、、、というところです。
一応行数とかバイト数を書いてみると
Lines bytes
Crystal 8 238
(Crystal 3 62) (FDPS ではない型使えたら)
C++ 8 236
C 15 404
Fortran 17 554
という感じで、Crystal での実装のソースのサイズは FDPS のネイティブな実
装をそのまま呼べるC++ のコードに比べても大差なく、他の言語に比べると圧
倒的に簡潔に書けることがわかります。
世の中にはコードが簡潔であっても意味不明な言語というのは太古の APL から
始まって無数にあるわけですが、Crystal は Ruby の文法を継承しているので
難しいことをしなければ読みやすいのではと思います。
とはいえ、Crystal には今のところ色々な問題もあり、特に大きな制限は
マルチスレッドでは動かないことです。もちろんそれ以前にまだ Ver 0.x であって何が起こるかわからないというのもありますがそれはおくとして。FDPS では、
OpenMP では動かないですが MPI では動く、という状態です。
ガーベージコレクタ周りがまだ色々あるようです。開発グループの
規模も小さいので(これは必ずしも悪いことではないですが)性能周りの改善には
時間がかかるかなという気もします。
なんで Julia じゃねーの?という意見もある気がしますが、現行の FDPS の
実装では必要な、C++ のメインプログラムからターゲット言語の(実効的な)メ
インプログラムを呼び出し、そこからさらにC API をターゲット言語のコール
バック関数付きで呼ぶ、ということができるほど私が Julia を理解してない
のでできてません。誰か C API から作って下さい。
素の Ruby や Python から呼べるインターフェースを生成することも
原理的にはできますが、性能のボトルネックが必ず発生するので
あえてしないほうがよい気がしています。今からフレームワークとか開発
するなら、 90年代のスクリプト言語である Python/Ruby ではなく、
2010年代の言語といえる Crystal/Nim/Julia/(あと他になんかありましたらなんでも)
で、コンパイラがあって自明な性能ボトルネックがないものを
使うほうが幸せになれるように思います。