simpleFoamでの周期境界条件の設定方法について

2,233 views
Skip to first unread message

Y.T.

unread,
Oct 18, 2012, 5:30:23 AM10/18/12
to open...@googlegroups.com
はじめて質問させていただきます。 Y.T.と申します。
OpenFoam初心者なので、初歩的な質問となってしまいますが申し訳ございません。

現在、ソルバーは[simpleFoam]で2次元で周期境界条件を用いて定常状態を求める簡単なケースで解析をしようとしております。

最初は、流入面と流出面にcyclicを使用しておりましたが、
以前の別の方からの投稿にありましたように圧力が減衰していってしまうので収束しませんでした。
この解決方法として以下の方法があるようでしたが、どれを試してもうまくいきませんでした。
①cyclicではなく、jumpCyclicを使用する。
②mappedを使用してu,k,epsilon,nutは流出面の値を用い、pに関してはzeroGradientを用いる。

以下、上記の解決方法に取り組んだ結果を書いていきますので、何かお気づきの点がございましたらご指摘いただけると非常に助かります。

①cyclicではなく、jumpCyclicを使用する。
=>まず境界条件としてjumpCyclicが使えない。
  typeにjumpCyclicを使用しておりましたが、そんなものは使えないとはじかれてしまいました。
  そこで少し探してみると、jumpCyclicはそのもととなるファイルが「fvPatchFields」の中にはありましたが、「fvPatchs」の中にはありませんでした。
  この違いなのか、もっと根本的に使い方を間違っているのか。C++は今まで使用したことがないので、内容が読めずに困っていました。
  web上にもjumpCyclicに関する記述はほとんどなく、手詰まりでした。


②mappedを使用してu,k,epsilon,nutは流出面の値を用い、pに関してはzeroGradientを用いる。
=>エラー文として
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time

Create mesh for time = 0

Reading field p

Reading field U



--> FOAM FATAL ERROR: 

    patch type 'genericPatch' not type 'mappedPatchBase'
    for patch fieldIn of field U in file "/home/****************************/2D_windtunnels/cyclic/case2/0/U"

    From function mappedFixedValueFvPatchField<Type>::mappedFixedValueFvPatchField
(
    const fvPatch& p,
    const DimensionedField<Type, volMesh>& iF,
    const dictionary& dict
)

    in file fields/fvPatchFields/derived/mappedFixedValue/mappedFixedValueFvPatchField.C at line 115.

FOAM exiting
////////////////////////////////////////////
と出てしまい、解析を回すまでいきません。
以下に各種設定ファイルを載せておきますのでお気づきの点があればご指摘いただけると大変嬉しいです。
ちなみに、流入面と流出面は平行でDX=0.24です。使用しているOpenFOAMのバージョンは2.1.1です。

FoamFile
{
    version     2.0;
    format      ascii;
    class       polyBoundaryMesh;
    location    "constant/polyMesh";
    object      boundary;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

6
(
    fieldIn
    {
        type            mapped;
        nFaces          201;
        startFace       77964;
        offset          ( 0.24 0 0 );
        samplePatch     fieldOut;
        sampleMode      nearestPatchFace;
        faces           ( ( 0 22 58 36 ) ( 22 29 65 58 ) );
    }
    fieldOut
    {
        type            patch;
        nFaces          201;
        startFace       78165;
    }
///////////////////////////////////////////////////////////////////////////////////////
FoamFile
{
    version     2.0;
    format      ascii;
    class       volVectorField;
    object      U;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 1 -1 0 0 0 0];

internalField   uniform (1 0 0);

boundaryField
{
    fieldIn
    {
        type                mapped;
        value               uniform (1 0 0);
        interpolationScheme cell;
        setAverage          true;
        average             (1 0 0);
    }

    fieldOut
    {
        type            zeroGradient;
    }
/////////////////////////////////////////////////////////////////
FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      epsilon;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 2 -3 0 0 0 0];

internalField   uniform 1e-06;

boundaryField
{
    fieldIn
    {
        type                mapped;
        value               uniform 1e-06;
        interpolationScheme cell;
        setAverage          false;
        average             1e-06;
    }
    fieldOut
    {
        type            cyclic;
    }
//////////////////////////////////////////////////
FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      k;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 2 -2 0 0 0 0];

internalField   uniform 1e-08;

boundaryField
{
    fieldIn
    {
        type                mapped;
        value               uniform 1e-08;
        interpolationScheme cell;
        setAverage          false;
        average             1e-08;
    }

    fieldOut
    {
        type            zeroGradient;
    }
//////////////////////////////////////////////
FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      nut;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 2 -1 0 0 0 0];

internalField   uniform 0;

boundaryField
{
    fieldIn
    {
        type                mapped;
        value               uniform 0;
        interpolationScheme cell;
        setAverage          false;
        average             0;
    }
    fieldOut
    {
        type            zeroGradient;
    }
//////////////////////////////////////////////////
FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      p;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 2 -2 0 0 0 0];

internalField   uniform 0;

boundaryField
{
    fieldIn
    {
        type            zeroGradient;
    }

    fieldOut
    {
        type            fixedValue;
        value           uniform 0;
    }
///////////////////////////////////////////////

以上です。
初歩的な質問なのに長くなってしまい大変申し訳ございませんが、ご助言いただけるととても嬉しいです。
どうぞよろしくお願いします。

Y.T.

unread,
Oct 19, 2012, 5:02:48 AM10/19/12
to open...@googlegroups.com
昨日このトピックを投稿したY.T.です。

②mappedを使用してu,k,epsilon,nutは流出面の値を用い、pに関してはzeroGradientを用いる。
に関してですが、
・boundaryの流入面のtypeが、mappedではなくmappedPatch
・epsilonの流出面のtypeがcyclicだった
等の理由で回りませんでしたが、それらを修正した結果解析は回るようになりました。
初歩的なミスで申し訳ありません。恥ずかしい限りです。

しかしながら、解析はいまだに収束は向かえないようです。
メッシュ等も確認してみながら、解析を続けていきたいと思います。

また、どなたかjumpCyclicの使用方法をご存じの方がいらっしゃったらご教授いただけると幸いです。
どうぞよろしくお願いいたします。


2012年10月18日木曜日 11時30分23秒 UTC+2 Y.T.:

ohbuchi

unread,
Oct 21, 2012, 7:54:14 AM10/21/12
to open...@googlegroups.com
こんばんは。
cyclicJumpの件ですが、これを直接使うことはせず、これから派生した
fixedJumpか、更に派生したfanを使うことになると思います。
fanの使用例は、pimpleFoam/TJunctionFanチュートリアルが参考になります。


2012年10月19日金曜日 18時02分48秒 UTC+9 Y.T.:

Y.T.

unread,
Oct 22, 2012, 7:41:48 AM10/22/12
to open...@googlegroups.com
ohbuchi様、

情報ありがとうございます。
cyclicJumpはそのまま使えないのですね。
私も過去の投稿を参考にfanを使用してみましたが、逆に風速がどんどん増加していき収束する気配がありません。
0/U,k,epsilon,nutはcyclicを使用し、pは以下のように設定してあります。
流入面と流出面をfanとして、圧力が一定になるようにしていると思うのですが、これはおかしいのでしょうか。
初歩的な質問を連投してしまい本当に申し訳ありませんが、何かお気づきの点がございましたら指摘していただけると嬉しいです。
どうぞよろしくお願い致します。


FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      p;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 2 -2 0 0 0 0];

internalField   uniform 0;

boundaryField
{
    fieldIn
    {
        type            fan;
        patchType       cyclic;
        jump            uniform 0;
        f               1(-0.1);
        value           uniform 0;
    }

    fieldOut
    {
        type            fan;
        patchType       cyclic;
        jump            uniform 0;
        f               1(-0.1);
        value           uniform 0;
    }



2012年10月21日日曜日 13時54分14秒 UTC+2 ohbuchi:

Masashi Imano

unread,
Oct 22, 2012, 8:24:41 AM10/22/12
to open...@googlegroups.com
今野です。

この件に関連して以前ohbuchiさんが過去のスレッドで以下のように回答されていました。

https://groups.google.com/forum/#!msg/openfoam/qRaaWl5k2l4/AiJ7Tfyd2xMJ

ここでfanのソースに以下の部分があるので、圧力jump量は強制的に0以上にされていまいます。

jump_ = max(jump_, scalar(0));

従って以下のように f = -0.1 と負であるのは意味がありません。

> f 1(-0.1);

なお、以下のスレッドでも書きましたが、Ver.2.1以降であれば、simpleFoamにて
平均圧力勾配を制御して風量を一定に制御することもできると思います。

https://groups.google.com/d/topic/openfoam/qRaaWl5k2l4/discussion

以上、ご参考まで。

2012年10月22日 20:41 Y.T. <yoshiter...@gmail.com>:
> --
> このメールは Google グループのグループ「OpenFOAM」の登録者に送られています。
> このディスカッションをウェブ上で閲覧するには、https://groups.google.com/d/msg/openfoam/-/c5ph5dhMZ8AJ
> にアクセスしてください。
>
> このグループに投稿するには、open...@googlegroups.com にメールを送信してください。
> このグループから退会するには、openfoam+u...@googlegroups.com にメールを送信してください。
> 詳細については、http://groups.google.com/group/openfoam?hl=ja からこのグループにアクセスしてください。



--
IMANO Masashi, Ph.D.

Y.T.

unread,
Oct 23, 2012, 12:03:43 PM10/23/12
to open...@googlegroups.com
今野様、

ご助言ありがとうございます。
現在国外で勉強中なのですが、周りにあまりOpenFoamを使用できる方がいないため大変励みになります。

>従って以下のように f = -0.1 と負であるのは意味がありません。 

そうだったのですか。
正の値をいれると望んでいる方向と逆向きに風が流れてしまったので、負としてしまいました。
しかしながら正の値を入れてもうまく動作・収束しませんでした。

今野さんがご紹介されていたsourcesPropertiesを用いて解析を回してみましたが、こちらもなかなか収束しませんでした。
流入面、流出面はすべてcyclicにし、constantに上記ファイルを入れて解析を回しました。
不躾ながら上記のファイルについて少し質問してもよろしいでしょうか。
上記ファイルに下記のような文があると思います。
      Ubar        (1 0 0);                 // desired average velocity
      gradPini    gradPini [0 1 -2 0 0] 0; // initial pressure gradient
私は上記のようにしましたが、Ubarはドメイン内の平均風速値で、
下段のgradPiniは最初の圧力勾配のみを指定しているだけで、
後はUbarの値になるように圧力勾配を付けているという考えでよろしいのでしょうか。

また、収束しないのはドメインのせいではないかという気もしてきました。
現在、street canyonの解析をしていますが、建物の幅とstreet canyonの幅が
等しいようなものが無限に続いているような状態を、周期境界条件を用いて解析しております。
途中の結果を見てみると、全体的にはそれらしい結果が出ていますが、
street canyonの周辺はLESの解析のように変動し続けていて、収束する気配がありません。
そこにくると残差も変動し続けています。

流入面と流出面の位置、緩和係数等、いろいろいじってみましたが、中々収束は迎えません。

こちらに関しましても何かご助言等いただければ大変嬉しいです。
どうぞよろしくお願いいたします。

Y.T.

2012年10月22日月曜日 14時24分53秒 UTC+2 Masashi Imano:
--
IMANO Masashi, Ph.D.

Y.T.

unread,
Nov 1, 2012, 7:04:28 AM11/1/12
to open...@googlegroups.com
Y.T.です。
少し進展したのでご報告と、いまだに問題もあるのでまたご相談させていただきたく投稿いたします。
現在sourcesPropertiesを用い、以前と同様のstreet canyonを周期境界条件で解いています。

以前は、kとepsilonのinternalFieldに0を入力してしまっていたため、うまく収束しなかったようです。
現在はある値を入力しています。そうすると以前のように結果が変動し続けるような結果ではなく、
ある程度までは安定したような分布が得られるようになりました。
しかしながら、あるところまでは残差は落ちるようになったのですが、それ以降残差が落ちなくなってしまいました。
特に圧力と乱流エネルギーです。
そのまま解析を続けると、風速の分布も変な分布になってしまい困ってしまいました。
緩和係数もいじってみましたが、効果はあまりありませんでした。

そのときの設定は、風速の平均値は1.0m/sになるようにsourcesPropertiesを用い、
kとepsilonのintenalFieldの初期値には双方ともに0.015を使用しました。


また、とりあえずメッシュの妥当性を確かめてみようと思い、乱流計算ではなくlaminarにして層流計算にしてみたのですが、
こちらもやはり残差があるところで変動し続けるという結果になってしまっています。
メッシュはもう十分細かく、アスペクト比もそれほど大きくはないと思うのですが。。。
一応以下にcheckMeshの結果も載せさせていただきます。

周期境界だと、収束しにくかったりするのでしょうか。
ある程度まで残差等が落ちたら、そこで切り上げてしまってもよろしいのでしょうか。
ご意見いただけると幸いです。どうぞよろしくお願いいたします。


Create time

Create polyMesh for time = 0

Time = 0

Mesh stats
    points:           182868
    internal points:  0
    faces:            363455
    internal faces:   180589
    cells:            90674
    boundary patches: 6
    point zones:      0
    face zones:       0
    cell zones:       0

Overall number of cells of each type:
    hexahedra:     90674
    prisms:        0
    wedges:        0
    pyramids:      0
    tet wedges:    0
    tetrahedra:    0
    polyhedra:     0

Checking topology...
    Boundary definition OK.
    Cell to face addressing OK.
    Point usage OK.
    Upper triangular ordering OK.
    Face vertices OK.
    Number of regions: 1 (OK).

Checking patch topology for multiply connected surfaces ...
    Patch               Faces    Points   Surface topology                  
    fieldIn             397      796      ok (non-closed singly connected)  
    fieldOut            397      796      ok (non-closed singly connected)  
    fixedWalls          144      292      ok (non-closed singly connected)  
    top                 192      386      ok (non-closed singly connected)  
    street              388      778      ok (non-closed singly connected)  
    frontAndBack        181348   182868   ok (non-closed singly connected)  

Checking geometry...
    Overall domain bounding box (0 -0.12 0) (0.24 0.6 0.02)
    Mesh (non-empty, non-wedge) directions (1 1 0)
    Mesh (non-empty) directions (1 1 0)
    All edges aligned with or perpendicular to non-empty directions.
    Boundary openness (-3.577e-19 2.1462e-19 5.27011e-18) OK.
    Max cell openness = 1.66531e-16 OK.
    Max aspect ratio = 2.5625 OK.
    Minumum face area = 3.6e-07. Maximum face area = 3.85177e-05.  Face area magnitudes OK.
    Min volume = 7.2e-09. Max volume = 4.81472e-08.  Total volume = 0.003168.  Cell volumes OK.
    Mesh non-orthogonality Max: 44.8504 average: 10.2797
    Non-orthogonality check OK.
    Face pyramids OK.
    Max skewness = 0.736999 OK.
    Coupled point location match (average 2.77556e-17) OK.

Mesh OK.

End

以上です。
どうぞよろしくお願いいたします。

Y.T.



2012年10月23日火曜日 18時03分43秒 UTC+2 Y.T.:

ohbuchi

unread,
Nov 1, 2012, 7:34:39 PM11/1/12
to open...@googlegroups.com
こんにちは。
まず、乱流計算に必ず定常解が存在する訳ではありません。領域内で剥離やせん断があると渦が生成され、
渦は本質的に不安定なので、揺動を続けます。
限定された範囲で揺動を続ける流れ場の場合、連続の式の残差が0にならないためpの残差が減らずに
ほぼ一定に落ち着きます。
Y.T.さんのケースでは乱流エネルギーが収束せず、計算を続けると変な流れパターンになるということが
問題だと思います。乱流モデルの設定が不適切である可能性があります。
変なパターンになったときのk、epsilon、nutの値をチェックして下さい。乱流モデルが原因でうまく計算できない
場合は乱流量が異常に大きくなっていることが多いです。

ご参考まで。


2012年11月1日木曜日 20時04分28秒 UTC+9 Y.T.:

Y.T.

unread,
Nov 2, 2012, 7:29:27 AM11/2/12
to open...@googlegroups.com
ohbuchi様、

ご意見ありがとうございます。大変励みになります。

>まず、乱流計算に必ず定常解が存在する訳ではありません。領域内で剥離やせん断があると渦が生成され、
>渦は本質的に不安定なので、揺動を続けます。
>限定された範囲で揺動を続ける流れ場の場合、連続の式の残差が0にならないためpの残差が減らずに
>ほぼ一定に落ち着きます。

そうだったのですか。そこまで複雑な流れ場の解析をしたことがなかったので、全て必ず定常解が存在するものと思っていました。
そのようにpの残差が一定で落ち着いたら、それはその結果を使用するものなのでしょうか。

また、計算を進め、残差をよく眺めてみましたが、乱流エネルギーのkは結構落ちるのですが、
やはりpだけなかなか残差が小さくならない印象でした。
paraView等で結果を確認してみてもnut、k、epsilon等に極端な以上は見られませんでした。

現在、fvSolutionのSIMPLEのresidualControlでは全ての変数を1e-06にしているのですが、
pが中々小さくならないため、計算し続け、変になってしまっているのかなとも思いました。

幾度にわたり初歩的な質問をしてしまい申し訳ありません。
どうぞよろしくお願いいたします。

Y.T.

2012年11月2日金曜日 0時34分39秒 UTC+1 ohbuchi:
Reply all
Reply to author
Forward
0 new messages