はじめに
いつもFotranで流体力学計算（初期値から始めてその時間発展）をしているものです。今回はそういう用途のプログラムのコーディングで僕が気をつけていることを紹介します。共同研究者に見せるために公開しますが、ちょっと内容が不足しているように思うので、加筆の予定です。
縦に見て、きれいに揃うようにする
コードをデバッグするとき、縦にコードをみていくことが多いです。その時変数の文字数がそろっているとデバッグしやすくなります。
以下ではn個の変数がある3次元のテーブル(i,j,k)の補間をしています。例えばテーブルからi方向についてiとi+1の間のri:1-riの位置にある値を取り出したいときを考えています。j,kについても同様です。こういうのは書き間違いが生じやすいのですが、以下のように書いておけば書き間違いに気がつきやすいです。
y(1:n) = (1.-ri)*( (1.-rj)*( (1.-rk)*x(1:n,k ,j ,i )
& + rk *x(1:n,k+1,j ,i ) )
& + rj *( (1.-rk)*x(1:n,k ,j+1,i )
& + rk *x(1:n,k+1,j+1,i ) ) )
& + ri *( (1.-rj)*( (1.-rk)*x(1:n,k ,j ,i+1)
& + rk *x(1:n,k+1,j ,i+1) )
& + rj *( (1.-rk)*x(1:n,k ,j+1,i+1)
& + rk *x(1:n,k+1,j+1,i+1) ) )
似た変数は同じ文字数にする
コードを縦に見ると良いですよという話を続けます。以下のコードのほうが、
va(i,j,k) = fa(i,j,k)*dt
vb(i,j,k) = fb(i,j,k)*dt
vc(i,j,k) = fc(i,j,k)*dt
以下のものより読みやすいでしょう。
vone(i,j,k) = fone(i,j,k)*dt
vtwo(i,j,k) = ftwo(i,j,k)*dt
vthree(i,j,k) = fthree(i,j,k)*dt
文字数が違う場合にはしょうがなく配列の引数や()、=の位置を合わせるとよいと思います。ただ、行の始まる位置（インデント）には意味があることが多いのに、それを無視することになるので次善の策です。
vone(i,j,k) = fone(i,j,k)*dt
vtwo(i,j,k) = ftwo(i,j,k)*dt
vthree(i,j,k) = fthree(i,j,k)*dt
変数をnamespaceにまとめる
物理の計算をしているとエネルギーをe、フラックスをfとかの変数で置きたくなりますが、いろんなエネルギー、いろんなフラックスがあるし、そうでない用途でe,fを使いたくなるときもあります。こういうよく使う変数の名前の被りをさけるためにnamespaceを使いましょう。 Fortranではmoduleです。
僕は流体力学の計算をよくしているのですが、そのようにある座標系の中で場が時間発展していく場合を想定します。その場合以下の3つを違うnamespaceにまとめると良いです。
- 場の変数（密度とか速度）
- 座標系（各点の位置とか体積）
- 物理定数等
場の変数と座標系を分けるのは，サブルーチンの中で場の変数をコピーしたような似た変数が現れることがあるので、そのときに座標系は使うけど場の変数は使わないことがあり得ます。例えば以下のようにRunge Kutta（の最初の部分）を考えます。
call func(rho,v,drho1,dv1)
rho1(:,:,:) = rho(:,:,:) + drho1(:,:,:) * dt
v1 (:,:,:) = v(:,:,:) + dv1(:,:,:) * dt
call func(rho1,v1,drho2,dv2)
このfuncの中では座標系は使いますが、場の変数はコピーしたものが使われます。そういう時のために場の変数と座標系を分けておいたほうが良いです。
メインループの中身を簡潔にする
メインプログラムはそのプログラムで何をしているの概要が分かるように書くのが良いです。あとで例を追加します。
参考文献
- Fortranをコーディングする際に気をつけていること ：Fortranのコードの書き方に関することはこの記事が参考になります。
- C++を使って数値計算をしている学部生へのアドバイス：C++の場合はこちらを参照のこと。