# Assessment of Risks Induced by Countermining Unexploded Large-Charge Historical Ordnance in a Shallow Water Environment-Part II: Modeling of Seismo-Acoustic Wave Propagation

Nathalie Favretto-Cristini, Fang Wang, Paul Cristini, Thierry Garlan, Olivier Morio, Diego Mercerat, Vadim Monteiller, Anne Deschamps, Éric Beucler

► To cite this version:

Nathalie Favretto-Cristini, Fang Wang, Paul Cristini, Thierry Garlan, Olivier Morio, et al. Assessment of Risks Induced by Countermining Unexploded Large-Charge Historical Ordnance in a Shallow Water Environment-Part II: Modeling of Seismo-Acoustic Wave Propagation. IEEE Journal of Oceanic Engineering, 2022, 47 (2), pp.374-398. 10.1109/JOE.2021.3111791. hal-03491720

HAL Id: hal-03491720

https://hal.science/hal-03491720

Submitted on 1 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>
# Assessment of Risks Induced by Countermining Unexploded Large-Charge Historical Ordnance in a Shallow Water Environment—Part II: Modeling of Seismo-Acoustic Wave Propagation

Nathalie Favretto-Cristini<sup>10)</sup>, Fang Wang, Paul Cristini<sup>10)</sup>, Member, IEEE, Thierry Garlan, Olivier Morio, E. Diego Mercerat, Vadim Monteiller, Anne Deschamps, and Éric Beucler

**Abstract**—The goal of this work presented in a two-companion paper is to pave the way for reliably assessing the risks of damage to buildings on the shore, induced by the detonation of large-charge historical ordnance (i.e., countermining) in variable shallow water environments. Here, we focus on the impact of the marine environment, more specifically the unconsolidated sedimentary layer, on detonation-induced seismo-acoustic wave propagation. We rely on a multidisciplinary cross-study including real data obtained within the framework of a countermining campaign, and numerical simulations of the seismo-acoustic propagation using a spectral-element method. We first develop a strategy relying on physical insights into the different kind of waves that can propagate in a coastal environment, to provide clues for a computational cost reduction. The geological surveys and the hydroacoustic measurements provide input data for the 3-D axisymmetric modeling of wave propagation. The numerical simulations, obtained for one specific source-receiver path with a variable sedimentary facies, are compared with the real seismic data induced by the detonation of a charge either on the seabed, or in the water column, and recorded on the coast. Numerical analysis sheds light on the strong interaction between surface waves and the sedimentary facies. The short-scale and deep

Manuscript received October 9, 2020; revised July 2, 2021; accepted September 5, 2021. This work was supported in part by the French National Research Agency (ANR), and in part by the DGA (French Ministry of Defense Procurement Agency) under Grant ANR-15-ASTR-0001 POSA. This work was granted access to the French HPC resources of TGCC under allocation gen 7165 and mam0305 and of CINES under allocation A0020407165 and A0030410305, both made by GENCI, and of the Aix-Marseille Supercomputing Mesocenter under allocations h025. (Corresponding author: Nathalie Favretto-Cristini.)

**Associate Editor:** N. R. Chapman.

Nathalie Favretto-Cristini, Paul Cristini, and Vadim Monteiller are with the Laboratory of Mechanics and Acoustics, Aix-Marseille Univ., CNRS, Centrale Marseille, LMA, 13013 Marseille, France (e-mail: favretto@lma.cnrs-mrs.fr; cristini@lma.cnrs-mrs.fr; monteiller@lma.cnrs-mrs.fr).

Fang Wang was with the Laboratory of Mechanics and Acoustics, Aix-Marseille Univ., CNRS, Centrale Marseille, LMA, 13013 Marseille, France. She is now with CGG, 91300 Massy, France (e-mail: alinewf.wang@gmail.com).

Thierry Garlan and Olivier Morio are with the Marine Geology Department, Shom/DOPS/STM/SEDIM, 29228 Brest Cedex 2, France (e-mail: thierry.garlan@shom.fr; olivier.morio@shom.fr).

E. Diego Mercerat is with Equipe Repsody, CEREMA Méditerranée, 06903 Sophia Antipolis, France (e-mail: diego.mercerat@cerema.fr).

Anne Deschamps is with the Observatoire de la Côte d'Azur, Géoazur, Université Côte d'Azur, CNRS, IRD, 06560 Sophia Antipolis, France (e-mail: anne.deschamps@univ-cotedazur.fr).

Éric Beucler is with the LPG - Université de Nantes, Université d'Angers, CNRS, UMR6112, 44322 Nantes Cedex 3, France (e-mail: eric.beucler@univ-nantes.fr).

Digital Object Identifier 10.1109/JOE.2021.3111791

sedimentary basins favor a local wave-amplitude amplification and a frequency shift toward the low-frequency domain. However, the seismo-acoustic waves are globally attenuated during their propagation because of intrinsic attenuation and geometrical spreading, which generally prevents any large damage to nearby buildings on the shore.

**Index Terms**—Numerical modeling, seismic risks, seismo-acoustic propagation, surface waves (SAWs), underwater explosion.

## I. INTRODUCTION

HISTORICAL unexploded ordnance (UXO) from World War II that is discovered almost every week close to the French coast, must be destroyed quickly after discovery to ensure the safety of divers and ships. The favored destruction method is countermining, i.e., the use of a high-order detonation conducted by exploding an additional donor charge placed adjacent to the UXO munition [1]. Depending on whether the UXO is safe to move, such countermining occurs at specific safe locations or at the location of discovery.

The risks for people in charge of the UXO countermining are well known by the mine warfare experts. In contrast, the possible consequences of underwater explosions on the marine environment and on the buildings located on the coast are more difficult to control because they are too complex to be reliably evaluated. Indeed, they depend mostly on the environment geology and the characteristics (weight and location) of the explosive charges and, hence, on the detonation-induced wave propagation. However, it would be useful to rely on one (or several) reliable metric(s) that would help develop a decision support tool for the risk assessment regarding inland infrastructures.

One of the main goals of our work, presented here in a two-companion paper, is then to pave the way for determining these metrics and for assessing in a reliable manner the seismic risks induced by the detonation of large-charge historical ordnance in variable shallow water environments. By “large charge” we mean here charges of between 80- and 680-kg TNT-equivalent weights, and by “shallow water” we mean water depth less than 50 m.

In Part I of this article [2], we suggested assessing the seismic risks through the seismic magnitude on the Richter scale induced
<figure><img src="image_2.png" /><figcaption>Fig. 1. Experiment site, namely the Rade d'Hyères in the Mediterranean Sea (south-eastern part of France). The two locations, labeled 3TY and 3TZ, are the locations where the UXO were detonated during the campaign of December 2018. The locations of the temporary seismological stations deployed along the coast are labeled PS01–PS19. The bathymetry is also indicated.</figcaption></figure>

by the underwater explosion (for more details, see [2]). This requires a proper understanding of the relationship between the charge weight and location, the propagation of the explosion-induced waves, and the marine environment beforehand. In other words, we have to understand how the seabed (in particular, the viscoelastic sedimentary layer with a varying thickness) and the water column with a varying depth influence the propagation of the seismo-acoustic waves that are generated by the UXO detonation and that reach the coast.

Part I of this article investigates this influence experimentally within the framework of a countermining campaign, called the POSA project, and conducted in December 2018 in the Mediterranean Sea in the Rade d'Hyères (south-eastern part of France, Fig. 1). Analysis of the data induced by the detonation of different large explosive charges (from 80- up to 680-kg TNT-equivalent weights) and recorded by acoustic and seismic recording systems, coupled to information provided by geological surveys, confirmed that the seismic magnitude globally increases with increasing charge weight, wherever the explosion occurs. However, the seismic magnitude we inferred from the seismic signals is half a unit smaller for the case of an explosion in the water column than for the case of an explosion on the seabed. Moreover, the seismic signals seem to be of greater amplitude about the sedimentary basin. However, these observations require a deeper interpretation, in particular regarding the influence of the sedimentary layer. Better understanding the physics of seismo-acoustic wave propagation along the variable shallow marine environment can help. This is the goal of this article (Part II).

Here, we rely on numerical modeling of the seismo-acoustic propagation in the Rade d'Hyères, developed using the *in situ* acoustical and geological measurements as input data. The cross-validation between the numerical simulations and the real seismic data can highlight the key parameters (with the associated uncertainties and limitations) that really impact the wave propagation induced by large-charge detonations. Particular attention is paid to the impact of the (acoustical and geometrical) properties of the sedimentary layer, since it is known to possibly

generate specific phenomena that alter the wave propagation characteristics [3]–[5]. This feedback may be very useful for the future development of a decision support tool for the seismic risk assessment regarding inland infrastructures.

Acoustic and seismo-acoustic wedge problems in ocean acoustics have been addressed in numerous papers over the past three decades. Different numerical methods have been proposed to solve these problems, including, e.g., normal mode and parabolic-approximation techniques, and finite-difference, finite-element, or boundary-element methods. Their own advantages and limitations regarding their capabilities, the accuracy of their results, and the computational cost have been discussed with the support of a set of well-defined benchmarks for range-dependent seismo-acoustic problems, e.g., in [6] and [7]. Despite still ongoing improvements, none of these methods is able to solve, with both a high accuracy and a low computational cost, real range-dependent ocean scenarios that generally include changes not only in water depth and sound-speed profile, but also in seabed properties with range and depth. The methods that can account for the variable geometry and the heterogeneity of the seabed (e.g., methods of finite-element type) can model accurately the real seismo-acoustic wave propagation, but generally at a very expensive computational cost.

Recently, a time-domain spectral-element method (SEM) (see [8]–[11]) has been shown to efficiently solve full-wave propagation problems in ocean acoustics [12], [13]. Beyond its capability of handling complex geometries and rheologies accurately as any finite-element technique, the time-domain SEM runs efficiently on very large computers, exhibiting scaling that is almost linear with respect to the number of CPUs or GPUs. This property can lead to a drastic reduction of the duration of numerical simulations compared to some more classical time-domain finite-element techniques.

The SEM is based upon a high-order piecewise polynomial approximation of the weak formulation of the wave equation. It combines the accuracy of the pseudospectral method with the flexibility of the finite-element method [14]. In this method, the wavefield is represented in terms of high-degree Lagrange interpolants, and integrals are computed based upon Gauss–Lobatto–Legendre quadrature. This combination leading to a diagonal mass matrix, which leads in turn to a fully explicit time-scheme, is very well suited to numerical simulations on parallel computers (e.g., [15]). The SEM is also particularly well suited to handling interface matching conditions and complex geometries, even geometries for which very distorted mesh elements can occur [16]. Very efficient perfectly matched layers (PML) can be used to limit the size of the studied domain (e.g., [17]), and thus, to reduce the required computational resources that may otherwise become prohibitive for large size of 3-D domains and high-frequency (HF) simulations. In addition, the effect of wave attenuation can be accurately taken into account, and it has been shown that the behavior of surface and interface waves is accurately modeled (e.g., [18]), which is particularly important in this study. More specifically, we use the SPECFEM software package (https://geodynamics.org/cig/software/specfem2d/).

In past years, the SEM has been successfully applied in geophysical exploration and in seismology (e.g., [19] and [20],
among others), and validated by marine seismic experiments carried out in controlled conditions at the laboratory scale [21], [22]. Recently, this method has been validated in underwater acoustics [23] as well. In particular, transmission loss maps for the classical 2-D upslope wedge problem with a fluid bottom, obtained by time-domain full-wave simulations (using the SPECFEM software package), were in perfect agreement with those reported in [24, Sec. 6.9.2] that were obtained by a split-step implementation of the Thomson-Chapman parabolic equation. In [23], an efficient procedure was derived to compute transmission losses and time dispersion maps for broadband signals and for elastic media as well. In the past decade, the time-domain SEM has been shown to accurately model T-wave generation and propagation in ocean [25], [26], and hence, enable deeper understanding of the seismic-to-hydroacoustic conversion at the sea bottom and the hydroacoustic-to-seismic conversion at the shore [27], [28]. The method has been successfully validated by real experiments [25]. Thanks to their accuracy, time-domain SEM simulations have also proven to be very useful for providing valuable sets of synthetic waveforms to complement the limited set of *in situ* recorded waveforms, and hence, to assist in the characterization of arrivals from explosion-generated hydroacoustic waves recorded at the T-stations [28]. It is worth noting that the configurations considered in these works are quite similar to our wedge problem.

Part II of this article mainly reports the results of the cross-validation between the numerical simulations of the explosion-induced seismo-acoustic wave propagation in the Rade d'Hyères and the real seismic data, with a focus on the interaction between waves and the sedimentary layer.

This article is organized as follows. First, we discuss in Section II, the specific physical and geometrical properties of the environment that are chosen to be the input data for the generation of the numerical model. Then, in Section III, we present the strategy for numerically modeling most realistically the seismo-acoustic wave propagation in the Rade d'Hyères. Indeed, as it will be shown, modeling the 3-D propagation in a large environment with very small wave velocities implies a huge computational cost that may be not yet affordable, even at the current largest supercomputers. Physical insights into the seismo-acoustic wave propagation can then provide valuable clues for a cost reduction. Section IV is devoted to the results obtained by a 3-D axisymmetric modeling of the seismo-acoustic wave propagation induced by an explosion on the sea bottom. Analysis of the wave propagation along one specific source-receiver path with a very interesting sedimentary facies sheds light on the contribution of the different types of waves that can be generated, and hence, highlights the impact of the properties of the sedimentary layers on the explosion-induced seismo-acoustic signals. Section V focuses on the influence of the charge location (on the seabed *versus* in the water column). Finally, Section VI concludes this article.

## II. INPUT DATA

Before the countermining campaign, acoustic and geological surveys were conducted in 2015 and 2016 to get information on

<figure><img src="image_3.png" /><figcaption>Fig. 2. Detailed 3-D map of the thickness of the sedimentary layer in the Rade d'Hyères and in its vicinity. The black dots correspond to the seismic stations deployed along the shoreline.</figcaption></figure>

the seabed geomorphology. Particular attention was paid to the characteristics of the sediments.

Ultrahigh-resolution sub-bottom profiling surveys (performed with the Kongsberg SPB120 and the iXblue Echoes 10 000 for the deep and shallow water parts of the Rade d'Hyères, respectively) provided a detailed 3-D map of the thickness of the sedimentary layer [29], evaluated every 20 m along a horizontal spatial grid with a thickness uncertainty of ± 1–2 m (see Fig. 2). The sediment thickness for the coastal area was evaluated from new core samples and from information already published in literature (e.g., [30] and [31]). Globally, the sedimentary layer is relatively thin (including at the detonation locations 3TY and 3TZ) and varies within the 1–5-m range. However, there is a sedimentary basin, close to the western part of the land coast, whose thickness ranges from 15 to 30 m.

Based on the geological and acoustical analyses of new core samples combined with archival results from previous core samples (e.g., [30] and [31]), the seabed is globally described as sandy sediments with a fine grain size, except locally and close to the western part of the study area where there is a mixture of fine sands and muds (see [2, Fig. 3]). It has to be pointed out here that the presence of muds that are absorbing and low-speed media, may have a strong effect on wave propagation, in particular on surface wave propagation [32].

The sediment properties have been measured either in laboratory conditions, or *in situ* in the Rade d'Hyères using a portative celerimeter [33] working in the frequency range 50–270 kHz. They have also been estimated from the database of the Shom institution composed of thousands *in situ* measurements performed in similar coastal areas. The density (*ρ*), P-wave velocity (*V*<sub>P</sub>) and attenuation (*α*<sub>P</sub>) in sediments are globally constant within the experiment site (see Table I). However, close to the western part of the study area they have greater values (see Table I and Fig. 3).

It is worth noting here that the intrinsic limitations of the tools we used for the acoustic measurements result in considering the sediment properties as constant with depth. However, unconsolidated sediments are known to exhibit significant vertical
<table><caption>TABLE I
RANGE OF THE MEASURED PROPERTIES OF THE MARINE ENVIRONMENT IN THE RADE D'HYÈRES THAT ARE CONSIDERED AS INPUT DATA FOR THE NUMERICAL MODELING</caption><thead><tr><th>Properties</th><th>Density<br>ρ<br>(kg/m³)</th><th>P-wave<br>velocity V<sub>p</sub><br>(m/s)</th><th>S-wave<br>velocity V<sub>s</sub><br>(m/s)</th><th>P-wave<br>atten. α<sub>p</sub><br>(dB/m/kHz)</th><th>S-wave<br>attenuation α<sub>s</sub><br>(dB/m/kHz)</th></tr></thead><tbody><tr><td>Media</td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Water</td><td>1 000</td><td>1 507 ± 2</td><td>-</td><td>-</td><td>-</td></tr><tr><td>Sedimentary layer</td><td>1550-2000<br>± 20</td><td>1625-1750<br>± 20</td><td><b>200</b></td><td>0.49-0.69</td><td><b>40</b></td></tr><tr><td>Rocky basement</td><td>2600-2650<br>± 30</td><td>4100-4450<br>± 50</td><td>2700-2910<br>± 200</td><td>-</td><td>-</td></tr></tbody></table>

Note that these properties vary within the study area, but they are assumed to be constant with depth. The associated measured uncertainties are also indicated. Note that the P-wave parameters of the sedimentary layer have been measured within the frequency range 50–270 kHz. The S-wave properties of the sediments (in bold) have been estimated (see text).

<figure><img src="image_4.png" /><figcaption>Fig. 3. Acoustical properties of the sediments in the Rade d'Hyères and in its vicinity: (a) P-wave velocity; (b) P-wave attenuation. The black dots correspond to the seismic stations deployed along the shoreline.</figcaption></figure>

gradients associated with the P- and S-wave velocities and attenuations (e.g., [34] and [36]).

In contrast to P-wave parameters, the S-wave parameters in the (unconsolidated) sediments could not be measured because we had no appropriate measurement system, e.g., like the one described in [37].

As a result, we chose to rely on the H/V spectral ratio method (i.e., the horizontal-to-vertical-components spectral ratio associated with surface waves), well-known in seismology (e.g., [38]), to estimate the range of possible values for the S-wave speed $V_s$, namely 180–210 m/s (see Appendix A). Hereafter, we assume that the shear-wave properties of the sedimentary basin are quite homogeneous, which makes sense from Fig. 2, and that $V_s = 200$ m/s. This value is in agreement with the values provided in the literature for marine sediments in shallow water conditions (e.g., [34], [39], and [40]).

Hereafter, we also assume that the S-wave attenuation $α_s$ is equal to 40 dB/m/kHz, and hence that it is linearly dependent on the frequency, which may not be the case for marine sediments in the low-frequency (LF) regime (see, e.g., [40] and [41]). The value of 40 dB/m/kHz is higher than those reported in the literature (e.g., [39]–[42]). However, this value is consistent with the fact that our HF measurements provide effective (and not only intrinsic) P-wave attenuations that are ten times higher than the intrinsic attenuations commonly found in the literature. It is worth noting that the thin layer of mud that overlies the sedimentary layer close to the coast (see [2, Fig. 3]), may also have an impact on the wave attenuation [37].

Fig. 4 summarizes information on the characteristics of the bedrock underlying the sedimentary layer. The nature of rocks is mainly extrapolated from nearby geological outcrops. The location in depth of the top of the rocky basement is also extrapolated from the measured thickness of the sedimentary layer. The rock properties are either estimated from geological maps and *in situ* observations, or estimated from measurements performed on pieces of rocks present in sedimentary core samples. Depending on the nature of rocks, the density, P- and S-wave velocities in rocks vary within the study area (see Table I). However, in the absence of information on their variation with depth, these properties are assumed to be constant with depth.

It is worth noting that no specific survey was conducted to get information on the sound-speed profile in the Rade d'Hyères. From a single measurement of the water temperature, and considering the salinity of the Mediterranean Sea, the sound speed was assumed constant within the Rade d'Hyères and equal to 1 507 m/s (see Table I).
<figure><img src="image_5.png" /></figure>

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

<figure><img src="image_7.png" /><figcaption>Fig. 4. Acoustical properties of the rocky basement in the Rade d'Hyères and in its vicinity. The black dots correspond to the seismic stations deployed along the shoreline.</figcaption></figure>

### III. STRATEGY FOR NUMERICALLY MODELING THE SEISMO-ACOUSTIC WAVE PROPAGATION

#### A. Numerical Issues

Ideally, we would like to reproduce all the physical phenomena underlying the seismo-acoustic propagation in the Rade d'Hyères (in particular the diffraction effects due to the bay curvature and to the presence of the islands). To do this, it would be necessary to carry out 3-D modeling:

1) in a domain of about $15 \times 15 \times 0.5 \text{ km}^3$;

2) in the frequency range 1–300 Hz (where 90% of the energy from the explosive sources is concentrated, see [2]);

3) with a fine mesh of minimum size about 0.5 m (if we consider that the lowest velocity in the environment is the S-wave speed of 200 m/s in the unconsolidated sediments), to capture the influence of the thinnest sediment thickness (namely, 1–2 m) on the wave behavior at high frequencies.

Considering the Courant–Friedrichs–Lewy condition, the numerical stability of the simulations would be ensured for a time step of $5 \times 10^{-5}$. Unfortunately, such 3-D simulations are simply out of reach at this time, even if we optimize the size of the mesh elements and use the most powerful super-computing center.

Here, the goal is then to define an appropriate strategy for a parsimonious, but efficient modeling of the seismo-acoustic wave propagation in the Rade d'Hyères. By parsimonious and efficient, we mean a modeling that preserves the most significant wave phenomena (and hence, an optimal modeling from the wave propagation viewpoint), while reducing as much as possible the computational cost.

Therefore, we choose to do the following.

1) Limit the study area (i.e., the entire Rade d'Hyères) to a corridor around one specific source-station path presenting a major interest for our purpose of understanding the influence of the sedimentary layer on the wave propagation, namely the 3TY-PS13 path of length 13 km (see Fig. 1). This leads to modeling the wave propagation in 2.5-D (i.e., 3-D axisymmetric), rather than in full 3-D, to limit computational costs. Considering only a part of the full physical domain results in the need for artificial boundaries in the simulations. To avoid spurious reflections from these boundaries, we use the PML technique of [17] implemented in the SPECFEM code.

2) Reduce the frequency range initially considered (i.e., 1–300 Hz) and focus only on the frequencies that may impact the infrastructures located on the coast, i.e., the frequency range 1–30 Hz. This frequency range corresponds mainly to that associated with the bubble signal (see Part I of this article).

#### B. Physical Insights

In addition, it is useful to identify the parts of the marine environment that have the greatest impact on waves, to properly adapt the domain mesh. To do this, let us first examine the following three main possible propagation regimes that can occur in a (very) shallow environment with a decreasing water depth like the Rade d'Hyères, namely:

1) the modal (acoustic) propagation in the water layer;

2) the propagation of waves/modes of interface/surface wave (SAW) types about the sedimentary layer;

3) the body (*P* and *S*) wave propagation in the elastic (sedimentary and rocky) layers (see Fig. 5).

These different types of waves may not propagate independently. Indeed, the thickness of the sedimentary layer together with the bathymetry profile that both vary from the source to the coast greatly influences the wave behavior. In shallow-water waveguides with a sloping bottom, it is well-known that the acoustic modes in water may interact with the interface waves supported by the ocean bottom [43], [44].
<figure><img src="image_8.png" /><figcaption>Fig. 5. Sketch of the three main propagation regimes (corresponding to different colors) that can occur in a (very) shallow water environment.</figcaption></figure>

Furthermore, depending on whether the source is placed on the sea bottom or in the water layer, the generation of the different types of propagation and the distribution of the associated wave energies may differ.

1) *Modal Acoustic Propagation:* The purpose here is not to recall in detail the basics of modal propagation in shallow water, as they can be found in reference books dealing with underwater acoustics (e.g., [24] and [45]). The aim is rather to identify the (finite) number of modes that can be generated by the explosive source located at 3TY, and their behavior during their propagation in the variable water layer toward the coast.

Indeed, in environments with a variable water depth the modes of order N have a cutoff frequency (or, equivalently, a cutoff height) below which they become evanescent and leak energy into sediments (see Fig. 5). For an ideal configuration, namely the Pekeris model [46], the cutoff frequency $f_c$ of the mode N is given by

$$f_c = \frac{V_w \left(N - \frac{1}{2}\right)}{2h \left[1 - \left(\frac{V_w}{V_{P_{sed}}}\right)^2\right]^{\frac{1}{2}}}. \quad (1)$$

Considering from our *in situ* measurements that the water depth at the location 3TY is $h \approx 46$ m, the sound speed in water is $V_w = 1507$ m/s, and that the P-wave velocity in sediments is $V_{P_{sed}} = 1625$ m/s, the cutoff frequency of the fundamental ($N = 1$) mode is then $f_c \approx 22$ Hz, while the cutoff frequency for mode 2 is $f_c \approx 65$ Hz.

In the frequency range of interest here, i.e., 1–30 Hz, only the fundamental mode can be generated by the explosive source and its spectral content is very narrow (between 23 and 30 Hz). During its propagation toward the coast, in a thinner and thinner water layer, this mode will lose progressively its different spectral propagative components when they reach their own cutoff height $h_c$ [easily calculable from (1)]. Since its spectrum is very narrow, the fundamental mode will not be very dispersive. We can then reasonably consider that this mode is more or less

harmonic, and that it will propagate toward the coast up to a water depth of about 33 m (cutoff height for the frequency 30 Hz). In other words, the modal acoustic propagation should no longer take place for a water depth less than 30 m.

Practically, this means that we could erect a “water wall” in the water domain mesh at the location of the seismic profile where the water depth is 30 m (corresponding approximately to the range 5500 m; see the red-dotted vertical line in Fig. 6). This would avoid specific issues related to the very elongated shape of the mesh elements in the water domain very close to the coast, while preserving the physics of acoustic propagation in the water layer.

It is worth noting here that, when the water depth $h_w$ is such that $h_w \le \lambda_w/10$ ($\lambda_w$ being the wavelength in water), the water layer has been shown to behave only like a load, and it is no longer necessary to mesh the area occupied by water [47]. Therefore, in our case, in theory it would be no longer necessary to mesh the water area below $h_w = 5$ m (see the green-dotted vertical line in Fig. 6). However, this remains to be refined according to the behavior of the SAWs.

2) *SAW Propagation:* The type of SAWs involved in the propagation from the source to the coast mainly depends on the combination of several parameters, including the presence/absence of the water layer and/or the sedimentary layer and their perception by the waves, depending on the ratio between their thickness and the wavelengths. To optimize the numerical modeling and the associated domain meshing, and considering the variable geometry of the marine environment in the Rade d’Hyères, together with the selected frequency content of the explosive sources, it may be useful here to look at the conditions of generation, evolvement, and transformation of the SAWs during their path from the source to the receiver, together with their associated characteristics. To do this, the wave propagation problem can be considered in terms of HF and LF components for each kind of waves. Potentially, the following five main types of SAWs may propagate in the Rade d’Hyères.

1) A Stoneley–Scholte wave at the water/sediments interface. This case implicitly considers that the thickness of both the water layer and the sedimentary layer is much greater than the SAW wavelength (i.e., mainly in the HF limit). The phase speed $V_{Sch}$ of the Stoneley–Scholte wave at the water/unconsolidated sediments interface is typically such that $V_{Sch} \approx 0.85 V_{S_{sed}} < V_w$. Therefore, here $V_{Sch} \approx 170$ m/s if we assume that the S-wave speed in sediments $V_{S_{sed}}$ is 200 m/s. The energy associated with the Stoneley–Scholte wave is partially trapped in an $1-\lambda_{Sch}$-thick waveguide around the interface, with a roughly equal distribution in water ($\sim 50\%$) and in sediments ($\sim 50\%$) (e.g., [48] and [49]). Since the unconsolidated sediments are viscoelastic, the SAW is dispersive and its energy may be strongly attenuated during its propagation path. It has to be pointed out here that, as the sediments are unconsolidated, no wave of the Leaky–Rayleigh type can propagate at the water/sediments interface [50].

2) A Stoneley–Scholte wave or a Leaky–Rayleigh wave (also named the pseudo-Rayleigh wave) at the water/rocky basement interface. This case may be fulfilled for a large
<figure><img src="image_9.png" /><figcaption>Fig. 6. Two-dimensional bathymetry and the depth profile of the line 3TY-PS13, obtained from measurements by a sub-bottom profiler. Note the difference in scale between the horizontal and the vertical axes. The charge detonation takes place at the offset 0. According to the experiment conditions, the explosive charge is located either on the sea bottom (~ 46 m deep), or in the water column (at ~ 11 m from the sea surface). Note that the changes in the media characteristics, in terms of geometry and physical properties, all along the source 3TY-receiver PS13 path are accounted for in the mesh definition. The vertical dashed lines indicate the location of the virtual “water walls” with a depth of 30 m at the range ~ 5500 m (red dotted line) and with a depth of 5 m at the range ~ 12260 m (green dotted line), respectively. Beyond these walls, it is no longer necessary to mesh the area occupied by the water to the coast to preserve the physics of the modal acoustic propagation (red wall), and the physics of both the modal acoustic and the SAW propagation (green wall), in the water layer.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Fig. 7. Mesh of the marine environment close to the coast. The green part corresponds to the water column, the blue part to the sedimentary layer and to a part of the rocky basement, and the yellow part to the rocky basement. The pink part illustrates the transition zone between the real marine environment (on the left side) and the PMLs located on the right side (shown in mustard color).</figcaption></figure>

water depth and a thin sedimentary layer (i.e., mainly in the LF limit). The Stoneley–Scholte wave at the water/rock interface propagates in a nondispersive way and faster than that for the previous case, since now $V_{Sch} \approx 0.9 V_w < V_{S \text{ rock}}$. Therefore, here $V_{Sch} \approx 1350$ m/s as the sound speed in water $V_w$ has been measured equal to about 1507 m/s. The energy associated with the Stoneley–Scholte wave is still concentrated at the interface and decays rapidly with distance from the interface. The penetration is still roughly equal to $1 \lambda_{Sch}$. Nevertheless, the energy is fairly concentrated in water (~ 80%) (e.g., [48]–[49]).

Besides the Stoneley–Scholte wave, a Leaky–Rayleigh wave can also propagate at the water/rock interface [51]. Its speed $V_{LR}$ is such that $V_w < V_{LR} \approx 0.9 V_{S \text{ rock}} < V_{S \text{ rock}}$ (e.g., [52]–[53]). So, typically, $V_{LR} \approx 2500$ m/s

here. The energy associated with the Leaky–Rayleigh wave is leaking in water during the SAW propagation, while decaying with distance from the interface in the rocky basement. The wave is then dispersive and its energy is attenuated during its propagation path.

3) Stoneley–Scholte modes in the water /unconsolidated sedimentary layer/rocky basement configuration. In the 1–30-Hz range, there are indeed frequencies for which the sedimentary layer is thin enough for waves to be influenced by the bedrock and, at the same time, thick enough to support modal propagation. It has to be mentioned that the water layer also may be seen as a finite-thickness layer or as a semi-infinite layer, depending on the frequencies.

In the HF limit (i.e., for $h_{sed}/\lambda_{Sch} \gg 1$), the speed $V$ of the fundamental mode tends toward the speed $V_{Sch}$ of the
<figure><img src="image_11.png" /><figcaption>Fig. 8. Signals induced by the detonation of a 400-kg TNT-equivalent charge located on the seabed at 3TY (corresponding to the explosion S6), and recorded at the station PS13 located on the coast. They represent the horizontal component of the velocity (colinear to the source-receiver direction), together with their spectrograms. (a) Simulated data. (b) Real data. Note that as the real source function was not accounted for in the numerical simulations, the amplitude of the simulated data is only indicative and cannot be compared with the real data amplitude. The amplitudes of the time signals have been normalized.</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>Fig. 9. Focus on the first 20 s of the time series shown in Fig. 8. Note that the amplitudes of the signals have been normalized.</figcaption></figure>
Stoneley–Scholte wave at the water/sediments interface, while in the LF limit (i.e., for $h_{sed}/\lambda_{Sch} \ll 1$) and close to the cutoff frequency, $V$ tends toward the S-wave speed $V_{Ssed}$ in the sediments [54]–[57], [42]. Below the low cutoff frequency, the sedimentary layer is no longer seen. Provided the water layer is thick enough, the wavespeed then tends toward either the speed $V_{Sch}$ of the Stoneley–Scholte wave, or the speed $V_{LR}$ of the Leaky–Rayleigh wave, propagating at the water/rock interface [38]. Otherwise, the wavespeed tends toward the speed $V_R$ of a Rayleigh wave propagating at the “free” surface of the bedrock [55].

The energy associated with the fundamental Stoneley–Scholte mode is localized very close to the water/sediments interface and is generally fully concentrated within the sedimentary layer [57]. In addition to the fundamental mode, higher Stoneley–Scholte modes may also appear. Fairly insensitive to the structure of the surficial sediments, they are more influenced by the deeper layers as their energy may deeply penetrate in the bedrock [32], [56].

4) A Rayleigh wave at the air/sediments or air/rock interface. In the very LF limit, and considering the ratio between the water depth and the wavelengths associated with the SAWs potentially involved, the water layer may indeed not be seen by these waves. We can then consider that SAWs of the Rayleigh type can propagate either at the “free” surface of the bedrock (when the sedimentary layer is thin enough), or at the “free” surface of the sediments (when they are thick enough). In both cases, a good approximation of the phase speed $V_R$ of the Rayleigh wave is given by [58]

$$V_R = \frac{0.862 + 1.14\nu}{1 + \nu} V_S \quad (2)$$

with

$$\nu = \frac{1}{2} \left( 1 - \frac{V_S^2}{V_P^2 - V_S^2} \right).$$

So, typically, $V_R \approx 190$ m/s or 2500 m/s here, following the S-wave speed under consideration. Moreover, the associated energy is almost fully concentrated (95%) in an $1-\lambda_R$-thick waveguide just below the interface [59]. It has to be pointed out that the Rayleigh wave propagating at the sediment surface is dispersive and its energy may be then strongly attenuated during its propagation path.

5) Rayleigh–Sezawa modes in the air/sedimentary layer/rocky basement configuration. As the sediments are unconsolidated and

$$\frac{V_{S \text{ sed}}}{V_{S \text{ rock}}} \ll \frac{1}{\sqrt{2}}$$

our case is similar to the “loading” case [60]. Therefore, regardless of the frequency, there is a fundamental Rayleigh mode whose speed $V_{R1}$ tends, in the limit of LF, toward the speed of a Rayleigh wave at the free surface of the bedrock (so typically, around 2500 m/s here), and in the limit of HF, toward the speed of a Rayleigh wave

at the free surface of the sediments (so typically, around 190 m/s here). The energy associated with this Rayleigh mode, also known as the $M_1$ wave [61], [62], is mainly trapped below the sediment surface and decreases exponentially with depth. The penetration depth in the bedrock decreases with increasing sedimentary layer thickness for a given frequency [60].

Depending on the considered frequency range, there can also be an infinite number of higher order modes, known as Sezawa modes [62], [60] with a low cutoff frequency. Their phase speed lies between $V_S$ rock (at the LF limit) and $V_S$ sed (at the HF limit). As illustrated in [62] for the mode 2, also known as the $M_2$ wave [61], [62], the energy associated with these modes has a fairly significant penetration depth in the bedrock in particular close to the cutoff frequencies, while for HF the energy is rather concentrated in the sedimentary layer with a penetration depth of one $\lambda$ in the basement. The sinusoidal nature of the displacement amplitudes of the Rayleigh–Sezawa modes in the sedimentary layer is much more pronounced for the higher Sezawa modes at large $h_{sed}/\lambda$ values [60].

3) Body Wave Propagation: Body ($P$ and $S$) waves can propagate in the elastic media such as the sediments and the rocky basement, especially if the explosive source is located on the seafloor. Since the distances between the source and the receivers on the coast are rather small (at least for LF wave components), the penetration of the LF components of the bulk waves in the rocky basement is expected to be shallow. Although $P$- and $S$-wave velocity profiles with depth could not be evaluated from measurements, our numerical modeling may be able to reproduce in a quite realistic way the direct LF propagation, related to this kind of waves, between the source and the receiver.

### C. Chosen Seismic Line Path and the Associated Mesh Generation

For most of the experiments, the explosive charge was detonated in the bay at the location labeled 3TY in Fig. 1 (for more details, see the companion paper Part I). The water depth at the location 3TY was measured approximately equal to 46 m.

For our purpose of investigating the influence of the sedimentary layer on the seismo-acoustic wave propagation, we choose to consider for the numerical modeling the seismic line that connects the source located at 3TY to the station located at PS13 on the western part of the Rade d'Hyères (see Fig. 6). Indeed, this line of approximately 13-km length passes through relatively thick sedimentary layers, namely one 7-m-thick reservoir around the offset 3500 m and another one with a thickness of up to 30 m that covers the offsets 8500–12 440 m, the offset 12 440 m being the location of the seismological station PS13.

In Fig. 6, it is worth noting the difference in scale between the horizontal and the vertical axes to prevent the first-sight illusion that the sea bottom is rough and strongly sloped.

Obtaining an equivalent pressure signal induced by large explosions at the fictitious boundary between the hyperacoustic near-field (nonlinear regime) and the linear acoustic far field (linear regime) [63], [64] is unfortunately far from being simple.
<figure><img src="image_13.png" /><figcaption>Fig. 10. (a) Location of the fictitious horizontal receiver array (colored dots) along the line 3TY-PS13; (b) simulated seismograms related to the horizontal component of the velocity; (c) associated spectra; and (d) focus on the first 10 s of the seismograms represented in (b). Note that, for the sake of a better visualization, the signals and the spectra obtained at the different offsets are represented in (c) and (d) at different amplitude scales, while the seismograms are represented in (b) at the same scale. In (d), the label “SS w/r” denotes a Stoneley–Scholte wave at the water/rock interface, the label “SS w/s” denotes a Stoneley–Scholte wave (or mode) at the water/sediments interface, and the label “Rayl r” a Rayleigh wave at the surface of the rocky basement, while the labels “P<sub>dir</sub> w”, “P<sub>dir</sub> s”, and “P<sub>dir</sub> r” denote the direct P-wave in the water, the sediments and the rock. The label “S<sub>dir</sub> r” denotes the direct S-wave in the rock.</figcaption></figure>
<figure><img src="image_14.png" /><figcaption>Fig. 10. (Continued.) (a) Location of the fictitious horizontal receiver array (colored dots) along the line 3TY-PS13; (b) simulated seismograms related to the horizontal component of the velocity; (c) associated spectra; and (d) focus on the first 10 s of the seismograms represented in (b). Note that, for the sake of a better visualization, the signals and the spectra obtained at the different offsets are represented in (c) and (d) at different amplitude scales, while the seismograms are represented in (b) at the same scale. In (d), the label “SS w/r” denotes a Stoneley–Scholte wave at the water/rock interface, the label “SS w/s” denotes a Stoneley–Scholte wave (or mode) at the water/sediments interface, and the label “Rayl r” a Rayleigh wave at the surface of the rocky basement, while the labels “P<sub>dir</sub> w”, “P<sub>dir</sub> s”, and “P<sub>dir</sub> r” denote the direct P-wave in the water, the sediments and the rock. The label “S<sub>dir</sub> r” denotes the direct S-wave in the rock.</figcaption></figure>

This is currently one of our undergoing works. Therefore, as input data for the simulations, we consider an omnidirectional source of Dirac delta function type. In this way, we can obtain the transfer function associated with the wave path between the source and the station under consideration, and then the simulated signals at the station by subsequently applying an appropriate lowpass filter (namely, an eighth-order Butterworth filter) in the 1–15-Hz range to the simulated data. Once the signal associated with the explosion of the different charge weights is perfectly known, it will be straightforward to obtain “more realistic” simulated signals at the station by convolving this source signal with the transfer function.

The creation of a nonstructured mesh made of quadrangles for a complex geometry is a difficult task. In this work, we use Cubit/Trelis software [65] to mesh the computational domain of interest. In the case of a nonstructured mesh, we need to keep in mind the subsequent computational cost and the accuracy by paying attention to the following three points.

1) The first point concerns the element size that must be small enough to accurately model the propagation of the highest frequencies.

2) The second point concerns the size of the different elements within one material that should be as equal as possible (depending on the geometry) to avoid too small elements. This is important because smaller elements require smaller time steps, and hence, a higher computational cost and more memory.

3) Finally, the creation of too distorted/elongated elements must be avoided, since this could result in a mesh of

poor quality, or even a mesh with elements with negative Jacobians, which would make the simulations unstable.

The main difficulties for our configuration are two-fold: the bathymetric profile near the coast, and the nature of the sediments that are unconsolidated.

The first difficulty is removed thanks to the above-considered physical insights related to SAWs. Indeed, we have previously shown (see Section III-B1) that since the modal acoustic propagation no longer takes place at a water depth of 30 m (corresponding to the range of ~5500 m along the line 3TY-PS13, see the red-dotted vertical line in Fig. 6), it is no longer necessary to mesh the area occupied by the water from this range to the coast. Nevertheless, in addition to the acoustic modes in water, SAWs also propagate close to the seabed. As a result, we have to be sure to preserve the integrity of the SAW propagation all along the line 3TY-PS13 as well. Therefore, considering the highest frequencies (i.e., 15–20 Hz) and the smallest wave velocity involved (i.e., $V_S$ sed = 200 m/s), together with the fact that, near the coast, the SAWs involved in the seismo-acoustic propagation are mainly of the Rayleigh(-Sezawa) type with their energy trapped in the sedimentary layer, we can finally choose to erect a “water wall” in the water domain mesh where the water depth is 5 m (corresponding to the range of ~12260 m along the line 3TY-PS13, see the green-dotted vertical line in Fig. 6).

The second difficulty we have to face in the numerical modeling concerns the presence of viscoelastic media with a high $V_P/V_S$ ratio (close to 10). This requires careful design of the mesh [66], since the regions with the smallest elements (corresponding to small values for $V_S$) drive the mesh structure.
<figure><img src="image_15.png" /><figcaption>Fig. 11. (a) Location of the fictitious vertical receiver array (colored dots) at the offset 3000 m from the source 3TY. (b) Simulated seismograms related to the horizontal component of the velocity. (c) Associated spectra.</figcaption></figure>

Additionally, since there is a long distance between the source and the receiver, and since small wavelengths are involved in the propagation, the number of time steps is high (1.8 10⁶). We have to rely on a time scheme with less numerical dispersion than the Newmark scheme that is classically used with the SEM. We choose to use a fourth-order low-dissipation low-dispersion Runge-Kutta scheme [67] to ensure the quality of our numerical results.

To use PMLs as absorbing layers at the right side of the computational domain, we need to add an artificial transition

layer at the end of the real domain to avoid unstabilities in the CPML behavior due to the strong velocity contrasts. As a consequence, the velocities inside this transition layer are changed continuously in the horizontal direction, so that at the right side of this transitional layer they are all identical (see Fig. 7). The velocities have then higher values, thus allowing the use of bigger elements.

The final mesh of the domain of interest is composed of 20 201 elements. The maximum and minimum edge lengths are 91.34 and 2.24 m, respectively. The time step for ensuring the
<figure><img src="image_16.png" /><figcaption>Fig. 12. (a) Location of the fictitious vertical receiver array (colored dots) at the offset 3500 m from the source 3TY. (b) Simulated seismograms related to the horizontal component of the velocity. (c) Associated spectra.</figcaption></figure>

stability of the time scheme is 7.5 10⁻⁵ s. We used 64 cores for the simulations. The computational cost was 3 h 10 min to simulate 135 s of wave propagation.

It is worth noting that in reality the mesh definition is composed of two steps. The first one that consists of accounting for the changes (with range and depth) in the

geometry of the environment, has been just described earlier. The second step that directly follows the first one, integrates a tomography model that describes the changes (with range and depth) in the properties of the environment, provided by the geological and acoustical analyses presented in Section II.
<figure><img src="image_17.png" /><figcaption>Fig. 13. (a) Location of the fictitious vertical receiver array (colored dots) at the offset 12440 m from the source 3TY (i.e., at the station PS13). (b) Simulated seismograms related to the horizontal component of the velocity. (c) Associated spectra.</figcaption></figure>

## IV. SEISMO-ACOUSTIC PROPAGATION INDUCED BY AN EXPLOSION ON THE SEA BOTTOM: NUMERICAL RESULTS AND DISCUSSION

We consider here the case where the detonation of the explosive charge occurs at the location 3TY ($h_w \approx 46$ m) on the seafloor. Since it is very difficult to put a source exactly at a fluid/solid interface in the SPECFEM code, the source is located more exactly at an arbitrary distance of 0.08 m below the water/sediments interface.

To identify the main limitations of the 3-D axisymmetric modeling, we first compare the numerical result of the seismo-acoustic propagation along the path 3TY-PS13 with the real data

recorded at the station PS13. We then simulate the signals at fictitious receivers distributed along the path between the source 3TY and the station PS13 to better understand the influence of the sedimentary layer on the seismo-acoustic wave propagation.

### A. *Comparison Between the Simulated and the Real Data*

The time signal representing the horizontal component of the velocity (that is colinear to the source-receiver direction) and obtained numerically, together with the associated spectrogram, is shown in Fig. 8(a). For comparison, Fig. 8(b) presents the real seismogram, and more specifically the East component of the
<figure><img src="image_18.png" /><figcaption>Fig. 14. Signals induced by the detonation of a 80-kg TNT-equivalent charge located at 3TY in the water column (explosion S8), and recorded at the station PS13. They represent the horizontal component of the velocity (colinear to the source-receiver direction), together with their spectrograms. (a) Simulated data. (b) Real data. Note that, as the real source function was not accounted for in the numerical simulations, the amplitude of the simulated data is only indicative and cannot be compared with the real data amplitude.</figcaption></figure>

velocity (also colinear to the source-receiver direction), induced by the detonation of a 400-kg TNT-equivalent charge (i.e., the explosion S6; see Table I in Part I of this article) and recorded at the station PS13 during the detonation campaign in December 2018. Fig. 8(b) also presents the associated spectrogram.

It has to be pointed out here that as the real source function is not accounted for in the numerical simulations, the amplitude of the simulated data is only indicative and cannot be compared with the real data amplitude. Therefore, the amplitudes of the time series have been normalized in Fig. 8.

We can see that despite the differences in the global wave-
forms and in the signal duration, the spectrograms associated
with the real and the simulated data are globally qualitatively
similar. In particular, both exhibit similar very LF components
(around 1.5 Hz) with significant amplitudes that are present for
the whole signal duration. In addition, the simulated spectrogram
exhibits a component spreading around 10–11 Hz that appears
more clearly in the real spectrogram at a sensibly lower fre-
quency (around 9 Hz). The LF components (around 1.5 Hz) are
likely associated with SAW propagation, while the “higher” fre-
quency components (around 9 Hz) are likely associated mainly
with body wave propagation (see [2]).

In contrast, the agreement between the simulated and the
real-time series is quite poor. This is particularly marked for the

waveforms from 15 to 50 s. Some of the corresponding events in
the real data likely correspond to SAWs that have evolved along
the propagation path according to the varying environment.
For instance, the events occurring at ~ 30 s may correspond
to a Stoneley-Scholte wave propagating over 9000 m at the
water/rock interface that subsequently evolved into a Rayleigh
wave propagating at the surface of the thick sedimentary basin to
the shore. Others likely correspond to scattering and multipaths
from the nearby islands and the coastline as well. These events
have been identified as coda in [2]. Unfortunately, all these
real features do not appear in the simulated results, probably
because of the assumptions made (namely, the 2.5-D model-
ing and the misestimation of the shear wave properties of the
sediments).

A focus on the first 15 s of the time series (see Fig. 9) shows that the simulated time series matches the real one slightly better, although the match is not perfect there either. Note that the amplitudes of the time signals have been normalized in Fig. 9. From the arrival times, we can identify both in the real and simulated time series the direct body waves in rock, followed by the Rayleigh wave at the rocky basement surface, and then the *P*-wave in sediments. However, the LF signal arriving at 15 s in the real data, and likely associated with SAWs, does not appear in the simulated data, probably because the value we chose for
<figure><img src="image_19.png" /><figcaption>Fig. 15. (a) Signals, and (b) their associated spectra, simulated for the case of the source on the seafloor (red curves) and for the case of the source located in the water column (blue curves), and obtained at the station PS13. Note that for the sake of better visualization the signal and spectrum amplitudes related to the source in the water column are enlarged by a factor of 40.</figcaption></figure>

the S-wave attenuation in the sediments greatly overestimates the actual one.

It is worth noting that, at this point, it is extremely difficult to associate, without doubt, each event with a specific wave. Indeed, the waves initially generated by the source are transformed along the propagation path. For instance, the arrival times between 15 and 20 s in the real data in Fig. 9 may correspond to two scenarios: either a Stoneley–Scholte wave that propagated at the water/rock interface over ~10 500 m, or a Rayleigh wave that propagated at the rocky basement surface over ~10 000 m, both having subsequently evolved into a Rayleigh wave that propagated at the surface of the sedimentary basin to the receiver. Some waves can even disappear, according to the varying geometry and the varying acoustical properties of the environment.

Only a full 3-D modeling of the wave propagation in the time domain that takes into account the complexity of the Rade d'Hyères could help. Indeed, the wavefronts at different propagation times could be monitored and the different types of waves together with their interaction during the propagation could be better identified. Unfortunately, as mentioned earlier, this modeling is out of reach at this time.

The differences observed between the simulated and the real data likely come from the fact that:

1) the numerical modeling, carried out in 2.5-D instead of 3-D, does not reproduce some important physical phenomena such as, for instance, 3-D scattering effects and multipaths from the nearby islands;

2) the uncertainties associated with the sedimentary layer thickness, although fairly small (±1–2 m), have a significant impact on the existence and the propagation

of SAWs, especially when the layer is thin, as discussed previously;

3) a misestimation of the S-wave attenuation in sediments can be dramatic, as shown by the absence of the real events arriving between 15 and 20 s in the simulated data in Fig. 9. Furthermore, neglecting the variations in the sediment and rock properties with depth does not allow accurate modeling of the SAW dispersion;

4) to a lesser extent, the uncertainties related to the explosion characteristics (location, emitted signal) may slightly modify the signals simulated at the station.

Nevertheless, despite all these uncertainties and the strong assumptions we made, the 2.5-D simulations can provide global results that are informative enough for the purpose of understanding the influence of the environment, including the sedimentary layer, in the generation of the large-amplitude LF events that may impact the infrastructures located on the shore.

## B. Analysis of the Wave Propagation Along the Source–Receiver Path

To better understand the seismo-acoustic wave propagation along the path between the source 3TY and the station PS13, we simulate the signals at 10 fictitious receivers distributed along the line path 3TY-PS13. The receivers are located in the sedimentary layer, more precisely at 1–2 m below the water/sediment interface at 100, 500, 1000, 3000, 3500, 5000, 7000, 9000, 11 000, and 12 440 m (station PS13) from the source, respectively (see Fig. 10). The seismograms, thus, obtained, together with their associated spectrum, are illustrated in Fig. 10. For the sake of a better visualization, the seismograms are shown at different amplitude scales.

To provide additional insight, we also simulate the seismo-acoustic wave propagation (time series and associated spectra) at fictitious receivers distributed along vertical arrays, above and below the sea bottom, and located at different offsets, namely 3000 (see Fig. 11), 3500 (see Fig. 12) and 12 440 m (i.e., at the location of the station PS13, Fig. 13), respectively.

It is worth noting that the seismograms in Figs. 10–13 are related to the horizontal component of the velocity (colinear to the source–receiver direction), while the seismograms related to the normal component are provided as supplements in Appendix B.

On Fig. 10(b), we note that for relatively small offsets near the source (i.e., for offsets ranging from 100 to 1000 m), the seismograms exhibit numerous signals with high amplitudes that, however, are attenuated during their propagation along the path to the station PS13. It is worth noting that amplitude amplifications can be observed close to the beginning of the sedimentary basins, namely for the offsets 3500 and 9000 m. This amplification, together with the emergence of LF components in the time series, is very local and more pronounced for the 7-m-thick basin (see Figs. 11 and 12). This likely results from the site effect induced by the particular shape of the short-scale sedimentary basin, as it is well known in engineering seismology (e.g., [3] and [4]).
<figure><img src="image_20.png" /><figcaption>Fig. 16. (a) Location of the fictitious horizontal receiver array (colored dots) along the line 3TY-PS13, (b) simulated seismograms, and (c) associated spectra obtained for the case of an explosion located in the water column at ~ 11 m below the sea surface. Note that for the sake of better visualization the signals and the spectra obtained at the different offsets are represented at different amplitude scales.</figcaption></figure>

For the near offsets, the frequency content of the time series globally spreads from 5 to 15 Hz with the main contribution in the 10–15-Hz range [see Fig. 10(c)]. However, for the large offsets where the sedimentary layer is thick, the frequency content of the seismograms shifts toward the LF domain (between 1 and 3 Hz).

Presenting the time series at different scales allows to emphasize the most significant events [see Fig. 10(d)]. If we pay attention to the first 10 s of the time series, we can identify some of the different types of waves discussed in Section III-B from the arrival times.

Quite logically, the first arrivals correspond to the direct *P* and *S* waves propagating in the rocky basement. Since rock is

an elastic solid medium, these body waves propagate quite well from the source to the receiver. However, their amplitude is small compared to that of the other waves.

These body waves are followed by a SAW of the Rayleigh type propagating at the surface of the rocky basement [see Fig. 10(d)]. Since its speed is slightly less than the *S*-wave speed (see Section III-B), the SAW can be clearly identified only for offsets larger than 3000 m. The water layer does not appear to affect the SAW propagation, meaning that the ratio of the water depth to the SAW wavelength is small enough for the water layer to be considered by the SAW as a load. As a result, the frequency content of the Rayleigh-type wave likely lies in the LF domain [see Figs. 10(d) and 13(b)]. Moreover,
<figure><img src="image_21.png" /><figcaption>Fig. 17. (a) H/V spectral ratio resulting from ambient noise data at the station PS15, and the best-fitting model emphasizing the natural frequency of the sedimentary basin in the Rade d'Hyères at the station site. (b) P- and S-wave velocity profiles as a function of depth, and the best-fitting models for the P- and S-wave velocities in the unconsolidated sediments.</figcaption></figure>

the thickness of the sedimentary layer $h_{sed}$ is much smaller than the SAW wavelength $\lambda_R$ over the 9000-m propagation path. The SAW then sees the sediments only as a coating. However, as the thickness of the sedimentary layer increases steadily from 7 to 30 m (i.e., for offsets ranging from 9000- to 12440-m station PS13), the impact of the sediments on the SAW propagation is clearly noticeable. Compared to the other events, the amplitude of the Rayleigh wave is greatly increased, while its frequency content is logically shifted toward the lower frequency domain to keep the $h_{sed}/\lambda_R$ ratio very small [see Fig. 10(d)]. A close examination of the signals and the associated spectra obtained at the vertical receiver array located at the offset 12440 m highlights the very LF components (around 1–2 Hz) of the Rayleigh wave [see Fig. 13(b)]. Its global amplitude, maximal just below the free surface, decreases with increasing depth in the bottom [see Fig. 13(b)].

A direct P-wave propagating in the sediments follows the Rayleigh wave in the time series [see Fig. 10(d)]. This body wave can be well monitored from the offset 1000 m to the offset 12440 m (station PS13), despite the attenuation in sediments.

This wave is likely fed by the direct P-wave that closely succeeds it in time and that propagates in water close to the sea bottom [see Fig. 10(d)]. Compared to the other waves, the amplitude of the P-wave in the water steadily increases when the water depth decreases. At the offset 7000 m, its amplitude is of the same magnitude as the amplitude of the P-wave in the sediments and as the amplitude of the Stoneley–Scholte wave that propagates at the water/rock interface and that closely succeeds it in time in the time series [see Fig. 10(d)].

It is worth specifying here that at first sight the successive P-wave arrivals and waveforms could resemble the sequence of waveforms of a Pekeris waveguide described, for instance, in [46], namely the sequence of the ground wave, the water wave and the Airy phase occurring in a two-layered liquid half-space. However, in theory in our configuration, below a water depth of 30 m (see Section III-B), and thus, for offsets greater than

6000 m, no acoustic mode is supposed to propagate in the water layer.

The Stoneley–Scholte wave propagating at the water/rock interface can be observed for offsets ranging from 500 m to 7000 m [see Fig. 10(d)]. This wave propagates as if the sedimentary layer did not exist. However, the presence of the 7-m-thick sedimentary basin within the range 3000–5000 m seems to indirectly affect the SAW (see Figs. 11 and 12).

Because their velocities have close values, the P-wave in the sediments, the P-wave in the water, and the Stoneley–Scholte wave at the water/rock interface interfere more and more as they approach the coast, more specifically for the offsets ranging from 9000 to 12440 m [see Fig. 10(d)]. Indeed, as the water depth decreases, the Stoneley–Scholte wave at the water/rock interface likely evolves gradually into a Rayleigh wave at the rocky basement surface with lower frequency components that do not see the thick sedimentary basin. As a result, the velocity of the evolved SAW wave steadily increases from ~ 1350 m/s (Stoneley–Scholte wave velocity) to ~ 2500 m/s (Rayleigh wave velocity) (see Section III-B). Moreover, the energy associated with the P-wave in the water is likely transmitted to the direct P-wave in the sedimentary layer when the water depth tends to zero. Consequently, the three waveforms associated with the P-wave in the water and in the sediments and with the SAW likely mix at far offsets close to the shore. This results in stretched waveforms such as those observed for instance between 4 and 6.5 s for the offset 9000 m and between 7 and 9 s for the offset 12440 m [see Fig. 10(d)].

If we focus on Fig. 13 (corresponding to the offset 12440 m), it is evident that the waveform that stretches from 7 to 9 s has higher frequency components (around 9–10 Hz) than the Rayleigh wave (around 1–2 Hz). Moreover, in contrast to the Rayleigh wave whose amplitude is maximal below the free surface, the maximal amplitude of the waveform is localized deeper within the sea bottom. The waveform likely corresponds to a Sezawa wave [43], [44].
<figure><img src="image_22.png" /><figcaption>Fig. 18. (a) Location of the fictitious horizontal receiver array (colored dots) along the line 3TY-PS13. (b) Simulated seismograms related to the normal component of the velocity. (c) Associated spectra. (d) Focus on the first 10 s of the seismograms represented in (b). Note that for the sake of better visualization the signals and the spectra obtained at the different offsets are represented in (c) and (d) at different amplitude scales, while the seismograms are represented in (b) at the same scale.</figcaption></figure>

The presence of the 30-m-deep sedimentary basin seems to
favor seismic events with LF contents and long durations. These
events are associated with SAWs that result from complex wave
evolvements during the propagation from the source to the shore.
If they are of high amplitude, this kind of events may be harmful
for the buildings on the shoreline. In our configuration, except
the slight amplification of the seismo-acoustic events observed

<figure><img src="image_23.png" /><figcaption>Fig. 19. (a) Location of the fictitious vertical receiver array (colored dots) at the offset 3000 m from the source 3TY. (b) Simulated seismograms related to the normal component of the velocity. (c) Associated spectra.</figcaption></figure>

locally at offsets close to 9000 m, no signal amplification has
been observed for larger offsets, since the thick sedimentary
basin is laterally well extended. The signal amplitudes then
decrease monotonically with the propagation distance, due to the
absorbing property of sediments. The LF signals, thus, become
weak enough to likely be harmless for the infrastructures on the
coast (see [2]).

It is worth noting that a SAW of the Stoneley–Scholte wave
type or the Rayleigh type should occur in the time series at
near offsets, or even at far offsets once waves evolve about the
thick sedimentary basin. On Fig. 10(d), for near offsets, we can
identify a SAW of the Stoneley–Scholte wave type from the
arrival times. But, unfortunately, the value we have chosen for
the S-wave attenuation is too high to allow the waves propagate
over far distances. The same conclusion holds for the S-wave in
the sediments.
<figure><img src="image_24.png" /></figure>

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

<figure><img src="image_26.png" /><figcaption>Fig. 20. (a) Location of the fictitious vertical receiver array (colored dots) at the offset 3500 m from the source 3TY. (b) Simulated seismograms related to the normal component of the velocity. (c) Associated spectra.</figcaption></figure>

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

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

<figure><img src="image_29.png" /><figcaption>Fig. 21. (a) Location of the fictitious vertical receiver array (colored dots) at the offset 12440 m from the source 3TY (i.e., at the station PS13). (b) Simulated seismograms related to the normal component of the velocity. (c) Associated spectra.</figcaption></figure>

## V. SEISMO-ACOUSTIC PROPAGATION INDUCED BY AN EXPLOSION IN THE WATER COLUMN: NUMERICAL RESULTS AND DISCUSSION

We now consider the case where the detonation of the explosive charge still occurs at the location 3TY but close to the sea surface (namely, about 11 m below the air/water interface).

It is worth specifying here that considering the expressions for the maximum bubble radius and for the vertical migration provided in [68], a charge with a weight of up to 80 kg TNT-equivalent is deep enough to generate at least the first bubble pulse. The *in situ* observations reported in [2] support these theoretical results.

The time signal representing the horizontal component of the velocity (that is colinear to the source-receiver direction) and obtained numerically, together with the associated spectrogram, is shown in Fig. 14(a). For comparison, Fig. 14(b) presents the real seismogram, and more specifically the East component of the velocity (also colinear to the source-receiver direction), induced by the detonation of a 80-kg TNT-equivalent charge

(explosion S8, [2, Table I]) and recorded at the station PS13 during the detonation campaign in December 2018. Fig. 14(b) also presents the associated spectrogram.

At first sight, we note significant differences in the global waveforms and in the signal duration, as well as in the spectrograms, associated with the real and the simulated data. In particular, the real data exhibit signals that have been seemingly produced initially by 3-D reverberation by the islands and the coast, and subsequently transformed into seismic events (e.g., signals arriving at ~ 90 s). Carried out in 2.5-D instead of 3-D, the numerical modeling is not able to reproduce such physical phenomena. In addition, it is worth recalling here the abovementioned uncertainties and the strong assumptions underlying the numerical modeling.

However, despite the significant differences, the main common feature that needs to be mentioned here is that the real and the simulated spectrograms do not exhibit the energetic LF component around 1.5 Hz observed for the case of the explosion located on the seabed (see Fig. 8 for comparison). This feature is
confirmed by comparing the signals and their associated spectra, obtained at the station PS13, and simulated for the case of the source on the seafloor and for the case of the source located in the water column (see Fig. 15). In particular, the seismic signal generated by the source located in the water column has a much lower amplitude and its spectral content is restricted to the highest frequencies (mainly above 9 Hz).

As shown in Fig. 16, the sedimentary basins have a much smaller impact on the seismo-acoustic wave propagation when the explosion occurs in the water column rather than on the seabed (see Fig. 10 for comparison). In particular, no significant signal amplification or frequency shift toward the LF domain are observed. This means that only a small amount of the hydroacoustic energy related to the source in the water column is transformed into seismic energy that subsequently propagates toward the coast with some loss due to the presence of sediments. Although this energy transformation is largely favoured by a soft sea bottom and less observed in the case of a hard seabed [27], and despite the fact that the lower the frequency, the more hydroacoustic energy is transformed into the seismic energy [28], this likely cannot lead to huge impact on inland infrastructures.

## VI. CONCLUSION

The main goal of this work, presented here in a two-companion paper, is to pave the way for assessing in a reliable manner the risk of building damage on the adjacent shore, induced by the detonation of large-charge historical ordnance (of between 80 and 680-kg TNT-equivalent weights) in variable shallow water environments with a water depth less than 50 m. While Part I of this article [2] focuses on the real acoustic and seismic data measured within the framework of a countermining campaign in the Rade d'Hyères (Mediterranean Sea, France), this article deals with the numerical modeling of the seismo-acoustic propagation induced by the charge explosion. Particular attention has been paid to the impact of the acoustical and geometrical properties of the sedimentary layer on wave propagation.

As modeling the 3-D propagation in a large marine environment with very low velocities implies a huge computational cost, we have first developed a strategy to provide clues for cost reduction. This strategy relies on physical insights into the different kind of waves that can propagate in a coastal water environment. Particular attention has been paid to the contribution of the surface/interface waves (SAWs), according to the shear-wave properties and the thickness of the sedimentary layer.

Considering information provided by the geological surveys and the hydroacoustic measurements (reported in [2]) as input data, a 3-D axisymmetric modeling of the seismo-acoustic wave propagation has been then conducted using a spectral-element method. The numerical simulations, obtained for one specific source-receiver path with a variable sedimentary facies including two basins of different shapes and thicknesses, have been compared with the real seismic data induced by a charge detonation located either on the seabed, or in the water column,

and recorded on the shore. Despite the uncertainties related to the environment properties and the assumptions underlying the numerical modeling, the simulations have provided global results that are informative enough for the purpose of understanding the influence of the environment on the seismo-acoustic wave propagation.

Analysis of the simulated data induced by the source on the seabed and obtained at different fictitious receivers along the propagation path has shed light on the strong interaction between the SAWs and the sedimentary facies, in particular close to the shore and at the beginning of the sedimentary basins. Indeed, as well known in seismology, short-scale and thick sedimentary basins favor a local wave-amplitude amplification and a frequency shift toward the LF domain. However, the seismo-acoustic waves are globally attenuated during their propagation from the source to the shore because of intrinsic attenuation of sediments and geometrical spreading.

For the case of a source in the water column, no strong signal amplification or frequency shift toward the LF domain has been observed despite the presence of the sedimentary basins. Only a small amount of the hydroacoustic energy related to the source in the water column is transformed into seismic energy that subsequently propagates with attenuation toward the coast.

The work presented in this article has to be viewed as a preliminary work for assessing the seismic risks induced by the large-charge detonation in shallow water environments, through the understanding of the seismo-acoustic wave propagation. More work definitely needs to be done for developing a reliable decision support tool.

The varying marine environment considered in this article induces a very complex seismic-acoustic wave propagation. Still associated with physical insights, additional numerical simulations of wave propagation conducted for various simpler cases of marine environment would be of great interest for better understanding how the different kinds of waves (i.e., the body waves, the acoustic modes, the SAWs) appear, interact each other, and evolve along the propagation path. This work is currently under investigation.

The acoustical and geometrical properties of the marine environment, including the water layer, the sedimentary layer, and the rocky basement, have to be assessed as accurately as possible, since they strongly influence the seismo-acoustic wave propagation. Careful attention must be paid to the thickness and to the weak shear rigidity of the sediments (including its variation with depth), as these parameters control the SAW propagation [69]. Besides more accurate direct measurements provided by specific tools, analysis of the Stoneley-Scholte wave propagation can be an effective alternative to refraction or reflection surveys for the estimation of the weak shear-wave velocity in the shallow subseafloor and its gradient with depth [25], [23], [70], [71].

In addition, the accurate knowledge of the response of the large-charge detonation, as well as the coupling between the charge detonation and the seabed, would be of great importance for modeling in a more reliable way the seismo-acoustic propagation induced by explosions.
# APPENDIX A
## ESTIMATION OF THE SHEAR-WAVE VELOCITY OF THE SEDIMENTS IN THE RADE D'HYÈRES USING AN H/V SPECTRAL RATIO ANALYSIS

Because we could not measure directly the S-wave velocity of the sediments in the Rade d'Hyères, we had to evaluate it.

We chose to rely on the H/V spectral ratio method (i.e., the horizontal-to-vertical-components spectral ratio), well-known in seismology to infer near-surface shear wave velocities (e.g., [38]). Indeed, in the case of large impedance contrasts, the frequency peaks of the H/V spectral ratio can be related to the resonant frequencies of the sedimentary layer. The fundamental frequency ($f_0$) is equal to $V_S / 4 h_{sed}$, where $h_{sed}$ is the layer thickness. Knowing $h_{sed}$ then allows to estimate $V_S$ under the assumption of an homogeneous sedimentary layer.

Here, we calculated the H/V curve using one hour of ambient noise data recorded at the PS15 station located in the Rade d'Hyères. The inversion was carried out using HV-Inv, a recent computer code for forward calculation and inversion of H/V spectral ratios of ambient noise based on the diffuse field assumption [72]. The inversion results for the station PS15 are shown in Fig. 17. Observing the first frequency peak $f_0$ around 1.5 Hz and considering that the uncertainty on the sedimentary layer thickness is about ± 1–2 m, the value for $V_S$ could be estimated within the range of 180–210 m/s.

The inversion of the H/V curves implies a strong tradeoff between the sedimentary layer thickness $h_{sed}$ and the S-wave velocity $V_S$ in sediments. However, considering one hour of ambient noise data recorded at the other stations on the sedimentary basin (e.g., stations PS12–PS14), instead of the station PS15, resulted in a similar range of values for $V_S$ (between 180 and 300 m/s).

It must be stressed that this analysis assumes that the shear-wave velocities in the sediments in the basin (offshore and onshore) are homogeneous, which in fact is a strong assumption. This assumption, however, is reasonable, given the nature and the total thickness of the sedimentary layer (less than 30–35 m deep), given the slight variations of the P-wave velocities in the sedimentary basin, and given the other assumptions and approximations made throughout this work.

# APPENDIX B
## SIMULATED SEISMOGRAMS RELATED TO THE NORMAL COMPONENT OF THE VELOCITY AND ASSOCIATED SPECTRA

For the sake of comparison, and as supplement materials, we provide here the simulated seismograms related to the normal component of the velocity (perpendicular to the source-receiver direction) and the associated spectra.

The signals are simulated at 10 fictitious receivers distributed along the path 3TY-PS13. The receivers are located in the sedimentary layer, more precisely at 1–2 m below the water/sediment interface at 100, 500, 1000, 3000, 3500, 5000, 7000, 9000, 11000, and 12440 m (station PS13) from the source, respectively (see Fig. 18).

Time series and their associated spectra are also simulated at fictitious receivers distributed along vertical arrays, above and

below the sea bottom, and located at different offsets, namely 3000 m (see Fig. 19), 3500 m (see Fig. 20), and 12440 m (i.e., at the location of the station PS13, Fig. 21), respectively.

## ACKNOWLEDGMENT

The authors would like to thank the Editor E. N. Ross Chapman and three anonymous reviewers for their constructive reviews and for their highly valuable suggestions that greatly improved this article. The authors also would like to thank the Prédacture de la Méditerranée, the Groupe des Plongeurs Démineurs (GPD) and the Force de Guerre des Mines for their management of the two countermining operations. Shom/GHOA, and in particular X. Mathias, P. Guyomard, E. Brenon, H. Gauduin, M. Gosselin, J.-P. Boivin, and F. Jourdin, are acknowledged for the hydrographic and geological surveys. D. Ambrois, X. Martin and J. Chèze (Géoazur), M. Pernoud, and Ph. Langlaude (Cerema), and M. Bonnin (LGP) are also acknowledged for the seismic surveys. This article is dedicated to our late colleague Dimitri Komatitsch (1970–2019).

## REFERENCES

[1] S.-H. Cheong, L. Wang, P. Lepper, and S. Robinson, Characterisation of Acoustic Fields Generated by UXO Removal – Phase 2. National Physical Laboratory Report AC19. Teddington, U.K.: National Physical Laboratory, 2020.

[2] N. Favretto-Cristini et al., "Assessment of risks induced by countermining unexploded large-charge historical ordnance in a shallow water environment—Part I: Real case study," *IEEE J. Ocean. Eng.*, early access, 2021, doi: 10.1109/JOE.2021.3111819.

[3] P. Y. Bard and M. Bouchon, "The two-dimensional resonance of sediment-filled valleys," *Bull. Seismol. Soc. Amer.*, vol. 75, no. 2, pp. 519–541, 1985.

[4] F. J. Sánchez-Sesma and F. Luzón, "Seismic response of three-dimensional alluvial valleys for incident P, S, and Rayleigh waves," *Bull. Seismol. Soc. Amer.*, vol. 85, no. 1, pp. 269–284, 1995.

[5] V. M. Cruz-Atienza et al., "Long duration of ground motion in the paradigmatic Valley of Mexico," *Sci. Rep.*, vol. 6, no. 1, 2016, Art. no. 38807. [Online]. Available: https://doi.org/10.1038/srep38807.

[6] F. B. Jensen and C. M. Ferla, "Numerical solutions of range-dependent benchmark problems in ocean acoustics," *J. Acoust. Soc. Amer.*, vol. 87, pp. 1499–1510, 1990.

[7] J. T. Goh, H. Schmidt, P. Gerstoft, and W. Song, "Benchmarks for validating range-dependent seismo-acoustic propagation codes," *IEEE J. Ocean. Eng.*, vol. 22, no. 2, pp. 226–236, Apr. 1997.

[8] D. Komatitsch and J.-P. Vilotte, "The spectral-element method: An efficient tool to simulate the seismic response of 2D and 3D geological structures," *Bull. Seismol. Soc. Amer.*, vol. 88, no. 2, pp. 368–392, 1998.

[9] D. Komatitsch and J. Tromp, "Introduction to the spectral-element method for three-dimensional seismic wave propagation," *Geophys. J. Int.*, vol. 139, no. 3, pp. 806–822, 1999.

[10] A. Fichtner, *Full Seismic Waveform Modelling and Inversion*. Berlin, Germany: Springer, 2010.

[11] D. Peter et al., "Forward and adjoint simulations of seismic wave propagation on fully unstructured hexahedral meshes," *Geophys. J. Int.*, vol. 186, pp. 721–739, 2011.

[12] P. Cristini and D. Komatitsch, "Some illustrative examples of the use of a spectral-element method in ocean acoustics," *J. Acoust. Soc. Amer.*, vol. 131, no. 3, pp. EL229–EL235, 2012.

[13] A. Bottero, P. Cristini, D. Komatitsch, and M. Asch, "An axisymmetric time-domain spectral-element method for full-wave simulations: Application to ocean acoustics," *J. Acoust. Soc. Amer.*, vol. 140, pp. 3520–3530, 2016.

[14] J. Tromp, D. Komatitsch, and Q. Liu, "Spectral-element and adjoint methods in seismology," *Commun. Comput. Phys.*, vol. 3, pp. 1–32, 2008.

[15] L. Carrington et al., "High-frequency simulations of global seismic wave propagation using SPECFEM3D GLOBE on 62 thousand processor cores," in *Proc. ACM/IEEE Conf. Supercomput.*, 2008, pp. 1–11.
[16] S. P. Oliveira and G. Seriani, "Effect of element distortion on the numerical dispersion of spectral-element methods," *Commun. Comput. Phys.*, vol. 9, no. 4, pp. 937-958, 2011.

[17] Z. Xie, R. Matzen, P. Cristini, D. Komatitsch, and R. Martin, "A perfectly matched layer for fluid-solid problems: Application to ocean-acoustics simulations with solid ocean bottoms," *J. Acoust. Soc. Amer.*, vol. 140, pp. 165-175, 2016.

[18] D. Komatitsch, J. Carcione, F. Cavallini, and N. Favretto-Cristini, "Elastic surface waves in crystals. Part 2: Cross-check of two full-wave modeling methods," *Ultrasonics*, vol. 51, pp. 878-889, 2011.

[19] D. Komatitsch, Q. Liu, J. Tromp, P. Süss, C. Stidham, and J. H. Shaw, "Simulations of ground motion in the Los Angeles basin based upon the spectral-element method," *Bull. Seismol. Soc. Amer.*, vol. 94, no. 1, pp. 187-206, 2004.

[20] S. Tsuboi, D. Komatitsch, C. Ji, and J. Tromp, "Computations of global seismic wave propagation in three-dimensional earth models," *Lecture Notes Comput. Sci.*, vol. 4759, pp. 434-443, 2008.

[21] B. Solymosi *et al.*, "How to adapt numerical simulation of wave propagation and ultrasonic laboratory experiments to be comparable—A case study for a complex topographic model," *Geophysics*, vol. 83, no. 4, pp. T195-T207, 2018.

[22] B. Solymosi, N. Favretto-Cristini, V. Monteiller, P. Cristini, B. Ursin, and D. Komatitsch, "Seismic surveying and imaging at the laboratory scale: A framework to cross-validate experiments and simulations for a salt-body environment," *Geophysics*, vol. 85, no. 3, pp. T123-T139, 2020.

[23] A. Bottero, P. Cristini, D. Komatitsch, and Q. Brissaud, "Broadband transmission losses and time dispersion maps from time-domain numerical simulations in ocean acoustics," *J. Acoust. Soc. Amer.*, vol. 144, no. 3, pp. EL222-EL228, 2018.

[24] F. B. Jensen, W. A. Kuperman, M. B. Porter, and H. Schmidt, *Computational Ocean Acoustics*, 2nd ed. Berlin, Germany: Springer, 2011, p. 794.

[25] G. Jamet, C. Guennou, L. Guillon, C. Mazoyer, and J.-Y. Royer, "T-wave generation and propagation: A comparison between data and spectral-element modeling," *J. Acoust. Soc. Amer.*, vol. 134, no. 4, pp. 3376-3385, 2013.

[26] J. Lecoulant, C. Guennou, L. Guillon, and J.-Y. Royer, "Three-dimensional modeling of earthquake generated acoustic waves in the ocean in simplified configurations," *J. Acoust. Soc. Amer.*, vol. 146, no. 3, pp. 2113-2123, 2019.

[27] A. Bottero, P. Cristini, and D. Komatitsch, "On the influence of slopes, source, seabed and water column properties on T-waves generation at shore," *Pure Appl. Geophys.*, vol. 177, no. 12, pp. 5695-5711, 2020.

[28] J. L. Stevens *et al.*, "Calculation of hydroacoustic propagation and conversion to seismic phases at T-stations," *Pure Appl. Geophys.*, vol. 178, pp. 2579-2609, 2020, doi: 10.1007/s00024-020-02556-3.

[29] O. Morio, T. Garlan, X. Mathias, P. Guyomard, Y. Le Faou, and H. Gauduin, "ANR POSA: Rapport descriptif du modèle 3D géologique multi-paramètres de la Rade d'Hyères et de ses abords (in French)," *Service of Hydrography and Oceanography of the Navy, Brest, France*, Final Rep. 35, 2017.

[30] J. Blanc, "Etude sédimentologie de la presqu'île de Giens et de ses abords (in French)," *Recueil des Travaux de la Station Mar. d'Endoume*, vol. 33, no. 20, pp. 35-52, 1960.

[31] A. Stépanian and I. Thinion, Evaluation Des Risques Naturels Littoraux—SCoT Provence Méditerranée—Evaluation Des Stocks Sédimentaires Littoraux—Prospection Géophysique (in French), Final Rep. BRGM/RP-64158-FR, 2014.

[32] H. B. Ali and L. D. Bibee, "The influence of a thin sediment layer on the dispersion of Scholte interface waves," *Acoust. Lett.*, vol. 17, no. 1, pp. 13-18, 1993.

[33] X. Demoulin, T. Garlan, L. Guillon, and P. Guyomard, "In situ acoustic measurements of water saturated beach sands," in *Proc. Inst. Acoust.*, 2015, vol. 37, pp. 40-47.

[34] J. H. Beebe and C. W. Holland, "The effect of unconsolidated sediment rigidity on low-frequency acoustic propagation," in *Proc. Ocean Seismo-Acoust.: Low-Freq. Underwater Acoust. (NATO Conf. Ser., Mar. Sci.)*, 1986, pp. 207-215.

[35] D. T. Smith, "Geotechnical characteristics of the sea bed related to seismo-acoustics," in *Proc.Ocean Seismo-Acoust.: Low-Freq. Underwater Acoust. (NATO Conf. Ser., Mar. Sci.)*, 1986, pp. 483-500.

[36] A. G. Soloway, P. H. Dahl, and R. I. Odom, "Modeling explosion generated Scholte waves in sandy sediments with power law dependent shear

wave speed," *J. Acoust. Soc. Amer.*, vol. 138, no. 4, pp. EL370-EL374, 2015.

[37] K. M. Lee *et al.*, "In situ measurements of sediment acoustic properties in Currituck Sound and comparison to models," *J. Acoust. Soc. Amer.*, vol. 140, no. 5, pp. 3593-3606, 2016.

[38] S. Foti, C. G. Lai, G. J. Rix, and C. Strobbia, *Surface Wave Methods For Near-Surface Site Characterization*. Boca Raton, FL, USA: CRC Press, 2015.

[39] S. T. McDaniel and J. H. Beebe, "Influence of semi-consolidated sediments on sound propagation in a coastal region," in *Proc. Bottom-Interacting Ocean Acoust. (NATO Conf. Ser., Mar. Sci.)*, 1980, pp. 493-506.

[40] R. D. Stoll, "Velocity dispersion in water-saturated granular sediment," *J. Acoust. Soc. Amer.*, vol. 111, no. 2, pp. 785-793, 2002.

[41] A. C. Kibblewhite, "Attenuation of sound in marine sediments: A review with emphasis on new low-frequency data," *J. Acoust. Soc. Amer.*, vol. 86, no. 2, pp. 716-738, 1989.

[42] D. Rauch, "On the role of bottom interface waves in ocean seismo-acoustics: A review," in *Proc. Ocean Seismo-Acoust.: Low-Freq. underwater Acoust. (NATO Conf. Ser., Mar. Sci.)*, 1986, pp. 623-641.

[43] O. A. Godin, "Surface-to-volume wave conversion in shallow water with a gently sloping bottom," *Acoust. Phys.*, vol. 53, no. 6, pp. 714-720, 2007.

[44] G. I. Dolgikh, S. S. Burdin, S. G. Dolgikh, V. V. Ovcharenko, V. A. Chupin, and S. V. Yakovenko, "Particulars of a transmitted acoustic signal at the shelf of decreasing depth," *J. Acoust. Soc. Amer.*, vol. 142, no. 4, pp. 1990-1996, 2017.

[45] B. Katsnelson, V. Petnikov, and J. F. Lynch, *Fundamentals of Shallow Water Acoustics (The Underwater)*. New York, NY, USA: Springer, 2012.

[46] C. L. Pekeris, "Theory of propagation of explosive sound in shallow water," *Memoirs—Geological Soc. Amer. 0072-1069*, vol. 27, pp. 1-117, 1948.

[47] D. Komatitsch and J. Tromp, "Spectral-element simulations of global seismic wave propagation—II. Three-dimensional models, oceans, rotation and self-gravitation," *Geophys. J. Int.*, vol. 150, pp. 303-318, 2002.

[48] N. Favretto-Anrès, "Theoretical study of the Stoneley-Scholte wave at the interface between an ideal fluid and a viscoelastic solid," *Acta Acustica, Acustica*, vol. 82, no. 6, pp. 829-838, 1996.

[49] N. Favretto-Anrès and G. Rabau, "Excitation of the Stoneley-Scholte wave at the boundary between an ideal fluid and a viscoelastic solid," *J. Sound Vib.*, vol. 203, no. 2, pp. 193-208, 1997.

[50] J. M. Carcione and H. B. Helle, "The physics and simulation of wave propagation at the ocean bottom," *Geophysics*, vol. 69, no. 3, pp. 825-839, 2004.

[51] W. L. Roever, T. F. Vining, and E. Strick, "Propagation of elastic wave motion from an impulsive source along a fluid/solid interface," *Philos. Trans. Roy. Soc. Lond. A*, vol. 251, pp. 455-523, 1959.

[52] E. Strick, "Propagation of elastic wave motion from an impulsive source along a fluid/solid interface: III. The pseudo-Rayleigh wave," *Philos. Trans. Roy. Soc. Lond. A*, vol. 251, pp. 488-523, 1959.

[53] J. R. Chamuel, "Laboratory studies on pulsed leaky Rayleigh wave components in a water layer over a solid bottom," in *Shear Waves in Marine Sediments*, J. M. Hovem, M. D. Richardson, and R. D. Stoll, eds. Dordrecht The Netherlands: Kluwer, 1991, pp. 59-66.

[54] D. Rauch, "Experimental and theoretical studies of seismic interface waves in coastal waters," in *Proc. Bottom-Interacting Ocean Acoust. (NATO Conf. Ser., Mar. Sci.)*, 1980, pp. 307-327.

[55] F. Schirmer, "Experimental determination of properties of the Scholte wave in the bottom of the North Sea," in *Proc. Bottom-Interacting Ocean Acoust. (NATO Conf. Ser., Mar. Sci.)*, 1980, pp. 285-298.

[56] G. Nolet and L. M. Dorman, "Waveform analysis of Scholte modes in ocean sediment layers," *Geophys. J. Int.*, vol. 125, pp. 385-396, 1996.

[57] N. Favretto-Anrès, "Use of the Stoneley-Scholte and Love waves for the acoustical characterization of marine sediments (in French)," Ph.D. dissertation, Université de la Méditerranée, Marseille, France, 1997.

[58] J. D. Achenbach, *Wave Propagation in Elastic Solids*. Amsterdam, The Netherlands: Elsevier, 1975.

[59] I. A. Viktorov, *Rayleigh and Lamb Waves: Physical Theory and Applications*. Berlin, Germany: Springer, 2013.

[60] G. W. Farnell and E. L. Adler, "Elastic wave propagation in thin layers," in *Physical Acoustics IX*, W.P. Mason and R. N. Thurston, eds. New York, NY, USA: Academic, 1972, pp. 35-127.

[61] K. Sezawa and K. Kanai, "The M₂ seismic waves," *Bull. Earthq. Res. Inst. (Tokyo)*, vol. 13, pp. 471-475, 1935.

[62] A. Kubotera, "Rayleigh and Sezawa waves generated by explosions," *J. Phys. Earth*, vol. 5, no. 1, pp. 33-41, 1957.
[63] K. S. Hunter and T. L. Geers, "Pressure and velocity fields produced by an underwater explosion," *J. Acoust. Soc. Amer.*, vol. 115, no. 4, pp. 1483-1496, 2004.

[64] D. Blair, "Seismic radiation from an explosive column," *Geophysics*, vol. 75, no. 1, pp. E55-E65, 2010.

[65] T. Blacker, "CUBIT mesh generation environment users manual," 1994. Accessed: Oct. 8, 2018. [Online]. Available: https://cubit.sandia.gov/public/13.2/help_manual/WebHelp/cubit_users_manual.html

[66] P. Moczo, J. Kristek, M. Galis, E. Chaljub, and V. Etienne, "3-D finite-difference, finite-element, discontinuous-Galerkin and spectral-element schemes analysed for their accuracy with respect to P-wave to S-wave speed ratio," *Geophys. J. Int.*, vol. 187, pp. 1645-1667, 2011.

[67] J. Berland, C. Bogey, and C. Bailly, "Low-dissipation and low-dispersion fourth-order Runge-Kutta algorithm," *Comput. Fluids*, vol. 35, pp. 1459-1463, 2006.

[68] W. D. Reid, "The Response of Surface Ships to Underwater Explosions," DSTO Aeronautical and Maritime Res. Lab., Fairbairn, ACT, Australia, Rep., 1996.

[69] O. A. Godin, "Shear waves and sound attenuation in underwater waveguides," *J. Acoust. Soc. Amer.*, vol. 149, no. 5, pp. 3586-3598, 2021.

[70] J. Ewing, J. A. Carter, G. H. Sutton, and N. Barstow, "Shallow water sediment properties derived from high-frequency shear and interface waves," *J. Geophys. Res.*, vol. 97, pp. 4739-4762, 1992.

[71] N. Favretto-Anrès and J.-P. Sessarego, "Identification of shear wave parameters of viscoelastic solids by laboratory measurements of Stoneley-Scholte waves," *Acta Acustica Acustica*, vol. 85, pp. 505-516, 1999.

[72] A. García-Jerez, J. Piña-Jerez, F. J. Sánchez-Sesma, F. Luzón, and M. Perton, "A computer code for forward calculation and inversion of the H/V spectral ratio under the diffuse field assumption," *Comput. Geosci.*, vol. 97, pp. 67-78, 2016.

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

Nathalie Favretto-Cristini received the M.S. degree in mechanics and the Ph.D. degree in underwater acoustics from the Université de la Méditerranée, Marseille, France, in 1993 and 1997, respectively.

From 1997 to 1999, then from 1999 to 2008, she was successively a Postdoctoral Fellow and a Researcher with the French National Centre for Scientific Research (CNRS), Laboratory of Modeling and Imaging in Geosciences, Pau, France. Since 2008, she has been a CNRS Researcher, Laboratory of Mechanics and Acoustics (LMA), Marseille, France. She was

the Deputy Director of this laboratory from 2012 to 2017. Her past and present research interests are focused on the interaction of waves with the interfaces in the sea bottom, and more specifically on seismo-acoustic wave propagation for underwater acoustics and seismic exploration applications.

Dr. Favretto-Cristini is a member of the American Society of Exploration Geophysicists (SEG), the European Association of Geoscientists and Engineers (EAGE), and the French Society of Acoustics (SFA). She was a Co-Chief Scientist for the POSA project.

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

Fang Wang received the Graduate degree from Ecole Centrale de Marseille, Marseille, France, in 2012 and the Ph.D. degree in geophysics from Mines ParisTech, Paris, France, in 2016.

Her thesis work concerned waveform inversion. From 2017 to 2019, she was a Postdoctoral Fellow with the Laboratory of Mechanics and Acoustics (LMA), Marseille, France, working on numerical simulation of seismo-acoustic wave propagation. Since 2020, she has been a Geophysicist for Compagnie Générale de Géophysique (CGG), Massy, France,

currently working on land seismic imaging projects.

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

Paul Cristini (Member, IEEE) received the Graduate degree from Ecole Nationale Supérieure de Physique de Marseille, Marseille, France, in 1986 and the M.S. and Ph.D. degrees in acoustics from the University of Aix-Marseille II, Marseille, France, in 1988 and 1991, respectively.

Since 1991, he has been a Researcher with the French National Centre for Scientific Research (CNRS). From 1991 to 2000, then from 2000 to 2008, he was with the Laboratory of Mechanics and Acoustics (LMA), Marseille, France, and with the Laboratory of Modeling and Imaging in Geosciences, Pau, France, respectively. Since 2008, he has been with LMA working on numerical modeling of wave propagation in underwater acoustics for various applications with emphasis on the use of High-Performance Computing.

Dr. Cristini is a member of the Acoustical Society of America (ASA) and of the French Society of Acoustics (SFA).

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

Thierry Garlan received the M.S. and Ph.D. degrees in geology from the University of Caen, Caen, France, in 1982 and 1985, respectively.

Since 1988, he has been an Expert Researcher in marine sedimentology and sediment dynamics modeling with the Service of Hydrography and Oceanography of the Navy (Shom), Brest, France. From 1988 to 1992, he was in charge of the development of nautical charts of atoll environments from satellite data. Since 2008, he has been the Head of the Department of Marine Geology, Shom. He has been the Manager

of many research projects concerned with sediment dynamics and funded by the French Ministry of Research or the French Ministry of Defence. He has set up and co-organized many international conferences devoted to marine and river dunes dynamics and to coastal dynamics.

Dr. Garlan is a member of two French Expert Committees focused on Solid Earth (Commission Spécialisée Terre Solide, CNRS/INSU and GT Terre Solide ALLENVI). He was the Lead Chief Scientist for the POSA project.

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

Olivier Morio received the M.S. degree in geophysics from the European Institute of Marine Sciences, Plouzané, France, in 2013 and the Ph.D. degree in marine geosciences from the University of South Brittany, Lorient, France, in 2017.

He was a Research Engineer with the Department of Marine Geology, Service of Hydrography and Oceanography of the Navy (Shom), Brest, France. He has been a part of many research projects regarding coastal dynamics and numerical modeling, especially on shallow water that were funded by the French

Ministry of Research or the French Ministry of Defence. Since 2018, he has been in charge of the onboard acoustic systems for Shom hydrographic fleet (Multibeam EchoSounders, Sub-bottom Profilers, and Acoustic Doppler Current Profilers). His current research works include seabed acoustic backscattering and classification.
<figure><img src="image_35.png" /></figure>

**E. Diego Mercerat** received the M.S. degree in applied geophysics from ITC, Delft, The Netherlands, in 2001, and the Ph.D. degree in civil engineering and geomechanics from the National Polytechnical Institute of Lorraine, Nancy, France, in 2007.

He was a Postdoctoral Fellow with the Institut de Physique du Globe de Paris (IPGP), Paris, France, and with the University of Nice-Sophia Antipolis, Nice, France, working on seismic modeling and tomography. Since 2011, he has been a Researcher in Seismology with the French Ministry of Ecology, Sustainable Development and Energy, with joint appointments at CEREMA Méditerranée and at the Géoazur Laboratory, University of Côte d'Azur, France. His research interests include seismic wave modeling and inversion, seismic response of civil engineering structures and geomorphological features, and seismic hazard assessment.

Dr. Mercerat is a member of the Scientific Committee of the RAP-RESIF consortium (French Accelerometric Network).

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

**Anne Deschamps** received the Ph.D. degree in geophysics from University Paris VI, Paris, France, in 1986.

From 1979 to 1992, she was a Researcher in Seismology with the Institut de Physique du Globe de Paris (IPGP), Paris, France. From 1992 to 2002, she was in charge of the regional seismic network in Provence-Alpes-Côte d'Azur region (south-eastern part of France) with the Géoazur Laboratory, Sophia Antipolis, France. She developed accelerometric and broadband seismological observation. She is currently an Emeritus Research Director with the French National Centre for Scientific Research (CNRS). Her research interests are mainly focused on the physics of earthquakes. They include the large earthquakes, seismicity, seismogenesis and tectonics in the Mediterranean area, and also different geometrical contexts like the large subduction of South America, the collision in Mongolia, and the rifting in Tanzania. She actively participated in the creation of the Corinth Rift Laboratory, namely the Single Transnational European Near Fault Observatory that was recently integrated in EPOS (a Research European Infrastructure). Her current research interests concern the seismological site effects and the impact of seismic waves on buildings.

Dr. Deschamps was a Co-Chief Scientist for the POSA Project.

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

**Vadim Monteiller** received the Ph.D. degree in geophysics from the University of Savoie, Chambéry, France, in 2005.

From 2008 to 2015, he was successively a Postdoctoral Fellow with the University Paul Sabatier, Toulouse, France, with the Laboratory of Mechanics and Acoustics (LMA), Marseille, France, and with the University of Nice-Sophia Antipolis, France. Since 2015, he has been a Research Engineer in Scientific Computing with LMA. His current research interests concern high-performance computing, computational seismology, and seismic imaging by full waveform inversion.

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

**Éric Beucler** received the M.S. degree in internal geophysics in 1998 and the Ph.D. degree in seismology in 2002, both from the University Paris VII—Institut de Physique du Globe de Paris (IPGP), Paris, France.

After a Postdoctoral Fellowship with the Department of Earth Sciences, University of Oxford, Oxford, U.K., he was an Assistant Professor in 2004 with the Laboratory of Planetology and Geodynamics, University of Nantes, Nantes, France. He became a Full Professor in 2018. His past and present research works concern global and regional seismology. He is deeply involved in the construction of the French permanent Broadband seismic array. He has developed nonlinear inverse procedures to infer the interior of planets and to detect and locate small magnitude seismic sources.

Dr. Beucler is the Co-I of the instrument SEIS (InSight NASA mission).