大気シミュレーション

目標

大気の流れを支配する法則の数値シミュレーションをやってみる。基礎となる法則

流体の運動：　Navie-Stokes方程式、連続の式
熱力学：　気体の状態方程式、熱伝導方程式

を離散化した差分方程式を計算機上で解いて、流体の運動を調べる。

以下の課題を一つずつ確認していく。

偏微分方程式を数値的に解いて、既知の問題に対して妥当な結果を得られることを確認する。
大気の物理量の可視化、流体運動の特徴の抽出。
特徴ある流体運動がおこる条件を調べる。

自前のプログラム

上の熱力学関係式と流体の運動方程式を差分方程式として実装したものを作ってみた。コードはhttps://github.com/tkono17/Airdog.gitにおいてあるので、そこからダウンロードして使う。

まず、作業用のディレクトリを作って以下の作業は全てそのディレクトリの下で行う。ここでは、作業ディレクトリをFluidという名前のディレクトリとして、以下の手順で示すcd Fluid は実際には、Fluidの絶対パスを指定する必要がある場合があるので自分の環境に応じて適切なコマンドを実行する。

mkdir Fluid;

Gitレポジトリからコードをダウンロードする

cd Fluid
git clone https://github.com/tkono17/Airdog.git

いろいろな作業を行うためにディレクトリを用意する

cd Fluid
mkdir build sw scripts run analysis

簡単に各ディレクトリの役割をまとめると、こんな感じである。

ディレクトリ	使用用途
Airdog	プログラムのソースコード
build	プログラムをコンパイルするための作業用
sw	コンパイルしたプログラムを保管するため
scripts	いろいろなスクリプト（環境設定等）
run	プログラムを実行するところ
analysis	生成したデータを解析するところ

環境設定

今後は上で示したディレクトリ構造のもとでプログラムの開発やデータ解析を行っていきたい。そのために、簡単な環境設定をしておくと便利である。 scripts/setup.shというディレクトリを作り、以下の内容を記述する。

export SWDIR=Fluid/sw
export PATH=${SWDIR}/bin:${PATH}
export LD_LIBRARY_PATH=${SWDIR}/lib:${LD_LIBRARY_PATH}

Airdogプログラムをコンパイルするためには、環境変数SWDIRでインストール先ディレクトリを設定するようになっている。このファイルを用意できたら、ターミナルを開いて作業を開始する前に毎回Fluidディレクトリから以下のコマンドを実行する。

source ./scripts/setup.sh

プログラムのコンパイル

cd Fluid/build
cmake ../Airdog
make
make install

makeを実行して100%処理が終了していれば、コンパイル作業は正常に終了したことになる。プログラムのソースコードを更新した場合は再度コンパイルが必要なため、cmake, makeをやり直す必要がある。 make installを実行すると、コンパイル済みのファイルがFluid/swにコピーされて利用可能になる。

プログラムの内容

有限の大きさの空間に格子を設定して、各格子点における空気の物理量を保持する。全ての格子点における物理量が系の状態を決定する。プログラムでは基礎方程式の時間発展を解いて、系の状態の時間発展を追跡する。データ量が膨大になるため、何回かの時間ステップを経過した後での状態の2次元断面や2次元面への射影をした結果として得られる2次元分布のみを記録するようにしている。

状態変数：　密度、圧力、温度、速度 (u, v, w)

境界条件

OpenFOAMプログラム

無料の流体解析プログラムである。 OpenFoamExample

ROOTによる解析

ここでは、ROOTによる解析のためのディレクトリ(analysis/)を用意して作業することにする。このディレクトリがまだ存在しない場合には、analysis.tar.gzからanalysisディレクトリを圧縮したファイルをダウンロードしてから以下のコマンドを実行する。

tar xvfz analysis.tar.gz

これを実行した後では、作業ディレクトリの下にanalysis/ができているはずである。

ls -l

でディレクトリが存在していることを確認する。作業はanalysis/で行うので、以下のコマンドでそこに移動する。

cd analysis

ROOTはデータ解析に有用な機能をC++のライブラリとして用意しているものなので、ROOTを用いた解析では基本的にC++のプログラムを書くようなものである。ただし、コンパイル可能なちゃんとしたプログラムを書かずに、ソースコードを逐次的に解釈して実行してくれる機能（インタープリタ）が用意されているので、プログラムの一部（通常、C++の関数）のみ書いて実行することができる。このようにインタープリタで実行するコードの一部をマクロと呼んでいる。マクロには実行すべき関数を記述すればよい。例えば以下のような内容：

#include <iostream>

void helloWorld() {
  std::cout << "Hello, World!" << std::endl;
}

マクロのファイル名は関数名と一致している必要がある。そのため、このファイルはhelloWorld.Cというファイルに記述されていないと正しく実行されない。マクロを実行するにはrootコマンドでROOTを立ち上げてから

root [0] .x helloWorld.C

".x"はROOT内部でのみ有効なコマンドで、その後に指定するファイル名で保存されたマクロを実行するという指示である。この例ではマクロファイルは全てmacros/というディレクトリの下に保存されている。

マクロ例

マクロ名	内容
mk_vortex.C	渦が作る速度場を生成してヒストグラムに保存する。work/vortex_ex1.root
pl_vec.C	速度場を矢印で可視化するためのマクロ。ついでに各点での速度の大きさをヒストグラムで重ねて表示している。
pl_vorticity.C	pl_vec.Cと同様だが、速度の大きさの代わりに渦度の大きさを表示する。

vortex_ex1.rootの中には、2種類の渦のデータがある。

1つ目は単一の渦がある設定。
2つ目は渦が2つあり、それぞれが逆回転している設定である。

それぞれのデータは、速度場を2次元ヒストグラムとして保存しているので、vx, vy用に2つのヒストグラムがある。（h1x, h1y; h2x, h2y）この例では、渦の中心から距離1以下では速度場はゼロで、中心からの距離1のところで速度50、そこからさらに離れるにつれて1/sqrt(r)で速度が小さくなっていくような状況でデータを生成している。

> root
root [0] .x mk_vortex.C
root [1] .x pl_vec.C

とやれば渦データを生成して、その後速度場をプロットさせることが可能である。結果はfigures/vec.pngに画像ファイルとして保存される。一度データを作った後は、mk_vortex.Cは実行する必要はない。結果はこんな感じになる。

AirSimulation