iLog

About

Formuraではじめる数値計算

2017.12.04 Physics

はじめに

こんにちは、ishiy(@ishiy1993)です。

この記事は CAMPHOR- Advent Calendar 2017 5日目の記事です。

今回は Formura を使いはじめるために必要なことを解説します。

Formuraの紹介

Formura はステンシル計算のためのDSL(Domain Specific Language)です。

ステンシル計算とは、配列として保持している状態を各点まわりの情報だけで更新していくような計算のことです。ライフゲームや偏微分方程式の陽解法がその例になります。多くの物理現象は偏微分方程式で記述でき、コンピュータを使って数値的に偏微分方程式を解くことは超新星爆発、重力波などの物理の探求や車、飛行機などの製品の開発で不可欠です。

しかし、これらの大規模な数値計算コードの開発は非常に労力のかかる仕事です。可能ならば、本質的な部分(どういう方程式をどう離散化して解くのか)だけを記述し残りの最適化や並列化は機械的に行ってくれるのが理想です。この理想の実現にFormuraは挑戦しています。

余談

Formuraは、@nushioさんによって開発されました。僕は今年の2月ごろ@nushioさんに連絡をとり、今後一緒に開発をしていくつもりでした。しかし今年の7月に@nushioさんが夭折してしまったため、今後僕を含め複数人でFormuraの開発を進めていくことになりました。彼の「物理屋さんがもっと物理に集中できる世界をつくっていこう」という意思を継いで開発していきます。

Formuraのインストール

Formuraで数値計算をするには、Formuraの処理系(formura)とMPIのライブラリが必要になります。 これらのセットアップはDockerを利用すれば一瞬で完了します。

$ docker pull formura/formura

以下では、ユーザーがdockerグループに属していることを仮定します。 またLinuxでのみ動作確認を行っています。

Formuraの使い方

以下では、Formuraの使い方と文法を説明します。 適宜、formura-sampleのコードを参考にしてください。

Formuraはスキームを記述する fmr ファイルと設定を記述する yaml ファイルからC/Fortranのライブラリを生成します。それぞれのファイルに書く必要があるものは次の節で解説します。 また、ライブラリに定義された関数を呼び出すコードも必要になります。(以下の diffusion-main.cpp)

% docker run -it --rm -v $PWD:/work -u $UID:$GID -w /work formura/formura bash
$ formura diffusion.fmr #diffusion.c, diffusion.h, diffusion_internal_*.cが生成される
$ mpicxx diffusion-main.cpp diffusion.c diffusion*.c -o diffusion
$ mpirun -n 1 ./diffusion

ここで mpirun-n オプションで指定しているのは並列度です。この値は yaml ファイルに書かれている並列度と一致している必要があります。

実行したときに、次のようなエラーが出る場合はファイルを出力するディレクトリが存在するか確認してみてください。

===================================================================================
=   BAD TERMINATION OF ONE OF YOUR APPLICATION PROCESSES
=   PID 46 RUNNING AT 997dc4cc3a53
=   EXIT CODE: 139
=   CLEANING UP REMAINING PROCESSES
=   YOU CAN IGNORE THE BELOW CLEANUP MESSAGES
===================================================================================
YOUR APPLICATION TERMINATED WITH THE EXIT STRING: Segmentation fault (signal 11)
This typically refers to a problem with your application.
Please see the FAQ page for debugging suggestions

(補足: 現在のFormuraは1次元、2次元の計算におけるMPIでの並列化にバグがあり、2並列以上ではうまく動きません。この問題は将来的には解決される予定です。)

Formuraの文法

dimensionとaxes

fmr ファイルの一番上にあるのが dimension 宣言と axes 宣言です。 これらはその名の通り対象とする次元と軸の名前を与えるものです。

関数定義とラムダ式

Formuraには2種類の関数定義の書き方があります。

begin function res = add1(x)
  res = x + 1
end function
begin function add1(x) returns res
  res = x + 1
end function

これらはどちらも1を足す関数 add1 の定義です。

またFormuraはラムダ式もサポートしているので、 add1 をもっと簡潔に表現できます。

add1 = fun(x) x + 1

init関数とstep関数

fmr ファイルでは init 関数と step 関数を定義する必要があります。

init関数では, 扱う状態に識別子と型を与えます。Formuraはいろいろな型を提供していますが、とりあえず、 double 型、Grid型、タプル型あたりを知っておけばよいと思います。 Grid型とタプル型については後の節で解説します。Formuraにおいて変数宣言は

<型名> :: <変数名>

という形で行います。

step 関数では、状態を更新する計算を記述します。 Formuraはステンシル計算のみを扱うので、配列のある点の値をどう更新するかだけを記述すればよいです。 サンプル(diffusion.fmr)では、3次元拡散方程式を解いているので参考にしてください。

if式

Formuraは if 式をサポートしており、次のように使用できます。

min = fun(a,b) if a < b then a else b

Grid

Grid型は配列のようなものです。

double [] :: xs

とすると double のGrid型を宣言できます。

Grid型は近傍の値を参照することができ、これを使うと差分式は次のように書けます。

d_xx = fun(q) (q[i+1] + q[i-1] - 2*q[i])/h/h

ここでiは注目点のインデックスを表しています。 よって q[i]は注目点の値、q[i+1]は注目点の右隣の格子点の値を意味します。

タプル

タプルはいくつかの値の組です。 Formuraにおいてタプルはベクトルのような振る舞いをします。 つまり、タプル同士の和 (a1,a2) + (b1,b2)(a1+b1,a2+b2) になります。 また、スカラー倍 c * (a1,a2)(c*a1,c*a2) となります。 この機能のおかげで、オイラー方程式を解くコードなどを簡潔に表現できます。

さらにタプルqsの1つ目の要素を取り出すには qs 0 とします。 もしくは

(q1,q2) = qs

のようにパターンマッチで取り出すこともできます。

設定ファイル

解像度や並列度などはyamlファイルで指定します。

initial_walls:
    x: [100]
intra_node_shape: [1000]
monitor_interval: 4
mpi_grid_shape: [1]
temporal_blocking_interval: 2

yamlファイルの例を上に示しています。これらのプロパティが必要になります。 intra_node_shape には、1ノードあたりの解像度を指定します。 mpi_grid_shape には、各軸ごとのノード数を指定します。合計のノード数(並列度)は mpi_grid_shape のリストの要素をすべて掛けあわせたものになります。 initial_wallstemporal_blocking_interval はTemporal Blocking関係のパラメータです。 Temporal Blockingは空間方向のブロッキングと時間方向のブロッキングを同時に行う高速化技法です。Formuraは、Temporal Blockingされたコードを生成します。inital_walls は空間方向のブロッキングの幅、temporal_blocking_interval は時間方向のブロッキングの段数を指定するものです。 monitor_interval は一度に進める時間ステップの数です。Temporal Blockingの関係で、monitor_interval は、temporal_blocking_interval2*N 倍である必要があります。 (N は自然数)

(補足: 現在Formuraがサポートしているのは台形型のTemporal Blockingです。そのため initial_wallstemporal_blocking_interval の間にも制約がありますが、ここでは書きません。)

Formuraが生成するもの

Formuraが生成したライブラリには次のようなものが定義されています。

Formura_Navigator 構造体は現在のタイムステップや各軸方向の添字の走る範囲などを管理するものです。また Formura_Init 関数はこの構造体を初期化するための関数です。 Formura_Forward 関数は、monitor_interval で指定した分だけ状態を更新する関数です。 また状態はグローバルな配列として保持されています。配列の名前はfmrファイル中のstep関数の引数の名前と同じです。

diffusion-main.cpp などの中ではこれらの関数を呼び出して計算を行います。 やることは、各種初期化(状態、Navigator、MPI)とファイルへの書き出しです。

まとめ

Formuraを使いはじめるのに必要なことを解説しました。(やや説明不足かもしれないですが…) Formuraは、まだまだ発展途上でありこれからも進化していきます。 興味のある方はぜひ使ってみてください。 また、formuraのexamplesにはより多くの例があるのでそちらものぞいてみてください。

参考ページ:

間違いなどありましたら、@ishiy1993に連絡していただけると助かります。

CAMPHOR- Advent Calendar 2017 の明日の担当は @marty1martie です。お楽しみに。

ではでは。