トップ > ライブラリ > 技術計算 > 流体解析

技術計算：流体解析

現在工事中です。

【注】数式表示にはMathJaxを利用しています。IE8以下では表示が遅くなる可能性があります。FireFox などIE8以外のブラウザを利用下さい。

流体解析概論

　流体解析（CFD：Computer Fluid Dynamics）は、流体の運動を計算機上でシミュレーションし、流体挙動を明らかにすることにより 流路形状の最適化、条件の最適化などを費用をかけることなく検討するために行われる。空洞実験など小スケールな実物を用いた試験では、ものの製作、実験の実施、 データの解析など時間がかかり、また労力・経費も必要とする。計算機上でモデル化し、機上検討することにより、大幅に時間と費用を低減することにメリットがある。

基礎方程式

　流体解析を実施するということは、少なくとも運動方程式（またはNavier-Stokes式）を解き、流体の速度場を求める解析を行うことをいいます。流体の速度を一様とみなし、 温度場、濃度場などを解く反応器シミュレーターとは区別される。反応器シミュレーションに関しては技術計算：反応器シミュレーションに詳述しています。
　ここで運動方程式とは運動量保存式であり、独立変数として少なくとも一次元以上の空間、時間をもち、従属変数として、 空間に沿った方向の流体速度をもつ。\((x,y,z)\)直交座標三次元空間の場合には、それぞれの方向の速度成分\((u,v,w)\)が従属変数となる。非定常の運動量保存則を直交座標系で表すと 次の(1)～(3)式となる。ここで、\(\rho\)は流体密度、\(\mu\)は流体粘度、\(g\)は重力加速度を示し、密度は解析領域内で一定とした。すなわち非圧縮性（液に近い）流体とした。

【注】圧縮性流体（たとえば気体）では密度は時間・場所により変動する。したがって密度を未知数とすれば、気体の状態方程式が支配方程式として現れる。

\[ \begin{align*} \rho \frac{\partial u}{\partial t}+ \rho \biggl( u \frac{\partial u}{\partial x} + v \frac{\partial u}{\partial y} + w \frac{\partial u}{\partial z} \biggr) & = -\frac{\partial p}{\partial x} + \mu \biggl( \frac{\partial u^2}{\partial x^2} + \frac{\partial u^2}{\partial y^2} + \frac{\partial u^2}{\partial z^2} \biggr) + g_x \rho \tag{1} \\ \rho \frac{\partial v}{\partial t}+ \rho \biggl( u \frac{\partial v}{\partial x} + v \frac{\partial v}{\partial y} + w \frac{\partial v}{\partial z} \biggr) & = -\frac{\partial p}{\partial y} + \mu \biggl( \frac{\partial v^2}{\partial x^2} + \frac{\partial v^2}{\partial y^2} + \frac{\partial v^2}{\partial z^2} \biggr) + g_y \rho \tag{2} \\ \rho \frac{\partial w}{\partial t}+ \rho \biggl( u \frac{\partial w}{\partial x} + v \frac{\partial w}{\partial y} + w \frac{\partial w}{\partial z} \biggr) & = -\frac{\partial p}{\partial z} + \mu \biggl( \frac{\partial w^2}{\partial x^2} + \frac{\partial w^2}{\partial y^2} + \frac{\partial w^2}{\partial z^2} \biggr) + g_z \rho \tag{3} \\ \end{align*} \]

　流体解析では速度場のみを解くことはなく、通常圧力場を同時に解く。圧力場は連続の式（または質量保存式という）を連立して解く。連続の式は次の (4)式となる。

\[ \begin{align*} \frac{\partial \rho}{\partial t}+ + \frac{\partial \rho u}{\partial x} + \frac{\partial \rho v}{\partial y} + \frac{\partial \rho w}{\partial z}=0 \tag{4} \end{align*} \]

　温度場は熱収支式を、濃度場は成分物質収支式を、同時に連立して解くことにより、温度・濃度の空間分布を解析結果として得ることができる。式の表示は省略する。 　気象関連の数値予報では、水蒸気の保存式を解き、空気中の水分の変化・降水量などを予測する。これは空気中の水分濃度に相当し、成分（水）の保存式を解くことに相当する。 　もっと複雑な二相系（気液、気固、液固）、三相系（気液固）のときは、それぞれの相で保存式となる支配方程式を構築することになる。相間で運動量移動、物質移動、熱移動が ある場合にはこれら移動速度の関係式をあらかじめ求め、モデル化しておかないと解くことができない。こうした複雑系は専門書に譲るとして、ここでは単相系について述べる。

　速度場、圧力場、温度場および濃度場を支配する基礎方程式（支配方程式ともいう）を、連立して解くことを流体解析と称する。 基礎方程式は連立の偏微分方程式を構成しており、解くためには必要な境界条件、初期条件を必要とする。

　上でも述べたように、流体解析とは少なくとも運動量保存則と連続の式とから速度・圧力場を解くことをいう。温度場、濃度場の 解析・計算は必須ではない。

市販ソフト

　2000年代に入ると、流体解析ソフトとして、数多くの優れたソフトが市販されるようになってきている。また数多くのオプションを備え、 いろいろなケースの問題が解けるようになってきている。 　たとえばFluent、Star CD、Phenicsなど世界的に有名なソフトがすでにあって、さまざまな工学的問題 に適用されている。こうした流体解析ソフトでは、前処理の格子形成、流体力学による計算、後処理（結果の表示）まで一連のソフト群 から構成されている。汎用の流体解析ソフトは、前処理、後処理を含め計算手法などの多くのオプションが用意されていて、初心者が使いこなすにはかなり期間のトレーニング を受ける必要がある。また当然のことながら、流体解析に関する基礎知識をもっていないと、誤った結果を導くことにもなる。

　筆者自身はこうした汎用の流体解析ソフトを用いた解析経験はないものの、三次元流体解析コードをゼロから開発し、流体場の解析を 実施した経験がある。有限要素法を用いて、高粘性流体からの揮発性物質の蒸発解析を実施した。また有限差分法を用いて、撹拌槽 の三次元流体解析、特殊円筒型気相反応器の三次元流体解析をコード開発から経験している。 したがって、市販の流体解析ソフトなどの結果を評価する上で、その得失を十分理解しており、 正当な評価・判断を下すことができると思っている。市販ツールには用意されてい ない機能を実現するため、コードを独自に開発することが有利な場合もある。

実行環境

　流体解析は、以前（20年くらい前）は大学の計算機センターにあるスーパーコンピュータと呼ばれる大型ベクトル計算機を用いて、高速演算を 行っていた。最近でも「地球シミュレーター」とか理化学研究所の「京」とかのマルチCPU計算機で行っているが、Intel XEONやマルチコアCPUのワークステーションや汎用パソコンでも 流体解析ソフトが稼働できる環境にある。最近では、2012年頃にIntel社のXEON PHIは60コアのCPUをもつ、ちょっと高級なパソコンで、マルチCPUの並列計算が 利用可能となっている。IntelのマルチコアCPU、Core i7やCore i5などでもマルチスレッドなど並列計算が可能な 計算環境をもっており、ひと世代前のスパコンの性能を、計算速度とメモリ容量の点において、凌駕するパソコンが現時点では利用できている。

　また解析結果を出力する方法も、20年前はプロッタによる二次元表示が主流であり、画像を表示する ディスプレイもCRTを使った表示で低精細の出力画面であった。しかし、現在では液晶ディスプレイに移り変わり、高精細となり、 表示用の描画ソフトもOpen-GLに代表される美しい出力が可能になっている。

　Open-GLのサンプルについては、計算機：データ処理・グラフィックスを参照されたい。

解析手法

　通常、運動方程式は非線形項があり、解析的に解くことができず、数値計算により、反復計算をして解く必要がある。 流体解析で支配方程式である、連立の偏微分方程式を解く主な方法として、次の解法がある

　ここでは、これら解法の得失を簡単に述べる。

　差分法は、解析領域を格子に分割し、各格子点の座標値、変数値を用いて、偏微分方程式を差分化する。これを離散化するという。 運動方程式には移流項が含まれ、非線形項となっている。差分化方程式は、オイラー法で代表されるような陽解法により、格子点を順番に 解くことを繰り返し、全系が式を満足するまで反復計算を行う方法がある。また陰的オイラー法（Backward Euler）に代表される陰解法があり、 差分化方程式を部分的に連立方程式と見なし、行列演算で解ベクトルをえ、これを全系が式を満足するまで反復する方法がある。 　差分法は、それほど計算機のメモリを消費することはなく、パソコンやワークステーションの計算機で解くことが可能である。 　したがって、多くの市販流体解析ソフトはこの差分法を採用している。

　有限要素法は、解析領域を有限個の要素（エレメント）に分割し、すべての要素の節点（格子点）の変数を未知数とする、巨大な正方行列を もつ連立一次方程式とし、これを解くことで格子点の変数値を得る。そのため、格子点の数や変数の数が大きくなると、巨大な正方行列を 計算機メモリ内に保存しておく必要があり、メモリを大量に必要とする。

　境界要素法は、解析領域の境界だけを要素に分割し、積分定理を利用し、境界要素の節点だけを変数とすることで、行列のサイズを小さくできる。 なお詳細は、コード作成の経験がないので専門書を参照されたい。

流体解析の一般的な作業手順

　ここでは、三次元の流体解析、特に非圧縮性流体で速度・圧力場を解く場合を例に取り上げ、一般的な作業の流れを簡単に紹介する。また、解法としては 有限差分法と有限要素法を対象として話を進める。

解析モデルの定式化

　流体解析を実施するに当たって、まず何のために実施するのかの目的を明確化する必要がある。 その目的を達成するためには、必要十分な支配方程式を構築するとともに、解析領域（時間・空間）を設定する必要がある。 同時に、解析領域の初期条件・境界条件を選定する必要がある。
　解析モデルとなる支配方程式を定式化する。多くは保存則を表す時空間を独立変数とする偏微分方程式となる。 このとき解析対象の三次元座標系を直交座標系、円筒座標系、極座標系かを選定する。
　解析モデルを構築するとき、解に大きな影響を及ぼす項を含める必要がある。また大きく影響しない項は省略するなどを行い、 支配方程式をできるだけ単純化しておく。これにより無駄な計算を省くことができる。 　一般座標系を用いて、保存則をベクトル表示すれば、上記の任意の座標系に、容易に変換適応させることができる。 一般座標系についてはここでは省略するが、別途解説する予定である。

　解析目的を満たすようモデルを構築するとき、熱収支、成分収支を考慮する必要があるのか。 また流れが乱流となっていたら、乱流モデルを考慮する必要がある。 単一相流れでなく、二相流・混相流となるのかなどで解析すべきモデルの支配方程式が増減する。
　運動量保存則、熱保存則などの支配方程式を見ると、非定常項、移流項、拡散項、重力項、摩擦項、 物理・化学的な変化項などのうち、どの項が支配的であるかの目安を立てておかないと、流体解析結果として全く予想外の結果になることがある。
　思わぬ項を無視したため、実際の現象をシミュレートできなかったという結果になりかねない。

前処理

　有限差分法・有限要素法どちらにしても、解析領域をメッシュに刻む。有限差分法では格子点（Grid Point）と呼ぶが、格子点に変数を定義する。 一方、有限要素法では、メッシュに刻んだ最小単位を要素（三角、矩形、立方体、ヘキサゴナルなど）と呼び、その頂点を節点と呼ぶ。この節点に変数を定義する。
　これら格子点、節点（以下格子点と呼ぶ）に割り振られた従属変数を未知数とし、連立の非線形方程式を解くことになる。
　格子点に変数をもつのでなく、三次元の最小単位体積に代表値をもたせる解析方法（有限体積法とよぶ）もある。これも差分法の変法である。

　格子の種類として、境界に合わせたメッシュ分割を行う境界適合格子を通常は採用する。 格子間隔は等間隔、非等間隔どちらでも構わないが、隣り合う格子幅は滑らかに変化するように刻む。 また、変数が急変する部分の格子幅を狭くすることを行う。
　格子点の総数は、解くべき支配方程式の未知数であり、利用できる計算機資源（メモリ、計算速度）に依存し、 ある程度有限の値に抑える必要がある。

　格子形成法には、直交格子を直接指定する方法や、ラプラス方程式を解くことで格子を解析的に得る方法が ある。単純な解析領域であれば直接指定で格子を作った方が理解しやすい。

数値解析

　有限差分法では、支配方程式と境界条件を含め、離散化し差分方程式に変換する。大抵は、非線形の連立方程式となっており、解析的に解くことはできず、 数値解析により、解くことが行われる。一方、有限要素法では、全変数を未知数とした、巨大な連立一次方程式を作り、これを解くことで解を得る。 　どちらも、連立一次方程式を繰り返し解くこととなり、膨大な計算量を費やすことになる。 そのため、スーパーコンピュータやマルチコアプロセッサによる並列計算を行うことになる。 　もちろん、時間がかかるけれどもシングルコアのパソコンでも解くことができる。 　一昔前のスーパーコンピュータはベクトル演算器を備え、ベクトル計算を瞬時に行うようコード開発 されていたが、昨今ではマルチCPUをもつスーパーコンピュータで、並列計算を行えるよう、 連立一次方程式の高速解法が開発されている。 　マルチコアのパソコンでも、MPIとかOpenMPといった並列計算用の規格があり、並列化による 計算が身近なものとなっている。 　非定常の計算では、Time Marching法により、時間ステップ毎の計算を行い、定常時の解や途中の解を得ることをしている。 　差分法の陽解法では、時間刻みΔtと空間刻みΔxとの間にCFL条件を満たすような制限がある。 　また、拡散項に比べ、移流項が支配的になると解が不安定化する。これを抑制するため風上差分というテクニックを利用する。

後処理

　流体解析の結果は、格子点値（Grid Point Value）としてとびとびの値が得られる。これを値として出力しても 理解することは難しい。そのため見やすくするため、昔はプロッタによる二次元的なワイヤーフレーム的な表示出力 を行っていた。またはディスプレイ上で同じく二次元的な出力表示を行っていた。 　最近では、OpenGLなどのグラフィック出力など美麗なディスプレー表示が可能となっている。

有限差分法と有限要素法

境界要素法、有限体積法

乱流と乱流モデル

は、ここで取り扱っています。

Literature Cited

先頭に戻る