プログラム全体は, 制御を行うメインプログラム (arare.f90) と, 個々の設 定と計算を行う複数のモジュールおよびサブルーチンから構成されている. モ ジュールとサブルーチンは 気象庁標準コーディングルール にしたがい,
のいずれかに分類されるように作成している.
複数のサブルーチンで用いられる定数の値を決めるモジュール. 定数 の値はそのモジュール内で確定し, モジュールを呼び出す親プログラ ム内ではその定数の値は変更されない.
親プログラム側では use 文でモジュールを呼び出すだけで定数の値が 確定する.
複数のサブルーチンで用いられる変数の型, 変数が配列であった場合 の大きさ, 変数の値を決めるモジュール. 定数の値はそのモジュール 内で確定し, モジュールを呼び出す親プログラム内ではその定数の値 は変更されない.
親プログラム側では use 文でモジュールを呼び出すだけではなく, モ ジュール内に含まれる変数の値を確定するサブルーチンを呼び出す必 要がある.
以下のモジュールが該当する.
同じような機能をもつサブルーチンをまとめたモジュール. このモジュー ル内で値を代入した変数は他のモジュールまたは外部サブルーチン内 でその値を変更されることがある.
以下のモジュールが該当する.
上記のモジュールには含まれないもの. 以下のサブルーチンが該当する.
プログラムの書法は基本的には 気象庁標準コーディングルール に準拠している. これに加え, 変数と変数配列に作用させる微分や平均操作を行う配列関数を, SPMODEL の書法に準じて記述することにより, ソースコードの可読性の向上をめざしている.
変数は
[格子点位置情報]_[変数名]_[時間方向添字]
と記述し, 配列関数は
[出力変数の格子点位置情報]_[演算名]_[入力変数の格子点位置情報]
と記述する. 変数名の具体例は変数名リストを参照.
その変数が配置されている格子点の位置を表す. x-z の 2 次元交互格子 (Arakawa C グリッド)の場合, 以下の接頭詞を付ける.
3 次元の場合は
とする.
基本的なルールは以下の通りである.
時間レベルを表す. 以下の添字を用いる.
A: 時刻 t+Δt (1 ステップ計算後:After ) B: 時刻 t (現在の時刻 :A の次なので B) C: 時刻 t-Δt (1 ステップ計算前: ) D: 時刻 t-2Δt (2 ステップ計算前: ) E: ...
座標 X, Lon
Y, Lat
Z, Rad, Press, Sigma
座標格子間隔 DelX, DelY, DelZ
暦(構造体にする?) Date
Date%Year, Date%Month, Date%Day,
Date%Hour, Date%Min, Date%Sec
時間 Time
時間格子間隔(長い) DelTLong
(短い) DelTShort
速度 VelX, VelLon
VelY, VelLat
VelZ, VelRad, VelPress, VelSigma
流線関数 Psi
速度ポテンシャル Chi
渦度 Vortx, RotVx
Vorty, RotVy
Vortz, RotVz, Zeta
水平発散 Div, DivVxVy
温度 Temp
仮温度 Vtemp
温位 Theta, Ptemp
仮温位 Vtheta, Vptemp
密度 Dens, Rho
圧力 Press
エクスナー関数 Exner
ジオポテンシャル Gph, Phi
トレーサー混合比 Q****
* 水蒸気 Qvap
* 雨水 Qrai
* 雲水 Qclw
* ダスト Qdust
降水量 Rain
円周率 Pi ステファン-ボルツマン定数 ?? 重力加速度 Grav 気体定数 GasR 定圧比熱 Cp 定積比熱 Cv 非熱比 Gamma 自転角速度 Omega 気体分子量 GasWt
ここではメインプログラム arare.f90 の処理のおおまかな流れを記す.
1 階微分演算を行う.
call ss_dx_ss(var) call ss_dz_ss(var) call ss_dx_fs(var) call ss_dz_fs(var) call ss_dx_sf(var) call ss_dz_sf(var) call fs_dx_ss(var) call fs_dz_ss(var) call fs_dx_fs(var) call fs_dz_fs(var) call fs_dx_sf(var) call fs_dz_sf(var) call sf_dx_ss(var) call sf_dz_ss(var) call sf_dx_fs(var) call sf_dz_fs(var) call sf_dx_sf(var) call sf_dz_sf(var)
平均操作を行う.
call ss_avr_ss(var) call ss_avr_fs(var) call ss_avr_sf(var) call fs_avr_ss(var) call fs_avr_fs(var) call fs_avr_sf(var) call sf_avr_ss(var) call sf_avr_fs(var) call sf_avr_sf(var)
設定ファイル名および出力ファイル名を設定する.
設定ファイル名はモジュール内のパラメータ cfgfile で指定する. 一方, 出力ファイル名は設定ファイルに記述したファイル名を 読み込むことで決定される.
call fileset_init
座標軸情報を設定する.
設定ファイルから格子点数, 領域幅の情報を取得する. これらの値から格子点間隔を計算する.
配列添字の範囲を定義する. 配列添字の範囲には物理的に意味をもつ領域と 計算に用いられる領域の 2 種類がある. 計算に用いられる領域の範囲は 物理的に意味のある領域に境界外の格子点数を加えることによって定義される.
各独立変数に対して, 格子点上および半格子ずれた点上の座標値を定義する.
call gridset_init
このモジュール内では境界条件の設定も行っている. 設定ファイルから x, z 方向の境界の種類(周期境界/固定境界/自由境界) を読み込み, 格子点ごとの境界条件(対称/反対称)に変換する.
この機能が gridset.f90 に含まれているのは開発経過の事情による ものであり, あえてこのモジュールに含めた特別な理由はない.
時間積分に必要なパラメータの設定をする.
設定ファイルから積分時間, 時間間隔の値を読み込み, ループ回数を計算する. あまりが出た場合は警告を出す.
クランクニコルソン係数の値が 0 以上 1 以下になっているかを チェックする. 条件を逸脱した場合, 警告を出してプログラムを終了させる.
call timeset_init
物理定数の値を設定する. 設定ファイルから惑星名を読み込んで その惑星のパラメータを設定する. 設定する物理定数は
である.
call physprm_init
2004/11/18 現在では地球の物理定数の値しか実装されていない.
netcdf ファイルに記述する情報を設定ファイルから読み込む
call expname_init
基本場の初期分布を決定する
call basicprm_init
将来的には初期分布をファイルから読み込むようにしてこのプログラムはなくなる.
行列計算に用いる線形計算ライブラリを決定する.
設定ファイルから線形計算ライブラリに SSL2 , LAPACK のどちらを用いるかを 読み込む.
SSL2 を選択した場合は 3 重対角行列の係数を引数として ssl2_linear に渡して計算を実行させる.
call linlib_init(N)
LAPACK にはまだ対応していないため LAPACK が選択された場合は メッセージを表示してプログラムを終了させる.
ssl2 を用いて線形計算を行うためのラッパー
call ssl2_lax_init(M) call ssl2_lax(A, B, sw) call ssl2_lsx_init(M) call ssl2_lsx(A, B, sw) call ssl2_lsix_init(M) call ssl2_lsix(A, B, sw) call ssl2_ltx_init(M) call ssl2_ltx(A, B, C, D, sw)
作業配列には private 属性を与えている.
音速から cfl 条件を計算してチェックする. 音速が負の場合は警告を表示する.
call cflchk
変数配列の大きさを割り付ける
call ArareAlloc
変数配列, および時間を初期化する. これらの配列および変数には全て 0 が 代入される.
call ArareArrayInit
擾乱の初期分布を与える.
call ArareDisturb
基本場を与える basicprm.f90 と同様, 将来的には初期分布はファイルから 読み込むようにすれば, このサブルーチンは必要なくなる.
境界条件を設定する.
call boundary(type, var)
境界条件のタイプは
C : 周期境界(ふらっくす格子点, ベクトル格子点共通) fS : フラックス格子点(ベクトル格子点)において対称 fA : フラックス格子点(ベクトル格子点)において反対称 sS : スカラー格子点において対称 sA : スカラー格子点において反対称
という記号を
(X 方向の境界条件)(Z 方向の境界条件)
と並べることによって表す. 例えば CfS ならX 方向は周期境界, Z 方向はフラックス格子点で対称な境界条件を意味する.
変数の格子点情報と境界条件の格子点情報が合致するように引数を与える. 例えば, 変数が ss_Exner の場合, 境界条件のタイプは CC, CsS, CsA などとなる.
速度の収束の値を計算する.
call VelDiv(fs_VelX, sf_VelZ, ss_VelDiv)
HE-VI 法で速度の水平成分を計算する.
call VelXHeVi(ss_Exner_B, ss_Div_B, fs_Fx_B, fs_VelX_B, fs_VelX_A)
HE-VI 法で速度の鉛直成分を計算する.
call VelZHeVi(ss_Exner_B, ss_Exner_A, ss_Div_B, sf_Fz_B, sf_VelZ_B, sf_VelZ_A)
HE-VI 法でエクスナー関数を計算する.
call ExnerHeVi(ss_VelDiv_B, sf_Fz_B, fs_VelX_A, sf_VelZ_B, ss_Exner_B, ss_Exner_A)
ソースコードは こちら よりダウンロードできる.
以下の環境での動作が確認されている.
事前に以下のソフトウェアをインストールしておく必要がある.
$ tar xvzf arare.tar.gz $ cd arare
詳細は実ファイルを参照. 必要なライブラリの格納ディレクトリの指定などを 変更する.
出力ファイル名, 空間解像度などの設定指定する NAMELIST ファイル arare.conf を編集する. NAMELIST 名はそれぞれが読み込まれるモジュールま たはサブルーチン名に対応している. 設定項目は以下の通り.
expname.f90 で読み込まれる. 実験名などを指定する.
physprm.f90 で読み込まれる. 惑星に固有の物理パラメータを読み込む.
fileset.f90 で読み込まれる. 出力される netCDF ファイル名を指定する. 格子点ごとに出力ファイルが異なることに注意.
gridset.f90 で読み込まれる. 格子点数, 計算領域の始点と終点, 境界条件を 指定する.
timeset.f90 で読み込まれる. 時間積分に関するパラメータを指定する.
linlib.f90 で読み込まれる. 行列計算ライブラリを指定する.
ソースディレクトリの直下で make する. $ make 実行ファイル arare が作成される.
$ ./arare
初期値を変えたい場合, arare.f90 の内部サブルーチン ((<ArareDisturb>)) 内を適宜編集する.
出力ファイルは netCDF ファイルなので, netCDF を解釈するソフトウェアを 用いることで結果の図示, 解析を行うことができる. たとえば dennou ruby 製品 の gave 等を使用するとよい.
Last Update : 2004/11/23 (北守 太一, 杉山 耕一朗, 小高 正嗣)