Theory

In a vertically heterogeneous, isotropic and elastic medium occupying a half-space, a $P-S_V$ wave travelling along X axis generates displacements along X and Z axis of the form

$$u_x = r_1 (k,z,\omega) e^{i(kx-\omega t)}$$

$$u_y =$$ 0 (3.24)

$$u_z = r_2 (k,z,\omega) e^{i(kx-\omega t)}$$

where $u_x$ , $u_y$ and $u_z$ are the radial, transversal and vertical components, respectively, $r_1$ and $r_2$ are the complex amplitudes (including phase shifts between components), $\omega$ is the angular frequency, and $V_R=\frac{\omega}{k}$ is the velocity of Rayleigh waves (m/s). A motion-stress vector for Rayleigh waves is defined in a similar way as for Love case (section 3.1.3).

$$r(k,z,\omega)=\left(\begin{array}{c}r_1(k,z,\omega) \ r_2(k,z,\omega) \ r^\sigma_1(k,z,\omega) \ r^\sigma_2(k,z,\omega) \end{array}\right)$$ (3.25)

where,

$$r^\sigma_1 = i \left( (\lambda+2\mu) \frac{d r_2}{dz} +k \lambda r_1 \right)$$ (3.26)

$$r^\sigma_2 = \mu \left( \frac{d r_1}{dz}-k r_2 \right)$$

$r^\sigma_1$ is the amplitude of the vertical compression stress $\sigma_{zz}$ and $r^\sigma_2$ is the amplitude of the radial shear stress $\sigma_{xz}$ . From equation of motion (3.9), the solution must satisfy the following differential equation:

$$\frac {d}{dz} \left(\begin{array}{c}r_1 \ r_2\ r^\sigma_1 \ r^... ... \left(\begin{array}{c}r_1 \ r_2\ r^\sigma_1 \ r^\sigma_2 \end{array}\right)$$ (3.27)

where $\rho$ is the density, $\lambda$ and $\mu$ are Lamé moduli. The $z$ dependencies of $\rho$ , $\lambda$ and $\mu$ have been dropped for simplicity. For surface waves, the boundary conditions require that:

$$r_1\rightarrow 0 \quad \textrm{and} \quad r_2\rightarrow 0 \quad \textrm{when} \quad z \rightarrow \infty$$ (3.28)

$$r^\sigma_1=r^\sigma_2=0 \quad \textrm{at the free surface} \quad z=z_0$$ (3.29)

Similarly to Love case, the equation of motion is reduced to an equation of the same form as equation (3.1). The solution for the motion-stress vector is given by equation (3.5). The constraint on the motion-stress at infinity is transformed into a radiation condition that no up-going waves are found in the bottom half-space. For $P-S_V$ plane waves, the amplitudes of downgoing ( $\grave P_n$ and $\grave S_n$ for P and S-waves, respectively) and up going ( $\acute P_n$ and $\acute S_n$ for P and S-waves, respectively) waves traveling across an homogeneous half space are function of the motion-stress vector at the top of the half space ($z=z_n$ ) (Aki and Richards (2002); Dunkin (1965)). The motion-stress vector is propagated to $z_0$ by the means of equation (3.5). The subscript $n$ is added to all parameters defined for layer $n$ .

$$\left(\begin{array}{c} \acute P_n \ \acute S_n \ \grave P_n \ ... ...ma_2(z_n) \end{array}\right)= T_n^{-1} G(z_n,z_{n-1}) \ldots{} G(z_1,z_0)r(z_0)$$ (3.30)

where,

$$T_n^{-1}= \frac{-V_{sn}^2}{2\mu_n \hat h_n \hat k_n \omega^ 2} \l... ..._n k \hat h_n \hat k_n & -i k \hat h_n & - \hat h_n \hat k_n \end{array}\right)$$ (3.31)

where, ${\hat h_n}^2 = k^2-\frac{\omega^2}{V_{pn}^2}$ , $V_{pn}$ is the velocity of P waves (m/s), ${\hat k_n}^2=k^2-\frac{\omega^2}{V_{sn}^2}$ , $V_{sn}$ is the velocity of S waves (m/s), $l_n=k^2+ {\hat k_n}^2$ , $\mu_n=\rho_n V_{sn}^2$ is the rigidity, and $\rho_n$ is the density (t/m$^3$ ).

Merging boundary conditions with equation (3.30),

$$\left(\begin{array}{c} 0 \ 0 \ \grave P_n \ \grave S_n \end{array}\right) = T_n^{-1} G_n \ldots{} G_1 r(z_0)$$

$$= \underbrace{\left(\begin{array}{cccc} r_{11} & r_{12} & r_{13} & ... ...R(z_0)} \left(\begin{array}{c}r_1(z_n) \ r_2(z_n) \ 0 \ 0 \end{array}\right)$$ (3.32)

This equation is always true when the sub-determinant $(r_{11}r_{22}-r_{12}r_{21})$ vanishes. Like in the Love case, the problem of finding the dispersion curves is thus reduced to a root search along the slowness or the velocity axis for a given frequency. However as stated by Dunkin (1965), the terms of the sub-determinant can become very large. Subtracting two large numbers results in a loss of significant digits, which implies the use of very high precision computations (e.g. 128 bit numbers or even more whereas computers are classically limited to 32 or 64 bits). Hence, Dunkin proposed an alternative way of propagating motion-stress vector by the means of the following theorem. If $P=A^{(0)} A^{(1)} \ldots{} A^{(n-1)} A^{(n)}$ then

$$p \left\vert \begin{array}{cc}i & j \ k & l\end{array}\right\ver... ...ht\vert a^{(n)}\left\vert \begin{array}{cc}u & v \ k & l\end{array}\right\vert$$ (3.33)

where $p \left\vert \begin{array}{cc}i & j \ k & l\end{array}\right\vert = p_{ik}p_{jl}-p_{il}p_{jk}$ is the second order sub-determinant of matrix P. The notation $P_{ijkl}$ is also used in appendix A.

In equation (3.33), the summation rules^3.1 apply for indices appearing two times like $m$ and $n$ . In this case, the summed pairs of indices are to be only distinct pairs of distinct indices^3.2 (by convention, $m<n$ , $o<p$ ,..., $s<t$ , $u<v$ ). It follows from equation (3.33) that:

$$r_{11}(z_0) r_{22}(z_0) -r_{12}(z_0) r_{21}(z_0) = r(z_0) \left\vert \begin{array}{cc}1 & 2 \ 1 & 2\end{array}\right\vert$$

$$= t_n^-1 \left\vert \begin{array}{cc}1 & 2 \ a & b\end{array}\righ... ...\ldots{} g_1 \left\vert \begin{array}{cc}e & f \ 1 & 2\end{array}\right\vert=0$$ (3.34)

With the condition on indices, the factor $t_n^-1 \left\vert \begin{array}{cc}1 & 2 \ a & b\end{array}\right\vert$ has 6 components: 12, 13, 14, 23, 24 and 34 (Dunkin (1965)). On the other hand, $g_n \left\vert \begin{array}{cc}a & b \ c & d\end{array}\right\vert$ has 6x6 components. Hence, like in the Love case, for a given frequency ( $\frac{\omega}{2\pi}$ ), only a few discrete values are possible for the velocity of the Rayleigh surface wave ( $V_R=\left[\frac{\omega}{k(\omega)}\right]_i$ ), corresponding to the dispersion curves of various modes. The dispersion curve is found by seeking the roots of $r(z_0) \left\vert \begin{array}{cc}1 & 2 \ 1 & 2\end{array}\right\vert$ .