113

日本ガスタービン学会誌

Vol.47

No.2

2019.3

フィルム冷却効率予測のためのベイジアンモデル較正Bayesian Model Calibration for Prediction of Film Cooling Effectiveness三坂 孝志＊1＊4

淺海 典男＊2

出田 武臣＊2MISAKA Takashi

ASAUMI NorioIDETA Takeomi

大林 茂＊3
OBAYASHI Shigeru

ABSTRACT
Film-cooling effectiveness of a slanted round hole in a jet-in-crossÀow con¿guration was predicted by Reynoldsaveraged Navier-Stokes (RANS) simulation with turbulence model parameters optimized based on measurement
data. Computationally ef¿cient approximate models based on Kriging and proper orthogonal decomposition (POD)
were employed for estimating the posterior density distributions of the RANS model parameters using a Markovchain Monte Carlo method, which forms a probabilistic parameter estimation framework based on measurement
data. Maximum a posteriori estimation of the parameters realized better prediction compared to that with the original
parameters. The effect of model parameters was investigated by correlating the change of the thermo-Àuid ¿eld and
POD bases along with estimated POD coef¿cients. In addition, the impact of measurement position and quantity forthe estimation of model parameters were assessed by using the parameter estimation framework.
キーワード：ガスタービン，フィルム冷却，RANS乱流モデル，ベイズパラメータ推定
Key words : Gas Turbine, Film Cooling, RANS Turbulence Model, Bayesian Parameter Estimation．緒言

うな高コストな非定常解析により冷却孔形状や冷却流

ガスタービンエンジンの熱効率向上やCO2削減を実現

条件の最適化を行うことが依然として難しいことから，するためにはタービン入口温度の上昇が必要であり，よフィルム冷却流れのRANS予測精度の向上が必要である。

り高温に耐えうる材料の開発に加えてタービン翼の冷却近年，データ駆動型の乱流モデリング手法として，が重要となる。フィルム冷却は最も効果的なタービン翼

RANS乱流モデルのパラメータおよびモデル式の不確実冷却技術の一つであり，フィルム冷却流れの熱・流体数性に対応するためのアプローチが提案されている。パラ

値解析によるフィルム冷却効率の予測は，冷却孔やそのメータの不確実性に関しては，データ同化などのベイズ配列を設計するための重要な技術である。

推定手法により，特定の流れ場に対して実験計測デー

フィルム冷却効率の予測精度は，一様流中に吹き出すタに基づきパラメータを推定する方法が提案されてい冷却流とその壁面との干渉の解析精度に依存する。こる⑵‑⑸。一方，モデル式の不確実性に関しては，レイノのような流れ場のレイノルズ平均ナビエ・ストーク

ルズ応力のブジネスク近似からのずれを直接数値計算

ス（Reynolds-Averaged Navier-Stokes, RANS）方程式に（DNS）の結果を用いて時間平均速度場の関数として機基づく数値流体解析は，一般に用いられる

械学習によりオフラインで求め，オンラインではRANS

方程式⑴のRANS乱流モデルでは難しいことが知られている 。

平均速度場を参照しつつ学習した上記関係式を用いて計フィルム冷却流れはラージエディシミュレーション

算する方法が提案されている⑹‑⑼。

（Large Eddy Simulation, LES）などの非定常解析によっベイズ推定手法は実験計測データを決定論的な数値解て精度良く予測できることは確認されているが，このよ析と融合するための統計的枠組みとして期待できるものの，乱流モデルパラメータを半ば機械的に推定しただけ

＊
＊
＊
＊

原稿受付 2018年 月16日
査読完了 2019年 月16日
東北大学 学際科学フロンティア研究所〒980-8578 仙台市青葉区荒巻字青葉 ‑㈱IHI〒235-8501 横浜市磯子区新中原町
東北大学 流体科学研究所〒980-8577 仙台市青葉区片平 ‑ ‑現在 国立研究開発法人 産業技術総合研究所
E-mail: takashi.misaka@aist.go.jp

では，そのパラメータの有効性や適用可能範囲に関して
解析条件毎に実験計測値による検証が必要となる。このような制限を緩和するアプローチの一つとして，パラメータの変化と流れ場への影響を関連づけることにより，推定されたパラメータによる解析結果の改善理由を把握する方法が考えられる。
本研究では，RANS乱流モデルのパラメータを冷却効率計測値により較正することでフィルム冷却流れの予測
ー 43 ー

Download service for the GTSJ member of ID , via 157.1.136.4, 2024/10/20.

114

フィルム冷却効率予測のためのベイジアンモデル較正

日本ガスタービン学会誌

精度を向上させることを目的としている。冷却効率の計に対応した近似流れ場を再構築し，任意の位置での冷測値と計算値との差で定義される尤度に基づくサンプリ

却効率を計算することができるため，そのコストが計測

ング計算をマルコフ連鎖モンテカルロ（Markov Chain

点数に依存しない。また，PODによって得られるPOD

Monte Carlo, MCMC）法により行い，パラメータの事

基底ベクトルがパラメータの変化に対する流れ場の変化

後確率分布を得る。高コストなサンプリング計算の負の特徴（主成分）を表すことから，パラメータ変化と流荷を低減するために，Kriging応答曲面や固有直交分解

れ場への影響を把握するのに役立つ。加えて，数値解析（Proper Orthogonal Decomposition, POD）に基づく次結果から擬似的な計測値を作成することで，上記の枠組

元縮約モデルを用いる。通常，MCMC法では数千からみを利用し，どのような計測がパラメータ推定に有効で数万回の尤度評価を必要とするが，本研究では数百程度あるかを検討する。計測位置・計測量を変えた多数回のサンプリング結果から構築された近似モデルにおいての評価は，POD次元縮約モデルを用いることで効率的パラメータ推定を行うことで，計算コストを百分の一程に行う。最終的に，最尤推定パラメータ値（Maximum

度に抑えることが期待できる。本研究では特にRANS乱

A Posteriori解, MAP解）を用いて再度数値解析を行い，流モデルのパラメータと流れ場を関連づけることを目

実験結果やオリジナルのパラメータ値による解析結果と的として，POD次元縮約モデルを用いた推定を検討す比較する。

る。最終的に，最尤推定パラメータ値（最大事後確率推定値）を用いて再度数値解析を行い，実験結果やオリジ
ナルのパラメータ値による解析結果を比較する。加えて，
実験計測値に基づくベイズ推定の枠組みを利用して，計測方法のパラメータ推定精度への影響評価を行う。

．RANS乱流モデルのベイジアンモデル較正
.

ベイジアンモデル較正の流れ

パラメータ推定の流れをFig. 1に示す。まず，修正するパラメータの選定と上下限値の設定を行う。モデルパ
ラメータに理論的な制約条件が存在する場合にはそれに基づき決定する。平板境界層などの特定の流れ場に対して，乱流モデルの輸送方程式から決まるパラメータ値は理論的なパラメータ値の目安となるが，対象とする流れ

Fig. 1 A flow chart of the data-driven probabilistic parameter
estimation procedure

場が大きく異なる場合には大胆な変更を許すというアプ.

ローチも可能である。本研究では実験計画法としてラテ⑽

Kriging応答曲面法

入力パラメータ（乱流モデルパラメータ）と評価値

ン超方格法 を用い，設定した上下限値の範囲内に効率（壁面の実験計測点で定義される冷却効率）の関係を得的にパラメータを生成する。
次に各パラメータに対して数値解析を行い，得られたるための労力が大きな問題に対して，その応答を近似す

流れ場から近似モデルを構築する。近似モデルを通してるために応答曲面法を利用することができる。応答曲面

マルコフ連鎖モンテカルロ法によるパラメータ推定を行法 は 近 年 特 に 数 値 流 体 解 析（Computational Fluid

うことで計算コストの増加を抑制する。近似モデルを作

Dynamics, CFD）を用いた最適化に利用され，その有成する前処理として，各ケースから得られる冷却効率

効性が確認されている⑽。応答曲面として，多項式法やとその実験値の二乗平均平方根誤差（RMSE）を評価し，動径基底関数（Radial Basis Function, RBF）法が知ら

RMSEの小さなケースを利用して近似モデルを作成する。

れているが，実際の応答分布への適合性の良さや新たな本研究では近似モデル作成に関して二通りのアプロー

サンプル点の追加指標が得られる利点から，Kriging法

チを検討する。一つ目は計算された各流れ場において実を用いる。Kriging法では， 個のパラメータからなるベ験計測値に対応する位置での冷却効率を計算し，各流れ

クトル

場の冷却効率とパラメータ値に関するKriging応答曲面

する。

を構築する方法である。このアプローチでは一つの冷却の関数y(p)を以下のように近似⑵

効率計測点に関して一つの応答曲面が必要であり，計測

ここでȝ(p)大域的定数モデルであり，与えられたすべて点数が増える毎に応答曲面構築のコストが増加する。
もう一方のアプローチは，パラメータを変化させたのサンプル点の平均値である。一方，İ(p)は局所的モデときの流れ場データセットに対してPODを適用するこルであり，点pにおけるȝからの偏差に相当する。任意とで次元縮約モデルを構築し，これを近似モデルとしの二点間pi, pjの相関行列の大きさと分散を最大化するて用いる方法である。このアプローチではパラメータ

ような尤度関数を考えると，ȝ(p)およびİ(p)は以下のよー 44 ー
Download service for the GTSJ member of ID , via 157.1.136.4, 2024/10/20.

Vol.47

No.2

2019.3115

フィルム冷却効率予測のためのベイジアンモデル較正

うに表される。

⑼⑶

ここで〈・〉は内積を表す。この係数Įiを乱流モデルパラ

ここでRは相関行列，fは応答ベクトル，1は単位ベクトメータの関数として補間することにより，非計算条件でル，そして，rは相関から得られるベクトルである。

の流れ場を再構築する（POD次元縮約モデル）。本研究.

では放射基底関数（Radial Basis Function, RBF）法に

固有直交分解に基づく次元縮約モデル

PODはデータに含まれる変動を最も良く表現する基
底ベクトル（POD基底ベクトル）を導出する手法であり，よりPOD係数を補間する。すなわち，モードiのPOD係
数は放射基底関数φの線形結合によって，

いわゆる主成分分析である。すなわち，⑽⑷
ここで，Uはデータの含まれる行列，Ȍは直交基底ベクと表される。ここでwi,jは重み係数，pは乱流モデルパラ

トル（POD基底ベクトル）である。CFD結果に対してメータをまとめたベクトル，rjはrj= p-pj

PODを適用する場合には，一定間隔で流れ場を保存し，る距離，pjはサンプルされたケースのパラメータベクト

その流れ変数ベクトルを並べることで以下のように行列Uを構成する（Snapshot POD）⑾。本研究では乱流モデルである。また，基底関数φ(rj)は

で定義され

⑾

ルパラメータを変更して得られる流れ場を以下のようにと定義される。ここでr0はパラメータであり，補間誤差

行列Xにまとめる。

が小さくなるように交差判定等によって定義する。
⑸.

マルコフ連鎖モンテカルロ法
マルコフ連鎖モンテカルロ（Markov Chain Monte

ここで，添え字nはCFDのデータ数（流れ場変数×格子

Carlo, MCMC）法では，あるパラメータ値における実点数），mは乱流モデルパラメータセット数に対応する

験値と計算値の近さを表す尤度Plを評価し，ある割合

T

流れ場の数である。これを利用して行列U Uの固有値問（遷移確率）でパラメータ値を更新することでパラメー

題

タの事後確率分布を得る方法である。
⑹⑿

を解くことで固有値Ȝiと固有ベクトルĳiを得る。これか

ここでȘexpおよびȘcalは冷却効率の実験値および計算値をN個並べたベクトルである。標準偏差σは計測誤差など

らPOD基底ベクトルは⑺と 定義される。 こ の と き， 流 れ 場 の 変 動 に 対 する各本推定で陽にパラメータ化されていない誤差を表現する。
よく用いられるメトロポリス法の処理手順は以下のようになる。

POD基底ベクトル（モード）の寄与率が，固有値Ȝiからわかる。すなわち，主要なPOD基底ベクトルは乱流モ① 初期位置（パラメータ）を設定:

デルパラメータを変化させたときにCFDから得られる② 新しい位置の計算:

流れ場の変動をよく表すように構成されることから，パラメータの変化が流れ場に与える影響を把握するのに適していると考えられる。POD基底ベクトルを用いると元の流れ場は以下のように再構成することができる。

③ 遷移確率の計算:

④ パラメータの更新:⑤ ②に戻って繰り返す⑻
ここでrand（）は[0,1]の乱数である。上記アルゴリズムで
ここで，uは分解される流れ場の変数，uは時間平均成は尤度比Pl (pnew,ınew2)/Pl (pnow,ınow2)が分で空間の関数，Ȍiは式⑺で求めたPOD基底ベクトルが更新されやすくなる。

に近いほどpnow

であり，そして，Įiは基底ベクトルの重ね合わせに用い一般にMCMC法によるパラメータの推定では尤度のる係数である。右辺第二項の和Pは流れ場の近似度に応評価回数が非常に多くなることから計算コストが大きい

じて設定する。パラメータを変化させたときの各流れ場

が，ここではKriging応答曲面およびPOD次元縮約モデに対応したPOD係数は以下のようにして求めることがルを通してMCMC法によるパラメータ探索を行うことできる。

でコストの問題を回避する。尤度の最大化によるパラメータの探索は，勾配法や遺伝的アルゴリズムによるー 45 ー

Download service for the GTSJ member of ID , via 157.1.136.4, 2024/10/20.

116

フィルム冷却効率予測のためのベイジアンモデル較正

日本ガスタービン学会誌

最適化計算と類似であるが，最適化計算においてはパラメータの最適値を求めているのに対して，MCMC法では最適値周辺を広く探索することにより，パラメータ値の確率分布を得ることができる。
.

サポートベクターマシン

近似モデルの構築においては，流れ場が振動するな
ど，実際の流れ場から大きく乖離しているような非物理的ケースを除くことで，近似モデルの精度を向上させると共にMCMCによるパラメータの推定範囲を適切に設定する。本研究ではRayらと同様にサポートベクター
マシン（Support Vector Machine, SVM）を用いた⑷,⑿。
まず，各ケースの冷却効率のRMSEを評価し，RMSE

Fig. 2 A computational domain for a jet-in-crossﬂow
conﬁgurationの小さなケースと大きなケースを分ける閾値を設定し，
RMSEの小さなケースを，大きなケースをとラベル付けする。本研究ではCFD解析を行った100ケースのう
ちRMSEの小さな50ケースを近似モデル構築に利用した。各モデルパラメータと上記ラベルの組に対して，二値SVMを適用することにより任意のパラメータ値入力に対する‑判定を行い，前述のMCMC推定において，と判定されたケースを棄却する。これはMCMC推定においてPl(p,ı2)→Pl(p,ı2)Pp(p)と置き換えることで実現することができる。ここでPp(p)はSVMの判定により以下のように定義される。
⒀
Fig. 3 A close-up view near a cooling hole exit of thecomputational mesh

．フィルム冷却流れの解析.

.

流体計算手法SST乱流モデルと推定パラメータの設定

流体解析にはANSYS社の流体解析ソルバー Fluent

本研究ではMenterによって提案されたSST k-ω乱流

16.2を 次元・単精度で利用した。Semi-Implicit Method

モデルのパラメータに関して検討を行った。SST k-ω乱for Pressure Linked Equations（SIMPLE）アルゴリズ

流モデルは内層と外層でモデル係数を変化させることにムによる圧縮流体の定常解析を行った。ナビエ・ストーよりモデルを切り替えているため，内層および外層にそクス方程式の運動量方程式，エネルギー式および乱流モ

れぞれモデルパラメータが存在する。ここではTable 1デル方程式にはセルベース最小二乗法による勾配を使っに示すパラメータを可変とし，計測データに基づき推た定することとした。パラメータの定義はMenter⒀ に従には

次精度MUSCL，密度には次精度風上差分，圧力次精度のスキームを利用している。時間積分はう。パラメータ変更後に流れ場の収束に必要な反復回数

次精度の陰解法である。乱流モデルにはSST k-ω乱流モの事前検討を行い，本研究ではすべてのケースで2000デルを利用した。Fig. 2およびにそれぞれ計算領域と回の反復を行った後に冷却効率を評価することとした。

計算格子（冷却孔付近）を示す。要素数は6,424,800，総Bergmannらはパラメータの可動範囲として，SST k-ω

節点数は6,586,780である。主流部の流入境界条件は流乱流モデルのモデルパラメータ値の半分および速20 m/s，温度298 K，流出部は圧力出口境界条件としを最小および最大値としている⒁。本研究ではより広範た。冷却流は温度188 K（密度比1.6），質量流量2.98x10‑3

囲に設定し，解析結果の妥当性や計算可能性をもとにパ

kg/sの境界条件（運動量比1.0）を与えた。主流部のス

ラメータの範囲を絞っていくことにした。比散逸率ωのパン方向には周期境界条件を適用している。冷却効果を

輸送方程式の生成項に含まれる係数γ1およびγ2は以下示す指標として本論文では式⒁で定義される冷却効率η

の式⒂から決定されている．本研究では生成項が負になを用いる。

らないという制約条件を与え，γ1またはγ2が負になっ⒁倍の値た場合にそのケースを除外した．また，本流れ場において，パラメータa1が0.25より小さい場合に流れ場のはくここでTmain, TjetおよびTwallはそれぞれ主流温度，冷却流

離・非定常化が著しくなることが事前検討から明らか

温度および冷却された壁面の温度である。

だったため，0.25以上に制限した。
ー 46 ー

Download service for the GTSJ member of ID , via 157.1.136.4, 2024/10/20.

Vol.47

No.2

2019.3

フィルム冷却効率予測のためのベイジアンモデル較正

⒂

117

ことから，Table 1のパラメータ範囲では冷却効率の小さなケースを表現できていないことがわかる。SVMに今回扱ったフィルム冷却流れの条件では，Table 1のよって非物理的なケースを除いた結果をFig. 4 ⒝に示す。

上下限値の設定でも計算が発散することはなかったが，ここではRMSEの小さな上位50ケースを残しており，冷流れが非定常になり定常流れ場に収束しないケースは散却効率が実験値から特に逸脱したケースが除かれている。

見された。

本研究ではFig. 4 ⒝に示すケースを用いて近似モデルを構築する。
Fig. 5にパラメータの事前分布を示す。横軸がパラTable 1 Parameter values considered in the estimation
β1

β2
σω1
σω2
σk1
σk2
a1
β*
Prt

Original
0.075
0.0828
0.5
0.856
0.85
1.0
0.31
0.09
0.85

Minimum
0.01

0.01
0.2
0.2
0.2
0.2
0.25
0.01
0.1

メータ値，縦軸が密度分布を示している。上記のようMaximum
0.5
0.5
10.0
10.0
10.0
10.0
1.5
0.5
1.5

に100ケースのうちRMSEの小さな50ケースとそれ以外の50ケースに分類して学習したSVMから得られた式⒀を用いて，MCMCよって推定したパラメータの確率分
布を事前分布として示している。Fig. 5の破線はSST乱流モデルの標準パラメータ値を示している。σk1および(a)

フィルム冷却流れでは主流の境界層厚さが冷却効率の評価結果に影響するため，実験条件⒂に境界層プロファ
イルを合わせるための流入条件の調整を行った。計算にはSST k-ω乱流モデルのモデルパラメータ値を用いた。
Table 2に境界層厚さの実験との比較を示す。この比較
では冷却流流量をゼロとし，冷却孔中心軸上のいくつかの点で排除厚さδ1，運動量厚さδ2および形状係数Hの評価を行った。流入境界条件として，x/d=‑2.0位置での境界層厚さを実験値に合わせるために，流入境界条件として厚みは0.21Dの層流境界層プロファイルを設定した。
kおよびωは主流乱れ度を実験条件に合わせて一様流入(b)とした。Table 2から下流では境界層厚さの発達が実験値と比較して小さいことがわかる。また，形状係数から実験，計算共に発達した乱流境界層となっていることが確認できる。
Table 2 Boundary layer thickness for the case without coolingﬂow in comparison with the experiment

/
‑2.0
6.5
15.0
23.0

δ1/
CFD Exp(15)
0.088
0.089
0.124
0.120
0.139
0.143
0.154
0.167

δ 2/
CFD

0.057
0.080
0.094
0.106

Exp(15)
0.059
0.083
0.098
0.115

CFD
1.55
1.56
1.49
1.46

Exp(15)
1.50
1.47
1.46
1.44

Fig. 4 Spreads of centerline ﬁlm cooling eﬀectiveness,⒜ simulated cases and ⒝ those after excluding
unphysical cases．パラメータのベイズ推定.

実験値による絞り込み

Fig. 4 ⒜に100ケースの解析から得られた冷却孔中心軸に沿った冷却効率の分布を示す。Table 1に示した範囲のパラメータ変化によって冷却効率が大きく変化する
ことがわかる。特に流れに沿って振動している分布は流れ場自体が振動しているケースである。Fig. 4 ⒜には同条件で得られた実験値を黒丸で示している ⒃,⒄。x/d=5以降では冷却効率が実験値よりも小さなケースが少ない
ー 47 ー
Download service for the GTSJ member of ID , via 157.1.136.4, 2024/10/20.

Fig. 5 Prior density distributions of model parameters

118

フィルム冷却効率予測のためのベイジアンモデル較正

日本ガスタービン学会誌

σk2を除くパラメータは標準値よりも大きな値にするこす。パラメータを機械的に振った結果では，流れ場の振とで冷却効率のRMSEが小さくなることがわかる。また，動など非物理的なケースも含まれるが，冷却効率の実験

分布のピークがパラメータの上下限値に張り付くことも値からの誤差が小さい推定ケースでは，SST k-ω乱流モ

ないため，設定したパラメータの上下限値が妥当である

デルのパラメータで解析されたに近い流れ場が得られてことが確かめられる。

いることが確認できる。

.

Kriging応答曲面法によるパラメータ推定結果

Fig. 6に冷却効率の計測値を用いてMCMC法によって

(a)

[m/s]

(b)

[m/s]

推定したモデルパラメータの事後分布を示す。尤度を評価するための実験値として，Fig. 4に丸印で示す中心軸
上およびスパン平均された冷却効率を用いた。Fig. 5ではRMSEが比較的小さなケースを用いて事前分布としていたため緩やかな密度分布となっているが，Fig. 6の事
後分布ではRMSEの特に小さくなるパラメータ値に対して尤度が高くなるため急峻な密度分布となっている。また，事前分布のピークが事後分布のピークと必ずしも一致していないことが確かめられる。Fig. 6には尤度関数で定義した標準偏差σの確率分布も示している。標準偏Fig. 7 Streamwise velocity distribution on a central slice,⒜ SST original, ⒝ SST estimated

差の値は0.21であり，式⑿に示す尤度関数の定義ではモデル式の不確実性や計測値の不確実性がこの標準偏差として表現されていることになる。σω1およびa1は多峰分布となっており，このような分布はアンサンブルカルマ(a)[deg C](b)[deg C]

ンフィルターや粒子数の少ない粒子フィルターでは推定の困難な分布である。Table 3にパラメータの密度分布のピーク値から得た推定値を示す。ここで得られた推定値は一般的なモデル係数と大きく異なっており，特に乱流エネルギーの散逸率に関わるパラメータβ*は倍以上の値になっている。このパラメータセットが乱流場予測に及ぼす影響については，慎重な検討が必要である。
Fig. 7およびにSST k-ω乱流モデルおよび推定されたパラメータで得られた主流速度および温度の分布を示
Fig. 8 Temperature distribution on a central slice, ⒜ SST
original, ⒝ SST estimatedFig. 9に冷却孔前後の主流方向速度分布の比較を示す。
推定されたモデルパラメータを用いた結果においても，冷却孔直前の境界層速度分布は大きく変わらないことが確かめられる。一方で，冷却孔直後の速度分布ではz/
D=1より小さい領域で減速されていることがわかる。こFig. 6 Posterior density distribution inferred by the measured
ﬁlm-cooling eﬀectiveness

β2
σω1
σω2
σk1
σk2
a1
β*
Prt

Original
0.075
0.0828
0.5
0.856
0.85
1.0
0.31
0.09
0.85

メータに対する推定パラメータの冷却効率の低下に対応していると考えられる。

Table 3 Estimated SST model parameters
β1の冷却流を含む流れの減速がSST k-ω乱流モデルのパラ推定されたパラメータを用いて冷却効率を評価した

Minimum
0.1939

結果をFig. 10に示す。中心線上の冷却効率はSSTモデル

0.27880.52762.63440.25540.22820.80560.3023の標準パラメータ値を用いた結果と比較して，推定したモデルパラメータを用いることで大きく実験値に近づく
ことがわかる。スパン平均した冷却効率に関しては，x/d=5より小さい位置では推定パラメータの結果が実験値に近く，それ以降では実験値の冷却効率の上昇を追従できていないことがわかる。一方で中心線の冷却効率は下1.0614

流側で計測値よりも大きくなっていることから，実験値ー 48 ー
Download service for the GTSJ member of ID , via 157.1.136.4, 2024/10/20.

Vol.47

No.2

2019.3

1.2

1.2

SST original
SST esmated

1.0

示すKriging応答曲面の結果と同様であった。MAPパラ

SST original
SST esmated

1.0

メータに対応するモード

u/U

0.6

までのPOD係数の確率分布をFig. 11に示す。また，確率が最大となるPOD

0.6

係数をTable 4に示す。これらはFig. 11の横軸の数値に0.4対応している。この0.2

捉えている。

0.4
0.2

から

0.8
0.8

u/U

119

フィルム冷却効率予測のためのベイジアンモデル較正

モードで流れ場の変動の約75%を

0.0

(a)

0.0
1.E-02

1.E+00

1.E+02-0.2
1.E-02z/D(b)
1.E+00

1.E+02z/D

Fig. 9 Streamwise velocity proﬁles at ⒜ x/D = -2.0, ⒝ x/D = 2.0.
1.0

(a)

0.9

SST esmated

0.8

SST original

0.7
centerlin η

Fig. 11 Posterior density distribution of POD coeﬃcients from
1st to 4th mode inferred by the experiment

Table 4 Maximum a posteriori estimates of POD coefficientsfrom 1st to 4th mode

Exp (Sinha)

0.6

Mode #
α1
α2
α3
α4

0.5
0.4
0.3
0.2
0.1

Estimated value
507.507‑1333.959
793.065‑655.168

0.0

0

5

10

15

20

Fig. 12に壁面温度の平均，

x/D
0.4

(b)

0.3

average η

次および

次モードを

示す。ここで，これらの主要PODモードは壁面温度の
SST esmated

主要な変化の仕方に対応している。これらの図では温度

SST originalの高低を濃淡で示している。壁面温度は平均および各Exp (Sinha)

モードの重ね合わせとなるので，次モードに関して冷却孔下流の中心軸上を考えると，冷却孔付近で壁面温度0.2を下げ，下流側で温度を上げるようなモードになっていることがわかる（

0.1

次モードの係数はTable 4に示す

ように正の値）。平均にこの

次モードが加わることで，Fig. 10 ⒜に示すように冷却効率が実験値に近づく。ま
0.0
0

5

10

15

た，Fig. 13では冷却流によって発生する縦渦をスパン

20

x/D

方向速度で示しているが，上記のように平均とFig. 10 Film-cooling eﬀectiveness, ⒜ centerline, ⒝ span-average

ドの足し合わせを考えることで，次モー

次モードが渦対の高さを低くするように働いていることがわかる。渦対の高によって最適化したパラメータを用いたとしても，SST
k-ω乱流モデルでは冷却流の広がりを表現するのが難し

さが低くなることで，冷却域がスパン方向に広がると考えられ，温度モードの結果とも矛盾しない。

いことがわかる。これはレイノルズ応力の線形近似が影響していると考えられ，非線形の構造式やレイノルズ応力の輸送方程式を解くレイノルズ応力モデルを用いるこ

averageとで改善する可能性がある⑹。
これらのパラメータの妥当性の評価は一般には難しい問題であると考えられる。次節では，パラメータの変化
が流れ場に与える影響を把握しつつパラメータ推定を行

mode 1
Lower temperatureう手法として，POD次元縮約モデルの可能性を検討する。
.

POD次元縮約モデルによる推定結果の解釈

mode 2

POD次元縮約モデルの作成には25次までのPOD基底を用いた（P=25）。これはパラメータ変更時の流れ場の変化の99%を近似するように選んだ。POD次元縮約モデ
ルによる乱流モデルパラメータの推定結果は，Fig. 6にー 49 ー
Download service for the GTSJ member of ID , via 157.1.136.4, 2024/10/20.

Fig. 12 Temperature average, 1st- and 2nd-modes on wall
surface.

120

フィルム冷却効率予測のためのベイジアンモデル較正

average

日本ガスタービン学会誌

mode 1wall
x/d=10sectionexit
hole

Reducevortex heightx/d=5
x/d=0 section
sectionFig. 13 Spanwise velocity average and 1st mode on several
cross-sections

Fig. 14 The measurement positions considered in the
numerical experiment．パラメータ推定における計測位置の影響
データ同化においては，一般にどの位置でどのような計測を行うかによって状態推定の精度が変わる。本研究(a)で行うパラメータ推定にも同様のことが言えるため，計(b)

測位置・計測物理量の違いがパラメータ推定にどのような影響を与えるかを検討した。
任意の実験計測データを得るのは現実的ではないことから，ここでは基準となるパラメータ値の計算結果を擬
似的な計測値として数値実験を行う。これにより，基準パラメータで得られた基準流れ場から任意の位置で物理量を抽出して，擬似的な計測値として扱うことができる。
計測位置はFig. 14に示すように，冷却孔出口面，x/D=0,5および10の断面，冷却孔下流の壁面，そして，Fig. 10で用いた実験値に相当する中心軸およびスパン平均のFig. 15 The impact of the measurement positions on the
reduction of estimation errors

物理量である。加えて，x/D=5および20の位置における壁面温度のスパン方向分布をx/d=5 spanおよびx/d=20
spanとしてFig. 15に示している。これらの計測位置にラメータ推定誤差が小さいことがわかる。一方で，下流おいて温度，よどみ圧力および流速（それぞれ，Fig.

断面での圧力の計測結果はパラメータ推定には向かない

15横軸のtemp，ptotおよびvelcに対応）を抽出し，疑似ことが確認できる。

的な計測値とした。

．結言

上記のようにして得られた基準流れ場，疑似計測値および基準パラメータを用いて，計測の評価は以下のよう

本研究では，RANS解析によるフィルム冷却効率の予にして行われる。すなわち，式⑿でȘexpを疑似計測値Șref

測精度向上を実現するために，RANS乱流モデルパラに置き換えた式⒃で表される尤度関数を最大化するパラメータを実験値に基づいて較正する手法を検討した。近メータpを推定する。

似モデルとしてKriging応答曲面法およびPOD次元縮⒃約モデルを利用することで尤度評価の計算コストを下げ，MCMC法によってパラメータの事後確率分布を得

このとき，Șrefは基準流れ場u(pref)から得られた疑似計測

た。予測されたパラメータの最大事後確率推定値を用い値，Șcalはサンプルされた流れ場u(pcal)から得られた疑た検証解析から，実験値に近いフィルム冷却効率が得ら

似計測値に対応する。

れることを確認した。ただし，推定値は一般的なモデル⒄

係数と大きく異なっており妥当性についてはさらなる検討が必要である。また，POD次元縮約モデルによるパ

ここでhは観測演算子である。このようにして推定されラメータ推定から，PODモードと温度場・流れ場の修たモデルパラメータpと基準流れ場に対応する基準パラ

正との関連づけを行った。そして，パラメータ推定の枠

メータprefのRMSEを式⒅のように求めることで計測方組みを利用して計測方法のパラメータ推定精度への影響

法の善し悪しを評価する。

評価を行った。一方で，近似モデルの精度はパラメータ⒅

推定の精度に直接影響するため，十分に精度検証を行う必要がある。
本研究はフィルム冷却流れにおける冷却効率の予測精Fig. 15 ⒜および⒝に計測位置を変えたときのパラメータβ1およびβ2の推定誤差Emを示す。相対的な比較度向上を目的としており，特にフィルム冷却流れの条件を濃淡で示しており，色の濃い部分で推定誤差が小さい（密度比，運動量比など）や冷却孔形状の変化に対してことになる。結果から流出孔で温度を計測した場合のパ

汎用的な予測性能を持たせることを目指している。今後

ー 50 ー
Download service for the GTSJ member of ID , via 157.1.136.4, 2024/10/20.

Vol.47

No.2

2019.3

フィルム冷却効率予測のためのベイジアンモデル較正

の課題としては，推定されたパラメータセットの適用範⑼

Wang, J.-X., Wu, J.-L. and Xiao, H., Physics-Informed
Machine Learning for Predictive Turbulence Modeling:
Using Data to Improve RANS Modeled Reynolds
Stresses, arXiv:1606.07987, 2016.

⑽

Forrester, A., Sobester, A. and Keane, A., Engineering
Design via Surrogate Modelling: A Practical Guide,（2008）, Wiley.

⑾

Sirovich, L., Turbulence and the Dynamics ofCoherent Structures I - Coherent Structures, Vol. 45

囲を明確にすることが挙げられる。

参考文献⑴

Harrison, K. L. and Bogard, D. G.,

Comparison of

RANS Turbulence Models for Prediction of Film
Cooling Performance, Proceedings of ASME Turbo
Expo 2008, GT2008-51423（2008）.
⑵⑶

Kato, H., Ishiko, K. and Yoshizawa, A., Optimization of
Parameter Values in the Turbulence Model Aided by
Data Assimilation, AIAA Journal, Vol. 54, No. 5（2016）,
pp. 1512-1523.
Guillasa, S., Gloverb, N. and Malki-Epshteinb, L.,
Bayesian Calibration of the Constants of the k-ε
Turbulence Model for a CFD Model of Street Canyon

⑷

⑸

⑹

Flow, Computer Methods in Applied Mechanics andEngineering, Vol. 279（2014）, pp. 536-553.
Ray, J., Lefantzi, S., Arunajatesan, S. and Dechant, L.,
Bayesian Parameter Estimation of a k-ɽ Model for
Accurate Jet-in-Crossﬂow Simulations, AIAA Journal,
Vol. 54, No. 8（2016）,pp. 2432-2448.
Duraisamy, K., Zhang, Z. J. and Singh, A. P., New
Approaches in Turbulence and Transition Modeling
Using Data-driven Techniques, AIAA Paper 20151284, 2015.
Ling, J. and Templeton, J., Evaluation of Machine
Learning Algorithms for Prediction of Regions ofHigh Reynolds Averaged Navier Stokes Uncertainty,

Physics of Fluids, Vol. 27（2015）,pp. 085103-1-22.
⑺ Ling, J., Ruiz, A., Lacaze, G. and Oefelein, J., Uncertainty
Analysis and Data-Driven Model Advances for a Jet-inCrossﬂow, ASME Journal of Turbomachinery, Vol. 139（2017）, 021008-1.
⑻ Xiao, H., Wu, J. L., Wang, J. X., Sun, R. and Roy, C. J.,
Quantifying and Reducing Model-Form Uncertainties
in Reynolds-Averaged Navier-Stokes Simulations: A⑿

Processing Letters, Vol. 9（1999）, pp. 293-300.
Menter, F. R., Two-Equation Eddy Viscosity⒁

Turbulence Models for Engineering Applications,
AIAA Journal, Vol. 32（1994）, pp. 1598-1605.
Bergmann, C., Ormiston, S. and Chatoorgoon, V.,
Sensitivity Studies of Shear Stress Transport
Turbulence Model Parameters on the Predictionof Seven-Rod Bundle Benchmark Experiments,
Transactions of ASME, Journal of Nuclear Engineeringand Radiation Science, Vol.2（2015）,pp. 011012-1-10.

⒂

Pietrzyk, J. R., Bogard, D. G. and Crawford, M. E.,
Effects of Density Ratio on the Hydrodynamics ofFilm Cooling, Transactions of the ASME, Journal of
Turbomachinery, Vol. 112（1990）, pp. 437-443.
⒃ Kohli, A. and Bogard, D. G., Adiabatic Film
Effectiveness, Thermal Fields, and Velocity Fields forFilm Cooling with Large Angle Injection, Transactionsof the ASME, Journal of Turbomachinery, Vol. 119⒄

115-136.

ー 51 ー

Powered by TCPDF (www.tcpdf.org)（1987）, pp. 561-571.
Suykens, J. A. K. and Vandewalle, J., Least Squares
Support Vector Machine Classifiers, Neural⒀

Data-Driven, Physics-Informed Bayesian Approach,
Journal of Computational Physics, Vol. 324（2016）, pp.

Download service for the GTSJ member of ID , via 157.1.136.4, 2024/10/20.

121（1997）, pp. 352-358.
Sinha, A. K., Bogard, D. G. and Crawford, M. E., Film
Cooling Effectiveness Downstream of a Single Row ofHoles with Variable Density Ratio, Transactions of the
ASME, Journal of Turbomachinery, Vol. 113（1991）,
pp. 442-449