Shear wave splitting

Shear wave splitting, also called seismic birefringence, is the phenomenon that occurs when a polarized shear wave enters an anisotropic medium. The incident shear wave splits into two polarized shear waves. Shear wave splitting is typically used as a tool for testing the anisotropy of an area of interest. These measurements reflect the degree of anisotropy and lead to a better understanding of the area's crack density and orientation or crystal alignment.^[1] As an analogy, one can think of the anisotropy of a particular area as a black box and the shear wave splitting measurements as a way of looking at what is in the box.

An incident shear wave may enter an anisotropic medium from an isotropic media by encountering a change in the preferred orientation or character of the medium. When a polarized shear wave enters a new, anisotropic medium, it splits into two shear waves. One of these shear waves will be faster than the other and oriented parallel to the cracks or crystals in the medium. The second wave will be slower than the first and sometimes orthogonal to both the first shear wave and the cracks or crystals in the media. The time delays observed between the slow and fast shear waves give information about the density of cracks in the medium. The orientation of the fast shear wave records the direction of the cracks in the medium.

When plotted using polarization diagrams, the arrival of split shear waves can be identified by the abrupt changes in direction of the particle motion.

In a homogeneous material that is weakly anisotropic, the incident shear wave will split into two quasi-shear waves with approximately orthogonal polarizations that reach the receiver at approximately the same time. In the deeper crust and upper mantle, the high frequency shear waves split completely into two separate shear waves with different polarizations and a time delay between them that may be up to a few seconds.^[2]

Hess^[3] (1964) made the first measurements of P wave azimuthal velocity variations in oceanic basins. This area was chosen for this study because oceanic basins are made of large, relatively uniform homogeneous rocks. Hess observed, from previous seismic velocity experiments with olivine crystals, that if the crystals had even a slight statistical orientation this would be extremely evident in the seismic velocities recorded using seismic refraction. This concept was tested using seismic refraction profiles from the Mendocino fracture zone. Hess found that the slow compressional waves propagated perpendicular to the plane of slip and the higher velocity component was parallel to it. He inferred that the structure of oceanic basins could be recorded quickly and understood better if these techniques were used.

Ando^[4] (1980) focused on identifying shear-wave anisotropy in the upper mantle. This study focused on shear wave splitting recorded near the Chubu Volcanic Area in Japan. Using newly implemented telemetric seismographic stations, they were able to record both P wave and S wave arrivals from earthquakes up to 260 km beneath the volcanic area. The depths of these earthquakes make this area ideal for studying the structure of the upper mantle. They noted the arrivals of two distinct shear waves with different polarizations (N-S, fast and E-W, slow) approximately 0.7 seconds apart. It was concluded that the splitting was not caused by the earthquake source but by the travel path of the waves on the way to the seismometers. Data from other nearby stations were used to constrain the source of the seismic anisotropy. He found the anisotropy to be consistent with the area directly below the volcanic area and was hypothesized to occur due to oriented crystals in a deep rooted magma chamber. If the magma chamber contained elliptical inclusions oriented approximately N-S, then the maximum velocity direction would also be N-S, accounting for the presence of seismic birefringence.

Crampin^[5] (1980) proposed the theory of earthquake prediction using shear wave splitting measurements. This theory is based on the fact that microcracks between the grains or crystals in rocks will open wider than normal at high stress levels. After the stress subsides, the microcracks will return to their original positions. This phenomenon of cracks opening and closing in response to changing stress conditions is called dilatancy. Because shear wave splitting signatures are dependent on both the orientation of the microcracks (perpendicular to the dominant stress direction) and the abundance of cracks, the signature will change over time to reflect the stress changes in the area. Once the signatures for an area are recognized, they may then be applied to predict nearby earthquakes with the same signatures.

Crampin^[6] (1981) first acknowledged the phenomenon of azimuthally-aligned shear wave splitting in the crust. He reviewed the current theory, updated equations to better understand shear-wave splitting, and presented a few new concepts. Crampin established that the solution to most anisotropic problems can be developed. If a corresponding solution for an isotropic case can be formulated, then the anisotropic case can be arrived at with more calculations. The correct identification of body and surface wave polarizations is the key to determining the degree of anisotropy. The modeling of many two-phase materials can be simplified by the use of anisotropic elastic-constants. These constants can be found by looking at recorded data. This has been observed in several areas worldwide.^[7]

The difference in the travel velocities of the two shear waves can be explained by comparing their polarizations with the dominant direction of anisotropy in the area. The interactions between the tiny particles that make up solids and liquids can be used as an analogue for the way a wave travels through a medium. Solids have very tightly bound particles that transmit energy very quickly and efficiently. In a liquid, the particles are much less tightly bound and it generally takes a longer time for the energy to be transmitted. This is because the particles have further to travel to transfer the energy from one to another. If a shear wave is polarized parallel to the cracks in this anisotropic medium, then it may look like the dark blue wave in Figure 1. This wave is acting on the particles like energy being transferred through a solid. It will have a high velocity because of the proximity of the grains to each other. If there is a shear wave that is polarized perpendicular to the liquid-filled cracks or elongated olivine crystals present in the medium, then it would act upon these particles like those that make up a liquid or gas. The energy would be transferred more slowly through the medium and the velocity would be slower than the first shear wave. The time delay between the shear wave arrivals depends on several factors including the degree of anisotropy and the distance the waves travel to the recording station. Media with wider, larger cracks will have a longer time delay than a media with small or even closed cracks. Shear wave splitting will continue to occur until the shear-wave velocity anisotropy reaches about 5.5%.^[7]

Mathematical Explanation (Ray theory)^[8]

The equation of motion in rectangular Cartesian coordinates can be written as

∂ ∂ x i [ c i j k l ∂ U k ∂ x l ] = ρ ∂ 2 U j ∂ t 2 {\displaystyle {\frac {\partial }{\partial x_{i}}}\left[c_{ijkl}{\frac {\partial U_{k}}{\partial x_{l}}}\right]=\rho {\frac {\partial ^{2}U_{j}}{\partial t^{2}}}}

1

where t is the time, ρ {\displaystyle \rho } is the density, U j {\displaystyle U_{j}} is the component of the displacement vector U, and c i j k l {\displaystyle c_{ijkl}} represents the elastic tensor.
A wavefront can be described by the equation

t = τ ( x i ) {\displaystyle t=\tau \left(x_{i}\right)}

2

The solution to (1) can be expressed as a ray series

U k ( x i , t ) = ∑ n = 0 ∞ U k ( n ) ( x i ) f n ( t − τ ( x i ) ) {\displaystyle U_{k}\left(x_{i},t\right)=\sum _{n=0}^{\infty }U_{k}^{\left(n\right)}\left(x_{i}\right)f_{n}\left(t-\tau \left(x_{i}\right)\right)}

3

where the function f n ( ϑ ) {\displaystyle f_{n}\left(\vartheta \right)} satisfies the relation

d f n + 1 ( ϑ ) / d ϑ = f n ( ϑ ) {\displaystyle df_{n+1}\left(\vartheta \right)/d\vartheta =f_{n}\left(\vartheta \right)}

4

Substitute (3) into (1),

N ( U ( n ) ) − M ( U ( n − 1 ) ) + L ( U ( n − 2 ) ) = 0 {\displaystyle N\left(U^{\left(n\right)}\right)-M\left(U^{\left({n-1}\right)}\right)+L\left(U^{\left({n-2}\right)}\right)=0}

5

where the vector operators N,M,L are given by the formula:

{ N j ( U ( n ) ) = Γ j k U k ( n ) − U j ( n ) M j ( U ( n ) ) = p i   a i j k l ∂ U k ( n ) ∂ x l + ρ − 1 ∂ ∂ x i ( ρ   a i j k l   p l U k ( n ) ) L j ( U ( n ) ) = ρ − 1 ∂ ∂ x i ( ρ   a i j k l ∂ U k ( n ) ∂ x l ) {\displaystyle {\begin{cases}N_{j}\left(U^{\left(n\right)}\right)=\Gamma _{jk}U_{k}^{\left(n\right)}-U_{j}^{\left(n\right)}\\M_{j}\left(U^{\left(n\right)}\right)=p_{i}~a_{ijkl}{\frac {\partial U_{k}^{\left(n\right)}}{\partial x_{l}}}+\rho ^{-1}{\frac {\partial }{\partial x_{i}}}\left(\rho ~a_{ijkl}~p_{l}U_{k}^{\left(n\right)}\right)\\L_{j}\left(U^{\left(n\right)}\right)=\rho ^{-1}{\frac {\partial }{\partial x_{i}}}\left(\rho ~a_{ijkl}{\frac {\partial U_{k}^{\left(n\right)}}{\partial x_{l}}}\right)\end{cases}}}

6

where

Γ j k = p i   p l   a i j k l , a i j k l = c i j k l / ρ , p i = ∂ τ ∂ x i {\displaystyle \Gamma _{jk}=p_{i}~p_{l}~a_{ijkl},\quad a_{ijkl}=c_{ijkl}/\rho ,\quad p_{i}={\frac {\partial \tau }{\partial x_{i}}}}

7

For the first order n = 0 {\displaystyle n=0} , so U ( − 1 ) = U ( − 2 ) = 0 {\displaystyle U^{\left(-1\right)}=U^{\left(-2\right)}=0} , and only the first component of the equation (5) is left.
Thus,

N j ( U ( 0 ) ) = Γ j k U k ( 0 ) − U j ( 0 ) == ( Γ j k − δ j k ) U k ( 0 ) == 0 {\displaystyle N_{j}\left(U^{\left(0\right)}\right)=\Gamma _{jk}U_{k}^{\left(0\right)}-U_{j}^{\left(0\right)}==\left(\Gamma _{jk}-\delta _{jk}\right)U_{k}^{\left(0\right)}==0}

8

To obtain the solution of (8), the eigenvalues and eigenvectors of matrix Γ j k {\displaystyle \Gamma _{jk}} are needed,

D e t ( Γ j k − G δ j k ) = 0 {\displaystyle Det\left(\Gamma _{jk}-G\delta _{jk}\right)=0}

9

which can be rewritten as

G 3 − P G 2 + Q G − R = 0 {\displaystyle G^{3}-PG^{2}+QG-R=0}

9

where the values P , Q {\displaystyle P,Q} and R {\displaystyle R} are the invariants of the symmetric matrix Γ j k {\displaystyle \Gamma _{jk}} .
The matrix Γ j k {\displaystyle \Gamma _{jk}} has three eigenvectors: g 1 ,   g 2 ,   g 3 {\displaystyle g_{1},~g_{2},~g_{3}} , which correspond to three eigenvalues of G 1 ,   G 2 , {\displaystyle G_{1},~G_{2},} and   G 3 {\displaystyle ~G_{3}} .

• For isotropic media, G 1 = α 2 p i p i {\displaystyle G_{1}=\alpha ^{2}p_{i}p_{i}} corresponds to the compressional wave and G 2 = G 3 = β 2 p i p i {\displaystyle G_{2}=G_{3}=\beta ^{2}p_{i}p_{i}} corresponds to the two shear waves traveling together.
• For anisotropic media, G 2 ≠ G 3 {\displaystyle G_{2}\neq G_{3}} , indicates that the two shear waves have split.

^[9] In an isotropic homogeneous medium, the shear wave function can be written as

u ( ω ) = A   w ( ω ) exp ⁡ [ − i ω T 0 ] ⋅ p ^ {\displaystyle u\left(\omega \right)=A~w\left(\omega \right)\exp \left[-i\omega T_{0}\right]\cdot {\hat {p}}}

10

where A is the complex amplitude, w ( ω ) {\displaystyle w\left(\omega \right)} is the wavelet function (the result of the Fourier transformed source time function), and p ^ {\displaystyle {\hat {p}}} is a real unit vector pointing in the displacement direction and contained in the plane orthogonal to the propagation direction.
The process of shear wave splitting can be represented as the application of the splitting operator to the shear wave function.

Γ = exp ⁡ [ i ω δ t / 2 ] f ^ f ^ + exp ⁡ [ − i ω δ t / 2 ] s ^ s ^ {\displaystyle \Gamma =\exp \left[i\omega \delta t/2\right]{\hat {f}}{\hat {f}}+\exp \left[-i\omega \delta t/2\right]{\hat {s}}{\hat {s}}}

11

where f ^ {\displaystyle {\hat {f}}} and s ^ {\displaystyle {\hat {s}}} are eigenvectors of the polarization matrix with eigenvalues corresponding to the two shear wave velocities.
The resulting split waveform is

u s ( ω ) = A   w ( ω ) exp ⁡ [ − i ω T 0 ] Γ ( ϕ , δ t ) ⋅ p ^ {\displaystyle u_{s}\left(\omega \right)=A~w\left(\omega \right)\exp \left[-i\omega T_{0}\right]\Gamma \left(\phi ,\delta t\right)\cdot {\hat {p}}}

12

Where δ t {\displaystyle \delta t} is the time delay between the slow and fast shear waves and ϕ {\displaystyle \phi } is the angle between the polarization of the incident shear wave p ^ {\displaystyle {\hat {p}}} and the polarization of the fast shear wave f ^ {\displaystyle {\hat {f}}} . These two parameters can be individually estimated from multiple component seismic recordings.

Shear wave splitting measurements have been used to explore earthquake prediction, and to map fracture networks created by high pressure fracturing of reservoirs.

According to Crampin^[5] shear wave splitting measurements can be used to monitor stress levels in the earth. It is well known that rocks near an earthquake-prone zone will exhibit dilatancy. Shear wave splitting is produced by seismic waves traveling through a medium with oriented cracks or crystals. The changes in shear wave splitting measurements over the time leading up to an impending earthquake can be studied to give insight to the timing and location of the earthquake. These phenomena may be observed many hundreds of kilometers from the epicenter.

The petroleum industry uses shear-wave splitting measurements to map the fractures throughout a hydrocarbon reservoir. To date, this is the best method to gain in situ information about the fracture network present in a hydrocarbon reservoir.^[10] The best production in a field is associated with an area where there are multiple small fractures that are open, allowing for constant flow of the hydrocarbons. Shear-wave splitting measurements are recorded and analyzed to obtain the degree of anisotropy throughout the reservoir. The area with the largest degree of anisotropy will generally be the best place to drill because it will contain the largest number of open fractures.^[11]

On October 27, 1998, during a four-year study of shear wave splitting in Iceland, Crampin and his coworkers recognized that time delays between split shear-waves were increasing at two seismic recording stations, BJA and SAU, in southwest Iceland. The following factors lead the group to recognize this as a possible precursor to an earthquake:^[12]

• The increase persisted for nearly 4 months.
• It had approximately the same duration and slope as a previously recorded magnitude 5.1 earthquake in Iceland.
• The time delay increase at station BJA started at about 4 m s / k m {\displaystyle 4ms/km} and escalated to approximately 10 m s / k m {\displaystyle 10ms/km} .
• 10 m s / k m {\displaystyle 10ms/km} was the inferred level of fracture for the previous earthquake.

These features suggested that the crust was approaching fracture criticality and that an earthquake was likely to occur in the near future. Based on this information, an alert was sent to the Iceland Meteorological Office (IMO) on October 27 and 29, warning of an approaching earthquake. On November 10, they sent another email specifying that an earthquake was likely to occur within the next 5 months. Three days later, on November 13, IMO reported a magnitude 5 earthquake near the BJA station. Crampin et al. suggests that this is the first scientifically, as opposed to precursory or statistically, predicted earthquake. They proved that variations of shear-wave splitting can be used to forecast earthquakes.

This technique was not successful again until 2008 due to the lack of appropriate source-geophone-earthquake geometry needed to evaluate changes in shear wave splitting signatures and time delays.^[7]

Volti and Crampin observed temporal increases in Band-1 time-delays for 5 months at approximately 240 kilometer depth in directions N,SW and W,SW before the 1996 Gjalp Eruption in Vatnajökull Icefield. This was the largest eruption in Iceland in several decades.

The pattern of increasing shear wave splitting time-delays is typical of the increase now seen before many earthquakes in Iceland and elsewhere. The time delays just before earthquakes characteristically decrease immediately following the eruption because the majority of the stress is released at that one time. The increase in normalized time-delays in volcanic eruptions does not decrease at the time of the eruption but gradually declines at about 2 m s / k m / y e a r {\displaystyle 2ms/km/year} over several. This decrease is approximately linear and there appeared to be no other significant magmatic disturbances during the period following the eruption.

More observations are needed to confirm whether the increase and decrease time delay pattern is universal for all volcanic eruptions or if each area is different. It is possible that different types of eruptions show different shear wave splitting behaviors.^[7][13]

Bokelmann and Harjes reported the effects on the shear waves of fluid injection at about 9 kilometer depth in the German Continental Deep Drilling Program (KTB) deep drilling site in southeast Germany. They observed shear-wave splitting from injection-induced events at a pilot well offset 190 meters form the KTB well. A borehole recorder at a depth of 4,000 meters was used to record the splitting measurements.^[14]

They found:

• Temporal variations in shear-wave splitting as a direct result of injection-induced events.
• That the initial ~1% shear wave splitting decreases by 2.5% in the next 12 hours following the injection.
• The largest decrease occurred within two hours after the injection.
• The splitting time to be very stable after the injection ceased.

No direct interpretation of the decrease is proposed but it is suggested that the decrease is associated with stress release by the induced events.

Shear-wave splitting measurements can provide the most accurate and in depth information about a particular region. However, there are limits that need to be accounted for when recording or analyzing shear wave splitting measurements. These include the sensitive nature of shear waves, that shear wave splitting varies with incidence and azimuth, and that shear waves may split multiple times throughout an anisotropic medium, possibly every time the orientation changes.^[15]

Shear wave splitting is very sensitive to fine changes in the pore pressure in the Earth's crust. In order to successfully detect the degree of anisotropy in a region there must be more several arrivals that are well distributed in time. Too few events cannot detect the change even if they are from similar waveforms.^[7] The Shear wave splitting varies with both incidence angle and propagation azimuth. Unless this data is viewed in polar projection, the 3-D nature is not reflected and may be misleading.^[7] Shear wave splitting may be caused by more than just one layer that is anisotropic and located anywhere between the source and the receiver station. The shear wave splitting measurements have extensive lateral resolution but very poor vertical resolution.^[16] The polarizations of shear waves vary throughout the rock mass. Therefore, the observed polarizations may be those of the near surface structure and are not necessarily representative of the structure of interest.^[17]

Due to the nature of split shear waves, when they are recorded in typical three-component seismograms, they write very complicated signatures. Polarizations and time delays are heavily scattered and vary greatly both in time and space. Because of the variation in signature, it is easy to misinterpret the arrivals and polarization of incoming shear waves.^[18] Below is an explanation of a few of the common misunderstandings associated with shear waves, further information can be found in Crampin and Peacock (2008).^[7]

• Polarizations of split shear waves are orthogonal.^[7][18]

Shear waves that propagate along the ray path at a group velocity have polarizations that are only orthogonal in a few specific directions. Polarizations of body waves are orthogonal in all phase velocity directions, however this type of propagation is generally very difficult to observe or record.

• Polarizations of split shear-waves are fixed, parallel to cracks, or normal to spreading centers.^[7][18]

Even when propagating through parallel cracks or perpendicular to spreading centers or parallel to cracks, the polarizations of shear waves will always vary in three dimensions with incidence and azimuth within the shear wave window.

• Crack anisotropy always decreases with depth as fluid filled cracks are closed by lithostatic pressure.^[7][18]

This statement only holds true if the fluid in the cracks is somehow removed. This may be accomplished via chemical absorption, drainage, or flow to the surface. However, these occur in relatively rare instances and there is evidence that supports the presence of fluids at depth. This includes data from the Kola deep well and the presence of high conductivity in the lower crust.

• Signal-to-noise ratios of shear-wave splitting above small earthquakes can be improved by stacking.^[7][18]

Stacking seismic data from a reflection survey is useful because it was collected with a predictable, controlled source. When the source is uncontrolled and unpredictable, stacking the data only degrades the signal. Because recorded shear wave time delays and polarizations vary in their incidence angle and azimuth of radio propagation, stacking these arrivals will degrade the signal and decrease the signal to noise ratio, resulting in a plot that is noisy and hard to interpret at best.^[7]

Our^[who?] understanding of shear wave splitting and how to best use the measurements is constantly improving. As our^[who?] knowledge improves in this area, there will invariably^[why?] be better ways of recording and interpreting these measurements and more opportunities^{[according to whom?]} to use the data. Currently,^[when?] it is being developed^{[by whom?]} for use in the petroleum industry and for predicting earthquakes and volcanic eruptions.

Shear wave splitting measurements have been used^{[by whom?]} successfully to predict several earthquakes. With better equipment and more densely spaced recording stations, we^[who?] have been able to study the signature variations of shear wave splitting over earthquakes in different regions. These signatures change over time to reflect the amount of stress present in an area. After several earthquakes have been recorded and studied, the signatures of shear wave splitting just before an earthquake occurs become well known and this can be used to predict future events. This same phenomenon can be seen before a volcanic eruption and it is inferred that they may be predicted in the same manner.^{[citation needed]}

The petroleum industry^[who?] has been using shear wave splitting measurements recorded above hydrocarbon reservoirs to gain invaluable information^[which?] about the reservoir for years. Equipment is constantly^[when?] being updated to reveal new images and more information.^[7]

• Alfred Wegener Institute for polar and Marine Research(AWI)(Germany)
• Shear-wave splitting in Matlab(France)
• A lot of interesting seismic images(ASU)
• Information on Solids, Liquids, and Gasses
• Shear-Wave-Splitting-in-Anisotropic-Media Collection (Goethe University Frankfurt, Germany)

You can download a MATLAB code and create a demonstration movie by yourself here on MathWorks website.