Microsoft Word - ...W_.....l.doc

Transcription

1 263 特集 オープンソースの大きな流れ OpenFOAM によるターボポンプの流れ解析 * 株式会社荏原製作所汎用水力技術開発室大渕真志 CFD Simulation of Turbo Pump with OpenFOAM Masashi OHBUCHI, EBARA Corporation 1 はじめに OpenFOAM 1) はオープンソース CFD ソフトウエアであり 近年大いに注目を集めている 80 年代末から英国 Imperial College で開発され 1999 年に設立された英国 Nabla 社によって FOAM(Field Operation And Manipulation) という名称で販売されていたが 2004 年にオープンソースとなった このため 研究コードからオープンソース化された他の CFD ソフトと比べて完成度が高いことが特徴となっている Fluent や StarCCM 等の商用コードが GUI ベースの統合解析環境を提供するのに対し OpenFOAM は有限体積法ライブラリと標準ソルバ及びユーティリティプログラム群という形で提供される このため ユーザは自分の問題を解くために 多くのユーティリティを使いこなす必要がある また メッシュや条件設定などは多数のテキストファイルで定義するため あまり使い易いとは言えない しかし 一般的な問題なら標準ソルバで十分に扱える完成度の高さと 一度内部構造を理解できれば比較的簡単にプログラムを拡張できる点がここまで注目を集めている所以であると思われる 特に カスタムプログラムの開発では 根幹部分を隠蔽している市販 CFD ソフトのユーザプログラミング機能に比べて遥かに柔軟性が高く 生産性も高い 例えば ラプラス方程式 T t 2 ( D T ) = 0 T (1) を解くプログラムは次の様に記述できる 2) * 神奈川県藤沢市本藤沢 [email protected] solve ( fvm::ddt(t)-fvm::laplacian(dt,t) ); ここに fvm::ddt() は式 (1) の時間微分項 fvm::laplacian() は拡散項に対応している 2) このプログラムは偏微分方程式を有限体積法で陰的に離散化してマトリクスに格納し 解を求めるまでを表している この様に OpenFOAM では離散化や境界条件の詳細に立ち入ることなく 極めてコンパクトに偏微分方程式の数値解法を記述することができる OpenFOAM は C++ 言語の機能 ( オブジェクト指向 ジェネリック プログラミング 演算子オーバーロードなど ) を駆使することでこの表記法を実現しており Fortran 等の構造化言語での手続き的プログラミングでは考えられない生産性の高さを実現している 3) 以下では OpenFOAM をターボポンプの流れ解析に適用する方法について解説する 2 ターボポンプの流れ解析に必要な機能ターボポンプは回転機械であるため 回転領域をいかに扱うかが課題となる 定常解析では 回転領域で遠心力とコリオリ力を付加する Multiple Reference Frame(MRF) という手法が用いられ 非定常解析では回転領域のメッシュを実際に回転させる Sliding Mesh という計算手法が用いられる 回転系と静止系のインターフェースでメッシュを一致させるには手間がかかるため 回転系 静止系のメッシュを別々に作成し メッシュ間のインターフェースで相互に補間しあう計算ができれば便利で

2 264 ある この不連続メッシュ インターフェース機能がターボポンプの実用解析では極めて重要である 3 OpenFOAM の対応状況 OpenFOAM には MRF やダイナミックメッシュの機能が以前から提供されていた しかし 最近まで不連続メッシュ インターフェース機能は存在しなかった 3.1 Extend 版について欧米のユーザが中心となって作られたコニュニティによって開発が進められている OpenFOAM の派生版 (Extend 版 ) がある 4) 2009 年にこの Extend 版で GGI(Generalized Grid Interface) と呼ばれる不連続メッシュインターフェース機能が実装された Extend プロジェクトのターボ機械向け機能について 5) は プロジェクト中心メンバ H.Jasak 氏による解説で詳しく論じられている Extend 版では GGI 以外にも正規版にない多くの機能や独自のソルバが提供されている 現在のバージョン 1.6-ext( 正規版の version1.7 相当 ) に含まれる独自機能には以下の様なものがある 豊富なダイナミックメッシュ機能 粘弾性流体ソルバ 液膜流ソルバ ブロックマトリックスのサポート Extend 版の開発はオープンに行われているため 自由に参加でき CFD-Online の Forum 6) で開発者にバグの報告や質問もできる また コミュニティには幾つかの SIG(Special Interest Group) が作られ Internet 上で成果が公開されている ターボ機械関連では TurboMachinery SIG があり 下記の成果を公開している 7) 性能評価ライブラリ(turboPerformance) CGNS コンバータ ガイドラインとサンプルデータこれらの成果を学ぶことで 実用的なターボポンプの流れ解析を実施できる様になる Extend 版の GGI には3つの種類がある Ggi 通常の GGI cyclicggi 周期境界に用いる GGI overlapggi ピッチの異なる動静翼境界用近日中に mixing plane 機能が正式にリリースされる予定である 5) 尚 overlapggi と連立方程式ソルバ GAMG は併用できない また overlapggi を適用する境界パッチの周方向エッジが放射状でないと極 めて不安定になるので注意されたい 3.2 正規版の対応正規版では最近まで不連続メッシュ インターフェース機能がなかったため ターボポンプ解析には Extend 版が必須という状況だった 2011 年 12 月にリリースされた OpenFOAM-2.1 で AMI(Arbitrary Mesh Interface) が導入され ようやく不連続メッシュ インターフェースを扱える様になった AMI の実装は概ね GGI と同様だが cyclicami という1 種類の境界パッチのみで Ggi と cyclicggi の機能をカバーする 現時点では overlapggi 機能のサポートはない 4 具体的方法と解析事例それでは OpenFOAM でターボポンプ流れ解析を行う具体的な手順を解説する 主な作業手順は下記の通りである (1) 作業用ディレクトリの作成 (2) メッシュの作成または変換 (3) 乱流モデル 物性値の定義 (4) 境界条件 初期条件の定義 (5) 回転領域の設定 (6) 並列計算の設定 (7) 計算条件 モニタリング設定 (8) 計算実行 (10) 可視化および後処理 OpenFOAM では問題ごとに作業用のディレクトリを作成する必要がある その構造は図 1に示す様にする 図 1 OpenFOAM の作業ディレクトリ構造

3 265 ここで system ディレクトリには計算時間や時間刻み 出力間隔などを定義する controldict や各フィールド変数の離散化スキームを定義する fvschemes 各方程式を解くソルバの設定や SIMPLE PISO などアルゴリズムのパラメータを定義する fvsolution などが含まれる 標準ユーティリティの多くは system/controldict の存在を確認して動作するものが多いため まずこのファイルを作成する必要がある 次にメッシュの作成であるが OpenFOAM のメッシュは polymesh と呼ばれる極めて柔軟な形式で定義される polymesh ではセルは任意の数のフェースで囲まれたものであり 各フェースは任意の数の節点で定義される polymesh によりテトラ ヘキサなどは勿論 歪み要素やハンギングノード # を持つ様なセル 最近普及したポリヘドラ要素なども問題なく表現することができる polymesh は constant/polymesh ディレクトリで定義され 下記のファイルで構成される points 節点データ face フェースデータ ( 節点 ID リスト ) owner 各フェースのオーナセル ID neighbour 各フェースの隣接セル ID boundary 境界パッチの定義図 2に示す様にフェースを定義する節点を右ネジのルールでフェースの向きを定義する 面積ベクトル Sf の内側をオーナセル 外側を隣接セルとする 従って 境界に設置されるフェースでは必ず隣接セルが存在しないことになる することで 6 面体セルを生成する 領域のエッジは 線分 円弧 スプラインで定義でき 領域分けさえ 適切に行えば品質の良いメッシュを生成できる し かし 設定に手間がかかるので複雑な実用的問題に はあまり向かない 一方 snappyhexmesh は形状デ ータを STL ファイルで与え blockmesh で作成し基 本メッシュが STL の 3 角パッチをカットした部分 を細分割し STL フェースにスナップさせることで 複雑なメッシュを自動生成するツールである 指定 した面にレイヤーメッシュを生成させることもでき る OpenFOAM-2.0 以降では STL 形状の特徴エッジ を保持することも可能になり 極めて実用的なメッ シュ生成ツールである 特徴エッジの抽出に利用す るツール (surfacefeatureextract) は 3 角パッチ間の なす角度で自動的に判定する 図 3 にインペラ形状 から抽出された特徴エッジの例を示す 図 3 抽出された特徴エッジ図 3 の様に特徴エッジは正しく認識されているのに snappyhexmesh でメッシュを作成すると正しくエッジを保持できない場合がある ( 図 4) 図 2 節点番号と面積ベクトルの関係 boundary ファイルでは各パッチの種類を定義する OpenFOAM に付属のメッシュ生成ユーティリティで作成したり メッシュ変換コンバータで変換されると自動的に boundary が生成されるが 境界パッチの種別は適当に割り当てられるため 自分で再定義する必要がある OpenFOAM には 2つのメッシュ作成ユーティリティ blockmesh と snappyhexmesh が含まれる blockmesh はブロック構造格子を作成するためのツールで 計算領域を6つの面で囲まれた領域に分割 図 4 特徴エッジが正しく再現されない例 snappyhexmeshdict の設定によってある程度改善されるが 形状が完全には再現できない場合があることに注意する必要がある やはり複雑な形状を扱う実用解析では市販のメッシュ生成ツール (ICEM や Pointwise など ) を使う方が格段に生産性が高い この場合には メッシュ変換ユーティリティを使って PolyMesh 形式に変換する ここでは fluent 形式のメッシュファイルか

4 266 ら polymesh 形式に変換する手順について説明する これには標準ユーティリティの fluentmeshtofoam を用いる >fluentmeshtofoam -writezones -writesets ***.msh ここで fluent メッシュファイル名を ***.msh とする -writezones オプションはメッシュが複数の領域からなる場合にメッシュ領域を cellzones というファイルに 境界パッチを facezones に出力するオプションである また -writesets オプションは領域や境界パッチを polymesh/sets ディレクトリにセルやフェースの集合ファイルとして出力するものである 尚 cellzones は単一セル領域のメッシュファイルからは生成されない 回転領域と静止領域を別の fluent メッシュファイルで作成し polymesh に変換してからマージ (mergemeshes ユーティリティを利用 ) すると cellzones は生成されない セル領域は回転領域の回転数を定義する際に必要となる マージされたメッシュから不連続な領域を抽出し cellzones ファイルを作成するには regioncellsets というユーティリティを使うと便利である このユーティリティは正規版 OpenFOAM には含まれず Extend 版 ( バージョン1.6) で提供されるものであるが 正規版のバージョン2.1でも問題なくコンパイルして利用することができる 境界条件を定義する前に 利用するソルバと乱流モデルを決めておく必要がある これにより扱うフィールド変数 (U p k epsilon など ) が決まる ソルバには 定常解析では MRFsimpleFoam を 非定常解析では pimpledymfoam を用いることが多い 乱流モデルには RAS LES DES など多くが選択可能であるが ターボポンプ用としては k-ω SST がよく利用されている 境界条件 初期条件の設定には下記のファイルを定義する constant/polymesh/boundary 0/U, p, k, omega, nut( 全フィールド変数 ) 図 1の time directories は時刻値を名前に持つディレクトリで 初期から計算する場合には 0 というディレクトリを作成し そこに各フィールド変数の初期 境界値を設定するファイルを定義する必要がある OpenFOAM の境界条件は図 5に示す様な階層構造になっている boundary ファイルで指定するのは Base type と呼ばれる境界条件で wall patch cyclic など (AMI や GGI も ) がある Primitive type や Derived type は Time directory の各フィールド変数ファイルで利用する 同じ境界パッチでも速度 U と圧力 p では設定する境界条件の種類が異なる 例えば 流入境界では U は fixedvalue で速度を規定するが p は zerogradient とする 一方流出境界では U は zerogradient p は fixedvalue にするなどである RANS 乱流モデルのフィールド変数 (k epsilon nut) には壁面で壁関数を用いるが それぞれ kqrwallfunction epsilonwallfunction など異なるタイプで設定すべきパラメータも異なる これらの組み合わせや設定例は TurboMachinery SIG のサンプルや OpenFOAM に付属の Tutorial を参考にすればよい 図 5 境界条件の構成次に回転領域の定義であるが 定常解析の場合には constant/mrfzones を 非定常解析の場合には constant/dynamicmeshdict を設定する MRFZones ファイルには 回転領域のセルゾーン名と回転軸 回転数を定義する 下記に例を示す 1 // MRF 領域数 ( rotor // MRF 領域のセルゾーン名 nonrotatingpatches (); origin origin [ ] (0 0 0); axis axis [ ] (0 0 1); omega omega [ ] ; ) 回転速度 omega の単位は [rad/s] である 一方 dynamicmeshdict では利用するダイナミックメッシュライブラリの種類 回転領域のゾーン名 回転軸 回転数などを記述するが 正規版と Extend 版で記述内容が異なる 例えば Extend 版では dynamicfvmeshlib "libtopochangerfvmesh.so"; dynamicfvmesh mixerggifvmesh;

5 267 mixerggifvmeshcoeffs coordinatesystem type cylindrical; origin (0 0 0); axis (0 0 1); direction (100); rpm ; // rpm slider moving ( IMP.IN IMP.OUT ); static ( SUCT.OUTLET CAST.INLET ); などと記述する 特徴は slider に回転側と静止側の GGI パッチ名を記述することと ダイナミックメッシュライブラリに mixerggifvmesh という Extend 版にのみ存在するものを指定する点である 尚 この mixerggifvmesh では回転領域は必ず rotor という単一のセルゾーンでなければならない 一方 正規版では dynamicfvmesh solidbodymotionfvmesh; motionsolverlibs ( "libfvmotionsolvers.so" ); solidbodymotionfvmeshcoeffs cellzone rotor; solidbodymotionfunction rotatingmotion; rotatingmotioncoeffs CofG (000); radialvelocity ( ); // deg/s ライブラリに solidbodymotionfvmesh を用い AMI パッチの記述が不要になっている また 角速度はベクトルで指定し単位は [deg/s] となっている 以上で 設定はほぼ完了である 作業フォルダに移動して ソルバ名を入力すれば計算が実行できるが 規模の大きい問題は並列化することで大幅に高速化できる 並列計算は次の様な手順で実行する system/decomposepardict を作成 decomposepar ユーティリティで領域分割 ソルバを並列実行して計算実施 計算結果を reconstructpar でマージする decomposepardict では並列分割数と分割方法を指定する OpenFOAM では単純に方向別で領域分割するだけでなく Metis 8) や Scotch 9) などのライブラリも利用できる Extend 版では GGI パッチを globalfacezones に指定する必要がある つまり Extend 版では GGIパッチを全計算ノードで共有する このためか Extend 版で GGI を使った計算の並列効率は正規版で AMI を使った場合より若干劣る様である ソルバの並列実行は -parallel オプションと並列数を指定すればよい OpenFOAM のインストールディレクトリ下 bin/tools/runfunction で定義されている runparallel という関数を使うと簡単である >. $WM_PROJECT_DIR/bin/tools/RunFunction >runparallel MRFSimpleFoam 4 最後の4は並列数である decomposepardict の設定と一致させる必要がある 計算が完了したら reconstructpar ユーティリティを実行すれば シリアル計算と同じように後処理を行うことができる ベクトル図 流線 等値面 表面および断面コンタ図といった一般的な可視化処理には オープンソースの可視化ツール paraview を用いることができる 10) ターボポンプの解析を実行したら ポンプ性能を評価する必要が出てくる Extend 版では turboperformance ライブラリを利用すれば計算実行と同時にヘッド 軸動力 効率が計算されファイルに記録される 利用方法は controldict の末尾に下記の記述を加えればよい functions ( turboperformance type turboperformance; functionobjectlibs ("libturboperformance.so"); inletpatches (inflow); // 流入パッチ名 outletpatches (outflow); // 流出パッチ名 patches ( vane hub shroud ); // 回転壁パッチ log true; rhoinf 998.0; // 密度 CofR (0 0 0); // 回転中心座標 omega ( ); // 角速度 rad/s );

6 268 一方 正規版では turboperformance ライブラリを利用できないので 自分で評価する必要がある まずヘッドの算出には 流入 流出パッチでの全圧平均値を求めて差を求める必要があるが パッチ平均を求めるには標準ユーティリティの patchaverage を用いる まず ptot ユーティリティで全圧を計算する >ptot これにより各 time directory に ptot ファイルが作成される 次に patchaverage ユーティリティでパッチ平均値を流入 流出パッチについて求める >patchaverage ptot inflow >patchaverage ptot outflow この差を計算して重力加速度で割ればヘッドが求められる 尚 OpenFOAM の非圧縮ソルバでは圧力は p/ρ で定義されている patchaverage は面積平均であるが 圧力はエネルギーなので正確には質量加重平均にすべきである 標準ユーティリティに質量加重平均を求めるものはないが patchaverage を元にカスタムユーティリティを作ることは比較的簡単であるし 非公式の情報サイト OpenFOAM Wiki 11) を主催する Bernhard Gschaider 氏が開発したライブラリ swak4foam 12) に含まれる simplefunctionobjects を使って質量加重平均を求めることができる 次に トルクは標準で提供される functionobject の forces を利用する functions forces type forces; functionobjectlibs ("libforces.so"); outputcontrol timestep; outputinterval 1; patches (hub shroud vane); // 壁面パッチ pname p; UName U; rhoname rhoinf; log true; rhoinf 998; // 密度 CofR (000); // モーメント中心 これにより forces ディレクトリが作成され 指定したパッチに作用する力とモーメントの各成分が出 力される トルクはモーメントの Z 成分 ( 回転軸が Z の場合 ) である 4.1 遠心ポンプ段落解析 遠心ポンプのインペラの 1 段落のみを取り出し周 期境界を適用して計算することは インペラ設計の 初期段階のスクリーニングとして頻繁に行われる 計算規模は小さいが 沢山のケースを計算する必要 があるため 市販の CFD ソフトウエアを利用する場 合にはライセンス数の制限を受けやすい問題でもあ る 我々は標準ソルバ SRFSimpleFoam MRFSimpleFoam を用いて Windows 上で動作する羽 根段落流れ専用解析ツールを開発した ( 図 6) 図 6 Windows 上での設計用解析ツール群この様に OpenFOAM はオープンソースであるた め 計算内容やファイル I/O を自由に拡張できるの で特定用途向けの専用開発ツールを開発するのに適 している 尚 Windows 版 OpenFOAM には市販品 も幾つかあるが Symscape 社が公開しているパッチ を使えば自分で作ることもできる 13) し バイナリ 版もインターネット上に公開されている 14) 4.2 斜流ポンプ定常一体解析 インペラと静止系流路 ( ガイドベーンやボリュー トケーシングなど ) を一体で計算することはポンプ 設計において性能予測のために不可欠である 計算 の負荷を考え 設計流量付近の運転条件では定常解 析を行うことが多い 図 7 に斜流ポンプの計算例を 示す 4.3 遠心ポンプ非定常解析 非設計点ではポンプ内部流れは非定常性が強まるため 定常解析では性能を正確に予測することが困難になる そのため計算負荷は大きくなるが Sliding Mesh を使った非定常解析を行う必要がでてくる 定常解析との設定方法の違いは ダイナミックメッシ

7 269 ュの定義くらいである 図 8に遠心ポンプの非定常解析例を示す 羽根目玉部から吸込ケーシング内へ渦線が伸びているのがわかる を差し引いたもの p_ρgh を計算している このため 壁面では buoyantpressure という境界条件を用いて p_ρgh の勾配ゼロを指定している 図 7 斜流ポンプ定常一体解析の例 図 9 水槽内に設置された斜流ポンプの流れ解析 5 まとめ 図 8 遠心ポンプの非定常解析例 4.4 水槽とポンプの一体解析自由水面を有する吸込水槽流れの解析は良く実施する問題であるが 自由水面解析そのものが時間がかかる非定常計算であるだけでなく 特性時間が長いため長時間の積分を要するため なかなかポンプと一体で計算することは難しい しかし OpenFOAM で提供されている MRFInterFoam ではポンプ部分に MRF を利用できるため自由水面の特性時間で非定常解析を行うことが可能である 図 9に自由水面を有する水槽内に設置した斜流ポンプの解析例を示す 定格よりもかなり少ない流量での運転のため羽根入口逆流が水槽内に新入してベルマウスと底面の間に強い渦を形成しているのがわかる 本計算では自由水面を扱うため空気相の物性値や重力加速度の設定が必要になる 圧力は静水圧 ρgh OpenFOAM の概要とターボポンプの流れ解析に適用する具体的な手順について簡単に説明した 商用の CFD ソフトウエアとの比較を行ったところ 良好な一致を確認している 15) ので OpenFOAM は安心して利用できるレベルにあると考えている OpenFOAM の開発は正規版 Extend 版ともに活発に続けられており 毎年機能拡張している Extend 版ではネイティブ Windows 版の開発 Mac OS への移植も行われており 利用のための敷居が低くなることが期待できる 市販ソフトウエアに比べると使い勝手の点で大きく見劣りすることは否めないが コミュニティの努力により徐々に障害が取り除かれつつある フリーソフトであるためライセンス料金がかかからないというコストメリットばかりが注目されがちだが 一番の良さはカスタマイズし易さにある OpenFOAM をカスタマイズして使い勝手の良い専用ソルバを作り 設計現場で CFD を普及させることも十分可能である この様に 内部構造を理解しながら使うため 理論や数値解法 プログラミング技術も含めた本物の利用技術を獲得する理想的なツールであると思われる ぜひ OpenFOAM を利用してみて頂きたい OpenFOAM は SGI 社の登録商標である Fluent は ANSYS 社の登録商標である StarCCM は CD-Adapco 社の登録商標である

8 270 ICEM CFD は ANSYS 社の登録商標である Pointwise は Pointwise 社の登録商標である # ハンキングノードはセルを細分割するときに生ずるフェースエッジを内分するノードのことを言う 引用文献 1) OpenCFD s Home Page, 2) OpenFOAM User Guide オープン CAE 学会で翻訳版を公開している ソフトウエアマニュア ル翻訳 3) OpenFOAM Programming Guide ob/master/doc/guides-a4/programmersguide.pdf オープン CAE 学会で翻訳版を公開している 4) OpenFOAM Extend Project s Home Page, 5) Hrvoje Jasak: OpenFOAM Turbo Tools: From general purpose CFD to turbomachinery simulations (2011), AJK2011-FED 6) OpenFOAM Forum on CFD Online 7) Turbomachinery SIG s Home Page ery 8) Metis library s Home Page, 9) Scotch library s Home Page, 10) Paraview 11) OpenFOAM Wiki 12) swak4foam swak4foam 13) Symscape Inc,OpenFOAM2.1.x on Windows 64-bit with MS MPI,(2012), ws-64-mpi 14) bluecfd-singlecore 15) 大渕真志 :OpenFOAM によるターボポンプ流れ解析 CFX との比較. オープン CAE シンポジウム 2010(2010), オープン CAE シンポジウム 2010/