数理解析研究所講究録
第 1800 巻 2012 年 161-170161

相互作用を考慮した
非球形気泡群の非線形振動解析
大分大工栗原 央流 (Eru Kurihara)

Faculty of Engineering, Oita University

1

はじめに
高速度カメラによる撮影技術の進歩により，周囲圧力の急激な変動にともなう微小な気泡群の力学的な挙動が実験的にあきらかにされつつある [1]. 複数の気泡からなる系の力学挙動は単一気泡の場
合とは異なり，気泡間の相互作用が気泡の並進運動や変形といった特徴的なふるまいに対して大きな影響をおよぼす．とくに気泡間距離が小さい場合や気泡が固体壁面などの境界近傍に存在するときは気泡間の相互作用が気泡の変形や崩壊時の挙動に対して支配的な役割を果たす．隣接する 2 気泡の変形にともなって対向噴流が形成されること，壁面近傍の気泡の崩壊の際に壁面へ向かって高速の噴流が形成されることは実験的によく知られており [2], このような高速の噴流がキャビテーションによる壁面の壊食に一定の寄与を果たしていると考えられている [3].
非球形気泡の崩壊にともなう高速噴流の形成は，衝撃波による体内結石の破砕 (Shock wave lithotripsy)に代表される超音波 (衝撃波) の医療への応用において病変組織の破壊に寄与する一方で，周辺組織の損傷の原因ともなる．このため噴流の形成を含む非球形気泡の有限振幅の振動に関する理論的な研究がさかんに行われているが，そのほとんどが境界要素法に代表される数値的な研究であり [4,5],
計算機資源の制約から大規模な気泡群のふるまいを予測することは困難である．このような気泡群に対しては気相と液相を平均化したモデル方程式による解析がなされているが，気泡崩壊にともなう高速噴流が毛細血管や細胞へ与える影響といった局所的な現象をとらえることはできない．
著者らは気泡境界の多重極展開と解析力学の手法を用いて 2 気泡間の相互作用を考慮した非球形気
泡の運動を記述する力学方程式を導出し，これにより気泡の体積振動，並進運動，形状変形といった基本的な力学挙動の解析を行った [6]. 本研究では，高次の形状変形モードが気泡の運動に与える影響の評価するとともに，従来の手法の 3 次元問題への拡張を試みた．

2問題の定式化

図 1 に相互作用するふたつの非球形気泡の模式図を示す．球面調和関数を用いて変形する気泡境界面を表現するために，それぞれの気泡に対して気泡の中心点を定義し，2 つの気泡の中心間距離 (分離距離) を $D$ とする．

気泡表面の位置はそれぞれの気泡に対して中心からのベクトル $r_{si}(i=1,2)$ によって表される．ふたつの気泡からなる系では，気泡

図 1:2 つの非球形気泡

162の中心点を結ぶ直線を対称軸 (X 軸) とすることにより軸対称の 2 次元問題としての扱いが可能となる．この場合，気泡境界を表す位置ベクトル，はルジャンドル関数$P_{n}(cos\theta_{i})$

を用いて

(1)

$r_{si}(t, \theta_{i})=[R_{\dot{\eta}}(t)P_{0}(cos\theta_{i})+\sum_{n=2}^{\infty}A_{ni}(t)P_{n}(cos\theta_{i})]e_{r^{}i}$に対応する気泡 の
ただし，$i=1,2$ . ここで，瑞は気泡 の代表半径， はルジャンドル関数
は軸から測った角形状変形モードの振幅である．また， は気泡の中心を示す位置ベクトル，$A_{ni}$$i$$i$$P_{n}$$\theta_{i}$$r_{i}$$X$

度， は半径方向の単位ベクトル， は気泡 の表面における外向き法線ベクトルである．ここで，$i$$n_{i}$$e_{ri}$

次の 2 つの無次元パラメータ$\epsilon$

と

$\mu$

を導入する．
$\epsilon=O(\frac{A_{ni}}{R_{\dot{\eta}}})\ll 1$(2)$\mu=\frac{R_{\dot{\eta}}}{D}\ll 1$(3)式 (2) と式 (3) は，それぞれ球の平均半径に対して形状変形の振幅が十分に小さいこと，気泡同士は十分に離れていることを表す．本研究では，これら 2 つのパラメータそれぞれについて 2 次の項までを考慮した解析を行う．

いま，気泡の周囲の液体は非圧縮かつ粘性を無視できるものとする．初期に周囲の液体は静止一様状態にあったとすれば，単一の気泡の周囲のながれは速度ポテンシャル によって次のように定義さ$\phi$

れる．

(4)$v=grad\phi$(5)$\phi_{i}=\sum_{n=0}^{\infty}\frac{a_{ni}}{r_{i}^{n+1}}P_{n}(\cos\theta_{i})$ここで， は速度ポテンシャル$\phi_{i}$$\phi$に対する気泡 による寄与であり，$i$$a_{ni}$は式 (5) における $n$ 次のモードの振動に対応する係数である．

運動学的境界条件は，気泡境界の移動速度$\dot{r}$$i$。

と気泡表面での流速 $v_{i}=grad\phi_{i}$ を用いて

th si. $n_{i}=v|_{r=r_{si}}\cdot n_{i}=grad\phi_{i}|_{r_{si}}$ . ni

(6)$\dot{r}_{si}=[\dot{R}+\sum_{n=2}^{\infty}A_{ni}P_{n}(\cos\theta_{i})]e_{r^{l}i}$(7)と表される．ただし，である． は$n_{i}$(8)$n_{i}= \frac{N_{i}}{|N_{i}|}$瓦$= \frac{\partial r_{s}i}{\partial\theta_{i}}\cross\frac{\partial r_{i}}{\partial\varphi_{i}}=r_{si}^{2}\sin\theta_{i}e_{\dot{n}}-r_{si}\sin\theta_{i}\frac{\partial r_{si}}{\partial\theta_{i}}e_{\theta i}$(9)で与えられる．ここに，腕は対称軸に対する周方向の角度である．境界条件 (6) に式 (5) を代入する
ことにより，単一の気泡に対する速度ポテンシャル$\Phi_{i}$の係数$a_{ni}$が決定される．2 気泡を含む液体の163
流れを記述する速度ポテンシャルは(10)$\phi=\sum_{i=1}^{2}[\sum_{n=0}^{\infty}\frac{a_{ni}}{r_{i}^{n+1}}P_{n}(\cos\theta_{i})]$によって与えられるが，これが式 (4),(6) を満足するように次の関係を用いて係数$a_{ni}$を修正する．

$\frac{P_{n}(\cos\theta_{2})}{r_{2}^{n+1}}=\frac{1}{n!D^{n+1}}\sum_{m=0}^{\infty}\frac{(-1)^{m}(m+n)!}{m!}(\frac{r_{1}}{D})^{m}P_{m}(\cos\theta_{1})$(11)$\frac{P_{n}(\cos\theta_{1})}{r_{2}^{n+1}}=\frac{(-1)^{n}}{n!D^{n+1}}\sum_{m=0}^{\infty}\frac{(m+n)!}{m!}(\frac{r_{2}}{D})^{m}P_{m}(\cos\theta_{2})$(12)これにより相互作用する 2 気泡のまわりの流れを記述する速度ポテンシャルの係数に対する 2 次の項までの近似で以下のように定まる．

$a_{ni}$が$\epsilon$

および

$\mu$

$a_{0i}=-h\dot{R}_{i}$$- \frac{2A_{2i}\dot{A}_{2i}}{5}-\frac{2A_{3i}\dot{A}_{3i}}{7}-\frac{2A_{4i}\dot{A}_{4i}}{9}-\frac{2A_{5i}\dot{A}_{5i}}{11}$$- \frac{A_{2i}^{2}\dot{h}}{5h}-\frac{A_{3i}^{2}\dot{\hslash}}{7R_{i}}-\frac{A_{4i}^{2}\dot{h}}{9\hslash}-\frac{A_{5i}^{2}\dot{R}_{\dot{\triangleleft}}}{11R_{i}}$(13)$a_{1i}=- \frac{RUi}{2}+\frac{3A_{2i}Ui}{10}$$- \frac{6A_{3i}\dot{A}_{2i}}{35}-\frac{27A_{2i^{\lrcorner}}4_{3i}}{70}-\frac{5A_{4i\lrcorner}4_{3i}}{42}-\frac{2A_{3i}\dot{A}_{4i}}{7}-\frac{A_{5i}A_{4i}}{11}-\frac{5A_{4i}\dot{A}_{5i}}{22}$$- \frac{183A_{2i}^{2}Ui}{350h}+\frac{18A_{4i}A_{2i}Ui}{35lb}-\frac{127A_{3i}^{2}Ui}{245h}-\frac{12A_{3i}A_{2i}\dot{R}}{35Pt}-\frac{5A_{4i}A_{3i}\dot{R}}{21\hslash}-\frac{2A_{5i}A_{4i}\dot{R}_{i}}{11\hslash}$$+ \frac{40A_{5i}A_{3i}Ui}{77h}-\frac{241A_{4i}^{2}Ui}{462R_{i}}+\frac{874A_{5i}^{2}Ui}{1573h}$$+ \frac{(e_{ij}e_{x})}{D^{2}}[-\frac{hR_{j}^{2}\dot{R}_{j}}{2}-\frac{183A_{2i}^{2}R_{j}^{2}\dot{R}_{j}}{350\hslash}+\frac{18A_{4i}A_{2i}R_{j}^{2}R_{j}}{35h}-\frac{127A_{3i}^{2}R_{j}^{2}R_{j}}{245h}$$+ \frac{40A_{5i}A_{3i}R_{j}^{2}R_{j}}{77R}-\frac{241A_{4i}^{2}R_{j}^{2}\dot{R}_{j}}{462h}+\frac{874A_{5i}^{2}R_{j}^{2}R_{j}}{1573h}+\frac{3A_{2i}R_{j}^{2}R_{j}}{10}$$- \frac{A_{2j}^{2}k\dot{R}_{j}}{10}-\frac{A_{3j}^{2}R_{\dot{4}}R_{j}}{14}-\frac{A_{4j}^{2}hR_{j}}{18}-\frac{A_{5j}^{2}\hslash\dot{R}_{j}}{22}$$- \frac{A_{2j}hR_{j}\dot{A}_{2j}}{5}-\frac{A_{3j}hR_{j}\dot{A}_{3j}}{7}-\frac{A_{4j}hR_{j^{A}}4_{4j}}{9}-\frac{A_{5j}P*R_{j}\dot{A}_{5j}}{11}]$

(14)

164

$a_{2i}=- \frac{\hslash_{4}\dot{4}_{2i}}{3}-\frac{2A_{2i}\dot{R}}{3}+\frac{3A_{3i}Ui}{7}$$- \frac{2A_{2i}\dot{A}_{2i}}{7}-\frac{4A_{4i}\dot{A}_{2i}}{63}-\frac{11A_{3i\lrcorner}\dot{4}_{3i}}{63}-\frac{25A_{5i}\dot{A}_{3i}}{693}$$- \frac{8A_{2i}\dot{A}_{4i}}{21}-\frac{260A_{4i}\dot{A}_{4i}}{2079}-\frac{200A_{3i}\dot{A}_{5i}}{693}-\frac{125A_{5i}\dot{A}_{5i}}{1287}$$- \frac{2A_{2i}^{2}\dot{\hslash}}{7\hslash}-\frac{20A_{4i}A_{2i}\dot{\hslash}}{63\hslash}-\frac{10A_{3i}^{2}\dot{\hslash}}{63\hslash}-\frac{50A_{5i}A_{3i}\dot{\hslash}}{231R_{i}}-\frac{20A_{4i}^{2}\dot{\hslash}}{189\hslash}-\frac{100A_{5i}^{2}\dot{\hslash}}{1287\hslash}$$- \frac{233A_{3i}A_{2i}Ui}{245\hslash}+\frac{10A_{\Re}A_{2i}Ui}{11h}-\frac{4430A_{4i}A_{3i}Ui}{4851\hslash}+\frac{6125A_{5i}A_{4i}Ui}{14157\hslash}$$+ \frac{(e_{ij}\cdot e_{x})}{D^{2}}[-\frac{233A_{3i}A_{2i}R_{j}^{2}\dot{R}_{j}}{245\hslash}+\frac{10A_{5i}A_{2i}R_{j}^{2}R_{j}}{11\hslash}$$- \frac{4430A_{4i}A_{3i}R_{j}^{2}\dot{R}_{j}}{4851\hslash}+\frac{6125A_{5i}A_{4i}R_{j}^{2}\dot{R}_{j}}{14157\hslash D^{2}}+\frac{3A_{3i}R_{j}^{2}\dot{R}_{j}}{7}]$(15)$a_{3i}=- \frac{\hslash\dot{A}_{3i}}{4}-\frac{A_{3i}\dot{\hslash}}{2}j+\frac{A_{4i}Ui}{2}-\frac{9A_{2i}Ui}{10}$$- \frac{A_{3i}\dot{A}_{2i}}{5}-\frac{17A_{2i}\dot{A}_{3i}}{60}-\frac{5A_{4i}\dot{A}_{3i}}{44}-\frac{2A_{3i}\dot{A}_{4i}}{11}-\frac{A_{5i4}\dot{4}_{4i}}{13}-\frac{25A_{2i}\dot{A}_{5i}}{66}-\frac{115A_{4i}\dot{A}_{5i}}{858}$$- \frac{17A_{3i}A_{2i}\dot{\hslash}}{30h}-\frac{10A_{5i}A_{2i}\dot{\hslash}}{33\hslash}-\frac{7A_{4i}A_{3i}\dot{\hslash}}{22\hslash}-\frac{7A_{5i}A_{4i}\dot{\hslash}}{33\hslash}$$- \frac{9A_{2i}^{2}Ui}{25\hslash}-\frac{56A_{4i}A_{2i}Ui}{55\hslash}-\frac{117A_{3i}^{2}Ui}{385\hslash}+\frac{585A_{5i}A_{3i}Ui}{847\hslash}-\frac{118A_{4i}^{2}Ui}{429\hslash}+\frac{9A_{5i}^{2}Ui}{26\hslash}$$+ \frac{(e_{ij}\cdot e_{x})}{D^{2}}[-\frac{9A_{2i}^{2}R_{i}^{2}R_{j}}{25\hslash}-\frac{56A_{4i}A_{2i}R_{j}^{2}\dot{R}_{j}}{55\hslash}-\frac{117A_{3i}^{2}R_{j}^{2}\dot{R}_{j}}{385h}+\frac{585A_{5i}A_{3i}R_{j}^{2}R_{j}}{847\hslash}$$- \frac{118A_{4i}^{2}R_{j}^{2}R_{j}}{429\hslash}+\frac{9A_{5i}^{2}R_{j}^{2}\dot{R}_{j}}{26\hslash}-\frac{9A_{2i}R_{j}^{2}R_{j}}{10}+\frac{A_{4i}R_{j}^{2}R_{j}}{2}]$

(16)

165

$a_{4i}=- \frac{h\dot{A}_{4i}}{5}-\frac{2A_{4i}\dot{h}}{5}-\frac{6A_{3i}Ui}{7}+\frac{6A_{5i}Ui}{11}$$- \frac{96A_{2i\lrcorner}\dot{4}_{2i}}{175}-\frac{12A_{4i\lrcorner}\dot{4}_{2i}}{77}-\frac{81A_{3i}\dot{A}_{3i}}{385}-\frac{81A_{5i}\dot{A}_{3i}}{1001}$$- \frac{108A_{2i\lrcorner}\dot{4}_{4i}}{385}-\frac{648A_{4i}\dot{A}_{4i}}{5005}-\frac{186A_{3i}\dot{A}_{5i}}{1001}-\frac{6A_{5i}\dot{A}_{5i}}{65}$$- \frac{138A_{2i}^{2}\dot{a}}{175\hslash}-\frac{216A_{4i}A_{2i}\dot{\hslash}}{385R_{\dot{4}}}-\frac{108A_{3i}^{2}\dot{\hslash}}{385h}-\frac{318A_{5i}A_{3i}\dot{a}}{1001\hslash}-\frac{162A_{4i}^{2}\dot{\hslash}}{1001\hslash}$$- \frac{6A_{5i}^{2}\dot{h}}{55h}-\frac{1524A_{3i}A_{2i}Ui}{2695h}+\frac{642A_{5i}A_{2i}Ui}{605h}-\frac{348A_{4i}A_{3i}Ui}{637\kappa}+\underline{480A_{5i}A_{4i}Ui}$

1573

.

$R_{\eta}$$+ \frac{(e_{ij}\cdot e_{x})}{D^{2}}[-\frac{1524A_{3i}A_{2i}R_{j}^{2}\dot{R}_{j}}{2695h}+\frac{642A_{5i}A_{2i}R_{j}^{2}R_{j}}{605h}-\frac{348A_{4i}A_{3i}R_{j}^{2}R_{j}}{637R_{i}}$(17)$+ \frac{480A_{5i}A_{4i}R_{j}^{2}\dot{R}_{j}}{1573h}-\frac{6A_{3i}R_{j}^{2}R_{j}}{7}+\frac{6A_{5i}R_{j}^{2}\dot{R}_{j}}{11}]$$a_{5i}=- \frac{\hslash\dot{A}_{5i}}{6}-\frac{A_{5i}\dot{k}}{3}-\frac{5A_{4i}Ui}{6}$$- \frac{10A_{3i}\dot{A}_{2i}}{21}-\frac{5A_{5i}\dot{A}_{2i}}{39}-\frac{65A_{2i}\dot{A}_{3i}}{126}-\frac{95A_{4i}\dot{A}_{3i}}{546}$$- \frac{58A_{3i}\dot{A}_{4i}}{273}-\frac{4A_{5i\lrcorner}\dot{4}_{4i}}{39}-\frac{5A_{2i}\dot{A}_{5i}}{18}-\frac{16A_{4i}\dot{A}_{5i}}{117}$$- \frac{9A_{2i}^{2}Ui}{7\hslash}-\frac{121A_{4i}A_{2i}Ui}{273\hslash}-\frac{120A_{3i}^{2}Ui}{637h}$$- \frac{95A_{3i}A_{2i}\dot{\hslash}}{63\hslash}-\frac{5A_{5i}A_{2i}\dot{\hslash}}{9h}-\frac{151A_{4i}A_{3i}\dot{h}}{273R}-\frac{38A_{5i}A_{4i}\dot{h}}{117h}$$+ \frac{558A_{5i}A_{3i}Ui}{1001h}-\frac{185A_{4i}^{2}Ui}{819R_{\dot{4}}}+\frac{900A_{5i}^{2}Ui}{2431h}$$+ \frac{(e_{ij}\cdot e_{x})}{D^{2}}[-\frac{9A_{2i}^{2}R_{j}^{2}R_{j}}{7h}-\frac{121A_{4i}A_{2i}R_{j}^{2}R_{j}}{273h}-\frac{120A_{3i}^{2}R_{j}^{2}R_{j}}{637R_{i}}$$+ \frac{558A_{5i}A_{3i}R_{j}^{2}\dot{R}_{j}}{1001R_{i}}-\frac{185A_{4i}^{2}R_{j}^{2}R_{j}}{819h}+\frac{900A_{5i}^{2}R_{j}^{2}R_{j}}{2431h}-\frac{5A_{4i}R_{j}^{2}R_{j}}{6}]$

(18)

は
ただし，
軸に沿った気泡の並進速度を表し， は 軸に沿った単位ベクトルである．
気泡の運動を記述する方程式は，各モードの振動の振幅を一般化座標とした Lagrange 方程式$U_{i}(=\dot{X}_{i})$$X$$e_{x}$(19)$\frac{d}{dt}(\frac{\partial \mathcal{L}}{\partial\dot{q}})=\frac{\partial \mathcal{L}}{\partial q}$

によって与えられる．Lagrangian$\mathcal{L}$は系の運動エネルギー$X$$\mathcal{K}$とポテンシャルエネルギー $V$ により以下のように定義される．
$\mathcal{L}(q,\dot{q})=\mathcal{K}(q,\dot{q})-\mathcal{V}(q)$

(20)

166
ここで， は一般化座標であり，本研究では島 $(t),$$q$$X_{i}(t),$ $A_{ni}(t)$に対応する．

気泡周囲の液体は非圧縮であるので，二つの気泡を含む系の全運動エネルギー$\mathcal{K}$は，気泡表面亀についての速度ポテンシャル (10) の積分により次のように求められる．
$\mathcal{K}=\frac{\rho}{2}\int_{V_{liq}}|\nabla\phi|^{2}dV_{liq}$$=- \frac{\rho}{2}\sum_{i=1}^{2}\int_{S_{i}}\phi(r_{si},\theta_{i})[\dot{r}_{si}+U_{i}P_{1}(\cos\theta_{i})]dS_{i}$(21)$=- \pi\rho\sum_{i=1}^{2}\int_{0}^{\pi}\phi(r_{si}, \theta_{i})[\dot{r}_{si}+U_{i}P_{1}(\cos\theta_{i})]r_{si}^{2}\sin\theta_{i}d\theta_{i}$ポテンシャルエネルギー$V_{\sigma i}$の和で与えられる．

$\mathcal{V}_{gi}$$\mathcal{V}$は，気泡内部の気体の圧縮性による$\mathcal{V}_{gi}$と気泡表面の界面張力による

は

(22)

$v_{\mathfrak{X}}=V_{i}(P_{0}+ \frac{P_{i}}{\kappa-1})=\frac{4}{3}\pi R_{ei}(P_{0}+\frac{P_{i}}{\kappa-1})$と表される．ここに は無限遠方での圧力， は気泡 の表面における圧力， は気泡内部の気体が次ののポリトロープ指数である． は気泡の体積であり，これを用いて非球形気泡の有効半径$P_{i}$$P_{0}$$i$$\kappa$$R_{ei}$$V_{i}$

ように定義される．
$R_{ei}=( \frac{3V_{i}}{4\pi})^{1/3}$$R_{\dot{\eta}}=R_{ei}- \frac{1}{R_{ei}}\sum_{n=2}^{\infty}\frac{A_{ni}^{2}}{2n+1}+O(\epsilon^{3})$

界面張力によるポテンシャルエネルギー$\mathcal{V}_{\sigma i}$

(23)
(24)

は有効半径 (23) を用いて以下で与えられる．

$\mathcal{V}_{\sigma i}=\sigma S_{i}=2\pi\sigma\int_{0^{r_{si}}}^{\pi}\sqrt{1+\frac{1}{r_{si}^{2}}(\frac{\partial r_{si}}{\partial\theta_{i}})^{2}}\sin\theta_{i}d\theta_{i}$$=4 \pi\sigma[R_{ei}^{2}+\frac{1}{2}\sum_{n=2}^{\infty}\frac{(n-1)(n+2)}{2n+1}A_{ni}^{2}]+O(\epsilon^{3})$(25)こうして得られた運動エネルギー (21) とポテンシャルエネルギー (22), (25) を用いて系の Lagrange

方程式 (19) を構成することで気泡の各モードごとの気泡の振幅が決定される．

3

2 つの非球形気泡のふるまい

ここでは気泡の振動モード $n=5$ までを考慮し，過去に得られた $n=3$ に対する結果 [6] と比較を行うことで高次の変形モードが気泡のふるまいに対して及ぼす影響について考察する．図 2 に初期
半径琉 $=5\mu m$ , 気泡間距離 $D/R_{4}=10$ の 2 気泡の自由振動における各振動モードの振幅の時間発展を示す $((a):n=3, (b):n=5)$ . すべての振幅は気泡の初期半径埼で規格化されており，横軸は

167

Normalized time

Normalized time

図 2: 気泡の形状振動の振幅の時間発展 $((a): n=3, (b): n=5)$.

Miraert 周波数で規格化されている．図 2(a),(b) それぞれについて，各モードの振幅を見てみるとどちらの場合も気泡の形状振動 $(n\geq 2)$ は時間的に発展するが，その振幅は体積振動 $(n=0(R$ モード $))$ に比較して十分に小さいことが確認できる．また，図 2(b) において高次 $(n=4,5)$ の変形モードの振幅は低次 $(n=2,3)$ のモードに比較してより小さな値にとどまっており，形状振動をともなう気泡の振動特性は定性的に 4 重極 $(n=2)$ と 8 重極 $(n=3)$ 振動によって決定されると言える．しかしながら 8 重極振動の振幅をそれぞれのケースで比較すると，$n=5$ まで展開した場合の振幅 (図中の太い破線) は $n=3$ までの展開による結果の約半分程度になっていることがわかる．これは，$n=3$よりも高次の項を考慮したために 8 重極振動のエネルギーが高次の振動モードへと伝達された結果と考えることができる．4 重極モード (図中の太い実線) がどちらのケースでもほぼ同じふるまいを示していることも形状振動モードの階層的な関係を示唆していると言える．
図 2 における時刻 $t=5.5$ (リバウンド直後) の気泡の形状を図 3 に示す．ここで，図 3(a),(b) は図$2(a),(b)$ にそれぞれ対応している．どちらの場合にも気泡崩壊にともなう高速噴流形成の要因となる
気泡背後のくぼみが観察される．気泡の位置 (並進速度に対応) や気泡前面の形状などにわずかな相違が見られるが，二つの気泡の定性的なふるまいについては展開の次数に関する依存性は少なく，気泡崩壊初期における噴流の方向やおおよその速度，気泡の位置などの予測に関しては $n=3$ までの展
開でも十分な予測が可能であると考えられる．

4

3 次元問題への拡張

ここまでの解析では 2 つの気泡からなる系を考えていたため，2 次元軸対称問題としての取扱いが可能であった．しかしながら，一般に 3 つ以上の気泡からなる系はすべての気泡がひとつの直線上にあるような例外を除いて本質的に 3 次元の問題である．したがって多数の気泡を含む気泡群のふるまいを記述するためには，これまでの理論を 3 次元へと拡張する必要がある．
図 4 に 3 次元空間における任意の気泡 に関する座標系を示す． は原点から気泡 の中心までの距離， と腕は気泡 の中心に対する原点からの方位角を表す．また，気泡 と との中心間距離を
$i$

$i$

$r_{i}$$\theta_{i}$

$i$

$i$

$j$

168

図 3: $t=5.5$ における 2 気泡の形状 $((a): n=3, (b): n=5)$$L_{ij}$とする．このとき気泡

$i$

の境界面

$r_{si}$は，球面調和関数巧 m$(\theta i, \varphi_{i})$

.

を用いて

(26)

$r_{si}=R_{\dot{\eta}}+ \sum_{l=2}^{\infty}\sum_{m=-l}^{l}A_{lmi}(t)Y_{lm}(\theta_{i}, \varphi_{i})$と表される．ただし，$A_{lmi}(t)$ンドル陪関数$P_{l}^{|m|}$は球面調和関数巧 m に対応した振動モードの振幅であり，$Y_{lm}$

はルジャ

を用いて

(27)

$Y_{lm}(\theta_{i}, \varphi_{i})=P_{l}^{|m|}(\cos\theta_{i})e^{im\varphi}$で定義される．2 気泡の場合と同様に速度ポテンシヤル をによって展開し，式 (6) と同様の境を一般化座標とするラグラ界条件を課すことで速度ポテンシャルの近似解を得る．これにより$\phi_{i}$$Y_{lm}$$A_{lmi}$に対応する極ンジュ方程式が構成される．2 気泡の場合とは異なり，3 次元多体問題の場合にはめて多くの非線形相互作用項が方程式中に現れる．実際，8 重極振動までを考慮した解析では，2 気泡の場合には 8 本の常微分方程式についてそれぞれが数十項程度の規模となるが，3 気泡の問題では数百項以上の非線形項を持つ 48 本の方程式と解かなくてはならない．したがって本稿では，3 次元問題の詳細については触れずに，単純な 3 気泡の運動の解析結果を示すにとどめる．
$A_{lmi}$図 5 は，$l\geq 3$までの展開によって得られた 3 つの非球形気泡の自由振動の様子である．3 つの気泡の初期半径馬は $7\mu m$ , 初期位置 $(xi/l$砺 $, y_{i}/RO, z_{i}/R_{0})$ はそれぞれ $(0,0,0),$ $(8,0,0),$ $(0,8,0)$ とした．
時刻 $t=0$ に振幅$R\sqrt{}$琉 $=1.4$ で体積振動を開始した．2 気泡の場合と同様に，時間の経過とともに高次の変形モードが励起されることが観察された．また，気泡の崩壊にともなう噴流の方向は気泡間相互作用の強さ，すなわち気泡間距離，振動の振幅，気泡の大きさなどに依存することが分かった．

5

まとめ

本研究では，相互作用する 2 つの非球形気泡の運動についての 2 次元軸対称問題に関して多重極展開の次数と気泡のふるまいの依存性を調べた．その結果，展開の打ち切りのために高次の振動モードの振幅については定量的な差が見られたが，気泡の定性的なふるまいに関しては 8 重極 $(n=3)$ まで

169

$x$

図 4: 気泡$i$の座標系

図 5: 相互作用する 3 つの非球形気泡の展開によって十分な近似が成り立つことが示された．このことは，8 重極振動が気泡崩壊にともなう高速噴流の形成に対して本質的な役割を担っていることを意味する．
また，3 次元問題については，本手法によって得られる力学方程式が極めて多数の項からなる大規

模な方程式系となってしまうため，得られた方程式自体を用いた数理物理的な議論を行うことは難し
い．したがって今後はより詳細な考察を行うことで気泡の運動に対して本質的な影響をもつ項を抽出し，それらの情報をもとに 3 次元の気泡群の運動を理解することが必要となる．

謝辞

本研究の一部は，財団法人マツダ財団ならびに科学研究費補助金 (23760142) の助成によるもので
ある．

参考文献
[1] Yu. A. Pishchalnikov, Oleg A. Sapozhnikov, Michael R. Bailey, James A. McAteer, James
C. Williams, Jr., Andrew P. Evan, Robin O. Cleveland, and Lawrence A. Crum. Interactionsof cavitation bubbles observed by high-speed imaging in shock wave lithotripsy. In A. A.
Atchley, V. W. Sparrow, and R. M. Keolian, editors, Innovations in Nonlinear Acoustics: 17th
International Symposium on Nonlinear Acoustics, pages 299-302, New York, 2006. American
Institute of Physics.

[2] Werner Lauterbom and Claus Dieter Ohl. The peculiar dynamuics of cavitation bubbles. Appl.
Sci. Res., 58:63-76, 1998.

170[3] Y. Tomita and A. Shima. Mechanisms of impulsive pressure generation and damage pit formation by bubble collapse. J. Fluid Mech., 169:535-564, 1986.
[4] J. R. Blake, B. B. Taib, and G. Doherty. ‘] ransient cavities near boundaries part 1. rigid
boundary. J. Fluid Mech., 170:479-497, 1986.
$]$[5] J. R. Blake, B. B. Taib, and G. Doherty. TYansient cavities near boumdaries part 2. fxee surface.
J. Fluid Mech., 181:197-212, 1987.
[6] E. Kurihara, T. A. Hay, Y. A. Ilinskii, E. A. Zabolotskaya, and M. F. Hamilton. Model for the
dynamics of two interacting axisymmetric spherical bubbles umdergoing small shape oscillation.
J. Acoust. Soc. Am., 130:3357-3369, 2011.