Feel＆Think

今回は、地下水流れによって物質が運搬されながら拡散する問題に、陽解法を適用してみます。まずは、物質移行の支配方程式を示します。

ここに、cは間隙水中の物質の濃度、nは空隙率です。また、v^cは分散と拡散による物質の移動速度（流束）であり、v^wは間隙中の水の速度です。ただし、濃度変化に伴う間隙水の密度変化は十分小さいことを前提とします。詳細は、「物質移行に関する支配方程式と陽解法」をご覧ください。

上式で、右辺第1項は分散と拡散により単位時間に単位体積中で生じた物質の出入り量の差を表しており、第2項は間隙水により運ばれることで同様に生じた物質の出入り量の差を表しています。また、左辺はこの物質量の変化が濃度の変化をもたらすことを表しています。なお、間隙水の平均流速vとv^wには、次の関係があります。

分散と拡散による物質の流束v^cは、Fickの法則に従うものとします。

Dは拡散係数であり、次のような形で記述されます。

ここに、α_Tとα_Lは横分散長（間隙内流速方向と垂直方向の分散長）と縦分散長（間隙内流速方向の分散長）、D_mは分子拡散係数、τは屈曲率、δijはクロネッカのデルタです。また、は間隙内流速ベクトルの長さであり、2次元であれば となります。

陽解法とするために、地下水流れの場合と同様に、まず式(3)にガラーキン法を適用し離散化します。このとき、一次の内挿関数を用い、要素重心で濃度が定義され物質の流束は節点で求まるものとします。

節点値について要素単位での合算が終わった状態では、次のような式が得られます。

ここに、{C_e^t}は節点における等価濃度勾配（重みの付いた濃度勾配）、{C_s^t}は境界面の濃度と等価な節点濃度勾配、Wは重みです。また、tは現在の時刻を表しており、後に示す時間差分では既知量として取り扱うことを示しています。

左辺第1項の係数[W]は対角マトリクスであることから、節点ごとの割り算によりそれぞれの節点における物質の流束を求めることができます。

式(1)において、移流項である右辺第2項には、

として、同様な手続きによりそれぞれの節点で求めた濃度勾配を用います。

式(1)は、要素の物質の流束を節点値より内挿して求め、時間微分は差分で近似することとします。ただし、地下水流れの解析より間隙内流速が得られていることを前提とします。

[N]は内挿関数からなるマトリクス、[∂]は空間に関する偏導関数マトリクスです。あるいは、勾配マトリクス[B]を用いると次式となります。

間隙内の流速は次式で求められます。

なお、間隙水の平均流速vtが既知の間隙水圧φ^tより求めることができることは、前回示したとおりです。

ここに、F_e^tは間隙水圧の勾配から求まる等価節点流速、F_s^tは境界面の水圧と等価な節点流速です。
したがって、次のようにして既知の値から次の時間の濃度を求めることができます。

これまで説明した陽解法による、物質移行の解析例をご紹介します。この解析モデルでは、左から右へ地下水が流れており、流速などは事前に地下水流れの解析で求めてあります。ここでは、濃度変化による地下水の密度変化が十分小さいと仮定しており、物質移行は地下水流れに影響を及ぼしません。したがって、地下水流れを別個に解いて、流速を求めておくことができます。

解析結果を動画-1に示します。濃度変化が地下水流れの方向に生じていくことがわかります。最終場面は、解析開始から2日後の状態です。分散がなければ、濃度が1に近い領域が流れの方向に170mほど拡大していくはずですが、解析結果も拡散を伴いながらほぼ同様の距離を移行していることがわかります。

同じモデルで縦分散長を100m、横分散長を10mとして計算した結果を、動画-2に示します。先の結果に比べて、分散が強くなっていることがわかります。

このように、陽解法によって物質移行の解析も可能であることがわかりました。次回は、同じ物質移行ですが濃度の変化により地下水の密度が変化する、いわゆる密度流の解析をご紹介します。