平面的流れ場・拡散場の高精度計算法の開発

九州大学学術情報リポジトリ

Kyushu University Institutional Repository

平面的流れ場・拡散場の高精度計算法の開発

矢野, 真一郎

Graduate School of Engineering, Kyushu University

https://doi.org/10.11501/3142491

出版情報：Kyushu University, 1998, 博士（工学）, 課程博士 バージョン：

権利関係：

平面的流れ場・拡散場の高精度計算法 の開発

平成 10 年 3月

矢 野 真 一 郎

圃色』

「平面的流れ場・拡散場の高精度計算法の開発 J

目 次

第 1 章 序 論

1.1本研究の背景 1.2既往の研究

1.2.1平面

2

次元乱流モデルについて 1.2.2移流分散係数について

1.3本研究の目的と本論文の構成 参考文献

2今ノ釘

AU

寸

ζ J

門/

第 2 章 平 面 的 流 れ 場 の 計 算 法 の 高 精 度 化

2.1はじめに 10 

2.2平 面

2

次元 k‑

ε

乱流モデルの誘導 11  2.2.1標 準 型

3

次元 k‑e乱流モデルについて 11  222Ui，k，Eの鉛直分布形 11  2.2.3渦動粘性係数vtの鉛直分布形 13  2.2.4平面

2

次 元

k‑

ε一方程式の修正 ，r  15  2.2.5補正係数の推定式の導出 15  2.2.6モデル定数 C

kとC

Eの推定式の導出 16 

2.3補正係数とモデル定数のチューニング 17  2.3.1補正の効果について 17  2.3.2底面摩擦による productionに対するモデル定数 C_k，C_Eについて 18  2.3.3修正された

ι ε

モデルを用いた平面的流れ場の数値シミュレー

ションについて 18 

2.4おわりに お

参考文献 24 

第 3章 平 面 的 潮 流 場 に お け る 粒 子 追 跡 計 算 法 の 高 精 度 化

3.1はじめに

3.2粒子追跡計算の高精度化

3.2.1粒子追跡計算 (Euler‑Lagrange法)の基礎式 3.2.2ランダム変動速度の評価

3.3数値実験について

3.3.1潮流シミュレーション 3.3.2拡散能と海水交換率の関係

3.3.3湾内水の滞留に対する拡散能の影響について

41  41  41  42  44  44  45  46 

圃色』

3.4おわりに 参考文献

第

4 章

閉 鎖 性 海 域 に お け る 平 面

2

次 元 拡 散 シ ミ ュ レ ー シ ョ ン の た め の 分 散 係 数 の 推 定 法

4.1はじめに

4.2平面

2

次元モデルのための分散係数の評価法 4.2.1分散係数のモデル化

4.2.2αの決定のための潮流・水質拡散シミュレーション 4.2.3計算結果について

4.3比 例 定 数 αの普遍関数表示化

4.3.1成層度パラメータに関する考察 4.3.2吹送流パラメータに関する考察 4.3.3 'l'と φによる比例定数 αの関数表示 4.4分散係数の異方性についての検討

4.4.1非定常性の導入について 4.4.2異方性の導入について 4.4.3異方性の修正について 4.5おわりに

参考文献

48  49 

8 9 9 9 1 2 2 3 4 4 4 5 6 7 8  

ζ J ζ J ζ J ζ J f O / o f O / O / O / O f o f o f o f o f o  

第 5

章 新 し い 分 散 係 数 の 評 価 法 を 採 用 し た 物 質 循 環 モ デ ル に よ る 閉 鎖 性 内 湾 の 水 質 予 測

5.1はじめに 82 

5.2博多湾における水質予測シミュレーションのための物質循環モデル 83  5.2.1分散係数の評価法について 83  5.2.2博多湾における既往の水質予測シミュレーションについて 83  5.2.3物質循環モデルについて 84  5.3博多湾における水質予測シミュレーションの計算結果 86  5.3.1比例定数 α，

s

の評価 86 

5.3.2計算結果について 87 

5.4おわりに 89 

参考文献 90 

第

6

章 結 論 99 

謝 辞 101 

第 1 章 序

1 . 1本研究の背景

三とト

員同

今世紀の急激な科学技術・産業の発展と人口の増加に伴う自然環境の悪化により，地 球規模から我々の身近な生活空間に関わる小さな規模までの様々な環境問題が，近年世 界的・全人類的な共通課題としてクローズアップされてきている.そのような状況の中 で，水工学 (hydraulicengineering)の分野においても自然水域の水環境保全を研究の主目 的とした環境水理学 (environmentalhydraulics)が中心的な存在になりつつある.

また，近年のコンピューターの演算性能の急速な発展と大容量化，並びに操作性の向 上と低価格化は 水 工 学 の 分 野 に お け る 学 術 的 研 究 の 方 法 論 を 大 き く 変 え た . す な わ ち，かつて理論的解析と実験的解析が主流であったものが， 1960年代から始まった流れ と物質輸送の解析をモデル化した支配方程式を基に数値シミュレーションを行う数値水 理学 (computationalhydraulics)と， 1990年代から始まった情報工学的要素も含めたより 総合的な数値水理学の発展形である hydroinformatics [Abbott (1991， 1994)]がその座を 取って代わろうとしている.

これらの水工学の最近の動向は実務レベルにおいても同様であり 1980年代以降必ず しも水工学・環境水理学の専門家ではない技術者が end‑userになり，パッケージ化もし くは製品化された数値解析コードを利用して水環境問題の解析を行うとわう Abbott(199  3)の分類するところの第 4世代モデリング (fourth‑generationmodelling)  (もしくは，人 工知能 (AI: artificial inteligence)や expertsystemを組み込んだ第 5世代モデリング (fifth‑ generation modelling) )が頻繁に行われている.

このような時流の中で水工学が実学である工学として求められていることのーっとし て，流動・拡散場解析のためのより汎用性・柔軟性の高い数値モデル(現象を表す物理 モデルや水質・生態系モデル 解析結果の可視化ツール等も含めて)の開発があると言 えよう.

水 工 学 が 対 象 と す る 水 域 は 実 験 室 レ ベ ル の 非 常 に 小 さ な も の か ら 河 川 ・ 湖 沼 ‑ estuary .海洋などの自然水域のように非常に大きなものまで様々である.それらの水域

には大小様々な時空間スケールの流れが混在しており このことが流れや物質拡散など の現象のモデリングを行う際に困難をもたらす原因と考えられている.そこで，ある時 空間スケールで表される(周)波数以下のスケールの現象を解析で得られる解により表 現し，それ以上の現象は陽的に表現せずに別の形式で(解析で得られる物理量で)モデ ル化することでこの困難を緩和するフィルタリング (filtering)という方法論(例えば，

Reynolds平均など)が考案され，広く適用されてきた.

水工学の中でも特に環境水理学が対象とする自然水域は 水平方向の空間スケールに 比べて鉛直方向のスケールが非常に小さい平面的な水域が一般的である.そのような平 面的な水域においては 流れ場・拡散場の特性量の平面的な分布が第一の関心の対象と

なってきた. 1990年代以降の計算機環境の整備に伴い ようやく 3次元的な現象につい

圃t

ても議論できるようになってきたが，社会基盤整備のための開発事業に対し事前に行わ れる環境アセスメントや水環境の改善策の検討等の際には，今後も平面的分布特性が第 一の関心の対象となり続けるものと思われる.

フ ィ ル タ リ ン グ の 一 種 と 考 え ら れ る 水 深 平 均 操 作 か ら 得 ら れ る 平 面

2

次 元 モ デ ル (depth‑averaged model)は，上述のような昨今の状況からも，またその利便性・経済性か らみても，今後も水工学・環境水理学・計算水理学・ hydroinformaticsの中において重要 な位置を占め続けるものと思われる.

1 . 2 既往の研究

水工学で取り扱う現象は一般的に，流体(通常は，水)の流れを支配する方程式と物 質の輸送方程式により記述される.それらの基礎方程式に対するフィルタリングには大 きく分けて表‑1.1に示したような 4つのものが考えられる [Vieira(1993)参照].フィ ルタリングにより発生する，基礎方程式中に含まれる各物理量の平均量からの偏差分に ついての相関(モーメント)を，平均操作された方程式から解として求められる平均量 により，モデル化する必要がある.

通常，最もよく用いられる偏差量同士の相関のモデリングとして，平均値の勾配で線 形近似する勾配拡散近似 (Boussinesq近似)があり，そのモデル係数である拡散能とし て表‑1.1 に示したような運動量・物質・熱などの各種物理量の拡散係数がある.拡散 係数や相関自体のモデリングのうちスケール

l

を除く残りの

3

つについては，表中に示 しているように様々なモデリングが試みられているものの普遍的なモデルの確立までに は未だ至っていない.

本論文で対象とする平面2次元モデルは，スケール 3のフィルタリングに相当してい るが，大きく分けて流動シミュレーションのための平面

2

次元乱流モデルと拡散シミュ レーションのための移流分散係数の

2

つ の モ デ ル に つ い て 高 精 度 化 が 必 要 と さ れ て い る.ここでは， これらのモデルに関する既存の研究を紹介する.

1.2.1平 面

2

次元乱流モデルについて

平面2次元乱流モデルとは， Reynolds応力，もしくは SGS(Sub‑Grid Scale)応力の水 深 平 均 値 を 勾 配 型 近 似 し た 場 合 の 見 掛 け の 渦 動 粘 性 係 数 を モ デ リ ン グ し た も の と 言 え る.一般的に適用されているものを以下に整理して示す.

a).0方程式モデか

混合距離理論を基に構築された

O

方程式乱流モデルは，渦動粘性係数¥を乱れの代表 長さスケール(混合距離)と代表流速スケールの積により表す.通常 平面

2

次元モデ ルにおいては，代表長さとして水深 hを代表流速として底面摩擦速度 ll. を採用し，次 式のように表されることが多い.

Vt=

α

^U* 

h 

︐ ^︑︐^{︑ ' ︐ ︐}

..  

z . 

/ 唱.EA •

' EK 

岡t

比例定数 αとしては，開水路等流の場合から α=0 (10‑¹)が標準値としてよく利用さ れている [Vieira(1993)Jが，乱流構造が底面摩擦のみに支配されていない複雑な流れ 場では

α

のチューニングが必要であり一般性に乏しい.

b).2方程式モデル

数ある乱流モデルの中でも最もよく利用されているものとして

k ‑ ε

乱流モデルがあげ られる.これは，その利便性とこれまでの実績が評価されてのことと考えられる.

平面2次元 k‑ε乱流モデルの流動シミュレーションへの適用は， Rastogi  and Rodi  (1978)  に よ り 流 入 の あ る 開 水 路 流 に 対 し て 初 め て 試 み ら れ て 成 功 を 収 め て 以 来， McGuirk and Rodi (1978)による横流入開水路流の解析， Chapman and Kuo (1985)による 急拡のある関水路流の解析， ASCE Task Committee (1989)にレヴューされている Rodiら の研究グループによる一連の解析， Younus and Chaudhry (1994)による跳水のある流れ場 等の解析に利用され一定の成功を収めている.

しかしながら，標準型として現在利用されている Rastogiand Rodiのモデルについて は適用限界を指摘するものも幾つかある.例えば， Chu and Babarutsi (1989) と Booij (1989)は，水平方向と鉛直方向の乱れの production項の修正の必要性を指摘している.

木村ら(1997)は標準型の平面2次元 k‑ε乱流モデルでは， mixing layerに発達する非定 常な大規模水平渦を再現出来ないことを指摘している.また， Rastogi and Rodiのモデル は彼ら自身が指摘しているように，厳密には乱れエネルギーとその散逸率の水深平均値 を評価できるモデルにはなっていないなど問題点が多い.

c). SGS (sub‑grid scale)モデル

LES (Large Eddy Simulation)で用いられる SGS(Sub‑Grid Scale)渦粘性係数の考え方を 平面

2

次元モデルに適用した研究が Madsenet al. (1988)と中辻ら(1992)により試みられ ている.それらによると， Grid Scale L1でフィルタリングすることにより得られた SGS 応力%を，次式で表される Smagorinskyモデルで評価している.

内 n ̲ ̲ 

f  d U i   ̲  dU

i ¥ 

‑τij= LV~ υ = Vol .::.....ーニ

+̲̲ 

^J  I 

.J  " . J   " ¥  

d x  

j 

d x  

i }  (1.2a) 

Ve 

= 

(1 

C

_s 

L 1   f  ( 

²

^むら

^{)112 (1.2b)} 

ここで， i，j=l，2，

c s

はSmagorinsky定数であり，フィルター幅 Aが慣性小領域にある 場合の理論値として C キ0.2が 得 ら れ て い る [ 数 値 流 体 力 学 編 集 委 員 会 (1995)J .  Madsen et al.と中辻らは，それぞれ潮流場におけるC5の値を実測結果を基に評価してい るが，現段階においてそれらの係数の一般性についての確認はされていない.また，局 所平衡性が保証されていない波数領域におけるスケール分離には物理的根拠が無いとい

う批判もある[灘岡・八木(1993)J .  d).その他のモデル

その他の平面

2

次元乱流モデルとしては，乱れの構造を水平方向と鉛直方向に分離し て，それぞれに起因する渦動粘性係数を独立にモデル化した Chuand Babaruttsi (1989)の

2‑component modelがある.また，乱流現象を水深以下のスケールを持つ

3

次元的な乱流 (SDS (Sub‑Depth scale)乱流)と平面2次元的な大規模渦とにスケール分離し， SDS乱 流成分を

1

方程式モデルとしてモデル化した灘岡・八木(1993)による SDS&2DHモデル がある. 2‑component modelは水平・鉛直方向の乱れの相互作用が考慮されていないこ と， SDS&2DHモデルは SDS乱流の長さスケールの評価が暖昧であることや水深以下 のスケールを持つ計算格子(LU， L1y )を用いる場合に対してモデリングの物理的整合性 がないことなどの問題点があると考えられる.

1.2.2移流分散係数について

移流拡散方程式を水深平均 もしくは断面平均することにより移流項から生じる移流 速度と濃度の偏差分の相関項は一般的に見かけ上の拡散効果を生じさせる.この効果を 移流分散 (dispersion)もしくはシアー拡散 (sheardiffusion)と呼び¥一般的に乱流拡散に 比べて非常に強い拡散効果を持つことが知られている.相関を勾配拡散近似する場合の 拡散能を移流分散係数 (dispersioncoefficient)と呼ぶ.移流分散の概念の起源は， Taylor  (1953， 1954)による有界な領域 (boundedregion)である円管路内の層流・乱流場における 理論と Corrsin(1953)による有界でない領域 (unboundedregion)における理論にあると言 われている[国司(1976)J . 

水深平均により得られる平面2次元モデルにおける移流分散係数のモデル化は， Elder  (1959)により与えられた対数分布則が成り立つ開水路乱流場における理論に始まったと 言え，現在に至るまで平面的流れ場における最も標準的な評価法として広く利用されて いる.また， Bowden (1965)は幾つかの流速分布形を持つ場における移流分散係数と単 一 周 期 を 持 つ 往 復 流 に お け る 一 周 期 平 均 の 移 流 分 散 係 数 を 評 価 し て い る . しかしなが ら，これらの評価式は理論的であるものの単純な場を仮定しているので 現実の自然水 域に適用する場合にはモデル係数の調整が必要とされる.

自然水域の中でも，特に閉鎖性内湾などの estuaryにおいては，成層や吹送流などの 影響を受けて移流分散係数の値は大きく変化する [Fischeret  a1.  (1979)が，現場実測と 室内実験の結果をまとめている

J

.速水・宇野木(1970)は瀬戸内海の年平均塩素量分布

を断面平均 1次元モデルで解析し，移流分散係数として 107cm2/sのオーダーが必要であ ることを示し，我が国における移流分散現象の研究の端緒を開いた.また，成層の影響 が顕著な感潮河川河口部において， Fischer (1972)と Dyer(1974)件断面 l次元モデルに おける移流 fluxについて流速と濃度を各種の平均量 (Reynolds平均・ー潮汐間平均・横 断方向平均・水深平均など)とその偏差分とに分離し 物質輸送機構の解明を試みた.

移流分散係数については，その後も特に自然水域において無数の研究(特に，地域的 な研究)が行われてきたが， Taylorが提出した理論的な解析以上の成果は上げられてい ない.それは，第一に個々の自然水域に固有の流動場の複雑さが障壁となっており，第 二に系統立てた研究が少なかったためと考えられる.よって 工学的な利用への便宜を 図るためにも，簡便かつより汎用的な分散係数の評価法の確立が望まれている.

1.3

本研究の目的と本論文の構成

上述のように 従 来 か ら 水 工 学 ・ 環 境 水 理 学 の 分 野 に お い て 頻 繁 に 利 用 さ れ ， こ れ か ら も な お 利 用 さ れ 続 け る も の と 期 待 さ れ て い る 平 面

2

次 元 モ デ ル に よ る 流 動 場 と 拡 散 場 の数値シミュレーションには，依然として改善の余地が多く残されている.

本 論 文 で は ， 平 面

2

次 元 モ デ ル に よ る 流 動 ・ 拡 散 シ ミ ュ レ ー シ ョ ン を よ り 信 頼 性 の 高 いものとすることを目的とし 平 面

2

次 元 乱 流 モ デ ル と 移 流 分 散 係 数 の 評 価 法 の 高 精 度 化 を 試 み る . ま た ， 実 用 計 算 と し て も 良 く 利 用 さ れ る 平 面

2

次元モデルによる粒子追跡、

計 算 に お け る 拡 散 効 果 の 採 り 入 れ 方 の モ デ ル 化 に も 取 り 組 み ， 具 体 的 な 適 用 例 と し て 拡 散が閉鎖性海域における海水交換率に与える影響について調べる.

これらにより 一 段 と 信 頼 性 の 増 し た 平 面

2

次 元 流 動 ・ 拡 散 モ デ ル が 開 発 で き る た め ， よ り 汎 用 性 の あ る 解 析 コ ー ド の 開 発 が 可 能 と な り ， 広 く 水 工 学 の 分 野 に お け る 実 務 計算への利用へ供されるようになることが期待される.

本論文は以下の様に構成されている.

第2章 で は ， 平 面2次 元 モ デ ル に よ る 流 動 シ ミ ュ レ ー シ ョ ン の 高 精 度 化 を は か る た め に，平面

2

次 元 k‑

ε

乱流モデルの改良を試みる.k‑

ε

モ デ ル に 含 ま れ て い る 各 物 理 量 に 開 水 路 乱 流 場 に お け る 半 経 験 的 な 鉛 直 分 布 形 を 仮 定 し ， 標 準 型 と な っ て い る Rastogi and Rodiの モ デ ル の 修 正 を 行 う . ま た ， 平 面

2

次 元 k‑e乱 流 モ デ ル に 特 有 の 底 面 か ら の production項 の モ デ ル 定 数 の 修 正 も 同 時 に 行 う . こ の 修 正 型 平 面

2

次元戸 ‑ε乱 流 モ デ ル を 平 面 的 流 れ 場 の 再 現 計 算 に 適 用 し ， 本 モ デ ル の チ ュ ー ニ ン グ を 行 い ， 妥 当 性 の 検 討 を 試みる.

第3章 で は ， 閉 鎖 性 内 湾 の 海 水 交 換 率 を 算 定 す る 際 に よ く 利 用 さ れ る 平 面2次 元 モ デ ルによる粒子追跡計算 (Euler‑Lagrange method)に お い て ， 拡 散 能 の 強 さ に 対 す る 海 水 交 換 率 の 応 答 性 を モ デ ル 湾 に お け る 数 値 実 験 に よ り 検 討 す る . 拡 散 能 と し て は 乱 流 拡 散 と 移 流 分 散 が 考 え ら れ る た め ， 先 ず そ れ ら に 対 し て 個 別 に ラ ン ダ ム 変 動 成 分 を 算 出 で き る よ う に 積 分 時 間 ス ケ ー ル の モ デ ル 化 を 行 う . そ し て 拡 散 係 数 の 強 さ と 異 方 性 の 大 き さ が 短 期 的 な 海 水 交 換 率 と 長 期 的 な 湾 内 水 の 滞 留 へ 与 え る 影 響 を 評 価 ・ 検 討 す る こ と を 試みる.

第

4

章 で は ， 閉 鎖 性 内 湾 に お け る 平 面

2

次 元 拡 散 シ ミ ュ レ ー シ ョ ン の た め の 拡 散 能 で あ る 分 散 係 数 に つ い て の モ デ ル 化 を 試 み る . 分 散 係 数 を 場 所 毎 の 潮 流 最 大 流 速 と 平 均 水 深 の 積 に 比 例 す る 形 式 で 表 し ， 各 湾 で 一 定 と お か れ た 比 例 係 数 を 博 多 湾 を は じ め と す る

6つ の 内 湾 に お け る 潮 流 ・ 拡 散 シ ミ ュ レ ー シ ョ ン を 行 う こ と よ り 同 定 す る . つ ぎ に ， 各 湾 に お け る 成 層 度 と 吹 送 流 が 分 散 に 与 え る 影 響 の 強 さ を 表 す

2

つ の パ ラ メ ー タ を 提 案 し ， そ れ ら に よ る 比 例 係 数 の 普 遍 関 数 表 示 を 試 み る . ま た 分 散 係 数 に 異 方 性 を 導 入 し，数値実験的に等方性の度合いの効果を検討する.

国L

ー圃圃』ーー

第5章では，第4章において提案された閉鎖性内湾における平面2次元拡散計算のた めの分散係数の評価法を適用して，博多湾における水質予測シミュレーションを行う.

博 多 湾 に お い て 同 定 さ れ た 生 化 学 的 モ デ ル パ ラ メ ー タ を 用 い ， 化 学 的 酸 素 要 求 量 COD と全リン濃度 T‑Pの循環を考慮した物質循環2層モデルを適用することにより，非保存 性物質の拡散シミュレーションに対する新しい分散係数の評価法の適用性を検討する.

第6章では，各章で、得られた知見を総括し結論とする.

圃L

一第 1 章 参 考 文 献 ‑

木 村 一 郎 ・ 細 田 尚 ・ 村 本 嘉 雄 ・ 楼 井 寿 之 (1997): "開水路横流入部における渦運動の数 値シミュレーション^H，水工学論文集，41，pp. 717‑722. 

国司秀明(1976): "海洋における物質の分散H，海洋物理学 II(寺本俊彦 編) ，東京大学 出版会，pp. 39‑61. 

数値流体力学編集委員会(編) (1995):乱流解析，東京大学出版会， 314p.

中辻啓二・狩野晋一・栗田秀明(1992): "SGS渦 動 粘 性 係 数 を 用 い た 大 阪 湾 潮 流 の 有 限 要素法解析"，水工学論文集，36， pp. 693‑696. 

灘岡和夫・八木宏(1993): "浅い水域の乱流場に関する数値計算モデルの開発と沿岸流 場への適用^H，土木学会論文集， II‑24， pp. 25・34.

速水領一郎、字野木早百(1970): "瀬戸内海における海水の交流と物質の拡散H，海岸工学 講演会論文集， 17， pp. 385‑393. 

Abbott， M. B. (1994):  "Hydroinformatics: A Copemican Revolution in Hydraulics"，よ Hydr. Res.， 32 (Extra Issue)， pp. 3‑13. 

Abbott， M. B. (1993):  "The Dynamic Environment: An Introduction"， CoastaJ， EstuariaJ and  Harbour Engineer's Reference Book (edited by M. B. Abbott and W. A. Price)， E & FN  SPON， pp. 3‑11. 

Abbott， M. B. (1991): Hydroinformatics: Information TechnoJogy and the Aq

. u

atic Environment，  Gower， Aldershot， UK. 

The ASCE Task Committee on Turbulence Models in Hydraulic Computations (1989): "Turbu‑ lence Modeling of Surface Water Flow and Transport"，よHydr.Engrg.， ASCE， 114 (9)， pp.  970‑1073. 

Booij， R. (1989):  "Depth‑averaged K‑ε'‑modelling"， Proc. of XXIIIrd Cong. IAHR， A， Ottawa，  pp.199・206.

Bowden， K. F. (1965): "Horizontal Mixing in the Sea Due to a Shearing Current"，よ FJuidMech.，  21 (2)， pp. 83‑95. 

Chapman， R. S. and Kuo， C. Y. (1985): "Application of the Two‑equation k‑εTurbulence Model  to  a Two‑dimensional， Steady， Free Surface Flow Problem with Separation"， Int. J. Num. 

Methods in Fluids， 5， pp. 257 ‑268. 

Chu， V. H. and Babarutsi， S. (1989):  "Modelling the  Turburent Mixing Layers in  a Shallow  Open‑Channel"， Proc. of XXllIrd Cong. IAHR， A， Ottawa， pp. 19ト198.

Corrsin， S. (1953): "Remarks on Turbulent Heat Transfer"， Proc. of 1st Iowa Sympo. on Thenη0‑

dynamics， Iowa State Univ. 

Dyer， K. R. (1974): "The Salt Balance in Stratified Estuaries"， Estuarine Coastal Mar. Sci.， 2， pp.  273‑281. 

Elder， J.  W. (1959): "The Dispersion of Marked Fluid in Turbulent Shear Flow"，よ FluidMech.，  5， pp. 554‑560. 

Fischer， H. B.， List， E. J.， Koh， R. C. Y.， Imberger， J. and Brooks， N. H. (1979): Mixing in Inland  and Coastal Waters， Academic Press. 

Fischer， H. B.  (1972):  "Mass transpo口mechanisms in  partially  stratified estuaries"，よ Fluid Mech.，  53 (4)， pp. 671‑687. 

Madsen， P.  A.， Rugbjerg， M. and Warren， 1.  R. (1988): "Subgrid Modelling in Depth Integrated  Flows"， Proc. of 21 st Int. Conf. Coastal Engineering， ASCE， 1， pp. 505‑511. 

McGuirk， J.  J. and Rodi， W. (1978): "A Depth‑averaged Mathematical Model for the Near Field  of Side Discharges into Open‑channel Flow"，よ FluidMech.， 86 (4)， pp. 761‑781. 

Rastogi， A. K. and Rodi， W. (1978): "Predictions of Heat and Mass Transfer in Open Channels"， 

よHydr.Div.， ASCE， 104 (HY3)， pp. 397‑420. 

Taylor， G. 1. (1953) : "Dispersion of Soluble Matter in Solvent Flowing Slowly through a Tube"，  Proc. Royal Society ofLondon， A219， pp. 186‑203. 

Taylor， G. 1. (1954) : "The Dispersion of Matter in Turbulent Flow through a Pipe"， Proc. Royal  Society of London， A223， pp. 446‑468. 

Vieira， J.  R. (1993): "Dispersive Processes in Two‑Dimensional Models"， Coastal， Estuarial and  Harbour Engineer's Reference Book (edited by M. B. Abbott and W. A. Price)， E & FN  SPON， pp. 179‑190. 

Younus， M. and Chaudhry， M. H. (1994):  "A depth‑averaged k‑εturbulence model for  the  computation of free‑surface flow"，よ Hydr.Res.， 32 (3)， pp. 415‑444. 

以上の他に本章を執筆するにあたり全体的に以下の文献を参考にした.

岩佐義朗(編著) (1995):数値水理学，丸善，214p.

玉井信行(1980):密度流の水理，技報堂出版，260p. 

Rodi， W. (1984): Turbulent Models and Their Application in  Hydraulics  A state‑of‑the‑art  review (2nd edition)， IAHR Monograph， BALKEMA， 104 p. 

表 ‑1.1流 れ 場 ・ 拡 散 場 の フ ィ ル タ リ ン グ の 分 類 [Vieira (1993)を 参 照 の 上 ， 加筆修 正 ]

ス ケ ー ル 物 理 過 程 フィルタリング 勾配型近似を用いた 分散(モーメント) ・拡散能

場合の拡散能 のモデリング

分子運動 Newtonの法則 動粘性係数

Fickの法則 分子拡散係数

2  乱 流 Reynolds平均 渦動粘性係数 0方程式モデル(混合距離理論)

、

。

(ensemble .時間平均) 渦動拡散係数 1方程式モデル

' 

2方程式モデル

(k ‑eモデル，k‑ωモデルなど)

Reynolds応力モデル，代数応力モデル 等 Grid Scaleに対する平均 渦動粘性係数 LES乱流モデル

3  水深・断面 水深平均・断面平均 移流分散係数 Taylor， Elderのモデル

流速分布 運動量補正係数 等

、. 

・・・・・・・ ・・・・・・・ ・・・ ・・ ・・ ・⁴・・

4  時空間的な モ デ ル 解 像 度 (L1x，  L1t  )  sub‑grid粘 性 係 数 SGSモデル 流速分布 に対する平均 sub‑grid分 散 係 数

第 2 章 平面的流れ場の計算法の高精度化

2 . 1 はじめに

内湾や沿岸域などにおける水質汚濁の予測や水環境の改善策の検討には，潮流場と汚 染物質の拡散場の数値シミュレーションが通常よく利用されている.潮流シミュレーシ ヨンには， しばしば水平方向の渦動粘性係数として時空間的に一定な値が仮定されてお り，経験的な値が採用されることが多い. しかしながら，潮流シミュレーションにおけ る水平方向の渦動粘性係数の推定精度はやはり重要であり，その点を指摘している研究 例も幾っかある.例えば， Oonishi (1 977)は内湾と外海の聞の海水交換に強く関連してい る潮汐残差流のパターンが，渦動粘性係数の強さに依存していることを報告している.

阿部ら(1995)も渦動粘性係数の値を系統的に変えた内湾における潮流シミュレーション を行い，渦動粘性係数と海水交換率の関連性が高いことを報告している.従って，信頼 性の高い潮流・水質拡散シミュレーションを行うためにも正確な渦動粘性係数の評価が 重要となっている.

水深が浅く水深方向に良く混合した(鉛直スケールに対して水平スケールの非常に大 きい)沿岸域における潮流・水質拡散シミュレーションに対して，平面

2

次元モデルは 3次元モデルに比べて利便性と経済性に優れている. Rastogi and Rodi (1978)は水深平均 された流れ場(開水路上流部に流入があり mixinglayerが発生する場)_{~の計算に初めて} 平面

2

次元 k‑ε乱流モデルを適用し成功を収めている.また， McGuirk and Rodi (1978)  は開水路中に横流入がある場に対して平面

2

次 元 k‑ε乱流モデルを適用し，良い結果を 得ている.その後も，平面

2

次元

ι ε

乱流モデルを適用した平面的流れ場の計算に関す る研究が数多く発表されている[例えば， Chapman and Kuo (1985)， Booij (1989， 1991)，  Chu and Babarutsi (1989)，小松ら (1993)，Younus and Chaudhry (1994)， J alil et a1. (1994)，木村 ら(1997)等].それらの研究では現在標準型となっている Rastogiand Rodiのモデルを 改良して適用している.特に， Booij (1989)は標準型の平面

2

次元 k‑εモデルにおいて底 面摩擦に起因する productionと水平方向の shearに起因する productionの比率に修正を 行って用いている.また， Chu and Babarutsi (1989)は水平方向の shearによる渦動粘性 係数を修正した平面

2

次元 k‑εモデルを用いて表し，また底面摩擦による渦動粘性係数 をO方程式モデルを用いて表して，これらを組み合わせた two‑componentモデルを提案 している.ところで，これらの修正モデルが一定の成功を収めているのは，標準型であ るRastogiand Rodiのモデルに適用限界があることを示唆している.

良く知られているように，水平方向流速 U

i，乱 れ エ ネ ル ギ ‑k，乱れエネルギーの散 逸率 ε，渦動粘性係数 V(などの各物理量が鉛直分布を持つために，

k

一方程式や ε一方 程式などの輸送方程式中の移流項に対する水深平均操作が，流下方向に付加的な拡がり の効果，いわゆるH移流分散効果"("Shear効果"とも呼ばれる)を生じさせる.また，複 数の物理量の積の形を持つ他の項についても，水深平均操作を行ってこれらの項を単に それぞれの水深平均量の積で表そうとすると，補正の必要性が生じる. しかしながら，

前述のRastogiand Rodiタイプの標準型の平面

2

次元 k‑

ε

乱流モデルにはこれらの補正は 考慮されてはいない.そのため，彼ら自身が言及しているように標準型の平面 2次元

ι

E乱流モデルでは厳密な意味での kと εの水深平均値を算出することはできない.

そこで，本章では標準型平面

2

次元 k‑

ε

乱流モデルの k一方程式と

ε

一方程式に対し て，水深平均操作に起因する誤差に対する補正係数の導入によりモデルの修正を試み，

より厳密な基礎方程式の導出を行った.具体的には，各物理量に対して開水路における fully‑developed turbulent flowでの半経験的な鉛直分布形を援用して仮定した.それらの 仮定をもとに補正係数を算出し，平面

2

次 元 k‑E乱流モデルの修正の必要性について考 察を行った.また，平面的な流れ場に対して修正された平面

2

次 元 k‑ε乱流モデルを適 用し，補正係数のチューニングと本モデルの妥当性の検討を試みた.

2 . 2 平 面 2 次元 k ‑ ε 乱 流 モ デ ル の 誘 導

2.2.1標 準 型

3

次元 k‑

ε

乱流モデルについて

標準型の

3

次元 k‑ε乱流モデルの基礎方程式は Einsteinの総和規約を適用すると次式 の様に表される.

(2.1) 

， 

p u 

p  + 

以 一

町 一

丸

q a

一 丸

弘 一

+  日 丸

以一 副 .k一方程式

(2.2) 

ぷ 一

k

ε 

司ノ

‑

c 

p  ε

一

k

Cし

C  +  ε

一町

︑

d

竹 一 て 乱

q a

一 以

ε

一 一

一向

︑

dてd

u  + 

従一 副

‑ ε

一方程式

(2.3) 

P 

= Vt (

担 i+ 担 L ¥ 担

L

¥ Oχj  OXi ) OXj  ^，  ここで，

(2.4)  Vf=CMK2. _，....

ε

^， 

，Uiはxi方向の流速 モデル定数の標準値 i.j=l，2，3，XIは空間座標 (X

1とX

2は水平方向，x

3は鉛直方向) 成分，¥は渦動粘性係数， Cμ， C

1e， C_2e， ^(J_k'  ^(J

e はモデル定数である.

[ Launder and Spalding (1974)参照]を表‑2.1 に示す.

222Ui'k，εの鉛直分布形

平面

2

次元 k‑ε乱流モデルに対して補正効果を導入するために，水平方向流速・乱れ

. . ̲  

エネルギー・乱れエネルギーの散逸率・渦動粘性係数については以下のように開水路等 流の fully‑developedturbulent flowにおける半経験的な鉛直分布形を援用して採用するこ

とにする.

a).水平方向流速の鉛直分布形について

先ず，水平方向の流速 Ui(i=l，2)には対数分布則を仮定する.対数分布則は相当砂 粒粗度 ksを代表長さにとることにより以下のような無次元形で表される.

立

L= 

lZn  { 

^̲L ¥ 

U * i  

1(  ¥ 

k s  

J 

=~ln(i)+A

^(2.5) 

ここで，

z

は底面から垂直上向きにとった座標， 1(は Karman定数(=0.4)，【人iはxi方 向の底面摩擦速度，仏 =(U$12+UJ)l

ペ

，C= 

z 

/h， hは水深， 55=K5/hである.

定数項 A については， Nikuradse (1933)により行われた円管流における一連の実験的 研究により， k+=U.ks/vに対して図‑2.1に示されているような関係が明らかにされ

ている.これを式で表すと以下のようになる.

‑完全粗面の場合(k + > 70 ) 

A  = A r =  8 . 5  

‑粗滑遷移領域の場合 (4‑‑5くk+< 70 ) 

A  = f ( 平)

‑水理学的滑面の場合 (k+<4‑‑5)

A 

= 

lZn ( 竺 ̲&)+Ac

1(  ¥  V  I  ~

. 

. 

(2.5a) 

(2.5b) 

(2.5c) 

ここで ，

v

は動粘性係数， A5=5.5で あ る . よ り 厳 密 な 分 布 を 要 求 す る な ら ば Coles (1956)の wake関数を導入する必要があるが，本研究ではモデルの簡単化のために敢えて 導入はしなかった.

b).乱れエネルギーと乱れエネルギ一散逸率の鉛直分布形について

乱れエネルギ ‑ kと乱れエネルギーの散逸率 εの鉛直分布形には 以 下 に 示 す 禰 津 (1977)による関水路乱流場における半経験式を採用した.

宅

= D

仰( ‑ 2  c )  

Uf 

回 E 仰( ‑ 3  c  ) 

U}

仔

(2.6) 

(2.7) 

ここで，D 

= 

4.78， E 

= 

9.76である.図‑2.2a)，b)に式 (2.6)と(2.7)による kと εの鉛直 分布形を示す.

式 (2.6)，(2.7)は壁面領域 (wallregion) 0く ご く 0.15を除く外部領域において妥当である ことが禰津・中川(1986)により明らかにされているが，本研究では壁面領域においても これらの分布形をそのまま拡張して用いることにした.

c).各物理量の水深平均値について

式 (2.5)，(2.6)， (2.7)を全水深

c = 

0 ‑‑1で積分することにより得られる，流速・乱れエ ネルギー・乱れエネルギーの散逸率の(無次元)水深平均値は次のように表される.

‑水平方向の流速

ω '  

1

一K

+  A 

一 一

1 7

らn 

l

一K

ι . Z

一 一

.  

^(2.8) 

， 

$ 

‑乱れエネルギー

よ‑偽=

2.07  Uf 

(2.9) 

‑乱れエネルギーの散逸率

ヰ=々= 9 . 8 4   U1 

(2.10) 

ここで， ψは流速係数である.

2.2.3渦 動 粘 性 係 数 ¥ の 鉛 直 分 布 形

良く知られているように，等流の乱流場において底面から水表面までの全領域で流速 の対数則が成り立つ場合には，渦動粘性係数¥は次式で表されるような放物線分布を持 つ

古=州 ^1‑ c )  

^(2.11 ) 

. . . ̲  

一方，

k

とεの分布形(式 (2.6)，(2.7) )を

ιE

乱流モデルによる¥の評価式 (2.4)に代 入すると，次のような

v t

の分布形が得られる.

̲̲̲̲}i̲ 

= 

Cμ

， ̲ n

² 

JC  exp (  ‑ c  ) 

U

^寧 h ^，‑ E  ^(2.12) 

図‑2.3a)に式 (2.11)と (2.12)により表される円の鉛直分布形を示す.これらの分布 形と Nezu& Rodi (1986)によりまとめられた開水路乱流場における実験結果図‑2.3b)

を比較すると，中間領域 ( intennediate region)  0.15 < c < 0.6においては式 (2.12)が， 自由水面領域 ( free‑surface region)  0.6 < c < 1.0と壁面領域 0<c<0.15では式 (2.11)  が比較的良く一致している.中間領域における両者の分布形の違いがあまり大きくない

ことから，渦動粘性係数の鉛直分布形として式 (2.11)を採用することとした.

式 (2.11)を c

= 

^{0 ‑‑ 1}の全範囲で積分することにより，渦動粘性係数の(無次元)水 深平均値は次式のようになる.

̲̲̲̲}i̲ = 

a " .   =

五二=0.0667 

U* h  ..  6  ^(2.13) 

Rastogi and Rodi (1978)は3次元 k‑ε乱流モデルとの analogyから，水深平均された渦 動粘性係数が式 (2.4)と同じ形式の関係式により評価できるものと仮定した. しかしな がら，その際のモデル定数 Cμ は

3

次元モデルと同じものを採用している.

そこで，次式で表すような補正係数

s

_v を導入することにより，

l Z

^{の評価式の修正を}

試みた.

五 = s V

_t C_J.1 

s

^{v は式 (2.9)，(2.10)， (2.13)}より計算され，次のようになる.

s V

^t 

= ^生

^t

^a

^e 

=  ^{1 . 7 1}  

Cμ

a f  

(2.14) 

(2.15) 

すなわち，修正された評価式 (2.14)とRastogiand Rodiモ デ ル に お け る 玄 の 評 価 式 の 比は約1.7ということになる.従って，鉛直構造として式 (2.11)の放物線分布を持つ渦 動粘性係数 vI の 水 深 平 均 値 五 は ， 式 ^(2.14)，(2.15)を用いることによってより精度良く 評価できるものと期待される.

2.2.4平面 2次元

k‑

，ε一方程式の修正

平面 2次元

ι ε

乱流モデルにおける克一，

E

一方程式は

3

次元の輸送方程式 (2.1)，(2.2)，  (2.3)を底面から水表面に渡って積分することにより得られる.次式に示すように補正が 必要と思われる各項にはそれぞれ補正係数を導入した.

(2.16) 

. k

一方程式

。 T T . d k  ̲  d { 

a. ̲̲̲ 

v

_t 

d k  

^¥.̲l  a 

一

_{d t}   + f J

kHA 

U

i ::'" 

=  ̲ 

^V 

^{I  f J}

^{kHD一 一 一}

^{I +  f J}

^{kHP P H}  

⁺ 

^{Ck一 一 一 ε}

I 'H~" ^•

d X i   d X i  ¥ 

I .U~~ (J"k 

dXi) 

I 

‑ 2

一方程式

一 ‑ ; ‑ ;

d " E   d 

( 予

;aE i 

i l  n 

竺 d t  

+ sI ε..，~~.. HA

U i

. 

  d ~C X i   d = 

，U 

X i  ¥ I 

^sEHD一 一εσ

dx  i ) 1   + 伽

1 " "  H PClε... ，

円 ¥

¥ 

k} 

^{P H  A A}  

+ 

CE ，..

笠 ‑ h 

2  I"~ -..~ s2εC2εζ 

k 

(2.17) 

(2.18) 

一 以 一 弘

J

一 均 一 以

+  一

川 一

一 門

例

H 

P 

ここで，

p

^凶刷'A'

s

m払HA'

ε

の

p r o d ω u c ω t i o ∞ n

に関連した経験的なモデル定数 ，i

ム

，j

= 

，

1

上

2

でで、ある.なお，標準型である

R a s t o g i  a n d  R o d i

モデルは，補正係数を全て1.

0

とし，モデル定数

C

k，

C

^t:を後で述べる式 (2.29)， (2.30)により評価したものである.

通常，移流項を水深平均することにより生じる移流分散効果は 勾 配 型 の 拡 散 項 の 形 式で表現され，その効果は移流分散係数として表されるが 本研究では補正係数の形式 で評価することとした.

2 . 2 . 5

補正係数の推定式の導出

k‑E乱流モデル中の物理量，すなわち Ui，k，ε，¥の任意の深さの点における値は，それ らの水深平均値(式 (2.8)，(2.9)， (2.10)， (2.13) )と仮定された鉛直分布形(式 (2.5)，(2.6)， 

なお，

ん( c 

)は物理量 φ(

c 

)を水深平均値 φ で無次元化した鉛直分布形である.それらを 3 次元の基礎式 (2.1)，(2.2)， (2.3)の中へ代入し，水深平均操作をすることにより各補正係 数は以下のように算出された.

φ ( c ) = φ f φ (  c 

^{)の形式で表される.}

を組み合わせることにより，

(2.7)， (2.11)  ) 

(2.19) 

skHA = 

ω φ f   " I ? / r ¥  (  ‑

1.

65  +  0 . 4 3 2   A' ) 

P

ε

^HA 

= 右(‑ ^{7.84+  1 . 0 1}  

^{A' ) 問 )}

伽 = o o z γ ^D =  ^{0 . 9 3 9}

^{仰 )}

P E m = 0 0 2 ! ? E = O

^{仰 向}

skHP 

= て し ( 0 . 4 4

ー

0 . 2 7 8A '   +  0 . 6 7 7  A  ，

2 )  (2.23) 

Qフ~

a V t  

ßêHP= 偽~

⁽ 

^{0 . 9 8 5   ‑ 0 . 4 4}  

^A' 

^{+  0 . 0 7 1 2}     1 (

A

，

2 )  (2.24)  G

ε

avtD

ψ2  

ここで， A' 

=  ‑

^(1/1() ln 

c .  

^+ A 

=

ψ+ 1/1(である.但し，

s

_2eの計算におし}てのみ積分の発 散を防ぐために， 5=0 55の領域において k= DU.² = const.と ε=ε(

え )

^{= const.を仮定}

した .

s

^{kHD と}

s

^e仰 を 除 い て ， 流 速 係 数 ψの 関 数 形 と し て 表 示 さ れ て い る . 三 は 式 (2.8) により ψと次式の様に関係づけられている.

ι ⁼ 

^exp [ 

^{1 ( (}  

^A ‑

^{ψ )   ‑ 1   ]} 

^(2.8a) 

2.2.6モ デ ル 定 数 CkとCeの推定式の導出

玉一方程式 (2.16)中の底面摩擦に起因する production項 P_kV= C_k U//hは，式 (2.3)の全 体の production項の内鉛直方向 shearによる production項のみを取り出して次式のよう

に積分することにより得られる.

れ

v= 

h 1 2   L 

^Vf { ( 

^{aa~l r  +  す ( } r dc} 

^(2.26a) 

領域

c

= 0 ‑‑

C s

^{において constantshear stress} :汽 (dUI/dX3)=U.f=const.と線形流速分 布 Ui/U.i=A5/55を仮定し，式 (2.5)，(2.11)を用いてこの積分を実行することにより，

モデル定数 Ckは次のように得られる.

C k   = ψ+ 会 =ψ

^(2.27) 

同様にして，

t

一方程式 (2.17)中の底面摩擦に起因する production項 Pu=CE UJ/h2  は，

P^EV 

= 吋州知

²

⁺  を ( } r

^dc 

に式 (2.5)，(2.6)， (2.7)， (2.11)を代入して計算することにより得られる.

これよりモデル定数 CEは次のようになる.

I 

1 

{ 1 仰(イ 1 (  ^{1  ‑ e  )} 

^de 

⁺ 

^A

^{仰(ゴ c s)} 

C

E

=

^{日 己}I

~ I 

^序

(1‑C)dC+ 

^{f 7}  

1

1 ( ん

^s

c ペ イ c s

2 . 3 補正係数とモデル定数のチューニング

2.3.1補正の効果について

， ， 

(2.26b) 

2.2.5において導出された各補正係数と流速係数の関係を図‑2.4に示す.通常，水理 学が対象とする場では流速係数 ψは 8....̲̲  25の間の値をとるといわれている[椿(1973) 参照]ため， ψ=6 (完全粗面の時，式 (2.8a)より生 =0に相当する) ....̲̲  25の間に対し て補正係数の計算を行っている.また こ こ で は 底 面 が 完 全 粗 面 に 相 当 す る 場 合 (A=  A̲ 

= 

^8.5) についてのみ図示した.

五 と 長 の 輸 送 方 程 式 中 の 移 流 分 散 効 果 は

s

_k削 と

s

f.仰の値から評価できる.図

‑ 2 . 4 a )

で，

s

_k削 は ほ ぼ lに近い値を示している.従って，平面 2次 元 モ デ ル に お け る 互 の 輸 送 に対しては移流分散効果はあまり重要でない.一方，

s

^{f.HAの値は1.0に比べて小さく，}

εの輸送に対してはその効果が大きいことが分かる.

同様に，式 (2.21)，(2.22)より拡散項への補正効果が，図‑2.4b)より水平方向のshear による production項への補正効果が，また図‑2.4c)よ り 五 の 散 逸 項 に 対 す る 補 正 効 果 が評価できる.これらの補正効果は総じて五に対しては重要であるが 正に対してはそ れ程重要で、はないことが分かった.

2 . 3 . 2

底面摩擦による

p r o d u c t i o n

に対するモデル定数 C

k，C_eについて

修正された平面

2

次元 k‑ε乱流モデルにおいて底面摩擦による k，e のproductionは，

式

( 2 . 2 7 )

，

( 2 . 2 8 )

に示したモデル定数 C_k，C_eにより摩擦速度 U.，水 深 hと関係づけられて いる.広く利用されている標準型の平面

2

次 元 k‑ε乱 流 モ デ ル で は，C_k， 

C

e に対して

Rastogi and Rodi (1978)により次の評価式が提案されている.

C k   =ψ 

C

_E 

=  3.6ψ3/2C2E c J / 2  

( 2 . 2 9 ) 

(2.30) 

Rastogi and Rodi は，C_kの値を等流場における bulkのエネルギーバランスから推定し た.また， Laufer (1951)の実験結果から渦動粘性係数の水深平均値を求め， CEの推定式 を算出した.なお， ASCE Task Cornmittee (1988)は同様の考え方から一般形として次式 を紹介している.

C ε = 

{e*(jt

ド ^/2ψ3βC 2 ε C J / 2

(2.30a) 

ここで，^(}r 

= ^玄

^{/I;は乱流 Prandtl/Schmidt数，e*} 

= 

^I; / U*h は無次元水深平均乱流拡散 係数である. e._{として，実験水路では} 0.15を，実河川では Fischeret al. (1979)より 0.6

を推奨値としている.

Rastogi and Rodiのモデルと本研究で修正されたモデルについてモデル定数の評価式を 比較すると，

C

^k についてはほぼ同じであるが Ce はかなり異なっている.完全粗面の場 合の流速係数 ψに対する両モデルの C^e の値を図‑2.5に示す .Ck に関する両者の一致 は，簡便な計算のために導入された仮定がbulkのエネルギーバランスに対しては適した ものであったことを示している. 一方， CEの違いは，定式化において用いられた実験結 果と導出方法の違いに起因するものである.

本修正モデル と Rastogiand Rodiモデルの比較により，厳密には正確とは言えないが 従 来 型 (Rastogi and Rodi型)の克一方程式の利用が可能であること， 一方 正 確 な 玉 と 長の算定のためには三一方程式の修正が不可欠であることが明らかになった.

2 . 3 . 3

修正された k‑

ε

モデルを用いた平面的流れ場の数値シミュレ‑ションについて 修正された平面2次元

ι ε

乱流モデルの有効性の確認のために， Booij(l991)により行 われた水理実験をもとに平面的流れ場における流れのシミュレーションを行った.

a).平面的流れ場の実験について

Booij (1991)は開水路の横に正方形の cavityを設置した場合について流れの実験を行 った.この流れでは開水路と cavityの接続部に生じる mixinglayer内の乱れに起因する 乱流拡散により運動量が cavity内に輸送される.その運動量に駆動される cavity内の circulationと 乱 れ エ ネ ル ギ 一 疋 の 分 布 の 測 定 を 行 っ て い る . 図‑2.6に実験に用いられ

た開水路の概形を，表‑2.2に実験条件を示す.

また， Booijに よ る 流 速 分 布 と 乱 れ エ ネ ル ギ ー だ の 分 布 の 測 定 結 果 を そ れ ぞ れ 図 ‑ 2.7a)， b)に示す.測定方法は，流速分布についてはフロートとプロペラ流速計を用いて いる.互については 電磁流速計により底面から

O

.4

h

の地点で水平方向 2成分のみを測 定 し て 疋 の 水 平 方 向 成 分 を 求 め ， 禰 津 (1977)による開水路乱流における

3

方向の乱れ 強度の比率から，疋の水平方向成分の和を1.2倍することにより算出している.

b).平面的流れ場の数値シミュレーションについて

本研究で修正された平面2次 元 k‑ε乱 流 モ デ ル の 適 用 性 を 検 証 す る た め に ， 平 面2次 元数値シミュレーションにより上述の Booijによる開水路実験の再現計算を試みた.

流れの基礎方程式は以下に示す水深平均された連続の式と運動方程式である.

(2.31a) 

2E+JL(h a t   axj‑ 可)

^= 0 

‑連続の式

‑運動方程式

a  ^{U i   h} 

^I 

a  ^{U i  U}

^j 

^{h  ̲  ̲} 

^L 

a 

^{I ̲ I  L  ¥ I} 

a  ^{h 京¥ ‑ ‑ ; ‑ ;   ， . ; 一 一

t  j = ‑ g  

h-(zb+h)+ ート~ij

^J ‑

c f   U i ゾ匂 E f i

d t   d X j  

^LJ 

d X i '  

‑ V   I 

d X j  ¥  P I 

‑J  ‑， 

(2.31 b) 

(2.31c) 

円 ︒

了

κ 2

一

3 一 同 一 丸

+  一 川

一 弘

J

一 % 一

一 円

ρ '

. Reynolds応力

ここで ，Zbは基準水平面から底面までの高さ ，gは 重 力 加 速 度 (

= 

9.8 m²/s)  ， pは水の 密度，

q j

はKroneckerのデルタ，今は底面の摩擦係数(ィマ=1 /ψ) ， i， j 

= 

1，2である.

乱流モデルとして平面

2

次元 k‑ε乱 流 モ デ ル ( 式 (2.16)，(2.17)， (2.14)) を用いて，運 動 方 程 式 中 の 渦 動 粘 性 係 数 百 と 乱 れ エ ネ ル ギ 一 疋 を 算 出 し た . 計 算 条 件 を 表‑2.3に 示す.流速係数が ψ=15.8の 場 合 の 各 補 正 係 数 と モ デ ル 定 数 の 値 を ， 式 (2.19)‑‑(2.25)，  (2.27)， (2.28)により算定し，表‑2.4に示した.ここで，これらの定数の算定には，流速 分布として仮定された対数分布則 (2.5)中の定数項 Aの値が必要となる.Booijの実験条 件から完全粗面を仮定すると，式 (2.8)より相当砂粒粗度 ksは約 2m mと推定される.

このとき k+は開水路部で約 63となり粗滑遷移領域にある. しかし 図‑2.1より推定 しでほぼ A

= 

A 

= 

8.5と見なせるものと判断し，底面は完全粗面として取り扱うことに した.また，今回の検討では底面の摩擦係数 C

fは一定値を採用しており その点からも 完全粗面が仮定されていることになる.

、

』

各基礎式の移流項に 2次精度の donnar‑cell法を，拡散項には中央差分法を適用して staggard格子系を用いて差分化し， ADI法[金子ら(1975)参照]により計算を行った.

計算領域は図‑2.6に示している通りである.

境界条件は，流入条件として図‑2.6に示す A‑A'線上において

h  =  h o

， 

U  = 

qo / 

h

， 

V  =  0 

， 

長。=偽

U*02，

ス̲，.

U

^牢03

句 ー uε

ヲ，

を与え，流出条件として図‑2.6に示す B‑B'線上において

a k   a E   a u h   a v h   o  d X   d X   d X   d X  

を与えた.ここで ，

U

， 

V

はx，y方向流速 ，

U .

_o = qo / 

h

ψ=cf l/2q

。

^/hである.

(2.32a)  (2.32b)  (2.32c) 

(2.32d) 

(2.32e) 

壁面における境界条件としては以下に示す壁関数法[Launder and Spalding (1974)参 照]を用いた.壁面上では，壁面に垂直方向の流速成分として

Unw = 

0 

(2.33a) 

を与えた.また 壁面から垂直方向に第 l番目の格子点上で 壁面と平行な方向の流速 成分と乱れエネルギー 並びにその散逸率として

万 一 ^p=E1ln i

^主旦叫、

r"  1(  ¥ 

1 0   J '  

τ^て一一2

k= と コ ι

べv /

C

'‑'

I

μ 

I  

E J ^{u * w}

³¹ 

1( x_nw 

(2.33b) 

(2.33c) 

(2.33d) 

を与えた.ここで，万二は壁面での摩擦速度，x

nw は壁面から第一格子点までの距離，l

。

、

^h

は壁面の粗さを示すパラメータであり， Booijにより 10= 0.33 X 10‑⁴ mが与えられた.

初期条件として，

0

方程式モデルを用いた予備計算より得られた流速分布を与え，収 束計算を行った.

c).数値シミュレーションの結果について

本研究で修正された平面 2次元 k‑ε乱流モデルによる Booijの水理実験の再現計算結 果を図‑2.8に示す.図‑2.8a)，b)は流速分布を，図‑2.8c)は 乱 れ エ ネ ル ギ 一 疋 の 分 布を，また図‑2.8d)は 乱 れ エ ネ ル ギ ー の 散 逸 率 五 の 分 布 を そ れ ぞ れ 示 し て い る . 計 算 結果の図‑2.8b)，c)と実測結果の図‑2.7a)，b)を比較すると明らかに異なっており，本 修正モデルには問題があることが分かった.

次に比較のために標準型の Rastogiand Rodiのモデルを用いて計算を行った.流速係 数が ψ=15.8の場合のモデル定数の値を，式 (2.29)，(2.30)より算定し表‑2.5に示す.

図‑2.9a)，b)に流速分布を，図‑2.9c)に 互 の 分 布 を ， ま た 図‑2.9d)に 五 の 分 布 を そ れぞれ示す.流速分布の計算結果図‑2.9b)と実験結果図‑2.7a)は， cavity内に発生す る circulationの形状や scaleが良く一致している.また，

ε

^，

t

の 分 布 の 計 算 結 果 図 ‑ 2.9c)， d)では， mixing layerにおける強い水平方向 shearにより生成された疋と五が，

cavity内に発生した circulationflowにより cavity内部に移流されながら散逸していく構 造が良く表現されている. しかし，実験結果図‑2.7b)と計算結果図‑2.9c)を比較する

と，

k

の 値 自 体 が 全 体 的 に 大 き く (2倍程度に)評価されており，標準型モデルの精度 はあまり良いものとは言えない.J alil  et al. (1994)により行われた開水路に cavityを設置 した流れ場の，標準型平面

2

次元 k‑ε乱流モデルによる計算においても'同様な結果が得 られており，これらのことから標準型モデルにはやはり特有の問題があるものと推測さ れる.

d).平面

2

次元 k‑ε乱流モデルの再修正

前述のように本章で修正された平面

2

次 元 k‑

ε

乱流モデルによる計算結果では実験結 果との良い一致は得られなかった.そこで，修正型の平面

2

次 元 k‑ε乱流モデルについ て再修正を試みることにする.

平面2次元 k‑ε乱流モデルには水平方向の shearによる乱れの productionと鉛直方向 の shearによる乱れの productionの両方がそれぞれ独自にモデル化されている.Booij  (1989)は等流場においてチューニングされた Rastogiand Rodiのモデルにおいて，これ

らの productionの weightの大きさに修正の必要性のあることを指摘した.また，流速分 布に対数分布則を仮定し，鉛直方向の渦動粘性係数には放物線分布を与えた時の水深平 均値(式 (2.13))を採用して，モデル定数 C_kをRastogiand Rodiの値に対して 0.1倍の 大きさに修正して用いている.そして，その値と Leanand Weare (1979)による水平方向 の渦動粘性係数の実験値をもとに CEを 0.023倍に修正している.この修正により

ε

，

E

に対して両方向の productionの weightが調整されたことになり，そのモデルを用いて Booij (1991)は開水路に cavityがある流れ場の再現計算を行った.その結果， cavity内の 流れや乱れエネルギー互の分布の計算結果は，共に良くその構造を再現していたが，計 算値はやはり過大に評価されていた.

Booijによる修正では，疋の鉛直 shearによる productionの強さを表すモデlレ定数 C