EARTHEN DAM MONITORING USING PASSIVE SEISMIC: A FORWARD MODELING STUDY

. Internal erosion caused by infiltration and piping is the main cause of earthen dams’ failures. Traditional methods of inspection are not representative and are unable to cover the entire structure. The seismic interferometry emerges as an alternative for continuous monitoring earthen dams from the detection of variations in seismic velocities caused by water saturation. The objective of this study is the geophysical investigation of the Paranoá’s dam (Brasília, DF, Brazil) in hypothetical infiltration scenarios. To obtain a geological model that represents the expected structure, we performed an acquisition with active sources in the region and applied inversion techniques to retrieve the seismic velocities. We simulated acquisitions in a simple geological model of two horizontal layers for better understanding the behavior of the retrieved seismograms. The results show that more recording time and a wide receiver coverage improves the signal-to-noise ratio. We simulated a layer saturation by varying the velocities and noticed that the method was able to detect such variations. We noticed alterations in the reflection hyperbola asymptotes and in the first arrival time through the seismograms. The method also verified velocity changes in a more complex geological model. Our results suggest that the traffic energy is sufficient to image the dam’s structure, even not satisfying the theoretical diffuse wave field condition. Furthermore, we implemented a technique to monitor velocity variations using the semblance calculation. The method detected variations in the structure Vp in the order of 10%.


INTRODUCTION
Dams and levees are structures built from different materials and have distinct purposes.There is a need to continuously monitor these structures due to the serious consequences caused by failures, which damage the environment and can result in loss of life.The leading cause of earthen dam failure is the internal erosion caused by water (Flores-Berrones et al., 2011).The main phenomena that cause erosion are emergence, piping, and saturation of the medium.Traditional inspection methods cannot cover the entire area of interest and do not detect early signs of internal erosion (Albuquerque et al., 2019).Different geophysical methods are already applied to supplement these inspections, such as electrical resistivity (Sjödahl et al., 2008), ground penetration radar (Martini et al., 2016), and seismics (Osazuwa and Chinedu, 2008).
Background noise is generally separated from the seismic data to enhance signals of earthquakes or other signals of interest.However, there are growing applications for the parts of the signal that are usually treated as noise.For example, multiple processing expands the illumination of the subsurface (Verschuur and Berkhout, 2015;Berkhout, 2017), and diffractions can be used to image structures below the resolution limit (Moser and Howard, 2008;Maciel and Biloti, 2020).
Seismic interferometry is a technique capable of recovering the Green's function between two receivers recording environmental noise (Wapenaar et al., 2010).The method has been widely applied to monitor different s tructures s uch a s volcanoes (Brenguier et al., 2008;Olivier et al., 2019;Yates et al., 2019), fault zones (Wegler and Sens-Schönfelder, 2007;Brenguier et al., 2019), underground mines (Olivier et al., 2015), oil and gas reservoirs (Bakulin et al., 2007), landslides (Mainsant et al., 2012;Hussain et al., 2019), groundwater levels (Clements and Denolle, 2018), and buildings (Mordret et al., 2017).Authors are already studying the application of the method for dam monitoring, as done by Planès et al. (2016) and Olivier et al. (2017).Both studies were able to detect velocity changes in the structures and suggest using the technique as a non-invasive way to monitor and prevent accidents caused by leaks and ruptures.
The scope of this work is the geophysical investigation of the structure of the Paranoá's dam (Brasília, DF, Brazil).We calculated Vp and Vs profiles from an acquisition with active sources and applied two inversion methods for retrieving the seismic velocities.We also evaluated the feasibility of using seismic interferometry for continuous monitoring of the dam using passive seismic data.Moreover, we developed synthetic passive seismic acquisition models based on the inversions obtained in our initial investigation.
After the cross-correlation of the passive data, we developed monitoring panels based on the semblance calculation to monitor the velocity changes.

METHODS
We performed an acquisition with active sources at the foot of the Paranoá's dam (Figure 1).The line has 475 meters, where the receivers were spaced 5 meters, and the sources, 10 meters.The receivers are 14 Hz vertical geophones with a sampling rate of 0.025 ms.We used a 10 kg sledgehammer as source, with a stack of 3 shots.We applied the seismic refraction tomography through the time-term inversion to obtain a velocity profile of the P wave.Additionally, we applied the MASW (Multichannel Analysis of Surface Waves) technique (Park et al., 1999;Xia et al., 1999) to obtain a Vs profile.The Paranoá's dam is an earthen dam, which has compacted landfills and rocks for support (Fell et al., 2005).The core is made of impervious material (compacted clay); the upstream and downstream slopes are composed of boulders and rockfill.Figure 2 illustrates the typical structure of an earthen dam.

Seismic Interferometry
It is possible to recover subsurface reflectors by crosscorrelating the ambient seismic noise (Draganov et al., 2009).To exemplify the theory, we consider a onedimensional situation, as done by Wapenaar et al. (2010).Figure 3a shows an impulse radiated by a source at x = x s and t = t 0 , propagating along the x-axis.There are two receivers: x a and x b .Figure 3b shows the response observed in the receiver at x a , denoted as G(x a , x s , t), where G represents the Green's function.The recorded response consists of an unitary impulse at t a = (x a − x s )/c, so that G(x a , x s , t) = δ(t − t a ), where δ(t) is the Dirac delta.
The method involves cross-correlating the recorded signals of both receivers.As the path of the wave emitted at x s has x a and x b in common, the travel time between x a and x b is canceled in the cross-correlation process, resulting in the travel time between x a and x b , i.e., t b − t a = (x b − x a )/c. Figure 3d shows the impulse registered by the receiver in x b if the source was placed in x a .
One does not need to know the information about the position x s and the velocity c.The analysis is done similarly if the source impulse occurs at an arbitrary t = t s instead of t = t 0 , since the ray path passes through x a and x b .The cross-correlation of the registered responses at x a and x b is defined as G(x b , x s , t) * G(x a , x s , t).The time reversal of the Green's function of x b turns the operation in a correlation: (1) Neglecting the rigour that the distributions manipulation require, we can substitute the Delta functions on the right side of the equation: (2) which is the Green's function G(x b , x a , t).This expression indicates that the cross-correlation of the impulses registered at two receivers is an impulse at x b if x a was the emitting source.

Forward Models
The theory of seismic interferometry is based on some hypotheses: the wave field needs to be diffuse and the sources cannot be correlated.Furthermore, the theory does not predict factors such as wave attenuation and geometric scattering.As a consequence, synthetic models for interferometry tests should simulate these conditions.We used the finite difference modeller fdelmodc developed by Thorbecke and Draganov (2011) to generate the geological models and forward modeling the passive seismic data used in this study.In addition to simulating long-term passive acquisitions by modeling noisy sources scattered in the subsurface, the program is open source and can be used by the scientific community to study the different applications of the seismic interferometry.To better understand how the method works, we initially generated simple models with two horizontal layers (Figure 4), inspired by the work of Calvet ( 2018), and we varied two different parameters: the total recording time and the receivers' offset.Additionally, we simulated a hypothetical saturation process of the first layer and performed a time-lapse.
The layer velocities are based on the Vp and Vs values shown in Figure 5.The boundaries of the model edges were defined as absorbents; thus, the second layer extends infinitely in-depth.However, it was limited to 50-meter depth due to computational limitations.We defined the source wavelet as the Ricker.
Based on the results obtained after the seismic acquisition with active sources in the dam, we developed a geological model of greater complexity to represent the dam's structure.The maximum frequency of the receivers was 30 Hz, the sampling (dt) was 0.009 seconds, and the model discretization (dx and dz) was 1 meter.We used the Gardner's relation (Gardner et al., 1974) to determine the density of the layers: where ρ is the density in kg/m 3 and Vp is the P-wave velocity in m/s.We propose a method for continuous monitoring earthen dams by assessing the variation of seismic velocities through a sequence of coherence panels.The method is derived from classical velocity analysis techniques used in seismic exploration.We developed the panels to assess whether there is a change in the structure velocity.We used the suvelan program, from the Seismic Unix package (Cohen and Stockwell, 2013).For each retrieved seismogram after the seismic interferometry, the program calculates the semblance, as presented by Neidell and Taner (1971).The coherence values are synthesized in velocity analysis panels, henceforth referred to as "monitoring panels", where the x-axis represents velocity (m/s), and the y-axis, time (s).The region with higher coherence values is expected to indicate the estimated velocity and thickness of each model's first layer.

Active Seismic Data
After the Seismic Refraction Tomography and MASW inversions, we generated the velocity profiles (Vp and Vs) shown in Figures 5a and 5b.It is possible to observe in both models a superficial low-velocity layer (Vp between 200 and 500 m/s and Vs between 100 and 350 m/s), marked by dark blue, between 0 and 10-meter depth.We interpret it as unconsolidated soil.Then there is a higher velocity zone, marked by light green (Vp between 800 and 2200 m/s and Vs between 600 and 850 m/s), interpreted as a region composed of compacted clay material from the dam landfill.We noticed the beginning of a third layer, marked by yellow (Vp between 2200 and 2500 m/s and Vs between 900 and 1100 m/s), close to the depth limit of the models.Furthermore, it is possible to observe higher velocity structures at the edges of the Vp profile, above 2500 m/s.Such features may represent the dam's abutments (rock foundation), but we did not include this velocity in the synthetic models because we do not have more evidence that they represent the rock foundation.

Passive Synthetic Data
The number of noisy sources was set as 100 in all models in this study.We chose the number of random noisy sources empirically from the results presented by Thorbecke and Draganov (2011).The authors emphasize that the contribution of more noisy sources in models does not alter the signal-to-noise ratio (S/N) for the retrieved reflections in subsurface.A more in-depth study on the contributions of varying the number of sources in models can be found at the referred paper.Figure 6 shows the positions of the random noisy sources in the two-layered model, marked as black dots.The sources are randomly distributed in the model to simulate the white noise naturally produced by the constant vibration of the subsurface.We also varied the distance between the receivers, i.e., the total number of geophones used in the acquisition.We fixed the recording time as 30 seconds.The offsets used in the modeling were 1, 2, 4, 8, and 20 meters.The results are shown in Figure 8. Until the case where the distance is 8 meters (Figure 8e), we did not observe any noticeable changes in the hyperbola resolution or in the noise amplitude.However, when analyzing the result generated by the distance of 20 meters (Figure 8f), we noticed a deterioration in the hyperbola resolution.It is still possible to identify the top of the reflection, but the rest of the data could be misinterpreted as noise.This result indicates that it is necessary to determine a wide coverage of receivers in an acquisition for seismic interferometry monitoring purposes.
Figure 8: Retrieved seismograms after seismic interferometry using different offsets.
To simulate a saturation process in the superficial layer, we gradually changed the physical properties and performed a continuous recording.According to Knight and Nolen-Hoeksema (1990), soil saturation contributes to an increase in Vp and a decrease in Vs.Therefore, we set the initial values of Vp and Vs to 400 and 200 m/s, respectively, and we varied them up to the values of 600 and 100 m/s.The properties of the second layer remained constant, as shown in the previous analysis.
The velocity values used in these simulations are based on the results of the Vp and Vs profiles generated after the inversions (Figure 5.) The Vp/Vs ratio is considered as a key for detecting alterations in the presence of fluids in the rock pores (Brantut and David, 2019).It can also be used to calculate the Poisson's ratio σ, a parameter used to detect the increasing of fluids in rocks and sediments.We present the expected ranges of Poisson's ratio for soils and rock in Table 1.Comparing it with the values presented in Table 2, we demonstrate that our time-lapse experiment simulates a hypothetical case of the first layer saturation.We split the data into 30-second files to run the correlations and generate the seismograms retrieved by seismic interferometry.The results are shown in Figure 9.It is possible to notice the change in the slope of the hyperbola asymptotes as the Vp increases, i.e., the angle of the asymptotes increases relatively to the normal.We also noticed a decrease in the time of the first arrival, i.e., the top of the hyperbola approaches the surface.The imaging was able to detect the changes caused by the simulation of a saturation process in the first layer of the model.

Paranoá's Dam Model
Based on the velocities obtained by the inversions (Figure 5), we developed a two-layered irregular geological model to simulate passive acquisitions in the Paranoá's dam.The most superficial layer represents the unconsolidated soil, followed by a layer of compacted clay from the landfill.The geometry and properties of the model are shown in Figure 10.The dam crest is a highway, so most of the noise is generated by vehicular traffic.Planès et al. (2017) assessed the application of passive monitoring in a levee near to a highway.The authors analyzed the spectral content of the noise and concluded that only the traffic contributed to the spectrum, restricted to the superficial portion of the subsurface.Thus, we changed the positions of the noisy sources in the model to evaluate if it would cause any difference in the correlations.We simulated an acquisition with the sources distributed along the entire model (Figure 11a) and another in which we restricted the sources in the first 5-meter depth (Figure 11b).The traffic-generated seismic energy is part of the anthropogenic noise spectrum (Nakata et al., 2019), with frequencies in the 10-15 Hz range (Diaz et al., 2020;Riahi and Gerstoft, 2015).It is important to investigate the nature of the sources because of the dissipative and high attenuation characteristics of the near surface (Lindsey et al., 2020).Hence, we set the peak frequency of the noisy sources in our models to 15 Hz to simulate the energy transmitted by vehicular traffic on the highway at the crest of the dam.The recording time was defined as 120 seconds, and the spacing between the receivers was set as 1 meter.We changed the top boundary to a free surface for better simulating field conditions.Thus, the energy is reflected when arriving at the receivers, which generates multiples in the seismograms.The retrieved reflection response is shown in Figure 12.The amplitudes in Figure 12a are stronger than in Figure 12b; regions with larger offsets were not well resolved.This result reinforces that one must take care when determining the offset to ensure that all targets of interest are retrieved.Nevertheless, it also indicates that the traffic energy is sufficient to image the dam's structure.The results obtained in our synthetic experiments indicate that the vehicular traffic energy spectrum could be used for imaging and monitoring earthen dams, even not satisfying the theoretical imaging condition that random sources should surround the receiver (Wapenaar et al., 2010).The maximum depth of investigation depends on the noise signature, such as spectral signature and intensity, which suggests a more in-depth study of the nature of the noise in the dam.These points will be the object of a further research study.
We applied the monitoring panels to verify velocity variations generated by a hypothetical case of infiltration and saturation of the first layer of the dam's model.We generated a panel for each retrieved seismogram after the time-lapse, as shown in Figure 14.The panels show the regions of higher coherence values marked by red, and lower coherence marked by blue.The uncertainty of the method is related to the size of the high coherence region; a more dispersed region leads to a less accurate panel.As the model's surface was defined as a free-surface during the modeling, the presence of multiples in the seismograms generates a sparse region of coherence.However, it is possible to observe a region of high coherence near 1 s depth and 400 m/s velocity in Figure 14a, as expected.As the velocity gradually increases, we noticed the displacement of the regions of higher coherence in the panels; even though the velocities indicated by the greater coherence regions in the panels are slightly higher than the velocities of each model, it is possible to notice an increase in the velocity of the first layer.These minor discrepancies might be adjusted with better parameter tuning and suitable pre-processing methods.The method allows identify-ing variations in the structure velocities in the order of 10% (40 m/s), that suggests its ability to be applied in passive monitoring of dams.Further research will focus on studying the limits and potentialities of the monitoring panels.

CONCLUSIONS
The seismic survey with active sources in the Paranoá's dam provided the main parameters for the synthetic models that approximate the expected structure.The simple two-layered horizontal model allows a better understanding of how the technique works and provides information that can be applied in the field; e.g., more recording time and a wide coverage of receivers improves the S/N in the retrieved seismograms.The method is sensitive to velocity variations in simple models and in the more complex model created for the dam.The time-lapse experiments detected velocity variations typically caused by water saturation, infiltration, and piping.We noticed changes in the reflection hyperbola asymptote slopes and in the first arrival times through the seismograms of both models.Our modelings showed that the traffic-generated energy may be used for seismic interferometry monitoring purposes, even not satisfying the theoretical condition of a diffuse wave field.The proposed monitoring panels are sensitive to velocity variations.For the present case, we show a sequence of variations in the order of 10% that are properly detected by our method.The precision is associated with the size of the coherence blur.The blur can become straighter using suitable pre-processing methods, improving the resolution of the method.Seismic interferometry is a viable method that can be applied in the monitoring of earthen dams in order to complement information on the structure and variations in physical properties.Future studies will focus on applying the method on passive seismic field data and discussing the maximum investigation depth and the resolution of velocity variation.

Figure 1 :
Figure 1: Location of the seismic survey line at the Paranoá's dam.Note that the dam's crest is a highway.

Figure 3 :
Figure 3: Example of 1D seismic interferometry.(a) The source at x s transmits a plane wave that travels along the x-axis.(b) Response observed by the receiver at x a , i.e., G(x a , x s , t).(c) Similar to the previous situation, but in the receiver at x b .(d) Crosscorrelation of the recorded responses at x a and x b , corresponding if the second receiver registered an impulse sent from the first.Modified from Wapenaar et al. (2010).

Figure 6 :
Figure 6: The sources (black dots) are positioned at random positions along the model.

Figure 7 :
Figure 7: Retrieved seismograms after seismic interferometry using different recording times.

Figure 9 :
Figure 9: Retrieved seismograms after seismic interferometry with variations in the first layer velocities.

Figure 10 :
Figure 10: Physical properties of the geological model of the Paranoá's dam.

Figure 11 :
Figure 11: Noisy sources (black dots) randomly distributed along the geological model of the Paranoá's dam.(a) The sources occupy the entire model; (b) The sources are restricted to the first 5-meter depth.

Figure 12 :
Figure 12: Retrieved seismograms of the dam's models.(a) The sources occupy the entire model; (b) The sources are restricted to the first 5-meter depth.

Figure 13 :
Figure 13: Retrieved seismograms after seismic interferometry with variations in the first and second layer velocities in the dam's model.

Figure 14 :
Figure 14: Resulting monitoring panels for the retrieved seismograms after the time-lapse.The expected velocity/time positions are marked as a red "x".

Table 2 :
First layer model physical properties used in the time-lapse simulation: P-and S-wave velocities (m/s), Poisson's ration and density (kg/m 3 ).