Project topics and materials
Get Complete Project Material File(s) Now! »
Equipartition of the energy
One common feature to all scattering systems is that they eventually reach to a diffusive state in which each normal mode is excited equally, independently of the initial conditions of the field. The amount of energy is proportional to the number of modes that can be sustained in the medium (Hill, 1986; Morse and Bolt, 1944), or in other words, to the density of states; for example, in an elastic medium there can be two transversal modes of propagation (s-waves), and one longitudinal (p-wave). Then, the energy partition between the two types of waves is (Weaver, 1982) where c and c are the velocities of P-waves and S-waves respectively. The density of states can also be calculated in another way: the Green’s function of a linear operator L can be defined as [z L(r)]G r; r0 = r r0 (2.12) where z is a complex variable, and L represents the wave equation operator. L has a set of eigenfunctions un and eigenvalues n that per definition obey the following relation L(r)un(r) = nun(r) (2.13)
Each eigenfunction represents the a normal mode in the system. It can be shown that the spectral representation of the Green’s function is an analytic function except for the points where z = n, which means that it diverges at the resonant modes of the system (Economou, 2006). This important property can be used to isolate and count the normal modes of the system with the Green’s function. As an example, the density of states for the Helmholtz equation can be expressed as (Sheng and Tiggelen, 2007)
2(!) 1 Im G !; r = r0
(!; r) = dkd! (2.14)
where ! represents the frequency, and k the wavenumber.
The equipartition of the energy was observed by Shapiro et al. (2000) and later by Hennino et al. (2001) who made a series of observations theoretical estimations of equipartition ratios in different combinations of seismic waves. Margerin, Campillo, and Van Tiggelen (2000) studied numerically the equipartition time in inhomogeneous elastic media using Monte Carlo simulations that included mode conversions and polarization.
Coda wave interferometry
The seismic wavefronts propagating through the earth get distorted by the interaction with heterogeneities that scatter a fraction of their energy in different directions. This process produces the characteristic coda wave segments after the main arrivals in the seismograms (Aki, 1969). Although complicated to analyze, the coda wave is not random and can be understood as the superposition of all the possible wave fronts that propagated through many different paths when going from the source to the receiver (Sato and Fehler, 2012). This implies that a complex seismogram generated by a seismic event can be reproduced if, by some mean, we repeat the same seismic source at the same position; this has been observed with two or more earthquakes that occur in the same location with very similar magnitudes (Geller and Mueller, 1980; Poupinet, Ellsworth, and Frechet, 1984; Beroza, Cole, and Ellsworth, 1995). However, although very similar, the received waveforms at the surface are not exactly equal: small variations are produced as a consequence of small changes of the properties of the medium that happen in the period of time between the two seismic events. This is the main principle behind some studies that use these alike earthquakes, called doublets, to track the changes in the velocity (Poupinet, Ellsworth, and Frechet, 1984; Ratdomopurbo and Poupinet, 1995) or in the attenuation (Beroza, Cole, and Ellsworth, 1995) of the crust. We also refer the reader to Sato (1988) for a review of temporal variations of attenuation detected with coda waves. Coda wave interferometry (CWI) is the use of the coda of a waveform to detect an locate all these changes in the medium.
CWI has been applied in experiments to monitor changes of ultrasound waves propagating through granite, (Snieder et al., 2002; Grêt, Snieder, and Scales, 2006), in a volcano using its natural mild repeating explosions(Grêt et al., 2005), in concrete induced by thermal vari-ations (Larose et al., 2006), in the crust as a consequence of earthquakes (Nishimura et al., 2000; Peng and Ben-Zion, 2006), and associated to structural changes preceding eruptions (Ratdomopurbo and Poupinet, 1995; Wegler et al., 2006).
From the assumption of superposition of the field, the travel-time perturbation h i pro-duced by a small perturbation of the propagation velocity between a source and a receiver, can be understood as the average of the travel time perturbations P generated between all the possible trajectories P (Pacheco and Snieder, 2005; Snieder, 2006) t = h i = PPP IP (2.15)
I the wave that has propagated along path where P represents the energy intensity ofP P . This result is very interesting as it establishes a relationship between the eventual phase perturbations of the coda waves and the energies arriving to a certain point which are the main physical quantity of the radiative transfer theory that assumes that the phases cancel each other at the receiver. This is fundamental in the conceptual development of the sensitivity kernels.
With coda wave interferometry it is possible to relate different type of perturbations to measurements made in the coda. These perturbations can be of different nature: velocity fluctuations, random displacements of punctual scatterers and perturbations in the position of the source (Snieder et al., 2002; Snieder, 2006); here we focus on the first one. A typical phase change in the waveform recorded at the surface, generated by a change in the velocity of a medium can be seen in figure 2.3. The difference between the waveforms obtained with and without the velocity change increases proportionally to the time that the wave has spent traveling in the medium. In this case, the relationship between the velocity variation and the phase delay at different lapse times is t=t = v=v. The quantity t=t is called the apparent velocity variation, also denoted » by some authors. The term apparent is added because when the velocity is changed in only one part of the medium, there is still a change in the phases at the surface, but the quantity t=t is not necessarily equal to the velocity pertubation v=v in the medium anymore.
MCSW method
In monitoring, a linear relation between the phase delay and the velocity change in the medium is assumed; therefore, the apparent velocity variation is obtained as the slope of the phase delay t for different lapse times t. One of the main methods to estimate this slope is the Moving Cross Spectral Window (MCSW) also known as the doublets method, as it was first proposed and used to analyze the phase differences between the coda generated by doublets earthquakes by Poupinet, Ellsworth, and Frechet (1984). The MCSW method is based on a well know property of the Fourier transform: let u(t) represent a signal and U(!) = F[u(t)] its Fourier transform. Then the Fourier transform of temporal shifted version of this signal u(t t) is (Pinsky, 2008) F [u(t t)] = U(!)e i! t (2.16)
Notice how this temporal shift affects different frequencies at different degrees; a compar-ison between the Fourier transform of the original and the shifted signals would show a linear relationship between their phase differences and the frequency. In this way, the t can be estimated for a small section of the coda of the two signals. Since this delay is different for different lapse times, this procedure is repeated over a window that is systematically moved along the coda. This procedure is schematized in figure 2.4. In the end, a linear regression is performed between the obtained values of t and their respective lapse times to obtained the apparent velocity estimation.
Seismic interferometry
The link between recordings of correlations of noise and the properties of the medium can be traced back to early works done in helioseismology (Duvall et al., 1993; Giles et al., 1997). Later, a formal connection of these cross-correlations with the retrieval of the Green’s func-tion of the medium was established (Weaver and Lobkis, 2004; Lobkis and Weaver, 2001; Wapenaar, 2004). The Green’s function represents the response that one part of the system (in this case, the displacements on the surface) would have as a consequence of an impulsive unit source in the medium. This relationship permits the use of tools of analysis in seismol-ogy that pass from monitoring to imaging, to regions where there is no presence of active sources. The reconstruction of the Green’s function using ambient seismic noise is called Seismic Interferometry (SI).
Let u (t; x1) and u (t; x2) denote time dependent displacements recordings of noise seismic field at two points over the surface. Then the noise cross correlation over a time interval [0; T ] can be written as CT ( ; x1; x2) = u (t; x1) u (t + ; x2) dt where represents the time lag between the two signals. In a diffusive, equipartitioned field it can be shown that (Snieder, 2004; Roux et al., 2005a)
@ CT ( ; x1; x2) = G ( ; x1; x2) G ( ; x1; x2) (2.18) as long as T is sufficiently large. The positive and negative parts of the Green’s function are related to the causal and anti-causal parts of the cross-correlation, shown in figure 2.5. Each of them represents two possible waves traveling in opposite directions from one station to another. The equipartitioning of the seismic field can be generated with a good distribution of seismic sources, or as a consequence of the scattering produced by heterogeneities in the earth. Several factors can affect the formation of a diffusive field like strong attenuation, an asymmetric distribution of noise sources or a limited observation time window of the coda (Paul et al., 2005; Stehly et al., 2008). in practice, the Green’s function is recovered through the average of many correlations obtained with a dense seismic network or over several days. Maybe the most visible feature of the cross-correlations of ambient seismic noise is the reconstruction of the travel times between several stations as is shown in the image 2.6. It has been shown, nonetheless, that even in situations where the Green’s function is not perfectly reconstructed, it is still possible to track changes in the medium using the cross-correlation (Hadziioannou et al., 2009).
Passive image interferometry
Sens-Schönfelder and Wegler (2006) showed that not only the ballistic surface waves and P-waves can be retrieved from correlations of continuous recordings of ambient seismic noise (Roux et al., 2005b) but also the scattered waves that form the coda. Coda Wave Interfer-ometry applied to the Green’s function coda recovered by ambient seismic noise methods is called Passive Image Interferometry (PSI).
PSI has been used to measure the variation of the seismic field velocity caused by a variety of phenomena like earthquakes (Wegler and Sens-Schönfelder, 2007; Brenguier et al., 2008a; Wegler et al., 2009; Nakata and Snieder, 2011), volcanic activity (Rivet, Brenguier, and Cappa, 2015; Brenguier et al., 2008b), the thermoelastic response of the soil (Meier, Shapiro, and 5000s) correspond to surface waves that have traveled between the two stations making a turn around the earth. Taken from Nakata, Gualtieri, and Fichtner (2019)
Brenguier, 2010), clay landslide (Mainsant et al., 2012), pressurized water injection (Hillers et al., 2015), internal erosion in earthen dams (Planès et al., 2016), the earth tides produced by the sun and the moon (Sens-Schönfelder and Eulenfeld, 2019), temperature variations due to periodic heating of the lunar surface by the sun (Sens-Schönfelder and Larose, 2008), and variations in the loading and pore pressure perturbations over and below the glacial cover (Mordret et al., 2016). Analysis of the ambient noise can also be used to predict earthquake ground motion (Prieto and Beroza, 2008) and study the seismic response of civil structures (Prieto et al., 2010; Nakata et al., 2013).
Many studies focus on the hydrological effects on the v=v: Sens-Schönfelder and We-gler (2006) estimated the underground water level using a model developed by Akasaka and Nakanishi (2000), to make a direct relation with the measured velocity variations in a vol-cano. Meier, Shapiro, and Brenguier (2010) analyze velocity variations within the Los Angeles basin and conclude that the seasonal variations are strongly influenced by groundwater level changes and thermo-elastic strain variations. Tsai (2011) proposed periodic models to recre-ate displacements and velocity changes from thermo-elastic stresses or hydrological loadings. Wang et al. (2017) found a direct relation between velocity variations and several hydrological and meteorological processes across Japan, mainly based on the pore pressure generated by the rainfall water through a diffusion process. Hillers, Campillo, and Ma (2014) also show the correlation between velocity changes and periodicity of precipitation events in Taiwan.
Sensitivity kernels
Observation and motivation
The phase fluctuations generated in the surface as a consequence of the presence of a bulk velocity perturbation do not necessarily have a linear relationship with the lapse time; their relationship depends on the position and size of the perturbation. The behavior of the ap-parent velocity variation t=t for two opposite cases is shown in figure 2.7. These results are obtained in a 2-D elastic medium with both the source and the receivers located at the surface. Both are examples of the apparent velocity perturbations measured on the surface when the velocity of a whole layer in a medium is slightly modified. At the left the t=t mea-surements are obtained when the layer is located close to the surface, and at the right when the perturbation is deeper in the medium. In the first one, the apparent velocity variation at early times is high as the perturbed layer is quickly sampled by most of the energy emitted by the source; however, it quickly decreases as the seismic field expands into deeper regions, making the perturbation relatively smaller in comparison with the total explored area. In the second case, the perturbation is only detected by the front of the propagating energy at the beginning and is later sampled progressively more with the ongoing expansion of the wave field. Notice that the magnitudes also change dramatically: in the latter case, t is one order of magnitude smaller than in the former case. This simple example shows that the difference of behavior of the apparent velocity variation opens the possibility of estimating the shape and location of the perturbation with the measurements made on the surface: this is the main objective of developing the travel-time sensitivity kernels.
The travel-time sensitivity kernel was first introduced by Pacheco and Snieder (2005) for the diffusive regime, and by Pacheco and Snieder (2006) for the single scattering regime. These sensitivity kernels relate a travel-time perturbation measured between a source and a receiver with all the possible velocity perturbations in the medium around them. For the first one, an analytical expression was obtained for the autocorrelation configuration (coincident source and receiver) for two and three dimensions and was compared against numerical simulations where a small perturbation is introduced. The travel-time sensitivity kernel for the single scattering regime was calculated following the approach of Sato (1977) to obtain the energy density of scattered waves assuming single isotropic. Later Planès et al. (2014) introduced the decorrelation sensitivity kernel that related a perturbation in the medium that alters the paths going from the source to the receiver, with the distortion that it generates to the waveform at the receiver. Later Margerin et al. (2016) made an alternative derivation of this kernel and extended the travel-time sensitivity kernel to include the direction of the energy propagation with the use of the specific intensity. In this work was proposed the idea that the difference between the travel-time and the decorrelation sensitivities depended on whether the perturbation had an active or passive role over the generation of new propagation paths between the source and the receiver. Applications like localizing perturbations in numerical simulations have been done using the sensitivity kernel. The decorrelation sensitivity kernel was successfully used in locating millimetric holes drilled in a concrete sample (Larose et al., 2010). The sensitivity kernels can be expressed as convolutions of intensity which makes the radiative transfer theory a natural tool for their estimation; Lesage, Reyes-Dávila, and Arámbula-Mendoza (2014), for example, makes use of a solution of the transport equations for 2-D to locate changes produced by volcanic activity. On the other hand, Obermann et al. (2016) makes use of the 3-D radiative transfer solution proposed by Paasschens (1997) to estimate the body wave sensitivity kernel.
The coupling between body and surface waves has remained one challenging factor to the development of the sensitivity kernels in a 3-D halfspace that is the usual setting for most seismic applications. Obermann et al. (2013a) and Obermann et al. (2016) approached this problem by expressing the sensitivity as a linear combination of two independent sensitivities, one for surface and another for body waves, with a controlling factor mediating between them that changes with time, and that is estimated by comparisons with full wavefield numerical simulations. These studies show the predominance of surface waves at early lapse times and of body waves at late lapse times. However, Wu et al. (2016) measured velocity variations from the phase of the Rayleigh waves obtained by seismic interferometry, and observed a progressive decrease of the velocity drop produced by the 2004 Parkfield earthquake (Brenguier et al., 2008a) at lower frequencies. Since the surface waves penetration increases at lower frequencies, they conclude that the velocity variation is mostly constrained in the surface. Furthermore, they concluded that the observations can be explained only by the surface wave dispersion as body wave sensitivity would produce both velocities decreases and increases at different points; this last affirmation was supported by travel-time sensitivity kernels that show that in a spherical geometry, direct shear waves may have alternating early or late arrivals depending on the position on the surface of the receiver (Stein and Wysession, 2009; Zhao, Jordan, and Chapman, 2000). Similar findings were reported by Nakata and Snieder (2011) that found that the magnitude (MW ) 9.0 Tohoku-Oki earthquake produced a shear wave velocity reduction limited to the first 100m of the crust, despite being extended to an area 1200km wide by comparing arrival times between stations located at the surface and in a borehole in ranges of 100m to 337m depth. This, however, contradicts Wang et al. (2019) that through an inversion made with seismic and long-period tiltmeters data, found that the same earthquake produced changes that reach 40km depth in the crust.
Basic theory of traveltime sensitivity kernels
We begin with a concise derivation of the travel-time sensitivity kernel for a diffusive medium following its intensity interpretation (Pacheco and Snieder, 2005), and later make a connection with the more general version based on the specific intensities (Margerin et al., 2016). Let us imagine a normalized intensity impulse generated at a source, that propagates through a medium. The intensity at lapse time t at a point r is then described by P (r; t). If there is no mechanism of intrinsic attenuation acting on the system, the total energy at time t W (V; t) = P (r; t)dV (r) (2.19) is always equal to 1 , that is, the total energy emitted by the source. The normalization of the intensity allows us to make two additional interpretations of this quantity: P (r; t) represents the probability that a particle, following a random walk, arrives or passes by the position r at time t (Roepstorff, 2012). P (r; t) can also be understood as the Green’s function of the diffusion equation at position r at time t. The probability that a particle that was emitted at the source at s arrives to the receiver at r at time t, can be written in terms of the transit of the phonon through a volume dV located at r0 at time t0 (illustrated in figure 2.8).
P r; s; t + t0 = P r; r0; t t0 P r0; s; t0 dV r0 (2.20)
This is the Chapman-Kolmogorov equation (Ross, 2014; Papoulis and Pillai, 2002; Roep-storff, 2012) which states that the probability of going from the source to the receiver can be written by using the intermediate point r0 at a time t0, as long as all of them are taken into account. Each of the terms involved in the last equation are transition probabilities; for example P r0; s; t0 = P r0; t0js; t = 0 (2.21) is the probability that particle transits the position r0 at time t0, under the condition that it was emitted from the source at position s and at time t = 0. An integration over all the possible travel-times of these equations leads to t = ZV Z t P (r; r ; t ) P (r ; s; t ) dV (r0) 0 t 00 0 dt0 (2.22) 0 P (r; s; t)
The term inside the square parenthesis is the sensitivity kernel K. One of the most remarkable features of the kernel is that it expresses the travel time as a spatial distribution over the volume around the source and the receiver. An observer that could follow all the phonons going from the source to the receiver, could measure the sensitivity as the average time spent inside each volume of the medium (Margerin et al., 2016); therefore, the sensitivity kernels show which parts of the medium are preferred by the particles when going from the source to the receiver. As such, the regions around the source, the receiver, and between them, are zones that we expect to have high sensitivity as they are probably very frequented by the phonons that travel between the two. A simple schematic explanation of the travel-time sensitivity can be seen in figure 2.9.
Table of contents :
1 Introduction
2 Concepts and methods
2.1 Radiative transfer theory
2.1.1 Radiative transfer equation
2.1.2 Diffusion approximation
2.2 Equipartition of the energy
2.3 Coda wave interferometry
2.3.1 MCSW method
2.4 Seismic interferometry
2.4.1 Passive image interferometry
2.5 Sensitivity kernels
2.5.1 Observation and motivation
2.5.2 Basic theory of traveltime sensitivity kernels
3 Separation of phenomena from v=v measurements
3.1 Abstract
3.2 Introduction
3.3 Data processing description
3.4 Procedure and results
3.4.1 Measured and modeled velocity variations
3.4.2 Analysis of geodetic data
3.5 Loading effect of the rainfall
3.6 Conclusions
3.7 Supplementary information
4 Coupling between surface and body waves
4.1 Abstract
4.2 Introduction
4.3 Scalar wave equation model with surface and body waves
4.3.1 Equation of motion
4.3.2 Eigenfunctions and Green’s function
4.3.3 Source radiation and density of states
4.4 Single scattering by a point scatterer
4.4.1 Scattering of a surface wave
4.4.2 Scattering of body waves
4.5 Equation of radiative transfer
4.5.1 Phenomenological derivation
4.5.2 Energy conservation and equipartition
4.5.3 Diffusion Approximation
4.6 Monte-Carlo Simulations
4.6.1 Overview of the method
4.6.2 Numerical results
4.7 Conclusions
4.8 Supplementary information
4.8.1 Variational formulation for mixed boundary conditions
4.8.2 Far-field expression of the Green’s function for scalar waves in a halfspace with mixed B.C
5 Sensitivity Kernels
5.1 Introduction
5.2 Scalar model
5.3 Penetration depth of the surface wave
5.4 Phase velocity variation for surface waves
5.5 Group velocity
5.6 Time densities
5.7 Sensitivity kernels
5.8 Monte Carlo simulations
5.8.1 General outline
5.9 Surface and body wave sensitivity
5.9.1 Time partition coefficient
5.10 Spatial scaling parameter
5.11 Conclusions
6 Recovery of velocity variations at depth
6.1 Introduction
6.1.1 The inversion problem
6.2 Setting of the seismic inverse problem
6.2.1 Construction of the operator G
6.2.2 Parameters election
6.3 Performance of the inversion in synthetic tests
6.3.1 Noise-free inversion: Effects of the downsampling process
6.3.2 Performance at different depths
6.3.3 Performance with different lapse times
6.3.4 Performance with different levels of noise
6.4 Conclusions
7 Conclusions and perspectives
