はじめに

本記事は多孔質体を含めた領域の流体解析を行うことを目的としています。

こちらに多孔質体ソルバーporousSimpleFoamの詳細がまとめられていますが、ここでは、pisoFoam等の基本ソルバーに多孔質体モデルを導入するやり方(fvOptionを利用）をまとめたいと思います。

環境

なお、本記事はFoundation版OpenFOAM v8、pisoFoam/laminar/porousBlockageを使用しています。

Porous model

多孔質モデルは、基本、運動方程式のソース項 $S_m$ に組み込まれます。

\frac{\partial \rho U}{\partial t} + (U \cdot \nabla) U = \nabla \sigma + S_m

Darcy-Forchheiman equation

Darcy-Forchheiman 多孔質モデルでは、ソース項は以下のように表されます。

S_m = - \left( \mu D + \frac{\rho |U|}{2} f  \right) U

ここで、D は Darcy coefficient [1/m$^2$]、f は Forchheimer coefficient [1/m]を表しています。

Power law

べき乗則モデルでは、

S_m = - \rho C_0 |U|^{(C_1 -1)/2}

のようにべき乗則で流動抵抗が表されます。ここで、$C_0$ は model linear coefficient、$C_1$ は model exponent coefficientであり、ユーザーが設定する必要があります。

その他モデル

その他にも、以下のようなモデルが用意されています。

設定方法

ここでは、先に述べたとおりpisoFoam/laminar/porousBlockageチュートリアルを見ていきます。

fvOptionを使った多孔質体モデルの実装では、

topoSetによる多孔質領域の設定
多孔質体の物性値の設定
の大きく２つを行う必要があります。

topoSet

system/topoSetDict

    // porousBlockage
    {
        name    porousBlockageCellSet;
        type    cellSet;
        action  new;
        source  boxToCell;
        sourceInfo
        {
            box (-0.5 -0.5 -1) (0.5 0.5 1);
        }
    }
    {
        name    porousBlockage;
        type    cellZoneSet;
        action  new;
        source  setToCellZone;
        sourceInfo
        {
            set porousBlockageCellSet;
        }
    }

fvOption

constant/fvOption

    type            explicitPorositySource;

    explicitPorositySourceCoeffs
    {
        selectionMode   cellZone;
        cellZone        porousBlockage;

        type            DarcyForchheimer;

        // D 100;  // Very little blockage
        // D 200;  // Some blockage but steady flow
        // D 500;  // Slight waviness in the far wake
        D 1000; // Fully shedding behavior

        d   ($D $D $D);
        f   (0 0 0);

        coordinateSystem
        {
            type    cartesian;
            origin  (0 0 0);
            coordinateRotation
            {
                type    axesRotation;
                e1  (1 0 0);
                e2  (0 1 0);
            }
        }
    }

このチュートリアルでは多孔質体が等方的な圧損を有している場合を表しています。

ハニカムなどの一方向のみ流れやすく、多方向には圧損が生じるような場合は、$d$、$f$の値を以下のように方向性を持たせて設定します。

    DarcyForchheimerCoeffs
    {
        d   d [0 -2 0 0 0 0 0] (7.65e7 3e10 3e10);
        f   f [0 -1 0 0 0 0 0] (255 1e5 1e5);

        coordinateSystem
        {
            type    cartesian;  // global co-ordinate system (redundant)
            origin  (0 0 0);    // redundant in this case
            coordinateRotation
            {
                type    axesRotation; // local Cartesian co-ordinates
                e1      (0.998 0.061 0);
                e2      (0.016 -1.000 0);
            }
        }
    }

なお、coordinateSystemは、多孔質体のlocal座標を設定するためのようです（どうしてlocal座標を設定するのか、その意図は指し図りきれていないです・・・）。local座標は３つのベクトル（$\bf{e_1}$, $\bf{e_2}$, $\bf{e_3}$）を使って設定しますが、入力ではうち２つを設定する必要があります。特にlocal座標を考える必要がない場合は、

coordinateRotation
{
    type    axesRotation; // local Cartesian co-ordinates
    e1      (1 0 0);
    e2      (0 1 0);
}

のようにx, yベクトルを設定します。なお、$\bf{e_3} = \bf{e_1} \times \bf{e_2}$ で自動的に規定されます。

porousBlockageの解析結果例

さいごに

以上、駆け足でfvOptionを使った多孔質体の設定方法を見てきました。以前は、自らでソース項をUEqn.Hなどに追加する必要がありましたが、fvOptionを使えばそのようなソースからのビルドをする必要がありません。
まだ、説明を端折っているところもありますので、順次、加筆していきたいと思います。