Simulating propagation of separated wave modes in general anisotropic media, Part I: qP-wave propagators

Plane-wave analysis of the elastic wave equation

Vector and component notations are used alternatively throughout the paper. The wave equation in general heterogeneous anisotropic media can be expressed as (Carcione, 2001),

$$\rho\frac{\partial^2\mathbf{u}}{\partial t^2} = [{\bigtriangledown}{\mathbf{C}{\bigtriangledown}^{T}}]\mathbf{u} + \mathbf{f},$$ (1)

where $\mathbf{u}=(u_x,u_y,u_z)^{T}$ is the particle displacement vector, $\mathbf{f}=(f_x,f_y,f_z)^{T}$ represents the force term, $\rho$ is the density, $\mathbf{C}$ is the matrix representing the stiffness tensor in a two-index notation called the “Voigt recipe”, and the symmetric gradient operator has the following matrix representation:

$$\tensor{\bigtriangledown} = \begin{pmatrix}\frac{\partial}{\parti... ...} & \frac{\partial}{\partial y} & \frac{\partial}{\partial x} &0 \end{pmatrix}.$$ (2)

Assuming that the material properties vary sufficiently slowly so that spatial derivatives of the stiffnesses can be ignored, equation 1 can be simplified as

$$\rho\frac{\partial^2\mathbf{u}}{\partial t^2} = \Gamma\mathbf{u} + \mathbf{f},$$ (3)

where $\Gamma$ is the $3\times3$ symmetric Christoffel differential-operator matrix, of which the elements are given for locally smooth media as follows (Auld, 1973),

$$\begin{split}\Gamma_{11} &= C_{11}\frac{\partial^2}{\partial ... ...C_{66})\frac{\partial^2}{{\partial x}{\partial y}}. \end{split}$$ (4)

For the most important types of seismic anisotropy such as transverse isotropy and orthorhombic anisotropy, some terms in equation 4 vanish because the corresponding stiffness coefficients become zeros.

Neglecting the source term, a plane-wave analysis of the elastic anisotropic wave equation yields the Christoffel equation,

$$\widetilde{\mathbf{\Gamma}}\widetilde{\mathbf{u}} = \rho{\omega}^2\widetilde{\mathbf{u}},$$ (5)

or

$$(\widetilde{\mathbf{\Gamma}} - \rho{\omega}^2\mathbf{I})\widetilde{\mathbf{u}} = \mathbf{0},$$ (6)

where $\omega$ is the frequency, $\widetilde{\mathbf{u}}=(\widetilde{u}_x,\widetilde{u}_y,\widetilde{u}_z)^{T}$ is the wavefield in Fourier domain, $\widetilde{\Gamma} = \widetilde{\mathbf{L}}\mathbf{C}\widetilde{\mathbf{L}}^{T}$ is the symmetric Christoffel matrix, $\mathbf{I}$ is a $3\times3$ identity matrix. To support the sign notation in equations 5 and 6, we remove the imaginary unit $i$ of the wavenumber-domain counterpart of the gradient operator $\tensor{\bigtriangledown}$ and thus express matrix $\widetilde{\mathbf{L}}$ as:

$$\widetilde{\mathbf{L}}= \begin{pmatrix}k_x & 0 &0 &0 & k_z & k_y \cr 0 & k_y & 0 & k_z &0 & k_x \cr 0 & 0 & k_z & k_y & k_x &0\end{pmatrix}.$$ (7)

Setting the determinant of $\widetilde{\mathbf{\Gamma}} - \rho{\omega}^2\mathbf{I}$ in equation 6 to zero gives the characteristic equation, and expanding that determinant gives the (angular) dispersion relation. For a given spatial direction specified by a wave vector $\mathbf{k} = (k_x, k_y, k_z)^{T}$ , the characteristic equation poses a standard $3\times3$ eigenvalue problem. The three eigenvalues correspond to the phase velocities of the qP-wave and two qS waves. Inserting one of the eigenvalues back into the Christoffel equation gives ratios of the components of $\mathbf{\widetilde{u}}$ , from which the polarization or displacement direction can be determined for the given wave mode. In general, these directions are neither parallel nor perpendicular to the wave vector, and depend on the local material parameters for the anisotropic medium. For a given wave vector or slowness direction, the polarization vectors of the three wave modes are always mutually orthogonal.

Applying an inverse Fourier transform to the dispersion relation yields a high-order PDE in time and space and contains mixed space and time derivatives. Setting the shear velocity along the axis of symmetry to zero while using Thomsen's parameter notation yields the pseudo-acoustic dispersion relation and wave equation in VTI media (Alkhalifah, 2000). Most published methods instead have used coupled PDEs (derived from the pseudo-acoustic dispersion relation) that are only second-order in time and eliminate the mixed space-time derivatives, e.g., Zhou et al. (2006b). Many kinematically equivalent coupled second-order systems can be generated from the dispersion relation by similarity transformations (Fowler et al., 2010). In the next section, we present a particular similarity transformation to the Christoffel equation in order to derive a minimal second-order coupled system, which is helpful for simulating propagation of separated qP-waves in anisotropic media.

Simulating propagation of separated wave modes in general anisotropic media, Part I: qP-wave propagators