成層流体中に発生する内部モード波の非線形挙動に関する研究

成層流体中に発生する内部モード波の非線形挙動に関する研究

Non-linear behavior of internal mode wave generated in stratified fluids

  土木工学専攻 

19

号  佐藤  学

SATO Gaku

１． はじめに

密度成層の流れに与える効果を積極的に利用したものの 一つとしてダム貯水池や火力及び原子力発電所の取水に 等に行われる選択取水（selective withdrawal)がある．ダム 貯水池ではプランクトンの異常増殖や洪水時における長 期濁水現象により放流水が濁る問題が生じる．また夏場 に下層の冷たい水を下流に放流すると稲作などの成長阻 害の要因の一つとなると考えられる．そのため任意の層 の水を必要量最適に確保することが選択取水の問題であ り，本研究の目的は非線形の干渉を考慮した取水量の定 量的な把握を目的としている．密度成層した貯水池から の取水に関する研究は実験的にも理論的にも多くある．

特に密度成層が線形の場合は

Yih¹⁾

による解析的研究や

Debler²⁾

の実験的研究により，取水口上部に淀み層(静止

領域)が現れる臨界値を示す臨界フルード数

Fr

が

1/π以

下であることを

Yih

が解析により求め，Debler は

2

次元 水槽の一方の端に

line sink

を設けてこれより一定流量を 取水する下層取水の実験より独自に臨界フルード数

Fr

≦

0.28

を得ている．さらに，内部フルード数が小さい時 の流れに関して，

Koh³⁾

の解析解があり，実験との良い一 致を得ている．また密度成層した貯水池からの取水時に は特異な波動が存在することが

Pao⁴⁾

らによって示され ている．この波動に着目した研究として山田ら

^{5) 6)}

の研究 がある．山田らは移流項を部分的に残した

Oseen

近似に 基づく理論展開を行い，ラプラス変換を用いた線形の解 析解を示している． 山田らは合わせて実験も行っており，

定性的な選択取水の流れの状況を示すとともに，線形理 論解の限界について言及している． 選択取水の問題は線 形理論解で現象を十分に説明できるが定量的な把握には 非線形での解析が必要であり，純解析的に解くには限界 がある．したがって数値計算に頼らざるを得ない現状が あり，定量的な把握には非線形性を考慮した数値計算を 行う必要がある．本研究は選択取水の問題として流体内 に発生するモード波の非線形挙動に着目し， 非線形性を

考慮した数値計算を行い，山田らにより示されている線 形理論解と数値計算結果との比較を行なった． あわせて，

非線形の現象をある時点の線形モードを初期条件として 取り扱い，線形解を順々に繰り込むことで非線形解を表 現する繰り込み理論を選択取水の非線形現象に適用し，

その実用性を考察した．

2．理論

山田らの研究は基礎方程式として，運動方程式 (1)(2)及 び密度の保存式として式(3)，連続式として式(4)を用いて いる．

式(1)(2)より，圧力

P

を消去すると式(5)となる．ここに

ω:∂u/∂x-∂v/∂y(渦度）

，

D/Dt=∂/∂t-U∂/∂x

である．

さらにラグランジュの流れ関数を用いると

u=∂ψ/∂y

，

v=-∂ψ/∂x

，

ω=

∇

²ψ

より式(5)は式(6)になる．

2 2

0 x y

g Dt D

∂ +∂

∂

= ∂ρ ω ρ

ω ₍₅₎

2 2

0 2

x y g Dt D

∂ + ∂∇

∂

= ∂

∇ ρ ν ψ

ρ

ψ ₍₆₎

貯水池

•プランクトンの異常増殖,洪水時 に生じる長期濁水現象による下流 域での浄水障害

•冷水放流に伴う農業への冷害

火力及び原子力発電所

•河床からの異物混入,表層からの 温度の高い水,クラゲ等 の海洋生物 混入による冷却水の冷却効率の低下

任意の層の水を目的に応じて選択的に取水することが重要

取水

濁水

特定の層からのみ取水 取水における問題点

選択取水

=0

∂ + ∂

∂

− ∂

∂

∂

v y U x t

ρ ρ ρ

=0

∂ + ∂

∂

∂ y v x u

2 2

0

1

y u x

p x

U u t u

∂ + ∂

∂

− ∂

∂ =

− ∂

∂

∂ ν

ρ

2 2

0 0

1

y g v y p x

U v t u

∂ + ∂

∂ −

− ∂

∂ =

− ∂

∂

∂ ν

ρ ρ ρ

(3)

(4) (1)

(2)

図-1 取水における問題点 

式(3)(6)から密度ρを消去し式(7)を得る．

密度分布は線形と仮定し，密度勾配を

ε

とすると式(8)と なる． ここに ε

=(ρ0-ρs)/ ρ0d

，

ρs

は水面の密度である．

式(8)がこの解析の基本式である．無次元化は長さの基準 としては水深

d，時間の基準は密度成層の勾配ε

から決 まる時間スケール

T(Brunt-Vaisara振動数の逆数)を用いる

と式(9)になる ．

ξ=x/d

，

η=y/d

，

τ=t/T(T=1/ ε_g)，

ただし，

φ=

ψ

T/d²

，レイノルズ数

Re= d²/νT

， 内部フル

ード数

Fr=UT/d

とする．式(9)に対する境界条件及び初期

条件としては全壁面で滑りを許すと,以下の様になる．

式(11)に式(10)の(ⅳ)を考慮してラプラス変換を行うと式

(12)となり，ϕ

を

y

方向にフーリエ変換を行い，その中の

任意の

1

項に関して式(12)の解の形を式(13)のように与 える．ここに

ϕ:ϕ

の偶関数，

An:境界条件より決まる定数，

k=nπ(n=1

，

2

，

3

，････

)である．

4 4 2

2 2

2 ( )

Re ) 1

( η

ϕ ξ ξ

ϕ ϕ

ξ ∂

∂

∂

− ∂

∂ = +∂

∂ ∇

−Fr ∂ p Fr

p

) exp(

sin cos )

( η ξ

ϕ kb k ak

Kp p Aⁿ

n = −

式(13)を式(12)に代入し，若干の計算の後 α を決める方程 式(14)を得る．

 

(p+Fr⋅k⋅a)²(a²−1)+a²=1/Re(p+Fr⋅k⋅a)

これから根

α

を求めることは困難である．よって，極限 の場合を考えることにする．

a)．Fr→０

，

Re→∞(移流項

無視，非粘性)この時(15)式の根は α

=

±

P/ P²+1

となり，

これを式(13)に代入すると式

(15)を得る．

式(15)のラプラス逆変換をとれば式(16)を得る．

ところで，式(15)の解はある種の波動現象を表しており，

この波動は固有関数（モード）によって決まる伝播速度

Cn

をもっており，

Cn

＝

1/n

π と示されている． C

n

を次 元をもつ形で表すと

d/(n

π

T)での上流に遡上し得るモー

ド(n)は

d/(n

π

T)

＞

q/d

を満たし，これを内部フルード数で 表すと

1/(n

π

)

＞

Fr

となる．つまり，モード(n)によって は波速が平均流速よりも小さくなり上流に遡上できない ものも存在する．すべてのモードが出揃うと以後は定常 な流れとなる．式(16)よりモードの違い(Fr 数の違い)に伴 う流速分布を図示すると図-2 になる．次々に到達するモ ードによって流速分布が鋭くなっていることが分かる．

2．繰り込み理論 

非線形現象を純解析的に解ける現象は稀であり，種々の 非線形現象は解析的には近似やオーダー比較などにより 現象が理解される．しかし現実には現象の定性的・定量 的な評価には数値計算に頼る現状がある．このような中 で，本研究では線形解析の範囲の中で非線形性の取り扱 いに着目し，波動の非線形性を表現するために非線形現 象に対して線形解の繰り込みにより非線形解を導く繰り 込み理論を提案する．繰り込み理論とはある時刻におけ る場の線形解を初期条件として，次の場を解き順々に線 形解を繰り込むことで非線形解を表現する理論である．

本研究では繰り込み理論を選択取水時に生じる非線形挙 動に適用した．本研究における繰り込み理論の考え方を

[ ]

0 , 0 )

(

1 0 ) 0

(

1 , 0 , 1 0 ) ( )

( ) 0 ( ) (

, 0 , 0 ) ( 1 , 0 , 0

, 0 ) (

2 2

<

=

∇

∂ =

∂∇

⎭⎬

⎫

⎩⎨

⎧

<

<

<

= <

=

≤

≤ +

−

=

=

∞

→

≤

≤

=

=

∞

<

≤

=

τ τ ϕ

ϕ

η η η

ξ η τ

η η ξ

ϕ

ξ ξ ϕ

η ξ ϕ

 

ⅳ

   

  

  

ⅲ

１

ⅱ

ⅰ

b Fr F b

H Fr F

2 2 2 2 2 2 2 2

y Dt

D g x

Dt D

∂

∇

= ∂

∂ + ∂

∇ψ ε ψ ν ψ (8)

1 exp

sin cos )

(

2+

−

=

p k pk

Kp kb p Aⁿ

n

η ξ ϕ

) , , ( )

0 η (τ ϕ ξ η τ

ϕ =− FrH +

⎩⎨

⎧

⎭⎬

− ⎫

−

= ∑

=

πη π π

η

ϕ n b n

Fr n

N

n

sin 1 cos

2

1 0

2 2 2 2

2 2

0 2 2 2

)

( Dt y

D y x x y x g Dt D

∂

∇ + ∂

∂

∂

∂

∂ +∂

∂

∂

∂

= ∂

∇ ψ ρ ψ ρ ν ψ

ρ

ψ (7)

(11)

(12)

(13)

(14)

(16)

) Re (

) 1

( ₂

2 2 2

2 2 2

η ϕ τ ξ ϕ ϕ ξ

τ ∂

∇

= ∂

∂ +∂

∂ ∇

− ∂

∂

∂

D

Fr D ⁽⁹⁾

(10)

(15)

図-2 山田の線形理論解から求めた Fr 数の違いに伴う流速分布

0 0.5 1

1.0

η u/u_max

Fr=0.159 Fr=0.055 Fr=0.03

0.2

(n=2) (n=6) (n=10)

図-3 に示す． 選択取水時にはフルード数によって決定さ れるモード波が存在することが示されている．取水開始 直後，貯水池内全域に鉛直方向に一様な流速分布をもつ ポテンシャル流れが発生する．その場を初期条件とし，

第

2

モードの場において分離流線で囲まれる取水層が形 成され，この層内での流速を

V1

，取水層厚を

δ

とする．

取水流量は常に一定と考えられるので第

1

モードの場と 第

2

モードの場で連続式を用いることで第

2

モードのフ ルード数を求める．さらにその内部フルード数を次の場 へ順々に繰り込むことによって式(17)を得ることができ る．ここで，山田らの研究より

n

次のモードにおける臨 界フルード数は理論的に

Fr＝1/2πと示されているので

式(18)を得る．その両辺に対数をとり，

n

を求める式(19) を得ることができる．導出された

n

を

n*=n+1

に繰り込 むことによって取水層厚の式を得る．

繰り込みの手順を以下にまとめた．

(1)

第

1

モードの流速，全水深，フルード数を初期条件 として第

2

モードの取水層内のフルード数を計算する．

(2)

同様に第

2

モードの諸条件を新たな初期条件として 第

3

モードのフルード数を計算する．

(3)

これを

n

次まで繰り込むことにより式

(17)を得る．

(4)

中層取水の場合の限界フルード数はFr

_n =1/2πなので

式(18)の関係を得ることができ，式

(18)の両辺に対数をと

り整理すると

n

に関する式(19)が得られる．

(5)

導出された

n

を

n*=n+1

に代入し， 取水幅を示す 式(20)の

n*に代入することでn

次モードでの取水幅を求 める．ここで

δ/d

は過去の実験より

δ/d=0.5

と理論的に導 かれているので，その値を代入する．

本理論におけるパラメータは粘性や非線形の効果を表現 する

δ/d

及び取水形式によって変化するパラメータ

a

で ある．パラメータ

a

は取水形式によって異なり，既往の 研究より求められている臨界フルード数

Fr

となる時の 取水層厚の比(d/h)が

1

となるように

a

を決定する．

粘性の効果や非線形の効果を表現するためには数値計算 の結果や実際の観測値などにより，a や

δ

を決定する．

Yih

の解析解や

Debler

の実験値と繰り込み理論から求め た解の整合性を確かめるためにパラメータを適当に選び それぞれ図示したのが図-4 である．Yih の線形解は粘性 を考慮していないので，理論より導出された取水幅比

δ/d=0.5

を用いることで，線形理論から導出された

Yih

の

解析解と一致することが分かる．また

Debler

が実験より 求めた近似式との比較では内部フルード数の大きい実験 値で，取水層厚比が一致している．これは

Debler

の実験 値が内部フルード数の大きい値から得られているためで あろう．

3．数値計算

3．1

基礎方程式  基礎方程式としてレイノルズ方程式及 び連続式，さらに移流拡散方程式を用いている．

3. 2

初期条件及び境界条件 

本研究の数値実験で設定した貯水池を図-5に示す． 2次 元取水を再現するため，下流端の河床から鉛直上方35m 地点に取水口を設定した．初期条件としては貯水池内全 域の各地点における流速を0[m/s]とする．(u=0[m/s]，

v=0[m/s]，w=0[m/s]:0≦x≦100[m]， 0≦y≦3[m]，0≦z≦

70[m])．また，貯水池内の密度分布としては表面密度ρ_s

＝

999.7042[kg/㎥]，底面密度ρ0

＝

1019.484[kg/㎥]とし，貯 1 0

2 =

∂ +∂

∂

∂

j j

s x

u t P ρc

) ( )

( _S

j j j

j x

D S Su x

x t S

∂

∂

∂

= ∂

∂ + ∂

∂

∂

3 ) } 2 {

1 ( ) (

T k

x v y u g x

x P x

u u t u

ij i j j i j i i j

j

i i − δ

∂ +∂

∂

∂

∂ + ∂

∂ +

− ∂

∂ = +∂

∂

∂ ν

ρ

0.20.40.60.81 0.20.40.60.81 0.20.40.60.81

繰り込み理論の考え方

① ② ③

第１モードを初期条件とし 第２モードのFr数を計算

第２モードを初期条件とし 第３モードのFr数を計算

V₀ V₁ V₂

d

δ

^d

δ2

第１モード 第２モード 第３モード

図-3 繰り込み理論の考え方 

図-4 繰り込み理論と既往の研究の比較 

0.1 0.2 0.3 0.4

InternalFr0 0.2

0.4 0.6 0.8 1 1.2

selective withdrawal breadth

Yih(1958) Debler(1959)

δ/d=0.50 δ/d=0.526

(21) (20)

(22)

n n

Fr d Fr

2

0 ⎟

⎠

⎜ ⎞

⎝

= ⎛

δ ^2π

2 1

0 ⎟⎟⎠ =

⎜⎜ ⎞

⎝

⎛d ⁿ

Fr δ

[ ] _⎟⎟

⎠

⎜⎜ ⎞

⎝

× ⎛

=

d Fr

n δ

π) /log 2 log(

2 1

0

(4) (5)

(6)

* n

a d ⎟

⎠

⎜ ⎞

⎝

⎛δ

(7)

(17) (18)

(19) (20)

水池表面から線形に増加する密度分布を仮定した．本論 文では内部フルード数を変えて行った計算のうち，内部 フルード数の小さい計算ケースを掲げている．それぞれ の計算条件は表-1に示す．Run1〜3は図-5に示す貯水池 での計算を行い，Run4，5では内部フルード数を大きく したため全水深を90mに変えて計算を行った．また，上 流端境界条件には密度成層を破壊しない流速で取水口か ら取水される流量と等量の流量を与えている．さらに全 ての境界をスリップ条件とした．

3．3

数値計算結果

Run1の条件の下，取水口から距離40m離れた地点での各

時刻(5分間隔)の流速分布を図-6に示す．次々に到達する モードによって最大流速が増していることが分かる．ま た線形解析で見られた取水層厚の減少はある時間をもっ て止まり，一定の値に収束することが分かる．したがっ て，線形解析では内部フルード数が小さい解析において 取水層厚が際限なく細くなったが，非線形を考慮するこ とで細くならないことを示した．

3．4

繰り込み理論解と比較

内部フルード数を変えて行ったRun1〜Run5の結果及び

Deblerの実験値を用いて繰り込み理論の妥当性を考える．

パラメーラδをδ=0.5，

0.5125，0.525に変化させた取水層

厚比と内部フルード数の関係を図-7に示す．数値計算で 得られた取水層厚比と実験値とが同じライン上に並び繰 り込み理論の妥当性を示した．これにより数値計算に頼 らずに非線形を考慮した取水層厚比が求まる．

4．まとめ

線形理論解で選択取水の現象を説明することは十分可能 であるが，内部フルード数が小さい実際の取水では非線

形の干渉や粘性の効果が大きく，定量的な解析には非線 形の領域まで解析を必要とする．現実には数値計算に頼 らざるを得ないが，本研究ではある時点の線形モードを 初期条件として取り扱い，線形解を順々に繰り込むこと で非線形解を表現する繰り込み理論を示し，過去の実験 値と数値計算からその手法の妥当性を示した．

＜参考文献＞

1)Yih，C.S., "On the Flow of a Stratified Fluid," Proceedings of the Third National Congress of Applied Mechanics, (1958).

2)Debler，W.R., "Stratified Flow into a line Sink," Proc.

Am.Soc. Civ. Engrs., J.of the Engineering Mechanics Division, Vol. 85, No.EM3,July 1959, pp. 51-65.3)Koh，T.W.  J.Fluid Mech., 24, 1966. 4)Pao，H.P. and Kao，T.W. : Dynamics of establishment of selective withdrawal from a line sink, Part 1,Theory, J. Fluid Mech. 65, 1974. 5)吉川秀夫,山田正,水谷

俊孝,第22回水理講演会論文集,pp95-100,1978年. 

6)山

田正,吉川秀夫,水谷俊孝,土木学会論文報告集 第280 号,pp 67-79,1978年12月. 

図-5  数値実験で設定した貯水池

 

100m 70m

ρ₀ ρ_s

取水口 35m 3m ｑ

y z x

100m 70m

ρ₀ ρ_s

取水口 35m 3m ｑ

y z x

–0.01 0 0.01 0.02

0 10 20 30 40 50 50 70

水深(m)

流速

(m/s)

t=5min t=10min t=15min t=20min t=5min t=10min t=15min t=20min

 

図-7 数値計算結果，Debler の実験値，繰り込み理論の比較 

0.5 1.0 1.5 2.0

Internal Fr0 0.1

0.2 0.3 0.4 0.5 0.6 0.7

selective withdrawal breadth

(×10^-2)

Debler(1959) 数値計算結果(Run1〜5)

δ/d=0.50 δ/d=0.5125 δ/d=0.525

表-1  数値計算条件 

7.5×10^-3 21.75

2.16×10^-4 2.79

Run5

5.0×10^-3 21.75

2.16×10^-4 1.863

Run4

2.71×10^-4 18.96

2.84×10^-4 0.07

Run3

2.17×10^-4 18.96

2.84×10^-4 0.056

Run2

1.55×10^-4 18.96

2.84×10^-4 0.04

Run1

Fr(UT/d) T(sec)

ε(1/m)  q(m³/sec)

7.5×10^-3 21.75

2.16×10^-4 2.79

Run5

5.0×10^-3 21.75

2.16×10^-4 1.863

Run4

2.71×10^-4 18.96

2.84×10^-4 0.07

Run3

2.17×10^-4 18.96

2.84×10^-4 0.056

Run2

1.55×10^-4 18.96

2.84×10^-4 0.04

Run1

Fr(UT/d) T(sec)

ε(1/m)  q(m³/sec)

–0.01 0 0.01 0.02

0 10 20 30 40 50 50 70

水深

(m)

流速

(m/s)

t=5min t=10min t=15min t=20min

–0.01 0 0.01 0.02

0 10 20 30 40 50 50 70

水深

(m)

流速

(m/s)

t=5min t=10min t=15min t=20min

図-6 数値計算より導いた各時刻の流速分布(x=40m 地点)