Re=0.1の二次元円柱の抗力係数について

[md]

いつもお世話になっております．

現在、二次元円柱の抗力係数をレイノルズ数が0.1の条件で計算しております．

以下が今回の計算に用いているファイルです．

https://github.com/MurakiDaisuke/cylinder_2D_Re0.1.git

円柱直径0.01m，一様流速0.1m/s，動粘度0.01でicoFoamで計算しています．

同じ大きさの領域を用いて動粘度を変更し，Re=1,10,100でも計算を行った際には

抗力係数の値は実験値に近い値が得られました．

しかし，Re=0.1の計算では抗力係数がかなり大きな値となってしまいました．

どこに問題があるのかわからない状態なので，何か知っている方がいましたら

ご教授頂けると幸いです．

そもそもOpenFOAMでこのようなかなりレイノルズ数が低い場合に正しい計算が行えるのでしょうか．

よろしくお願いいたします．

[k k]

こんにちは。kです。

自分も球体周りの抗力係数が文献値に合うか検証したことがありますが、レイノルズ数の小さい領域（Re=0.01）くらいまでそこそこ合います。
こちらをご参考ください。
GitHub – kamakiri1225/sphere_forceCoeffs

こちらは、すべての領域で定常解析をしています(simpleFoam)が、高レイノルズ数の抗力係数が落ち込む部分（ドラッグクライシス？）以外は、定常にしても非定常にしてもそこそこ合いました。ドラッグクライシスも合わせようとすると、LESで相当数のメッシュが必要と論文で見かけたので諦めました。

モデルを見ましたが、代表長さやレイノルズ数の計算はあっているとして気になる点を示しておきます。

・使用しているのはOpenFOAMv2206かと思います。extrudeMeshDictだけOpenFOAM 9のバージョンになっているので、Allrunを実行するとここでエラーになります。

・自分が検証した限り抗力係数はかなり物体まわりのメッシュにこだわる必要がありましたので、境界層は入れた方が良いです。
snappyHexMeshで「addLayers       true;」としてする必要があります。

・以下余談です。
「relativeSizes false;」となっていますが、これだと境界層の厚みの数値をexpansionRatioとfinalLayerThicknessで与える必要があり、メッシュサイズや形状の大きさが変わると境界層ができているかどうか気付きにくいので、以下の方が良いかと思います。
以下でもだいたい同じものができると思います。

relativeSizes true;

finalLayerThickness 0.3;

minThickness 0.03;

これでまわりのメッシュに合わせて相対的にメッシュが作られます。

・sidewallは自然流出のような境界条件にしていますが、slip条件の方が良いと思います。
正直これが抗力係数に影響するとは思えませんが。

あとは代表長さが間違っていないかは再度確認してみた方が良いかと思います。
ここを間違っていると偶然合うこともありますが、流速条件を変えると合わなくなります。

2023年5月15日月曜日 15:49:14 UTC+9 md:

[md]

k様

ご返信ありがとうございます。

メッシュ作成の方は頂いたアドバイスを参考に現在行っている途中ですが、

いくつか行ったことと疑問に感じたことを以下に示します。

・sidewallの境界条件をslipに変更して計算したところ、

抗力係数が104(自然流出条件での計算時87)となり実験値の60からさらに大きく離れた結果となりました。

メッシュの方は特に変更していないのでこの結果だけ見ると境界条件は自然流出条件の方が良いように感じました。

・次のように広いレイノルズ数で球の抗力係数を計算していると思います。

https://takun-physics.net/13482/

この際レイノルズ数を変更するために、速度を変更しているように見えました。私の計算では、レイノルズ数を合わせるために速度ではなく動粘度を変更して計算を行っています。速度を変更すると代表長さ/速度で表される無次元時間が変化し時間刻みを変える必要があると考えています。そうすると変更するパラメータが増えるため動粘度だけを変更する方が簡単に感じたのですが、速度を変えることと動粘度を変えることで計算に違いはありますでしょうか。なぜ速度を変更したのか教えていただきたいです。

（そもそもsimpleFoamだと時間刻みはあまり問題にならないかもしれませんが、私の場合はicoFoamを使用しているので時間刻みも変える必要があるからこのように感じるのでしょうか。）

・上記の質問とも関連するのですが、レイノルズ数によってメッシュをどのように変更していますでしょうか。レイノルズ数が低くなるにつれてメッシュを細かくしすぎると計算が発散してしまうと思っています。

あまり分かっていないことも多いため、変なことを言ってしまっていたら申し訳ないですが、

よろしくお願いいたします。

2023年5月16日火曜日 19:55:16 UTC+9 k k:

[k k]

こんにちは。

・外部の境界条件を変えてそこまで抗力係数に影響出るものかなと自分の手元でもシミュレーションしてみましたが、層流の流れなのにかなり乱れた流れになっていますが、そちらでもそのようになっていますか？

　外部の境界条件が影響ない程度には広くとっているように見えますが、どこか設定がおかしいような気がしました。

・メッシュの決め方については以下の記事のようにReの大きさから最小渦スケールを見積もって、そのサイズを再現するためには4メッシュ（変化をとらえれるくらいのサイズ）と私は見積もっています。

　層流の場合はメッシュで対応できますが、乱流域になるとそうはいかないため結局最小渦スケールよりは粗めのメッシュになるので、乱流モデルを導入しています。

　https://www.cradle.co.jp/media/column/a428

    ちなみに私はRe=1000くらいで乱流モデルを使う使わないを切り替えてパラスタしました。

・球体周りの解析はボールの速度を変えた時の変化を見たかったので流速を変えましたが、相似則から考えて粘性でレイノルズ数を合わせても良いと思っています。

　もしdelta tの調整が面倒でしたらクーラン数の自動調整してはどうでしょうか。

　クーラン数が大きいのでΔtを調整する必要がありそうです。

　クーラン数固定のΔtの自動調整はcontrolDIctで↓このような感じだったと思います。

adjustTimeStep  yes;

maxCo           0.9;

icoFoamで使えるかは知りません。

使えなかったらdelta tを小さくしてみてください。

計算は発散はしていないですが、挙動がおかしいような気がしますので、残差と連続式の誤差も合わせて見た方が良いかもしれません。

既知かもしれませんが参考記事です。

http://opencae.gifu-nct.ac.jp/pukiwiki/index.php?plugin=attach&refer=%C2%E8%A3%B2%A3%B4%B2%F3%CA%D9%B6%AF%B2%F1%A1%A7H250810&openfile=2Dflow_around_cylinder_with_pisoFoam_rev.pdf

2023年5月18日木曜日 17:50:31 UTC+9 md:

[k k]

連投で申し訳ありません。

Δtを小さくしたらおっしゃっている通り抗力係数が少し大きめの値になりました。

改善策は提示できないですが、自分なら以下を試していきます。

・レイノルズ数をふって文献の曲線になるか試す(ふった結果合わないのがRe=0.1ということかもしれませんが）

・メッシュを色々と試す

・離散化スキームを試す

・定常なので定常ソルバで検証する（非定常解析よりはかなりはやく検証できます）

別で気になるところですが、流入条件を流速と圧力の両方を固定値にすると安定しないので以下の方が良いです。

    inlet

    {

        // type            fixedValue;

        // value           uniform 0;

        type            zeroGradient;

    }

2023年5月18日木曜日 21:27:08 UTC+9 k k:

[k k]

連投ですみません。

実験値の点を補間した線上の値がRe=0.1でCd=60というだけで、Re=0.1の実験データがあるのかというのは疑問に思いました。

そのへん知識がなく申し訳ありません。

色々検証して何かわかりましたら、是非報告いただければと思います。

2023年5月18日木曜日 22:17:16 UTC+9 k k:

[md]

k 様

ご返信ありがとうございます。

・境界条件について

以下の左の写真がsidewallの境界条件をinletOutletにした時、右がslipとして時になります。

乱れた流れとまではいきませんが、境界条件によって流速分布がかなり違うように見えます。

また、流入条件を流速と圧力両方固定値にしていたのを、圧力：zerogaradientとして計算してみましたが

抗力係数は少し大きくなり、さらに理想値との差が大きくなってしまいました。

これは計算が安定しないときに変えた方がよいということでしょうか。

・メッシュサイズについて

メッシュサイズ1辺の長さはRe=0.1の場合、(0.1)^(3/4)=0.18となり円柱の直径0.01mを代表長さLとするとL/0.18=0.06 mとなるので

最低でも0.06mよりメッシュを細かくするれば良いという認識であっているでしょうか。

レイノルズ数がかなり小さくなっていくと、メッシュの最低サイズは大きくなっていくと思いますが逆にどの程度の細かさで計算を行うべきかが難しいと感じています。

メッシュを細かくするとdtを小さくとらないとクーラン数が大きくなってしまいますし、計算時間が長くなってしまうと考えています。

・Re=0.1での抗力係数について

仰るとおり、今のところRe=0.1での抗力係数の実験値は見つけられていません。そのためOseen-Lambの式を用いて計算したCdの値(58.4くらい)と比較しています。

Re=1, 10, 100でも計算をしましたが、その際は 以下の文献にある実験値と比較しています。

 D.J.Tritton, Experomental on the flow past a circular cylinder at low Reynolds numbers,J.Fluid Mech. 6 (1959) pp547-567.

メッシュなど検討すべき箇所はまだ多いので、また何か分かりましたらご報告させていただきます。

よろしくお願いいたします。

2023年5月18日木曜日 22:23:54 UTC+9 k k:

[k k]

こんばんは。

＞ これは計算が安定しないときに変えた方がよいということでしょうか。

境界で流速と圧力の固定値にした場合、両方を満たすように解くことが難しいため、方程式が成り立たなくなるためです。

音速を超えるような場合は情報が左から右へ順番に伝わる（右へ伝わることがない）として、両方固定値にすることはあるようですが、非圧縮の場合はどちらかを固定値にして片方は微分値0にしないと成り立たないからです。

イメージですが、境界の隣の点は方程式を解いた結果で値が決まり例えばp=100となったとします。しかし、境界が無理やり固定値0にされていると、∇pが大きくなっていつか発散してしまいます。

離散の数値なのでどちらも固定値でも変動がなかったらあまり影響ないかもしれませんが。

チュートリアルなどで確認してみてください。

$FOAM_TUTORIALS/incompressible/simpleFoam/motorBike

定常解析にして比較してみて、もし同じような結果なら定常解析で検討した方が良いかもしれません。

お役に立てず申し訳ありません。

よろしくお願いいたします。

2023年5月18日木曜日 23:37:34 UTC+9 md:

[k k]

✕(右へ伝わることはない)

○(左へ伝わることはない)

2023年5月19日(金) 午前0:04 k k <kusuko...@gmail.com>:

--
このメールは Google グループのグループ「OpenFOAM」に登録しているユーザーに送られています。
このグループから退会し、グループからのメールの配信を停止するには openfoam+u...@googlegroups.com にメールを送信してください。
このディスカッションをウェブ上で閲覧するには https://groups.google.com/d/msgid/openfoam/db56f2fc-6e5d-46af-8720-abd9a9797949n%40googlegroups.com にアクセスしてください。

[md]

k様

ご返信ありがとうございます。

拘束がきつすぎて計算がうまくできない可能性があるということなのですね。

色々と丁寧にご教授いただきありがとうございます。

また何か気づいた点がありましたらよろしくお願いいたします。

2023年5月19日金曜日 10:15:53 UTC+9 k k:

[md]

皆様、お世話になっております．

二次元円柱の抗力係数について検討を行いましたので共有いたします．

今回は境界条件と領域のサイズについて検討しました．

・境界条件について
以前の計算では今回の資料のcase1の条件で計算をしていましたが
速度と圧力の条件にどちらも固定値を設定するのは良くないということで条件を変更しました．

ただし、抗力係数の値としては以前の条件の方が実際の値に近い値となっているため
以前の条件の方が正しい境界条件といえるのでしょうか．

・領域サイズについて
Re=0.1の時の計算結果が合わなかったため、領域サイズをさらに大きくしてみました．

結果としては、特にレイノルズ数が低いときに抗力係数に大きな変化がありました．
以前計算が合っていなかったRe=0.1においてもかなり近い値にはなってきましたが

ここまで領域を大きくしないといけないのは不思議な感じがします．

何か気づいた点がありましたら、ご教授いただけると幸いです．

2023年5月19日金曜日 13:04:49 UTC+9 md: