- Email: [email protected]
Perturbation method for the second-order nonlinear effect of focused acoustic field around a scatterer in an ideal fluid
Ultrasonics 54 (2014) 576–585
Contents lists available at ScienceDirect
Ultrasonics journal homepage: www.elsevier.com/locate/ultras
Perturbation method for the second-order nonlinear effect of focused acoustic field around a scatterer in an ideal fluid Gang Liu, Pahala Gedara Jayathilake, Boo Cheong Khoo ⇑ Department of Mechanical Engineering, National University of Singapore, Kent Ridge Crescent, Singapore 119260, Singapore
a r t i c l e
i n f o
Article history: Received 18 March 2013 Received in revised form 23 July 2013 Accepted 18 August 2013 Available online 11 September 2013 Keywords: Perturbation method Westervelt equation Mach number Nonlinear acoustic wave
a b s t r a c t Two nonlinear models are proposed to investigate the focused acoustic waves that the nonlinear effects will be important inside the liquid around the scatterer. Firstly, the one dimensional solutions for the widely used Westervelt equation with different coordinates are obtained based on the perturbation method with the second order nonlinear terms. Then, by introducing the small parameter (Mach number), a dimensionless formulation and asymptotic perturbation expansion via the compressible potential flow theory is applied. This model permits the decoupling between the velocity potential and enthalpy to second order, with the first potential solutions satisfying the linear wave equation (Helmholtz equation), whereas the second order solutions are associated with the linear non-homogeneous equation. Based on the model, the local nonlinear effects of focused acoustic waves on certain volume are studied in which the findings may have important implications for bubble cavitation/initiation via focused ultrasound called HIFU (High Intensity Focused Ultrasound). The calculated results show that for the domain encompassing less than ten times the radius away from the center of the scatterer, the non-linear effect exerts a significant influence on the focused high intensity acoustic wave. Moreover, at the comparatively higher frequencies, for the model of spherical wave, a lower Mach number may result in stronger nonlinear effects. Ó 2013 Elsevier B.V. All rights reserved.
Introduction Applications of high intensity sonic energy in the industry have increased during the past decades [1]. The applications of intense sound are based on nonlinear effects produced by finite amplitude pressure variations [2,3]. It is well known that linear acoustic theory is based on the assumption of small amplitude waves and linear constitutive theory of the fluid medium [4]. Although these assumptions hold for many vibro–acoustic interactions, they become inapplicable in sound fields with high amplitude pressure [5]. Unlike the linear acoustic wave equation, the nonlinear counterpart can handle waves with large finite amplitudes and allow accurate modeling of nonlinear constitutive models in the fluid. Interesting phenomena undetected via linear acoustics can be observed in these studies. For example, waveform distortion, formation of shock waves, increased absorption, nonlinear interaction (as opposed to superposition) are seen when two sound waves are mixed, amplitude dependent directivity of acoustic beams, cavitation and sonoluminescence [6].
⇑ Corresponding author. Address: Singapore-MIT Alliance, 4 Engineering Drive 3, National University of Singapore, Singapore 117576, Singapore. E-mail addresses: [email protected], [email protected] (B.C. Khoo). 0041-624X/$ - see front matter Ó 2013 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.ultras.2013.08.011
The widely used linear wave equation used to describe the small signal (infinitesimal amplitude) propagation is the Helmholtz equation [7] given as below:
r2 H
1 @2H ¼ 0: c21 @t2
ð1Þ
Here, H can be pressure, velocity or velocity potential in the field of fluid mechanics, or displacement and displacement potential for solid mechanics [8]; c1 is the small signal sound speed since there are no shear waves inside the ideal fluid. For the linear wave propagation in elastic solids, the in-plane longitudinal and transverse wave [9] as well as anti-plane transverse wave [10] can co-exist and governed by the same Helmholtz equation. If we neglect the time term in Eq. (1), the Laplace equation is also widely employed for the incompressible ideal fluid [11,12] and static linear solid mechanics [13]. To the best of authors’ knowledge, there are various models to simulate the nonlinear characteristic of the acoustic wave propagating through the fluid. For instance, the one-dimensional Burgers equation has been found to be an excellent approximation of the conservation equations for plane progressive waves of finite amplitude in a thermoviscous fluid [14]. An effective model that combined effects of diffraction, absorption, and nonlinearity in directional sound beams (i.e. radiated from sources with
577
G. Liu et al. / Ultrasonics 54 (2014) 576–585
assumed a solution for the acoustic pressure consisting of two terms in the form of p = p1 + p2, where p1 represents the first order approximation and p2 the second order correction. Here p2 is much smaller than p1 (p2 p1). Neglecting the terms of third and higher orders, we can obtain the approximate equation as:
r2 ðp1 þ p2 Þ
1 @ 2 ðp1 þ p2 Þ b @ 2 p21 þ ¼ 0: 2 2 c1 q1 c41 @t2 @t
ð4Þ
According to Eq. (4), the second order solution for the acoustic wave can be written as:
p ¼ P1 expðixtÞ þ P2 expð2ixtÞ:
Fig. 1. Schematic description of the model for multiple acoustic waves focused around a scatterer inside an ideal fluid.
dimensions that are large compared with the characteristic wavelength) are taken into account by the Khokhlov–Zabolotskaya– Kuznestov (KZK) parabolic wave equation [15,16]. Several other models have been developed in response to specific needs. For example, the Westervelt equation is an almost incidental product of Westervelt’s discovery of the parametric arrays [17–19]. In the past, much computational work have been done on nonlinear effects in beams based on the KZK equations [20–23]. The one dimensional case of Westervelt equation is studied extensively by finite element method [24]. However, when acoustic fields of more complex conditions are considered, advanced numerical calculation related to finite-differential becomes necessary [25,26]. In this work, our interest is to develop the analytical solution of the multiple harmonic acoustic waves focused on the volume near a scatterer, where the second order nonlinear effects dominate (see Fig. 1). Our study has important implications for further work on bubble/nucleation cavitation by HIFU (High Intensity Focused Ultrasound) and others. (Our approach based on asymptotic perturbation expansion of small parameter and non-dimensionalisation carried out w.r.t the dimension of the scatterer, however, is not applicable in the absence of the scatterer.) Second order nonlinear solution for Westervelt equation On the problems of multiple acoustic waves focusing on intensive volume, one effective model governing the dynamics is the Westervelt equation as given below:
r2 p
1 @2p b @ 2 p2 d @3p þ þ ¼ 0: 2 2 c21 @t c41 @t 3 q1 c41 @t
ð2Þ
Here, the third term is the nonlinear term and the forth term stands for the thermoviscous effect. q1 is the static medium density, c1 is the wave speed of the acoustic wave in the undisturbed liquid, b = (c + 1)/2 is the coefficient of nonlinearity, and c is defined as the ratio of specific heats of the gas. Here, d = m(# + (c 1)/Pr) is the diffusivity associated with sound absorption, where m = l/q1 is the kinematic viscosity, # = 4/3 + lB/l is the viscosity number, lB is the bulk viscosity coefficient, and Pr is Prandtl number. If we consider the case of ideal fluids, the forth term can be set equal to zero, and the reduced Westervelt equation is given as:
r2 p
1 @2p b @ 2 p2 þ ¼ 0: c21 @t2 q1 c41 @t 2
ð3Þ
To solve the second-order nonlinear equation, the successive approximation method can be applied. In this method, it is
The proposed analytical solutions of the one dimensional Eq. (4) for the planar, cylindrical, and spherical waves are obtained as shown below. It is pointed out, however, that the coefficients on the right hand side of the initial governing equation for the second-order nonlinear term shown as Eq. (10) in [24] are deemed incorrect. In particular, Eq. (10) in Ref. [24] is shown with the term 2 2k b=ðq1 c21 Þ, which we reckoned should be reflected as 2 4k b=ðq1 c21 Þ after the separation of the variables. The details for our derivations were not shown here, and the results with different expressions and coefficients are given below as:
pp ¼ P 0 eikx eixt þ pc
¼ P0
qffiffiffiffiffiffi 2
iðkrp=4Þ
pkr e ð2Þ
bP20 ð1 þ 2ikxÞe2ikx e2ixt ; 2q1 c21 4ibP20 2 1 c1
eixt þ pq
P0 H0 ðkrÞ eixt þ
2ibkrP20 q1 c21
ð2Þ
ð5Þ
eip=2 e2ikr e2ixt 2
ð6Þ
ðH0 ðkrÞÞ e2ixt ;
! 2 1 P 0 ikr ixt ibkP 0 ln r X ð1Þn ðn 2Þ! e2ikr e2ixt : ps ¼ e e þ þ n1 r q1 c21 r ð4ikÞ rn n¼2 ð7Þ It is clear that the first term in Eqs. (5)–(7) are the general solutions for the one dimensional linear acoustic wave with time harmonic condition, and the second terms are the approximate solutions for the second-order nonlinearity which satisfy the Sommerfield radiation condition at the far field. Here, P0 stands for harmonic excitation of the wave amplitude, and k ¼ x=c1 ¼ 2p=k, where x is the circular frequency of the incident wave, and k stands for the wave length. It should be noted that the calculated results plotted for the second-order solutions for the cylindrical and spherical waves differ from Ref. [24] on two counts. As discussed above, it is due to the different governing equations for the second-order nonlinear term and the boundary conditions prescribed in [24] is the Neumann type as it enforces the normal derivative of pressure at the boundary of the rigid scatterer. However, for our simulation, we imposed the Dirichlet boundary conditions. See Section 4 below for further discussion on the observed difference. Perturbative method with a small parameter for the nonlinear acoustic wave Mathematical formulation of the nonlinear acoustic wave From the analytical solutions in Section 2 for the one dimensional Westervelt equation, it is clear that one can only get the pressure for the fluid domain without the explicit results for the velocity and velocity potential. Below, in order to give a more clearer description about the nonlinear effects of the (nonlinear) acoustic wave scattering by inclusion, we propose a mathematical formulation and a different model for the nonlinear acoustic wave that goes beyond the traditional models based on the Burgers
578
G. Liu et al. / Ultrasonics 54 (2014) 576–585
equation, KZK equation, and Westervelt equation. The compressible potential flow theory is invoked together with the dimensionless formulation and asymptotic perturbation expansion referring to the small parameter Mach number. It will be shown that our approach allows the decoupling of the potential and enthalpy terms to the second order, and which accord much flexibility for the calculation of the other physical quantities inside the fluid domain. For our model, we shall take the liquid as inviscid and compressible. As such, the liquid flow is governed by the equation of mass conservation:
~ @u 1 ~ ru ~ ¼ rp: þu @t q
ð9Þ
In many fluid dynamic problems, significant liquid compressibility due to high speed motion occurs when thermal effects in the liquid are unimportant. We therefore assume that thermal effects in the liquid are insignificant. With this assumption, the liquid state is defined by a single thermodynamic variable. The sound speed c and the enthalpy h of the liquid can be defined as follows [27]:
h¼
Z
p
dp
q
p1
ð10Þ
;
where the reference pressure p1 is the pressure in the undisturbed liquid (hydrostatic pressure). By assuming that the flow is irrotational, we may introduce a ~ ¼ rw. By applying velocity potential w such that u _ q_ rp=rq, and u ~ ¼ rw, where the dot ‘‘’’ stands c2 ¼ dp=dq p= for the differential w.r.t time, we can re-write Eq. (8) as:
r2 w þ
1 @h þ rw rh ¼ 0: 2 c @t
ð11Þ
Similarly for Eq. (9), the integration leads to the unsteady Bernoulli equation:
@w 1 þ jrwj2 ¼ h h1 : @t 2
ð12Þ
Here, h1 may be set equal to zero since the enthalpy is referenced to the undisturbed fluid at infinity. To find the expressions for the speed of sound c and enthalpy h, we use the Tait model equation of state for water to relate the pressure and density as follows [28,29]:
p¼
p1 þ B
qn1
qn B:
ð13Þ
The values B = 3049.13 bars and n = 7.15 give an excellent fit based on experiments in water for up to 105 bars [30]. Here, p1 and q1 are the pressure and density in the undisturbed liquid (hydrostatic model). From Eq. (13) together with c2 = dp/dq, we get:
c2 ¼
nðp þ BÞ
q
¼
n
q1
p p1
q1
2 1 p p1 þ ...; 2c21 q1
ð16Þ
1 1 pp ¼ ðn 1Þ 4 1 þ . . . c2 c21 c1 q1
ð17Þ
ð8Þ
and the momentum conservation equation:
dp ; dq
h¼
Non-dimensional formulation of the governing equations
@q ~ Þ ¼ 0; þ r ðqu @t
c2 ¼
Next, we expand the enthalpy h and the sound speed c around the equilibrium pressure p1 using a Taylor series expansion as follows [28]:
ðp1 þ BÞ1=n ðp þ BÞðn1Þ=n
ð14Þ
!
ðn1Þ=n pþB 1 ; p1 þ B
r¼
Rs
e
r ;
c ¼ c1 c ;
t¼
Rs t ; U
w ¼ Rs Uw ;
h ¼ U 2 h ;
p ¼ p1 þ q1 U 2 p ;
ð18Þ ð19Þ
where the Mach number e = U/c1 (is also the small parameter to be used for the perturbation analysis). It is evident that e – 0, or otherwise the dimensionless formulation in Eq. (18) will have a singularity. From the above, we have: 2 @2w 2 U @ w ¼ e ; @r 2 Rs @r 2
@h U 2 @h ¼e ; @r Rs @r
@w @w ¼ eU ; @r @r
@w @w ¼ U2 ; @t @t
@h U 3 @h ¼ @t Rs @t
ð20Þ
ð21Þ
Next, by substituting Eqs. (20), (21), and (19) into (11) and (12), respectively, we obtain the governing equations in terms of dimensionless variables:
@ 2 w 1 @h 2 @w @h ¼ 0; þ þ e c2 @t @r2 @r @r
ð22Þ
2 @w 1 2 @w þ e þ h ¼ 0: @t 2 @r
ð23Þ
We stress that the dimensionless variable r⁄ also stands for the space vector dependent on the dimensional coordinate system which is planar, cylindrical or spherical. Next, from Eq. (15), the dimensionless speed of sound of the liquid is:
c2 ¼ 1 þ e2 ðn 1Þh ;
ð24Þ
and the dimensionless enthalpy of the liquid in Eq. (15) is given by:
h ¼
1
e2
1 n1
! ðn1Þ=n pþB 1 : p1 þ B
ð25Þ
Substituting p = p1 + q1U2p⁄ and h = U2h⁄ into Eqs. (16) and (17), and referring to the Taylor series expansion, we can get.
h ¼ p
1 2 2 e p þ oðe2 Þ; 2
1 ¼ 1 ðn 1Þe2 p þ oðe2 Þ: c2
and
c2 c21 c2 h¼ ¼ 1 n1 n1
We would like to introduce two typical scales Rs and U, where Rs denotes the radius of the scatterer and U (p1/q1)1/2 stands for the liquid velocity near the scatterer. The dimensionless quantities indicated by asterisk are given as below:
ð15Þ
ð26Þ
ð27Þ
By substituting Eqs. (26) and (27) into Eq. (22), we can obtain:
@ w @h @h 2 @w @h þ oðe2 Þ ¼ 0: þ þ e ðn 1Þh @r2 @t @r @r @t 2
where c21 ¼ nðpq1 þBÞ is the wave speed of the acoustic wave in the 1 undisturbed liquid.
ð28Þ
G. Liu et al. / Ultrasonics 54 (2014) 576–585
Next, the perturbation expansions for the potential w⁄ and enthalpy h⁄ are defined as follows:
w ðr ; t Þ ¼ w0 ðr ; t Þ þ e2 w1 ðr ; t Þ þ oðe4 Þ;
ð29Þ
h ðr ; t Þ ¼ h0 ðr ; t Þ þ e2 h1 ðr ; t Þ þ oðe4 Þ;
ð30Þ
where these are then introduced into Eqs. (23) and (28). We can get the separate governing equations for the first and second order w.r.t. e2 as below:
@w0 þ h0 ¼ 0; @t 2 @w1 1 @w0 þ þ h1 þ oðe4 Þ ¼ 0; 2 @r @t
ð31Þ ð32Þ
2
@ w0 @h0 þ ¼ 0; @r 2 @t
ð33Þ
@ 2 w1 @h1 @w0 @h0 @h0 þ þ ðn 1Þh0 þ oðe4 Þ ¼ 0: @r 2 @t @r @r @t
ð34Þ
@ 2 w0 ¼ 0; @t2
r2 w1
@ 2 w1 1 @ @ jr w0 j2 r w0 ðr w0 Þ 2 @t @t @t 2
ðn 1Þ
@w0 @ 2 w0 þ oðe4 Þ ¼ 0: @t @t 2
ð35Þ
ð36Þ
where the Hamilton operator r⁄ is defined in terms of r⁄. As elucidated before, the parameter n is assumed to reasonably match the experimental pressure-density relation for water up to 10 [5] bars at n = 7.15 [30]. Analytical solution for the one-dimensionless equation Analytical solution for the plane wave The analytical solution for the dimensional second-order nonlinear plane wave is given below as (details are provided in the Appendix A):
w0 ¼ A cosðkrÞ eixt ;
ð37Þ
w1 ¼ A cosð2krÞ e2ixt ixA2 nþ1 e2ixt ; ðn þ 1Þkr sinð2krÞ þ cosð2krÞ þ n 3 4 8U 2 ð38Þ where A, k, and x are the amplitude, wave number, and frequency, respectively. Analytical solution for the cylindrical wave The analytical solution for the dimensional second-order nonlinear cylindrical wave is given below as: ð2Þ
w0 ¼ A H0 ðkrÞ eixt ;
Analytical solution for the spherical wave The analytical solution for the dimensional second-order nonlinear spherical wave is given below as:
w0 ¼
w1 ¼
ARs cosðkrÞ eixt ; er
ð41Þ
ARs ixA2 R2 cosð2krÞ e2ixt 2 s er 8U e2 r2 ! krðn þ 1ÞðCið4krÞ þ 2Cið2krÞ þ ln Rers Þsinð2krÞ e2ixt : krðn þ 1ÞðSið4krÞ þ 2Sið2krÞÞcosð2krÞ 8cos2 ðkrÞ ð42Þ
By substituting Eqs. (37) and (38), Eqs. (39) and (40), Eqs. (41) and (42) into w = w0 + e2w1 + O(e4), respectively, we obtain the corresponding asymptotic expression for the velocity potential up to and including O(e2). Results and discussions
By combining Eqs. (31) and (33), Eqs. (32) and (34), we can obtain the fundamental equation for the second order nonlinear acoustic wave (dimensionless formulation) as:
r2 w0
579
ð39Þ
ipA2 x3 e2 ð2Þ w1 ¼ AH0 ð2krÞe2ixt þ 2U 4 0 1 R 2 2 ð2Þ ð2Þ J 0 ð2krÞ ððn1ÞðH0 ðkrÞÞ 2ðH1 ðkrÞÞ ÞY 0 ð2krÞrdr @ A e2ixt : R 2 2 ð2Þ ð2Þ Y 0 ð2krÞ ððn1ÞðH0 ðkrÞÞ 2ðH1 ðkrÞÞ ÞJ0 ð2krÞrdr ð40Þ
The obtained analytical solutions are analyzed with regard to the nonlinear effects around the volume inside the liquid where the multiple incident waves are focused. It is important to emphasize that the initial interactions among the multiple linear harmonic sources are not taken into account since each source is not having a sufficiently large intensity to induce nonlinear deformation. If one assume the scatterer to be the equilibrium bubble, the Rayleigh–Plesset equation [31] or modified Rayleigh model [32] as well as the relationship between the pressure in the fluid and the gas pressure inside the bubble [33] (where the nonlinear effects have not been taken into account previously), can be employed to describe the bubble cavitation/deformation. Alternatively, we may use the general dimensional equation dr/dt = rw together with Eqs. (41) and (42) to construct and express the nonlinear variation of the bubble radius as well. In our paper, we shall assume the radius of the scatterer Rs = 1. Since we know that the nonlinear acoustic wave effect will only dominate the fluid region and its vicinity outside the scatterer, the distance defined as x or r should larger than the radius of the scatterer (Rs = 1); in the subsequent analysis the effective distance Re will range from 1 to 100. By using the analytical solutions for the one dimensional Westervelt equation, the pressure distribution around the scatterer for the fundamental and the second order harmonic of planar, cylindrical and spherical coordinates in water are calculated. Here, in our calculation, the wave speed c1 is assumed to be 1500 m/s, water density q1 is 1000 kg/m3, nonlinear coefficient b is 3.5, initial pressure P0 is 4 105 Pa, and frequency is f = 15, 000 Hz [24]. (For example, one can interpret the presence of four incident waves at pressure intensity of Pi = 105 Pa each. As such, there are the multiple incident waves focusing with linear superposition giving rise to P0 = 4 105 Pa at the vicinity of the scatterer.) The calculated results for the pressure ratio between the second order harmonic to fundamental (first order term) are displayed in Fig. 2. Fig. 2a shows that the ratio between the second harmonic to the first order harmonic for plane wave increases linearly with increasing distance away from the focused point. For the second harmonics of spherical and cylindrical waves as shown in Fig. 2b and c, respectively, the increase with distance from the focused point is no longer linear. The mentioned increase for the spherical wave depicts a higher rate of change with distance compared to the cylindrical wave; both of them depict a lower magnitude and rate of increase compared to the planar wave. As mentioned in Section 2, the differences observed for Fig. 2b and c applicable to the cylindrical and spherical waves, respectively, vis-à-vis the corresponding results in [24] are attributed to the different governing
580
G. Liu et al. / Ultrasonics 54 (2014) 576–585
Fig. 3 presents a more explicit comparison between the fundamental and the second order harmonic pressure amplitudes in planar, cylindrical and spherical coordinates. According to the definition of the pressure level (dB) set at 20 log10(p/preference), we can plot the pressure level distribution w.r.t. the distance away from the focused point (the reference pressure is assumed to be 106 Pa, and the range is from 1.0 to 100 which is different from Ref. [24] found in the range from 0 to 1.0). The calculated pressure amplitude variation shown in Fig. 3, however, can be compared to the results obtained by [24]. From Fig. 3, it is observed that the second order harmonic distribution for the planar wave increases with distance and becomes larger than the first harmonic at about 25 units. The second order harmonic distribution of cylindrical wave seems to exhibit minor variation. This is in line with the approximate expanded solution for the second order of the cylindrical 4ibP 2
wave taken to be pq c02 eip=2 e2ikr e2ixt as shown in Eq. (6). It is 1 1 well known that the absolute value of this complex number has no variation with respect to the increase of the distance r. Lastly, it is observed that the second order for the harmonic spherical wave initially increases slightly and then decreases with the increment of the distance away from the focused point. The summation of the fundamental and second order terms are shown in Fig. 4. The numerical results of the linearized governing equations are separately calculated by Finite Difference method to check on the accuracy and consistency of the analytical results. It is evident that our analytical results concur well with the numerical results. It is also observed that the amplitude for the cylindrical and spherical waves decrease monotonically. In particular, the spherical wave decreases at a faster rate than the cylindrical wave since the pressure for the spherical wave varies with the approximate factor r1 (see Eq. (7)), whereas the cylindrical waves with 1 the factor of r2 (see Eq. (6)). On the contrary, the pressure amplitude increases with distance away from the focused point for the planar wave. This can be attributed to the absolute value of the second order pressure planar wave increasing with respect to the distance away from the focused point (see Eq. (5)). It is clear that within the distance of less than 10, the decrease of the aggregate pressure amplitudes of cylindrical and spherical wave appear to be very fast and the tendency to moderate beyond this point for the region between 10 and 100. Being so, one can suggest that the domain encompassing less than ten times of the radius of the scatterer, the non-linear effect exerts a significant influence on the focused high intensity acoustic wave. Next, we shall discuss the nonlinear effect for the focused acoustic wave with the small parameter Mach number by
Fig. 2. The ratio of the pressure second harmonic to the fundamental term vs. the variation of the distance away from the focused point.
equation for the second order nonlinear term and different boundary conditions imposed. It is clear that although the general/broad variation of increasing of magnitude of second harmonic to fundamental ratio with distance from the focal point is similar, the observed trend and rate of increase are different (for our case the boundary conditions of Dirichlet type are imposed). For the second harmonic to fundamental ratio versus distance from focused point observed between our Fig. 2a and Fig. 1 of [24] applicable for the plane wave, the difference in magnitude can also be attributed to the different governing equation for the second order nonlinear term employed.
Fig. 3. The comparison between the analytical results of planar, cylindrical and spherical wave including the fundamental and the second harmonic.
G. Liu et al. / Ultrasonics 54 (2014) 576–585
Fig. 4. The variation of pressure amplitude vs. the distance away from the focused point for the planar, cylindrical and spherical wave.
compressible potential flow theory. We will assume the amplitude of velocity potential A = 1, and radius of the scatterer Rs = 1. As interpreted earlier, x⁄r⁄ = kr is defined as the dimensionless wave number and the velocity potential is only dependent on the variable U in the frequency domain and the small parameter e = U/c1 (c1 is a constant at 1500 m/s). For the numerical calculations, we shall take the dimensionless wave number as kRs = 1.0, 3.0 and
581
5.0 and the small parameter Mach number e as 0.1, 0.3 and 0.5, so that U = ec1 becomes 150, 450 and 750, respectively. Fig. 5 shows the variation of the second order term of plane wave w.r.t. the increase of the distances from the focused point for different parameters. It is evident that with the increase of the wave number k (shown as the horizontal rows), we can observe that the general tendency of the amplitude for the second order term to increase and the variation takes on shorter periodic time. For the fixed low wave number of k = 1.0 (first column in Fig. 5), the increase of the Mach number (e = 0.1, 0.3 and 0.5) can bring about comparatively large amplitude for the second order term. However, the maximum amplitude for each second term appears to be constant without any variation w.r.t. the increase of distance from the focused point. This phenomena can be interpreted as that the low frequency stands for quasi-static model, and thus the variation of the distances have lesser influence on the second-order nonlinear term. At the high frequency of k = 5.0 (third column in Fig. 5), the maximum amplitudes tend to increase w.r.t. the position away from the focused point. Fig. 6 illustrates the variation of the nonlinear effect for the cylindrical wave at different wave number, Mach number and distance away from the focused point. The general tendency is that with the increase of the distance away from the focused point, this results in ever smaller amplitude for the second order term. For the fixed wave number k, the increase of the Mach number will bring about a larger amplitude for the second order term. With the fixed higher Mach number (e = 0.5), the increment of the wave number can also weaken the amplitude for the second order term. Fig. 7
Fig. 5. The variation of the second order term of plane wave vs. the variation of wave number k, Mach number e and the distance away from the focused point.
582
G. Liu et al. / Ultrasonics 54 (2014) 576–585
Fig. 6. The variation of the second order term for cylindrical wave vs. the different wave number, Mach number and the distance away from the focused point. Fig. 7. The variation of the second order term for spherical wave vs. the different wave number, Mach number and the distance away from the focused point.
shows the variation of the nonlinear effect for the spherical wave for the different wave number, Mach number and distance away from the focused point. Similar to that observed for the cylindrical wave, the general tendency is that the increase of the distance away from the focused point will result in smaller amplitude for the second order term. At the fixed low frequency of k = 1.0, the increase of the Mach number (e = 0.1, 0.3 and 0.5) will bring about a larger amplitude for the second order term. Conversely, at the higher frequencies of k = 3.0 and 5.0, a lower Mach number results
in stronger nonlinear effects and this trend is very different from the cylindrical wave. Finally, Fig. 8 gives the total velocity potential (the summation of the first order term and second order term) at k = 2.0 and e = 0.3. We can observe the variation of the velocity potential for both the cylindrical wave and spherical wave appears to be consistent with the conclusion from the Westervelt model for the pressure distribution. The general tendency of the velocity potential will decrease
583
G. Liu et al. / Ultrasonics 54 (2014) 576–585
@ 2 w0 @ 2 w0 ¼ 0; @r 2 @t2
ðA:1Þ
2 @ 2 w1 @ 2 w1 1 @ @w0 @w @ 2 w0 @w ¼ þ 0 þ ðn 1Þ 0 2 2 2 @t @r @r @r @r @t @t @t
@ 2 w0 : @t2
ðA:2Þ
Here, we would like to point out that r⁄ can be directly replaced by the Cartesian coordinate x⁄ for the plane wave. We propose the one dimensionless solution for the linear acoustic wave in Eq. (A.1) as below:
w0 ¼ A cosðk r =eÞ eix t A cosðx r Þ eix t ;
ðA:3Þ
where A⁄, k⁄, x⁄ are the dimensionless amplitude, wave number, frequency, respectively, and defined as:
Fig. 8. The velocity potential distribution near the scatterer (the summation of the first order term and the second order term) at k = 2.0 and e = 0.3.
w.r.t. the increase of the distance away from the focused point while the amplitude of the velocity potential for the planar wave seems to be periodic constant. For the spherical wave, since our analytical solution is expressed by cosine/sine function (or cosine/sine integral function) rather than the exponential function, a characteristic feature of this function is periodic rather than just a monotonic curve as from the Westervelt equation. Conclusions In this study, we adopted two nonlinear models to investigate the multiple incident acoustic waves focused on certain domain where the nonlinear effect is not negligible in the vicinity of the scatterer. The general solutions for the one dimensional Westervelt equation with different coordinate systems (plannar, cylindrical and spherical) were analytically obtained based on the perturbation method with keeping only the second order nonlinear terms. By introducing the small parameter (Mach number), we applied the compressible potential flow theory and proposed a dimensionless formulation and asymptotic perturbation expansion for the velocity potential and enthalpy which is different from the existing fractional nonlinear acoustic models (e.g. the Burgers equation, KZK equation and Westervelt equation). Our analytical solutions and numerical calculations have shown the general tendency of the velocity potential and pressure to decrease w.r.t. the increase of the distance away from the focused point. At least, within the region which is about 10 times the radius of the scatterer, the nonlinear effect exerts a significant influence on the distribution of the pressure and velocity potential. It is also interesting to note that at high frequencies, lower Mach numbers appear to bring out even stronger nonlinear effects for the spherical wave. Our approach with a small parameter for the cylindrical and spherical waves could serve as an effective analytical model to simulate the focused nonlinear acoustic near the scatterer in an ideal fluid and be applied to study bubble cavitation dynamic associated with HIFU in our future work. Appendix A Analytical. solution for the plane wave For the one dimensional (1D) plane wave, the governing equations are (the higher order small terms will be omitted subsequently):
A ¼
A ; Rs U
k ¼ Rs k;
x ¼
Rs x; U
r ¼
e Rs
r;
t ¼
U t: Rs
ðA:4Þ
It is evident that Eq. (A.3) is a solution of Eq. (A.1). By substituting Eq. (A.3) into Eq. (A.2) and assuming that w1 ¼ /1 ðr Þ e2ix t , we obtain the simplified second-order non-homogeneous differential equation as:
@ 2 /1 2 2 þ 4x2 /1 ¼ 2iA x3 sin ðx r Þ iðn 1ÞA2 x3 @r2 cos2 ðx r Þ:
ðA:5Þ
The general analytical solution for this second-order nonhomogeneous Eq. (A.5) is:
w1 ¼ C 1 cosð2x r Þ e2ix t þ C 2 sinð2x r Þ e2ix t ix A2 ððn þ 1Þx r sinð2x r Þ 8 nþ1 cosð2x r Þ þ n 3Þ e2ix t : þ 4
ðA:6Þ
where C1 and C2 are unknown coefficients to be determined by the imposed boundary conditions. For simplicity in the discussion below, the coefficients for the second-order nonlinear terms can be reasonably assumed as C1 = A⁄, C2 = 0 without loss of generality. Analytical. solution for the cylindrical wave For the one dimensional (1D) cylindrical wave, the governing equations are:
@ 2 w0 1 @w0 @ 2 w0 þ ¼ 0; r @r @r 2 @t 2
ðA:7Þ
2 @ 2 w1 1 @w1 @ 2 w1 1 @ @w0 @w0 @ 2 w0 þ ¼ þ r @r 2 @t @r @r 2 @r @r @t @t 2 þ ðn 1Þ
@w0 @ 2 w0 : @t @t 2
ðA:8Þ
The proposed solution for the linear acoustic wave (second kind of Hankel function with zero order) is assumed as: ð2Þ
w0 ¼ A H0 ðx r Þ eix t :
ðA:9Þ
By substituting Eq. (A.9) into Eq. (A.8) and assuming that w1 ¼ /1 ðr Þ e2ix t , we obtain the simplified second-order nonhomogeneous differential equation as: 2 @ 2 /1 1 @/1 2 ð2Þ þ þ 4x2 /1 ¼ 2iA x3 ðH1 ðx r ÞÞ r @r @r2 2
iðn 1ÞA2 x3 ðH0 ðx r ÞÞ ; ð2Þ
ðA:10Þ
584
G. Liu et al. / Ultrasonics 54 (2014) 576–585
and the analytical solution is (the softwave named MAPLE will be used here): 2ix t
w1 ¼ C 1 J0 ð2x r Þe 0 B @
J0 ð2x r Þ
R
þC 2 Y 0 ð2x r Þe 2
ð2Þ ððn1ÞðH0 ð
Y 0 ð2x r Þ
x r ÞÞ
R
2ix t
x r ÞÞ
ð2Þ
w0 ¼ A H0 ðkrÞ eixt ;
ðA:19Þ
3 2 A
ipx þ 2
1
2
ð2Þ 2ðH1 ð 2
ð2Þ ððn1ÞðH0 ð
Therefore, the dimensional 1D cylindrical wave is given as:
x r ÞÞ ÞY 0 ð2x r Þr dr C Ae2ix t :
ð2Þ 2ðH1 ð
2
x r ÞÞ ÞJ0 ð2x r Þr dr
ðA:11Þ Here, J0() is the Bessel function of order zero, Y0() is the second ð2Þ kind of Bessel function of order zero, and H0 ðÞ ¼ J 0 ðÞ iY 0 ðÞ is second kind Hankel function with order zero. For ease of discussion and without loss of generality, one can reasonably assume that C1 = A⁄ and C2 = iA⁄.
ipA2 x3 e2 ð2Þ w1 ¼ A H0 ð2krÞ e2ixt þ 2U 4 0 1 R 2 2 ð2Þ ð2Þ J 0 ð2krÞ ððn 1ÞðH0 ðkrÞÞ 2ðH1 ðkrÞÞ Þ Y 0 ð2krÞ rdr A e2ixt ; @ R 2 2 ð2Þ ð2Þ Y 0 ð2krÞ ððn 1ÞðH0 ðkrÞÞ 2ðH1 ðkrÞÞ Þ J 0 ð2krÞ rdr ðA:20Þ
and the dimensional 1D spherical wave as:
w0 ¼
w1 ¼
Analytical. solution for the spherical wave For the one dimensional (1D) spherical wave, the governing equations are:
ARs cosðkrÞ eixt ; er
ðA:21Þ
ARs ixA2 R2 cosð2krÞ e2ixt 2 s er 0 8U e2 r2 1 krðn þ 1Þ Cið4krÞ þ 2Cið2krÞ þ ln Rers sinð2krÞ A e2ixt : @ krðn þ 1ÞðSið4krÞ þ 2Sið2krÞÞ cosð2krÞ 8cos2 ðkrÞ ðA:22Þ
@ 2 w0 2 @w0 @ 2 w0 þ ¼0 r @r @r 2 @t 2
ðA:12Þ
2 @ 2 w1 2 @w1 @ 2 w1 1 @ @w0 @w @ 2 w0 þ ¼ þ 0 2 r @r 2 @t @r 2 @r @r @r @t @t @w @ 2 w0 þ ðn 1Þ 0 : @t @t 2
ðA:13Þ
A cosðx r Þ eix t : r
ðA:14Þ
By substituting Eq. (A.14) into Eq. (A.13) and assuming that w1 ¼ /1 ðr Þ e2ix t , we have:
@ 2 /1 2 @/1 þ þ 4x2 /1 r @r @r 2 2
1 ¼ 2iA x ðr 2 cosðx r Þ þ r x sinðx r ÞÞ
2
2 iðn 1ÞA2 x3 r 2 cos ðx r Þ;
ðA:15Þ
in which the analytical solution is: w1 ¼
C1 C2 ix A2 cosð2x r Þe2ix t þ sinð2x r Þe2ix t r r 8r2
x r ðnþ1ÞðCið4x r Þþ2Cið2x r Þþlnr Þsinð2x r Þ x r ðnþ1ÞðSið4x r Þþ2Sið2x r ÞÞcosð2x r Þ8cos2 ðx r Þ
! e2ix t :
ðA:16Þ Here Ci and Si stand for the Cosine Integral and Sine Integral functions, respectively. As before, C1 and C2 are unknown coefficients to be solved based on the imposed boundary conditions. For simplicity, we shall assume C1 = A⁄ and C2 = 0. Since we have assumed r ¼ Res r , t ¼ RUs t , w ¼ Rs Uw , h ¼ U 2 h , c ¼ c1 c , p ¼ p1 þ q1 U 2 p , A ¼ RAs U, x ¼ RUs x, t ¼ RUs t, r ¼ Res r,
e ¼ cU1 , k ¼ cx1 , the waves above can be transformed accordingly to the dimensional 1D wave as:
w0 ¼ A cosðkrÞ eixt ; w1 ¼ A cosð2krÞ e2ixt
total
dimensional
velocity
potential
is
then
References
Again, the solution for the linear acoustic wave is assumed as:
w0 ¼
The
w = w0 + e2w1 + O(e4).
ðA:17Þ ixA2 2
ððn þ 1Þkr sinð2krÞ þ
8U cosð2krÞ þ n 3Þ e2ixt ;
nþ1 4 ðA:18Þ
[1] M.F. Hamilton, D.T. Blackstock, Nonlinear Acoustics, Academic Press, New York, 1998. [2] S.W. Ohl, A. Shrestha, B.C. Khoo, A. Kishen, Characterizing bubble dynamics created by high-intensity focused ultrasound for the delivery of antibacterial nanoparticles into a dental hard tissue, Proceedings of the Institution of Mechanical Engineers Part H-Journal of Engineering in Medicine 224 (2010) 1285–1296. [3] S.W. Fong, E. Klaseboer, B.C. Khoo, Interaction of microbubbles with high intensity pulsed ultrasound, Journal of the Acoustical Society of America 123 (2008) 1784–1793. [4] G.B. Whitham, Linear and Nonlinear Waves, John Wiley & Sons, New York, 1974 (pp. 639). [5] T. Walsh, M. Torres, Finite element methods for nonlinear acoustics in fluids, Journal of Computational Acoustics 15 (2007) 353–375. [6] M.J. Crocker, Handbook of Acoustics, John Wiley & Sons, New York, 1998 (pp. 1469). [7] G. Liu, B. Ji, H.T. Chen, D.K. Liu, Antiplane harmonic elastodynamic stress analysis of an infinite wedge with a circular cavity, ASME Journal of Applied Mechanics 76 (2009) 061008. [8] G. Liu, H.T. Chen, D.K. Liu, B.C. Khoo, Surface motion of a half-space with triangular and semicircular hills under incident SH waves, Bulletin of the Seismological Society of America 100 (2010) 1306–1319. [9] D.K. Liu, B.Z. Gai, G.Y. Tao, Applications of the method of complex functions to dynamic stress concentrations, Wave Motion 4 (1982) 293–304. [10] G. Liu, B. Ji, D.K. Liu, Analytical solution for ground motion of a half space with a semi-cylindrical canyon and a beeline crack, Proceedings of the Royal Society A 464 (2008) 1905–1921. [11] Q.X. Wang, K.S. Yeo, B.C. Khoo, K.Y. Lam, Strong interaction between a buoyancy bubble and a free surface, Theoretical and Computational Fluid Dynamics 8 (1996) 73–88. [12] A.A. Doinikov, S.T. Zavtrak, On the mutual interaction of two gas bubbles in a sound field, Physics of Fluids 7 (1995) 1923–1930. [13] R. Ballarini, A rigid line inclusion at a bimaterial interface, Engineering Fracture Mechanics 37 (1990) 1–5. [14] D.T. Blackstock, Generalized Burgers equation for plane waves, Journal of the Acoustical Society of America 77 (1985) 2050–2053. [15] X.Z. Liu, J.L. Li, X.F. Gong, D. Zhang, Nonlinear absorption in biological tissue for high intensity focused ultrasound, Ultrasonics 44 (2006) e27–e30. [16] O.V. Bessonova, V.A. Khokhlova, Spatial structure of high intensity focused ultrasound beams of various geometry, Ultrasonics 17 (2009) 45–49. [17] M. Sun, D. Zhang, X.F. Gong, A novel approach for description of nonlinear field radiated from a concave source with wide aperture angle, Ultrasonics 44 (2006) e1435–e1438. [18] B.N. Kim, S.W. Yoon, Nonlinear parameter estimation in water-saturated sandy sediment with difference frequency acoustic wave, Ultrasonics 49 (2009) 438– 445. [19] G.V. Norton, R.D. Purrington, The Westervelt equation with viscous attenuation versus a causal propagation operator: A numerical comparison, Journal of Sound and Vibration 327 (2009) 163–172. [20] J.N. Tjotta, S. Tjotta, E. Vefring, Propagation and interaction of two collinear finite amplitude sound beams, Journal of the Acoustical Society of America 88 (1990) 2859–2870.
G. Liu et al. / Ultrasonics 54 (2014) 576–585 [21] Y.S. Lee, M.F. Hamilton, Time domain modeling of pulsed finite amplitude sound beams, Journal of the Acoustical Society of America 97 (1995) 906–917. [22] T. Kamakura, M. Akiyama, K. Aoki, A higher-order parabolic equation for describing nonlinear propagation of ultrasound beams, Acoustic Science and Technology 25 (2004) 163–165. [23] T. Kamakura, Two model equations for describing nonlinear sound beams, Japanese Journal of Applied Physics 43 (2004) 2808–2812. [24] C.C. Pozuelo, B. Dubus, J.A. Gallego-Juarez, Finite element analysis of the nonlinear propagation of high intensity acoustic waves, Journal of the Acoustical Society of America 106 (1999) 91–101. [25] C. Vanhille, C.C. Pozuelo, Numerical model for nonlinear standing waves and weak shocks in thermoviscous fluids, Journal of the Acoustical Society of America 109 (2001) 2660–2667. [26] C. Vanhille, C.C. Pozuelo, A numerical formulation for nonlinear ultrasonic waves propagation in fluids, Ultrasonics 42 (2004) 1123–1128. [27] Q.X. Wang, J.R. Blake, Non-spherical bubble dynamics in a compressible liquid. Part 1. Travelling acoustic wave, Journal of Fluid Mechanics 659 (2010) 191– 224.
585
[28] A. Prosperetti, A. Lezzi, Bubble dynamics in a compressible liquid. Part 1. First order theory, Journal of Fluid Mechanics 168 (1986) 457–478. [29] A. Lezzi, A. Prosperetti, Bubble dynamics in a compressible liquid. Part 2. Second order theory, Journal of Fluid Mechanics 185 (1987) 289– 321. [30] S. Fujikawa, T. Akamatsu, Effects of the non-equilibrium condensation of vapour on the pressure wave produced by the collapse of a bubble in a liquid, Journal of Fluid Mechanics 97 (1980) 481–512. [31] J.P. Franc, J.M. Michel, Fundamentals of Cavitation, Kluwer Academic Publishers, Boston, 2004 (pp. 305). [32] S.W. Gong, S.W. Ohl, E. Klaseboer, B.C. Khoo, Scaling law for bubbles induced by different external sources: Theoretical and experimental study, Physical Review E 81 (2010) 056317. [33] T.J. Matula, P.R. Hilmo, B.D. Storey, A.J. Szeri, Radial response of individual bubbles subjected to shock wave lithotripsy pulses in vitro, Physics of Fluids 14 (2002) 913–921.
Contents lists available at ScienceDirect
Ultrasonics journal homepage: www.elsevier.com/locate/ultras
Perturbation method for the second-order nonlinear effect of focused acoustic field around a scatterer in an ideal fluid Gang Liu, Pahala Gedara Jayathilake, Boo Cheong Khoo ⇑ Department of Mechanical Engineering, National University of Singapore, Kent Ridge Crescent, Singapore 119260, Singapore
a r t i c l e
i n f o
Article history: Received 18 March 2013 Received in revised form 23 July 2013 Accepted 18 August 2013 Available online 11 September 2013 Keywords: Perturbation method Westervelt equation Mach number Nonlinear acoustic wave
a b s t r a c t Two nonlinear models are proposed to investigate the focused acoustic waves that the nonlinear effects will be important inside the liquid around the scatterer. Firstly, the one dimensional solutions for the widely used Westervelt equation with different coordinates are obtained based on the perturbation method with the second order nonlinear terms. Then, by introducing the small parameter (Mach number), a dimensionless formulation and asymptotic perturbation expansion via the compressible potential flow theory is applied. This model permits the decoupling between the velocity potential and enthalpy to second order, with the first potential solutions satisfying the linear wave equation (Helmholtz equation), whereas the second order solutions are associated with the linear non-homogeneous equation. Based on the model, the local nonlinear effects of focused acoustic waves on certain volume are studied in which the findings may have important implications for bubble cavitation/initiation via focused ultrasound called HIFU (High Intensity Focused Ultrasound). The calculated results show that for the domain encompassing less than ten times the radius away from the center of the scatterer, the non-linear effect exerts a significant influence on the focused high intensity acoustic wave. Moreover, at the comparatively higher frequencies, for the model of spherical wave, a lower Mach number may result in stronger nonlinear effects. Ó 2013 Elsevier B.V. All rights reserved.
Introduction Applications of high intensity sonic energy in the industry have increased during the past decades [1]. The applications of intense sound are based on nonlinear effects produced by finite amplitude pressure variations [2,3]. It is well known that linear acoustic theory is based on the assumption of small amplitude waves and linear constitutive theory of the fluid medium [4]. Although these assumptions hold for many vibro–acoustic interactions, they become inapplicable in sound fields with high amplitude pressure [5]. Unlike the linear acoustic wave equation, the nonlinear counterpart can handle waves with large finite amplitudes and allow accurate modeling of nonlinear constitutive models in the fluid. Interesting phenomena undetected via linear acoustics can be observed in these studies. For example, waveform distortion, formation of shock waves, increased absorption, nonlinear interaction (as opposed to superposition) are seen when two sound waves are mixed, amplitude dependent directivity of acoustic beams, cavitation and sonoluminescence [6].
⇑ Corresponding author. Address: Singapore-MIT Alliance, 4 Engineering Drive 3, National University of Singapore, Singapore 117576, Singapore. E-mail addresses: [email protected], [email protected] (B.C. Khoo). 0041-624X/$ - see front matter Ó 2013 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.ultras.2013.08.011
The widely used linear wave equation used to describe the small signal (infinitesimal amplitude) propagation is the Helmholtz equation [7] given as below:
r2 H
1 @2H ¼ 0: c21 @t2
ð1Þ
Here, H can be pressure, velocity or velocity potential in the field of fluid mechanics, or displacement and displacement potential for solid mechanics [8]; c1 is the small signal sound speed since there are no shear waves inside the ideal fluid. For the linear wave propagation in elastic solids, the in-plane longitudinal and transverse wave [9] as well as anti-plane transverse wave [10] can co-exist and governed by the same Helmholtz equation. If we neglect the time term in Eq. (1), the Laplace equation is also widely employed for the incompressible ideal fluid [11,12] and static linear solid mechanics [13]. To the best of authors’ knowledge, there are various models to simulate the nonlinear characteristic of the acoustic wave propagating through the fluid. For instance, the one-dimensional Burgers equation has been found to be an excellent approximation of the conservation equations for plane progressive waves of finite amplitude in a thermoviscous fluid [14]. An effective model that combined effects of diffraction, absorption, and nonlinearity in directional sound beams (i.e. radiated from sources with
577
G. Liu et al. / Ultrasonics 54 (2014) 576–585
assumed a solution for the acoustic pressure consisting of two terms in the form of p = p1 + p2, where p1 represents the first order approximation and p2 the second order correction. Here p2 is much smaller than p1 (p2 p1). Neglecting the terms of third and higher orders, we can obtain the approximate equation as:
r2 ðp1 þ p2 Þ
1 @ 2 ðp1 þ p2 Þ b @ 2 p21 þ ¼ 0: 2 2 c1 q1 c41 @t2 @t
ð4Þ
According to Eq. (4), the second order solution for the acoustic wave can be written as:
p ¼ P1 expðixtÞ þ P2 expð2ixtÞ:
Fig. 1. Schematic description of the model for multiple acoustic waves focused around a scatterer inside an ideal fluid.
dimensions that are large compared with the characteristic wavelength) are taken into account by the Khokhlov–Zabolotskaya– Kuznestov (KZK) parabolic wave equation [15,16]. Several other models have been developed in response to specific needs. For example, the Westervelt equation is an almost incidental product of Westervelt’s discovery of the parametric arrays [17–19]. In the past, much computational work have been done on nonlinear effects in beams based on the KZK equations [20–23]. The one dimensional case of Westervelt equation is studied extensively by finite element method [24]. However, when acoustic fields of more complex conditions are considered, advanced numerical calculation related to finite-differential becomes necessary [25,26]. In this work, our interest is to develop the analytical solution of the multiple harmonic acoustic waves focused on the volume near a scatterer, where the second order nonlinear effects dominate (see Fig. 1). Our study has important implications for further work on bubble/nucleation cavitation by HIFU (High Intensity Focused Ultrasound) and others. (Our approach based on asymptotic perturbation expansion of small parameter and non-dimensionalisation carried out w.r.t the dimension of the scatterer, however, is not applicable in the absence of the scatterer.) Second order nonlinear solution for Westervelt equation On the problems of multiple acoustic waves focusing on intensive volume, one effective model governing the dynamics is the Westervelt equation as given below:
r2 p
1 @2p b @ 2 p2 d @3p þ þ ¼ 0: 2 2 c21 @t c41 @t 3 q1 c41 @t
ð2Þ
Here, the third term is the nonlinear term and the forth term stands for the thermoviscous effect. q1 is the static medium density, c1 is the wave speed of the acoustic wave in the undisturbed liquid, b = (c + 1)/2 is the coefficient of nonlinearity, and c is defined as the ratio of specific heats of the gas. Here, d = m(# + (c 1)/Pr) is the diffusivity associated with sound absorption, where m = l/q1 is the kinematic viscosity, # = 4/3 + lB/l is the viscosity number, lB is the bulk viscosity coefficient, and Pr is Prandtl number. If we consider the case of ideal fluids, the forth term can be set equal to zero, and the reduced Westervelt equation is given as:
r2 p
1 @2p b @ 2 p2 þ ¼ 0: c21 @t2 q1 c41 @t 2
ð3Þ
To solve the second-order nonlinear equation, the successive approximation method can be applied. In this method, it is
The proposed analytical solutions of the one dimensional Eq. (4) for the planar, cylindrical, and spherical waves are obtained as shown below. It is pointed out, however, that the coefficients on the right hand side of the initial governing equation for the second-order nonlinear term shown as Eq. (10) in [24] are deemed incorrect. In particular, Eq. (10) in Ref. [24] is shown with the term 2 2k b=ðq1 c21 Þ, which we reckoned should be reflected as 2 4k b=ðq1 c21 Þ after the separation of the variables. The details for our derivations were not shown here, and the results with different expressions and coefficients are given below as:
pp ¼ P 0 eikx eixt þ pc
¼ P0
qffiffiffiffiffiffi 2
iðkrp=4Þ
pkr e ð2Þ
bP20 ð1 þ 2ikxÞe2ikx e2ixt ; 2q1 c21 4ibP20 2 1 c1
eixt þ pq
P0 H0 ðkrÞ eixt þ
2ibkrP20 q1 c21
ð2Þ
ð5Þ
eip=2 e2ikr e2ixt 2
ð6Þ
ðH0 ðkrÞÞ e2ixt ;
! 2 1 P 0 ikr ixt ibkP 0 ln r X ð1Þn ðn 2Þ! e2ikr e2ixt : ps ¼ e e þ þ n1 r q1 c21 r ð4ikÞ rn n¼2 ð7Þ It is clear that the first term in Eqs. (5)–(7) are the general solutions for the one dimensional linear acoustic wave with time harmonic condition, and the second terms are the approximate solutions for the second-order nonlinearity which satisfy the Sommerfield radiation condition at the far field. Here, P0 stands for harmonic excitation of the wave amplitude, and k ¼ x=c1 ¼ 2p=k, where x is the circular frequency of the incident wave, and k stands for the wave length. It should be noted that the calculated results plotted for the second-order solutions for the cylindrical and spherical waves differ from Ref. [24] on two counts. As discussed above, it is due to the different governing equations for the second-order nonlinear term and the boundary conditions prescribed in [24] is the Neumann type as it enforces the normal derivative of pressure at the boundary of the rigid scatterer. However, for our simulation, we imposed the Dirichlet boundary conditions. See Section 4 below for further discussion on the observed difference. Perturbative method with a small parameter for the nonlinear acoustic wave Mathematical formulation of the nonlinear acoustic wave From the analytical solutions in Section 2 for the one dimensional Westervelt equation, it is clear that one can only get the pressure for the fluid domain without the explicit results for the velocity and velocity potential. Below, in order to give a more clearer description about the nonlinear effects of the (nonlinear) acoustic wave scattering by inclusion, we propose a mathematical formulation and a different model for the nonlinear acoustic wave that goes beyond the traditional models based on the Burgers
578
G. Liu et al. / Ultrasonics 54 (2014) 576–585
equation, KZK equation, and Westervelt equation. The compressible potential flow theory is invoked together with the dimensionless formulation and asymptotic perturbation expansion referring to the small parameter Mach number. It will be shown that our approach allows the decoupling of the potential and enthalpy terms to the second order, and which accord much flexibility for the calculation of the other physical quantities inside the fluid domain. For our model, we shall take the liquid as inviscid and compressible. As such, the liquid flow is governed by the equation of mass conservation:
~ @u 1 ~ ru ~ ¼ rp: þu @t q
ð9Þ
In many fluid dynamic problems, significant liquid compressibility due to high speed motion occurs when thermal effects in the liquid are unimportant. We therefore assume that thermal effects in the liquid are insignificant. With this assumption, the liquid state is defined by a single thermodynamic variable. The sound speed c and the enthalpy h of the liquid can be defined as follows [27]:
h¼
Z
p
dp
q
p1
ð10Þ
;
where the reference pressure p1 is the pressure in the undisturbed liquid (hydrostatic pressure). By assuming that the flow is irrotational, we may introduce a ~ ¼ rw. By applying velocity potential w such that u _ q_ rp=rq, and u ~ ¼ rw, where the dot ‘‘’’ stands c2 ¼ dp=dq p= for the differential w.r.t time, we can re-write Eq. (8) as:
r2 w þ
1 @h þ rw rh ¼ 0: 2 c @t
ð11Þ
Similarly for Eq. (9), the integration leads to the unsteady Bernoulli equation:
@w 1 þ jrwj2 ¼ h h1 : @t 2
ð12Þ
Here, h1 may be set equal to zero since the enthalpy is referenced to the undisturbed fluid at infinity. To find the expressions for the speed of sound c and enthalpy h, we use the Tait model equation of state for water to relate the pressure and density as follows [28,29]:
p¼
p1 þ B
qn1
qn B:
ð13Þ
The values B = 3049.13 bars and n = 7.15 give an excellent fit based on experiments in water for up to 105 bars [30]. Here, p1 and q1 are the pressure and density in the undisturbed liquid (hydrostatic model). From Eq. (13) together with c2 = dp/dq, we get:
c2 ¼
nðp þ BÞ
q
¼
n
q1
p p1
q1
2 1 p p1 þ ...; 2c21 q1
ð16Þ
1 1 pp ¼ ðn 1Þ 4 1 þ . . . c2 c21 c1 q1
ð17Þ
ð8Þ
and the momentum conservation equation:
dp ; dq
h¼
Non-dimensional formulation of the governing equations
@q ~ Þ ¼ 0; þ r ðqu @t
c2 ¼
Next, we expand the enthalpy h and the sound speed c around the equilibrium pressure p1 using a Taylor series expansion as follows [28]:
ðp1 þ BÞ1=n ðp þ BÞðn1Þ=n
ð14Þ
!
ðn1Þ=n pþB 1 ; p1 þ B
r¼
Rs
e
r ;
c ¼ c1 c ;
t¼
Rs t ; U
w ¼ Rs Uw ;
h ¼ U 2 h ;
p ¼ p1 þ q1 U 2 p ;
ð18Þ ð19Þ
where the Mach number e = U/c1 (is also the small parameter to be used for the perturbation analysis). It is evident that e – 0, or otherwise the dimensionless formulation in Eq. (18) will have a singularity. From the above, we have: 2 @2w 2 U @ w ¼ e ; @r 2 Rs @r 2
@h U 2 @h ¼e ; @r Rs @r
@w @w ¼ eU ; @r @r
@w @w ¼ U2 ; @t @t
@h U 3 @h ¼ @t Rs @t
ð20Þ
ð21Þ
Next, by substituting Eqs. (20), (21), and (19) into (11) and (12), respectively, we obtain the governing equations in terms of dimensionless variables:
@ 2 w 1 @h 2 @w @h ¼ 0; þ þ e c2 @t @r2 @r @r
ð22Þ
2 @w 1 2 @w þ e þ h ¼ 0: @t 2 @r
ð23Þ
We stress that the dimensionless variable r⁄ also stands for the space vector dependent on the dimensional coordinate system which is planar, cylindrical or spherical. Next, from Eq. (15), the dimensionless speed of sound of the liquid is:
c2 ¼ 1 þ e2 ðn 1Þh ;
ð24Þ
and the dimensionless enthalpy of the liquid in Eq. (15) is given by:
h ¼
1
e2
1 n1
! ðn1Þ=n pþB 1 : p1 þ B
ð25Þ
Substituting p = p1 + q1U2p⁄ and h = U2h⁄ into Eqs. (16) and (17), and referring to the Taylor series expansion, we can get.
h ¼ p
1 2 2 e p þ oðe2 Þ; 2
1 ¼ 1 ðn 1Þe2 p þ oðe2 Þ: c2
and
c2 c21 c2 h¼ ¼ 1 n1 n1
We would like to introduce two typical scales Rs and U, where Rs denotes the radius of the scatterer and U (p1/q1)1/2 stands for the liquid velocity near the scatterer. The dimensionless quantities indicated by asterisk are given as below:
ð15Þ
ð26Þ
ð27Þ
By substituting Eqs. (26) and (27) into Eq. (22), we can obtain:
@ w @h @h 2 @w @h þ oðe2 Þ ¼ 0: þ þ e ðn 1Þh @r2 @t @r @r @t 2
where c21 ¼ nðpq1 þBÞ is the wave speed of the acoustic wave in the 1 undisturbed liquid.
ð28Þ
G. Liu et al. / Ultrasonics 54 (2014) 576–585
Next, the perturbation expansions for the potential w⁄ and enthalpy h⁄ are defined as follows:
w ðr ; t Þ ¼ w0 ðr ; t Þ þ e2 w1 ðr ; t Þ þ oðe4 Þ;
ð29Þ
h ðr ; t Þ ¼ h0 ðr ; t Þ þ e2 h1 ðr ; t Þ þ oðe4 Þ;
ð30Þ
where these are then introduced into Eqs. (23) and (28). We can get the separate governing equations for the first and second order w.r.t. e2 as below:
@w0 þ h0 ¼ 0; @t 2 @w1 1 @w0 þ þ h1 þ oðe4 Þ ¼ 0; 2 @r @t
ð31Þ ð32Þ
2
@ w0 @h0 þ ¼ 0; @r 2 @t
ð33Þ
@ 2 w1 @h1 @w0 @h0 @h0 þ þ ðn 1Þh0 þ oðe4 Þ ¼ 0: @r 2 @t @r @r @t
ð34Þ
@ 2 w0 ¼ 0; @t2
r2 w1
@ 2 w1 1 @ @ jr w0 j2 r w0 ðr w0 Þ 2 @t @t @t 2
ðn 1Þ
@w0 @ 2 w0 þ oðe4 Þ ¼ 0: @t @t 2
ð35Þ
ð36Þ
where the Hamilton operator r⁄ is defined in terms of r⁄. As elucidated before, the parameter n is assumed to reasonably match the experimental pressure-density relation for water up to 10 [5] bars at n = 7.15 [30]. Analytical solution for the one-dimensionless equation Analytical solution for the plane wave The analytical solution for the dimensional second-order nonlinear plane wave is given below as (details are provided in the Appendix A):
w0 ¼ A cosðkrÞ eixt ;
ð37Þ
w1 ¼ A cosð2krÞ e2ixt ixA2 nþ1 e2ixt ; ðn þ 1Þkr sinð2krÞ þ cosð2krÞ þ n 3 4 8U 2 ð38Þ where A, k, and x are the amplitude, wave number, and frequency, respectively. Analytical solution for the cylindrical wave The analytical solution for the dimensional second-order nonlinear cylindrical wave is given below as: ð2Þ
w0 ¼ A H0 ðkrÞ eixt ;
Analytical solution for the spherical wave The analytical solution for the dimensional second-order nonlinear spherical wave is given below as:
w0 ¼
w1 ¼
ARs cosðkrÞ eixt ; er
ð41Þ
ARs ixA2 R2 cosð2krÞ e2ixt 2 s er 8U e2 r2 ! krðn þ 1ÞðCið4krÞ þ 2Cið2krÞ þ ln Rers Þsinð2krÞ e2ixt : krðn þ 1ÞðSið4krÞ þ 2Sið2krÞÞcosð2krÞ 8cos2 ðkrÞ ð42Þ
By substituting Eqs. (37) and (38), Eqs. (39) and (40), Eqs. (41) and (42) into w = w0 + e2w1 + O(e4), respectively, we obtain the corresponding asymptotic expression for the velocity potential up to and including O(e2). Results and discussions
By combining Eqs. (31) and (33), Eqs. (32) and (34), we can obtain the fundamental equation for the second order nonlinear acoustic wave (dimensionless formulation) as:
r2 w0
579
ð39Þ
ipA2 x3 e2 ð2Þ w1 ¼ AH0 ð2krÞe2ixt þ 2U 4 0 1 R 2 2 ð2Þ ð2Þ J 0 ð2krÞ ððn1ÞðH0 ðkrÞÞ 2ðH1 ðkrÞÞ ÞY 0 ð2krÞrdr @ A e2ixt : R 2 2 ð2Þ ð2Þ Y 0 ð2krÞ ððn1ÞðH0 ðkrÞÞ 2ðH1 ðkrÞÞ ÞJ0 ð2krÞrdr ð40Þ
The obtained analytical solutions are analyzed with regard to the nonlinear effects around the volume inside the liquid where the multiple incident waves are focused. It is important to emphasize that the initial interactions among the multiple linear harmonic sources are not taken into account since each source is not having a sufficiently large intensity to induce nonlinear deformation. If one assume the scatterer to be the equilibrium bubble, the Rayleigh–Plesset equation [31] or modified Rayleigh model [32] as well as the relationship between the pressure in the fluid and the gas pressure inside the bubble [33] (where the nonlinear effects have not been taken into account previously), can be employed to describe the bubble cavitation/deformation. Alternatively, we may use the general dimensional equation dr/dt = rw together with Eqs. (41) and (42) to construct and express the nonlinear variation of the bubble radius as well. In our paper, we shall assume the radius of the scatterer Rs = 1. Since we know that the nonlinear acoustic wave effect will only dominate the fluid region and its vicinity outside the scatterer, the distance defined as x or r should larger than the radius of the scatterer (Rs = 1); in the subsequent analysis the effective distance Re will range from 1 to 100. By using the analytical solutions for the one dimensional Westervelt equation, the pressure distribution around the scatterer for the fundamental and the second order harmonic of planar, cylindrical and spherical coordinates in water are calculated. Here, in our calculation, the wave speed c1 is assumed to be 1500 m/s, water density q1 is 1000 kg/m3, nonlinear coefficient b is 3.5, initial pressure P0 is 4 105 Pa, and frequency is f = 15, 000 Hz [24]. (For example, one can interpret the presence of four incident waves at pressure intensity of Pi = 105 Pa each. As such, there are the multiple incident waves focusing with linear superposition giving rise to P0 = 4 105 Pa at the vicinity of the scatterer.) The calculated results for the pressure ratio between the second order harmonic to fundamental (first order term) are displayed in Fig. 2. Fig. 2a shows that the ratio between the second harmonic to the first order harmonic for plane wave increases linearly with increasing distance away from the focused point. For the second harmonics of spherical and cylindrical waves as shown in Fig. 2b and c, respectively, the increase with distance from the focused point is no longer linear. The mentioned increase for the spherical wave depicts a higher rate of change with distance compared to the cylindrical wave; both of them depict a lower magnitude and rate of increase compared to the planar wave. As mentioned in Section 2, the differences observed for Fig. 2b and c applicable to the cylindrical and spherical waves, respectively, vis-à-vis the corresponding results in [24] are attributed to the different governing
580
G. Liu et al. / Ultrasonics 54 (2014) 576–585
Fig. 3 presents a more explicit comparison between the fundamental and the second order harmonic pressure amplitudes in planar, cylindrical and spherical coordinates. According to the definition of the pressure level (dB) set at 20 log10(p/preference), we can plot the pressure level distribution w.r.t. the distance away from the focused point (the reference pressure is assumed to be 106 Pa, and the range is from 1.0 to 100 which is different from Ref. [24] found in the range from 0 to 1.0). The calculated pressure amplitude variation shown in Fig. 3, however, can be compared to the results obtained by [24]. From Fig. 3, it is observed that the second order harmonic distribution for the planar wave increases with distance and becomes larger than the first harmonic at about 25 units. The second order harmonic distribution of cylindrical wave seems to exhibit minor variation. This is in line with the approximate expanded solution for the second order of the cylindrical 4ibP 2
wave taken to be pq c02 eip=2 e2ikr e2ixt as shown in Eq. (6). It is 1 1 well known that the absolute value of this complex number has no variation with respect to the increase of the distance r. Lastly, it is observed that the second order for the harmonic spherical wave initially increases slightly and then decreases with the increment of the distance away from the focused point. The summation of the fundamental and second order terms are shown in Fig. 4. The numerical results of the linearized governing equations are separately calculated by Finite Difference method to check on the accuracy and consistency of the analytical results. It is evident that our analytical results concur well with the numerical results. It is also observed that the amplitude for the cylindrical and spherical waves decrease monotonically. In particular, the spherical wave decreases at a faster rate than the cylindrical wave since the pressure for the spherical wave varies with the approximate factor r1 (see Eq. (7)), whereas the cylindrical waves with 1 the factor of r2 (see Eq. (6)). On the contrary, the pressure amplitude increases with distance away from the focused point for the planar wave. This can be attributed to the absolute value of the second order pressure planar wave increasing with respect to the distance away from the focused point (see Eq. (5)). It is clear that within the distance of less than 10, the decrease of the aggregate pressure amplitudes of cylindrical and spherical wave appear to be very fast and the tendency to moderate beyond this point for the region between 10 and 100. Being so, one can suggest that the domain encompassing less than ten times of the radius of the scatterer, the non-linear effect exerts a significant influence on the focused high intensity acoustic wave. Next, we shall discuss the nonlinear effect for the focused acoustic wave with the small parameter Mach number by
Fig. 2. The ratio of the pressure second harmonic to the fundamental term vs. the variation of the distance away from the focused point.
equation for the second order nonlinear term and different boundary conditions imposed. It is clear that although the general/broad variation of increasing of magnitude of second harmonic to fundamental ratio with distance from the focal point is similar, the observed trend and rate of increase are different (for our case the boundary conditions of Dirichlet type are imposed). For the second harmonic to fundamental ratio versus distance from focused point observed between our Fig. 2a and Fig. 1 of [24] applicable for the plane wave, the difference in magnitude can also be attributed to the different governing equation for the second order nonlinear term employed.
Fig. 3. The comparison between the analytical results of planar, cylindrical and spherical wave including the fundamental and the second harmonic.
G. Liu et al. / Ultrasonics 54 (2014) 576–585
Fig. 4. The variation of pressure amplitude vs. the distance away from the focused point for the planar, cylindrical and spherical wave.
compressible potential flow theory. We will assume the amplitude of velocity potential A = 1, and radius of the scatterer Rs = 1. As interpreted earlier, x⁄r⁄ = kr is defined as the dimensionless wave number and the velocity potential is only dependent on the variable U in the frequency domain and the small parameter e = U/c1 (c1 is a constant at 1500 m/s). For the numerical calculations, we shall take the dimensionless wave number as kRs = 1.0, 3.0 and
581
5.0 and the small parameter Mach number e as 0.1, 0.3 and 0.5, so that U = ec1 becomes 150, 450 and 750, respectively. Fig. 5 shows the variation of the second order term of plane wave w.r.t. the increase of the distances from the focused point for different parameters. It is evident that with the increase of the wave number k (shown as the horizontal rows), we can observe that the general tendency of the amplitude for the second order term to increase and the variation takes on shorter periodic time. For the fixed low wave number of k = 1.0 (first column in Fig. 5), the increase of the Mach number (e = 0.1, 0.3 and 0.5) can bring about comparatively large amplitude for the second order term. However, the maximum amplitude for each second term appears to be constant without any variation w.r.t. the increase of distance from the focused point. This phenomena can be interpreted as that the low frequency stands for quasi-static model, and thus the variation of the distances have lesser influence on the second-order nonlinear term. At the high frequency of k = 5.0 (third column in Fig. 5), the maximum amplitudes tend to increase w.r.t. the position away from the focused point. Fig. 6 illustrates the variation of the nonlinear effect for the cylindrical wave at different wave number, Mach number and distance away from the focused point. The general tendency is that with the increase of the distance away from the focused point, this results in ever smaller amplitude for the second order term. For the fixed wave number k, the increase of the Mach number will bring about a larger amplitude for the second order term. With the fixed higher Mach number (e = 0.5), the increment of the wave number can also weaken the amplitude for the second order term. Fig. 7
Fig. 5. The variation of the second order term of plane wave vs. the variation of wave number k, Mach number e and the distance away from the focused point.
582
G. Liu et al. / Ultrasonics 54 (2014) 576–585
Fig. 6. The variation of the second order term for cylindrical wave vs. the different wave number, Mach number and the distance away from the focused point. Fig. 7. The variation of the second order term for spherical wave vs. the different wave number, Mach number and the distance away from the focused point.
shows the variation of the nonlinear effect for the spherical wave for the different wave number, Mach number and distance away from the focused point. Similar to that observed for the cylindrical wave, the general tendency is that the increase of the distance away from the focused point will result in smaller amplitude for the second order term. At the fixed low frequency of k = 1.0, the increase of the Mach number (e = 0.1, 0.3 and 0.5) will bring about a larger amplitude for the second order term. Conversely, at the higher frequencies of k = 3.0 and 5.0, a lower Mach number results
in stronger nonlinear effects and this trend is very different from the cylindrical wave. Finally, Fig. 8 gives the total velocity potential (the summation of the first order term and second order term) at k = 2.0 and e = 0.3. We can observe the variation of the velocity potential for both the cylindrical wave and spherical wave appears to be consistent with the conclusion from the Westervelt model for the pressure distribution. The general tendency of the velocity potential will decrease
583
G. Liu et al. / Ultrasonics 54 (2014) 576–585
@ 2 w0 @ 2 w0 ¼ 0; @r 2 @t2
ðA:1Þ
2 @ 2 w1 @ 2 w1 1 @ @w0 @w @ 2 w0 @w ¼ þ 0 þ ðn 1Þ 0 2 2 2 @t @r @r @r @r @t @t @t
@ 2 w0 : @t2
ðA:2Þ
Here, we would like to point out that r⁄ can be directly replaced by the Cartesian coordinate x⁄ for the plane wave. We propose the one dimensionless solution for the linear acoustic wave in Eq. (A.1) as below:
w0 ¼ A cosðk r =eÞ eix t A cosðx r Þ eix t ;
ðA:3Þ
where A⁄, k⁄, x⁄ are the dimensionless amplitude, wave number, frequency, respectively, and defined as:
Fig. 8. The velocity potential distribution near the scatterer (the summation of the first order term and the second order term) at k = 2.0 and e = 0.3.
w.r.t. the increase of the distance away from the focused point while the amplitude of the velocity potential for the planar wave seems to be periodic constant. For the spherical wave, since our analytical solution is expressed by cosine/sine function (or cosine/sine integral function) rather than the exponential function, a characteristic feature of this function is periodic rather than just a monotonic curve as from the Westervelt equation. Conclusions In this study, we adopted two nonlinear models to investigate the multiple incident acoustic waves focused on certain domain where the nonlinear effect is not negligible in the vicinity of the scatterer. The general solutions for the one dimensional Westervelt equation with different coordinate systems (plannar, cylindrical and spherical) were analytically obtained based on the perturbation method with keeping only the second order nonlinear terms. By introducing the small parameter (Mach number), we applied the compressible potential flow theory and proposed a dimensionless formulation and asymptotic perturbation expansion for the velocity potential and enthalpy which is different from the existing fractional nonlinear acoustic models (e.g. the Burgers equation, KZK equation and Westervelt equation). Our analytical solutions and numerical calculations have shown the general tendency of the velocity potential and pressure to decrease w.r.t. the increase of the distance away from the focused point. At least, within the region which is about 10 times the radius of the scatterer, the nonlinear effect exerts a significant influence on the distribution of the pressure and velocity potential. It is also interesting to note that at high frequencies, lower Mach numbers appear to bring out even stronger nonlinear effects for the spherical wave. Our approach with a small parameter for the cylindrical and spherical waves could serve as an effective analytical model to simulate the focused nonlinear acoustic near the scatterer in an ideal fluid and be applied to study bubble cavitation dynamic associated with HIFU in our future work. Appendix A Analytical. solution for the plane wave For the one dimensional (1D) plane wave, the governing equations are (the higher order small terms will be omitted subsequently):
A ¼
A ; Rs U
k ¼ Rs k;
x ¼
Rs x; U
r ¼
e Rs
r;
t ¼
U t: Rs
ðA:4Þ
It is evident that Eq. (A.3) is a solution of Eq. (A.1). By substituting Eq. (A.3) into Eq. (A.2) and assuming that w1 ¼ /1 ðr Þ e2ix t , we obtain the simplified second-order non-homogeneous differential equation as:
@ 2 /1 2 2 þ 4x2 /1 ¼ 2iA x3 sin ðx r Þ iðn 1ÞA2 x3 @r2 cos2 ðx r Þ:
ðA:5Þ
The general analytical solution for this second-order nonhomogeneous Eq. (A.5) is:
w1 ¼ C 1 cosð2x r Þ e2ix t þ C 2 sinð2x r Þ e2ix t ix A2 ððn þ 1Þx r sinð2x r Þ 8 nþ1 cosð2x r Þ þ n 3Þ e2ix t : þ 4
ðA:6Þ
where C1 and C2 are unknown coefficients to be determined by the imposed boundary conditions. For simplicity in the discussion below, the coefficients for the second-order nonlinear terms can be reasonably assumed as C1 = A⁄, C2 = 0 without loss of generality. Analytical. solution for the cylindrical wave For the one dimensional (1D) cylindrical wave, the governing equations are:
@ 2 w0 1 @w0 @ 2 w0 þ ¼ 0; r @r @r 2 @t 2
ðA:7Þ
2 @ 2 w1 1 @w1 @ 2 w1 1 @ @w0 @w0 @ 2 w0 þ ¼ þ r @r 2 @t @r @r 2 @r @r @t @t 2 þ ðn 1Þ
@w0 @ 2 w0 : @t @t 2
ðA:8Þ
The proposed solution for the linear acoustic wave (second kind of Hankel function with zero order) is assumed as: ð2Þ
w0 ¼ A H0 ðx r Þ eix t :
ðA:9Þ
By substituting Eq. (A.9) into Eq. (A.8) and assuming that w1 ¼ /1 ðr Þ e2ix t , we obtain the simplified second-order nonhomogeneous differential equation as: 2 @ 2 /1 1 @/1 2 ð2Þ þ þ 4x2 /1 ¼ 2iA x3 ðH1 ðx r ÞÞ r @r @r2 2
iðn 1ÞA2 x3 ðH0 ðx r ÞÞ ; ð2Þ
ðA:10Þ
584
G. Liu et al. / Ultrasonics 54 (2014) 576–585
and the analytical solution is (the softwave named MAPLE will be used here): 2ix t
w1 ¼ C 1 J0 ð2x r Þe 0 B @
J0 ð2x r Þ
R
þC 2 Y 0 ð2x r Þe 2
ð2Þ ððn1ÞðH0 ð
Y 0 ð2x r Þ
x r ÞÞ
R
2ix t
x r ÞÞ
ð2Þ
w0 ¼ A H0 ðkrÞ eixt ;
ðA:19Þ
3 2 A
ipx þ 2
1
2
ð2Þ 2ðH1 ð 2
ð2Þ ððn1ÞðH0 ð
Therefore, the dimensional 1D cylindrical wave is given as:
x r ÞÞ ÞY 0 ð2x r Þr dr C Ae2ix t :
ð2Þ 2ðH1 ð
2
x r ÞÞ ÞJ0 ð2x r Þr dr
ðA:11Þ Here, J0() is the Bessel function of order zero, Y0() is the second ð2Þ kind of Bessel function of order zero, and H0 ðÞ ¼ J 0 ðÞ iY 0 ðÞ is second kind Hankel function with order zero. For ease of discussion and without loss of generality, one can reasonably assume that C1 = A⁄ and C2 = iA⁄.
ipA2 x3 e2 ð2Þ w1 ¼ A H0 ð2krÞ e2ixt þ 2U 4 0 1 R 2 2 ð2Þ ð2Þ J 0 ð2krÞ ððn 1ÞðH0 ðkrÞÞ 2ðH1 ðkrÞÞ Þ Y 0 ð2krÞ rdr A e2ixt ; @ R 2 2 ð2Þ ð2Þ Y 0 ð2krÞ ððn 1ÞðH0 ðkrÞÞ 2ðH1 ðkrÞÞ Þ J 0 ð2krÞ rdr ðA:20Þ
and the dimensional 1D spherical wave as:
w0 ¼
w1 ¼
Analytical. solution for the spherical wave For the one dimensional (1D) spherical wave, the governing equations are:
ARs cosðkrÞ eixt ; er
ðA:21Þ
ARs ixA2 R2 cosð2krÞ e2ixt 2 s er 0 8U e2 r2 1 krðn þ 1Þ Cið4krÞ þ 2Cið2krÞ þ ln Rers sinð2krÞ A e2ixt : @ krðn þ 1ÞðSið4krÞ þ 2Sið2krÞÞ cosð2krÞ 8cos2 ðkrÞ ðA:22Þ
@ 2 w0 2 @w0 @ 2 w0 þ ¼0 r @r @r 2 @t 2
ðA:12Þ
2 @ 2 w1 2 @w1 @ 2 w1 1 @ @w0 @w @ 2 w0 þ ¼ þ 0 2 r @r 2 @t @r 2 @r @r @r @t @t @w @ 2 w0 þ ðn 1Þ 0 : @t @t 2
ðA:13Þ
A cosðx r Þ eix t : r
ðA:14Þ
By substituting Eq. (A.14) into Eq. (A.13) and assuming that w1 ¼ /1 ðr Þ e2ix t , we have:
@ 2 /1 2 @/1 þ þ 4x2 /1 r @r @r 2 2
1 ¼ 2iA x ðr 2 cosðx r Þ þ r x sinðx r ÞÞ
2
2 iðn 1ÞA2 x3 r 2 cos ðx r Þ;
ðA:15Þ
in which the analytical solution is: w1 ¼
C1 C2 ix A2 cosð2x r Þe2ix t þ sinð2x r Þe2ix t r r 8r2
x r ðnþ1ÞðCið4x r Þþ2Cið2x r Þþlnr Þsinð2x r Þ x r ðnþ1ÞðSið4x r Þþ2Sið2x r ÞÞcosð2x r Þ8cos2 ðx r Þ
! e2ix t :
ðA:16Þ Here Ci and Si stand for the Cosine Integral and Sine Integral functions, respectively. As before, C1 and C2 are unknown coefficients to be solved based on the imposed boundary conditions. For simplicity, we shall assume C1 = A⁄ and C2 = 0. Since we have assumed r ¼ Res r , t ¼ RUs t , w ¼ Rs Uw , h ¼ U 2 h , c ¼ c1 c , p ¼ p1 þ q1 U 2 p , A ¼ RAs U, x ¼ RUs x, t ¼ RUs t, r ¼ Res r,
e ¼ cU1 , k ¼ cx1 , the waves above can be transformed accordingly to the dimensional 1D wave as:
w0 ¼ A cosðkrÞ eixt ; w1 ¼ A cosð2krÞ e2ixt
total
dimensional
velocity
potential
is
then
References
Again, the solution for the linear acoustic wave is assumed as:
w0 ¼
The
w = w0 + e2w1 + O(e4).
ðA:17Þ ixA2 2
ððn þ 1Þkr sinð2krÞ þ
8U cosð2krÞ þ n 3Þ e2ixt ;
nþ1 4 ðA:18Þ
[1] M.F. Hamilton, D.T. Blackstock, Nonlinear Acoustics, Academic Press, New York, 1998. [2] S.W. Ohl, A. Shrestha, B.C. Khoo, A. Kishen, Characterizing bubble dynamics created by high-intensity focused ultrasound for the delivery of antibacterial nanoparticles into a dental hard tissue, Proceedings of the Institution of Mechanical Engineers Part H-Journal of Engineering in Medicine 224 (2010) 1285–1296. [3] S.W. Fong, E. Klaseboer, B.C. Khoo, Interaction of microbubbles with high intensity pulsed ultrasound, Journal of the Acoustical Society of America 123 (2008) 1784–1793. [4] G.B. Whitham, Linear and Nonlinear Waves, John Wiley & Sons, New York, 1974 (pp. 639). [5] T. Walsh, M. Torres, Finite element methods for nonlinear acoustics in fluids, Journal of Computational Acoustics 15 (2007) 353–375. [6] M.J. Crocker, Handbook of Acoustics, John Wiley & Sons, New York, 1998 (pp. 1469). [7] G. Liu, B. Ji, H.T. Chen, D.K. Liu, Antiplane harmonic elastodynamic stress analysis of an infinite wedge with a circular cavity, ASME Journal of Applied Mechanics 76 (2009) 061008. [8] G. Liu, H.T. Chen, D.K. Liu, B.C. Khoo, Surface motion of a half-space with triangular and semicircular hills under incident SH waves, Bulletin of the Seismological Society of America 100 (2010) 1306–1319. [9] D.K. Liu, B.Z. Gai, G.Y. Tao, Applications of the method of complex functions to dynamic stress concentrations, Wave Motion 4 (1982) 293–304. [10] G. Liu, B. Ji, D.K. Liu, Analytical solution for ground motion of a half space with a semi-cylindrical canyon and a beeline crack, Proceedings of the Royal Society A 464 (2008) 1905–1921. [11] Q.X. Wang, K.S. Yeo, B.C. Khoo, K.Y. Lam, Strong interaction between a buoyancy bubble and a free surface, Theoretical and Computational Fluid Dynamics 8 (1996) 73–88. [12] A.A. Doinikov, S.T. Zavtrak, On the mutual interaction of two gas bubbles in a sound field, Physics of Fluids 7 (1995) 1923–1930. [13] R. Ballarini, A rigid line inclusion at a bimaterial interface, Engineering Fracture Mechanics 37 (1990) 1–5. [14] D.T. Blackstock, Generalized Burgers equation for plane waves, Journal of the Acoustical Society of America 77 (1985) 2050–2053. [15] X.Z. Liu, J.L. Li, X.F. Gong, D. Zhang, Nonlinear absorption in biological tissue for high intensity focused ultrasound, Ultrasonics 44 (2006) e27–e30. [16] O.V. Bessonova, V.A. Khokhlova, Spatial structure of high intensity focused ultrasound beams of various geometry, Ultrasonics 17 (2009) 45–49. [17] M. Sun, D. Zhang, X.F. Gong, A novel approach for description of nonlinear field radiated from a concave source with wide aperture angle, Ultrasonics 44 (2006) e1435–e1438. [18] B.N. Kim, S.W. Yoon, Nonlinear parameter estimation in water-saturated sandy sediment with difference frequency acoustic wave, Ultrasonics 49 (2009) 438– 445. [19] G.V. Norton, R.D. Purrington, The Westervelt equation with viscous attenuation versus a causal propagation operator: A numerical comparison, Journal of Sound and Vibration 327 (2009) 163–172. [20] J.N. Tjotta, S. Tjotta, E. Vefring, Propagation and interaction of two collinear finite amplitude sound beams, Journal of the Acoustical Society of America 88 (1990) 2859–2870.
G. Liu et al. / Ultrasonics 54 (2014) 576–585 [21] Y.S. Lee, M.F. Hamilton, Time domain modeling of pulsed finite amplitude sound beams, Journal of the Acoustical Society of America 97 (1995) 906–917. [22] T. Kamakura, M. Akiyama, K. Aoki, A higher-order parabolic equation for describing nonlinear propagation of ultrasound beams, Acoustic Science and Technology 25 (2004) 163–165. [23] T. Kamakura, Two model equations for describing nonlinear sound beams, Japanese Journal of Applied Physics 43 (2004) 2808–2812. [24] C.C. Pozuelo, B. Dubus, J.A. Gallego-Juarez, Finite element analysis of the nonlinear propagation of high intensity acoustic waves, Journal of the Acoustical Society of America 106 (1999) 91–101. [25] C. Vanhille, C.C. Pozuelo, Numerical model for nonlinear standing waves and weak shocks in thermoviscous fluids, Journal of the Acoustical Society of America 109 (2001) 2660–2667. [26] C. Vanhille, C.C. Pozuelo, A numerical formulation for nonlinear ultrasonic waves propagation in fluids, Ultrasonics 42 (2004) 1123–1128. [27] Q.X. Wang, J.R. Blake, Non-spherical bubble dynamics in a compressible liquid. Part 1. Travelling acoustic wave, Journal of Fluid Mechanics 659 (2010) 191– 224.
585
[28] A. Prosperetti, A. Lezzi, Bubble dynamics in a compressible liquid. Part 1. First order theory, Journal of Fluid Mechanics 168 (1986) 457–478. [29] A. Lezzi, A. Prosperetti, Bubble dynamics in a compressible liquid. Part 2. Second order theory, Journal of Fluid Mechanics 185 (1987) 289– 321. [30] S. Fujikawa, T. Akamatsu, Effects of the non-equilibrium condensation of vapour on the pressure wave produced by the collapse of a bubble in a liquid, Journal of Fluid Mechanics 97 (1980) 481–512. [31] J.P. Franc, J.M. Michel, Fundamentals of Cavitation, Kluwer Academic Publishers, Boston, 2004 (pp. 305). [32] S.W. Gong, S.W. Ohl, E. Klaseboer, B.C. Khoo, Scaling law for bubbles induced by different external sources: Theoretical and experimental study, Physical Review E 81 (2010) 056317. [33] T.J. Matula, P.R. Hilmo, B.D. Storey, A.J. Szeri, Radial response of individual bubbles subjected to shock wave lithotripsy pulses in vitro, Physics of Fluids 14 (2002) 913–921.