第 23 回数値流体力学シンポジウム
C9-1

URANS モデルによる三角形粗度を側壁に有する水路流れの数値解析
URANS computations on channel flows with triangular roughness on side walls○ 河内 友一, 中国電・エネルギア総研, 〒739-0046 東広島市鏡山 3-9-1 ,E-mail: 260113@pnet.energia.co.jp
細田 尚，京大・工・都市社会工学専攻, 〒615-8540 京都市西京区京都大学桂 C1-3,
E-mail: hosoda.takashi.4w@kyoto-u.ac.jp
木村 一郎，北大・工・環境ﾌｨｰﾙﾄﾞ工学部門，〒060-8628 北海道札幌市北区北 13 条西 8 丁目
E-mail:i-kimu2@eng.hokudai.ac.jp
木村 敏章，京大・工・都市社会工学専攻, 〒615-8540 京都市西京区京都大学桂 C1-3
Turbulent flow structures in open channels with various roughnesses are important topics in hydropower engineering.
In this study, applicability of numerical mode with URANS type turbulence models for simulating the channel flows
with triangular roughness on side walls..We carried out the experiments and showed the comparisons withexperimental data and calculated results were generally adequate. The numerical result showed that non-linear k-εmodels could reproduce the velocity profile ,formation of cross - sectional secondary currents in the channel and themaximum-velocity-dip phenomena. Furthermore, we visualized vortex pipe from the calculated result andinvestigated the three-dimensional unsteady flow structures around side walls.

そこで，本研究では，側壁に三角形粗度を連続的に配置した水路を対象として URANS による数値解析(1)を実施し，実験結果と比較した．乱流モデルは木村らによる乱れ非負条件の非線形 k-εモデル(2)と標準型線形 k-εモデルの両方を実施し，その結果を比較した．

１．はじめに温室効果ガス削減への関心の高まりや原油価格の乱高下を背景としたエネルギーセキュリティの観点から，再生可能な純国産エネルギーである水力発電所はその重要性を大きく増している．
しかし，国内において，経済性に有利な水力発電所の開発地点は数少なくなっており，その開発にあたっては大幅な建設費のコストダウンが求められている．
水力発電所の建設費の大きな部分を占める導水路トンネル工
事では，通水能力を向上させるためにコンクリート覆工を行う．
このコンクリート覆工を省略すると大幅なコストダウンが期待できるが，この場合には Fig.1 に示すように水路壁面の凹凸が大きく不均一であることから，通常水路設計で経験的に用いられているマニングの粗度係数の値を評価することが難しい．
このような壁面に凹凸がある水路の抵抗則を正確に評価するための手法の一つとして，実用性にすぐれた URANS による数値解析手法があげられる．現地における抵抗則の予測のためには様々な粗度形状の高レイノルズ流れについて数値計算を行う必要があるが，このような場合において URANS は DNS や LES と比較して計算負荷が小さいという利点がある．

２．対象とした流れ
今回実験および解析の対象とした流れは Fig.2 に示すような両側の側壁に凹凸を模擬した人工粗度を連続的に配置した開水路流れである．掘削したトンネルの壁面は Fig.1 に示すように角錐形状が多いことから，側壁の凹凸形状には二等辺三角形を採用した．
なお，水路底面には粗度を設置していない．二等辺三角形の形状
は Fig.2 に示すように底辺 4cm，高さ 2cm とした．その他水理諸
元は Table.1 のとおりである．
３．水理実験および数値解析手法の概要(1) 水理実験
水理実験は長さ 40m，幅 60cm の二次元水路で実施した．この水路区間の 10m 間に三角形粗度を設置した．水路および三角形粗度は鋼製である．

Flow

4cm2cmSide wall
Fig .1 Headrace tunnel without concrete lining

Fig .2 Sketch of channel and triangular roughness on side walls
1Copyright © 2010 by JSFM

第 23 回数値流体力学シンポジウム
C9-1
Table.2 Calculation conditions
移流項差分
2 次精度 TVD-MUSCLE
時間積分
2 次精度 Adams-Bashforh
0.0005sec
計算時間間隔Table.1 Hydraulic conditions
諸元値0.220m
水深 h0.6m
水路幅（谷部） B
0.08m3/s
流量 Q0.61m/s
断面平均流速（谷部）Vave1/500
水路勾配 S0.41
フルード数 Fr
レイノルズ数 Re
1.3×105

計算手順は以下のとおりである．まず静止水の状態から計算を開始し，下流端を自由流出状態として上流端の境界条件を 0.1 秒で所定の流量まで増加させた．所定の流量に達した後，数秒間計算を行い，水路全体の断面流量の誤差が充分小さくなった時点で，周期境界条件による計算を開始した．

水路の縦断方向の水位をポイントゲージで計測し，三角形粗度の谷部の断面におけるX方向流速分布を微小プロペラ流速計で計測した．流速はサンプリング周期 20Hz，サンプリング数 2000 個で計測し，平均流速を求めた．

４．解析結果と実験結果の比較
実験，解析（非線形 k-εモデル，線形 k-εモデル）の谷部断面の時間平均 X 方向流速コンターを Fig.4 に示す．解析結果には YZ方向流速ベクトルをあわせて示している．なお，解析結果は周期境界条件の解析開始から 30～60sec の 30sec 間の値をサンプリング周期 20Hz で平均した流速値である．

(2) 数値解析
数値解析手法は木村らによる 3 次元移動一般曲線座標系で流速ベクトルの反変成分で表示した基礎式による手法(1)を用いた．乱流モデルには木村らによる乱れ強さ非負条件を考慮した非線形 kεモデル(2)を用いた．なお，比較のために標準型の線形 k-εモデルによる解析も実施した．実際の流れは開水路であり，木村らの手法は自由水面を有する流れも再現できるが，ひとまず本研究では，簡単のために水面をリジッド・リッド条件とした．水面近傍の鉛直方向の乱れの減衰を考慮するため，減衰関数を水面付近の渦動粘性係数に乗じている．
水面の実験を実施した水路は 10m と長区間であることから，上下流端には周期境界条件を用い，水路の長さは三角形粗度 30 個分（1.2m）とした．水路は半断面を計算対象とし，水路中央はフリースリップ条件を設定した．なお壁面条件には壁法則を用いており，滑面として計算した．
計算格子はＸ方向 300×Y 方向 30×Z 方向 40 分割とし，粗度付
近は Fig.3 のように分割している．その他の解析条件についてはTable.2 のとおりである．

Experimental result

Water surface

Side wall

Center-line

Caluculated result (Non-Linear k-εmodel)(a) Cross section of calculation gridCalculated result (Standard k-εmodel)
Fig .4 Comparison of the calculated velocities with experimental data(b) Plan view of calculation grid around triangular roughness
Fig .3 Calculation grid2

Copyright © 2010 by JSFM

第 23 回数値流体力学シンポジウム
C9-1
中山ら(5)により，水面変動が発生するような流れでは，水面が穏やかな流れと比較して鉛直方向乱れ強度が増大することが指摘されている．今回の解析ではこれらの自由水面が大きく変動している場合の現象を考慮されていないことが，2 次流の再現性に課題が残る要因の一つと考えられる．また，三角形粗度谷部にははく離が発生しているが，このような箇所に壁法則を使用していることの妥当性についても今後の検討が必要である．
このように流速分布の再現精度に課題は残るものの，側壁に三角形粗度を設けた開水路流れにおける底面セルおよび側壁セルの形成，2 次流による velocity dip 現象などが非線形 k-εモデルにより概ね再現できることが示された．

実験では最大流速部分が水面より下に生じるいわゆる Velocitydip が現れている．非線形 k-εモデルによる解析結果においても同様に水面よりやや下に最大流速部が生じている．非線形 k-εモデルのベクトル図では粗度側の水面近傍および底部に明確な 2 次流が形成されている．非線形 k-εモデルによる解析により，実験の流速分布の傾向が概ね再現できていることがわかる．ただし，最大流速部が実験結果と比較してやや下側であり，さらに高速部分が側岸側に寄っている．一方，線形 k-εモデルにおいて 2 次流はほとんど発生しておらず，最大流速も水面付近に発生しており実験の再現性が悪い．Fig.5 に谷部の横断方向の流速分布の実験と解析結果の比較を示す．非線形 k-εの結果はとくに水面付近の再現性にすぐれているが，Z=0.5H 付近では高速部が側壁側により過ぎていることがわかる．
非線形 k-εモデルによる解析結果の 2 次流セル図と 2 次流の流速の大きさを最大流速で除した値のコンターを Fig.6 に示す．
Tominaga ら(3)や今本ら(4)の長方形断面開水路における実験で報告されている底面セルおよび側壁セルが再現されている．ただし，2次流の最大値は，主方向流速最大値の 10%程度と比較的大きい値である．水路断面のアスペクト比や粗度の形状が異なるので単純な比較はできないが，Tominaga ら(3)による側壁に相当粗度 1.2cmの粗度をつけた開水路流れの実験では，2 次流の大きさは主方向流速最大値の 2%程度である．非線形 k-εによる解析結果における高速部が，下部および側壁側によっているのは，これらの 2 次流が実際よりもやや大きく評価されていることによると推測される．
実験では三角形粗度付近の自由水面が大きく変動していたが，今回の解析では水面にリジッド・リッド条件を用いている．
Z=0.89h1.2

EXP.

1.0EXP.

1.0

NL

0.8

NL0.8LL

0.60.6U(m/s)

U(m/s)

1.20.4

Z=0.27h

0.4
0.2

0.2
0.0

0.0-0.2-0.2-0.02

0.03

0.08

0.13

0.18

0.23

-0.02

0.28

0.03

0.08

0.13

EXP.

0.8NLLU(m/s)

U(m/s)

1.0

0.4
0.2

1.2

EXP.

1.0

NL

0.8

L

0.18

0.23

0.28

0.18

0.23

0.28

0.18

0.23

0.28

0.6
0.4
0.0

0.03

0.08

0.13

0.18

0.23

0.28

-0.02

0.03

0.08

0.13

Y(m)

Y(m)

Z=0.64h

1.2

EXP.

NL

1.0NL

L

0.8

L

1.2

EXP.

1.0
0.8
0.6

U(m/s)

U(m/s)

0.28

-0.2

-0.02

0.4

Z=0.09h

0.6
0.4
0.2

0.2

0.0

0.0

-0.2

-0.2
-0.02

0.03

0.08

0.13

0.18

0.23

-0.02

0.28

0.03

0.08

0.13

Z=0.46h

1.2EXP.

NL

1.0NLL0.8LEXP.

1.00.8U(m/s)

1.2FlowY(m)Y(m)

U(m/s)

0.23

Z=0.14h0.2

0.0-0.20.60.4

Z=0.05h

0.6
0.4
0.2

0.2

0.00.0-0.2

-0.2
-0.02

0.18

Y(m)

Y(m)
Z=0.82h

1.20.6

５．流れの３次元構造
次に，非線形 k-εモデルによる解析結果を用いて，三角形粗度を有する開水路流れの 3 次元構造について検討を加える．
ここでは，三角形粗度付近における複雑な渦構造を可視化する方法として，中山ら(6)と同様に速度勾配テンソルの第 2 不変量による方法を用いた．
解析開始から 35.05sec と 36.15sec における第 2 不変量の等値面の鳥瞰図を Fig.7 に示す．特徴的なのは，水面壁側に発生している渦管構造である．
三角形粗度の壁面から水路の中央方向に伸び，渦の回転は流れ方向に対して反時計回りである．渦管は間欠的に発生しており，水面壁側の低速部を形成していることがわかる．
また，底部にも同様の渦管が発生しているが，こちらの管は流れ方向に対して時計回りである．Z=0.5H 付近には，水面・底部ほど管状ではなく水路中央への伸びは小さいが，やや弱い渦が塊状に発生している．

0.03

0.08

0.13

0.18

0.23

0.28

-0.02

0.03

0.08

Y(m)

0.13(a) t=35.05sec

Y(m)EXP:experimental results NL:Non-linear k-εmodel L:linear k-εmodel

Fig .5 Comparison velocity profiles

Flow(a) t=36.15sec
Fig .6 Secondary current vector and cells(calculated result of non-linear k-εmodel)

Fig .7 Vortex tubes (calculated result of non-linear k-εmodel)
3Copyright © 2010 by JSFM

第 23 回数値流体力学シンポジウム
C9-1
今後は自由水面を考慮に入れた URANS モデルによる解析の実施および，LES による解析結果と URANS モデルとの比較を行い，複雑な粗度形状の水路における流れ構造の解明とその水路を対象とした場合の URANS モデルの再現精度の向上を図りたい．

つぎに，非線形 k-εモデルによる解析結果における水面近傍の粗度山部周辺 5 点の主流速の時系列を Fig.8 に示す．とくに水面付近・粗度付近では流速の変動が大きく，間欠的に先ほどの 3 次元的な渦の発生による低速部が現れていることがわかる．この 5点における約 40 秒間における流速変動量の卓越周波数は 3～7secであった．

参考文献(1) 木村，細田, 音田, "完全スタガード移動一般座標系における三次元開水路乱流シミュレータの開発," 四日市大学環境情報論集, Vol.5 , No1･2 , (2002), pp. 145-170.
(2) 木村，細田, "乱れ強さ非負条件を考慮した非線形 k-εモデ
ルによる立方体周辺の流れの三次元解析,"水工学論文集,
Vol.44 , (2000), pp. 599-604.
(3) Tominaga, A., Nezu, I., Ezaki,K.,Nakagawa,H., “Three-dimensionalrbulent structure in straight open channel flows,” Journal ofhydraulic research, Vol.27,No.1, (1989), pp. 149-173.
(4) 今本，石垣，塩野，”開水路側壁近傍流れの 3 次元構造について”,京都大学防災研究所年報，第 38 号，B-2，(1995)，pp.529-538.
(5) 中山，禰津，”開水路流れにおけるバーストの水面近傍での挙動及び乱れ構造との関連性について”,土木学会論文集，No.635,/Ⅱ-49,(1999),pp.31-40.
(6) 中山，焦，”速度勾配テンソル不変量解析を応用した取水槽
流れの 3 次元渦構造特性の検索”,水工学論文集，第 45 巻，(2001)，pp.451-456.

６．おわりに今回，側壁に三角形粗度を有する開水路流れを対象としてURANS モデルにより解析を実施した．得られた結果は以下のとおりである．
(1) 非線形 k-εモデルにより，側壁に三角形粗度を設けた開水路流れにおける底面セルおよび側壁セルの形成，2 次流によるVelocity dip 現象などが概ね再現できることがわかった．ただし，水面付近の 2 次流が過大に評価されるなど流速の再現精度には課題が残る．
(2) 非線形 k-εモデルによる結果から，三角形粗度を有する開水路流れの 3 次元の渦構造が推定できた．水面粗度山部周辺の流速変動は 3～7sec 周期で発生していた．
point.A

0.8Distance from edge of roughness =0.028m Z=0.214m Tp=7.4sec

U(m/s)

0.6
0.4
0.2
030

35

40

45

point.B

0.8

50 time(sec) 55

60

65

70

Distance from edge of roughness=0.061m Z=0.198m Tp=7.4sec

U(m/s)

0.6
0.4
0.2
030

3540

50 time(sec) 55

45
point.C

0.860

65

70

Distance from edge of roughness =0.028m Z=0.198m Tp=5.3sec

U(m/s)

0.6
0.4
0.2
030

3540

45
point.D

0.8

50 time(sec) 55

60

65

70

Distance from edge of roughness =0.014m Z=0.198m Tp=5.3sec

U(m/s)

0.60.40.2
030

354050 time(sec) 55

45
point.E

0.860

6570

Distance from edge of roughness =0.028m Z=0.178m Tp=3.1sec

U(m/s)

0.6
0.4
0.2
030354050 time(sec) 55

45606570(a) Time series of mean velocity around roughness(calculated result of non-linear k-εmodel)AB

DCE(b) Location of pointFig .8 Time series of mean velocity around roughness

4Copyright © 2010 by JSFM