プロペラまわり流れの 数値計算法に関する研究 船舶海洋工学科　船舶コース 80392　八巻陽彦 船舶海洋工学科　船舶コースの八巻です。プロペラ周りの流れの数値計算法に関する研究 について発表致します。 2003/02/19（水）

研究背景 従来の舶用プロペラの設計法 ポテンシャル理論に基づく 理論計算と模型試験 問題点 粘性の影響 非線型性 ナビエ・ストークス方程式に基づく数値計算による設計が必要であるが、計算を行う際の前処理の困難さと精度の問題からまだCFDを用いての計算は一般的ではない。 舶用プロペラは従来、ポテンシャル理論に基づく理論計算と模型試験により設計されてきました。しかし、粘性や非線形性といった問題点が有ります。そこで、ナビエ・ストークス方程式に基づく数値計算による設計が必要になってくるのですが、精度の問題と計算を行う際の前処理の困難さからまだCFDを用いての計算は一般的ではありません。 2003/02/19（水）

研究目的 シャフト・ボスを含んだプロペラ幾何形状が複雑であり、 CADで完全に三次元形状を再現することが困難 舶用プロペラ周り流れの数値シミュレーションを行うための計算格子生成法を検討する。 特に、前処理の手順の確立、及び 格子解像度とプロペラ性能の予測精度との関係を検討する。 そこで今度の研究の目的です。 シャフトやボスを含んだプロペラ幾何形状は複雑で、CADで完全に三次元形状を再現することには困難があるんですが、これを何とかクリヤーしよう、というのと、 舶用プロペラ周り流れの数値シミュレーションを行うための計算格子生成法を検討。という目標です。 特に、前処理の手順の確立、及び 格子解像度とプロペラ性能の予測精度との関係を検討します。 2003/02/19（水）

従来の研究との比較 プロペラ周り流れへのCFDの適用例として以下が有る。 宇都ら（１９９２） 舩野（２００１） しかし、構造格子や部分的構造格子を用いているため、 格子生成が困難であった。 そこで、本研究では、複雑形状まわりの格子生成が 容易である完全四面体格子を採用した。 さらに、解像度が必要なところに格子を集めるために、 「解適合格子法」を導入した。 かつてCFDでプロペラ周り流れを計算した例には、1992年の宇都らや2001年のふねのの研究などがあるんですけれども、構造格子や部分的構造格子を用いているので、格子生成が困難でした。 そこで本研究では、複雑な形状まわりの格子生成が容易である完全四面体格子を採用しました。さらに、解像度が必要なところに格子を集めるために、「解適合格子法」を導入しました。 2003/02/19（水）

解適合格子とは 計算結果を利用して、解像度向上の必要な場所に充分細かい格子が配置されるように、分割し直す格子。 翼端渦などは、細かく分割すべき場所を事前に予測しにくいので、解適合格子は特に有効。 解適合格子というのは、計算結果を利用しまして、解像度を上げたい場所に充分細かい格子が配置されるように、分割し直す格子のことです。 翼端渦などは、細かく分割すべき場所を事前に予測しにくいので、解適合格子はとても有効です。 2003/02/19（水）

当発表の流れ １．検証実験 ２．数値計算 ３．まとめ この発表の流れは、以下の様になります。 ３．まとめ　 この発表の流れは、以下の様になります。 　数値計算結果の評価の基準とするための舶用プロペラまわり流れの実験、 　検証実験と同じ舶用プロペラのまわり流れを、標準的条件と改良した条件とで数値計算、 それぞれの結果と実験値と比較し、その改良法を検討、そして、まとめです。 2003/02/19（水）

検証実験 キャビテーションタンネルにて プロペラ単独性能試験 キャビテーションの観察 主流 主流 2003/02/19（水） まずは検証実験ということで、キャビテーションタンネルにて　プロペラ単独性能試験　と　キャビテーションの観察　を行いました。 この写真では、水は左から右へ流れています。 主流 2003/02/19（水）

プロペラ模型の主要目 翼数 4 直径D 0.2143 (m) ボス比 0.1800 ピッチ比(0.7R) 0.8493 展開面積比 0.6000 スキューバック 0.01257 (m) レーキ 10.00(deg.) 実験に使ったプロペラは、ご覧のようなコンベンショナルプロペラです。のちにご説明する数値計算では、前処理の段階でこの形状を計算機上に再現することになります。 ConventionalPropeller 2003/02/19（水）

実験の条件 一様流中でプロペラを回転させ、圧力一定の下、流速を変化させて推力・トルクを計測 プロペラ回転数n=25[rps] プロペラレイノルズ数Rn= =1.1×10６ タンネル減圧により、 σn= =1.0まで7通り それぞれのにσnについて、流速を系統的に変化 流速は、前進係数Jが0.05ずつ増すように調節 一様流中でプロペラを回転させまして、圧力一定の下流速を変化させて推力・トルクを計測しました。回転数は毎秒25回転、プロペラレイノルズ数は1.1×10６ です。圧力はタンネル減圧によってノン・キャビテーション状態から σn＝1.0まで7通りに変化させまして、それぞれのσnついて、流速を系統的即ち前進係数Jが0.05ずつ増すように変化させました。 2003/02/19（水）

実験結果(プロペラ単独性能) 実験の結果、プロペラ単独性能曲線はこのようになりました。黒い曲線が非キャビテーション状態のグラフで、他のカラーの曲線はそれぞれのσｎの値での曲線です。 2003/02/19（水）

sn=1.5のとき J=0.209 J=0.399 J=0.609 キャビテーションが発生しているが、推力低下はわずかである。 例としてσｎ＝1.5の場合を見てみます。キャビテーションは発生していますが、推力はあまり下がっていない事が分かります。 J=0.209 J=0.399 J=0.609 キャビテーションが発生しているが、推力低下はわずかである。 2003/02/19（水）

実験結果のまとめ キャビテーションは推力に影響を及ぼす事がわかるが、キャビテーションを入れた計算は現在、時間がかかる上に計算も合っていない 本研究では、非キャビテーション状態のみを数値計算する。 実験結果のまとめです。キャビテーションが推力に及ぼす影響はあるにはありますが、キャビテーションを入れた計算は時間がかかる上に計算も合っていないという現状ですし、実験してみたら僅かだった、ということで、本研究では、非キャビテーション状態のみを数値計算します。 2003/02/19（水）

数値計算 Pro/Engineer 2001 使用ソフトウェア Gambit 2.0 Fluent 6.0 形状モデリング 　（3次元CADソフト） 計算領域の定義と 格子生成 Gambit 2.0 （CFDプリプロセッサ） 数値計算 Fluent 6.0 （汎用熱流体解析ソルバー） それでは、数値計算の説明に移ります。 使ったソフトはご覧の通りです。 Pro/Engineerで形を作って、Gambitで格子を切って、 Fluentで計算します。 2003/02/19（水）

CADに読み込むために直交座標系へ変換する 一般にプロペラ形状は、迎角に沿った螺旋座標系で表されている（オフセットデータ） CADに読み込むために直交座標系へ変換する 普通、プロペラの形というのは、螺旋座標系で表されています。ですので、これをまず直交座標に変換します。 2003/02/19（水）

半径rにおける螺旋座標系(s,n) n 2πa(ｒ) rθ 2πｒ s x 2003/02/19（水） 半径 r における螺旋座標系(s,n)というのは、この図のように、翼後方にｓ、翼上方にｎを録る座標系です。円筒座標系(X、θ)から見た｢傾き」2πｒ分の2πa(ｒ)　つまりｒ分のa(ｒ)となります。 s x 2003/02/19（水）

座標変換 螺旋座標系 ⇒ 円筒座標系 ⇒ 直交座標系 2003/02/19（水） 螺旋座標系　⇒　円筒座標系　⇒　直交座標系 そして座標変換の式をまとめるとこのようになります。この(X、θ)というのは基本的に(s,n)に回転行列をかけたもので、ｘにはレーキがプラスされていてθは（ｒ分の1）倍されています。円筒から直交はお馴染みの式です。 2003/02/19（水）

CAD内で完成した三次元形状 翼・ボス・シャフトを別個に作成 立体の演算により一体化 STEP形式で出力する 以上のような手続で翼と、それにボス・シャフトを別個に作成しまして、 それらを立体の演算により一体化しまして、 STEP形式で出力します。 2003/02/19（水）

計算領域と境界条件の定義 4分の1対称領域 3R 2003/02/19（水） 計算領域は、半径方向に3R,上流下流方向とも5Rずつの円筒形を取りました。Rはプロペラの半径です。流れパターンの同じ領域が軸中心に90度ずつ４つ並ぶと考えられますので、そのうちの一つだけをとって4分の1円筒を計算領域としました。また、実際はプロペラが回るんですけれども、計算机上では座標が回転するという設定にしました。この領域を格子分割します。 3R 2003/02/19（水）

格子分割(翼面) 最小格子の大きさ (翼端) 2003/02/19（水） まずは翼面を格子分割した様子です。翼根や翼端は細かくしました。最も細かい翼端の格子は、一辺の長さが130分のR程度の四面体です。 2003/02/19（水）

格子分割(全体) 最大格子の大きさ (流入面) 2003/02/19（水） そして、領域全体を格子分割した様子です。プロペラから遠い場所ほど粗く分割していまして、その大きさは一辺が5分のR程度です。 2003/02/19（水）

Fluent推奨の計算条件(計算1) 四面体、非構造格子に分割 セル数 約20万 アルゴリズム SIMPLE 乱流モデル 標準k-εモデル 反復回数 数百回 圧力の補完 Standard 離散化のスキーム 1次上流差分 それではいよいよ計算に入るわけですが、まずは｢計算1｣と題しまして、標準といいますか、 Fluentの推奨するデフォルト的な条件で計算します。条件の内容はご覧の表の通りです。 2003/02/19（水）

計算結果(計算1) 誤差大 （J=0.6で推力12%、トルク47%の誤差） 2003/02/19（水）

圧力分布(計算1) 圧力分布はご覧のようになりました。 2003/02/19（水）

計算条件の改良(計算2) 圧力勾配による解適合格子を生成 計算条件を以下のように変更 計算1 計算2 セル数 約20万 約40万 乱流モデル 標準k-ε RNGk-ε 反復回数 数百回 数千回 運動量離散化スキーム 1次上流差分 QUICK (2次精度) 次に、先程の計算1とは条件を変えて計算してみました。まず圧力勾配によって解適合格子を生成しまして、更に計算条件も表の一番右の｢計算2｣のように変更しました。セル数は解適合によって倍の40万、乱流モデルはRNGk-ε、反復回数は収束条件を厳しくして数千回、運動量離散化スキームは2次精度のQUICKです。 2003/02/19（水）

計算結果(計算2) 誤差の低減 に成功 （J=0.6で 推力の誤差12%⇒3% トルクの誤差47%⇒19%） 2003/02/19（水） 　　推力の誤差12%⇒3% 　トルクの誤差47%⇒19%） 計算した結果を先程の｢計算1｣のグラフと重ねますと、ご覧のグラフの赤いプロットのようになりました。J=0.6で見てみますと、推力が12%から3%へ、トルクが47%から19%へと、誤差を減らす事に成功したことが分かります。 2003/02/19（水）

翼端渦(J＝0.6) 計算1 計算2 翼端渦の解像精度向上 2003/02/19（水） 計算1 計算2 最後は翼端渦です。計算1に比べて計算2の方が、翼端で渦を巻いている様子が分かりますので、翼端渦の解像精度が向上したと言えます。 翼端渦の解像精度向上 2003/02/19（水）

まとめ① CFD計算精度の検証・向上の為の検証データを実験により取得 オフセットデータで与えられるプロペラ形状を3次元ソリッドモデルに変換する手順を明確化 まとめです。まず、CFD計算精度の検証・向上の為の検証データを、実験により取得しました。 つぎに、オフセットデータで与えられるプロペラ形状を3次元ソリッドモデルに変換する手順を明確化しました。そして、 2003/02/19（水）

まとめ② Fluentの推奨する条件を用い、Gambitで格子生成、Fluentにより数値計算 推力係数 12%の誤差 推力係数 12%の誤差 トルク係数 47%の誤差 Fluentの推奨する条件を用いて格子生成・数値計算したところ、推力係数は12%、トルク係数は47%の誤差を持つ結果を得ました。 2003/02/19（水）

まとめ③ 次に、以下のように改良して計算 格子の解適合分割 乱流モデルをRNGk-εに変更 運動量方程式の対流項のスキームをQUICKに変更 　次に、以下のように改良して計算 格子の解適合分割 乱流モデルをRNGk-εに変更 運動量方程式の対流項のスキームをQUICKに変更 推力係数 誤差12%⇒3% トルク係数 誤差47%⇒19% 　そこで、格子の解適合分割、乱流モデルをRNGk-εに変更、運動量方程式の対流項のスキームをQUICKに変更といった改良を施して計算してみると、推力係数の誤差は12%から3%に、トルク係数の誤差は47%から19%に減少しました。 2003/02/19（水）

今後の課題 誤差を半減できたものの、トルク係数の誤差がまだ19%。格子解像度の不足が原因と考えられる。 格子解像度の不足している場所は、 翼面近傍の境界層（薄い三角柱の格子が必要） 下流側に放出された渦 解適合を解が収束するまで行う 今後の課題です。 誤差を半減できたものの、トルク係数の誤差がまだ19%ありました。格子解像度の不足が原因と考えられます。 格子解像度の不足している場所は、やはり翼面近傍の境界層と、下流側に放出された渦と考えられます。 特に境界層では、薄い三角柱の格子が必要と思われます。 それから、今回1回しか出来なかった解適合を、解が収束するまで行うということも重要だと思います。 2003/02/19（水）