オプティカルフローの適用による高速衝突噴流の画 像解析に関する研究
土黒, 聖斗
九州大学総合理工学府環境エネルギー工学専攻
http://hdl.handle.net/2324/4372241
出版情報：Kyushu University, 2020, 修士, 修士 バージョン：
権利関係：
令 和 ２ 年 度
修 士 論 文
オプティカルフローの適用による 高速衝突噴流の画像解析に関する研究
九州大学大学院総合理工学府 環境エネルギー工学専攻 エネルギー流体科学研究室
土黒 聖斗
指導教員
安養寺 正之
准教授提出年月日 令和３年 2 月 8 日
目次
第1章 背景と目的 ... 1
1.1 高速噴流現象 ... 1
1.2 高速流体現象における既存の計測技術 ... 6
1.3 オプティカルフローを利用した画像解析による流体解析の可能性 ... 7
1.4 研究目的 ... 9
第2章 シャドウグラフ画像に対するオプティカルフローを適用した画像解析原理 ... 10
2.1 シャドウグラフ法及び投影面での運動方程式 ... 10
2.2 オプティカルフロー ... 13
2.3 オプティカルフローのシャドウグラフ画像への適用 ... 15
第3章 実験装置および方法 ... 17
3.1 実験装置 ... 17
3.1.1 全体システム ... 17
3.1.2 二次元先細ノズル ... 19
3.1.3 衝突噴流システム ... 20
3.2 シャドウグラフ光学系 ... 21
3.3 実験条件 ... 24
3.4 画像解析手法 ... 25
3.4.1 解析領域とマスキング処理 ... 25
3.4.2 解析パラメータと解析条件 ... 27
第4章 結果および考察 ... 29
4.1 シャドウグラフ法による可視化結果 ... 29
4.2 解析パラメータの影響 ... 33
4.2.1 ラグランジュ未定定数の影響 ... 33
4.2.2 時空間解像度の影響 ... 35
4.2.3 画像間隔時間の影響 ... 41
4.3 オプティカルフローの適用結果 ... 43
4.3.1 オプティカルフローの定量的妥当性の検証 ... 43
4.3.2 非定常解析結果 ... 44
第5章 結論 ... 54
付録 ... 56
参考文献 ... 57
謝辞 ... 59
第 1 章 背景と目的
1.1 高速噴流現象 (i) 現象論
噴流は，速度をもった流体が小孔（スリット，ノズル）から空間中に噴出する現象である．
噴流現象は，噴流と周囲の流体やノズル形状，すなわち，それら流体の種類（気体，液体，あ るいは固気・気液・固液混相流など），速度やノズル・流路の形状などにより，多種・多様・
多岐にわたる流動形態を示す．また，噴流現象は，その本質が乱流現象を含む自由及び壁面剪 断流れ（流れの剥離，渦の生成などを伴う流れ）で，条件によっては超音速の圧縮性流れや振 動を伴う流れ，気体と液体あるいは固体と気体が混じりあう混相流れを含み，古くから工学 的・工業的に非常に重要で興味深い事象の一つである[1]．
噴流は，前述したように，ノズル形状，噴出速度，作動流体，周囲の状況などによりその流 動形態が大きく異なることから自由噴流，壁面噴流，衝突噴流の 3 つに大別される．ここで は，非圧縮性の各噴流現象における基本的な流動特性および圧縮性を伴う超音速の噴流構造 について整理する．
●自由噴流
二次元ノズルあるいは三次元軸対称ノズルから噴流を無限に広い同一流体の静止空間中に噴 出させた場合，噴流は周囲流体との間の大きな速度勾配と流体の有する粘性の作用により，
周囲の流体を巻き込み，噴流幅を広げながら，また，噴流自身は速度を減少させながら下流方 向に流れていく．自由噴流の概略をFig. 1.1に示す．ここで，x, yはそれぞれ，噴流の軸方向，
及びそれに直角な方向への座標軸，u, vはそれぞれ，x, y方向への平均速度成分である．いま，
幅Dの二次元ノズル，あるいは三次元円形ノズルから一様な速度分布u0で噴出する自由噴流 の平均速度分布はFig. 1.1に示すように，その分布形からノズル近傍とそれより下流の２つの 領域に分けることができる．速度が減少しないポテンシャルコア領域（Potential core region）
を含むノズル近傍を初期領域（Initial region）または遷移領域（Transition region）と呼び，ポテ ンシャルコアよりも下流を発達領域(Developed region)と呼ぶ．この時，遷移領域では剪断層の 不安定性に起因して紙面に対して垂直なスパン方向に軸をもつケルビン・ヘルムホルツ渦（以
下，K-H渦とする）が形成される．発達領域ではその渦同士の合体や融合を経て，2次元的な スパン方向の不安定性により主流方向に軸をもつ縦渦が発達する．また，ポテンシャルコア の終点付近では，噴流固有の不安定モード（preferred mode）がある[2]．
●壁面噴流
壁面噴流は，固体境界（壁面）に沿って噴出される噴流であり，境界層制御，高温壁の膜冷 却などに利用される．Fig. 1.2に二次元壁面噴流と近接平板との間（H / D ＜ 1，H：ノズルか ら壁面までの距離）から噴出する放射状壁面噴流の概略図を示す．一般的にH / D ＞ 1の場 合を衝突噴流と呼ぶ．いずれの場合も平板壁面上に壁面剪断層を，噴流の外側に自由剪断層 を形成する．
●衝突噴流
噴流を壁面に衝突させることにより，熱伝達あるいは乱流混合を促進する手法は，加熱・冷 却・乾燥を必要とする工学分野で広く利用されており，その性能は噴流の渦構造に関連して いることがわかっている．Fig. 1.3に二次元ノズル，あるいは三次元円形ノズルから出た噴流 が平板に垂直に衝突する衝突噴流の概略図を示す．衝突噴流は，ノズル・壁面間距離 Hが十 分に大きい場合には壁面の影響を受けない自由噴流領域（Free jet region），及び衝突噴流領域
（Impinging jet region）と壁面に衝突したあと，噴流が平板に沿って流れる壁面噴流領域（Wall
jet region）の3つの領域に大別される[3]．ここで，自由噴流領域及び壁面噴流領域については
前述とおおむね同様である．衝突噴流領域では，噴流が衝突する動圧的な効果により形成さ れた圧力勾配によって，衝突壁面近傍では噴流速度が急激に減少する[4]．また，壁面近傍での 流れの方向は壁面と垂直な方向から水平方向へと急激に変化するため，衝突噴流領域では非 回転な変形速度場で，乱れの非等方性も著しく強くなる．
Fig. 1.1 自由噴流の概略図[1]
(a) 二次元壁面噴流
(b) 放射状壁面噴流（H / D ＜ 1）
Fig. 1. 2 壁面噴流の概略図[1]
Fig. 1.3 衝突噴流の概略図[1]
⚫ 超音速における不足膨張流れ
Fig. 1.4 (a)は，超音速自由噴流における不足膨張流れの概略図を示す[5]．ノズル出口下流では，
噴流と周囲の空気の境界面に混合領域（Mixing region）があり，次に凸状のコア境界を示す．
また，噴流は，中心線に沿って超音速噴流の流れが伝播し，ノズル形状に起因した斜め衝撃波
（Oblique shock）の発生や，マッハ衝撃波ディスク（Mack shock disk）を形成する．斜め衝撃 波は衝撃波の連続反射である反射衝撃波（Reflected shock）をもたらし，マッハディスクの下
流側では，温度，圧力及び密度が急上昇し亜音速流れとなる．また，亜音速と超音速の境界面 は滑り面（Slip line）を形成する．
Fig. 1.4 (b)は，超音速衝突噴流における衝突噴流領域での流れの概略図を示す[6]．Fig. 1.4 (b)
に示すように，衝突噴流領域では，平面衝撃波（Plate shock）やテール衝撃波（Tail shock）な どの衝撃波構造が見られる．噴流が壁面に近づくと，平面衝撃波による速度の大幅な減少が 発生し，結果として、衝突噴流領域あるいはよどみ領域（Stagnation region）では、再循環流れ
（Recirculating flow）が現れる．
(a) 自由噴流[5]
(b) 衝突噴流[6]
Fig. 1. 4 不足膨張噴流の概略図
(ii) 工学的利用
噴流現象を伴う工学的アプリケーションには，飛翔体の推進装置[7]，衝突噴流による高温壁
の冷却[8][9]，高速水噴流によるジェットカッティング[10]などがあり，多岐にわたって応用され
ている．
これら噴流を伴う工学的アプリケーションの発展には，ノズルの形状，ノズル出口速度，噴 流と周囲の流体の種類や周囲の流路形状により多種多様な流動形態をとる噴流の挙動を明ら かにするために，噴流の速度や乱れ分布，圧力分布などの空間的，時間的な様相を把握する必 要がある．そのために，これらのアプリケーションを対象として，実験的計測手法の確立がな されてきた．次節では，従来の代表的な流体計測手法の紹介とそれらが抱える課題について 提示する．
1.2 高速流体現象における既存の計測技術
従来，流れ場の代表的な計測手法に用いられてきたものとして，ピトー管[11]や熱線流速計
[11]，レーザードップラー流速計が挙げられる．ピトー管は，管に設けられた全圧孔と静圧孔 から圧力差を計測することにより間接的に流速を計測する手法である．廉価であり，使用す る流体に対する制約が少なく，古くから利用されてきた．しかし，ピトー管の精度は使用する 圧力センサーに依存することや，時間分解能が低いといった問題がある．また，熱線流速計 は，流れの中に置かれた熱線の電気抵抗の変化によって発生するブリッジ回路における電位 差を取得し，流速と電圧の関係から流速を計測する手法であり，時間分解能が高く乱流現象 や渦の周波数解析など非定常計測に優れている．しかし，これら接触型の計測手法は流れの 中にセンサーを挿入するために，流れ場が乱され本来持っている特性を乱す恐れがある．特 に，衝撃波や膨張波を伴う高速流れにおいて，接触型は計測値に顕著な誤差を生じさせるこ とが考えられる．また，レーザードップラー流速計は，流れ場中にある固体あるいは液体粒子 をトレーサーとして混入させ，それら粒子の散乱光を観測して流速を計測する手法である．
非接触で，時間分解能に優れているが，衝撃波により流速がステップ上に減少する流れ場で はトレーサー粒子の追従性に問題がある．更に，これらは，点計測であることから流れ場全体 の分布情報を得るにはトラバースさせるなどして計測する必要があるが，いずれも乱流の微 細構造を捉えるほどの高い空間分解能は得られない．
そこで，1990 年代後半から，光学的可視化技術と画像解析技術を組み合わせた非接触流体 計測手法として，PIV（Particle Image velocimetry）[12]やPTV（Particle Tracking Velocimetry）[13]
が登場した．PIVは，流体中にトレーサー粒子を混入し，そこにレーザーシートを照射するこ とによってシート面上のトレーサー粒子の動きをカメラで撮影する．更に，撮影された連続 画像に対して，ある微小時間間隔 Δt におけるトレーサー粒子の変位 Δx を相関法によって処 理することで，計測面における流体の局所速度ベクトルを算出することができる．また，近年 では，ハイスピードカメラやレーザー光源などのハードウェア面の性能向上により高い時間 解像度と空間分解能による計測が可能となってきた．一方で，トレーサー粒子が入り込みに くい剥離剪断層内部や圧縮性流れにおける衝撃波や膨張波前後の速度が不連続な領域などト レーサーの追従が困難な領域では可視化が困難であるといったトレーサー粒子に起因した課 題が残る．Table 1-1 には，高速流体現象における代表的な既存の計測技術の特徴（利点・欠 点）について簡単にまとめている．
Table 1-1 代表的な流体計測手法の特徴
計測手法 利点 欠点
ピトー管 ⚫ 廉価で運用上の制約が少ない ⚫ 時間分解能が低い 熱線流速計 ⚫ 時間分解能が高い ⚫ 接触型で流れを乱す レーザードップラー流速計 ⚫ 非接触で時間分解能が高い ⚫ 流れに対する粒子の追従性
PIV
⚫ 非接触で流体の2D,3Dの計測 が可能
⚫ 流れに対する粒子の追従性
⚫ 計測機器が高価
1.3 オプティカルフローを利用した画像解析による流体解析の可能性
前節では，これまでに確立された流体計測手法について述べた．特に，PIVやPTVといっ た 2 次元あるいは 3 次元での定量的データが取得できる計測手法は，様々な流動形態を示す 噴流現象を理解する上で有用である．しかし，PIVやPTVなどの光学的可視化技術と画像解 析を用いた非接触流体計測手法が確立される以前から，特に，衝撃波や膨張波を伴うような 圧縮性流れにおいてはシャドウグラフ法やシュリーレン法といった光学的可視化手法が定性 的ではあるものの密度変化の大きい高速流体現象の理解に有効な手法として確立されている．
これら，様々な光学的可視化手法に対して Liu らはオプティカルフロー法[14-17]を用いて定量 的な画像解析のためのアルゴリズムを考案した[18]．このアルゴリズムの特徴は，撮影画像内 においてピクセルレベルでの解析が可能であり，取得画像分解能と同等の空間分解能で画像 上での速度ベクトルの算出が可能である点が挙げられる．実際，Liuらは，このアルゴリズム
を衝突噴流の粒子画像に適用して，相互相関法を用いた PIV 解析結果との比較をおこなった
[19]（Fig. 1.5）．Fig. 1.5の，ベクトル・コンター図は粒子画像に対するオプティカルフロー法 と相関法それぞれの速度ベクトル及び渦度場の結果を，グラフは，x, y方向それぞれの速度プ ロットを表す．Fig. 1.5より，オプティカルフロー法を用いた解析結果は渦と境界層間の複雑 な相互作用が働く領域において，相互相関法による PIV 解析結果よりも理論解に近い流れ場 を抽出できることを明らかにした．しかし，Liu らの検証では流速が4.5 m/sの低速流れに対 する適用に留まり，圧縮性を伴うような非定常で高速な流れにおいては，その適用の可否は 不明である．また，シャドウグラフ法やシュリーレン法のようなトレーサーを用いない光学 的可視化手法に対してオプティカルフロー法を実際に適用した例はなく，この方法により非 定常な流れの計測が可能かどうかについては学術的興味が非常に大きい．
Fig. 1.5 衝突噴流（流速=4.5 m/s）に対するオプティカルフロー法と相関法の比較[15]
1.4 研究目的
本研究では，亜音速から遷音速領域における衝突噴流現象を対象として，高速で非定常な 流れ場に対するトレーサーを用いない光学的可視化手法においてオプティカルフローを用い た画像解析手法の非定常解析への適用可能性について調べる．
第 2 章 シャドウグラフ画像に対するオプティ カルフローを適用した画像解析原理
2.1 シャドウグラフ法及び投影面での運動方程式
画像平面上における 3 次元空間座標系の輸送方程式や連続方程式をもとに Liu ら[18]は様々 な流れの可視化手法における投影面での運動方程式を導いた．ここでは，その中でも流れの 密度の 2 次微分を可視化するシャドウグラフ法について，その可視化原理とオプティカルフ ロー法で用いる物理モデルに基づいた投影面での運動方程式の導出について示す．
(i) シャドウグラフ法の原理[20]
シャドウグラフ法は媒質の密度変化による屈折率の変化を利用した光学的可視化手法であ る．シャドウグラフの原理をFig. 2.1 に示す．Fig. 2.1が示すようにシャドウグラフ法は，気 体あるいは液体の密度変化による光の影を，直接投影面に投影して観察することから直接投 影法とも呼ばれる．Fig. 2.1より，点光源からの光は，密度が一様でない測定部を通過すると き曲げられ，投影面に明るさの濃淡をもつ像を作る．ここでは，簡単のため，密度変化は2次 元的で紙面に垂直な x 方向には変化しないものとする．点 A を通過する光線 γAの屈折角は，
屈折率 n がグラッドストーン・デイルの関係より n = 1+Kρ（K：グラッドストーンデイル定 数，ρ：媒質の密度）で表される気体に対しては，
𝜃𝑦 = ∫ 𝐾 𝑛
𝜕𝜌
𝜕𝑦
𝑑 0
𝑑𝑧 (2.1)
密度はy方向のみに変化するものとすると，光線γAの全屈折角は，
𝜃𝑦= 𝐾𝑑 𝑛
𝑑𝜌
𝑑𝑦≅ 𝐾𝑑𝑑𝜌
𝑑𝑦 (2.2)
となる．もし密度勾配dρ/dyがy方向に一定ならば，すべての光は一様に曲げられ，投影面上 での像の明るさIは変わらない．したがって，明るさの変化ΔIは密度勾配の変化によって現 れる．すなわち，
∆𝐼 ∝ 𝑑𝜃𝑦
𝑑𝑦 ∝ 𝑑2𝜌
𝑑𝑦2 (2.3)
密度がx，y方向ともに変化する場合は，
∆𝐼 ∝ 𝜕2𝜌
𝜕𝑥2+ 𝜕2𝜌
𝜕𝑦2 (2.4)
の関係が成立する．式（2.4）から，シャドウグラフ法による投影面上における明るさの変化 は，対象とする媒質の密度の2次微分つまり密度勾配の変化に比例することがわかる．
Fig. 2.1 シャドウグラフ原理[2.2]
(ii) 投影面での運動方程式の導出
三次元空間における流体媒体から画像平面への投影法についてFig. 2.2 に示す．Fig. 2.2よ り，x = (x1, x2)は画像座標，X = (X1, X2, X3)はオブジェクト座標で，X2は主流方向，X3は光線に 沿っており，計測面に対して垂直である．計測面はX3 = Γ1 (X1, X2)，X3 = Γ2 (X1, X2) で与えら れる．シャドウグラフ法の原理より，シャドウグラフ画像の輝度 I は流体密度 ρ の 2 次微分 として示される[21]．すなわち，
𝐼 − 𝐼𝑇
𝐼𝑇 = 𝐶 ∫ 𝛁122
Γ2 Γ1
𝜌𝑑𝑋3 (2.5)
ここで，ITは流れなしの画像における輝度，Cはカメラの焦点距離や光源からの放射輝度など の光学系に関わる係数[22]，∇122 はラプラス演算子であり次式で表される．
𝛁122 = 𝜕2
𝜕𝑋12+ 𝜕2
𝜕𝑋22 (2.6)
計測面が光学窓のような場合，X3 = Γ1 = const.，X3 = Γ2 = const.が成り立つ．そして，式（2.5） において時間tで偏微分を行い，式（2.7）に示す連続の式を用いると，式（2.8）を得る．
∂𝜌
∂𝑡 + 𝛁 ∙ (𝜌𝑼) = 0
(2.7)
−1 𝐶
𝜕(𝐼 𝐼⁄ 𝑇− 1)
𝜕𝑡 = 𝛁122[𝛁12∙ ∫ 𝜌𝑼12𝑑𝑋3
Γ2 Γ1
] (2.8)
ここで，U12 = (U1, U2) は流体の座標平面 (X1, X2) に投影された2次元速度成分である．また，
式（2.5）より流体密度における積分のポアソン方程式が得られる．
𝛁122∫ 𝜌 𝑑𝑋3= 𝐶−1
Γ2 Γ1
(𝐼 𝐼⁄𝑇− 1) (2.9)
𝑔 = 𝐶 ∫ 𝜌𝑑𝑋ΓΓ2 3
1 と置き，ポアソン方程式𝛁122𝑔 = 𝐼 𝐼⁄ 𝑇− 1を適切な境界条件で解くことにより，
密度の積分値を得ることができる．したがって，式（2.8）よりシャドウグラフ法における投 影面での運動方程式が得られる．
𝜕𝑔
𝜕𝑡 + 𝐶 ∙ 𝛁12∙ ∫ 𝜌𝑼12𝑑𝑋3
Γ2 Γ1
= 0 (2.10)
ここで，流体密度ρで重み付けされた加重平均速度〈𝑼12〉𝜌を定義すると，
〈𝑼12〉𝜌=∫ 𝜌𝑼12d𝑋3 𝛤2
𝛤1
∫ 𝜌 d𝑋3 𝛤2 𝛤1
(2.11)
式（2.10）は次のように表すことができる．
∂𝑔
∂𝑡+ 𝛁12∙ (𝑔〈𝑼12〉𝜌) = 0 (2.12)
Fig. 2.2 流体から画像平面への投影[18]
2.2 オプティカルフロー
オプティカルフローとは，2次元空間での運動を捉えた連続画像により，物体の画像上にお ける瞬間的な速度または変位を推定する画像解析手法である． 速度ベクトルは画像上の各ピ クセルの輝度値の空間的・時間的変化から求める．
まずは，オプティカルフローの基礎方程式について述べるために，Fig. 2.3のような画像平 面における物体の輝度の移動について考える．Fig. 2.3 に示すように，画像中のある点 (x, y) の時刻tにおける物体の輝度はI (x, y, t) であり，微小時間δt後に物体が (x+δx, y+δy) に移動 したとすると，その点の輝度はI (x + δx, y + δy, t + δt) と表すことができる．この時，物体の 輝度は画像間で変化しないと仮定すると，以下の関係式が成り立つ．
𝐼(𝑥, 𝑦, 𝑡) = 𝐼(𝑥 + 𝛿𝑥, 𝑦 + 𝛿𝑦, 𝑡 + 𝛿𝑡) (2.13)
これは，画像間で移動した物体は位置や時間変化があったとしても輝度値の変化がないこ とを表す．ここで，輝度の時間変化と位置変化が無視できるほど微小であるととして，式（2.13）
の右辺をテイラー展開すると，
𝐼(𝑥 + 𝛿𝑥, 𝑦 + 𝛿𝑦, 𝑡 + 𝛿𝑡) = 𝐼(𝑥, 𝑦, 𝑡) +𝜕𝐼
𝜕𝑥𝛿𝑥 +𝜕𝐼
𝜕𝑦𝛿𝑦 +𝜕𝐼
𝜕𝑡𝛿𝑡 + 𝒪(𝛿2) (2.14)
となる．式（2.14）に式（2.13）を代入し，高次の項の影響は十分小さいとして無視すると，
𝜕𝐼
𝜕𝑥𝛿𝑥 +𝜕𝐼
𝜕𝑦𝛿𝑦 +𝜕𝐼
𝜕𝑡𝛿𝑡 = 0 (2.15)
が得られる．ここで，式（2.15）の両辺をδtで割り，δt → 0とすると，
𝜕𝐼
𝜕𝑥
𝜕𝑥
𝜕𝑡+𝜕𝐼
𝜕𝑦
𝜕𝑦
𝜕𝑡+𝜕𝐼
𝜕𝑡= 0 (2.16)
となり，この変位の時間変化は速度に相当するもので，
𝑢𝜕𝐼
𝜕𝑥+ 𝑣𝜕𝐼
𝜕𝑦+𝜕𝐼
𝜕𝑡 = 0 (2.17)
𝑑𝐼
𝑑𝑡+ 𝒖 ∙ ∇𝐼 = 0
(2.18)
が得られる．ここで，𝜵𝐼 = (𝜕𝐼/𝜕𝑥, 𝜕𝐼/𝜕𝑦)，𝜕𝐼/𝜕𝑡 = 𝑑𝐼/𝑑𝑡であり，𝒖 = (𝑢, 𝑣) = (𝜕𝑥/𝜕𝑡, 𝜕𝑦/
𝜕𝑡)が，オプティカルフローの速度を表す．式（2.18）はオプティカルフローにおける勾配法
[23]の輝度拘束方程式となる．式（2.18）を解くことによって画像間の輝度の差から速度を算出 することができるが，この式は窓問題にあたり，式を解くには拘束条件が足りない．故に，オ プティカルフローを求めるために拘束条件の異なる解析手法がこれまでに多く提案されてき
た[17, 23-27]．次の2.3節では，実際にシャドウグラフ画像の解析を行うにあたり実流体とオプテ
ィカルフローの流れを結びつけるためにHorn-Schunck法を修正したLiuら[18]のアルゴリズム について示す．
Fig. 2.3 画像平面上での物体の輝度の移動に伴うオプティカルフローの概略図
2.3 オプティカルフローのシャドウグラフ画像への適用
先ほど導出した式（2.12）と式（2.18）を比較すると，シャドウグラフ法における輝度Iの関 数であるgが発光量Iと，流体密度ρで重みづけされた加重平均速度〈𝑼12〉𝜌がオプティカルフ ローuに対応する形となっていることがわかる．
ここで，Horn ら[23]の輝度拘束方程式（2.18）は時系列画像からの視覚的な動きのみを決定 することに注意する．式（2.12）と式（2.18）を比較すると，∇ ∙ 𝒖 = 0の場合のみ，式（2.12）
はHorn-Schunck法と同形式となる．しかし，一般的に，∇ ∙ 𝒖 ≠ 0であるが，二次元流れ（U3
= 0）を仮定した場合，𝛁12∙ 〈𝑼12〉𝜌 = 𝛁12∙ 𝑼12 = 0故に，∇ ∙ 𝒖 = 0となるため，Horn-Schunck法 と同様に，オプティカルフローの計算手法を本計測に適応することができる．
また，式（2.18）は窓問題にあたり，オプティカルフローの速度を求めるためには追加の拘 束条件が必要となる．本研究では，空間的に速度ベクトルは滑らかに変化するというHornら によって提案された拘束条件[23]を導入した．これは，オプティカルフローの空間勾配の 2 乗 和，すなわち，
(𝜕𝑢
𝜕𝑥)
2
+ (𝜕𝑢
𝜕𝑦)
2
+ (𝜕𝑣
𝜕𝑥)
2
+ (𝜕𝑣
𝜕𝑦)
2
(2.19)
となる．式（2.12）と式（2.19）の2つの拘束条件が満たされるような汎関数J (u) を式（2.20）
のように定義し，この汎関数の最小化問題を解くことで，速度の2成分を決定していく．
𝐽(𝒖) = ∫ [∂𝑔
∂𝑡+ 𝛁 ∙ (𝑔〈𝑼12〉𝜌)]
2
𝑑𝑥𝑑𝑦 + 𝛼
𝛺
∫ (|𝛁𝑢|2+ |𝛁𝑣|2) 𝑑𝑥𝑑𝑦
𝛺
(2.20)
ここで，Ωは計算を行う画像領域である．αはラグランジュ未定定数を示しており，これは最 適解を得るために 2 つの制約の強さを調整するパラメータである．このパラメータを決定す る根拠は，解析対象となる流れ場の物理的性質に依存する．また，通常この汎関数の最小化問 題を解くには，オイラー・ラグランジュ方程式を導入する．J (u) を最小化するには，任意の 滑らかさ関数𝒗 = (𝑣1, 𝑣2)，および微小量𝑝を与え，ノイマン条件𝜕〈𝑼12 〉𝜌/𝜕𝑛 = 0が画像領域𝜕𝛺 に課せられるグリーンの定理を使用して，𝑑𝐽(〈𝑼12〉𝜌+ 𝑝𝒗) 𝑑𝑝⁄ |
𝑝=0を解くと，オイラー・ラグ ランジュ方程式は以下のようになる．
𝑔∇ [∂𝑔
∂𝑡 + 𝛁 ∙ (𝑔〈𝑼12〉𝜌)] + 𝛼∇2〈𝑼12〉𝜌= 0 (2.20)
また，式（2.20）を解くためにヤコビ反復法を使用する．Horn-Schunck法による方程式の解は，
解の収束を高速化するために式（2.20）の初期近似として使用する．
第 3 章 実験装置および方法
3.1 実験装置
3.1.1 全体システム
本研究に用いた実験装置全体の模式図をFig. 3.1に示す．作動流体には空気を用い，最大圧
力3.3 MPa，内容積5.0 m3の高圧ガスタンク（Fig. 3.2）に最大吐出圧力2.94 MPa，最大吐出容
量0.45 m3/minの圧縮機（田邊空気機械製作所，HC-54A）(Fig. 3.3)を用いて圧縮空気を充填す
る．タンク内の圧縮空気は高圧配管及び手動弁を介してよどみ室を通過後，ノズルから大気 に開放される．供給圧力は手動弁の開閉度によって調整する．ノズル両端には，アクリル製の
300mm四方の側壁が設置されている．側壁は光学計測及び静圧計測に応じて光学ガラスまた
は静圧孔付きの側壁に交換可能になっている．本研究では，光学計測のみを行い側壁の静圧 計測については同システムにおける原田ら[28]の計測結果を参照した．
Fig. 3.4には計測システムのブロック線図を示す．よどみ室内の圧力は圧力センサ（JTEKT，
PMS-5M-2 1M）を用いて計測し，増幅器（JTEKT，DC AMP AA6210）及びフィルター（NF MULTIFUNCTION FILTER 3611）を介してDAQ（National Instruments，NI cDAQ-9178，NI-9215） とオシロスコープ（KEYSIGHT，EDUX1002G）に取り込む．また，高精度ポータブル大気圧
計（GE Druck，DPI740）を用いて大気圧を計測する．通風中は，よどみ圧（p0）と大気圧（pb）
を用いてノズル作動圧力比（NPR : Nozzle Pressure Ratio = p0/pb）をリアルタイムでモニタリン グしながら手動弁の開閉度を調整し，一度目標のNPRをオーバーシュートさせたのちに緩や かに目標の NPR へと調圧する．オシロスコープは目標の NPR に相当する圧力センサの電圧 値をトリガーとしてシングルモードで設定し，目標のNPRに達すると自動的にハイスピード カメラでの撮影が開始される．
Fig. 3.1 実験装置全体模式図
Fig. 3.2 高圧ガスタンク
Fig. 3.3 圧縮機
Fig. 3.4 計測システムブロック線図
3.1.2 二次元先細ノズル
本研究では，亜音速及び遷音速領域における二次元の衝突噴流を生成するために，使用す るノズルは二次元先細ノズルとする（Fig. 3.5）．Fig. 3.5 (a)より，ノズル出口中心を原点とし，
主流方向に x軸を，壁面に沿う方に y 軸を，スパン方向に z軸を設定した．ノズル入口高さ
60 mm，出口高さD = 8 mm，ノズル入口から出口までの距離は75 mm，ノズルリップ厚さは
7.5 mmである．また，流れの二次元性を仮定するためにノズルのスパン長は80 mmで，ノズ
ル出口のアスペクト比は10となっている．Fig. 3.5 (b)より，ノズル形状は (x, y) = (-65, 30)，
(x, y) = (0, 4)と変曲点となる (x, y) = (-32.5, 17) を通る3次関数で近似している．ノズルは光造
形技術によってエポキシ系樹脂で製作された．
-40.0 -30.0 -20.0 -10.0 0.0 10.0 20.0 30.0 40.0
-80 -60 -40 -20 0 20
y [mm]
x [mm]
xx y
(a) ノズル寸法 (単位：mm) (b) ノズル曲線形状
( c) 実物写真
Fig. 3.5 二次元先細ノズル
3.1.3 衝突噴流システム
衝突噴流システムの模式図をFig. 3.6に示す．ノズルから衝突壁までの距離はLiuら[19]の衝 突噴流試験を参考に34 mm とした．本実験で用いた衝突壁は300 mm×80 mm×15 mmの鉄製 プレートの上流側に同形状の厚さ5 mmのアルミプレートを設け，下流側には厚さ50 mmの リブを設けた[28]．アルミプレートの表面粗さは並仕上げで衝突噴流時の壁面の表面粗さが流 れに及ぼす影響を軽減させている．また，リブ構造により壁面の強度は十分に確保され，遷音 速領域での衝突噴流における壁面のたわみは非常に小さい．
Fig. 3.6 衝突噴流システム模式図
3.2 シャドウグラフ光学系
本研究に使用したシャドウグラフ光学系のセットアップの一例を Fig. 3.7 に示す．光源に
は，250 Wメタルハライドランプ（協和株式会社，MID-25FC）（Fig. 3.8，Table 3-1）及び50
Wレーザー光源（CAVITAR，CAVILUX Smart UHS）（Fig. 3.9，Table 3-2）を用いた．メタル ハライドランプは，ハイスピードカメラ用調光機能付きのフリッカーレスで，太陽光に近い
高輝度（6750 K）を有し，集光レンズ及びピンホールを用いて点光源として使用する．レーザ
ー光源は，波長が645 nmの赤色の短波長レーザーで，低コヒーレント，超短パルス，超高繰 り返しで撮影が可能であり，ファイバー径が 1 mm のため，そのまま点光源として利用でき る．レーザー光源の場合は，輝度調整のために ND フィルターをファイバーの目の前に設置 し，ハイスピードカメラとの同期撮影にて使用する．凹面鏡には直径150 mm，焦点距離1000 mm のものを 2 つ用いた．また，連続画像の撮影には，Photron 社製のハイスピードカメラ
（Photron，FASTCAM SA-Z : 以後，SA-Zと省略する）（Fig. 3.10, Table 3-3）及びnac社製の ハイスピードカメラ（nac，MEMRECAM ACS-1 M60 : 以後，ACS-1と省略する）（Fig. 3.11, Table 3-4）のそれぞれにズームレンズ（Nikon，ED AF NIKKOR 80-200 mm）を取り付けて撮 影を行った．
（a）模式図
（b）写真
Fig.3.7 シャドウグラフ光学系セットアップ
Fig. 3.8 メタルハライドランプ Fig. 3.9 レーザー光源
Tabele 3-1 メタルハライドランプ主仕様 Table 3-2 レーザー光源主仕様
製品名 MID-25FC（協和株式会社） 製品名 CAIVLUX Smart UHS
出力 250W（100%～50%連続調光） 波長 645 ±10 nm以内
色温度 6750 K 出力 50 W±10%以内
パルス幅 10 ns～150 ns
Fig. 3.10 FASTACAM SA-Z Fig. 3.11 MEMRECAM ACS-1 M60
Tabele 3-3 FASATCAM SA-Z主仕様 Table 3-4 MEMRECAM ACS-1 M60主仕様 撮像素子 1024×1024 pixel CMOS 撮像素子 114万画素 CMOS
ISO感度 カラーISO 20000 モノクロISO 50000
ISO感度 カラーISO 5000~40000 モノクロISO 25000~200000 最短露光時間 159 nsec 最短露光時間 400 nsec
撮影速度 [fps]
× 解像度 [pixel]
20000 / 1024×1024 100000 / 640×280 2100000 / 128×8
撮影速度 [fps]
× 解像度 [pixel]
100000 / 1280×896 1000000 / 1280×32
3.3 実験条件
シャドウグラフ法による可視化計測の実験条件をTable 3-5，3-6，3-7に示す．気流条件は，
NPR = 1.5，1.7，1.9の亜音速及び遷音速領域を対象として実験をおこなった．それぞれのNPR
におけるノズル出口流速ueとマッハ数Meの関係は，NPR = 1.5のときue = 255 m/s，Me = 0.784，
NPR = 1.7のときue = 290 m/s，Me = 0.905，NPR = 1.9のときue = 316 m/s，Me = 1.00である．
各NPRにおけるノズル出口とマッハ数の関係は付録にも記載している．Table 3-5より，SA- Zとメタルハライドランプの場合，フレームレートは200 kfps，300 kfpsの2条件で，解像度 はノズル出口から壁面までのピクセル数を固定した状態で，それぞれのフレームレート時に おける最高解像度で撮影した．露光時間は，1.0 μsecとした．Table 3-6より，ACS-1とメタル ハライドランプの場合，フレームレートは200 kfps，300 kfpsの2条件で，解像度は，それぞ れのフレームレートにおける最高解像度とし，長手方向を y 軸に配置し，短手方向をノズル 出口から壁面までに配置した．露光時間は，Table 3-5 の実験条件時の流れ場全体の輝度値と 同等になるよう0.4 μsecとした．Table 3-7より，SA-Zとレーザー光源の場合，フレームレー
トは300 kfpsの1条件で，解像度はTable 3-5の実験条件時と同様とした．露光時間は，透過
率15 %のNDフィルターを使用した状態で，Table 3-5の実験条件時の流れ場全体の輝度値と
同等となるよう0.25 μsecとした．また，レーザー光源のパルス幅は10 nsecとし，立上りエッ ジによるハイスピードカメラとの同期撮影をおこなった．
Table 3-5 FASTCAM SA-Z & MID-25FCにおける可視化計測実験条件
気流条件
ハイスピードカメラの設定 フレームレート
[kfps]
解像度 [pixel]
露光時間 [μsec]
NPR = 1.5, 1.7, 1.9
200 256×232
1.0
300 256×128
Table 3-6 MEMRECAM ACS-1 M60 & MID-25FCにおける可視化計測実験条件
気流条件
ハイスピードカメラの設定 フレームレート
[kfps]
解像度 [pixel]
露光時間 [μsec]
NPR = 1.5, 1.7, 1.9
200 384×1280
0.4
300 192×1280
Table 3-7 FASTCAM SA-Z & CAVILUX Smart UHSにおける可視化計測実験条件
気流条件
ハイスピードカメラの設定 光源の設定 フレームレート
[kfps]
解像度 [pixel]
露光時間 [μsec]
パルス幅 [nsec]
NPR = 1.5, 1.7, 1.9 300 256×128 0.25 10
3.4 画像解析手法
本研究では，シャドウグラフ法を用いて撮影した画像の輝度変化を連続画像で捉えること により画像上における流体の速度ベクトルを算出する．ここでは，画像処理手法について述 べる．
3.4.1 解析領域とマスキング処理
例としてFig. 3.12には，Table 3-5及びTable 3-6の実験条件のフレームレート300kfps時に
おける生画像及び赤線で解析領域を示す．本研究では，画像解析における前処理として以下 の領域を取り除く処理（マスキング処理）をおこなった．
1）自由噴流領域において主流速である噴流内部の領域
2）対象である衝突噴流現象とは関係のない輝度変化の小さい領域
1）では，噴流内部の衝撃波が連続画像の時間進行に伴い主流方向とは逆向きに振動する様相 が可視化計測から捉えられるため，オプティカルフローでは解析が困難である．また，2）で は，誤ベクトルが算出されやすく，誤ベクトルによって流れ場全体の速度ベクトルの相対表 示が適切に行われないことがある．また，原田ら[28]によれば画像境界部において誤ベクトル が生じやすいため，画像境界部と1）2）のマスキング処理を含めた解析領域をFig. 3.13に示 す．Fig. 3.13では，青線で囲まれた領域の外をマスキングによって排除している．
(a) FASTCAM SA-Z&MID-25FC
(b) MEMRECAM ACS-1 M60&MID-25FC
Fig. 3.12 生画像及び解析領域
(a) FASTCAM SA-Z&MID-25FC
(b) MEMRECAM ACS-1 M60&MID-25FC
Fig. 3.13 マスキング領域
3.4.2 解析パラメータと解析条件
画像解析では，基本的に以下のパラメータの設定を行う．
⚫ 画像間隔時間（dt）；2枚の画像ペアの時間間隔
⚫ 画像オーバーラップ（tol）：画像ペア同士の時間間隔
⚫ ラグランジュ未定乗数（λ1，λ2）：初期推定（λ1），本推定（λ2）
⚫ 本推定イタレーション回数（ni2）：本推定における計算反復回数
本解析のフローチャートを Fig. 3.14示す．初めに，撮影された連続画像に対して先頭の画
像の時刻tをt = 0，その時の画像を初期画像I(tint)とする．ここからdtだけ時間が経過した時
の画像をI(tint+dt)とし，この 2枚の画像ペア I(tint)，I(tint+dt)から画像上における速度ベクトル
のスナップショット解を算出する．このスナップショット解を瞬間場として，連続画像に対 して本解析を適用して得られた時系列の瞬間場の結果を非定常解析結果とする．
画像オーバーラップ（tol）は画像ペア同士の時間間隔を示しており，tol ≠ 0の場合，1組目 の画像ペアに変化はないが2組目の画像ペアはI(tint+dt+ tol)とI(tint+2dt+ tol)の2枚の画像から 画像上における速度ベクトルを算出する．
Fig. 3.14 画像解析フローチャート
第 4 章 結果および考察
4.1 シャドウグラフ法による可視化結果
Fig.4.1 4.2, 4.3は，各実験条件におけるカメラのフレームレートが300 kfps時のNPR =1.5,
1.7, 1.9のシャドウグラフ法による可視化結果をそれぞれ示している．図の縦軸及び横軸は，
Fig. 3.6の衝突噴流システムの座標系に基づいて，ノズル出口高さDで正規化している．
Fig. 4.1, 4.2, 4.3より，シャドウグラフ法では，流れ場の密度変化率に対する輝度変化が自由
噴流領域や衝突噴流領域，壁面噴流領域で顕著に表れている．特に，自由噴流領域では，噴流 と周囲流体との境界面に存在する混合領域において千鳥状の輝度変化の様相が確認できるが，
その流れ構造については可視化画像から判断することは難しい．また，各実験条件において 遷音速領域にあたるNPR = 1.7 (Me = 0.905), 1.9 (Me = 1.00)では，噴流が壁面に衝突する直前に できる平面衝撃波やノズル出口形状に起因した斜め衝撃波とその反射衝撃波が噴流中で振動 する様相が確認できる[5]．特に，NPR = 1.9では衝撃波と壁面の干渉によって発生する音波が 周囲流体に伝搬する様相が顕著に表れる．これら，衝撃波の振動や音波に起因した輝度変化 は実際の流れとは異なるため，画像上の輝度変化を捉えるオプティカルフロー解析では誤ベ クトルとなる．そこで，本研究では，NPR = 1.5におけるシャドウグラフ可視化画像に対して オプティカルフロー解析を適用する．
また，Fig. 4.1とFig. 4.3 は，Photron社製のハイスピードカメラに対してメタルハライドラ
ンプとレーザー光源をそれぞれ用いて計測した結果である．Fig. 4.3は，Fig. 4.1と比較してよ り鮮明に流れ場の輝度変化が計測できている．これは，メタルハライドランプはフリッカー レスの連続光であり，瞬間場として捉えられた画像上の輝度情報は，カメラの露光時間（シャ ッタースピード）1.0 μsecに依存する．対して，レーザー光源はカメラのシャッタスピードと 同期したパルス幅10 nsecのパルス光であり，画像上の輝度情報はメタルハライドランプの場
合よりも1 / 100短い時間であるパルス幅に依存する．したがって，レーザー光源のほうがよ
り瞬間場における輝度の変化を鮮明に捉えることに成功している．これら，光源の違うシャ ドウグラフ法における可視化画像に対するオプティカルフロー法を用いた画像解析結果への 影響も後述で比較していく．
次節では，NPR = 1.5のシャドウグラフ可視化結果を中心に，オプティカルフローにおける
解析パラメータが解析結果に及ぼす影響について検討する．
(a) NPR = 1.5
(b) NPR = 1.7
(c) NPR = 1.9
Fig. 4.1 シャドウグラフ可視化結果（SA-Z & MID-25FC）
(a) NPR = 1.5
(b) NPR = 1.7
(c) NPR = 1.9
Fig. 4.2 シャドウグラフ可視化結果（ACS-1& MID-25FC）
(a) NPR = 1.5
(b) NPR = 1.7
(c) NPR = 1.9
Fig. 4.3 シャドウグラフ可視化結果（SA-Z & CAVILUX）
4.2 解析パラメータの影響
本節では，画像解析結果に影響を及ぼす解析パラメータとして，可視化画像の時空間解像 度，スナップショット解に用いる 2 枚の画像ペアの時間間隔 dt，最適解を得るための制約式 の各項の重み付けを決定するラグランジュ未定定数λの違いについてそれぞれ調査した．
4.2.1 ラグランジュ未定定数の影響
本項では，解析パラメータのうち，ラグランジュ未定定数が解析結果に及ぼす影響につい て調査する．ここで，本研究で用いるLiuらによって提案されたアルゴリズムでは，初期推定
としてHorn-Schunck推定量を算出し，そのHorn-Schunck解を初期近似としてLiuらのアルゴ
リズムによって本推定をおこなう．Horn-Schunck解は初期近似として使用されるため，Horn-
Schunck推定量のラグランジュ定数（λ1）の選択はここでは特にオプティカルフローの結果に
影響を与えないことがわかっている[18]．本研究では，本推定のラグランジュ未定定数（λ2）が 解析結果に及ぼす影響を調査するため初期推定とあわせて，(λ1，λ2) = (2, 200)，(200, 20000)の 2条件で比較した．メタルハライドランプとSA-Zのフレームレート300 kfpsで撮影したNPR
= 1.5のシャドウグラフ可視化画像に対して，ラグランジュ未定定数(λ1，λ2) = (2, 200)と(λ1，λ2)
= (200, 20000)をそれぞれ適用した画像解析結果をFig. 4.4，4.5に示す．ここで，t = 0はどち
らも同じシャドウグラフ可視化画像を基準として，(a) ~ (c)は時系列可視化画像に対する速度 ベクトル，(d) ~ (f)は流線を，それぞれシャドウグラフ可視化画像と重ね合わせた結果を示す．
Fig. 4.4では，自由噴流剪断層の揺れや剪断層の不安定性に起因したK-H渦が時系列の解析
結果において衝突噴流領域及び壁面噴流領域をへて移流していく様相が確認できる．しかし，
Fig. 4.5では，流れ場全体が平滑化され剪断層の揺れやK-H渦を瞬間場で捉えることができて
いない．ラグランジュ未定定数が大きい場合，変動は平滑化され，非定常な流れ場には適さな い．本研究で対象とする，高速で非定常な衝突噴流現象においてラグランジュ未定定数は，
(λ1，λ2) = (2, 200)が適切であると判断し，以後の解析結果においては同様とする．
(a) t = 0 μsec (a) t = 0 μsec
(b) t = 70 μsec (b) t = 70 μsec
(c) t =140 μsec (c) t =140 μsec
(d) t = 0 μsec (d) t = 0 μsec
(e) t = 70 μsec (e) t = 70 μsec
(f) t =140 μsec (f) t =140 μsec
Fig. 4.4 ，(λ1，λ2) = (2, 200)，NPR = 1.5
（SA-Z & MID-25FC）
Fig. 4.5 ，(λ1，λ2) = (200, 20000)，NPR = 1.5
（SA-Z & MID-25FC）
4.2.2 時空間解像度の影響
本項では，カメラのフレームレートと解像度に依存する時空間解像度が解析結果に及ぼす 影響を調査した．Tabel 4-1はSA-ZとACS-1の各フレームレートでの撮影時における1 pixel あたりの大きさを示している．ここで，1 pixel あたりの大きさはノズルから壁面までの距離 に割り当てられたピクセル数から算出している．ACS-1 の画角設定では長手方向の解像度で
ある1280 pixel はフレームレートに関係なく一定であり，短手方向がフレームレートに依存す
る．本実験では，剪断層の不安定性に起因した渦の移流の様相を捉えるために長手方向をy方 向に割り当てている．そのために，各フレームレートにおける1 pixelあたりの大きさは各フ レームレートにおける解像度との関係より，200 kfpsではACS-１のほうが優位であるが，300 kfpsではSA-Zのほうが優位な結果となった．そこで，より時間解像度の高い300 kfpsの場合 において空間解像度の違いによる解析結果について比較する．
Fig. 4.6，4.7は，それぞれメタルハライドランプと SA-Z及びACS-1の組合わせにおけるフ
レームレート300 kfps，NPR = 1.5での速度ベクトルと流線の時系列解析結果をシャドウグラ
フ可視化結果と重ね合わせたものを示す．ここで，t = 0はそれぞれ任意に選択したものであ る． Fig. 4.6 (a)，(d)より，(x / D, y / D) = (2.0, 1.0)で自由噴流剪断層の不安定性に起因したK- H渦が確認できる．K-H渦の大きさはx方向基準で1D程度で，Fig. 4.6 (b)，(c)，(e)，(f)にお いてK-H渦が衝突噴流領域で曲げられながら壁面噴流領域へと移流する様相が時系列画像か ら瞬間場として捉えられている．また，渦の回転速度は，渦が壁面に向かうにつれて渦の大き さが発達するとともに増加しているのが確認できる．Fig. 4.7 (a)，(d)より，(x / D, y / D) = (1.6, 2.6)及び(x / D, y / D) = (-1.5, 2.7)でほぼ左右対称に現れるK-H渦や，(x / D, y / D) = (5.6, 3.4)の 壁面近傍で境界層流れと壁面近傍の遅い流れに起因した二次渦が観測された．しかし，Fig. 4.7
(b)，(e)では，K-H渦とみられる大きさの渦の移流を捉えることができていない．また，Fig. 4.7
(c)，(f)では，K-H渦の左右対称性が崩れ放出されている様相が確認できる．これは，NPR = 1.5
が亜音速領域と遷音速領域の境界にあたり，噴流内部の衝撃波の振動に起因して噴流の左右 対称性が崩れることによる噴流剪断層の揺れによるものである．しかし，本研究で対称とす る高速衝突噴流現象の非定常解析においてACS-1とメタルハライドランプの組合わせは時系 列画像に対してK-H渦を瞬間場で捉えることができていないため，本研究では，SA-Zにおけ るフレームレート300 kfps及び解像度256×256の実験条件で撮影された時系列のシャドウグ ラフ可視化結果に対してオプティカルフロー法を適用する．
Table 4-1 各フレームレートにおける1ピクセルあたりの大きさ [mm]
SA-Z ACS-1
フレームレート [kfps]
200 0.136 0.117
300 0.136 0.189
(a) t = 0 μsec
(b) t = 70 μsec
(c) t =140 μsec
(d) t = 0 μsec
(e) t = 70 μsec
(f) t =140 μsec
Fig. 4.6 NPR = 1.5（SA-Z & MID-25FC）
(a) t = 0 μsec
(b) t = 70 μsec
(c) t = 140 μsec
(d) t = 0 μsec
(e) t = 70 μsec
(f) t = 140 μsec
Fig. 4.7 NPR = 1.5（ACS-1 & MID-25FC）
4.2.3 画像間隔時間の影響
本項では，スナップショット解を算出するための2枚の画像間隔時間dtが解析結果に及ぼ す影響を調査した．オプティカルフロー解析では，現象に対して十分な時間解像度で時系列 連続画像を取得し，適切な画像間隔時間のもとで，速度ベクトルを推定することが望ましい．
4.2,2 項の Fig. 4.6 で示した SA-Z とメタルハライドの組合わせにおけるフレームレート 300
kfps，NPR = 1.5の速度ベクトルと流線の結果は，最小画像間隔時間であるdt = 1（1/300000 sec
に相当）における結果である．Fig. 4.8は，Fig. 4.6に用いた同一のシャドウグラフ可視化結果 に対して画像間隔時間dt = 2の場合の速度ベクトルと流線の時系列解析結果を示す．Fig. 4.8 より，自由噴流剪断層の不安定性に起因したK-H渦が移流する様相が確認できるが， Fig.4.8
(d)より，その渦構造を示す流線が不安定であり，Fig. 4.6 (d)と比較して明確に捉えられていな
い．また，自由噴流領域の混合領域における速度ベクトルが主流方向の流れと一致しておら ず流れの連続性が保証されていない．また，Fig. 4.8 (e)，(f)よりノズル出口直後の剪断層にお いて噴流内部からの不連続な流れが確認できる．したがって，本研究で対象とする高速衝突 噴流現象におけるオプティカルフロー法を用いた画像解析では，フレームレート 300 kfps で 撮影された時系列のシャドウグラフ可視化画像に対して画像間隔時間 dt = 1 の場合に非定常 解析が成立すると考える．
よって，以後のオプティカルフロー解析において画像間隔時間はdt = 1を適用する．
(a) t = 0 μsec
(b) t = 70 μsec
(c) t =140 μsec
(d) t = 0 μsec
(e) t = 70 μsec
(f) t =140 μsec
Fig. 4.8 dt = 2, NPR = 1.5（SA-Z & MID-25FC）
4.3 オプティカルフローの適用結果
本節では，4.2節の解析パラメータが画像解析結果に及ぼす影響を踏まえ，NPR = 1.5の時 系列のシャドウグラフ可視化画像にオプティカルフロー法による非定常解析を適用した結果 にてついて高速衝突噴流現象における流れ場と比較していく．
4.3.1 オプティカルフローの定量的妥当性の検証
Fig. 4.9は，原田によって実施された(x / D, y / D) = (1.9, 0.8)でのNPR = 1.5における壁圧計
測から得られたパワースペクトル密度（Power Spectral Density : PSD）[28]と同位置におけるシ ャドウグラフ画像の輝度変動の PSD を比較したものである．Fig. 4.7 の横軸は，周波数
（Frequency）[Hz]，左の縦軸は壁圧の PSD [kPa2/Hz]，右の縦軸は輝度のカウント値の PSD
[count2/Hz]を示す．Fig. 4.9より，壁圧計測と可視化画像に基づく渦放出の周波数のピークは，
それぞれ7162 Hzと6800 Hzであり，両者には良好な相関関係があり，NPR = 1.5におけるシ ャドウグラフ時系列可視化画像に対する画像解析の非定常性に関する定量的妥当性を示して いる．
1.8 2.0 2.2 2.4 2.6 2.8
0 20 40 60 80 100
0 2000 4000 6000 8000 10000
Pressure Brightness
PSD(Pressure) [kPa2 /Hz] PSD(count) [count2 /Hz]
Frequency [Hz]
6800 [Hz]
7162 [Hz]
Fig. 4.9 パワースペクトル密度（NPR = 1.5）[29]
4.3.2 非定常解析結果
Fig. 4.10，4.11は，SA-Zとメタルハライドランプを用いてフレームレート300 kfpsで撮影
したシャドウグラフ時系列可視化画像に対する速度ベクトルと流線の解析結果を示す．Fig.
4.10，4.11の(a) ~ (g)より，自由噴流剪断層の不安定性に起因したK-H渦の移流を時系列に対
する瞬間場から確認できる．これらは，非定常解析における流れ場の連続性が成立している といえる．しかし，Fig. 4.10，4.11の(c) ~(f)で見られるx / D = 2.5 ~3.0付近の解析領域下部の 境界における時計回りの渦は，主流流れと噴流内部の平面衝撃波がノズルに向かって振動す る際の輝度変化を捉えた誤ベクトルであると考えられ，マスキング処理による境界付近の流 れ場の妥当性については実際の現象と比較しながら十分に議論する必要がある．
Fig. 4.12，4.13は，SA-Zとレーザー光源を用いてフレームレート300 kfpsで撮影したシャ
ドウグラフ時系列可視化画像に対する速度ベクトルと流線の解析結果を示す．Fig.4.12，4.13 より，(a)の自由噴流領域における大きさ 0.5D 程度の K-H 渦が時系列とともに渦の大きさが 発達しながら移流し，(d)では1D 程度まで発達したK-H渦が衝突噴流領域で壁面方向へと曲 げられ移流していく様相が確認できる．また，Fig. 4.12，4.13の(d)，(e)では，(x / D, y / D) =
(4.0, 0.5)の壁面近傍で時計回りの渦が確認できる．これは，衝突噴流領域内部のよどみ域にお ける再循環流れ[6]を捉えたものである．これらの結果から，レーザー光源を用いてシャドウグ ラフ可視化画像中における輝度変化が鮮明になったことによって，K-H 渦を瞬間場の時系列 として捉えることができ，更に，可視化画像からは判断の難しい再循環流れの抽出が可能と なった．
(a) t = 0 μsec
(b) t = 23.33 μsec
(c) t =46.66 μsec
(d) t = 70 μsec
(e) t = 93.33 μsec
(f) t = 116.66 μsec
(g) t = 140 μsec
Fig. 4.10 NPR = 1.5（SA-Z & MID-25FC）
(a) t = 0 μsec
(b) t = 23.33 μsec
(c) t =46.66 μsec
(d) t = 70 μsec
(e) t = 93.33 μsec
(f) t = 116.66 μsec
(g) t = 140 μsec
Fig. 4.11 NPR = 1.5（SA-Z & MID-25FC）
(a) t = 0 μsec
(b) t = 23.33 μsec
(c) t =70 μsec
(d) t = 130 μsec
(e) t = 156.66 μsec
Fig. 4.12 NPR = 1.5（SA-Z & CAVILUX）
(a) t = 0 μsec
(b) t = 23.33 μsec
(c) t =70 μsec
(d) t = 130 μsec
(e) t = 156.66 μsec
Fig. 4.13 NPR = 1.5（SA-Z & CAVILUX）
第 5 章 結論
本研究では，高速で非定常な流れ場に対してトレーサーを用いない光学的可視化手法にオ プティカルフロー法を適用する画像解析手法の非定常解析への適用可能性を調査するため，
亜音速から遷音速における衝突噴流現象に対してシャドウグラフ法を用いた可視化計測をお こない，得られた時系列可視化画像に対してオプティカルフローによる非定常解析をおこな った．
本研究で得られた結果を以下に示す．
シャドウグラフ可視化計測
⚫ NPR = 1.5，1.7，1.9の亜音速から遷音速領域における衝突噴流現象に対してシャドウグ
ラフ法による可視化計測を行った．可視化計測より，自由噴流剪断層の不安定性に起因し た噴流の揺れや，NPR = 1.7，1.9では，超音速流れで見られる壁面前方での平面衝撃波や ノズル形状に起因した斜め衝撃波とその反射衝撃波が振動する様相が捉えられた．また，
NPR = 1.9では，衝撃波と壁面の干渉によって発生する音波が流れ場全体に伝播する様相
が観測された．
⚫ SA-ZとACS-1の2台のハイスピードカメラを用いてそれぞれシャドウグラフ可視化計
測をおこなった結果，SA-Z ではノズル中心から片側半分の領域においてオプティカルフ ローを用いた非定常解析の成立性を左右する高時空間解像度での撮影に成功した．また，
ACS-1では，高解像度による衝突噴流全体の計測により，噴流剪断層の不安定性に起因し
た K-H 渦の対称性が噴流内部の衝撃波の振動などによって崩れ左右非対称に放出される 様相が観測された．
⚫ SA-Z を用いてフリッカーレスの連続光であるメタルハライドランプと短パルスレーザ ー光源それぞれでシャドウグラフ可視化計測をおこなった結果，短パルスによる高繰り返 し発光が可能なレーザー光源は，カメラの露光時間に依存するメタルハライドランプの場
合よりも1/100短い時間で現象を捉えられることから，連続光よりも鮮明なシャドウグラ
フ可視化画像を撮影することに成功した．
