数理解析研究所講究録 1430 巻 2005 年 209-217

209

2 流体モデルに基づく気泡流中の分散波動
LINEAR DISPERSIVE WAVES IN BUBBLY FLOWS
北大工 矢野 猛, 江頭 竜, 藤川重雄

Takeru Yano, Ryu Egashira and Shigeo Fujikawa
Hokkaido University, Sapporo 060-8628, Japan
Email: yano\copyright mech-me.eng.hokudai.ac.jp

ABSTRACT
One-dimensional linear dispersive wave motion inwater flows containing a number of small spherical air
bubbles is analytically studied on the basis of a set of averaged equations recently derived by the present authors.
The set of equations consists of the conservation laws foreach phase and the equation of motion of bubble wall.
In addition to the volume-averaged pressures in the liquid and gas phases, the surface-averaged liquid pressure at
the gas-liquid interface is incorporated as the third pressure. The compressibility ofwater is taken into account aswell as that of gas in bubbles, and the virtual mass force
is included, although the Reynolds stress, viscosity, heat
conductivity, and the phase change across the interface are
disregarded. The results of analysis are summarized asfollows: (i) the waves in uniform bubbly flows are composed of three types of propagation modes, namely, thefast mode, slow mode, and convection mode (void wave).
This means that the uniform bubbly flow is linearly stable under one-dimensional perturbations. (ii) In the flows
where the flow velocities in gas and liquid phases are different and the other quantities are uniform, the existence
and stability ofthe waves ofthe three types depends on themodel of the virtual mass force. For a model used in thepresent analysis, the three modes stably exist for all real
wave numbers, and hence the bubbly flow is stable. (iii) In
the limit of infinitesimal void fraction, the explicit representation of the elementary solution is obtained. (iv) The
instability does not appear in the range where the present
averaged equations are applicable, ifthe virtual mass force
is suitably chosen.
1. INTRODUCTION
The characteristics of wave motions in bubbly liquids are considerably different ffom those in a pure liquid, as shown by van Wijngaarden (1968), Caflisch etat. (1985), Pauchon and Baneljee (1986), Pauchon andBanerjee (1988), Commander and Prosperetti (1989), Lahey (1991), Park et al. (1998), Gavrilyuk and Saurel(2002), and many other researchers. The appearance ofdispersion is common there. The considerable decreaseof the speed of sound of long wavelength is well know
as one of the most remarkable features of the waves in the
bubbly liquids. The propagation ofthe so-called void wave
and the onset of its instability have also been main topics.
The dynamics of bubble oscillation, however, has not always been included into the analysis. In this paper, weshall analytically study the linear wave motions in a bubbly flow, where a number of small spherical gas bubbles
are contained and the bubble oscillations can be excitedby the wave motions.
$\mathrm{n}$Recently, the present authors have analyzed linear

wave motions in a bubbly liquid uniform and at rest on
the basis of a set of averaged equations derived by our-selves (Egashira et al. 2004), The set ofequations consistsof mass and momentum conservation laws for each phaseand the Keller equation of motion of bubble wall (Keller
and Kolodner 1956). According to Jones and Prosperetti(1985), the surface-averaged liquid pressure at the bubble wall is introduced as a dependent variable as well as
the volume-averaged liquid and gas pressures. The compressibility of fluid in the liquid phase is taken into account while the Reynolds stress, viscosity, heat conductivity, and phase change across the bubble wall are neglected. On the basis of these equations, we have studied
the one-dimensional linear wave motions for the case thatthe bubbly liquid is uniform and at rest. We have demonstrated that a fast mode, which is induced by the compressibility of fluid in liquid phase, appears and sometimes behaves like a precursor of the classical pressure wave (slowmode). In the following, we shall extend the previous
study to the linear wave motions in steady bubbly flows.

210
equation (Keller and Kolodner 1956),

2. GOVERNING EQUATIONS
We shall analyze the bubbly water flows containing a
number of small spherical air bubbles on the basis of theaveraged equations. The system of goveming equations iscomposed of the mass and momentum conservation laws
for gas and liquid phases, the Keller equabon for the oscillation of spherical bubble, the equations of state for gas
and liquil and so on (Egashira et al. 2004). For the onedimensional flows the mass and momentum conservation
laws can be written as follows (see, for example, Drew andPassman 1999):$\frac{\partial}{\partial t}(\alpha p_{G})+\frac{\partial}{\partial x}(\alpha\rho_{G}u_{G})=0$,$\frac{\partial}{\partial t}[(1-\alpha)\rho_{L}]+\frac{\partial}{\partial x}[(1-\alpha)\rho_{L}u_{L}]=0$

(1)

,

(2)

,

(3)

$\frac{\partial}{\partial t}(\alpha\rho_{G}u_{G})+\frac{\partial}{\partial x}(\alpha\rho_{G}u_{G}^{2})+\alpha\frac{\partial p_{G}}{\partial x}=F$$(1- \frac{1}{c_{L0}}\frac{D_{G}R}{Dt})R\frac{D_{G}^{2}R}{Dt^{2}}+\frac{3}{2}(1-\frac{1}{3c_{L0}}\frac{D_{G}R}{Dt})(\frac{D_{G}R}{Dt})^{2}$$=(1+ \frac{1}{c_{L0}}\frac{D_{G}R}{Dt})\frac{P}{\rho_{L0}}+\frac{R}{\rho_{L0}c_{L0}}\frac{D_{G}}{Dt}(p+P)$,(7)

are, respectively,a dwater),
unperturbed
anindensity
and
sound
of
the speed
the Tait equation of state for water,($R$ is the radius of bubble,$\mathrm{n}$$c_{L0}$$\rho_{L0}$$(p_{L}+B)\rho_{L0^{n}}=(p_{L0}+B)\rho_{L}^{n}$,(8)

is the water pressure in the unperturbed state, $B=$
, and $n=7.15$ ), the isothermal equation of statefor gas,(9),
are, respectively, the density and pressure in(andside the bubble in the unperturbed state), the conservatton
law of mass for the gas inside each bubble,{$p_{L0}$$304.9\mathrm{M}\mathrm{P}\mathrm{a}$$p_{G}\rho_{G0}=p_{G0}\rho_{G}$$p_{G0}$$\rho_{G\mathrm{O}}$$\frac{\partial}{\partial t}[(1-\alpha)\rho_{L}u_{L}]+\frac{\partial}{\partial x}[(1-\alpha)\rho_{L}u_{L}^{2}]$$\rho_{G}R^{3}=\rho_{G0}R_{0}^{3}$$+(1- \alpha)\frac{\partial p_{L}}{\partial x}+P\frac{\partial\alpha}{\partial x}=-F$

,

(4)

$G$

$R_{3}$

$F$

$p$

$u$

$\rho$(10)is the radius of the bubble in the unperturbed state),and the pressure ballance at the interface,$\langle$

where cy is the volume fraction of gas phase $(0<\alpha<1)$ ,is the density, is the velocity, is the pressure, represents the interfacial momentum transport, and the suband denote the volume-averaged quantities
scripts
in gas and liquid phases, respectively. In addition to theand , the pressure avervolume-averaged pressuresaged on the gas-liquid interface, $P$ , is introduced in thesystem of equations (see, Jones and Prosperetti 1985).
In the previous paper (Egashira et at. 2004), as theinterfacial momentum transport $F$, we have employed thevirtual mass force,,$p_{G}-(p_{L}+P)= \frac{2\sigma}{R}$

,

(11)

$L$

$p_{G}$$p_{L}$$F=- \beta[\frac{D_{G}}{Dt}(\alpha\rho_{L}u_{G})-\frac{D_{L}}{Dt}(\alpha\rho_{L}u_{L})]$,(5)

with $\beta=1/2$ used by Kameda and Matsumoto (1996).
How ever, the operator $u_{G}(D_{G}/Dt)-u_{L}(D_{L}/Dt)$ in Eq. (5)is not invariant under the Galilean hansformation, $\tau=t$ ,$\xi=x-Vt$ , and $v=u-V$ . We therefore introduce a different model,,(6)

which is invariant under the Galilean transformation. The, and are not given, and wevalues of coefficients
simply postulate Eq. (6) here. The importance of such invariance has been discussed in detail by Drew et al. (1979).
The previous analysis (Egashira et al. 2004) has been confined to the linear waves in quiescent bubbly Iiquids, andhence the results are not affected by the unfavorable character of Eq. (5).
The conservation laws ofmass and momentum for the, are augmented by the Keller
bubbly liquid, Eqs.
$(1\succ\prec 4)$$x=Lx^{*}$,$t= \frac{L}{\mathrm{q}_{0}}t^{*}$$\frac{u_{G}}{c_{L0}}=u_{G0}^{*}+\epsilon u_{G}’$$\alpha=\mathrm{q}^{*}(1+\epsilon\alpha’)$

,

,

,

$\frac{R}{R_{0}}=1+\epsilon R’$$\frac{u_{L}}{\mathrm{c}_{L0}}=u_{L0}^{*}+\epsilon u_{L}^{l}$,$\frac{p_{L}}{\rho_{L0}c_{L0^{2}}}=p_{L0}^{*}(1+\epsilon p_{L}’),$$\frac{p_{G}}{\rho_{L0}c_{L0^{2}}}=p_{G0}^{*}(1+\epsilon p_{G}’)$$\frac{\rho_{G}}{\rho_{L0}}=\rho_{G0}^{*}(1+\epsilon\rho_{G}^{J})$

,

,

,

$\frac{P}{\rho_{L0}\mathrm{c}_{L0^{2}}}=\epsilon P’$$\frac{\rho_{L}}{\rho_{L0}}=1+\epsilon\rho_{L}’$$\ovalbox{\tt\small REJECT}$

(12)

,

,

is a characteristic length of an initialis a measure of the amplitude of the perturbation, the prime indicates the normalized perturbation of
the order of unity. The quantities with are the nondimensionalized ones and the subscript 0 denotes the value
in the unperturbed steady flow;
are nondimensional constants, which are not .
andessarily of the order of unity. In the following analysis,
only the nondimensional variables are used, and therefore
be omitted.
the prime and the asterisk
As a result, we have the following set of the linearized
equations: the mass conservation law in gas phase,$L$$\mathrm{p}\mathrm{e}\Phi \mathrm{r}\mathrm{b}\mathrm{a}-$$\epsilon(\ll 1)$$- \beta_{2\beta_{L}}(u_{G}-u_{L})\frac{D_{G}\alpha}{Dt}-\beta_{3}\alpha(u_{G}-u_{L})\frac{D_{G}\rho_{L}}{Dt}$$\mathcal{B}_{2}$earized as follovs

where
tion,

$F=- \beta_{1}\alpha\rho_{L}(\frac{D_{G}u_{G}}{Dt}-\frac{\mathrm{R}^{u_{L}}}{Dt})$$\beta_{1},$where a is the surface tension.
AI1 the variables are nondimensionalized and lin-$\beta_{3}$$*$$u_{G0}^{*},$$u_{L0}^{*},$$\alpha_{0}^{*},$$p_{L0}^{*},$$p_{G0}^{*},$$\mathrm{n}\mathrm{e}\mathrm{c}$$\rho_{G0}^{*}$$\mathrm{w}\mathrm{i}\mathrm{U}$$\frac{\partial}{\partial t}(\alpha-3R)+\frac{\partial}{\partial_{X}}[u_{G0}(\alpha-3R)+u_{G}]=0$

,

(13)

211
represents a typical inter-bubble distance. The void ffaction in the unperturbed state % therefore is sufficiently
small compared with unity, i.e.,the mass conservation law in liquid phase,$\frac{\partial}{\partial t}[\alpha \mathit{0}\alpha+(\%-1)p_{L\mathrm{f}}p_{L}]$(14)$\alpha_{0}=\frac{4}{3}\pi L^{3}\delta^{3}\%\ll 1$$+ \frac{\partial}{\partial x}[\alpha \mathit{0}u_{L0}\alpha+(\circ_{\mathrm{b}}-1)(p_{L0}u_{Li^{p_{L}+u_{L})]}}=0,$

the momentum conservation law in gas phase,

.

(22}

has toThe wavelength of the perturbation
be sufficiently large compared with , and at least it is as. This leads to another restriction,
large as$\lambda=2\pi L/l\sigma$$R_{0}$$\frac{\partial}{\partial t}[(u_{G0}-u_{L0})(\beta_{2}\alpha+\beta_{3}p_{L0}p_{L})+(\beta_{1}+\rho_{G0})u_{G}-\beta_{1}u_{L}]$$\sqrt[3]{V}$$+ \frac{\partial}{\partial x}[u_{G0}(u_{G0}-u_{L0})(\beta_{2}\alpha+\beta_{3}p_{L\mathrm{f}}p_{L})-3p_{G0}R$$\frac{R_{0}}{\lambda}=\delta\kappa\ll 2\pi$$+(\beta_{1}+\rho_{G0})u_{G0}u_{G}-\beta_{1}u_{L0}u_{L}]=0$

,

.

(23)

(15)

the momentum conservation law in liquid phase,

3.1 Waves in uniform bubbly flow
Firstly, we shall consider the case that the gas and liquid phases move with the same speed in the same directionin the initial unperturbed state, namely $u_{G0}=u_{L0}=U$ .
Then, the dispersion relation can be written as$\frac{\partial}{\partial t}[\mathrm{q})(u_{G0}-u_{L0})(\beta_{2}\alpha+\beta_{3}p_{L0}p_{L})$$+\beta_{1}\alpha_{0}(u_{G}-u_{L})+(\mathrm{q}_{1}-1)u_{L}]$$+ \frac{\partial}{\partial x}[o_{b}u_{G0}(u_{G0}-u_{L0})(\beta_{2}\alpha+\beta_{3}p_{L0}p_{L})+(\%-1)p_{L0}p_{L}$$f(\kappa,\omega)=\{2\beta_{1}[\kappa^{2}-(\omega-\kappa U)^{2}][\delta^{2}(\omega-\kappa U)^{2}-\omega_{B}^{2}]$,$+\beta_{1}\alpha_{0}(u_{G0}u_{G}-u_{L0}u_{L})+(\mathrm{q}_{)}-1)u_{L0}u_{L}]=0$(16)$+2\mathrm{q}_{1}\beta_{1}(\omega-\kappa U)^{2}[3-2\omega_{B}^{2}+2\delta^{2}(\omega-\kappa U)^{2}]$and the linearized Keller equation,$-2\mathrm{q}^{2}\beta_{1}(\omega-\kappa U)^{2}[3-\omega_{B}^{2}+\delta^{2}(\omega-\kappa U)^{2}]$

,

(17)

$\frac{\partial R}{\partial t}+u_{G\mathit{0}}\frac{\partial R}{\partial_{X}}=S$,$+3\mathit{0}_{\mathrm{b}^{2}}\kappa^{2}(\omega_{B}^{2}-p_{L0})\}(\omega-\kappa U)^{2}=0$(18)$\frac{\partial S}{\partial t}+u_{G0}\frac{\partial S}{\partial x}=\frac{(2\hat{\sigma}-3p_{G0})R-p_{L0}p_{L}}{\delta^{2}}$, (24)

is the nondimensionalized eigenwhere
frequency of a bubble in the initial unperturbed state. The
linear dispersion relation is independent of the coefficients
here andin Eq. (6) Note that we neglectandis nottheresultsmall
andishereafter sinceaffected substantially.
In the case of $\beta 1=1/2$ , for an arbitrary real in all realistic conditions in bubbly flows of water and air bubbles,as the solutions of Eq. (24),
we have five real$\omega_{B}=\sqrt{3p_{L0}+4\hat{\sigma}}$

where$\delta=\frac{R}{L}$$-\mathrm{q}_{)}\kappa^{2}[(3+\beta_{1})\omega_{B}^{2}-3(1+\beta_{1})p_{L0}+2\beta_{1}\delta^{2}(\omega-\kappa U)^{2}]$,,$p_{G0}=p_{L0}+2\hat{\sigma}$$\hat{\sigma}=\frac{\sigma}{\rho_{L0}c_{L0}^{2}\mathrm{f}\mathrm{l}\}}$

(19)

.

$\beta_{2}$In the equation of motion of bubble wall (18), we neglectsome terms responsible for the wave attenuation due toacoustic radiation during the bubble oscillation.

$\beta_{3}$$\rho_{G0}$$(\cong 10^{-3})$$\rho_{G0}$$\kappa$$\omega’ \mathrm{s}$3. LINEAR DISPERSIVE WAVES
The linear analysis is started by substituting a plane)] etc., into thewave solution, i.e., $cx=a\exp$ [ (it(13}\prec18), where
equationssystem of the linearized
is a wavenumber ( is a normalized
is a constant,is an angular &equency ( iswavenumber), anda normalized angular frequency). The dispersion relation$f(\kappa,\omega)=0$ is obtained as the eigenequation (or secular
is a polyequation) of the resulting systenb where
representaexplicit
.
Its
in
and
order
of
sixth
nomialtion is too long and we may proceed the discussion without
showing it.
Here, we shall remark that the description on the basisof the averaged equations ( $1\succ(11)$ , and hence Eqs. (13)(18), may be subject to the restrictions,$\mathrm{i}$$\omega_{1}\pm=\kappa U\pm\sqrt\frac{\overline A-\sqrt{A^{2}-4B}}{2\delta^{2}}$

(25)

,

(26)

$-\kappa x$$a$$\kappa/L$,$\omega_{2\pm}=\kappa U\pm\sqrt{\frac{A+\sqrt{A^{2}-4B}}{2\delta^{2}}}$$\kappa$$\omega_{3}=\kappa U$

where$f(\kappa,(\mathit{0})$$/\sigma$$R_{0}\ll\sqrt[3]{V}<<L$

and

(27)

$\omega$$\omega \mathrm{q}_{0}/L$A$= \frac{\delta^{2}\kappa^{2}}{1-\mathcal{H}}+^{2}\alpha J_{B}+\frac{3\%}{1-\%}$

,

(28)

$\omega$$R_{0}\ll n_{\mathrm{D}}^{-1/3}\ll\sqrt[3]{V}$

,

(20)

or

$B=[u;_{B}^{2}+ \frac{3\alpha \mathit{0}(3-2\%)p_{G0}}{1-\%}]\frac{\delta^{2}\kappa^{2}}{1-\%}$.

When $U=0$, Eqs. (25) and (26) reduce to the results inthe bubbly liquid uniform and at rest studied in the Preandvious paper (Egashira et al. 2004), whereare, respectively, referred to as the slow mode and the fast
may be called a convecmode. The third mode w3
are real for an arbifive
the
fact
that
tion mode. The
hary real means that the one-dimensional steady flow of$u_{G0}=u_{L0}=U$ described by the present model is linearly$\omega_{1}\pm$$=\kappa U$$\delta^{3}\ll\frac{V}{L^{3}}<<1$

and$L^{3}\delta^{3}\%\ll 1\ll Vn_{0}$

,

(21)

$\omega’ \mathrm{s}$$\kappa$

where $V$ is the volume for averaging, is the number density of bubble in the unperturbed steady state, and$n_{0}$$\eta_{\}}^{-1/3}$(29)

stable under the one-dimensional perturbations.

$\omega 2\pm$212
then we have$f( \kappa,\omega)\cong\frac{1}{2}(\omega-\kappa U)^{2}\frac{\partial^{2}f}{\partial\omega^{2}}|_{0}+\frac{1}{2}(u_{G0}-U)^{2}\frac{\partial^{2}f}{\partial u_{G0}^{2}}|_{0}$. (32)$+( \omega-\kappa U)(u_{G0}-U)\frac{\partial^{2}f}{\partial\omega\partial u_{G0}}|_{0}=0$

Clearly, Eq. (32), the quadratic equation for , has tworeal roots, if and only if an inequality$\omega$$( \frac{\partial^{2}f}{\partial\omega\partial u_{G0}}|_{0})^{2}\geqq(\frac{\partial^{2}f}{\partial\omega^{2}}|_{0})(\frac{\partial^{2}f}{\partial u_{G0}^{2}}|_{0})$

,

(33)

holds. The partial derivatives in Eq. (33) can be given ffom
) as
the complete form of$f_{\backslash }^{/}\kappa,\omega$normalized wave number$\frac{\partial^{2}f}{\partial\omega^{2}}|_{0}=\frac{2\kappa^{2}p_{L0}}{\delta^{2}}[\beta_{1}(\omega_{B}^{2}+2\mathrm{q}_{\mathrm{J}}\hat{\sigma})+(\%-\alpha 0^{2})(_{\mathrm{t}}4^{2}+2\hat{\sigma})]$,$\kappa$

Fig 1. The dispersion relation of fast and slow modes inthe bubbly quiescent liquid: $u_{G0}=u_{L0}=U=0,$ $\delta=0.01$ ,
.
and

(34)
$\frac{\partial^{2}f}{\partial u_{G0}^{2}}|_{0}=\frac{2\kappa^{4}p_{L0}}{\delta^{2}}(1-\%)(\beta_{1}-\beta_{2})(\omega_{B}^{2}+2\mathrm{q}_{\mathrm{J}}\hat{\sigma}),$(35)$R0=1\mathrm{m}\mathrm{m}$$\frac{\partial^{2}f}{\partial\omega\partial u_{G0}}|_{0}=\frac{\kappa^{3}p_{L0}}{\delta^{2}}(\mathrm{q}_{\}}-1\rangle(2\beta_{1}-\beta_{2})(\omega_{B}^{2}+2\%\hat{\sigma})$Figure 1 shows the dispersion relation of the fast and. (2004).
slow modes for $U=0$ studied by Egashira$U\neq
0$
by
given
be
can
for
The corresponding curves.
adding$e\ell d$$\kappa U$

.

(36)

Note that the above three derivatives are independent ofand$\beta_{3}$$U$.

, the inequality (33) holds
In the simplest casefor arbitrary real and $U$ , and the solution of (32) for$\beta_{1}=\beta_{2}$$\beta_{1}=$$\kappa$$\beta_{2}$is given as

3.2 Waves in non-uniform flow$\kappa U$

In the case that $u_{G0}\neq u_{L0}$, the angular frequency

cannot be expressed explicitly in terms of simple functionsof . However, the stability of convection mode when theis sufficiently smalJ cananddifference between
readily be examined in the following way.
and $u_{G0}=u_{L0}=$
aroundExpanding
gives
,
and
to
with respect$\kappa$$u_{G0}$$\omega,$$u_{G0}$$\kappa[U-\frac{c}{D}(u_{G0}-u_{L0}\rangle]$

(37)

,

where

$u_{L0}$$C=\beta_{1}$$\omega=\kappa U$$f(\kappa,\omega)$$U$,$\omega=\{$$\omega$(1-%)$(\omega_{B}^{2}+2\%\hat{\sigma})$

,

(38)

$D=\beta_{1}\omega_{B}^{2}+\mathrm{q}_{1}[\omega_{B}^{2}+2(1+\beta_{1})\hat{\sigma}]-3\%^{2}(p_{L0}+2\hat{\sigma})(3^{\cdot}9)$$u_{L0}$The same as Eq. (37) with Eqs. (38) and (39) can beobtained even if we put $u_{G0}=U$ in Eq. (30). As a result, it is found that the convection mode when the differis sufficiently small is the staand
ence between
ble propagation mode for arbitrary real , if the inequality(33) holds. In other words, the bubbly flows with different
gas- and liquid-phase velocities are stable if the differencein flow velocities are sufficiently small. The third coeffidoes not affect the above
cients of virtual rnass force
results.
Simple analytical approach as shown above cannot beapplied to the casc that the difference in the flow velocities
is not small. It is, however, easy to examine the dispersion
relation numerically. The dispersion relations for the six
propagation modes are numerically evaluated and plottedin Fig. 2. As shown in the figure, all the six modes seemto be stable for moderately small in the entire range of, even if the difference in the flow velocities is not verysmall. Further numerical examination suggests that the$\omega$$f( \kappa,\omega)=\frac{1}{2}(\omega-\kappa U)^{2}\frac{\partial^{2}f}{\partial\omega^{2}}|_{0}+\frac{1}{2}(u_{G0}-U)^{2}\frac{\partial^{2}f}{\partial u_{c\mathrm{o}^{2}}}|_{0}$$u_{G0}$$+ \frac{1}{2}(u_{L\mathrm{f}\mathrm{J}}-U)^{2}\frac{\partial^{2}f}{\partial u_{L\tilde{0}}^{\gamma}}|_{0}+(\omega-\kappa U)\langle u_{G0}-U)\frac{\partial^{2}f}{\partial\omega\partial u_{G0}}|_{0}$$u_{L0}$$\kappa$$+( \omega-\kappa U)(u_{L\mathit{0}}-U)\frac{\partial^{2}f}{\partial\omega\partial u_{L0}}|_{0}$$\beta 3$$+(u_{G0}-U)(u_{L0}-U) \frac{\partial^{2}f}{\partial u_{G0}\partial u_{L0}}|_{0}+O(3)$

,

(30)

where the partial derivative with the subscript 0 denotes
and $u_{G0}=u_{L0}=U$ , and wethe value evaluated athave used the fact that$\omega=\kappa U$.

$f( \kappa,\kappa U)|_{u_{G0}=u_{L0}=U}=\frac{\partial f}{\partial\omega}|_{0}=\frac{\partial f}{\partial u_{G0}}|_{0}=\frac{\partial f}{\partial u_{L0}}|_{0}=0$$\alpha_{0}$(31)Without loss of generality, we can put$u_{G\mathrm{O}}\neq u_{L0}=U$, arid$\kappa$

213
nent, and cause some interesting nonlinear phenomenon.

The nonlinear problem will be studied in a ffiture work.

3.3 The limit of infinitesimal void fraction, the linearized mass and moIn the limit of
mentum conservation equations for liquid phase, Eqs. (14)and (16), are simplified as$\alpha_{0}arrow 0$

,

(40)

$\frac{\partial}{\partial t}(p_{L\mathrm{f}}p_{L})+\frac{\partial}{\partial x}(p_{L0}u_{L0}p_{L}+u_{L})=0$$\frac{\partial u_{L}}{\partial t}+\frac{\partial}{\partial x}(p_{L0}p_{L}+u_{L0}u_{L})=0$

,

(41)

which are the same as those for an inviscid compressible
fluid with the Tait equation of state 8). The dispersion
relation for the system (13), (40), (15), (41), (17), and (18)
can then be given as$\langle$

normalized wave number$\kappa$$(\omega-\kappa u_{G0})[\delta^{2}(\omega-\kappa u_{G0})^{2}-\omega_{B}^{2}][(\omega-\kappa u_{L0})^{2}-\kappa^{2}]$$\mathrm{x}[\beta_{2}\kappa(u_{G0}-u_{L0})+\beta_{1}(\omega-\kappa u_{G0})]=0$

.

(42)

disappears in the dispersion relation. From
Again,
Eq. (42), we immediately have$\beta_{3}$$\omega_{1}[perp]"=\kappa u_{G0}\pm\frac{\omega_{B}}{\mathit{5}}$$\omega_{2\pm}=\kappa u_{L0}\pm\kappa$

,

(43)

,

(44)

$\kappa u_{G0}$$\omega 3\pm=\{$$\kappa[u_{G0}-\frac{\beta_{2}}{\beta_{1}}(u_{G0}-u_{L0})]$

(45)

,

represents the ffee oscillation ofwhere the slow modeis theand the fast mode
bubbles convected withacoustic wave in the liquid moving with . The dispersion relation for each mode is shown in Fig. 3. As can beare almost equal
seen, the three modes in the limit
modes arefastandconvection. The
to those for.
almost dispersionless at
The expressions of the system of equations and the setof angular frequencies are now quite simple, and therefore
we can explicitly write down the linear dispersive waves
in the form suitable for the initial-value problem,$\omega 1\pm$$\omega z\pm$$u_{G0}$$u_{L0}$$\alpha_{0}arrow \mathrm{O}$normalized wave number$\alpha_{0}=10^{-5}$$\kappa$$\alpha_{0}=10^{-5}$

Fig 2. The dispersion relations of the six propagation
modes in the bubbly flows with different gas- and liquidphase velocities. The open circle indicates the maximumof the fast mode. (a): $=0.01$ , $\alpha_{0}=0.0002,$ $u_{G0}=0.4$ ,) $:\delta=0.01,$ %= 0.0002, $u_{G0}=0.5,$ $u_{L0}=$$0.2$ . The virtual mass coefficients are set as
and $R_{0}=0$ .lmm in both cases.
$\delta$$u_{L0}=0.3;\{\mathrm{b}$$\beta_{1}=\beta_{2}=$$q(x,t)= \sum_{n=1}^{3}f_{-\infty}^{\infty}a_{n}(\kappa)\hat{q}_{n+}(\kappa)\exp[\mathrm{i}(\omega_{n+}t-\kappa x)]d\kappa$$\beta_{3}=\frac{1}{2}$$+ \sum_{n=1}^{3}\int_{-\infty}^{\infty}b_{n}(\kappa)\hat{q}_{n-}(\kappa)\exp[i(\omega_{n-}t-\kappa x)]d\kappa$

angular frequency of convection mode does not becomecomplex in the wavenumber range shown in Eq. (23).
Incidentally, one can see two important features inFig. 2: one is the intersection of the slow mode and convection mode, which suggests the existence of interactionof these two modes in the nonlinear problem. The othercurve of fast mode propagatis the maximum in theing in the opposite direction to the bubbly flow, denotedby an open circle. Since the maximum means the zerogroup velocity, this may lead to a kind of localization ofwavenumber compothe wave energy of

,

(46)

where

$a_{n}(\kappa)$

and $b_{n}(\kappa)(n=1,2,3)$ are the expansion co-efficients, $q(x,t)$ is a column vector$\omega-\kappa$$\mathrm{c}\mathrm{o}\mathrm{r}\tau \mathrm{e}\mathrm{s}\mathrm{p}\mathrm{o}\mathrm{n}\mathrm{d}\mathrm{i}\mathrm{n}\mathrm{g}$$q(x,t)=\ovalbox{\tt\small REJECT}^{R(x,t)}u_{G}(x,t)\ovalbox{\tt\small REJECT} p_{L}(x,t)u_{L}(x,t)\alpha(x,t)S(x_{2}t)$

,

(47)

214
bubbly flows, is independent of the liquid-phase variables
and , whereas the fast mode (49) affects all the vanables in gas and liquid phases. The convection mode (50}consists of only and components, and hence xt may
deserve to be called the void wave.
Once an initial condition $q(x,0)$ is given, we haveffom Bq. (46)$p_{L}$$u_{L}$$\alpha$$u_{G}$$\sum_{n=1}^{3}[a_{n}(\kappa)\hat{q}_{n+}(\kappa)+b_{n}(\kappa)\hat{q}_{n-}(\kappa)]$$= \frac{1}{2\pi}\int_{-\infty}^{\infty}q(x,0)e^{i\kappa x}dx$. (S2)For the demonstration, we shall consider two simple ex-

normalized wave number
Fig 3.

$u_{G0}=0.5,$ $u_{L0}=0.2,$ $\delta=0.01,$

The dashed curves are those for$\alpha_{0}arrow 0$$\alpha_{0}=$

:

.

and
0.00001.

$\beta_{1}=\beta_{9}\sim$$\beta 1=\beta 2=\beta 3$$\kappa$

The dispersion relation in the limit ofamples: (i) $p_{L}(x,0)\neq 0$ and the other initial values are allzero; (ii) $\alpha(x,0)\neq 0$ and the other initial values are all.
zero. Furthermore, we assume
In the first example, the expansion coefficients are
given as$R\mathrm{o}=0.\mathrm{l}\mathrm{m}\mathrm{m}$’$a_{1}( \kappa)=\frac{\mathrm{i}p_{L0}p_{L}^{\mathrm{A}}(\kappa)}{\delta\omega_{B}^{2}[(\delta\kappa V+\omega_{B})^{2}-\delta^{2}\kappa^{2}]}$

,

(53)

,

(54)

$b_{1}( \kappa)=\frac{\mathrm{i}p_{L0}\hat{p}_{L}(\kappa)}{\delta\omega_{B}^{2}[\delta^{2}\kappa^{2}-(\delta\kappa V-\omega_{B})^{2}]}$and$\hat{q}_{n}\pm(\kappa)(n=1,2,3)$

are expressed as follows:,$\hat{q}_{1}\pm(\kappa)=\ovalbox{\tt\small REJECT}^{2\mathrm{i}\delta\omega_{B}(\mp\beta_{1}\omega_{B}-\beta_{2}\delta\kappa V)}3\mathrm{i}\overline{\delta}[\mp 2\beta_{1_{B}}\omega_{0}^{222}\pm\delta\kappa^{arrow Q]}\ovalbox{\tt\small REJECT} 3\mathrm{i}\delta\omega_{\epsilon_{0}}[\delta\kappa Q\pm 2\beta_{2}\omega_{B}V]2\omega_{B}^{2}(\beta_{1}\omega_{B}\pm\beta_{2}\delta\kappa V)$(48)$a_{2}( \kappa)=\frac{\mathrm{a}\hat{p}_{L}(\kappa)}{(1-V)[\omega_{B}^{2}-\delta^{2}\kappa^{2}(1-V)^{2}]}$

,

(55)

,

(56)

$b_{2}( \kappa)=\frac{\mathrm{i}\hat{p}_{L}(\kappa)}{(1+V)[\omega_{B}^{2}-\delta^{2}\kappa^{2}(1+V)^{2}]}$$a_{3}( \kappa)=\frac{2p_{L0}[3(\omega_{B}^{2}-p_{L0})+2\beta_{1}\omega_{B}^{2}V^{2}]\hat{p}_{L}(\kappa)}{\omega_{B}^{2}(V^{2}-1)}$$b_{3}(\kappa)=-4p_{L0}\hat{p}_{L}(\kappa)$

,

(57)
(58)

,

is the Fourier transform of $p_{L}(x,0)-$ The secwhere
ond example immediately gives(59)$a_{1}(\kappa)=a_{2}(\kappa)=b_{1}(\kappa)=b_{2}(\kappa)=b_{3}(\kappa)=0$ ,$\hat{p}_{L}(\kappa)$$\hat{q}2\pm(\kappa)=\ovalbox{\tt\small REJECT}^{2\mathrm{i}_{\mathrm{P}_{L0}}(V\mp 1)Y_{\mp}}-ip_{L0}(V\mp 1)[Z_{\mp}+6\beta_{2}(V\mp 1)V]\ovalbox{\tt\small REJECT} \mathrm{i}p_{L0}[Z_{\mp}+6\beta_{1}(V\mp 1)^{2}]\pm 2\mathrm{i}p_{L0\mp}T(V\mp 1)Y_{\mp}2\kappa p_{L0(\mp}2\mathrm{i}T_{\mp}(V\mp 1)Y_{\mp}V1)^{2}Y_{\mp}$, (49)$a_{3}(\kappa)=\hat{\alpha}(\kappa)$

where

,

$\hat{q}_{3+}(\kappa)=\ovalbox{\tt\small REJECT}_{0}^{0}10\ovalbox{\tt\small REJECT} 00$$\hat{q}_{3-}(\kappa)=\ovalbox{\tt\small REJECT}^{\beta_{1}}-\beta_{2}V)0000$

’

(50)

$\hat{\alpha}(\kappa)$

(60}

,

is the Fourier transform of$\alpha(x,0)$.

3.4 Numerical examples
In this subsection, we shall present typicai examples
obtained by solving the system of the linearized equations(13}\prec 18) for some initial conditions with a finite differshow the evolution ofwave proence method. Figures
files of pressure, bubble radius, and void ffiction. In thesefigures, wave profiles are plotted from $t=0$ to $t=50$, andthe bold dashed line indicates the locus of wave front.
The case that only the pressure perturbation is imposed at $t=0$ is shown in Figs, 4 and 5, The pressureis the fast mode, and the amplitudes of
profile in Fig.
the waves of the bubble radius and void ffaction are small, and asin this case. Figure 5 is the case for
are
mode
and
convection
mode
3,
the
fast
shown in Fig.
almost dispersionless, although the convection mode is not$\angle\vdash 7$$V=u_{G0}-u_{L0}$,$Q=\omega_{B}^{2}-p_{L0}(=2p_{L0}+4\hat{\sigma})$,,$T_{\mp}=\delta^{2}\kappa^{2}(V\%1)^{2}-\omega_{B}^{2}$$Z_{\mp}=-3Q+2T_{\mp}[\beta_{1}+\beta_{3}(V\mp1)V].\}$(51),$Y_{\mp}=\mp\beta_{1}+(\beta_{1}-\mathcal{B}_{2})V$$4(\mathrm{a})$may be regarded as eigenvectors asThe vectors.
sociated with the eigenvalues
From Eqs. ($48\vdash(50)$ , one can see that the slow mode(48), which corresponds to the classical pressure wave in$\hat{q}n\pm(\kappa)$$\omega_{n\pm}(\kappa\rangle$$\mathrm{q}=10^{-5}$clearly seen.

215
Figures 6 and 7 are the case that the perturbation ofvoid fraction only is applied at $t=0$ . In the both figures,a convection mode appears and its amplitude is large compared with the other modes. In particular the void wave isshown in Fig. 7.
strongly excited in the case ofsubscripts$B$

0: unperturbed steady state$\alpha_{0}=10^{-5}$

CONCLUSIONSREFERENCES

The detailed analysis of the linear dispersive wave inbubbly flows has been carried out based on the averaged
equations composed of conservation laws and equation of
motion of bubble wall. The result suggests that ali thethree types of waves may be stable in the wavenumber
range where the present set of averaged equations is applicable.
Gavrilyuk and Saurel (2002) have also obtained a similar stability result in a different method of analysis. Their
model equation contains an analogue of Rayleigh-Plesset
equation. The bubble dynamics may have some important
effect on the wave motion in the bubbly flows.

NOMENCLATURE

$c$

$L$

$n$

$p$

: speed of sound: characteristic length: number density of bubble: volume-averaged pressure: surface-averaged pressure: six dimensional column vector: radius ofbubble: oscillation velocity ofbubble walI
: time: volume-averaged velocity

$P$

$q$

$R$

$S$

$t$

$u$

$V$

$x$

: volume for averaging: space coordinate: Bubble: Gas and liquid phases$L,G$$\mathrm{m}/\mathrm{s}$$\mathrm{m}$$\mathrm{m}^{-3}$

Pa
Pa

Caflisch, R. E., Miksis, M. J., Papanicolaou, G. C.,
and Ting, L. (1985) Effective equation for wave propagation in bubbly liquids, J. Fluid Mech. 153, 259.
W. and Prosperetti, A. (1989) Linear
Commanderpressure waves in bubbly liquids: Comparison between
theory and experiments, J. Acoust. Soc. Am. 85, 732.
R., Yano, T., and Fujikawa, S. (2004) Linear
wave propagation of fast and slow modes in mixtures of. $34,317$ .
liquid and gas bubbles, Fluid $Dyn$.
L.,
Drew, D., Cheng,and Lahey, R. T., Jr. (1979) The
analysis of virtual mass effects in two-phase flow, Int. .
Multiphase Flow 5, 233.
Drew, D. A. and Passman, S. L. (1999) Theory ofMul$\mathrm{K}_{-}$$\mathrm{E}\mathrm{g}\mathrm{a}\mathrm{s}\mathrm{h}\dot{\alpha}*$${\rm Res}$$J$

ticomponent Fluids, Springer, New York.
Gavrilyuk, S. and Saurel, R. (2002) Mathematical andnumerical modeling of two-phase compressible flows withmicro-inertia, J. Comp. Phys. 175, 326.
Jones, A. V. and Prosperetti, A. (1985) On the suitability offirst order differential models for tvvo phase flow
prediction, Int. J. Multiphase Flow 11, 133,
Kameda, M. and Matsumoto, Y. (1996) Shock waves
in a liquid containing small gas bubbles, Phys. Fluids 8,322.

$\mathrm{m}$$\mathrm{m}/\mathrm{s}$$\mathrm{s}$$\mathrm{m}f\mathrm{s}$$\mathrm{m}^{3}$

Keller, J. B. and Kolodner I. I. (1956) Damping ofunderwater explosion bubble oscillations, J. Appl. Phys.
27 1152.
Lahey, R. T. Jr., (1991) Void wave propagation Phenomena in two-phase flow, AIChE J. 37, 123.
Park, J.-W., Drew, D. A., and Lahey, R. T., Jr.
(1998) The analysis of void wave propagation in adiabatic monodispersed bubbly two-phase flows using an
two-fluid model, Int. J. Multiphase
Flow 24, 1205,, S. (1986) Interphase moPauchon, C. andinthe averaged multifield
effects
interaction
mentummodel. Part : Void propagation in bubbly flows, Int. .
Mulliphase Flow 12, 559.
, S. (1988) Interphase moPauchon, C. andmentum interaction effects in the averaged multifieldmodel. Part : Kinematic waves and interfacial drag inbubbly flows, Int. J. Multiphase Flow 14, 253.
van Wijngaarden, L. (1968) On the equations of motion for mixtures of liquid and gas bubbles, J. Fluid Mech.
33, 465.
$\mathrm{e}\mathrm{n}\mathrm{s}\mathrm{e}\mathrm{m}\mathrm{b}\mathrm{l}\mathrm{e}rightarrow \mathrm{a}\mathrm{v}\mathrm{e}\mathrm{r}\mathrm{a}\mathrm{g}\mathrm{e}\mathrm{d}$$\mathrm{m}$$\mathrm{E}\mathrm{a}\mathrm{n}\mathrm{e}\iota \mathrm{j}\mathrm{e}\mathrm{e}$Greek symbols$\alpha$$\beta$$\delta$$\epsilon$$\kappa$$\rho$: void fraction (gas phase fraction): coefficient of interfacial momentum transport: ratio of to: sufficiently small parameter: wavenumber$R_{0}$$\mathrm{B}\mathrm{a}\mathrm{n}\mathrm{e}\dot{\eta}\mathrm{e}\mathrm{e}$$—-$$1\mathrm{I}$$L$: volume-averaged density

ci: angular frequency$J$$\mathrm{I}$$\mathrm{m}^{-1}$$\mathrm{k}\mathrm{g}/\mathrm{m}^{3}\mathrm{s}^{-1}$216

nondimensional space coordinate

nondimensional space coordinate

Fig 4. Wave profiles in the case of$\alpha_{0}=$0.0002, $u_{G0}=$
and $R0=$$0.5,$ $u_{L0}=0.2,$ $\delta=0.01,$ $\beta_{1}=\beta_{2}=\beta_{3}=0.5$ ,:( pressure, ) bubble radius, (c) void fraction.
The initial condition is $p_{L}=\exp(-x^{2})$ and $a=R=u_{L}=$$u_{G}=S=0$ at $t=0$ .
$0.\mathrm{l}\mathrm{m}\mathrm{m}$$\mathrm{a}\rangle$$\{\mathrm{b}$Fig 5. Wave profiles in the case of

0.00001, $u_{G0}=$
and $R_{0}=$pressure, (b) bubble radius, (c) void fraction.

$0.5,$ $u_{L0}=0.2,$ $\delta=0.01,$$0.\mathrm{l}\mathrm{m}\mathrm{n}:(\mathrm{a})$$\alpha_{0}=$$\beta_{1}=\beta_{2}=\beta_{3}=0.5,$

The initial condition is $p_{L}=\exp(-x^{2})$ and oc $=R=u_{L}=$$u_{G}=S=0$ at $t=0$ .

217

nondimensional space coordinate

nondimensional space coordinate

Wave profiles in the case of $\alpha_{0}=0.0002,$ $u_{G0}=$$0.5,$ $u_{L0}=0.2,6=0.01,$ $\beta_{1}=\beta_{2}=\beta_{3}=0.5$ , and $R0=$
pressure, { bubble radius, (c) void ffiction.
The initial condition is $\alpha=\exp(-x^{2})$ and $R=p_{L}=u_{L}=$$u_{G}=S=0$ at $t=0$ .

Fig 6.

$0.\mathrm{l}\mathrm{m}\mathrm{m}:(\mathrm{a})$$\mathrm{b})$Fig 7. Wave profiles in the case of

0.00001, $u_{G0}=$, and $R0=$.
pressure, (b) bubble radius, (c) void ffiction.
Ihe initial condition is $\alpha=\exp(-x^{2})$ and $R=p_{L}=u_{L}=$$u_{G}=S=0$ at $t=0$ .
$0.5,$ $u_{L0}=0.2,$ $\delta=0.01,$$0$$\mathrm{l}\mathrm{m}\mathrm{m}:(\mathrm{a})$$\alpha_{0}=$$\beta_{1}=\beta_{2}=\beta_{3}=0.5$