ここまではAthena++と共に配布されているパラメータファイルを使用してきました。ここでは解像度や出力ファイル、初期条件のパラメータなどを変更するためにパラメータファイルの編集方法について説明します。
前項で計算したOrszag-Tang Vortexで用いた解像度は決して高いものではありません。数値シミュレーションにおいて、興味ある現象・構造を分解するためには十分に高い分解能を設定する必要があります。そこでまずは分解能を二倍にして違いを見てみましょう。
1. パラメータファイルを編集する
パラメータファイルは単純なテキストファイルなので、任意のテキストエディタで編集することができます。編集する前に自分用のパラメータファイルを別名でコピーし、新しい計算用ディレクトリを作成することをお勧めします。
> cd /work/[username]/ > mkdir orszag_tang_high > cd orszag_tang_high > cp ~/athena/bin/athena . > cp ~/athena/inputs/mhd/athinput.orszag_tang ./athinput.orszag_tang_high
ではこのファイルをテキストエディタで開いてみましょう。ファイルは<...>というブロックで区切られていて、各項目は変数名 = 値(#以下はコメントとして無視される)という書式になっています。解像度は<mesh>ブロック内にありますので見てみましょう。
<mesh> nx1 = 200 # Number of zones in X1-direction x1min = -0.5 # minimum value of X1 x1max = 0.5 # maximum value of X1 ix1_bc = periodic # inner-X1 boundary flag ox1_bc = periodic # outer-X1 boundary flag nx2 = 200 # Number of zones in X2-direction x2min = -0.5 # minimum value of X2 x2max = 0.5 # maximum value of X2 ix2_bc = periodic # inner-X2 boundary flag ox2_bc = periodic # outer-X2 boundary flag
これはx1, x2方向ともに、0から1の範囲を200格子点で覆うという設定になっています。また境界条件として周期境界条件が指定されています。分解能を上げるためにはnx1, nx2を変更します。Athena++のインプットファイルはタブ区切りをサポートしていないので、空白は必ずホワイトスペースを使用してください。また、XC50を始めとするスパコンはLinux系のOSを使用していますので、改行コードはLFです。Windowsを使用している方はファイルを保存する際に注意してください(改行コードが何のことかわからない人は、まずは検索してください)。
<mesh> nx1 = 400 # Number of zones in X1-direction ... nx2 = 400 # Number of zones in X2-direction
2. 新しいファイルで再度計算を実行する
既にOrszag-Tang Vortex用にコードがコンパイルされていれば、コードを再度コンフィギュレーション・コンパイルする必要はありません。そのまま新しいファイルを使って計算を実行しましょう。実行すると以下のように出力されるはずです。
RootGrid = 1 x 1 x 1 MeshBlock 0, rank = 0, lx1 = 0, lx2 = 0, lx3 = 0, level = 0 is=2 ie=401 x1min=0 x1max=1 js=2 je=401 x2min=0 x2max=1 ks=0 ke=0 x3min=-0.5 x3max=0.5 Setup complete, entering main loop... cycle=0 time=0.00000000000000e+00 dt=4.51495485504003e-04 ... cycle=2774 time=1.00000000000000e+00 dt=3.40143584065472e-04 Terminating on time limit time=1.00000000000000e+00 cycle=2774 tlim=1.00000000000000e+00 nlim=-1 cpu time used = 8.07699951171875e+02 zone-cycles/cpu_second = 5.49511000000000e+05
ここで、時間刻みdtと、計算にかかった時間cpu time usedに注目してください。時間刻みがおよそ半分になっているのは、Athena++で採用している陽解法では時間刻みがCFL条件(Δt < Δx/c, cは情報が伝わる最大速度)で制限されているためです。分解能を2倍にしたことで1ステップ当たりの計算コストは4倍になる上に、必要な時間ステップが2倍になるため、合計で約8倍計算コストがかかることになります(3次元なら16倍!)。このように高解像度計算には大きなコストが必要になりますので、自分が解こうと思っている問題に応じて適切に分解能を設定する(あるいは自分が使える計算資源に応じて適切に「解ける」問題を設定する)ことが重要になります。
3. 可視化と比較
計算したファイルを前回と同様に可視化して結果を比較してみましょう。
解像度を上げたことによって細かい構造がより正確に表現されていることが確認できると思います。時間的に余裕があれば1次元の衝撃波管問題でも同様に解像度に対してどのように解が振る舞うか試してみましょう。
<output2> file_type = vtk # VTK data dump variable = prim # variables to be output id = prim # file ID dt = 0.01 # time increment between outputs <output3> file_type = vtk # VTK data dump variable = cons # variables to be output id = cons # file ID dt = 0.01 # time increment between outputs
これで計算を実行すると"OrszagTang.block0.prim.?????.vtk"と"OrszagTang.block0.cons.?????.vtk"という名前のファイルが出力されます。これらを可視化してみましょう。見たい物理量に合わせてこれらを使い分けることが重要です。また、ファイルサイズを節約するために特定の変数だけ出力することや、特定の平面でのスライスのみ出力することも可能です。詳しくはマニュアルを参照してください。