# Passive retrieval of Rayleigh waves in disordered elastic media.

Éric Larose, A. Derode, D. Clorennec, Ludovic Margerin, Michel Campillo

► To cite this version:

Éric Larose, A. Derode, D. Clorennec, Ludovic Margerin, Michel Campillo. Passive retrieval of Rayleigh waves in disordered elastic media. Physical Review E: Statistical, Nonlinear, and Soft Matter Physics [2001-2015], 2005, 7204 (4) Part 2, pp.046607. 10.1103/PhysRevE.72.046607. hal-00109948

HAL Id: hal-00109948

https://hal.science/hal-00109948v1

Submitted on 22 Nov 2021

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 CC BY 4.0 - Attribution - International License
Passive retrieval of Rayleigh waves in disordered elastic media

Eric Larose,<sup>1,2,*</sup> Arnaud Derode,<sup>2</sup> Dominique Clorennec,<sup>2</sup> Ludovic Margerin,<sup>1</sup> and Michel Campillo<sup>1</sup>

<sup>1</sup>Laboratoire de Géophysique Interne et Tectonophysique (LGIT), Université Joseph Fourier, CNRS UMR 5559, Grenoble, France

<sup>2</sup>Laboratoire Ondes Acoustiques (LOA), Université Paris 7, CNRS UMR 7587, ESPCI, Paris, France

(Received 19 April 2005; published 14 October 2005)

When averaged over sources or disorder, cross correlation of diffuse fields yields the Green’s function between two passive sensors. This technique is applied to elastic ultrasonic waves in an open scattering slab mimicking seismic waves in the Earth’s crust. It appears that the Rayleigh wave reconstruction depends on the scattering properties of the elastic slab. Special attention is paid to the specific role of bulk to Rayleigh wave coupling, which may result in unexpected phenomena, such as a persistent time asymmetry in the diffuse regime.

DOI: 10.1103/PhysRevE.72.046607

PACS number(s): 43.20.+g, 46.40.Cd, 91.30.Fn

## I. INTRODUCTION

Whatever the type of waves involved, knowing the Green’s function of a heterogeneous medium is the key to many essential applications, such as imaging, communication, or detection. In the last twenty years or so, mesoscopic physics has intensely studied wave phenomena in strongly disordered media: weak and strong localization, radiative transfer, diffusion approximation, etc. [1]. However, the *exact* Green’s function of a complex medium is not easily tractable, and usually theoreticians study ensemble-averaged quantities such as statistical correlations of intensities or wave fields. Moreover, from an experimental point of view, it is not always possible to measure the Green’s functions of a complex medium because it requires controllable arrays of sources and receivers that do not perturb the medium. In the context of laboratory ultrasound, this is (nearly) routine; but in other fields of wave physics, such as seismology, distant sources are not controllable. In that respect, a good deal of publications have followed recent works by Weaver. He proposed to cross-correlate the diffuse wave fields obtained at two passive sensors and showed that it yields the elastic Green’s function between the two receivers, as if one of the receivers was a source [2,3]. That correlations performed on passive sensors should yield the wave travel times is not that new. This principle was applied to helioseismology in the 1990s, where it provided tomographic images of the Sun’s interior [4,5]. Beyond travel time reconstruction, Weaver’s experimental retrieval of exact Green’s functions was a real breakthrough.

Weaver’s work has been followed by different contributions. Some theoretical works are based on the ergodic approximation where time and ensemble average coincide. This approach is very useful for estimating the role of scattering in the correlation asymmetry [6,7]. When several sources are available, another possibility is to average correlations over sources without moving the receivers. This was done in underwater acoustics [8]. Sabra *et al.* [9] showed the possibility of recovering the entire Green’s function of the sea wave-

guide (especially the late contributions of multiples) and proposed a model for estimating the signal-to-noise ratio of the correlations. Derode *et al.* [10] proposed to interpret the Green’s function reconstruction in terms of a time-reversal analogy, and showed that all the benefits of time-reversal devices in multiple scattering media could be fruitfully applied to “passive imaging” from correlations [11,12]. This idea is based on the mathematical principle of the representation theorem [13], equivalently referred to as the Helmholtz-Kirchhoff theorem [14]: if a perfect series of receivers has been sensing the wave field for ages, then one can mathematically have access to the wave field anytime and anywhere in the area enclosed by the sensors. Very recently, this principle was fully developed and applied to elastic waves in open media by Wapenaar [15].

From the very start of Weaver’s work, the correlation principle was applied to seismic waves [17–19]. The most energetic part of the Green’s function between two seismic stations (Rayleigh and Love wave trains) was retrieved from passive correlations of either coda waves or records of seismic noise. As we mentioned earlier, obtaining impulse responses without a controllable source is of high interest in seismology since it gives the possibility to simulate very energetic and punctual earthquakes everywhere around the Earth, and to do imaging without a controllable source (passive imaging). Results are quite encouraging, although several theoretical problems remain unsolved. They concern issues very specific to the Earth. In particular, elastic sources (earthquakes) are naturally never arranged in such a way that they would perfectly surround a couple of seismic stations. In addition, for distant earthquakes, the propagation directions of incident bulk waves are mostly vertical, and in a vertically layered medium they would not couple with Rayleigh waves. Why, then, do we observe a Rayleigh wave train in the correlation of waves generated by bulk sources? The Green’s function retrieval is linked to equipartition, which is due to wave scattering. What therefore, is the influence of scattering within the Earth’s crust in that process?

Here we present experimental results obtained in the lab with laser-induced and laser-detected ultrasonic waves propagating in a heterogeneous slab mimicking the Earth’s crust. Many previous experimental articles applied the passive imaging technique to acoustic waves. Here we propose

<footnote>*Electronic address: eric.larose@ujf-grenoble.fr</footnote>
<figure><img src="image_2.png" /><figcaption>FIG. 1. $10^6 \rightarrow 1$ scaled analogy between the Earth's crust and our ultrasonic experiment. In the crust, the scattering mean-free path was estimated in Mexico [16] at 1 Hz: $\ell^*=30$ km. In the aluminum waveguide we numerically found $\ell^*=5.5$ mm at 1 MHz.</figcaption></figure>

to investigate the emergence of the Green's function in the correlations of elastic waves propagating in a solid heterogeneous medium with a free surface. Because field experiments are tedious and the natural environment is mostly unknown (especially the scattering properties), we propose to build an Earth's crust model at the scale $1/10^6$ (presented in Fig. 1). Waves will be sensed at ultrasonic frequencies at the free surface of an elastic open medium. In our experiment, surface waves are not initially excited. This is different from the work of Malcolm *et al.* [7], wherein ultrasonic Rayleigh waves were generated on the same surface they were measured. In addition, they used a finite cylindrical medium with possibly round-trip wave trains, whereas our experiment is conducted in a nearly open medium. In our configuration, we would not expect the Rayleigh wave to be reconstructed in the correlation, except if scattering is present and mode conversion between Rayleigh and bulk waves occurs. To verify this assumption and study the role of mode conversion, we used two samples of identical dimensions, one with scatterers the other without.

The next section describes the experimental setup and the propagation medium. Section III presents a short theoretical study of the scattering properties (scattering cross sections $\sigma$ of a single scatterer and the transport mean-free paths $\ell_p^*$ and $\ell_s^*$ of the heterogeneous sample). This study is supported by experimental measurements. In Sec. IV, the passive imaging technique is applied to records acquired at the free surface. Time and frequency analyses are proposed, and a brief discussion on the time symmetry of the correlations concludes the article.

## II. EXPERIMENTAL SETUP

The experimental setup is depicted in Fig. 2. It was designed to mimic the propagation and scattering of elastic waves through the Earth's crust, at ultrasonic frequencies (0.8–3.2 MHz). A duraluminum slab whose dimensions are roughly $10^{-6}$ those of the crust was used. Fifty-four cylindrical

<figure><img src="image_3.png" /><figcaption>FIG. 2. Experimental setup: The propagation medium is an aluminum block drilled with 54 vertical cylindrical holes (diameter: 4 mm). A Q-switched Nd:YAG laser shoots on the bottom side 100 pulses (24 ns duration and 280 MW/cm$^2$ intensity each). On the top side a heterodyne interferometer senses the vertical displacement of the aluminum/air interface. This measure is repeated at seven different locations ($X_0=0$ mm→$X_6=60$ mm) for each source position. Plasticine was stuck on the edges and the four lateral sides to mimic absorbing boundary conditions and to avoid the generation of Rayleigh waves by mode conversion at the edges.</figcaption></figure>

holes (radius $a=2$ mm) were drilled at random along direction $y$ so that the waves propagating through the slab could undergo multiple scattering. The two-dimensional (2-D) spatial Fourier transform of the hole positions was calculated and is almost perfectly flat in the range of ultrasonic wavelengths involved in the experiment, which confirms the absence of spatial correlation between the holes. The density of scatterers was $n=0.0105$ mm$^{-2}$. Ideally, the slab should have had infinite dimensions along $x$ and $y$. To approach this condition, we stuck a thick layer of dense plasticine on the lateral sides and edges of the duraluminum sample. This was aimed at creating absorbing boundary conditions and avoiding the generation of Rayleigh waves by mode conversion at the edges. The energy decay time $\tau_{abs}$ was initially found to be 23 000 $\mu$s. With the plasticine it decreased to 120 $\mu$s (see next section for a detailed discussion on the absorption time). We therefore simulated an open slab in the $x$ and $y$ directions with a free surface at the top. Rigorously, the Earth's crust is a waveguide that partially leaks energy through the Moho to the underlying mantle. To perfectly match this feature, we should have placed an infinite medium with a different impedance at the bottom side of the aluminum slab. Yet we think that our results and conclusions do not suffer too much from this omission: the central point of our setup is to mimic an elastic and scattering medium, infinite in the horizontal directions with a free surface on the top.

The source we employed to simulate earthquakes was a Q-switched Nd:YAG laser that shot 24 ns pulses at the bottom side of the slab (each pulse energy: 9 mJ). Two regimes of elastic wave generation are possible with a laser source [20]. When the surfacic intensity of the pulse is weak (<15 MW/cm$^2$), the surface is heated up locally and its dilation creates Rayleigh waves (thermoelastic regime). At higher intensity, the laser evaporates part of the metal. In this ablation regime, both Rayleigh and bulk waves are gener-
<figure><img src="image_4.png" /><figcaption>FIG. 3. Directivity pattern (linear scale) of the laser source. Rayleigh waves are not taken into account since they are absorbed at the edges. The experiments were carried out in the ablation regime, where bulk waves are much more energetic than in the thermoelastic regime. Compressional wave directivity: thick line; shear wave directivity: dotted line.</figcaption></figure>

ated. Theoretical radiation patterns are displayed in Fig. 3 for compressional (P) and shear (S) waves. Rigorously, such a source is not truly reproducible since the laser impact can damage the surface. To make the experiment as reproducible as possible while staying in the ablation regime, the shot intensity was no more than 280 MW/cm². A 1 ms record was acquired 100 times without any observable change. In the following experiments, for a satisfactory signal-to-noise ratio, each impulse response was averaged over 100 consecutive shots.

As to the detection of the free surface motion, it was achieved with a contactless and quasipunctual device: a heterodyne optical interferometer developed by Royer *et al.* [21] that has the advantage of a very broadband response (20 kHz–45 MHz) and a sensitivity of $10^{-4}$ Å/√Hz. It was mounted perpendicularly to the slab and thus provided us with the absolute vertical component of the free surface displacement (top side), with a fine spatial resolution; the size of the laser spot was ~100 µm, whereas the typical elastic wavelengths here are ranging between 1 and 10 mm. This is similar to seismology, where sensors are nearly punctual compared to the wavelengths considered (several kilometers at 1 Hz). However, seismic sensors usually provide time records of the three components of the displacement field. Here the interferometer measured only the vertical movements of the free surface.

In an elastic body, three different kinds of wave polarization are possible. Compressional (longitudinal) waves are analogous to acoustic waves in fluids (velocity $v_P = 6.32$ mm/µs in duraluminum). Shear (transverse) waves have two possible polarizations (velocity $v_S = 3.13$ mm/µs): one we call SV (vertical) in the x-z plane (see Fig. 1) and one SH (horizontal) in the x-y plane. SH waves have no contri-

<figure><img src="image_5.png" /><figcaption>FIG. 4. Top: typical wave form obtained through the scattering slab around 1 MHz. The sensitivity of the heterodyne interferometer is $10^{-4}$ Å/√Hz corresponding the optimal level of optical reflection on the sensed surface. After averaging over 100 records we reached the precision of $10^{-2}$ Å. Bottom: intensity averaged over several source/sensor positions.</figcaption></figure>

bution in the z direction and therefore will not be detected by the interferometer. In addition to bulk waves, surface waves exist but here only Rayleigh waves ($v_R = 2.9$ mm/µs) will be taken into account since the others cannot be detected (no vertical displacement). The shortest wavelength in the aluminum slab is 0.9 mm, which is much greater than the duraluminum alloy grain size. Since the orientation of the grains is random, we consider the alloy to be isotropic for elastic waves in the frequency band of interest. Scattering at the grain edge is presumably also negligible compared to scattering by the void cylinders.

The overall translational symmetry along y of both the free surface and the cylindrical scatterers avoid any coupling between SH mode and the other SV and P modes. Therefore SH waves will not be considered in our article and SV waves will be referred to as S waves. The wave propagation in our experiment will be treated as two-dimensional (and quasi-one-dimensional for the surface Rayleigh waves). Waves initially propagating in the y direction are rapidly absorbed by the plasticine and lost.

The laser source and the laser interferometer could be translated independently: 35 sources and 7 sensing positions were used during the experiment, providing us with a set of 35×7 impulse responses. A typical wave form is depicted in Fig. 4. Around 1 MHz it lasts nearly 800 µs, and shows a long diffusive decay comparable to the seismic coda. Due to the strong scattering on the cylindrical cavities, no top-bottom reflection was observed in the data. This confirms the highly diffusive nature of the propagation in the scattering slab. The relevant scattering properties are discussed and evaluated in the next section.

### III. WAVE SCATTERING AND TRANSPORT PROPERTIES

#### A. Scattering cross section of an empty cylinder

In order to evaluate the amount of scattering and mode conversion, we calculated the differential cross section
<figure><img src="image_6.png" /><figcaption>FIG. 5. Differential scattering cross sections of a cylindrical cavity calculated for a compressional (P) or shear (S) incident plane wave around 1.2 MHz (top) and 2.4 MHz (bottom). The scattered wave is either compressional or transverse. Each pattern is normalized by its maximum.</figcaption></figure>

$\partial\sigma/\partial\theta$ and the total scattering cross section $\sigma$ of a cylindrical void in an elastic medium excited by a compressional or shear plane wave. A brief description of the calculation is given in the Appendix. For a detailed derivation we refer to [22–24]. The differential scattering cross section gives the angular distribution of the scattered surfacic intensity, normalized by the incident surfacic intensity. The total elastic cross section is $\sigma = \int (\partial\sigma/\partial\theta)d\theta$. In 2-D, it has the dimensions of a length. This corresponds to the scattering strength of an object at a given frequency. In an elastic medium, mode conversion can occur and different cross sections must be considered. In the case of an incident compressional wave, they are $\sigma_{PP}$, $\sigma_{PS}$, and $\sigma_{P} = \sigma_{PP} + \sigma_{PS}$, respectively, for the P-to-P, P-to-S, and total P elastic cross sections. We also calculated the elastic cross sections for an incident shear wave (S): $\sigma_{SP}$, $\sigma_{SS}$, and $\sigma_{S} = \sigma_{SP} + \sigma_{SS}$.

The differential cross sections plotted in Fig. 5 have been computed at 1.2 MHz and 2.4 MHz frequencies. The elastic scattering sections are plotted in Fig. 6 for frequencies ranging from 0.1 to 200 MHz. In average in the frequency band of interest (0.8–3.2 MHz), we obtained $\sigma_{P}=9.2$ mm. This value is comparable to measurements by White [25]. The same calculations were conducted for an incident shear

<figure><img src="image_7.png" /><figcaption>FIG. 6. Elastic ($\sigma$) scattering cross sections calculated for shear (S) and compressional (P) plane waves impinging on a cylindrical void with diameter 4 mm. Between 0.8 and 3.2 MHz, scattering is stronger for shear waves. At high frequencies, the elastic cross sections tend to the limit of twice the geometrical diameter.</figcaption></figure>

wave; we found an average of $\sigma_{S}=12$ mm in the 0.8–3.2 MHz frequency band.

## B. Transport properties

When the elastic wave propagates through the aluminum slab drilled with holes, it undergoes multiple scattering. Let $\varphi(t)[\varphi_{0}(t)]$ be the vertical displacements sensed at the free surface through the scattering (homogeneous) medium. Classically, this field is split into two contributions: the coherent and the incoherent part. The coherent wave is the ensemble-averaged field $\langle\varphi(t)\rangle$ (averaged over disorder configurations; here, the cylinder’s positions). We underline the difference between the coherent and the ballistic wave (i.e., the first arrival). For a detailed discussion about scattering effect on coherent and ballistic waves, see [26]. Away from resonances, the coherent wave can roughly be thought of as an attenuated version of the direct wave front $\varphi_{0}(t)$. When there is no intrinsic dissipation, the energy of the coherent wave decays with the slab thickness $H$ as $e^{-H\ell}$, where $\ell$ is the elastic mean-free path. Assuming a dilute set of scatterers, the elastic mean-free path is simply related to the scatterers density $n$ and their elastic cross section $\sigma$:

$$ \ell = \frac{1}{n\sigma}. $$

From the theoretical scattering cross section calculated above, we find $\ell_{P}=10.3$ mm. In order to measure the mean-free path experimentally, we used two aluminum slabs of exactly the same dimensions. The first served as a reference and provided measurements of $\varphi_{0}(t)$ for different source-sensor positions. The second one was drilled with holes. By translating the source-receiver device along the slab, we achieved something very similar to a configurational averaging and measured the energy of the coherent wave $\langle\varphi(t)\rangle^{2}$. Between 0.8 and 3.2 MHz, we obtained $\ell_{P}=9\pm2.5$ mm from these experiments.

The intensity of the incoherent part was also studied. The time evolution of the averaged incoherent intensity $I(t) = \langle\varphi(t)^{2}\rangle$ is governed by another parameter: the transport mean-free path $\ell^{*}$. In an elastic body, transport quantities have been theoretically defined by [27]

$$ \ell_{P}^{*} = \frac{1}{n} \frac{\sigma_{S} - \sigma_{SS}^{*} + \sigma_{PS}^{*}}{(\sigma_{P} - \sigma_{PP}^{*})(\sigma_{S} - \sigma_{SS}^{*}) - \sigma_{PS}^{*}\sigma_{SP}^{*}}, \quad (1) $$

$$ \ell_{S}^{*} = \frac{1}{n} \frac{\sigma_{P} - \sigma_{PP}^{*} + \sigma_{SP}^{*}}{(\sigma_{P} - \sigma_{PP}^{*})(\sigma_{S} - \sigma_{SS}^{*}) - \sigma_{PS}^{*}\sigma_{SP}^{*}}, \quad (2) $$

with $\sigma^{*} = \int (\partial\sigma/\partial\theta)\cos(\theta)d\theta$. It was evaluated numerically: $\ell_{P}^{*} \approx \ell_{S}^{*} = 5.5$ mm. In an experiment, this parameter is very hard to measure with a reasonable precision. The coherent backscattering effect [28–31] (also referred to as weak localization) does give a direct estimation of the transport mean-free path, but our experimental configuration did not allow this special measurement since we could not place a laser sensor in the vicinity of the laser source. Yet we checked that the experimental intensity decay $I(t)$ gives an order of mag-
<figure><img src="image_8.png" /><figcaption>FIG. 7. Averaged intensity $I(t)$ and theoretical fit for scalar wave diffusion in a 2-D semi-infinite slab with thickness $H=24$ mm. $\ell^*=5.5$ mm and $\tau_{abs}=120$ µs.</figcaption></figure>

nitude for $\ell^*$ that is consistent with the theoretical value.

For the sake of simplicity, we propose a 2-D scalar wave model for $I(t)$ [32], under the diffusion approximation. In an infinite slab of thickness $H$ with perfect reflections on both sides, the averaged transmitted intensity reads

$$I(X,t) = I_0 \left\{ \frac{1}{2H\sqrt{\pi Dt}} + \sum_{n=1}^{\infty} \frac{(-1)^n}{H\sqrt{\pi Dt}} e^{-n^2 \pi^2 Dt/H^2} \right\} e^{-X^2/4Dt - t/\tau_{abs}},$$

with $\tau_{abs}$ the absorption time (taking into account the intrinsic absorption in the aluminum and the lateral leaking due to the plasticine) and $D$ the diffusion constant. $X$ is the lateral distance between source and receiver. This formula is obtained using a modal decomposition of the diffusion equation in the $z$ direction. The intensity $I(t)$ is a mix of compressional and shear waves, each mode traveling with its own parameters (velocity, $\ell^*$, diffusion constant) and inter-changing their energy through scattering events. In our experiment, $\ell_S^*$ and $\ell_P^*$ are of the same order. The diffusion constant was approximated by $D=[1/(1+2v_P^2/v_S^2)](v_P\ell_P^* + 2v_S(v_P/v_S)^2\ell_S^*) \approx 40$ mm$^2/\mu$s. This assumption is valid after a couple of mean-free times, when the equipartition regime [33] is set. Equipartition means that the density of compressional and transverse modes equilibrates. Considering the specific velocities of each mode [34], we infer that 80% of the energy is transported by S waves, and only 19% by P waves (and an additional 1% for surface waves). Hence our best fit (Fig. 7) of the intensity decay in the coda gives $\tau_{abs} = 120$ µs $\pm$ 10% and $\ell^*=5-20$ mm.

## IV. TWO-POINT CORRELATION OF DIFFUSE FIELDS

In this section we focus on the experimental reconstruction of the direct Green's function from "passive" correlations. The main idea is to correlate diffuse fields sensed at two different locations on the top side when a source generates bulk waves at the bottom. Since we record the vertical component of the surface displacements, the two-point correlations should simulate a vertical source at the surface that generates mainly surface waves. Indeed, the experimental

correlations we obtained reveal a wave packet that travels at the speed of a Rayleigh wave. We insist that our sources do not generate surface waves at the top side of the slab. Moreover, if a surface wave happened to be generated anywhere, it would be completely absorbed by the plasticine. Under these conditions no Rayleigh wave should travel on the top surface, and no Rayleigh wave should be passively retrieved by correlations. Why, then, should passive imaging give rise to a Rayleigh wave train in our experiment?

We propose first to examine the role of scatterers for the emergence of the direct Rayleigh wave front in the correlations. To that end, we separately correlated coda records obtained through two different aluminum slabs: the first drilled with holes, the second without. Each impulse response was lasting $\approx 800$ µs before reaching the noise level (see Fig. 4). We underline that these record lengths are far from the Heisenberg time (break time) at which the modes of the aluminum block would be resolved (here $T_H \approx 10^6$ µs) and correlations would naturally converge to the Green's function. This modal approach is irrelevant to our experiment. The records were correlated and averaged over the 35 available sources. For the scattering slab, this reads

$$\langle C_{ij}(\tau) \rangle = \sum_{S=1}^{35} \int_{t=0}^{t=600} \mu s \varphi(S, X_i, t) \varphi(S, X_j, t + \tau) dt,$$

where $X_i$ and $X_j$ are the sensors' points (running from $X_0=0$ mm to $X_6=60$ mm along the array). For the homogeneous slab,

$$\langle C_{ij}^0(\tau) \rangle = \sum_{S=1}^{35} \int_{t=0}^{t=600} \mu s \varphi_0(S, X_i, t) \varphi_0(S, X_j, t + \tau) dt.$$

To enhance the signal-to-noise ratio, each correlation is time symmetrized $[=\langle C(+\tau) \rangle + \langle C(-\tau) \rangle]$ and normalized by its maximum. Results are displayed in Figs. 8 and 9. A propagating wave front (traveling at the Rayleigh wave velocity $v_R=2.9$ mm/µs) is clearly visible in the presence of scatterers, whereas it does not appear in the homogeneous slab. We also summed the six normalized propagating peaks after having delayed each signal according to the Rayleigh wave travel time. The summation is displayed in the enclosed box on each figure. In the scattering slab, its amplitude nearly corresponds to the coherent addition of six pulses. In the homogeneous device, the amplitude of the summation is $\approx 2.5$ (incoherent addition of six uncorrelated fields). We conclude on the necessity of mode coupling due to scattering for the Rayleigh wave to emerge from the passive correlations of diffuse fields generated by bulk wave sources. This is especially relevant for applications to seismology. Therefore, it appears once again [11] that the role of scattering is crucial in passive imaging. Firstly, because of multiple scattering, at late times the equipartition regime can be attained whatever the source/receiver positions. Secondly, because of mode conversions due to the scatterers, a Rayleigh wave emerges from the passive correlations even though no Rayleigh wave was generated by the sources. Note that the bandwidth in the upper band record is a little
<figure><img src="image_9.png" /></figure>

FIG. 8. Green's function reconstruction for different pairs of receivers. Correlations are averaged over the 35 sources and filtered in the 0.8–1.6 MHz frequency band. On the left (a) are displayed the seven cross-correlations $\langle C_{ij}(\tau) \rangle$ in the diffusive aluminum plate. On the right (b) the seven cross-correlations $\langle C_{ij}^0(\tau) \rangle$ are obtained in an equivalent aluminum block without any hole, where the mode conversion possibly occurring at the edges was suppressed by the plasticine. The insets show the summation of the six wave fronts from $X_0=10$ mm to $X_6=60$ mm after a time reduction based on the Rayleigh wave speed (dotted line: $v_R=2.9$ mm/µs).

wider than in the lower band. This was done to compensate for the coda shortening ($<$600 µs) so that the product $T\Delta f$ is kept constant.

We can go a little further and catch a glimpse of the time-symmetry properties. We performed the correlations in two consecutive time windows: from 0 to 45 µs and from 45 to 600 µs, and did not time-symmetrize the correlations (Figs. 10 and 11). 45 µs is twice the time after which the diffuse energy spreads homogeneously along the array of receivers (length $L=30$ mm): $L^2/4D \approx 22$ µs in an open 2-D scattering medium. The time series were filtered in two fre-

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

FIG. 9. Same figure as Fig. 8, except that correlations are filtered in the 1.6–3.2 MHz frequency band. The signal-to-noise ratio in (a) is increased compared to the results obtained at lower frequencies.

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

FIG. 10. Evolution of the asymmetry in the reconstructed Green's function in the low-frequency regime (0.8–1.6 MHz), for early (a) and late (b) times. At early times, the waves are exciting scatterers on one side of the receiving array more than the other. The energy flux is clearly going from the source to the rest of the medium. After 45 µs, the field is equipartitioned and the scattering halo fills the array, leading to a more symmetric correlation.

quency bands: 0.8–1.6 MHz and 1.6–3.2 MHz.

In the first time window, from 0 to 45 µs, correlations are asymmetric in time. This means that the *causal* part ($\tau>0$) and the *acausal* part ($\tau<0$) of the correlations are different (see left parts of Figs. 10 and 11). In the *causal* part, a Rayleigh wave front is clearly visible, whereas noise is dominating the *acausal* part. This is due to the preferential direction of Rayleigh wave propagation (waves traveling from $X_0$ to $X_3$ in our experiment). There is a net flux of energy from $X_0$

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

FIG. 11. Evolution of the asymmetry in the reconstructed Green's function in the high-frequency regime (1.6–3.2 MHz), for early (a) and late (b) time windows. At early times, the waves are exciting scatterers preferably on one side of the receiving array. After 45 µs, bulk waves are expected to be equipartitioned, but not Rayleigh waves. Mode conversion to Rayleigh waves occurs mainly at one scattering cavity, leading to an anisotropic energy flux, and an asymmetric correlation.
to $X_1, X_2$, and $X_3$ (distances 10, 20, and 30 mm in Figs. 10 and 11). This flux is due to the uneven distribution of sources in comparison to the receiver couples $X_0-X_1$, $X_0-X_2$, and $X_0-X_3$: most of the sources are on the $X_0$ side. At the early times of the coda (from 0 to 45 µs), the diffusion regime is not yet attained.

Later in the coda, from 45 to 600 µs, the wave field in the bulk of the aluminum slab is very likely to be equipartitioned. In the low-frequency band (0.8–1.6 MHz), the time symmetry of the correlation is indeed restored [6,36] (see right part of Fig. 10): Rayleigh waves travel in all directions. Nevertheless and surprisingly, the asymmetry persists in the high-frequency band (from 1.6 to 3.2 MHz, see right part of Fig. 11). To interpret this observation, we have carefully studied the location of the scattering sources around the array. In the high-frequency band, the Rayleigh wavelength is ~1.5 mm. The generation of Rayleigh waves by scattering necessarily occurs in the first half-wavelength beneath the free surface [35]. In our scattering slab, one hole was nearly showing on the surface (position $X<0$); another one was 0.61 mm beneath (position $X>60$ mm), the others being located much deeper. Rayleigh wave trains are generated mainly by the hole nearest to the surface, then propagate along the array of receivers (from $X_0$ to $X_6$ point). These waves are almost not perturbed (attenuated) until they reach the edges of the slab and the absorbing plasticine. They contribute to a very clear propagating pulse in the positive part of the correlation. The weak coupling due to the deeper hole ($z=-0.61$ mm) on the $X_6$ side contributes to a smaller pulse propagating from $X_6$ to $X_0$ in the negative part of the correlations.

This interpretation is in agreement with observations in the low-frequency band (0.8–1.6 MHz), where the average Rayleigh wavelength is 3 mm (Fig. 10). At least five holes are present in the first half-wavelength and should cause significant scattering of Rayleigh waves and coupling between surface and bulk waves. This time, the holes are evenly distributed along the sensor array. In the late coda, correlations $X_0 \times X_1$, $X_0 \times X_2$, and $X_0 \times X_3$ are nearly symmetric. The time symmetry is obtained thanks to scattering by a symmetric distribution of scatterers. Under these conditions, a global equipartition among bulk and surface waves is guaranteed. Furthermore the reconstructed surface wave is strongly attenuated along its path because it senses many scatterers along the array. In addition, these scatterers contribute signals around $\tau=0$ in the correlations, which degrade the reconstruction. We think this interpretation explains why the symmetric wave front is much more noisy at low frequencies (Fig. 8), where lots of cavities are encountered in the first half-wavelength, than at higher frequencies, where the Rayleigh wave front can propagate freely (Fig. 10).

Finally, we comment on the possible misidentification of the waves that are reconstructed in our experiment. Indeed, if many scatterers are present within a wavelength, the overall wave velocity may be different (effective medium). In the heterogeneous plate, a reduced-speed shear wave might propagate with the same wave speed as a Rayleigh wave in the bulk aluminum plate (i.e., without cavities). In our experiment, the hole interspacing is 10 mm on average, which is larger than the largest ultrasonic wavelength. Thus, it is

reasonable to assume that the shear waves propagate in the aluminum plate with the same wave speed as in the bulk. The measured wave velocity is 2.9 mm/µs at all frequencies (from 0.8 to 3.2 MHz), and indeed corresponds to the velocity of the Rayleigh wave.

## V. CONCLUSION

In this paper were presented laboratory experiments of elastic wave propagation in heterogeneous media at ultrasonic frequencies. An aluminum slab was made quasi-infinite by the use of absorbing boundaries to mimic the Earth's crust, and scattering was obtained by drilling cylindrical cavities. A relatively simple theoretical model for wave scattering properties was proposed. Wave field generation and detection were achieved using contactless and quasi-punctual laser devices. The cross-correlation of diffuse fields was performed, allowing us to retrieve passively the Rayleigh wave between two sensors only when scattering was present. Without scatterers, and in the case of bulk wave generation and surface detection, no Rayleigh wave was reconstructed. This illustrates the role of scattering and mode conversion in the Green's function passive reconstruction.

Analyses for different time windows and frequency bands were conducted. Previous works in acoustics [37] and seismology [17] observed that the Rayleigh wave reconstruction was harder with increasing frequency. Here we found that the Rayleigh wave reconstruction was more efficient in the high-frequency band. In addition, we observed that even in the late coda, where waves are expected to be equipartitioned, asymmetry in the correlations may remain. Both observations are due to the very specific coupling between bulk and Rayleigh waves, which occurs if scatterers are present in the first wavelength beneath the free surface. We emphasize that equipartition of bulk waves does not always mean equipartition of surface waves. On the one hand, our experiments show the need of scattering to passively retrieve the impulse response between two sensors; on the other hand, they show that scattering occurring between the sensors degrades the reconstructed Rayleigh waves. The trade-off between the two effects should be further investigated. Although our experiment was designed for seismological applications, results should be applicable to other fields of wave physics in which both surface and bulk waves are present.

## ACKNOWLEDGMENTS

The authors wish to thank Richard Weaver and Julien de Rosny for fruitful discussions, and Xavier Jacob and Samir Guerbaoui for experimental help. This work was supported by the Groupement de Recherche CNRS "Imagerie, Communication et Désordre" (GdR IMCODE 2253) and the CNRS program "DyETI."

## APPENDIX A: CALCULATION OF THE SCATTERING CROSS SECTION OF A CYLINDRICAL CAVITY

What follows is a brief calculation of the wave field scattered by a cylindrical cavity insonified by a plane compressional wave. Three displacement potentials are relevant: $\Phi_i$
is the displacement potential of the incident P wave, $\Phi_s$ is the scattered P wave potential, and $\Psi_s$ is the scattered S wave potential. They can be expanded as

$$ \Phi_i(r,t) = \sum_{m=0}^{\infty} \varepsilon^m i^m J_m(k_p r) \cos(m\theta) e^{-i\omega t}, $$

$$ \Phi_s(r,t) = \sum_{m=0}^{\infty} A_m H_m^{(1)}(k_p r) \cos(m\theta) e^{-i\omega t}, $$

$$ \Psi_s(r,t) = \sum_{m=1}^{\infty} B_m H_m^{(1)}(k_S r) \sin(m\theta) e^{-i\omega t}, $$

where $J_m$ and $H_m^{(1)}$ are, respectively, the Bessel and Hankel functions both of first kind and of order $m$, and $\varepsilon_m$ is the Neumann factor. Taking into account the null traction condition at the surface of the cylinder allows the calculation of the $A_m$ and $B_m$ coefficients. Those coefficients are given in [22–24]. The scattering cross sections are

$$ \sigma_{P \to P} = \frac{2}{k_P} \left[ 2|A_0|^2 + \sum_{m=1}^{\infty} |A_m|^2 \right], $$

$$ \sigma_{P \to S} = \frac{2}{k_S} \left[ \sum_{m=1}^{\infty} |B_m|^2 \right], $$

$$ \sigma_P = \sigma_{P \to P} + \sigma_{P \to S}. $$

The corresponding differential scattering cross-sections are given by

$$ \frac{\partial \sigma_{P \to P}}{\partial \theta}(\theta) = \frac{4}{k_P} \left| \sum_{m=0}^{\infty} A_m i^{-m} \cos(m\theta) \right|^2, $$

$$ \frac{\partial \sigma_{P \to S}}{\partial \theta}(\theta) = \frac{4}{k_S} \left| \sum_{m=1}^{\infty} B_m i^{-m} \sin(m\theta) \right|^2. $$

[1] M. C. W. van Rossum and T. M. Nieuwenhuizen, Rev. Mod. Phys. 71, 313 (1999).

[2] R. L. Weaver and O. I. Lobkis, Phys. Rev. Lett. 87, 134301 (2001).

[3] O. I. Lobkis and R. L. Weaver, J. Acoust. Soc. Am. 110, 3011 (2001).

[4] T. Duvall, S. Jefferies, J. Harvey, and M. Pomerantz, Nature (London) 362, 430 (1993).

[5] J. Rickett and J. Claerbout, Sol. Phys. 192, 203 (2000).

[6] B. A. van Tiggelen, Phys. Rev. Lett. 91, 243904 (2003).

[7] A. E. Malcolm, J. A. Scales, and B. A. van Tiggelen, Phys. Rev. E 70, 015601(R) (2004).

[8] P. Roux, W. Kuperman, and the NPAL Group, J. Acoust. Soc. Am. 116, 1995 (2004).

[9] K. G. Sabra, P. Roux, and W. A. Kuperman, J. Acoust. Soc. Am. 117, 164 (2005).

[10] A. Derode, E. Larose, M. Tanter, J. de Rosny, A. Tourin, M. Campillo, and M. Fink, J. Acoust. Soc. Am. 113, 2973 (2003).

[11] A. Derode, E. Larose, M. Campillo, and M. Fink, Appl. Phys. Lett. 83, 3054 (2003).

[12] E. Larose, A. Derode, M. Campillo, and M. Fink, J. Appl. Phys. 95, 8393 (2004).

[13] K. Aki and P. G. Richards, Quantitative Seismology, 2d ed. (University Science Books, Sausalito, CA, 2002).

[14] D. Cassereau and M. Fink, IEEE Trans. Ultrason. Ferroelectr. Freq. Control 39, 579 (1992).

[15] K. Wapenaar, Phys. Rev. Lett. 93, 254301 (2004).

[16] L. Margerin, M. Campillo, and B. A. van Tiggelen, Geophys. J. Int. 134, 596 (1998).

[17] M. Campillo and A. Paul, Science 299, 547 (2003).

[18] N. M. Shapiro and M. Campillo, Geophys. Res. Lett. 31, L07614 (2004).

[19] N. M. Shapiro, M. Campillo, L. Stehly, and M. H. Ritzwoller, Science 307, 1615 (2005).

[20] C. Scruby and L. Drain, Laser Ultrasonics—Techniques and Application (Adam Hilger, London, 1990).

[21] D. Royer and E. Dieulesaint, Proc.-IEEE Ultrason. Symp. p. 527 (1986).

[22] Y.-H. Pao and C.-C. Mow, Diffraction of Elastic Waves and Dynamic Stress Concentration (Crane Russak, New York, 1995).

[23] J. J. Faran, J. Acoust. Soc. Am. 23, 405 (1951).

[24] Y. Liu, R.-S. Wu, and C. F. Ying, Geophys. J. Int. 142, 439 (2000).

[25] R. M. White, J. Acoust. Soc. Am. 30, 771 (1958).

[26] A. Derode, A. Tourin, and M. Fink, Phys. Rev. E 64, 036605 (2001).

[27] J. Turner, Bull. Seismol. Soc. Am. 88, 276 (1998).

[28] P.-E. Wolf and G. Maret, Phys. Rev. Lett. 55, 2696 (1985).

[29] M. P. Van Albada and A. Lagendijk, Phys. Rev. Lett. 55, 2692 (1985).

[30] A. Tourin, A. Derode, P. Roux, B. A. van Tiggelen, and M. Fink, Phys. Rev. Lett. 79, 3637 (1997).

[31] E. Larose, L. Margerin, B. A. van Tiggelen, and M. Campillo, Phys. Rev. Lett. 93, 048501 (2004).

[32] N. P. Trégourès and B. A. van Tiggelen, Phys. Rev. E 66, 036601 (2002).

[33] R. Hennino, N. Trégourès, N. M. Shapiro, L. Margerin, M. Campillo, B. A. van Tiggelen, and R. L. Weaver, Phys. Rev. Lett. 86, 3447 (2001).

[34] R. L. Weaver, J. Acoust. Soc. Am. 71, 1608 (1982).

[35] T. Maeda, H. Sato, and T. Nishimura, in Proceedings of the Workshop on Probing Earth Media Having Small-Scale Heterogeneities, Sendai, 2004, p. 90.

[36] A. Paul, M. Campillo, L. Margerin, E. Larose, and A. Derode, J. Geophys. Res. 110, B08302 (2005).

[37] R. L. Weaver and O. I. Lobkis, J. Acoust. Soc. Am. 117, 3432 (2005).