まずは

できそうなことからやってみます.2020/04~

数値流体【導入】

もうそろそろ本題の流体の数値計算に触れていきたいと思います.
備忘録的な記事になるので,見苦しさもあるかもしれませんが間違えとかあったらご指摘ください.
しばらくは1次元 or 2次元の河川の水理計算(流れと河床変動を合わせた計算)を目指します.

流体の数値シミュレーションの流れとしては,

   ①物理モデルの選択 → ②離散化 → ③数値計算 → ④可視化

となります.

物理モデルの選択では対象とする流れから,NS eq.をどう扱うのか(そのまま or 条件を加えて近似),流れの基礎式となる偏微分方程式を決定します.

離散化では①の偏微分方程式離散化する手法を選択し,計算格子を設定することで,実際に解くべき代数方程式を導きます.

数値計算では②に導いた代数方程式を解く数値計算法を決めて(反復法や消去法など),プログラミング言語を用いて記述,計算結果を出力します.

可視化では③で出力された数値データをグラフやアニメを使って,可視化します.

つまり,流体力学はもちろん数値計算方法やプログラミングなど多岐にわたる勉強が必要になります... 基本的には,離散化手法に差分法,プログラミング言語fortranを使う予定です.

今回は水理学とのつながりを意識しつつ,数値流体(流体力学)の導入を示していきたいと思います.

流体の基礎式

基礎式は主に運動方程式(Equation of motion:EOM)と質量保存則(Law of conservation of mass)になります.流体の場合,運動方程式ナビエ-ストークス方程式(Navier-Stokes equation:NS eq.),質量保存則は連続の式(Continuity equation)と呼ばれます.

ナビエ-ストークス方程式

表面力として垂直応力とせん断応力が作用していると考えた流体微小要素にNewtonの第二法則  {\bf F}=m{\bf a} を適用することで誘導されます.次式が非圧縮流れにおける基礎式になります.
なお,左辺は加速度をラグランジェ的観測(空間座標が時間に対して従属関係)ではなくオイラー的観測(時間と空間が独立)で考え,テイラー展開を用いて近似することで示される流体の加速度表示になります.

 
\displaystyle \frac{\partial{\boldsymbol v}}{\partial t} + ({\boldsymbol v} \bullet {\nabla}) {\boldsymbol v} = {\boldsymbol K} - \frac{1}{\rho}{\nabla p} + {\nu}{\nabla}^2{\boldsymbol v}

ここに, {\boldsymbol v} は流速ベクトル,tは時間, {\boldsymbol K}は外力ベクトル, {\rho}は流体の密度,pは圧力, {\nu}は動粘性係数です.
各項は一般に次のように呼ばれ,対象とする流れによってそれぞれ扱いが変わります.

左辺

  • 第1項:局所加速度項(Variation term)
     時間的な(局所的な)速度の変化
  • 第2項:移流加速度項(移流項Advection termや対流項Convection termともいう)
     空間的な速度の変化,難しさの原因である非線形項...

右辺

  • 第1項:質量力項(Mass force term,外力項ともいう)
  • 第2項:圧力項(Pressure term)
  • 第3項:拡散項(Diffusion term)

ちなみに,左辺は実質微分の記号を用いて次式で表すことができます.

 
\displaystyle \frac{ \mathrm{D}{\boldsymbol v} }{ \mathrm{D}t } \equiv \frac{\partial{\boldsymbol v}}{\partial t} + ({\boldsymbol v} \bullet {\nabla}) {\boldsymbol v}

また,数値解析を行う上では無次元化することでスケールを考慮しなくて済むため,それぞれ代表物理量を用いて無次元化すると以下のように示されます.無次元量には上付きの*を付けています.

 
\displaystyle \frac{ \mathrm{D}{\boldsymbol v^*} }{ \mathrm{D}t^* } = {\boldsymbol K^*} - {\nabla^* p^*} + \frac{1}{R_e}{\nabla^*}^2{\boldsymbol v^*}

ここに, {R_e}レイノルズ数であり,外力項が無視できる場合,流れはレイノルズ数のみで決定するため流体力学水理学において重要な無次元パラメータになります.
次項で少し掘り下げます.

無次元パラメータ

レイノルズ数 {R_e},Reynolds Numberは慣性力と粘性力の比を表しています.つまり,この値が大きいほど粘性のない理想流体に近づきます. Reynoldsのpipe flowでの層流乱流の臨界レイノルズ数に関する実験はどこかで見たことがあるんじゃないかと思います.


{R_e} = \displaystyle \frac{VD}{\nu}

Vは代表速度,Dは代表長さです.
例えば,矩形断面開水路であればVを断面平均流速,Dを水深とすると,VD=q :単位幅流量より単位幅流量qと動粘性係数 {\nu}の比で求めることができます.

また,その他無次元パラメータとしては常流射流を決定するフルード数 {F_r}(Froude Number,代表流速と長波の伝播速度の比)や水理では水を扱うためほぼ考慮することはありませんが,圧縮流れを分類するマッハ数 MMach Number,代表流速と音速の比)があります.

フルード数  {F_r}重力場において(外力ベクトル  {\boldsymbol K} が重力ベクトル  {\boldsymbol g} であるとき),基礎式を無次元化すると出てくる項になります.


{F_r} = \displaystyle \frac{V}{\sqrt{gD}}

水理学的な立場としては射流(Supercritical Flow, {F_r} > 1)から常流(Subcritical (Normal) Flow, 1 > {F_r})へ遷移する際に跳水(Hydraulic Jump)という大きくエネルギーを減衰させる現象が起きます.
いずれ跳水についても記事を書きたいと思います.

連続の式

連続の式は流体における質量保存則になります.
流体中の微小要素(コントロール・ボリュームCV)を考え,CVにおける流体の質量の時間変化単位時間にCV表面を流入出する流体質量が等しいことから次式が誘導されます.

 
\displaystyle \frac{ \mathrm{D}{\rho}}{ \mathrm{D}t} + {\rho} ({\nabla} \bullet {\boldsymbol v})= 0

左辺の第1項は流体密度の時空間的な変化を表し,この項を無視できる場合を非圧縮流れといいます. すなわち,非圧縮流れとは単に圧力pによって密度 {\rho}が変化しない   \displaystyle \frac{\partial \rho}{\partial p} = 0 だけでなく,流体が圧力pや時間,場所によっても変化しない 
\displaystyle \frac{ \mathrm{D}{\rho}}{ \mathrm{D}t} = 0 のことを指します.

このあたりは日野先生の明解水理学の補遺に記述がありました.


おわりに

数式打つのが慣れなくて大変です.

正直,今回少し示した流体力学の導入部分(弾性体としての流体の扱いやNS eq.の誘導,流体の状態方程式など)は理解しきれていない部分が多くあります. のちのちこういう基礎部分の理解が応用につながると思うので,改めて詰めていきたいと思います.

次回から少しずつ流体の数値シミュレーションへ踏み込んでいきます.

参考図書,webサイト

数値流体力学ハンドブック

数値流体力学ハンドブック

  • 発売日: 2003/03/01
  • メディア: 単行本

流体力学 (JSMEテキストシリーズ)

流体力学 (JSMEテキストシリーズ)

明解水理学

明解水理学