Seismic Anisotropy in the Rift of the Reykjanes Peninsula, SW Iceland, Calculated Using a New Tomographic Method

Abstract

A new algorithm and computer code for seismic tomography in anisotropic inhomogeneous media was developed. The new tomographic approach is a generalization of classical isotropic seismic tomography which introduces spatially and directionally varying slowness. The velocity model was considered as a stack of homogeneous blocks in contact, with each block individually parameterized using background velocities of P- and S-waves and a set of 21 anisotropy parameters. The inverse problem was solved sequentially in five steps, using the velocity model from the previous step as the starting model for the subsequent step. These steps form a chain with increasing complexity: (1) isotropic homogeneous model; (2) isotropic velocity model with vertical velocity gradient; (3) 3-D inhomogeneous isotropic velocity model; (4) 3-D inhomogeneous model with uniform anisotropy; (5) 3-D inhomogeneous generally anisotropic model. The new algorithm was applied to real bulletin data of 18 seismic stations deployed in SW Iceland and operated favourably for the monitoring of local swarm-like seismicity. Next, the resolution, robustness and accuracy of the inversion were discussed using real and synthetic data. Real data inversion revealed a predominantly depth-dependent isotropic velocity background and additional general 3-D anisotropy. The parameterization of the medium was too flexible to allow for a reliable interpretation of the anisotropy inside the elementary blocks and a cluster analysis was applied to stabilize the inversion results. Three important clusters were identified as a result. The orientation of the anisotropy (fast and slow P-wave propagation directions) of two clusters coincided with the strike of the documented faults. The orientation of the anisotropy in the third cluster was interpreted as a consequence of the fluid dynamics around Kleifarvatn Lake. The P-wave anisotropy strength reached a value of $\pm$ 5–8%.

Appendices

Appendix A. Relation Between Anisotropy Parameters and Elements of the Voigt Matrix

The anisotropy parameters, denoted with Greek symbols, correspond to elements A${_{ij}}$ of the density normalized elastic moduli expressed in Voigt notation as follows:

$$\begin{aligned}&\epsilon _x=\frac{A_{11}-\alpha _0^2}{2\alpha _0^2},\nonumber \\&\quad \epsilon _y=\frac{A_{22}-\alpha _0^2}{2\alpha _0^2},\nonumber \\&\quad \epsilon _z=\frac{A_{33}-\alpha _0^2}{2\alpha _0^2},\nonumber \\&\quad \chi _x=\frac{A_{14}+2A_{56}}{\alpha _0^2},\nonumber \\&\quad \chi _y=\frac{A_{25}+2A_{46}}{\alpha _0^2},\nonumber \\&\quad \chi _z=\frac{A_{36}+2A_{45}}{\alpha _0^2},\nonumber \\&\quad \eta _x=\frac{2(A_{23}+2A_{44})-A_{22}-A_{33}}{2\alpha _0^2},\nonumber \\&\quad \eta _y=\frac{2(A_{13}+2A_{55})-A_{33}-A_{11}}{2\alpha _0^2},\nonumber \\&\quad \eta _z=\frac{2(A_{12}+2A_{66})-A_{11}-A_{22}}{2\alpha _0^2},\nonumber \\&\quad \xi _{15}=\frac{A_{25}+2A_{46}-A_{15}}{\alpha _0^2},\nonumber \\&\quad \xi _{16}=\frac{A_{36}+2A_{45}-A_{16}}{\alpha _0^2},\nonumber \\&\quad \xi _{24}=\frac{A_{14}+2A_{56}-A_{24}}{\alpha _0^2},\nonumber \\&\quad \xi _{26}=\frac{A_{36}+2A_{45}-A_{26}}{\alpha _0^2},\nonumber \\&\quad \xi _{34}=\frac{A_{14}+2A_{56}-A_{34}}{\alpha _0^2},\nonumber \\&\quad \xi _{35}=\frac{A_{25}+2A_{46}-A_{35}}{\alpha _0^2},\nonumber \\&\quad \gamma _x=\frac{A_{44}-\beta ^2}{2\beta _0^2},\nonumber \\&\quad \gamma _y=\frac{A_{55}-\beta ^2}{2\beta _0^2},\nonumber \\&\quad \gamma _z=\frac{A_{66}-\beta ^2}{2\beta _0^2},\nonumber \\&\quad \epsilon _{45}=\frac{A_{45}}{\beta _0^2},\nonumber \\&\quad \epsilon _{46}=\frac{A_{46}}{\beta _0^2},\nonumber \\&\quad \epsilon _{56}=\frac{A_{56}}{\beta _0^2} \end{aligned}$$

(A-1)

Note that there is some leeway in the selection of background velocities $\alpha$₀ and $\beta$₀. The transform (A-1) is therefore not unique. The fixed matrix A and different selections of background velocities generate different anisotropy parameters. However, the direction-dependent squared-velocity formula shown below is independent of the selection of $\alpha$₀ and $\beta$₀ for the fixed Voigt matrix A and a direction-dependent slowness formula should be calculated for such a combination of $\alpha$₀ and $\beta$₀ which produces a minimum norm of m (see below). If $\alpha$₀ and $\beta$₀ are fixed in this manner, both the transform (A-1) and inverse are unique.

Appendix B. Transformation of Anisotropy Parameters to Arbitrary Values of the Reference Isotropic Velocity

Let us reformulate the mapping in Eq. (A-1) between anisotropy parameters and elements of the Voigt matrix using vector/matrix notation. First the two parts of the vector m depending on $\alpha$₀ and $\beta$₀ will be separated:

$$\begin{aligned}&{\varvec{m}}=\left[ \begin{matrix}{\varvec{m}}_{\alpha }\\ {\varvec{m}}_{\beta }\end{matrix}\right] \nonumber \\&\quad {\varvec{m}}_{\alpha }=\left( \epsilon _x,\epsilon _y,\epsilon _z,\chi _x,\chi _y,\chi _z,\eta _x,\eta _y,\eta _z,\xi _{15,}\xi _{16,}\xi _{24,}\xi _{26,}\xi _{34,}\xi _{35,}\right) ^T\nonumber \\&\quad \varvec{m}_{\beta }=\left( \gamma _x,\gamma _e,\gamma _z,\epsilon _{45,}\epsilon _{46,}\epsilon _{56}\right) ^T. \end{aligned}$$

(A-2)

The transformation between the two parts ${\mathbf{m}}_\alpha$ and ${\mathbf{m}}_\beta$ and Voigt matrix elements can then be expressed as:

$$\begin{aligned} {\varvec{m}}_{\alpha }={\varvec{b}}_{\alpha }+\frac{1}{\alpha _0^2}\varvec{Q_{\alpha }}{\varvec{A}}_{\alpha },{\varvec{m}}_{\beta }={\varvec{b}}_{\beta }+\frac{1}{\beta _0^2}\varvec{Q_{\beta }}{\varvec{A}}_{\beta }. \end{aligned}$$

(A-3)

where

$$\begin{aligned}&\varvec{b}_{\alpha }= \left[ \begin{array}{lllllllllllllll} -\frac{1}{2}&-\frac{1}{2}&-\frac{1}{2}&0&0&0&0&0&0&0&0&0&0&0&0 \end{array} \right] ^T \nonumber \\&\varvec{Q}_{\alpha }= \left[ \begin{array}{lllllllllllllll} \frac{1}{2} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 \\ 0 &{} \frac{1}{2} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 \\ 0 &{} 0 &{} \frac{1}{2} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 \\ 0 &{} 0 &{} 0 &{} 1 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 \\ 0 &{} 0 &{} 0 &{} 0 &{} 1 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 \\ 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 1 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 \\ 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \frac{1}{2} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 \\ 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \frac{1}{2} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 \\ 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} \frac{1}{2} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 \\ 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 1 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 \\ 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 1 &{} 0 &{} 0 &{} 0 &{} 0 \\ 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 1 &{} 0 &{} 0 &{} 0 \\ 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 1 &{} 0 &{} 0 \\ 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 1 &{} 0 \\ 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 1 \end{array} \right] \nonumber \\&\varvec{A}_{\alpha }= \left[ \begin{array}{l} A_{11} \\ A_{22} \\ A_{33} \\ A_{14}+2A_{56} \\ A_{25}+2A_{46} \\ A_{36}+2A_{45} \\ 2(A_{23}+2A_{44})-A_{22}-A_{33} \\ 2(A_{13}+2A_{455})-A_{33}-A_{11} \\ 2(A_{12}+2A_{66})-A_{11}-A_{22} \\ A_{25}+2A_{46}-A_{15} \\ A_{36}+2A_{45}-A_{16} \\ A_{14}+2A_{56}-A_{24} \\ A_{36}+2A_{45}-A_{26} \\ A_{14}+2A_{56}-A_{34} \\ A_{25}+2A_{46}-A_{35} \end{array} \right] \end{aligned}$$

(A-4a)

as for ${\mathbf{m}}_\alpha$ and

$$\begin{aligned}&{\varvec{b}}_{\beta }= \left[ \begin{array}{llllll} -\frac{1}{2}&-\frac{1}{2}&-\frac{1}{2}&0&0&0 \end{array} \right] ^T\nonumber \\&{\varvec{Q}}_{\beta }= \left[ \begin{array}{llllll} \frac{1}{2} &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 \\ 0 &{} \frac{1}{2} &{} 0 &{} 0 &{} 0 &{} 0 \\ 0 &{} 0 &{} \frac{1}{2} &{} 0 &{} 0 &{} 0 \\ 0 &{} 0 &{} 0 &{} 1 &{} 0 &{} 0 \\ 0 &{} 0 &{} 0 &{} 0 &{} 1 &{} 0 \\ 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 1 \end{array} \right] ,\nonumber \\&{\varvec{A}}_{\beta }= \left[ \begin{array}{l} A_{44} \\ A_{55} \\ A_{66} \\ A_{45} \\ A_{46} \\ A_{56} \end{array} \right] . \end{aligned}$$

(A-4b)

as for ${\mathbf{m}}_\beta$ .

Let us first discuss the transformation of the {${\mathbf{m}}_\alpha$ , $\alpha$₀} part of m to another value of reference velocity $\alpha$₀. The index $\alpha$ in ${\mathbf{m}}_\alpha$ will be omitted from now on in order to simplify indexing. We shall consider the pair {m₀, $\alpha$₀} $\Leftrightarrow$ {m₁, $\alpha$₁} as two ${\mathbf{m}}_\alpha$ representations of the same Voigt matrix. Using Eq. (A-3) gives

$$\begin{aligned}{\varvec{m}}_0&={\varvec{b}}_{\alpha }+\frac{1}{\alpha _0^2}{\varvec{Q}}_{\alpha }{\varvec{A}}_{\alpha }\rightarrow {\varvec{A}}_{\alpha }=\alpha _0^2{\varvec{Q}}_{\alpha }^{-1}({\varvec{m}}_0-{\varvec{b}}_{\alpha })\nonumber \\ {\varvec{m}}_1&={\varvec{b}}_{\alpha }+\frac{1}{\alpha _1^2}{\varvec{Q}}_{\alpha }{\varvec{A}}_{\alpha }={\varvec{b}}_{\alpha } +\frac{\alpha _0^2}{\alpha _1^2}{\varvec{Q}}_{\alpha }{\varvec{Q}}_{\alpha }^{-1}({\varvec{m}}_0-{\varvec{b}}_{\alpha }) \\&={\varvec{b}}_{\alpha }+\frac{\alpha _0 ^2}{\alpha _1^2}({\varvec{m}}_0-{\varvec{b}}_{\alpha }). \end{aligned}$$

(A-5)

Once we have already {m₀, $\alpha$₀} and want to use other reference velocity $\alpha$₁, Eq. (A-5) gives us a rule how to calculate corresponding m₁. We may select $\alpha$₁ arbitrarily, suitable choice is such $\alpha$₁ which makes the norm of ${\mathbf{m}}_\alpha$ minimum, i.e. giving the most isotropic approximation of the actual anisotropic medium. Such value of $\alpha$${_{min}}$ corresponding to m${_{min}}$ can be determined as follows:

$$\begin{aligned}&{\varvec{m}}_{\textit{min}}= {\varvec{b}}_{\alpha }+\frac{\alpha _0^2}{\alpha _{\textit{min}}^2}({\varvec{m}}_0-{\varvec{b}}_{\alpha }), {\varvec{m}}_{\textit{min}}^T\varvec{m}_{\textit{min}}= \textit{min}(\alpha _{\textit{min}})\quad {\varvec{m}}_{\textit{min}}^T\varvec{m}_{\textit{min}}=\left( \frac{\alpha _0^2}{\alpha _{\textit{min}}^2}\right) ^2\left( {\varvec{b}}_{\alpha }^T{\varvec{b}}_{\alpha }+{\varvec{m}}_0^T{\varvec{m}}_0-2{\varvec{m}}_0^T{\varvec{b}}_{\alpha }\right)+\frac{\alpha _0^2}{\alpha _{\textit{min}}^2}\left( 2{\varvec{b}}_{\alpha }^T\varvec{m}_0-2{\varvec{b}}_{\alpha }^T{\varvec{b}}_{\alpha }\right) +{\varvec{b}}_{\alpha }^T\varvec{b}_{\alpha }=\Phi (\alpha _{\textit{min}})\quad \frac{d\Phi }{d\alpha _{\textit{min}}}=-4\frac{\alpha _0^2}{\alpha _{\textit{min}}^3}({\varvec{b}}_{\alpha }^T{\varvec{m}}_0-{\varvec{b}}_{\alpha }^T\varvec{b}_{\alpha } +\frac{\alpha _0^2}{\alpha _{\textit{min}}^2}({\varvec{b}}_{\alpha }^T{\varvec{b}}_{\alpha }\nonumber \\&+{\varvec{m}}_0^T{\varvec{m}}_0-2{\varvec{m}}_0^T{\varvec{b}}_{\alpha }))=0. \end{aligned}$$

(A-6)

Solving the last equation for $\alpha$${_{min}}$ results in

$$\begin{aligned} \alpha _{\textit{min}}^2=\alpha _0^2\frac{{\varvec{b}}_{\alpha }^T{\varvec{b}}_{\alpha }+\varvec{m}_0^T{\varvec{m}}_0-2{\varvec{m}}_0^T{\varvec{b}}_{\alpha }}{{\varvec{b}}_{\alpha }^T{\varvec{b}}_{\alpha }-{\varvec{b}}_{\alpha }^T{\varvec{m}}_0} \text {.} \end{aligned}$$

(A-7)

Once we have any specification of ${\mathbf{m}}_\alpha$ given as {m₀, $\alpha$₀}, we can calculate $\alpha$_min using Eq. (A-7) and recalculate m${_{min}}$ using Eq. (A-6). The calculation of the optimum $\beta$${_{min}}$ and m${_{min}}$ for the ${\mathbf{m}}_\beta$ -part of m is analogous.

Appendix C. Specification of direction-defining vectors in velocity formulas

Let us review the formulae for squared direction-dependent phase velocities $\alpha$² and $\beta$²:

$$\begin{aligned} \alpha ^2=\alpha _0^2+\alpha _0^2{\varvec{g}}^P{\varvec{m}},\beta ^2=\beta _0^2+(\beta _0^2{\varvec{g}}^{\textit{SA}}-\alpha _0^2{\varvec{g}}^{\textit{SB}}){\varvec{m}} . \end{aligned}$$

(A-8)

Separating the corresponding terms at individual elements of the vector m, one produces:

$$\begin{aligned}&{\varvec{g}}^P=2[n_1^{2,}n_2^{2,}n_3^{2,}\nonumber \\&\quad 2(n_2n_{3,}n_3n_{1,}n_1n_2),n_2^2n_3^{2,}n_1^2n_3^{2,}n_1^2n_2^{2,}\nonumber \\&\quad -2(n_1^3n_{3,}n_1^3n_{2,}n_2^3n_{3,}n_2^3n_{1,}n_3^3n_{2,}n_3^3n_1),\nonumber \\&\quad 0,0,0,0,0,0]\nonumber \\&\quad \varvec{g}^{\textit{SA}}=[ 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,n_2^2+n_3^{2,}n_1^2\nonumber \\&\quad +n_3^{2,}n_1^2+n_2^{2,}n_1n_{2,}n_1n_{3,}n_2n_3]\nonumber \\&\quad \varvec{g}^{\textit{SB}}=[ 0,0,0,0,0,0,n_2^2n_3^{2,}n_1^2n_3^{2,}n_1^2n_2^{2,}\nonumber \\&\quad n_1n_3(1-2n_1^2),n_1n_2(1-2n_1^2),n_2n_3(1-2n_2^2),\nonumber \\&n_1n_2(1-2n_2^2),n_2n_3(1-2n_3^2),n_1n_3(1-2n_3^2),\nonumber \\&\quad 0,0,0,0,0,0] \end{aligned}$$

(A-9)

Note that in case of P-waves, only the first 15 anisotropy parameters are important, while in the case of S-waves, only the last 15 anisotropy parameters play a role. The direction vector n = (n₁, n₂, n₃)^T.

Appendix D. Representative density-normalized stiffness tensors corresponding to clusters A-C in the Voigt matrix notation

The stiffness density-normalized matrices expressed in km²/s² units are evaluated for the background velocities $\alpha$₀ and $\beta$₀ given as the mean background velocities for all boxes connected to a particular cluster. The coordinate system has an x-axis oriented from west to east, y axis from south to north and z axis upwards. Equation (A-10) is for cluster A (red boxes in Fig. 15), Eq. (A-11) is for cluster B (blue boxes) and Eq. (A-12) is for cluster C (green boxes).

$$\begin{aligned} {\varvec{A}}_A= \left[ \begin{matrix} 32.29 &{} 11.91 &{} 5.22 &{} -2.89 &{} 1.85 &{} 2.33 \\ 11.91 &{} 29.26 &{} 10.09 &{} -1.93 &{} 0.52 &{} 0.90 \\ 5.22 &{} 10.09 &{} 27.64 &{} -0.94 &{} 0.75 &{} 0.55 \\ -2.89 &{} -1.93 &{} -0.94 &{} 9.43 &{} 0.11 &{} 0.02 \\ 1.85 &{} 0.52 &{} 0.75 &{} 0.11 &{} 11.06 &{} 0.83 \\ 2.33 &{} 0.90 &{} 0.55 &{} 0.02 &{} 0.83 &{} 9.44 \end{matrix} \right] \text {,} \end{aligned}$$

(A-10)

with an efficient P-wave anisotropy strength of $\pm$5.5%

$$\begin{aligned} {\varvec{A}}_B= \left[ \begin{matrix} 30.63 &{} 11.08 &{} 6.99 &{} -2.21 &{} 1.48 &{} 1.46 \\ 11.08 &{} 29.29 &{} 11.46 &{} -0.63 &{} 2.67 &{} 0.20 \\ 6.99 &{} 11.46 &{} 23.60 &{} -1.87 &{} -0.56 &{} 1.24 \\ -2.21 &{} -0.63 &{} -1.87 &{} 8.20 &{} -0.21 &{} -1.41 \\ 1.48 &{} 2.67 &{} -0.56 &{} -0.21 &{} 9.12 &{} 0.44 \\ 1.46 &{} 0.20 &{} 1.24 &{} -1.41 &{} 0.44 &{} 9.59 \end{matrix} \right] \text {,} \end{aligned}$$

(A-11)

with an efficient P-wave anisotropy strength of $\pm$7.3%

$$\begin{aligned} {\varvec{A}}_C= \left[ \begin{matrix} 35.70 &{} 12.90 &{} 4.27 &{} -3.69 &{} 0.05 &{} 1.78 \\ 12.90 &{} 30.66 &{} 13.07 &{} -3.75 &{} -0.02 &{} -1.98 \\ 4.27 &{} 13.07 &{} 29.49 &{} -2.31 &{} 0.45 &{} 0.18 \\ -3.69 &{} -3.75 &{} -2.31 &{} 8.25 &{} -0.55 &{} -0.00 \\ 0.05 &{} -0.02 &{} 0.45 &{} -0.55 &{} 12.37 &{} 0.97 \\ 1.78 &{} -1.98 &{} 0.18 &{} -0.00 &{} 0.97 &{} 10.17 \end{matrix} \right] \text {,} \end{aligned}$$

(A-12)

with an efficient P-wave anisotropy strength of $\pm$8.1%