熱力学計算コード oboro



|back to doc| |oboro とは| |平衡組成の計算| |エントロピー計算|
|dS=0 の探索| |大気科学的物理量| |oboro の構成|
[English]



oboro とは

oboro は上昇する気塊のとりうる対流平衡状態をギブス自由エネルギー最小化法を用いて計算する:

  1. 温度圧力を固定し, 元素数保存の条件のもとギブス自由エネルギーが最小化となる平衡組成を求める.

  2. 当該温度圧力における平衡組成を用いてエントロピーを計算する.

  3. 断熱線 dS = 0 をエントロピーの保存される温度圧力を探索することで求める.

  4. 断熱線 dS = 0 に沿った種々の大気科学的な物理量の計算
本計算では大気の擬断熱過程を仮定している. この過程では凝縮物は気塊から取り除かれるが, 取り除かれた凝縮物も含めた全エントロピーが保存する.



平衡組成の計算

平衡組成はギブス自由エネルギー最小化法を用いる. これは温度・圧力, 元素組成を与えた際の熱平衡状態を, ギブス自由エネルギーの最小化された状態として求める方法である. 数値計算的には, RAND 法という最適化法を利用する.

最適化法(Newton 法)のイメージ図. 適当な組成のまわりで 2 次関数に近似しその極値を反復的に求めることで, ポテンシャルの最も低い状態へ収束させる.

最適化法を用いる際に問題となるのは, 計算の初期値として与える初期組成である. この組成が平衡組成に十分に近くないと, ギブス自由エネルギーのグローバルな最小値を与える組成に収束せず, 局所的な最小値を与える組成に収束してしまう場合がある. 特に 2 相の境界(例えば水と氷)あたりでは, 2 相のギブス自由エネルギーがほとんど変わらないため, うまい初期組成から計算を開始しないと平衡状態に収束しないことがよくある.

適当な初期組成から計算を開始しないと, global minimum ではなく local minimum に収束することがある.


エントロピーの計算

平衡組成から解析的にエントロピーを計算することができる.

この式を用いることで, 任意の温度・圧力・平衡組成でのエントロピーが得られる.



dS = 0 の探索

擬湿潤断熱線に沿った温度, 圧力, 組成を以下の手順で求める. 通常, 観測から大気上端の温度圧力が与えられるので, 以下の手順で得られた断熱線のうち観測された温度圧力を通るものを探すこととなる.

断熱曲線の計算方法.

Step 1: 初期温度での平衡組成を計算し, エントロピー $S_0$ を求める(左上).

Step 2: 圧力を $dp$ だけ変化させる. 温度を変化させた時のエントロピーを 順次計算し, 前のステップでのエントロピー $S_0$ と一致する温度を求める (右上).

Step 3: ($T_{1}, p_{1}$) において凝縮が生じた場合, 次のステップにおいて保存 させるエントロピーは ($T_{1}, p_{1}$) での気相のエントロピー ${S_{4}}_{\rm gas}$ (左下).

Step 4: Step 1 から Step 3 において得られた温度・圧力を結んで 断熱線を引く(右下).



上図の Step 2) において, 温度を変化させながら前のステップでのエントロピーに収束させる必要がある. エントロピーを収束させる手法として, 以下のような「挟み込み法」を用いる (前のステップでのエントロピーを S_{pre} とする).

Step1:
適当な温度幅 dT の間隔でエントロピーを計算し, S_{pre} をはさむ温度を探す.
Step2:
ステップ 1 で得られた 2 点の温度の中点でのエントロピー S_{now} を反復的に計算し, S_{pre} と比較する. S_{now} と S_{pre} の変化率がマシンイプシロン程度になった場合には, S_{now} と S_{pre} は一致したものとみなす.


大気科学的な物理量の計算

断熱線として得られた鉛直 1 次元の温度, 平衡組成を用いて, 大気運動特性の支配パラメタの一つである「静的安定度」および「雲密度」を計算する.

静的安定度は地球と同様に積雲対流によって決まると仮定し, 静的安定度を湿潤断熱減率と乾燥断熱減率との差から評価する. 地球では積雲対流によって基本場の温度構造はほぼ湿潤断熱的な構造となる. なぜならば下降流の温度分布は放射冷却のために乾燥断熱減率に従った変化をせず, 上昇流域での湿潤断熱的な温度分布にほぼ一致するからである. そのような温度構造を持つ大気での静的安定度は, 乾燥気塊の乾燥断熱減率と基本場の湿潤断熱減率との差から決まる.

Air parcel is radiatively cooled while it descends. Temperature profile is roughly given as moist adiabat. Static stability is given as the difference between moist adiabat and dry adiabat



oboro の構成

oboro の構成は以下の通りである. モデルの I/O は gtool4 netCDF 規約に則った netCDF ファイルを利用している. また, 各々のメインルーチンで呼び出すサブルーチン群は src ディレクトリ以下に格納されている.

ObrAdbt
本熱力学計算コードの核となるメインルーチン. DataFile.dat に記載された系の設定と, ギブス自由エネルギー最小化法とを用いて, エントロピー一定の断熱状態を求める. 出力ファイルは gtool4 netCDF 形式となっている.

ObrMap
DataFile.dat に記載された系の設定と, ギブス自由エネルギー最小化法とを用いて, 温度圧力空間での平衡組成のマップを作成するメインルーチン.

ObrAtmos
雲密度, 静的安定度といった大気科学的な物理量を計算するメインルーチン. ObrAdbt の出力ファイルを読み込む必要がある. 出力ファイルは gtool4 netCDF 形式となっている.

ObrCheck
DataFile.dat の内容を基に, 物性値を計算し表示するためのメインルーチン. このメインルーチンの表示を見ることで, 物性値が正しく計算されているかがチェックできる.

ObrDebug
エントロピーの収束状況をチェックするためのメインルーチン. ObrAdbt の出力ファイルから必要な情報を取り出し, その圧力においてエントロピーの収束状況を逐次プロットする. エントロピーの収束が悪い圧力において何が起きているかをチェックするのに用いる.

ObrSvap
ObrAdbt の出力を用いて蒸気圧と飽和蒸気圧を同時にプロットするためのメインルーチン. 蒸気圧が飽和蒸気圧と一致するか否かで平衡状態にあるかチェックする (溶液の場合は多少飽和蒸気圧からずれることもある).

ObrData
計算に必要な物性値の書かれた netCDF ファイルを作成するためのメインルーチン. 質問にしたがって基準状態でのエントロピーとエンタルピー, 比熱, 飽和蒸気圧等を入力すると, それらの書かれた netCDF ファイルを生成する.