# Observation of Equipartition of Seismic Waves

R. Hennino, N. Trégourès, N. M. Shapiro, L. Margerin, M. Campillo, B. A. van Tiggelen, R. L. Weaver

► To cite this version:

R. Hennino, N. Trégourès, N. M. Shapiro, L. Margerin, M. Campillo, et al.. Observation of Equipartition of Seismic Waves. Physical Review Letters, 2001, 86, pp.3447-3450. 10.1103/Phys-RevLett.86.3447. insu-03606663

HAL Id: insu-03606663

https://insu.hal.science/insu-03606663

Submitted on 16 Mar 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Observation of Equipartition of Seismic Waves

R. Hennino,¹ N. Trégourès,² N. M. Shapiro,³ L. Margerin,¹,⁴ M. Campillo,¹ B. A. van Tiggelen,²* and R. L. Weaver⁵

¹Laboratoire de Géophysique Interne et Tectonophysique, Observatoire de Grenoble, Université Joseph Fourier, B.P. 53, 38041 Grenoble Cedex, France

²Laboratoire de Physique et Modélisation des Milieux Condensés, CNRS/Université Joseph Fourier, B.P. 166, 38042 Grenoble Cedex 09, France

³Instituto de Geofísica, Universidad Nacional Autónoma de México, CU 04510 México, D.F., Mexico

⁴Department of Geosciences, Princeton University, Guyot Hall, Princeton, New Jersey 08544

⁵Department of Theoretical and Applied Mechanics, University of Illinois at Urbana-Champaign, Urbana, Illinois 61801
(Received 29 September 2000)

Equipartition is a first principle in wave transport, based on the tendency of multiple scattering to homogenize phase space. We report observations of this principle for seismic waves created by earthquakes in Mexico. We find qualitative agreement with an equipartition model that accounts for mode conversions at the Earth's surface.

DOI: 10.1103/PhysRevLett.86.3447

PACS numbers: 91.30.Fn, 46.40.Cd, 62.30.+d, 91.30.Dk

Scattering of classical waves in controlled laboratory experiments has extended our knowledge of wave propagation considerably [1–3]. One new challenge is to apply this knowledge to waves in natural environments, whose properties are complex and not well known, for example, medical imaging [4], remote sensing [5], and seismology.

The seismic case is particularly difficult because several specific complications, such as near-field detection and mode conversions at the Earth's surface, rule out the application of methods developed in optics and acoustics. In addition, seismic records are intrinsically time resolved. Sheng et al. [6] considered the possibility of Anderson localization in time-dependent signals that have propagated through a layered Earth crust.

The seismic coda refers to the pronounced exponential time tail observed in the seismograms of regional earthquakes in the frequency band 1–10 Hz, reaching the level of seismic noise after sometimes more than 10 times the travel time of direct waves [7]. One important observational fact is that the time decay coefficient of the energy in the seismic coda, the so-called coda Q factor, is a regional constant [8], dependent on frequency, but independent of distance, depth, and magnitude of the seismic source.

Many studies have made an attempt to relate physical properties of the Earth's lithosphere to the observed regional seismic coda. The first pioneering studies by Aki and Chouet [7,8], followed by others [9,10], interpreted seismic coda as elastic waves that are singly scattered from inhomogeneities in the Earth's crust. It was shown that a combination of single scattering and dissipation can account for some aspects of the seismic coda [11]. Recent numerical studies suggested multiple scattering as the origin of seismic coda [12–16]. Our own study estimated mean free paths $\ell \approx 20–70$ km in the frequency band 1–10 Hz [16], comparable to the thickness H of the Earth's crust. For a shear wave velocity $\beta \approx 3.5$ km/s in the crust this implies a scattering mean free time of only 15 s, much smaller than the observed duration (180 s) of the coda.

A direct and *model-independent* confirmation of the persistence of multiple scattering is important for two reasons. First, it is important for future studies in extracting information correctly from the seismic coda, i.e., to do the inverse problem. Second, it would facilitate the application of many powerful mesoscopic tools developed in acoustics and optics [2,3]. However, several complications arise. One “elementary” test for multiple scattering is the randomization of the propagation direction of the waves [2]. Unfortunately, seismic observations measure local elastic displacements $u(r, t)$, and not “specific intensities,” i.e., angularly resolved fluxes, as in optics. Another elementary test would be an observation of the randomization of polarization [5] as indicating multiple scattering, since single scattering should be highly polarized for a polarized source such as an earthquake. The essential problem in seismic studies is that transverse S and longitudinal P waves cannot be separated in a ground displacement at one point. Only recently, Shapiro et al. [17] presented an approach to separate P from S waves using a small-aperture array. A second complication is the presence of the Earth’s free surface that coherently mixes transverse SV waves with longitudinal P waves [18].

In this Letter, observational data of seismic coda are confronted with a *first principle* in multiple scattering to resolve the controversy “single or multiple scattering in the seismic coda.” We have investigated the principle of *equipartition*, first discussed by Weaver [19,20] for elastic waves in a solid body, more recently by Papanicolaou et al. [21] and Turner [22] using a transport equation, and finally for light in liquid crystals [23]. It is based on the fact that multiple scattering tends to homogenize phase space: the spectral energy density around a specified frequency, originally distributed in phase space in a fashion that largely depends on the nature of the source, eventually becomes uniform. Equipartition has the remarkable property to be independent of the unknown fluctuations that cause the scattering. It would imply any energy ratio
to become *time independent*, and *independent* of the magnitude, polarization, and distance of the seismic source. In single scattering from a polarized source these energy ratios would vary significantly with time and among different sources.

The total elastic wave energy density is given by [24]

$$W = \frac{\rho}{2} (\partial_t \mathbf{u})^2 + \left(\frac{\lambda}{2} + \mu\right) (\mathrm{div} \, \mathbf{u})^2 + \frac{\mu}{2} (\mathrm{curl} \, \mathbf{u})^2 + I, \quad (1)$$

with $\mathbf{u}(\mathbf{r}, t)$ the local displacement vector, and $\rho$ the mass density. Deep-reflection and surface-wave studies have shown that the Lamé coefficients $\lambda$ and $\mu$ of the strain-stress tensor are roughly constant and equal in the Earth's crust, with fluctuations less than 5%. The different terms in Eq. (1) represent kinetic ($K$), compressional ($P$), and shear ($S$) energy density. The last term ($I$) is an interference term involving cross terms of the kind $\partial_i u_j \partial_k u_l$ whose averages vanish except near boundaries. We have chosen to study $S/P$, $K/(S + P)$, and $I/(S + P)$, and the ratio $H^2/V^2$ of kinetic energies for elastic displacements in different directions for which no spatial derivatives have to be carried out.

The most direct way to facilitate the evaluation of the three partial derivatives would be to install seismic receivers close to each other and at different depths. However, installing numerous seismic receivers, especially at depth, is extremely expensive. An array of four receivers at the corners of a 50 m side square [17] was temporarily set up close to the city of Chilpancingo (Mexico). The high seismic rate in Mexico made it possible to record, during the three months of the experiment, a series of local earthquakes with magnitudes between 4 and 5, and a spread of roughly 300 km in epicentral distance. The sensors are CMG-40T seismometers connected to Reftek digitizers, whose absolute time was provided by GPS satellites. The three vertical derivatives $\partial_z u_i$ were deduced from the measured two horizontal partial derivatives by imposing the stress-free condition $\lambda(\mathrm{div} \, \mathbf{u})\delta_{iz} + \mu(\partial_i u_z + \partial_z u_i) = 0$ for $i = x, y, z$ at the local Earth's surface.

We selected 12 earthquakes that exhibit a pronounced seismic coda in the frequency band 1-3 Hz. Figure 1 shows the observed seismograms of one selected event of magnitude 4.3 at a distance of 35 km. They have been averaged over time windows of $\Delta t = 15$ s separated by 0.3 s. The first arrivals of $S$ and $P$ waves, as well as the exponential coda at long times, can easily be identified. In the regime of seismic coda all ratios fluctuate around a constant value until the signal-to-noise level is reached. The horizontal lines in Figs. 1c and 1d locate the time average. We emphasize that in this regime, the local elastic energy density decays by 4 orders in magnitude.

We interpret the observed time-independent energy ratios in the coda as a sign of equipartition. The time-

<figure><img src="image_2.png" /><figcaption>FIG. 1. Observed seismogram, bandpassed between 1 and 3 Hz, for event 11 at an epicentral distance of 35 km from the detection array and a magnitude of 4.3. (a) Linear plot of the bandpassed displacement measured as a function of time. (b) Semilogarithmic plot of the energy density. A distinction is made between kinetic energy ($K$, dashed line), shear energy ($S$), and compressional energy ($P$). (c) Linear plot of the energy ratio $S/P$. (d) Linear plot of the energy ratios $K/(S + P)$ (solid line) and $H^2/V^2$ (dashed line). The horizontal lines denote the estimated time average.</figcaption></figure>

dependent fluctuations around the mean for one event are attributed to the Gaussian fluctuations that ought to be uncorrelated in time and of the order of $1/\sqrt{\Delta t \Delta f/2}$ for a bandwidth $\Delta f$. In spite of the large variations of the sources in seismic magnitude and distance, Fig. 2 shows that the measured energy ratios are roughly the same for all sources. Especially the uncertainties in the ratios of potential energies (Table I) are small. The nonzero value of $\langle I \rangle$ indicates the presence of mode conversions at the Earth's surface.

We now attempt to understand the observed ratios by assuming that the coda is an equipartitioned random field. If $\mathbf{u}(\mathbf{r}, t)$ is the local, time-dependent displacement vector, it can be expanded into the eigenfunctions $\mathbf{u}_n$ of the elastic medium with eigenfrequencies $\omega_n$,

$$\mathbf{u}(\mathbf{r}, t) = \sum_n \varepsilon_n e^{-i\omega_n t} \mathbf{u}_n(\mathbf{r}). \quad (2)$$

In the presence of disorder, all modes get mixed and the $\varepsilon_n$ become time-dependent random variables. For small
<figure><img src="image_3.png" /><figcaption>FIG. 2. Observed time-averaged energy ratios $S/P$ (top left), $K/(S + P)$ (top right), $H^2/V^2$ (bottom left), and $I/(P + S)$ for all 12 events. The shadowed bar denotes the mean value with the standard deviation. The error bars denote the observed time-dependent fluctuations.</figcaption></figure>

disorder the eigenfrequencies $\omega_n$ are not significantly altered. The simplest way of expressing equipartition is that the average $\langle\epsilon_n\rangle = 0$ and that $\langle\epsilon_n\epsilon_m^*\rangle = \sigma^2(t)\delta_{nm}$. This means that the mode amplitudes are independent random variables with random phase but with equal variance [25]. This statement implies an "equipartition" of the elastic energy (1) among the different modes. The time dependence $\sigma^2(t)$ cancels in any energy ratio.

In a homogeneous unbounded medium the modes are plane waves with either $S$ or $P$ polarization. In that case $\langle I\rangle = 0$, $\langle K\rangle/\langle S + P\rangle = 1$, and $\langle S\rangle/\langle P\rangle = 2(\alpha/\beta)^3$ [19,21,26], in terms of the $P$ velocity $\alpha = \sqrt{(\lambda + 2\mu)/\rho}$ and the $S$ velocity $\beta = \sqrt{\mu/\rho}$, including a factor of 2 to acknowledge the twofold degeneration of the transverse $S$ waves. This ratio is recognized as the ratio of the density of states (DOS) for both modes. For $\lambda = \mu$ we get $\langle S\rangle/\langle P\rangle = 10.39$. The observed energy ratios do not agree with these quantitative predictions.

A more realistic model should take into account the mode conversions at the Earth's free surface, and at the interface between crust and mantle (the Moho) at roughly 30 km depth. The equipartition process within a distance

of a few shear wavelengths ($\lambda_s \approx 2$ km) from the surface is dominated by the first. We have calculated the exact elastic eigenmodes of a homogeneous plate, bounded by two free surfaces [25,27]. For a thickness $H > 10\lambda_s$, the equipartition values near the surface, calculated from Eq. (2), are independent of $H$, with surface values that agreed with previous work for the elastic half-space [25]. Figure 3 shows how the ratios $\langle S\rangle/\langle P\rangle$ and $\langle K\rangle/\langle P + S\rangle$ depend on depth if equipartition between all modes is assumed. The oscillations in Fig. 3, associated with the local density of states (LDOS) near the boundary, are well known in optics [28,29]. Our present study has access to the surface values only (Table I). The agreement between model and observations is remarkable. Table I also shows the theoretical values when either only surface waves or only bulk waves are considered. In both cases they fall outside the estimated error bar of the observational values. This leads to a picture of the coda consisting of an equipartitioned mixture of both Rayleigh and bulk waves, which is expected when the absorption time of Rayleigh waves exceeds the time to mode convert into bulk waves.

The observed value $H^2/V^2 \approx 2.5$ agrees less well (Table I). One possible source of discrepancy could be the presence of a "true" effective surface different from both the local geographical free surface and the local surface of the Earth, due to fluctuations in the slope of the Earth's surface. To investigate this option, the full kinetic covariance tensor $\langle\frac{1}{2}\rho\partial_t u_i\partial_t u_j\rangle/\langle K\rangle$ was calculated from the data. Its off-diagonal elements are less than 5% $\pm$ 5%, i.e., consistent with zero and thus with our equipartition model. However, the observed eigenvalues $\langle X^2\rangle: \langle Y^2\rangle: \langle Z^2\rangle = 0.27:0.45:0.28$ ($Y$ is north) disagree with the theoretical diagonal elements $\langle X^2\rangle: \langle Y^2\rangle: \langle Z^2\rangle = 0.32:0.32:0.36$. No frame can be found that makes the entire covariance matrix consistent with our equipartition model. Some other source of discrepancy must exist.

Alternatively, the discrepancy could be attributed to excess absorption of surface waves. It is readily shown that if the absorption time of a surface wave is less than the time to mode convert to bulk waves, its contribution to the steady state partition will be diminished. Table I shows that a small deficit of Rayleigh waves could bring $H^2/V^2$ into accord with theory without seriously

<table><caption>TABLE I. Observed averaged energy ratios (in the frequency range 1–3 Hz, with standard deviation) compared to the theory near the free surface, with $\lambda = \mu$ independent of depth. The third and fourth columns are calculated by considering both bulk and Rayleigh waves. The last columns are calculated by considering only surface Rayleigh waves or only bulk waves.</caption><thead><tr><th>Energy ratio</th><th>Data z = 0</th><th>Theory z = 0</th><th>Theory z = &infin;</th><th>Theory Rayleigh only z = 0</th><th>Theory Bulk only z = 0</th></tr></thead><tbody><tr><td>S/P</td><td>7.30 &pm; 0.72</td><td>7.19</td><td>10.39</td><td>6.460</td><td>9.76</td></tr><tr><td>K/(S + P)</td><td>0.65 &pm; 0.08</td><td>0.534</td><td>1</td><td>0.268</td><td>1.19</td></tr><tr><td>I/(S + P)</td><td>-0.62 &pm; 0.03</td><td>-0.167</td><td>0</td><td>-1.464</td><td>-0.336</td></tr><tr><td>H<sup>2</sup>/V<sup>2</sup></td><td>2.56 &pm; 0.36</td><td>1.774</td><td>2</td><td>0.464</td><td>4.49</td></tr><tr><td>X<sup>2</sup>/Y<sup>2</sup></td><td>0.60 &pm; 0.20</td><td>1</td><td>1</td><td>1</td><td>1</td></tr></tbody></table>
<figure><img src="image_4.png" /><figcaption>FIG. 3. Theoretical prediction for the depth dependence of the average energy ratios S/P (normalized to bulk value), K/(P + S) and I/(P + S), assuming that $\lambda = \mu$.</figcaption></figure>

compromising the others. In addition, an anisotropic absorption of Rayleigh waves might even explain the anomalous value of $\dot{X}^2/\dot{Y}^2 = 0.6$. In view of the geographical location of the array, a mountain area stretched in the direction E-W, it is quite possible that the component Y is privileged over the X direction.

In conclusion, we have observed a first principle of multiple scattering—energy equipartition—in seismograms recorded in Mexico. This strongly supports the interpretation of seismic coda as a genuine multiple scattering process, and excludes single scattering as an alternative explanation.

This work was supported by the program *Intérieur de la Terre* of the INSU/CNRS, the *Groupement de Recherches 1847 PRIMA*, the Mexican CONACYT Project No. J32308-T, and the NSF Grant No. 9988645. We would like to thank R. Maynard, B. Chouet, J. Scales, C. Lacombe, and S. K. Singh for many helpful discussions. We are grateful to E. Guevera, G. Espitia, V. Kostoglakov, and J. Pacheco for their help in the setup of the experiment.

<footnote>*Corresponding author.</footnote>

<footnote>[1] *Wave Propagation in Complex Media*, edited by G. Papanicolaou (Springer-Verlag, New York, 1998).</footnote>

<footnote>[2] *New Aspects of Electromagnetic and Acoustic Wave Diffusion*, edited by POAN Research Group (Springer-Verlag, Heidelberg, 1998).</footnote>

<footnote>[3] *Wave Diffusion in Complex Media*, edited by J. P. Fouque (Kluwer, Dordrecht, 1999).</footnote>

<footnote>[4] A. Yodh and B. Chance, Phys. Today **10**, No. 3, 34 (1995).</footnote>

<footnote>[5] *Light Scattering by Non-Spherical Particles*, edited by M. I. Mishchenko, J. W. Hovenier, and L. D. Travis (Academic, San Diego, 2000).</footnote>

<footnote>[6] P. Sheng, Z. Q. Zhang, B. White, and G. Papanicolaou, Phys. Rev. Lett. **57**, 1000 (1986); B. White, P. Sheng, Z. Q. Zhang, and G. Papanicolaou, Phys. Rev. Lett. **59**, 1918 (1987).</footnote>

<footnote>[7] K. Aki, J. Geophys. Res. **74**, 615 (1969).</footnote>

<footnote>[8] K. Aki and B. Chouet, J. Geophys. Res. **80**, 3322 (1975).</footnote>

<footnote>[9] H. Sato, J. Phys. Earth **25**, 27 (1977).</footnote>

<footnote>[10] M. Herraiz and A. F. Espinoza, Pure Appl. Geophys. **125**, 1569 (1997).</footnote>

<footnote>[11] M. Hoshiba, Phys. Earth Planet. Inter. **104**, 15 (1997).</footnote>

<footnote>[12] I. R. Abubakirov and A. A. Gusev, Phys. Earth Planet. Inter. **64**, 52 (1990).</footnote>

<footnote>[13] M. Hoshiba, Phys. Earth Planet. Inter. **67**, 123 (1991).</footnote>

<footnote>[14] H. Sato, Geophys. J. Int. **112**, 141 (1993).</footnote>

<footnote>[15] R. S. Wu and K. Aki, Pure Appl. Geophys. **128**, 49 (1988).</footnote>

<footnote>[16] L. Margerin, M. Campillo, N. M. Shapiro, and B. A. van Tiggelen, Geophys. J. Int. **138**, 343 (1999).</footnote>

<footnote>[17] N. Shapiro, M. Campillo, L. Margerin, S. K. Singh, V. Kostoglakov, and J. Pacheco, Bull. Seismol. Soc. Am. **90**, 655 (2000); M. Campillo et al., in *Wave Diffusion in Complex Media* (Ref. [3]).</footnote>

<footnote>[18] K. Aki and P. G. Richards, *Quantitative Seismology: Theory and Methods* (Freeman, San Francisco, 1980); T. Lay and T. C. Wallace, *Modern Global Seismology* (Academic Press, San Diego, 1995).</footnote>

<footnote>[19] R. L. Weaver, J. Acoust. Soc. Am. **71**, 1608 (1982).</footnote>

<footnote>[20] R. L. Weaver, J. Mech. Phys. Solids. **38**, 55–86 (1990).</footnote>

<footnote>[21] G. C. Papanicolaou, L. V. Ryzhik, and J. B. Keller, Bull. Seismol. Soc. Am. **86**, 1107 (1996).</footnote>

<footnote>[22] J. A. Turner, Bull. Seismol. Soc. Am. **88**, 276 (1998).</footnote>

<footnote>[23] B. A. van Tiggelen, R. Maynard, and A. Heiderich, Phys. Rev. Lett. **77**, 639 (1996); H. Stark and T. C. Lubensky, Phys. Rev. Lett. **77**, 2229 (1996).</footnote>

<footnote>[24] P. M. Morse and H. Feshbach, *Methods of Theoretical Physics* (McGraw-Hill, New York, 1953).</footnote>

<footnote>[25] R. L. Weaver, J. Acoust. Soc. Am. **78**, 131 (1985).</footnote>

<footnote>[26] L. Margerin, M. Campillo, and B. A. van Tiggelen, J. Geophys. Res. **105**, 7873 (2000).</footnote>

<footnote>[27] R. L. Weaver, J. Sound Vib. **94**, 319 (1984).</footnote>

<footnote>[28] H. Khosravi and R. Loudon, Proc. R. Soc. London A **433**, 337 (1991).</footnote>

<footnote>[29] E. Snoeks, A. Lagendijk, and A. Polman, Phys. Rev. Lett. **74**, 2459 (1995).</footnote>