第 27 回数値流体力学シンポジウム	 
A05-3

剥離制御を目的とした二次元バンプ周りの層流剥離・乱流遷移現象の基礎解析
Laminar separation and turbulent transition phenomenon around a 2D bump geometry
for separation control○ 焼野 藍子, ISAS/JAXA, 神奈川県相模原市中央区由野台 3-1-1, E-mail: yakeno@flab.isas.jaxa.jp
河合 宗司, ISAS/JAXA, 神奈川県相模原市中央区由野台 3-1-1, E-mail: kawai@flab.isas.jaxa.jp
野々村 拓, ISAS/JAXA, 神奈川県相模原市中央区由野台 3-1-1, E-mail: kawai@flab.isas.jaxa.jp
藤井 孝藏, ISAS/JAXA, 神奈川県相模原市中央区由野台 3-1-1, E-mail: fujii@flab.isas.jaxa.jp
Aiko Yakeno, ISAS/JAXA, 3-1-1 Yoshinodai, Chuo-ku, Sagamihara-shi, Kanagawa, 252-5210
Soshi Kawai, ISAS/JAXA, 3-1-1 Yoshinodai, Chuo-ku, Sagamihara-shi, Kanagawa, 252-5210
Taku Nonomura, ISAS/JAXA, 3-1-1 Yoshinodai, Chuo-ku, Sagamihara-shi, Kanagawa, 252-5210
Kozo Fujii, ISAS/JAXA, 3-1-1 Yoshinodai, Chuo-ku, Sagamihara-shi, Kanagawa, 252-5210
We investigated laminar separation and turbulent transition phenomenon around two-dimensional bump geometry forthe purpose of building a simple way to determine the optimum control-condition of separation at the wing leadingedge. As a result, we found that there was a difference in the size and occurrence position of the two-dimensionalvortex. In Reh = 4000, two-dimensional vortex exists in the bump wake rather than the separating position. On theother hand, in Reh = 16000, 64000, two-dimensional vortex is generated in the vicinity of the separating position, notin the wake. Velocity shear rate of the base flow is considered to be one of the parameters of inflection point
instability. The size and location of two-dimensional vortex are different in Reynolds numbers, because the baseflows are different.
１．背景	 	  航空機翼や回転機器翼での剥離は機器性能を大幅に低減させるため問題である．そのため，剥離を抑制する制御手法が検討されて来た．従来，吹き出し吸い込みなどのマイクロデバイスを用いて前縁部で周期的な振動を付与することにより，剥離を抑制できることが実証されていた(1-11 など)．しかしながら，翼弦長と主流速度による F+のパラメータで整理した制御駆動の最適周期は文献によって異なり，レイノルズ数や翼形状が変化した場合に，アクチュエータの制御駆動条件を決めるのが難しい．そこで本研
究では，単純形状で層流剥離・乱流遷移現象の詳細を調べ，翼前縁部で流れ場を支配する物理パラメータを探索し，F+に代わる制御パラメータに関する考察を行う．それにより，剥離制御の最適な制御駆動条件を決定するための簡便な手法を構築することを目的とする．

(4)(5)	  熱流束 q に関しても同様で，フーリエの法則から式(6)となる．
ここで κ，Prl (=0.72)，a はそれぞれ熱伝導率，プラントル数，音速である．

(6)	  内部エネルギーe は理想気体を仮定し式(6)で計算される．比熱比は γ (=1.4)である．

２．計算手法	 	  計算ソルバーとして ISAS/JAXA で開発された圧縮性流体解析
ソルバーLANS3D(12)(13)を用いる．支配方程式は，密度 ρ，速度ui，圧力 p，内部エネルギーe，マッハ数 M を用いて，以下に示す圧縮性ナビエ・ストークス方程式である．

(7)
各変数は，外部流の速度 Uinf，バンプ高さ h，音速 a，外部流の密度 ρinf よりそれぞれ以下のように無次元化してある．

(1)
(2)
(8)
式(2)，式(3)中のレイノルズ数は以下のように定義される．
(3)	  応力テンソル τ は歪み速度テンソルの線形関数として，温度 Tに依存し，ストークスの仮定より右辺第二項の係数 λ は-2/3µl となるので，式(4)となる．動的粘性 µl (T) は空気に関するサザーラン
ドの法則から式(5)として算出する．T1 = 198.6 R とする．

(9)	  本計算は高解像度の非定常計算のため，移流項，粘性項，メトリック，ヤコビアンの空間差分は 6 時精度コンパクト差分スキーム(14)を用いている．また，移流項の数値振動を抑える 10 次精度
1Copyright © 2013 by JSFM

第 27 回数値流体力学シンポジウム	 
A05-3定とした．流出部は自由流出条件として，流出部から内側の 2 点目の値を流出部の値とした(式(12))．また，流出部の逆流を防ぐた
め，流出部の流れ方向速度は絶対値に設定している．解析対称の流れ場には影響が無い．の三重対角フィルター(15)を用い，フィルター係数は 0.42 とした．
時間積分には三段階の TDV Runge Kutta 法を用いている．クーラン数が最大で1.0 以下となるよう時間刻み幅はΔt = 0.0005 とした．
マッハ数は 0.2 としている．	  計算領域の模式図を図 1 に示す．単純形状として二次元の半円柱を考え，壁面との結合部を滑らかにした式(10)の数式で表す形
状(4)を用いる．
z1rr
J0 (A)I0 (A ) I0 (A)J0 (A )=h6.04844
aa (10)
ここで，J0 と I0 はそれぞれ 0 次元のベッセル関数と変形ベッセル関数を表す．パラメータ A は 3.1926 である．	  境界条件は，片側を壁面として滑り無し条件を課した．壁面は1 点目と 2 点目に 2 次精度`陽的差分を適用しており，数値振動を抑えるため圧力は 2 点目の値を与える(式(11))．
wall =

out =

out 1

uout = |uout 1 |vout = vout 1
wout = wout 1

eout = Pout 1 /(1) + 0.5( uk uk )out 1(12)	  二次元形状の三次元計算のため，スパン方向は周期境界条件を課した．

ρinf uinf Pinf

z*/h
ρinf

wall+1

uwall = 0.0

20h ~Constant

Pinfvwall = 0.0
wwall = 0.0

Free

Constant

ewall = Pwall+1 /(in1) + 0.5( uk uk )wall(11)	  流入部は密度，ブラジウス境界層の速度分布，圧力を一定として与えた．ブラジウス境界層は 99％境界層厚さδin をバンプ高さ
h の 1/4 とした．壁面からの遠方境界も，密度，速度，圧力で一h0.0

2.0

4.0

6.0

8.0

10.0

12.0

14.0

16.0

x*/h

Fig.1 Schematic of the computational domain

z*/h
4.0
3.0
2.5
2.01.5
1.0
0.5
0.0

1.0

2.0

3.0

4.0

5.0

6.0

7.0

8.0
x*/h

9.0

10.0

11.0

12.0

13.0

14.0

15.0

0.0
16.0

Fig. 2 Computation grid (Grid 2) and instantaneous flow field; every 5 meshes are plotted in each direction and the contour represents streamwise velocity ufrom -1.0(blue) to 1.0 (red) at Reh = 4000. Flow is coming from left to right.
Table. 1 Grid resolution of computation cases
Grid resolution for δno.

GridRehReδdt

Grid number
Nx, Nz, Ny(span)

114000183.60.0005

224000183.60.0005

3116000

593.6

0.00054216000

593.6

0.000553640001280

0.0005

Grid resolution for viscous scale

Δxmax / δΔzmax / δΔymax / δ

Δx +

Δz +

Δy +

527, 157, 127

0.159

0.01550.209

0.523~8.91

0.0512~0.872

0.690~11.7

757, 215, 167

0.09900.01370.150

0.366~5.83

0.0507~0.806

0.551~8.76

527, 157, 127

0.387

0.03780.510

3.00~26.5

0.294 ~2.60

3.96~34.9

757, 215, 167

0.09900.01370.149

0.733~11.7

0.101~1.61

1.10~17.5

1051, 249, 251

0.09980.01620.132

1.13~28.8

0.183~4.66

1.49~37.9

2Copyright © 2013 by JSFM

第 27 回数値流体力学シンポジウム	 
A05-3(a)(b)Fig. 4 Inatantaneous flow field. Iso-surface of Q and contour of streamwise velocity u from -1.0(blue) to 1.0(red). (a) Reh = 4000, (b) Reh = 16000	  計算格子は，解像度を変えた 3 種類を検討する．図 2 に中間解像度の格子(Grid 2)と流れ場との比較を示す．図中，グリッド線は 5 点おきに間引いてある．等値線は瞬時の流れ方向速度 u を表す．解析
対称の領域は格子を詰めて配置し，壁面からの遠方境界は，解析対称の領域から十分離れている．	  バンプ高さ h，外部流の速度 uinf を基本とするレイノルズ数 Reh(=ρinfuinfh/µl) を 500，4000，16000，64000 とする．その時，バンプの頂点付近で，接線方向速度の最大となる壁からの距離を，最小境界層厚さ δ とすると，
図 2 に示すように Reh の増加に伴い δ は減少し，一定値に近づいていく．	  表 1 に計算を行った数種類のケースについて，その時のレイノルズ数と使用した 3 種類の格子の格子数と，δ と粘性スケールで無次元化した格子解像度を示す．また，Reh に加えて δ を基本に再計算したReδ も示しており，
Reh = 500, 16000, 64000 はそれぞれReδ = 183.6,593.6, 1280 に相当する．統計量の算出に用いた格子は Reh = 500,4000 は Grid 1，Reh = 16000 は Grid 2，Reh = 64000 は Grid 3 である．

1in64

外部流速度 uinf，壁面せん断応力 τw を用いてそれぞれ以下の式で表せる．

Cf == 0.25
2D Computation
3D by Grid 22642310410Reh105(1)(2)	  バンプ周りの Cf ，Cp，形状係数 H をそれぞれ図 5 に示す．Cf より剥離点位置は Reh = 4000 では x = 5.90，Reh = 16000 では x = 5.94 となる．
Cp より圧力の極小となる位置は Reh = 4000, 16000 いずれも x =5.75 となる．平板の層流境界層で圧力の極小となる位置では形状係数 H = 2.6，剥離位置では H = 3.5 と見積もられる．本計算の結果で，剥離点位置は H = 3.5 となる位置とほぼ一致している．一方，圧力の極小となる位置は，実際には H = 2.6 となる位置よりも 0.05 上流にずれている．これは，流れがバンプに流入する手前でぶつかり，圧力が増加する影響を受けるためと考えられる．

0.1

0.01
210w

1/2 u2inf
P PinfCp =1/2 u2inf

610Fig.3 Minimum boundary layer thickness δ around a bump
３．結果	 
３−１	  瞬時流れ場の可視化	  レイノルズ数Reh = 4000, 16000の計算結果の可視化を図4に示す．
等値面は速度勾配テンソルの第二不変量 Q として，カラーコンターは流れ方向の瞬時速度 u である．可視化に用いた Q のしきい値は，Reh = 4000 では 0.1，Reh = 16000 では 1.0 であり，Q の値を上げることで，より細かいスケールの渦運動を可視化することになる．Reh =4000 では，比較的大きいスパン方向の二次元渦が周期的に発生している様子が示されているのに対して，Reh = 16000 では明確ではない．
３−２	  形状係数	  圧力係数 Cp，摩擦係数 Cf は，静圧 Ps，外部流の静圧 Pinf，密度ρ，3３−２	  剥離現象	  剥離位置付近での時間平均した流れ方向速度分布 U(x, z)，流れ方向速度の高さ方向剪断 ∂U(x, z)/∂z，高さ方向二階微分値 ∂2U(x, z)/∂z2を，図 6, 7, 8 に Reh = 4000, 16000, 64000 についてそれぞれ示す．剥離点はいずれも x = 5.9 付近である．図 7 より，流れ方向速度の高さ方向剪断 ∂U(x, z)/∂z は，Reh = 4000 では弱い最大値をとるが，Reh =16000, 64000 ではより大きな値となり，最大値がほぼ同じであることが分かる．また，流れ方向速度の高さ方向の二階微分値 ∂2U(x,z)/∂z2 がゼロとなる高さは変曲点であり，変曲点型不安定性が起こることが古くから知られている．図 8 より，いずれの Reh でも変曲点がおよそ x = 5.7 より壁面付近に見いだされ，後流になるに従い，変曲点は Reh が高い場合から順に壁面から遠ざかる様子が分かる．	  変曲点型不安定性が起こる場合は，まず始めにスパン方向に一様	 な二次元渦が成長し，基本流の流速で移流することが一般的に知られている．図 9 は，剥離位置付近での瞬時の圧力分布と二階微分値∂2U(x, z)/∂z2 を重ねて(変曲点位置を点線で示している)，Reh = 4000,16000, 64000 のそれぞれについて示したものである．Reh = 4000 の場
合は，二次元渦によって生じる圧力変動が見いだされない．
一方で，Reh = 16000, 64000 ではほぼ同じ大きさの二次元渦が同じ周期で発生している．
３−３	  乱流遷移現象	  バンプ後流の乱流遷移位置付近の時間平均した流れ方向速度分布
U(x, z)と，乱流エネルギーの分布を図 10，11 にそれぞれ示す．Reh =4000 より Reh = 16000 の場合でより早く壁面に付着し，剥離渦は小さい．
Copyright © 2013 by JSFM

0.051.20
Reh = 4000
Reh = 16000

0.04
0.03

z0.01

1.00

Cf

0.02

0.00
-0.01
5.6

5.7

5.8

5.9x6.0

6.1

0.06.2

謝辞	 	  本研究は「HPCI 戦略プログラム分野 4 次世代ものづくり課題 1『輸送機器・流体機器の流体制御による革新的高効率化・低騒音化に関する研究開発』」の枠組みとして，理化学研究所計算科学研究機構スーパーコンピュータ「京」の使用に依るものである．ここに記して謝意を表する．	 

1.20
Reh = 4000
Reh = 16000-0.2

第 27 回数値流体力学シンポジウム	 
A05-3(1) レイノルズ数により，発生する二次元渦の発生位置と大きさに違いがある．Reh = 4000 では，剥離位置付近ではなくバンプ後流において二次元渦が存在する．一方，Reh = 16000, 64000 では剥離位置付近で二次元渦が発生し，バンプ後流には発生しない．
(2) 基本流の速度剪断 ∂U/∂z は変曲点不安定性を支配するパラメータの一つであることが示唆される．
(3) レイノルズ数による二次元渦の発生位置と大きさの違いは，基本流の分布が異なるため，変曲点不安定性の違いが原因の一つと考えられる．その詳細は別途明らかにする予定である．

-0.4

Cpz-0.6

1.00

-0.8
-1.0
5.6

5.7

5.8

5.9
(H = 3.5) 5.91

3.5

6.0x5.97
Reh = 4000
Reh = 16000(H = 2.6) 5.81

0.06

1.200.04
Reh = 4000
Reh = 16000
Reh = 64000

0.020.00
5.6

5.87

5.7

1.00

2.0
1.5
5.6

0.086.2zH

3.0
2.5

6.1

5.7

5.8

5.9x6.0

6.1

5.85.9x6.0

6.1

6.2

Fig.6 U(x, z)

6.2
0.08

Fig.5 Cp, Cf and shape factor H

0.06

４．変曲点不安定性の一般的なパラメータに関する考察	 	  剥離位置，ならびに乱流遷移付近では変曲点型不安定性と見られる二次元渦が観察される．ところで，変曲点型不安定性を支配する因子については種々の議論がある．
例えば双曲線を仮定した場合に，基本流の最大剪断速度(∂U/∂z)と，最大速度差 ΔU であることが提案され整理されている(16)．しかしながらこのスケーリング則は，実際の流体の性質を表す，慣性力と粘性力の比としてのレイノルズ数の大小に対して矛盾がある．本研究の流れ場は，剥離位置では最大速度差ΔUは一定であり，基本流の最大剪断速度 ∂U/∂z が変化する．
∂U/∂z がほぼ同じ Reh = 16000 と 64000 の場合は二次元渦の大きさと放出周波数が同じになることから，∂U/∂z は変曲点不安定性を支配するパラメータの一つであることが示唆され，超音速境界層の乱流遷移に関する坂上(2010)(17)らの主張と関連すると思われる．一方乱流遷移付近では，Reh = 4000 では剥離位置より大きなスケールの二次元渦が観察される一方で，Reh = 16000 と 64000 では見られない．
この時，∂U/∂z は Reh = 4000 でより大きい．	 
５．結論	  剥離制御の最適な制御駆動条件を決定するための簡便な手法の構築を目的として，翼前縁部を模擬した二次元バンプ周りでの層流剥離・乱流遷移現象の詳細を調べた．その結果，以下のことが明らかになった．

0.04
Reh = 4000
Reh = 16000
Reh = 64000

0.02
0.00
5.6

5.7

5.8

5.9

x6.0

6.1

6.2

Fig.7 ∂U(x, z)/∂z

0.080.06
0.040.02

Reh = 4000
Reh = 16000
Reh = 64000

0.00
5.65.75.8

5.9x

6.0

6.1

6.2

Fig.8 ∂2U(x, z)/∂z2
4Copyright © 2013 by JSFM

第 27 回数値流体力学シンポジウム	 
A05-3
Fluid Mech, Vol. 16 (1984), pp. 365-424.
(4) M. R. Visbal, D. P. Rizzetta and J. Mathew, “Large-Eddy Simulationof Flow Past a 3-D Bump,” AIAA Paper, Vol. 917 (2007), pp. 1-18.
(5) M. Kiya, M. Shimizu and O. Mochizuki, “Sinusoidal forcing of aturbulent separation bubble,” Journal of Fluid Mechanics, Vol. 342(1997), pp. 119-139.
(6) A. Roshko, “On the development of turbulent wakes from vortex
streets,” Technicsl report of National Advisory Committee forAeronautics (1954), pp. 1-25.
(7) L. W. Sigurdson, “The structure and control of a turbulent reattaching
flow,” Journal of Fluid Mechanics, Vol. 298 (1995), pp. 139-165.
(8) B. S. Stratford, “The prediction of separation of the turbulent boundarylayer,” Journal of Fluid Mechanics, Vol. 5 (1959), pp. 1-16.
(9) Nonomura, T., Aono, H., Sato, M., Yakeno, A., Okada, K., Abe, Y.
and Fujii, K.,“Control mechanism of plasma actuator for separated
flow around NACA0015 at Reynoldsnumber63,000 -Separation
bubble related mechanisms -,” AIAAPaper, (2013), 2013-0853.
(10) Sato, M., Okada, K., Nonomura, T., Aono, H., Yakeno, A., Asada, K,Abe, Y. and Fujii, K.,“Massive parametric study by LES on
separated-flowcontrol around airfoil using DBD plasma actuator
atReynolds number 63,000” AIAAPaper,2013-2750 (2013).
(11) 青野，岡田，野々村，佐藤，焼野，藤井，”中レイノルズ数での DBD プラズマアクチュエータを用いた翼剥離制御におけるバースト発振周波数効果”, 日本機械学会流体工学部門論文集	 (2013), G0405.
(12) K. Fujii and S. Obayashi, “High-resolution Upwind Scheme for
Vortical-flow Simulations”, Journal of Aircraft, Vol. 26 (1989), pp.
1123-1129.
(13) H. Aono, T. Nonomura, N. Iizuka, T. Ohsako, T. Inari, Y. Hashimoto,R. Takagi, and K. Fujii, “Scalar Tuning of a Fluid Solver using
Compact Scheme on K Computer”, Proceedings of Parallel
Computational Fluid Dynamics Conference, Hunan China, (2013).
(14) S. K. Lele, “Compact Finite Difference Scheme with Spectral-likeResolution”, Journal of Computational Physics, Vol. 103, pp. 16-22.
(15) Gaitonde, D. V. and Visbal, R. M., “Pade Type Higher-Order
Boundary Filters for the Navier-Stokes Equations,” AIAA Journal,Vol. 38, No. 11 (2000), pp. 2103-2112.
(16) 松原雅春，”変曲点型速度分布を持つ剪断流の線形不安定性に対するスケーリング則について”，計算数理工学論文集，Vol. 7,
No. 1 (1997).
(17) 坂上昇史，“超音速境界層の乱流遷移に関する研究”，大阪府立大学博士論文 (2010)

Fig.9 ∂2U(x, z)/∂z2 and pressure contour around a bump.
a) Reh = 4000, (b) Reh = 16000 and (c)Reh = 64000.	 
参考文献	 
(1) D. Greenblatt and I. J. Wygnanski, "The control of flow separation byperiodic excitation," Progress in Aerospace Sciences, Vol. 36 (2000),
pp. 487-545.
(2) M. Gad-el-Hak and D. M. Bushnell, “Separation Control: Review,”
Journal of Fluid Engineering , Vol. 113 (1991), pp. 5-30.
(3) C.-M. Ho and P. Huerre, “Perturbed free shear layers,” Ann. Rev.

5

Copyright © 2013 by JSFM

第 27 回数値流体力学シンポジウム	 
A05-3

Fig.10

Fig.11U(x, z)Kinematic energy

6

Copyright © 2013 by JSFM