数理解析研究所講究録 1271 巻 2002 年 145-154145

落下する液体平面シートの振る舞い
阪大・基礎工菅健大郎 (Kentarou Kan)

阪大・基礎工吉永隆夫 (Takao Yoshinaga)

Faculty of Engineering Science, Osaka University

1

はじめに
液体シートの振る舞いはシート面上の表面張 ] 波の安定性に大きく依存することはよく知ら$f$れている. このような安定性の問題は流体$f\mathrm{J}$

学での代表的な問題\emptyset --$\vee\supset$であるばかりでなく》ポリマーシートの製造工程や化学反応等における遮蔽膜への応用においても重要である.
大変形する液体シートの解析では, シート表面が自由境界であるので境界条件が本質的に非線形となり, 解析的取り扱いは一般に困難である. かし, シートが薄い場合,\mbox{\boldmath $\zeta$}薄膜近似’ を用いて近似的にシートの運動を記述できることが知られてぃる. $[1][2]$ この近似ではシート内部での$\text{し}$諸量の変化を中心面上での変化に置き換えることにょり比較的簡単であるが強い非線形の発展方程式を導くことができる.

非定常な撹乱の伝播に関しては, 線形解析にょり一 な液体シートには二っの撹乱モードが存在することが知られている. -っは図 1.2 のようにシートの中心線に関して対称に厚みが変化することによって起こる 4 対称モード’, もうーっは図 1.3 のようにシートの厚みは一定で中心線の変位により起こる \mbox{\boldmath $\zeta$}反対称モード’ である. このような撹乱モードの伝播は重 f] 場で落下す146るシートにおいても, 重力の影響がそれほど大きくなく , 弱非線形性を考える$|$あることが示されている. さらに, このような条件 T ではシートに印加された減衰するのに対し, 反対称モード撹乱は増幅することが示されて 1 る.[4]

本研究では, 薄膜近似により得られた非線形方程式を重力場で落下する平 iる. 特}シートを伝播する大振幅の撹乱がどのように発展するかを調べ$\vee \mathrm{C}1^{\mathrm{a}}$め平面シートの厚みは下流に行くに従い薄くなり, その速度は増加する. こ状態での撹乱の伝播の問題は一般に解析的にも数値的にも取り扱いが難. しくが破断するような大振幅の変形に対する解析は行われていないようである.

2問題の定式化
図 2.1 に示すような液体シートを考える. 中心線は $y=\eta(x, t)$ で与えられている. 流体の$x,y$方向の速度成分をそれぞれ

ト厚みの半分を$a$で, 速度$x,a,\eta$[ま($t_{0}/u_{0}$で, 圧力は$a_{0}$とする.

れている. また,$u$$\sigma/a_{0}$$\sigma$と$u$$\mathrm{t},\mathrm{t}’$とし, シー

ただし, 以 T ではは$u_{0}$

で, 時間

$t$

はでそれぞれ無次元化さは流体の表面張力係数,$\rho$

は流体の密度, は重力加速度をであり, 流体は非$g$

粘性非圧縮と仮定する.
シート表面 $y=\pm a$ での圧力と表面張力による力の釣り合いより,$p \pm=\mp\frac{(\partial^{2}(_{l1}\pm a)/\partial x^{2})}{[1+(\partial(\eta\pm a)/\partial x)^{2}]^{\frac{3}{2}}}.$

ここで$\mathcal{P}+$と$p_{-}$

はそれぞれシート右側面と左側面にかかる圧力である. 重力トの基礎方程式は以 T のように示される

$j$

連続の式
$. \frac{\partial a}{\partial t}+\frac{\partial}{\partial x}.$(au)=0,

運動学的条件$\mathit{0}=\frac{\partial\eta}{\partial t}+u\frac{\partial\eta}{\partial x}$

運動方程式

.

,

147
$. \frac{\partial \mathrm{c}\iota}{cJt}+u\frac{\partial u}{\partial_{\mathrm{t}}x^{1}}=\frac{1}{\mathrm{L}\mathrm{L}e},(-.\frac{\dot{\mathrm{r}}fP}{\partial\iota I}.\cdot+\underline{.\frac{1}{)}}\frac{\triangle P}{c\iota}.‘\frac{\partial\uparrow 1}{\partial\alpha^{\backslash }})+\frac{1}{2F\iota^{2}}.$’.

$. \frac{\partial_{1’}}{\partial t}+1l.\cdot\frac{\partial_{1’}}{\partial a}$

ただし,は$\nu 1^{r}\prime e$$\backslash _{l}\mathrm{h}^{\Gamma}\mathrm{e}\mathrm{b}\mathrm{e}\mathrm{r}$$=-. \frac{1}{2\mathrm{T}’\mathrm{T},e}r\frac{\Delta\prime P}{\zeta\ell}$3(2.5)

数, $Fr$ は Froude 数であり以下のように定義される:$|/Ve= \frac{u_{0}^{2}}{(\sigma/\rho a_{0})},$

また.P.

(2.4)

$= \frac{1}{2}(p_{+}+p_{-}),$$\triangle P=p_{+}-p_{-}$$Fr= \frac{\iota\iota_{0}}{\sqrt{\underline{?}c\iota_{0\supset C}}}$.

$(^{\underline{\eta}}.6)$.

定常解
定常状態では液体シートの主流部分は対称モードとみなせるので $\eta=0$ であり, 定常状態では$/\partial’t=0$となるので式 (2.3) より $u=0$ となる. さら [こ式 (2.1) より

(3.1)

,

$p \pm=-\cdot\frac{(\partial^{2}a/\partial x^{2})}{[1+(\partial a/\partial\backslash \tau)^{2}]^{\frac{3}{}}}\underline,$となるので,$\triangle\prime P=0$,(3.2)$P=- \frac{(\partial^{2}a/\partial\backslash x^{2})}{[1+(r9a/\partial x)^{2}]^{\underline{\frac{3}{9}}}}.$となる. また, 式 (2.2) より

(a )
$u$

$/\partial x$

,

(3.3)

$–0$ となり, $a(0)=u(0)=1$ であるので以下のよう [こおくことができる:$a\tau\iota=1$

.

(3.4)

さらに式 (2.4) は$u \frac{\partial \mathrm{c}\iota}{\partial\tau}.\cdot=-\frac{1}{We}\frac{\partial P}{\partial x}+\frac{1}{2Fr^{2}}$

,

(3.5)

となり, 式 (3.3),(3.4) より.

$\}+‘\frac{1}{\underline{)}F\uparrow\backslash 2}$$- \frac{1}{a^{3}}\frac{\mathrm{r}Ja}{\partial x}=\frac{1}{\nu 7’\acute{e}}\frac{\partial}{\mathrm{r}J\alpha}.\{\frac{(\partial^{2}a/\partial x^{2})}{[1+(\partial a/\partial\alpha\cdot)^{2}]\backslash \overline{2}\wedge}$,(3.6)のような定常状態での液体シートの半厚み $a(x)$ を記述する方程式が得られる.
特に, 式 (3.6) で$\mathrm{I}/Vearrow\infty$での定常解を考えると(3.7)$- \frac{1}{a^{3}}\frac{\partial a}{\partial x}=.\frac{1}{2Fr^{2}}$となり, $a(0)=1$ より$a= \frac{1}{\sqrt{x/Fr^{2}+1}}$.

(3.8)

有限の垣 に対しては式 (3.6) を数値的に解く必要がある. 計算はルンゲークッタ法にょり, 初$\gamma e_{d}$期値 $a’(0),a”(0)$ に対してシートがどのように振舞うかを調べる. なお、定常状態でのシートの形状は対称なので, 以下の図ではシートの半厚みを正ととる.

$a$のみを描く. また, 重力加速度の方向は右向き

148
3.1

初期値の効果$\mathfrak{s},\mathrm{T}k=10.F=\sqrt{10}$とし, $a’(0)$ と $a”(0)$ の値によるシートの振る舞いに対する影響を調べるた

めに. $a’(0),a”(0)$ を変化させたときの様子を図 3.1\sim 図 3.3 に示す.
1$\approx 0.60.8$$\mathrm{q}$0.40\sear ow へ\Lambda ‘W.
5104$\cdot$

15\phi2025$\infty$

35旬45

父

図 3.2: $a’(0)=0.1.a”(0)=0.1$

図 3.1: $a’(0)=0,a”(0)=0$$\mathrm{q}$$X$

図 3.3: $a’(0)=-0.049821,a”(0)=0.007387$

図に示されているように $a’(0),a”(0)$ の値によってシート表面に現れる現象が異なる. さらに,
図 3.3 に見られるように, $a’(0),a”(0)$ の特定の値に対してはシート表面に振動が現れない場合があることを注意しておく. このようなシート表面に振動が現れない特定の条件を図 3.4, 図 3.5 に. の関数として示す. なお, $Wearrow\infty$ の場合は式 (3.8) 上り$F_{7}^{2}$$a’(0)=–1/(^{\underline{\mathrm{Q}}}F^{2},\backslash ),a’’(0)=3/(4Fr^{4})$が解析的に得られる.

$\hat{\check{\grave{.}}0\approx}$$\hat{\check{.}0\approx}$

Fr-,

Fr-,

図 3.4: 振動が現れない場合の $a’(0)$ の値

4

図 3.5: 振動が現れない場合の $a^{u}(0)$ の値

定常解の時間安定性
前節において初期値の与え方によって振動が現れない解が存在することを示した. そこで本節では振動する場合と, しない場合の二種類の解を初期条件として, 非定常解を数値的に調べることにより安定性を見る. ただし, パラメーターには $We=10,Fr=\sqrt{10}$ を用いている.

(1) 振動しない定常解の場合 $(a’(0)=-0.049821,c\iota’’(0)=0.007387)$

149$0.81.\backslash \backslash _{\sim}\backslash$0.40.6$\mathrm{t}_{3}$$0.81\backslash \backslash \backslash _{\backslash }$$-\sim_{\sim}\sim-\sim-$

5

科

10

——

15 20

2530

-\sim \sim \sim$\Im$$-\sim_{--\sim-}$$0.40.6$

35

40

45 50

05

$\mathrm{j}\mathrm{Q}$

15$x$$1_{\backslash }$

25

30

35

40 45 50$X$

図 4.1: $t=0$

図 4.

20

図 4.2: .

$t$$=\mathit{2}0$

図 4.2 は $t=0.20$ におけるシート表面の様子を示してぃる. この結果より振動しない’

解は時間的に変化せず安定であることが分かる.
(2) 振動する定常解の場合 $(a’(0)=0,a”(0)=0)$1
0.8$\backslash -\backslash arrow\backslash -\backslash -$$\mathrm{q}$$\backslash -\backslash \wedge\sim-^{-}\vee$$\mathrm{Q}$$0.40.6^{\cdot}\backslash .\backslash .\backslash \wedge.’\vee\backslash \vee\sim=\wedge\cdot.\backslash \sim\sim\backslash \backslash \wedge-."\backslash \sim-\sim--$

0510

1520

25

30

35 40

45

5科

$x$

$x$

図 4..3: $t=0$

図 4.4: $t=3$$\mathrm{q}$$X$

図 4.5: $t=20$

図 4.3\sim 図 4.5 } $\mathrm{h}t=0,3,20$ におけるシート表面の様子を示してぃる. 振動する定常解は時間がたつと下流部から振動部分がなくなる. この意味で振動する定常解は不安定である.
また, 図 4.2 と図 4.5 のシート形状がほぼ一致することから, 振動する定常解は不安定であり,時間が経つと振動しない定常解へと移行する. このことがら, 実際の定常状態でのシート形状には振動部は現れないと考えられる.

5

撹乱を与えた場合のシートの振る舞い
時間的に安定な定常解に撹乱を加えた場合のシートの振る舞いを調べる. 数値解析では対称モードの撹乱に対して中心線$\eta$は 0 のままにし, 半厚み (l($x$、ffi こ $x=0$ で以下のように撹乱を加えるとする :$c\iota(0, t)=1+\tilde{a}_{0}(t)$ここで$\tilde{a}_{0}$は. $=0$ で加えられた撹乱部であり, 適当な$\iota \mathrm{t}’$

.

(5.1)

$\omega,\omega_{1}(\omega>\omega_{1})$$\tilde{c\iota}_{0}=\frac{\overline{a}_{0}}{\cosh\{\omega_{1}(t-t_{0})\}}\cos\{\omega(t-t_{\mathrm{U}})\}$

,こ対して

(5.2)

150とする. 一方, 反対称モードの場合は $.r=0$ で半厚み

$a$は 1 のままにし, 中心線$\eta(.r.$, ffi こ以下のような撹乱を加える:$\mathrm{t}l(0.t)=\tilde{|7}\mathrm{o}(t)$ただし, 対称モードの場合と同様、適当$fx\omega,\omega_{1}$(5.3),に対して$7 \tilde{|}0=\frac{\overline{|l}0}{\cosh\{\omega_{1}(t-t_{0})\}}\cos\{\omega(t-t_{0})\}$.

(5.4)とする.
図 5.1(a),(b) にそれぞれ対称撹乱と反対称撹乱が加えられた様子を摸式的に示す. (a) 図では定常状態での半厚みを
は $\eta=0$ なので以下では$0.8,\omega=.5)$$\eta=?\tilde{|}$$c-\iota$とし、そこからのずれを$\tilde{a}_{0}$又はとする. 一方 (b) 図では中心線は定常状態でとなる.

$\mathrm{M}^{r}\prime e=10\backslash F,.2=10$

を$\tilde{a}$$\uparrow\tilde{|}0$の定常解に, 図 5.2 に示すような形状の撹乱$(\overline{c\iota}_{0}=\overline{\eta}_{0}=0.1,\omega_{1}=$として印加した場合のシートの様子を調べる.
0.1
0.05$|\epsilon\dot{e}$0-0 刀 6-0.1

05

図 5.2: 撹乱1015

20$(\overline{a}_{0}=\overline{\eta}0=0.1.\omega_{1}=0.8_{\backslash }\omega=5)$(1) 対称モード
図 5.3,5.4 こ対称モード撹乱として$\tilde{a}_{0}$を加えた場合の $t=10,20,.30$ でのシートの表面 \eta +a(実線) と中心線 7/(破線) を示す. 時間がたつにしたがいシート表面を撹乱が下流に伝播する様子がわかる. また, 中心線$\uparrow 7$

は時間が経過しても $\eta=0$ のままであることから, 反対称モード撹乱は発生していないことがわかる.

151$=\approx+\mathrm{Q}-$$arrowarrow\frac{+\mathrm{Q}}{--}$

$X$

$X$

図 5.3: シートの様子 $(t=10)$

図 5.4: シートの様子 $(t=20)$$\hat{\approx}\mathrm{P}+\mathrm{q}$$X$

図 5..5: シートの様子 $(t=30)$

次に図 $5.6\sim 5.8$ に $t=10_{\backslash }20,30$ での撹乱の発展する様子を示す. 図 $5.6\sim 5.8$ 上り下流に撹乱が伝播すると, 撹乱が減衰していくことがわかる. 図 5.9 は各場所での撹乱振幅の最大値の変化を示したものである. これからも伝播とともに撹乱が減衰してぃくことがはっきりとゎかる. この結果は弱非線形解析の結果と定性的に一致してぃる.
0. 1
0.05
00. 1
0.05
0-0.05-0. 1$\ovalbox{\tt\small REJECT}--$-O.05-0.1020

4060

1 科科

80

$\sim$

020

60$x$

80

100$X$

図 5.6: 撹乱の様子 $(t=10)$

図 5.7: 撹乱の様子 $(t=20)$

科.1

0.05

0.05
0-0.05-0. 1$.\mathrm{w}^{\eta_{4}}\ovalbox{\tt\small REJECT}$$1\mathrm{q}$

14020

科

60800.025
科

1 科科

1020304050$x$

60

70

80$X$

図 5. : 撹乱の様子 $(t=30)$

図 5.9:$\mathrm{S}$$\mathrm{x}$に対する撹乱振幅の最大値の変化(2) 反対称モード
図 5.10\sim 5.12 に反対称モード撹乱として示す (実線 $:\eta\pm a$ , 破線$:\prime l$).

$\tilde{7|}0$を加えた場合の $t=10_{\backslash }20,30$ でのシートの様子をこの場合も時間がたっにしたがいシート表面を撹乱が下流に伝播していく様子がわかる.
1$+\hat{7}arrowarrow\approx- 0.50.50|.\acute{|}^{\vee}\backslash \backslash ’\backslash ------------------------------------------------------$$+?\approx\approx-$-1010

20304050$X$

図 5.10: シートの様子 $(t=10)$60$X$

図 5.11: シートの様子 $(t=20)$

152
10.5
0$\mathrm{P}$$–\sim\cdot\backslash \vee\backslash \sim^{\prime\prime\backslash ,r_{\backslash }}\backslash$$—————–\backslash \cdot\backslash .’\backslash ’----\backslash -\backslash ----’.’\backslash \cdot\backslash \backslash \backslash \prime\prime\sim\prime\prime\backslash \cdot’-\backslash \prime\prime\backslash \backslash -\backslash \backslash -’\backslash -’\backslash ------$$.\backslash .\backslash -_{^{\backslash }\vee’}\sqrt\vee^{\wedge}\wedge\backslash$$\mathrm{B}\approx 4.5- 1$

010

2040305060$X$

図 5.12: シートの様子 $(t=30)$

図 5.13,5.15,5.17 に $t=10,20,30$ でのシート中心線 の変化を, 図 5.14,5.16,5.18 !こ $t=10,20,30$$\eta$

でのシートの半厚み$a$に現れる撹乱$\tilde{a}$の変化を示す. 図 5.14,5.16,5.18 より厚みが変化することから対称撹乱が時間が経つとともに現れてくることがわかる. つまり反対称モード撹乱のみを印加しても相互作用により対称モード撹乱が発生する事が示されている. また図 5.19 は各場所での反対称モード撹乱の振幅の最大値を示したものであり, $x=25$ 上り下流で撹乱が増幅していくことがはつきりとわかる. 弱非線形解析の結果では対称モードの発生は予測できないが, 反対称モードの増幅に関しては本解析と定性的に一致する.

$\mathrm{P}$$X$

図.5.14: 対称撹乱の様子 $(t=10)$

図 5.13: 反対称撹乱の様子 $(t=10)$

0.3$0.0\epsilon$

0.2$0.04$

0$\mathrm{A}_{\mathrm{V}\mathfrak{h}uparrow\int^{1}}^{\cdot}\Downarrow|\mathit{1}$$\mathrm{P}-.0.20.10.10$-0.04-0.302040$-0.0\epsilon_{0}$$60$

図 5.15: 反対称撹乱の様子 $(t=20)$

204060

図 5.16: 対称撹乱の様子 $(t=20)$

0.30.20.1$\mathrm{P}$

0$- 0.\iota$$- 0.2$-0.30

204060

図 5.17: 反対称撹乱の様子 $(t=30)$

図 5.18: 対称撹乱の様子 $(t=30)$

153
0.235
0.23$\mathrm{r}$

0.225
0.22

102030405060$x$

図 5.19: 各場所での撹乱振幅の最大値

5.1

シートの破断について

ここでは,We $=100,Fr^{2}=10$ の場合に対して定常解に図 5.2 で示すような撹乱を$\tilde{\eta}_{0}$として印加した場合のシートの様子を調べる. 図 5.20\sim 図 5.22 はそれぞれ $t=10,20,30$ での反対称撹乱
成分$7l$を示している. 図 5.17 と図 5.22 より垣とがわかる.

$r_{e}$が大きいほど撹乱振幅がより増幅されてぃるこ

これは弱非線形解析による結果とは逆であり, 対称モード撹乱が発生する事にょって起こるのではないかと考えられる. 図 5.23 は $t=10,20,30$ で発生した対称撹乱る.

鯢舛い討

このように時間がたつにっれて対称撹乱がかなり大きく発生してぃることがゎがる. 特に,対称撹乱のみ印加した場合には撹乱は減衰するが, この反対称撹乱により生成された対称撹乱の場合はほとんど減衰しないことに注意してお$\text{く}$.

$\mathrm{P}$$\mathrm{P}$

図5. 20: . $=10$

図 5.21: $t=20$$t$$\cdot$0.6

1.5
1$0.\mathrm{s}$0.5$\mathrm{P}$00-0.3-0.5-1-1.5$\prime^{-\prime_{1}^{\backslash }}\mathrm{V}..|||^{\mathrm{I}}.\acute{\mathfrak{s}}_{\mathfrak{i}^{\mathrm{I}}}\mathrm{I}_{1}^{1}\mathrm{I}_{1,\oint_{1}^{1}}^{\acute{|}}1\}’|’||\downarrow 1|||||\}_{1}||\mathrm{I},’$-0.6

0204060

図.5.22: $t=30$80

100010203040

図.5.23: 反対称撹乱により発生した

対称撹乱$\check{a}$の変化

図 5.24 は $t=30$ でのシート形状を描いたものである. $x=30$ 付近でシートがかなり薄くなっていることがわかる. この後 $t=39\sim 40$ でシートは破断してしまう. はじめに述べたように, これまで行われてきた弱非線形解析によると対称モードと反対称モードはカップリングしないため, 反対称モード撹乱のみを加えた場合には対称モード撹乱は発生しない. そのため反対称モー154
ド撹乱の場合, 局所的なシート厚みは定常でのそれと同じであり, シートが撹乱により破断することはない. しかし、実際には図 5.24 のように反対称モード撹乱を加えた場合には対称モード撹乱が発生することによりシート形状に薄くなる部分が現れ破断に至る可能性があると考えられる.

10$+\mathrm{I}\mathrm{P}8$-1-20102040305060$x$

図 5.24: シートの様子

6$(\mathrm{t}=30)$

結言
これまで得られた結果をまとめると以下のようになる:・落下する平面液体シートの定常解を調べた結果, 初期値により空間的に振動する解としない解が得られることがわかった. 特に 4 振動しない’ 解は特定の初期値に対してのみ得られることがわかった.

・振動する定常解と, しない定常解の時間安定性を調べた結果, 振動しない解が安定であり,振動する解は時間とともに振動がなくなっていくことがわかった.
・振動しない定常解に, 上流端で撹乱を加えた場合, 以下の結果が得られた

1. 対称撹乱のみを加えた場合には反対称モード撹乱は発生しない.
2. 反対称撹乱のみを加えた場合には対称モード撹乱が引き起こされる.
3. 対称モード撹乱は下流にいくにしたがい撹乱振幅は小さくなる.
4. 反対称モード撹乱では撹乱振幅が増幅する場合もある. また. 誘起された対称撹乱によってシートが破断する可能性がある.

参考文献[1] C.P.Lee&T.G.VVang: Phys. Fluids 29, 2076(1986).
[2]:J.Phys. Soc. .Jpn 70$\mathrm{T}.\mathrm{Y}’ \mathrm{o}\mathrm{s}\mathrm{h}\mathrm{i}\mathrm{n}\mathrm{a}\mathrm{g}\mathrm{a}\ \mathrm{I}_{\dot{\overline{\mathrm{C}}}}]_{\backslash }’\mathrm{o}\mathrm{t}\mathrm{a}.\mathrm{n}\mathrm{i}$$(^{\underline{7}}001),$

$70$

.

[3] Y.Inoue:J.Ph.vs. Soc. Jpn 38 (1975), 598.

[4] 中村康輔: 重力場で鉛直に落下する流体シート上の表面張力波(大阪大学基礎工学部修士論文,1993).