A conservative numerical scheme for modeling nonlinear acoustic propagation in thermoviscous homogeneous media

doi:10.1016/j.jcp.2018.02.005

Journal of Computational Physics

Volume 363, 15 June 2018, Pages 200-230

https://doi.org/10.1016/j.jcp.2018.02.005 Get rights and content

Highlights

•
A conservative nonlinear thermoviscous acoustic system model is proposed.
•
Strongly nonlinear acoustic waves are studied with high-order WENO schemes. A focused ultrasound device has been studied to physically validate the model.
•
Simulations using GPU accelerators were performed to illustrate the method in 2-d.

Abstract

A nonlinear system of partial differential equations capable of describing the nonlinear propagation and attenuation of finite amplitude perturbations in thermoviscous media is presented. This system constitutes a full nonlinear wave model that has been formulated in the conservation form. Initially, this model is investigated analytically in the inviscid limit where it has been found that the resulting flux function fulfills the Lax–Wendroff theorem, and the scheme can match the solutions of the Westervelt and Burgers equations numerically. Here, high-order numerical descriptions of strongly nonlinear wave propagations become of great interest. For that matter we consider finite difference formulations of the weighted essentially non-oscillatory (WENO) schemes associated with explicit strong stability preserving Runge–Kutta (SSP-RK) time integration methods. Although this strategy is known to be computationally demanding, it is found to be effective when implemented to be solved in graphical processing units (GPUs). As we consider wave propagations in unbounded domains, perfectly matching layers (PML) have been also considered in this work. The proposed system model is validated and illustrated by using one- and two-dimensional benchmark test cases proposed in the literature for nonlinear acoustic propagation in homogeneous thermoviscous media.

Introduction

The numerical description of high amplitude ultrasound propagation is a challenging and computationally intense task [1] for two main reasons. Specifically, the large spatial domain size typically required for acoustic applications and the indispensable need of high-order schemes to minimize the amount of numerical dissipation introduced by the numerical scheme.

The numerical modeling of nonlinear acoustics propagation plays an important role in the development of biomedical and engineering applications, for instance, in high intensity focused ultrasound (HIFU) for thermal ablation treatment [2], shock wave lithotripsy [3] and diagnostic ultrasound imaging [4]. For these applications the linear model is inadequate as the nonlinear effect in acoustic propagation plays the major role.

In the literature, many acoustic models have been proposed (see the recent review [5]). Arguably, the two most popular approaches for the direct modeling of nonlinear acoustic propagation are the (single sided) wave model of Khokhlov–Zabolotskaya–Kuznetsov (KZK) and the full (double sided) wave model of Westervelt. Likewise, the most notorious numerical approaches include the direct discretization with a finite difference time domain (FDTD) scheme [6], [7], [8], the dispersion relation preserving (DPR) method [9] and the k-space pseudo-spectral methods [10]. More recently, Albin and collaborators [11], [12] have presented viscous acoustic descriptions based on the Navier–Stokes (N–S) equations using the Fourier Continuation (FC) method. The FC method is an explicit spectral scheme for modeling the transient behavior of nonlinear hyperbolic conservation systems in long-term propagations at a reasonable cost.

All the previously mentioned schemes are shown to perform well as long as the solutions are smooth. However, in the presence of discontinuous solutions, these methods need to rely on numerical artifacts such as frequency filters to damp higher-order harmonics and avoid polluting the solution with Gibbs oscillations. By using those artifacts, the above numerical approaches increase their computational cost and/or make their implementations more involving.

Arguably, the most successful shock-capturing methods are those developed for systems of hyperbolic equations, more specifically, for conservation laws [13], [14]. This realization is clearly reflected in the work of Christov and collaborators [15], [16]. In their study of conservative acoustic formulations, shock structures were obtained for many one-dimensional weakly-nonlinear acoustic models using MUSCL schemes.

The acoustic explorations with second-order acoustic wave models and the full N–S system in the literature have led us to consider that an acoustic model formulated as a hyperbolic system of equations can yield simpler and artifact-free numerical descriptions of linear and nonlinear acoustics with shocks propagations in homogeneous and heterogeneous thermoviscous media. This work represents an initial step towards this goal. Like the classical thermoviscous wave models, our exploration concentrates on hyperbolic formulations to describe directly acoustic propagations as pressure wave fluctuations. For that purpose, we start by considering the finite amplitude model of Aanonsen [17] called second-order wave equation in [18]. Our investigation revealed that the second-order wave equation can lead to many suitable quasilinear hyperbolic systems. Although this was expected, it was not foreseen it could lead to a fully nonlinear formulation in conservation form. In the present work, this conservative formulation is numerically investigated for acoustic propagation in homogeneous media. In order to describe strongly nonlinear propagation, the methodology considered is based on the fifth- and seventh-order finite difference WENO reconstructions for the spatial evolution and is associated with low-storage Runge–Kutta schemes for the time integration. Acoustic wave propagation in open domains is only considered, therefore the implementation of perfect matching layers (PML) has been also taken into account to avoid unphysical reflections from the boundary of our discrete domain.

This work is organized as follows. The second-order thermoviscous acoustic wave models are reviewed in § 2. In § 3, a thermoviscous acoustic formulation as a conservative hyperbolic system of equations and its simplifications are presented. The numerical implementation of the model in homogeneous media is introduced for single and multiple dimensions in § 4. In § 5, the proposed models are validated and illustrated using one- and two-dimensional benchmark tests proposed in the literature. Lastly, concluding remarks and extensions of this work are addressed in § 6.

Section snippets

Nonlinear acoustic models

A model to describe the propagation of acoustics fluctuations through homogeneous media, under the assumptions of Newtonian, compressible, irrotational, heat conducting and viscous flow, might be obtained from the equations of fluid dynamics, namely: the equation of mass conservation $\frac{\partial ρ}{\partial t} + \nabla \cdot (ρ υ) = 0,$ the equation of motion $ρ (\frac{\partial υ}{\partial t} + (υ \cdot \nabla) υ) = - \nabla p + μ \nabla^{2} υ + (μ_{B} + \frac{1}{3} μ) \nabla (\nabla \cdot υ),$ the equation of heat transfer $ρ T (\frac{\partial s}{\partial t} + (υ \cdot \nabla) s) = \frac{μ}{2} {(\frac{\partial υ_{i}}{\partial x_{j}} + \frac{\partial υ_{j}}{\partial x_{i}} - \frac{2}{3} δ_{i j} \frac{\partial υ_{k}}{\partial x_{k}})}^{2} + μ_{B} {(\nabla \cdot υ)}^{2} + κ \nabla^{2} T,$ and the equation of state (EoS) $p = p (ρ, s),$ where $υ$

Thermoviscous hyperbolic acoustic system

Let us start our development by considering equations in (12a), (12b) as an intermediate step. The temporal derivative of the acoustic density $\frac{\partial ρ_{a}}{\partial t}$ using the EoS can be written in the following form [18] $\frac{\partial ρ_{a}}{\partial t} = \frac{1}{c_{0}^{2}} \frac{\partial p_{a}}{\partial t} - \frac{B / A}{2 ρ_{0} c_{0}^{4}} \frac{\partial p_{a}^{2}}{\partial t} - \frac{κ}{ρ_{0} c_{0}^{2}} (\frac{1}{c_{v}} - \frac{1}{c_{p}}) \nabla^{2} p_{a},$ so that substituting this last result in (12a), (12b) yields $\begin{matrix} ρ_{0} \frac{\partial υ_{a}}{\partial t} + \nabla p_{a} = & (μ_{B} + \frac{4}{3} μ) \nabla^{2} υ_{a} - \nabla L \\ \frac{1}{c_{0}^{2}} \frac{\partial p_{a}}{\partial t} + ρ_{0} \nabla \cdot υ_{a} = & \frac{β}{ρ_{0} c_{0}^{4}} \frac{\partial p_{a}^{2}}{\partial t} + \frac{κ}{ρ_{0} c_{0}^{2}} (\frac{1}{c_{v}} - \frac{1}{c_{p}}) \frac{\partial^{2} p_{a}}{\partial t^{2}} + \frac{1}{c_{0}^{2}} \frac{\partial L}{\partial t} . \end{matrix}$ To write this system in a hyperbolic system, it is necessary to write the nonlinear β-term

Numerical implementation

Describing numerically acoustic propagations with the NAS equations requires a high-resolution solver capable of dealing with discontinuous solutions. Moreover, as we target strongly nonlinear propagation, the number of points per wavelength required can yield a large computational domain. To obtain solutions in reasonable time, an explicit algorithm implemented using a high-performance computational (HPC) solution is desirable.

Given the numerical challenges posed by the NAS equations, a finite

Numerical tests

The following examples are all related to the numerical simulations of nonlinear acoustic propagations that result from a high intensity focused ultrasound (HIFU) test. This test can be performed as a two-dimensional problem with axial symmetry since the domain is assumed to be homogeneous. Here, the essential building blocks, namely a one- and two-dimensional Cartesian and axisymmetric formulation of the NAS equations, and extensions of this formulations were presented in § 4.

Here, we are

Concluding remarks

In this work, a conservative thermoviscous formulation to describe nonlinear acoustic propagation has been obtained. Moreover, it is found that at the core of this formulation there is a fully nonlinear hyperbolic system, denominated as the nonlinear acoustic system of equations.

Like the traditional second-order acoustic formulations, the present formulation solves directly for the acoustic pressure field. However, the present hyperbolic formulation also requires solving for the velocity

Acknowledgements

This research was supported by the Ministry of Science and Technology of Taiwan, R.O.C., under the grant MOST-105-2221-E-400-005 and by National Health Research Institutes's project BN-106-PP-08.

References (36)

P. Jordan
A survey of weakly-nonlinear acoustic models: 1910–2009
Mech. Res. Commun.
(2016)
N. Albin et al.
A spectral FC solver for the compressible Navier–Stokes equations in general domains, I: Explicit time-stepping
J. Comput. Phys.
(2011)
I. Christov et al.
Nonlinear acoustic propagation in homentropic perfect gases: a numerical study
Phys. Lett. A
(2006)
J.M. Burgers
A mathematical model illustrating the theory of turbulence
Adv. Appl. Mech.
(1948)
J.-P. Berenger
A perfectly matched layer for the absorption of electromagnetic waves
J. Comput. Phys.
(1994)
N.A. Petersson et al.
Discretizing singular point sources in hyperbolic wave propagation problems
J. Comput. Phys.
(2016)
C.-W. Shu et al.
Efficient implementation of essentially non-oscillatory shock-capturing schemes
J. Comput. Phys.
(1988)
G.-S. Jiang et al.
Efficient implementation of weighted ENO schemes
J. Comput. Phys.
(1996)
D.S. Balsara et al.
Monotonicity preserving weighted essentially non-oscillatory schemes with increasingly high order of accuracy
J. Comput. Phys.
(2000)
M.A. Diaz et al.
An efficient direct solver for rarefied gas flows with arbitrary statistics
J. Comput. Phys.
(2016)

M. Solovchuk et al.

Multiphysics modeling of liver tumor ablation by high intensity focused ultrasound

Commun. Comput. Phys.

(2015)

J. Kennedy et al.

High intensity focused ultrasound: surgery of the future?

Br. J. Radiol.

(2003)

J.E. Lingeman et al.

Shock wave lithotripsy: advances in technology and technique

Nat. Rev. Urol.

(2009)

T.L. Szabo

Diagnostic Ultrasound Imaging: Inside Out

(2004)

I.M. Hallaj et al.

FDTD simulation of finite-amplitude pressure and temperature fields for biomedical ultrasound

J. Acoust. Soc. Am.

(1999)

A. Karamalis et al.

Fast ultrasound image simulation using the Westervelt equation

A.A. Haigh et al.

Ultrasound simulation on the cell broadband engine using the Westervelt equation

S. Ginter et al.

Full-wave modeling of therapeutic ultrasound: nonlinear ultrasound propagation in ideal fluids

J. Acoust. Soc. Am.

(2002)

Cited by (0)

View full text

A conservative numerical scheme for modeling nonlinear acoustic propagation in thermoviscous homogeneous media

Highlights

Abstract

Introduction

Section snippets

Nonlinear acoustic models

Thermoviscous hyperbolic acoustic system

Numerical implementation

Numerical tests

Concluding remarks

Acknowledgements

Mech. Res. Commun.

J. Comput. Phys.

Phys. Lett. A

Adv. Appl. Mech.

J. Comput. Phys.

J. Comput. Phys.

J. Comput. Phys.

J. Comput. Phys.

J. Comput. Phys.

J. Comput. Phys.

Multiphysics modeling of liver tumor ablation by high intensity focused ultrasound

Commun. Comput. Phys.

High intensity focused ultrasound: surgery of the future?

Br. J. Radiol.

Shock wave lithotripsy: advances in technology and technique

Nat. Rev. Urol.

Diagnostic Ultrasound Imaging: Inside Out

FDTD simulation of finite-amplitude pressure and temperature fields for biomedical ultrasound

J. Acoust. Soc. Am.

Fast ultrasound image simulation using the Westervelt equation

Ultrasound simulation on the cell broadband engine using the Westervelt equation

Full-wave modeling of therapeutic ultrasound: nonlinear ultrasound propagation in ideal fluids

J. Acoust. Soc. Am.