Antoine Guillemot
email: antoine.guillemot@univ-grenoble-alpes.fr

N Bontemps

E Larose

D Teodor

S Faller

L Baillet

S Garambois

E Thibert

O Gagliardini

Christian Vincent

Investigating Subglacial Water-filled Cavities by Spectral Analysis of Ambient Seismic Noise: Results on the Polythermal Tête-Rousse Glacier (Mont Blanc, France)

published or not. The documents may come

Introduction

Subglacial cavities may trap a considerable liquid water content. Within glaciers meltwater reservoirs can be created by a number of factors, such as topographical traps, heterogeneous manuscript submitted to Geophysical Research Letters thermal flux regimes, variability of seasonal meteorological conditions and, consequently, hydraulic pressure gradients [START_REF] Björnsson | Hydrological characteristics of the drainage system beneath a surging glacier[END_REF][START_REF] Vincent | Origin of the outburst flood from Glacier de Tête Rousse in 1892 (Mont Blanc area, France)[END_REF]. An unexpected and instantaneous release of water from intraglacial or subglacial cavities [START_REF] Haeberli | Glacier Changes Following the Little Ice Age -A Survey of the International Data Basis and Its Perspectives[END_REF][START_REF] Vincent | Origin of the outburst flood from Glacier de Tête Rousse in 1892 (Mont Blanc area, France)[END_REF], known as a glacial outburst flood, often causes avalanches of ice and debris flows [START_REF] Breien | Erosion and morphology of a debris flow caused by a glacial lake outburst flood[END_REF][START_REF] Cenderelli | Flow hydraulics and geomorphic effects of glacial-lake outburst floods in the Mount Everest region, Nepal[END_REF][START_REF] Kattelmann | Glacial Lake Outburst Floods in the Nepal Himalaya: A Manageable Hazard?[END_REF][START_REF] Kershaw | Geomorphic and sedimentological signature of a two-phase outburst flood from moraine-dammed Queen Bess Lake, British Columbia, Canada[END_REF]Vilimek et al., 2005), leading to devastating damages in densely populated mountain areas (e.g., [START_REF] Haeberli | Frequency and Characteristics of Glacier Floods in the Swiss Alps[END_REF]. Excavation and drainage processes to evacuate subglacial water-filled cavities remain under-studied [START_REF] Fountain | Water flow through temperate glaciers[END_REF][START_REF] Boon | The role of hydrologically-driven ice fracture in drainage system evolution on an Arctic glacier[END_REF][START_REF] Roberts | JÖkulhlaups: A reassessment of floodwater flow through glaciers[END_REF][START_REF] Vincent | Origin of the outburst flood from Glacier de Tête Rousse in 1892 (Mont Blanc area, France)[END_REF].

Several studies have attempted to identify intraglacial cavities early during their formation [START_REF] Haeberli | Prevention of outburst floods from periglacial lakes at Grubengletscher, Valais, Swiss Alps[END_REF][START_REF] Haeberli | A surge-type movement at Ghiacciaio del Belvedere and a developing slope instability in the east face of Monte Rosa, Macugnaga, Italian Alps[END_REF][START_REF] Vincent | Origin of the outburst flood from Glacier de Tête Rousse in 1892 (Mont Blanc area, France)[END_REF][START_REF] Werder | Hazard assessment investigations in connection with the formation of a lake on the tongue of Unterer Grindelwaldgletscher, Bernese Alps, Switzerland[END_REF][START_REF] Gagliardini | Estimating the risk of glacier cavity collapse during artificial drainage: The case of Tête Rousse Glacier[END_REF] with the aim of preventing and mitigating the related hazards. Both active and passive geophysical methods have been deployed to investigate the glacier-bedrock interface and to characterize the physics of glaciers (e.g., [START_REF] Petrenko | Physics of ice[END_REF]. Some commonly used active methods are Ground Penetrating Radar (GPR), seismic refraction, and borehole measurements. In addition, passive seismic techniques can provide a reliable estimation of the seismic properties linked to glacier morphology by measuring ambient seismic noise. By using efficient survey logistics, these passive seismic measurements come at a contained cost [START_REF] Sergeant | On the Green's function emergence from interferometry of seismic wave fields generated in high-melt glaciers: implications for passive imaging and monitoring[END_REF]. Ambient noise studies have also been successfully used to reveal water flow activity (e.g., [START_REF] Björnsson | Hydrological characteristics of the drainage system beneath a surging glacier[END_REF][START_REF] Harper | Vertical extension of the subglacial drainage system into basal crevasses[END_REF]Helmstetter et al., 2015 a, b;[START_REF] Gimbert | A Multi-Physics Experiment with a Temporary Dense Seismic Array on the Argentière Glacier, French Alps: The RESOLVE Project[END_REF]Nanni et al., 2021[START_REF] Nanni | Dynamic Imaging of Glacier Structures at High-Resolution Using Source Localization With a Dense Seismic Array[END_REF].

The horizontal-to-vertical-spectral ratio (HVSR) technique, first applied to ambient noise recordings by [START_REF] Nogoshi | On the amplitude characteristics of ambient noise (part 2)[END_REF] and [START_REF] Nakamura | A METHOD FOR DYNAMIC CHARACTERISTICS ESTIMATION OF SUBSURFACE USING MICROTREMOR ON THE GROUND SURFACE[END_REF], is based on the analysis of H/V frequency content, obtained by dividing the vectorial sum of the horizontal components by the vertical components of the Fourier amplitude for the recorded wavefield. When impedance contrasts exist at interfaces (e.g., a soft layer over a hard bedrock), incident body waves (P-waves, horizontal and vertical S-waves) generate surface (Rayleigh and Love) waves that remain trapped within the shallow layer. When the wavelengths of the incident waves are comparable to the layer's thickness, resonance effects occur. Consequently, the resulting surface waves may have a much larger amplitude than the incident waves (e.g., [START_REF] Picotti | Estimation of glacier thicknesses and basal properties using the horizontal-to-vertical component spectral ratio (HVSR) technique from passive seismic data[END_REF]. In most cases where a 1D soft layer is underlying a stiffer one, the shape of HVSR curve is primarily determined by the fundamental-mode Rayleigh wave [START_REF] Konno | Ground-motion characteristics estimated from spectral ratio between horizontal and vertical components of microtremor[END_REF]. Under reasonable assumptions, the manuscript submitted to Geophysical Research Letters first HVSR natural peak gives a reliable estimation of the S-wave resonance frequency of a horizontally layered subsurface [START_REF] Lachetl | Numerical and Theoretical Investigations on the Possibilities and Limitations of Nakamura's Technique[END_REF]. Other spectral peaks correspond to a superposition of higher Rayleigh-wave, Love-wave, and S-wave resonance modes [START_REF] Bonnefoy-Claudet | Effects of Love Waves on Microtremor H/V Ratio[END_REF].

The HVSR technique is commonly used to retrieve the site amplification function and its resonance frequency in shallow geological layers (e.g. [START_REF] Lachetl | Numerical and Theoretical Investigations on the Possibilities and Limitations of Nakamura's Technique[END_REF]Ibs-von Seht & Wohlenberg, 1999;[START_REF] Fäh | A theoretical investigation of average H/V ratios[END_REF][START_REF] Fäh | Analysis of three-component ambient vibration array measurements[END_REF][START_REF] Bonnefoy-Claudet | The nature of noise wavefield and its applications for site effects studies[END_REF][START_REF] Tuan | On the relationship of peaks and troughs of the ellipticity (H/V) of Rayleigh waves and the transmission response of single layer over half-space models[END_REF][START_REF] Panzera | Lusi hydrothermal structure inferred through ambient vibration measurements[END_REF][START_REF] Panzera | Reconstructing a 3D model from geophysical data for local amplification modelling: The study case of the upper Rhone valley, Switzerland[END_REF] and Alpine valleys [START_REF] Panzera | Reconstructing a 3D model from geophysical data for local amplification modelling: The study case of the upper Rhone valley, Switzerland[END_REF]). It is also deployed to map soft sediment layers and estimate glacier thickness (e.g., [START_REF] Barnaba | The buried shape of an alpine valley from gravity surveys, seismic and ambient noise analysis[END_REF][START_REF] Picotti | Estimation of glacier thicknesses and basal properties using the horizontal-to-vertical component spectral ratio (HVSR) technique from passive seismic data[END_REF][START_REF] Preiswerk | Effects of geometry on the seismic wavefield of Alpine glaciers[END_REF][START_REF] Molnar | A review of the microtremor horizontal-to-vertical spectral ratio (MHVSR) method[END_REF][START_REF] Seivane | On the use of the microtremor HVSR for tracking velocity changes: a case study in Campo de Dalías basin (SE Spain)[END_REF]. [START_REF] Picotti | Estimation of glacier thicknesses and basal properties using the horizontal-to-vertical component spectral ratio (HVSR) technique from passive seismic data[END_REF] applied the HVSR method to five glaciers in the Alpine chain and the west Antarctic ice stream. Their results show how well the resonance frequency in the HVSR spectra correlates with the local bedrock depth over a wide range of ice thicknessesfrom a few tens of meters to more than 800 meters. The authors validated their results using other methods, such as GPR, geoelectrical, and active seismic techniques. The limitations of the HVSR method and its reliability for glacier characterization have been addressed in various studies, as reviewed by [START_REF] Picotti | Estimation of glacier thicknesses and basal properties using the horizontal-to-vertical component spectral ratio (HVSR) technique from passive seismic data[END_REF] and [START_REF] Molnar | A review of the microtremor horizontal-to-vertical spectral ratio (MHVSR) method[END_REF]. Both groups identify coupling of the sensors with the glacier and the basal impedance contrast as the main controlling factors. [START_REF] Preiswerk | Effects of geometry on the seismic wavefield of Alpine glaciers[END_REF] investigated how the glacier's geometry influences the seismic wavefield and reported an azimuth dependence of the HVSR curves' shape. They relate such directivity to higher-order resonance effects caused by the glacier's morphology (like bedrock irregularities). On glaciers, [START_REF] Carcione | Effect of soil and bedrock anelasticity on the S-wave amplification function[END_REF] demonstrated that attenuation and bedrock deformability affect both frequency and amplitude of HVSR peaks. Moreover, [START_REF] Saenger | A hydrocarbon microtremor survey over a gas field: Identification of seismic attributes | SEG Technical Program Expanded Abstracts[END_REF][START_REF] Saenger | A passive seismic survey over a gas field: Analysis of low-frequency anomalies[END_REF] used the HVSR method to extract hydrocarbon maps, focusing on HVSR troughs (or peaks of the inverse -Vertical-to-Horizontalspectral ratio, noted VHSR) to detect impedance anomalies. Specifically, they reported a peak in the VHSR (typically, values above 1: i.e., vertical energy higher than horizontal energy) at a specific frequency detected by seismic sensors placed above hydrocarbon reservoirs. Similarly, [START_REF] Antunes | Insights into the dynamics of the Nirano Mud Volcano through seismic characterization of drumbeat signals and V/H analysis[END_REF] noticed that VHSR peaks provides information about seismic energy anomalies generated by fluids in reservoirs since the wavefield is mainly polarized in the vertical direction. They identified a link between the location of the VHSR peak amplitude and the presence of mud and gas inside a volcanic conduit; they also noticed a spatial correlation between manuscript submitted to Geophysical Research Letters high VHSR anomalies and low resistivity anomalies in geoelectric profiles, corresponding to dome-like reservoirs and associated with high-salinity fluids.

Here, we analyze spectral ratios from ambient seismic noise measured for the Tête-Rousse Glacier (Mont Blanc massif) with the aim of locating potential cavities or areas containing melting water within the glacier. We first present the study area and the seismic data experiment in the field, we next detail the results of spectral analysis obtained from our dense, temporary seismic array. We interpret horizontal polarization of seismic noise around 9 Hz detected by HVSR peaks as the fundamental mode for the whole glacier. Vertical polarization anomalies at higher frequencies detected by HVSR troughs (i.e, more intuitive, VHSR peaks) are considered to reflect reduced basal impedance due to the presence of a cavity or a local anomalous basal response in some areas of the glacier. Results obtained in the field are supported by a proof-of-concept simulation based on the 3D finite-element method.

Experimental site and setup

The Tête-Rousse Glacier is located in the Mont Blanc massif in the French Alps (Figure 1a), along one of the classic ascension routes (Figure 1b), making it a very popular site among mountaineers. Fed by avalanches in its upper part, this glacier is approximately 450 m long and 150 m wide, with a maximum depth of 80 m, and its surface topography extends from 3100 to 3300 m elevation. There is at least one subglacial water-filled cavity under the central part of the glacier [START_REF] Vincent | Detection of a subglacial lake in Glacier de Tête Rousse (Mont Blanc area, France)[END_REF], known as the main cavity of this glacier.

In 1892 (July 11), an outburst flood originated from the Tête-Rousse Glacier following collapse of an intraglacial cavity. Around 200 000 m 3 of water and ice were released, dragging about 800 000 m 3 of sediments (boulders, soil, and mud) into the village of Saint Gervais -Le Fayet. The catastrophic event caused 175 fatalities and damaged multiple large infrastructures [START_REF] Vincent | Origin of the outburst flood from Glacier de Tête Rousse in 1892 (Mont Blanc area, France)[END_REF]. The mechanism controlling formation and growth of the subglacial cavity has been clearly identified and is probably related to the glacier's thermal regime. According to [START_REF] Vincent | Origin of the outburst flood from Glacier de Tête Rousse in 1892 (Mont Blanc area, France)[END_REF], the cold tongue acts as a thermal barrier to water flow. At the time of the outburst flood, factors combined to cause the collapse of the terminal tonguedue to a change in thermal regime (from cold to temperate) or a positive anomaly of the mass balance.

Between 2009 and 2010, data were acquired on the Tête-Rousse Glacier from dense GPR and Surface Nuclear Magnetic Resonance (SNMR) acquisitions. Their analysis revealed a new manuscript submitted to Geophysical Research Letters subglacial water-filled reservoir measuring 55 000 m 3 (+/-10 000 m 3 ) located at a depth of 40 m inside the glacier [START_REF] Vincent | Detection of a subglacial lake in Glacier de Tête Rousse (Mont Blanc area, France)[END_REF][START_REF] Legchenko | Monitoring water accumulation in a glacier using magnetic resonance imaging[END_REF]. To prevent another potential outburst flood, the reservoir has been drained artificially (downhole pumps) several timesin 2010, 2011, and 2012to prevent the cavity filling up again. In 2013, a large crevasse almost reaching the surface was detected, corresponding to the former cavity roof that collapsed during a pump in 2012. The other few crevasses present on the glacier seem to be the main channel via which water runs into the cavity [START_REF] Garambois | Ground-penetrating radar and surface nuclear magnetic resonance monitoring of an englacial water-filled cavity in the polythermal glacier of Tête Rousse[END_REF].

From SNMR surveys conducted in 2015 and 2018, no increase in the water content in the main cavity was detected. However, liquid water was detected at several points spread over the upper part of the glacier. These signals could be linked to a water-filled crevasse located upstream from the cavity together with a diffusive zone, revealed using GPR by [START_REF] Garambois | Ground-penetrating radar and surface nuclear magnetic resonance monitoring of an englacial water-filled cavity in the polythermal glacier of Tête Rousse[END_REF] (see Figure S6). These observations suggest a zone where basal partial melting occurs (Gilbert et al., 2012). The glacier is thus in an unsteady state that requires to monitor both its water content and the cavity's geometry. The seismic acquisition presented in this paper was performed in 2022, from May 4 to 22. We deployed 60 SmartSolo® three-component (3C) 5Hz seismic nodes, in a grid pattern with approximately 30 m by 25 m spacing (Figure 1c). In addition, three 3C Trilium Compact PostHole ® 20s seismometers were installed (named C1 to C3, red dots in Figure 1c), forming a transversal manuscript submitted to Geophysical Research Letters profile, hence crossing the area above the expected main subglacial cavity (see Supporting Information S1 and Figure S1 for further experimental details).

Data analysis and experimental results

We first performed a statistical analysis in the frequency domain for all the signals recorded by the seismic stations over the twenty days of experiment. Results show that the horizontal spectral components have a higher amplitude than the vertical component between 5 Hz and 30 Hz. This is the case over almost the entire seismic array, as expected for glaciers [START_REF] Preiswerk | Effects of geometry on the seismic wavefield of Alpine glaciers[END_REF]. HVSR peaks (see Fig. 2a) are expected to reveal the glacier's vibration modes, as reported previously [START_REF] Picotti | Estimation of glacier thicknesses and basal properties using the horizontal-to-vertical component spectral ratio (HVSR) technique from passive seismic data[END_REF][START_REF] Preiswerk | Effects of geometry on the seismic wavefield of Alpine glaciers[END_REF]. Additional examples of power density spectrum (PSDs) plots for several sensors can be found in Supporting Information S2 (Fig. S2 andS3). Interestingly, the amplitude of the vertical component is higher than the amplitude of the horizontal component for some seismic sensors at specific frequencies, including sensor C2 located in the central part of the glacier above the expected main cavity (Fig. 2b). Regions where the vertical component dominated (VHSR peaks above 1) are indicative of unusual impedance contrasts at depth, suggesting for instance the presence of cavities or ice/water-related contrasts.

The first HVSR peak recorded by the C2 sensor and most of the nodes is at about 9 Hz (Figure 2a), whereas the first VHSR peak is at about 15 Hz (Figure 2b). We then analyzed how the seismic noise is polarized at these frequencies (see Supporting Information S3 for details). We derived the polarization attributes from singular value decompositions by applying Flinn's method, using the seismic post-processing package Obspy. Our first-pass analysis, between 8 Hz and 10 Hz, encompasses the first mode observed in the HVSR peak (see Figure S4). This analysis reveals high incidence values (often between 75° and 90°), thus horizontal motion of the particle located in C2 in this frequency band. The corresponding azimuth values mostly ranges between 70° and 120° (E-W direction), indicating a longitudinally-oriented horizontal particle motion, in agreement with the directivity of the HVSR peak near 9 Hz (Figure 2a) and the geometry of the glacier. HVSR peak frequency for all other sensors is close to 9 Hz (see Figure S5), with small variability mainly controlled by bedrock depth [START_REF] Picotti | Estimation of glacier thicknesses and basal properties using the horizontal-to-vertical component spectral ratio (HVSR) technique from passive seismic data[END_REF]. The azimuth's temporal variability is low, and is probably related to the natural variability of the wavefield propagating in the glacier.

manuscript submitted to Geophysical Research Letters

We then looked at the polarization attributes of the data encompassing the first VHSR peak (or HVSR trough) observed -between 14 Hz and 16 Hz (Figure 2a). The seismic signal exhibits a mostly vertical motion of the particle located in C2 (highly variable but with a low incidence value, often between 10° and 50°). The corresponding azimuth values covered all directions, with consistently vertical polarization. The polarization analysis thus confirms our first observations from spectral analysis: the first HVSR peak around 9 Hz corresponds to a horizontal vibration mode along the glacier's main longitudinal axis, whereas the peak around 15 Hz revealed by VHSR peak is mostly in vertical mode. Polarization was stable throughout the duration of recording, adding support to our interpretation that these spectral peaks reflect natural vibration modes of the Tête-Rousse Glacier. These conclusions were further corroborated by the 3D mechanical modeling developed below.

To better locate areas of vertical resonance all over the glacier, we computed the VHSR for each station over the whole time interval, from May 5 to May 21, using a modification of the hvsrpy Python tool (Vantassel et al., 2022). We then picked the maximum amplitude averaged over the azimuth within the 2 Hz to 30 Hz frequency range, so as to focus only on first natural modes. Finally, using the ordinary-kriging method (based on a Matern model and a Matheron estimator), we interpolated the VHSR peak amplitude over the whole area covered by the seismic network (Figure 3). This analysis reveals a VHSR below 0.8 for most of the lower part of the glacier, whereas above the fourth line of nodes (see Figure 1c), the VHSR is slightly over 0.8. Moreover, two areas of the glacier were associated with VHSR peaks significantly exceeding 1 (around sensors C2 and 4.1). These areas are consistent with the location of the main subglacial cavity identified in 2010.

In the following section, we describe how we built a simple mechanical 3D model of the glacier to reproduce these resonance modes and their polarization (vertical and horizontal). Figure 3: VHSR peak amplitude for the Tête-Rousse Glacier, computed with time-and azimuth-averaging over the whole period of experiment and interpolated by kriging. The cavity outline is deduced from [START_REF] Gagliardini | Estimating the risk of glacier cavity collapse during artificial drainage: The case of Tête Rousse Glacier[END_REF].

Numerical validation of results

We performed direct numerical modeling based on finite-element analysis (FEA) to visualize vibrating modes in the Tête-Rousse Glacier and to reconstruct the resonance frequencies along with their modal polarization. We used COMSOL Multiphysics software (Comsol, 2023) for two purposes: 1) to perform modal analysis (MA) of the glacier, solving eigenfrequencies of the 3D object to reproduce the fundamental mode of the whole glacier, and 2) to determine the manuscript submitted to Geophysical Research Letters frequency response function (FRF) of a 2D axisymmetric model of a local area of the glacier, retrieving the vertical displacements related to the presence of a cavity.

3D modeling of the horizontal glacier mode

MA computes numerical eigenfrequencies of a solid structure and requires wellconstrained structural knowledge of the vibrating object. By combining all GPR and borehole measurements, [START_REF] Garambois | Ground-penetrating radar and surface nuclear magnetic resonance monitoring of an englacial water-filled cavity in the polythermal glacier of Tête Rousse[END_REF] produced a well-resolved estimation of the ice-bedrock interface's topography (see Figure S6). For our model, we simplified the 3D geometry of the glacier to an ellipsoid-shaped ice block with dimensions matching the field observations, lying on a rigid bedrock. Our modeled solid is composed of 94 000 quadratic tetrahedral finite elements, with a maximum dimension of 15 m. We applied consistent boundary conditions to mimic the configuration of the site: embedding the sides of the stable bedrock far from the glacier body, and leaving the upper surface free. The mechanical parameters of both glacier and bedrock were estimated from the literature for ice and fractured moraine [START_REF] Sergeant | On the Green's function emergence from interferometry of seismic wave fields generated in high-melt glaciers: implications for passive imaging and monitoring[END_REF]Table S1). For the sake of simplicity, the properties of both glacier and bedrock were assumed to be homogeneous.

We found the first eigenfrequency, corresponding to the fundamental mode, to be at 9.42 Hz (Figure 4a), close to the experimental HVSR peak measured at most nodesat around 9 Hz (Figure 2a). This first mode is a flexural mode polarized along the longest axis of the glacier (the longitudinal direction). The polarization of this first numerical mode is consistent with the observed directivity of the 9 Hz HVSR peak measured by most of the seismic sensors around the 100° azimuth. We noted that observed HVSR peak frequencies show a spatial variability over the glacier (Figure S5 in Supplementary Materials), mainly depending on the local bedrock interface depth and crevasses. We actually measured local "pseudo-modes" of compartments of the glacier, rather than a global resonance of the whole glacier, due to the attenuation on the glacier that could not be addressed by the model. Nevertheless, this 3D model suggests a longitudinal resonance around 9 Hz, similar to local horizontal modes detected by HVSR peaks.

Overall, this MA validates to the first order the mechanical parameters set of the ice and the bedrock and our interpretation of horizontal modes detected by HVSR. Then, a 2D modeling is needed to accurately interpretate local vertical resonance modes detected by VHSR peaks. manuscript submitted to Geophysical Research Letters

2D modeling of the cavity mode

We then focused on the presence of a subglacial cavity to assess its effects on spectral analysis. For simplicity, we used a 2D geometry with vertical axial symmetry, focusing on the central part of the glacier. We assumed the simulated subglacial water-filled cavity to be elliptical, with realistic dimensions. We only constrained the bedrock depth at 75 m and the location of the top of the cavity at around 35 m depth, in line with GPR measurements obtained in 2022 (Figure 4b). Our 2D solid model was composed of 12 000 quadratic tetrahedral finite elements, with a maximum size of 3 m side. A mesh refinement provided greater than 1% precision for all the mechanical fields computed below 30 Hz (minimum of 20 elements per wavelength).

Acoustoelastic coupling between liquid water and ice was implemented by imposing continuous normal acceleration at the cavity interface. We simulated ambient seismic sources by applying a linear harmonic loading between 5 Hz and 20 Hz along a horizontal line crossing the geometry.

After computing the FRF, we retrieved the corresponding frequency-dependent response displacement in the vertical direction at several points on the surface (closely reproducing the positions at which seismic recordings were made by a vertical seismometer).

We compared the vertical displacements (along the Z-axis) obtained with and without the subglacial cavity in Figure 4c. For a simulated sensor placed right above it, vertical vibrations has a clearly larger amplitude when the cavity is present. In these conditions, the peak amplitude was two-fold that measured without the cavity, and has a slightly different frequency.

Discussion

The aim of this study was to analyze ambient noise data from the Tête-Rousse Glacier to verify whether it could be used to detect and monitor subglacial cavities. Findings were confirmed by modeling. The results presented in Figure 3 identify two hotspots with high (above 1.1) VHSR amplitude peaks in the central part of the glacier. The first hotspot (centered around the C2 sensor) corresponds to the center of the main cavity first detected in 2010, and confirmed in 2018 using GPR and by borehole drilling (Figure S6). The second one was measured only by node 4.1, which corresponds to the area above the 2010 cavity, but with a more limited lateral extension. This sensor was also the closest to a transversal crevasse identified on the glacier's right bank, known as the former roof of the cavity that collapsed in 2012 [START_REF] Garambois | Ground-penetrating radar and surface nuclear magnetic resonance monitoring of an englacial water-filled cavity in the polythermal glacier of Tête Rousse[END_REF]. Note that this crevasse cannot explain the VHSR amplitude peak, since we expect an exacerbated horizontal manuscript submitted to Geophysical Research Letters energy increase at the border of the subvertical crack, due to shear vibration modes. However, this crevasse was previously shown to feed the cavity with surface water [START_REF] Garambois | Ground-penetrating radar and surface nuclear magnetic resonance monitoring of an englacial water-filled cavity in the polythermal glacier of Tête Rousse[END_REF]. These recent observations may then result from a reshaping of the central cavity into two (or more) independent smaller cavities. Additional geophysical observations will be required to confirm this hypothesis. Downward from the cavity line, the VHSR amplitude is always lower than 0.8 (Figure 3). This area corresponds to the cold part of the glacier [START_REF] Gilbert | The influence of snow cover thickness on the thermal regime of Tête Rousse Glacier (Mont Blanc range, 3200 m a.s.l.): Consequences for outburst flood hazards and glacier response to climate change[END_REF][START_REF] Garambois | Ground-penetrating radar and surface nuclear magnetic resonance monitoring of an englacial water-filled cavity in the polythermal glacier of Tête Rousse[END_REF], where the basal ice temperature is always much lower than that of the melting point, so that the ice layer sticks to the underlying bedrock and no basal melting water is expected.

Upward from the cavity line, the VHSR amplitude generally stands between 0.8 and 1.1. This area is very well spatially correlated with the temperate area of the glacier [START_REF] Gilbert | The influence of snow cover thickness on the thermal regime of Tête Rousse Glacier (Mont Blanc range, 3200 m a.s.l.): Consequences for outburst flood hazards and glacier response to climate change[END_REF][START_REF] Garambois | Ground-penetrating radar and surface nuclear magnetic resonance monitoring of an englacial water-filled cavity in the polythermal glacier of Tête Rousse[END_REF], where the basal ice temperature remains at the fusion point, allowing liquid water to emerge at the ice-bedrock interface, and some detachments to occur locally.

Moreover, within the temperate ice, a large amount of liquid water was detected by SNMR (Figure S6). No large cavity containing liquid water has been revealed by GPR up to now, but small diffuse cavities filled with liquid water, or wet patches, may exist within this volume, generating slight vertical vibrations. These results confirm the general basal conditions of the glacier, as already suggested by GPR and SNMR surveys.

The interpretation of VHSR amplitude has to be performed carefully, mainly due to the temporal variability of the seismic noise, the presence of body and surface waves in the wavefield, and local 2D or 3D resonance effects. In this study, threshold values delineating melting areas have been set

with the help of prior information on Tête-Rousse Glacier. Hence, it is worth noting than VHSR amplitude values are site-dependent, and conclusions based on this proxy should not yet be definitive without additional information.

It should also be noted that this seismic experiment was probably oversized: a shorter-lived, less intensive survey could have provided similar results. Statistical analysis of the temporal variability of the measurements will be the subject of further work, but at this point we estimate that only a few days of ambient noise data is required to conduct this type of experiment on glaciers. This duration might be decreased to a few hours through optimization of the data processing procedure.

Among other advanced techniques, the extraction of Rayleigh wave ellipticity by DOP-E method [START_REF] Berbellini | Constraining S-wave velocity using Rayleigh wave ellipticity from polarization analysis of seismic noise[END_REF][START_REF] Jones | Constraints on the cryohydrological warming of firn and ice in Greenland from Rayleigh wave ellipticity data[END_REF] and its inversion may significantly improve the manuscript submitted to Geophysical Research Letters identification of englacial features (see Supporting Information S6 and Figure S7 for further details).

Conclusion

The spectral analysis performed on continuous seismic records acquired at the surface of the Tête-Rousse Glacier exhibits both a HVSR frequency peak and a polarization that corresponds to the first longitudinal mode of the glacier, confirming the elastic parameters of ice and bedrock. The VHSR exhibits a clear positive anomaly indicating the presence of a subglacial cavity (areas where VHSR > 1.1) at a location that is consistent with the known position of the cavity. However, the VHSR mapping suggests that the cavity's geometry had been modified since previous observations, dating from 2010 and 2018. In addition, the upper part of the glacier, known to be composed of temperate ice, shows VHSR values close to 1, indicating a possible role for diffuse liquid water content and/or local detachment of the basal ice from the bedrock. These results indicate that ambient seismic noise-based spectral analysis of data collected by sensors at the surface can be used to map the presence of potential cavities in glaciers. Computing VHSR amplitude peaks is very straightforward, making this method accessible and rapid. But the level at which VHSR can detect cavities depends on the impedance contrast between the ice and the underlying layer (fractured or intact bedrock, sediments), and will need further investigation to be applicable on glaciers without prior information on their subglacial conditions.

  manuscript submitted to Geophysical Research Letters

Figure 1

 1 Figure 1: a) Mont Blanc massif in the Alps, the red star indicates the position of the Tête-Rousse Glacier. b) outline of the glacier and where it crosses the Mont Blanc Normal route. c) 2022 seismic sensor layout.

  manuscript submitted to Geophysical Research Letters

Figure 2

 2 Figure 2: a) Azimuth-dependent HVSR amplitude from seismometer C2 with time-averaging over the whole period of the experiment in May 2022. b) HVSR amplitude from some nodes along a longitudinal line at the center of the glacier, with time-and azimuth-averaging mean amplitude (dots) and standard deviation (errorbars) over the whole period of the experiment in May 2022.

  manuscript submitted to Geophysical Research Letters

Figure 4 :

 4 Figure 4: Fundamental mode of the glacier, based on the MA of the mechanical model applying simple elliptic geometry adapted from a). The main direction of displacement is depicted by the red arrow. (b) Vertical displacement of the 2D axisymmetric mechanical model at 12.8 Hz, showing the location of points 1, 2, and 3 at the glacier's surface. (c) Vertical displacement occurring at points 1, 2, and 3, recorded between 5 and 20 Hz with the cavity (continuous lines) and without the cavity (dashed lines), as outputs of FRF.

Acknowledgments

The authors wish to thank Sisprobe for renting their SmartSolo nodes and for assistance with data management. The three Centaur & Trillium PH sensors were from the CNRS-INSU-SISMOB pool. The authors gratefully acknowledge experimental help from Simon Bonin, Baptiste Camus, and Agnès Helmstetter. This work was partially funded by the French Ministère de la Transition écologique et de la Cohésion des territoires (Ministry for Ecological Transition and Territorial Cohesion) through its DGPR ROGP program, by the French national research agency (ANR) through its "Plan de Relance" program, and by the Engineering-Consultancy Bureau Géolithe. Open Research manuscript submitted to Geophysical Research Letters Data from seismic stations C1, C2, and C3 are available on RESIF (https://seismology.resif.fr/networks/#/YZ__2022).

Seismic data from nodes are temporarily available on ScienceDB through this link (https://www.scidb.cn/s/faUZf2), before being available permanently on RESIF.


Antoine Benady
email: antoine.benady@ens-paris-saclay.fr

Emmanuel Baranger
email: emmanuel.baranger@ens-paris-saclay.fr

Ludovic Chamoin
email: ludovic.chamoin@ens-paris-saclay.fr

l'erreur en relation de comportement modifiée

published or not. The documents may come

Apprentissage non-supervisé de lois de comportement nonlinéaires avec réseau de neurones thermodynamiquement consistent par minimisation de

Introduction

Depuis quelques années, une tendance visant à représenter des lois de comportements par réseau de neurones émerge. Initialement proposées sans intégrer de physique dans les réseaux [START_REF] Ghaboussi | Neural network material modelling[END_REF], ces approches ce sont petit à petit développées de manière à imposer des contraintes dans l'architecture [START_REF] Amos | Input convex neural networks[END_REF][START_REF] Benady | NN-mCRE : a modified Constitutive Relation Error framework for unsupervised learning of nonlinear state laws with physics-augmented Neural Networks[END_REF], ou à faciliter l'entrainement par l'utilisation de fonction coût ayant un sens physique [START_REF] Raissi | Physics-informed neural networks : A deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations[END_REF]. Plus spécifiquement, de nombreux travaux abordent la question de la représentation de modèles plastiques à l'aide de réseaux neuronaux (voir par exemple [START_REF] Bonatti | On the importance of self-consistency in recurrent neural network models representing elasto-plastic solids[END_REF][START_REF] Gorji | On the potential of recurrent neural networks for modeling path dependent plasticity[END_REF][START_REF] Wu | A recurrent neural networkaccelerated multi-scale model for elasto-plastic heterogeneous materials subjected to random cyclic and non-proportional loading paths[END_REF]). L'essentiel des méthodes proposant d'entrainer des réseaux de neurones pour représenter des lois de comportement utilise un apprentissage dit supervisé, c'est-à-dire faisant appel à des jeux de données de couples déformation-contrainte. Ces méthodes visent souvent à accélérer la résolution de loi d'évolution, ou sont utilisées comme technique d'homogénéisation.

La méthode proposée ici, est au contraire une méthode dite non-supervisée qui ne fait appel qu'à des informations mesurables et vise à améliorer la modelisation du comportement matériau. Ce comportement est ici représenté par des réseaux de neurones qui garantissent par construction le respect de la thermodynamique, permettant ainsi de corriger le biais de modèle. L'entrainement du réseau consiste en la minimisation de l'erreur en relation de comportement modifiée (mCRE), une fonctionnelle composée d'un terme d'erreur de modèle et d'un terme d'écart aux mesures, largement étudiée dans le cas d'identification de paramètres [START_REF] Ladevèze | Updating of finite element models using vibration tests[END_REF].

Une première partie de ce document vise à définir le cadre thermodynamique utilisé (Section 2). La Section 3 présente le cadre mCRE pour le problème d'apprentissage. Dans la Section 4 les spécifités liées aux réseaux de neurones sont détaillées (garanties thermodynamiques et réglage automatique des hyperparamètres). Enfin, un cas d'application de la méthode est présenté en Section 5.

Cadre thermodynamique

Cette section vise à introduire le cadre thermodynamique utilisé dans cet article. Dans cette formulation, la loi de comportement est formulée en termes d'équation d'état et de loi d'evolution par l'intermédiaire de deux potentiels ψ et φ, comme c'est le cas dans le cadre des matériaux standards généralisés. Cette section formule des contraintes sur les fonctions ψ et φ qui permettent de respecter automatiquement les principes de la thermodynamique. Tout d'abord, le premier principe s'écrit localement :

ρ ė = σ : ε + r -∇ • q (1)
avec e l'énergie interne spécifique, σ le tenseur des contraintes, ε le tenseur des déformations linéarisé, r une source volumique de chaleur et q le flux de chaleur. Le second principe s'écrit quant à lui :

ρ ṡ + ∇ • ( q θ ) - r θ ≥ 0 (2)
avec s la densité d'entropie spécifique. Après avoir introduit l'énergie libre de Helmholtz ψ = ρ(e-θs), l'inégalité de Clausius-Duhem est obtenue en combinant les équations (1) et (2) :

σ : ε -ψ - q • ∇θ θ ≥ 0 (pour un processus isotherme) (3)
Pour les problèmes de plasticité, le tenseur de déformation totale est écrit avec la décomposition suivante : ε = ε e + ε p , où ε e est la déformation élastique et ε p la déformation plastique.

Des variables internes X k (k = 1, . . . , K), regroupées dans le vecteur X, sont introduites pour représenter d'autres phénomènes (comme l'écrouissage, par exemple). Les équations d'état sont dérivées de l'énergie libre de Helmholtz utilisée comme potentiel thermodynamique ψ(ε, ε p , X).

Pour une transformation élastique réversible avec une température homogène, on a σ = ∂ψ ∂εe correspondant à la loi de Hooke. Par analogie, les forces thermodynamiques Y k regroupées dans Y sont définies de manière à être associées aux variables internes X avec Y = ∂ψ ∂X . L'utilisation d'une fonction convexe ψ satisfait automatiquement l'inégalité de Clausius-Duhem qui peut être réécrite comme suit :

(σ - ∂ψ ∂ε e ) : εe + σ : εp - K k=1 ∂ψ ∂X k Ẋk - q • ∇θ θ ≥ 0 (4)
Un potentiel de dissipation convexe, non négatif et nul à l'origine φ(ε p , Ẋ) peut être postulé pour définir des lois d'évolution garantissant la positivité de la dissipation mécanique D = σ : εp -Y : Ẋ. Dans la suite de cet article, les potentiels ψ et φ sont utilisés :

-dans la définition de l'erreur en relation de comportement modifiée (Section 3) ; -pour la formulation de loi de comportement sous la forme de réseau de neurones. Les réseaux de neurones représentant les fonctions ψ et φ sont choisies pour garantir par construction les propriétés énoncées, comme cela sera détaillé dans la Section 4.

3 L'erreur en relation de comportement modifiée comme fonction coût du problème d'apprentissage -Admissibilité cinématique, qui définit l'espace U ad des champs de déplacement satisfaisant les conditions aux limites de Dirichlet :

u |∂Ω 1 = u d (5)
-Admissibilité statique, qui définit l'espace S ad des champs de contrainte satisfaisant l'équilibre :

Ω σ : ε(v) = Ω f v d • v + ∂Ω 2 f s d • v ∀v ∈ U 0 (6)
où U 0 est l'espace d'admissibilité cinématique avec des conditions de Dirichlet homogènes.

Le tenseur de déformation est décomposé en les parties élastique ε e et plastique ε p telles que ε = ε e + ε p -Comportement :

-équations d'état :

s = ∂ψ ∂e e (7)
-lois d'évolution : Les potentiels ψ et φ sont paramétrés avec les paramètres p. Ces paramètres, qui sont traditionnellement des paramètres de modèles pré-établies, sont ici les poids et biais de réseaux de neurones. La minimisation de la mCRE vise à identifier les paramètres p des relations constitutives qui s'ajustent le mieux aux données expérimentales (des observations partielles et bruitées de déformation ε obs ou de déplacement u obs ). Ici, les conditions de Dirichlet et de Neumann sont supposées connues (bien que le cadre s'étende aux conditions limites peu fiables).

s = ∂φ ∂ ėp (8)

L'erreur en relation de comportement modifiée

Le concept d'erreur en relation de comportement (CRE) a été introduit dans les années 70 dans le contexte de la vérification des modèles éléments finis [START_REF] Ladevèze | Error estimate procedure in the finite element method and applications[END_REF]. L'idée générale de cette mesure d'erreur est de diviser les équations de la Section 3.1 en fonction de la fiabilité de l'information. La relation de comportement est considérée comme l'équation la moins fiable et étant donné un couple admissible (û, σ) ∈ (U ad × S ad ), la CRE est définie comme suit pour les problèmes indépendants du temps :

E 2 CRE (û, σ) = Ω ψ(û) + ψ * (σ) -σ : ε(û) (9)
La notation ψ * désigne la transformation de Legendre-Fenchel de ψ. Dans le cas où ψ est un potentiel quadratique décrivant la loi de Hooke en élasticité linéaire (ψ(ε) = 1 2 ε :

K : ε et ψ * (σ) = 1 2 σ : K -1 : σ), la relation de comportement est σ = K : ε(u) et la CRE s'écrit E 2 CRE (û, σ) = Ω (σ -Kε(û)) : K -1 : (σ -Kε(û)) = ||σ -Kε(û)|| 2 K -1 (10)
où ||•|| K -1 est la norme énergétique sur les champs de contrainte. À partir de [START_REF] Ladevèze | Constitutive relation errors for f.e. analysis considering (visco-) plasticity and damage[END_REF], il apparait que la CRE est un indicateur, pour un couple admissible, de la non-satisfaction de la relation de comportement.

Plus tard, toujours dans le contexte de l'estimation de l'erreur a posteriori, le concept de CRE a été adapté aux problèmes non-linéaires dépendants du temps. En présence de problèmes dépendants du temps, la CRE s'écrit pour une solution admissible ŝ = (ε e , ε p , σ, X, Y) comme suit :

E 2 CRE (ŝ) = T 0 Ω η ψ (ε e , σ, X, Y) + T 0 t 0 Ω η φ (ε p , σ, Ẋ, Y) (11)
avec,

η ψ (ε e , σ, X, Y) = ψ(ε e , X) + ψ * (σ, Y) -σ : ε e -X : Y (12) η φ (ε p , σ, Ẋ, Y) = φ(ε p , -Ẋ) + φ * (σ, Y) -σ : εp + Ẋ : Y (13)
Le concept CRE a ensuite été adapté dans les années 90 pour aborder les problèmes inverses [START_REF] Maia | Error localization for updating finite element models using frequency-response-functions[END_REF] ; les paramètres de la loi de comportement sont cherchés de manière à minimiser le CRE. Dans la première idée de la méthode, les observations étaient directement incluses dans la définition de l'espace admissible, ce qui conduisait à une faible robustesse au bruit de mesure. Pour améliorer la robustesse au bruit, la procédure modifiée (mCRE) [START_REF] Ladevèze | Application of a posteriori error estimation for structural model updating[END_REF] n'impose pas les observations dans l'espace admissible. Seules les informations fiables sont imposées par construction, tandis que les informations peu fiables sont relâchées. Récemment, le concept mCRE a été étendu au cas des lois d'évolution dans [START_REF] Marchand | Parameter identification and model updating in the context of nonlinear mechanical behaviors using a unified formulation of the modified constitutive relation error concept[END_REF]. Le présent article s'appuie sur le cadre mCRE développé dans [START_REF] Marchand | Parameter identification and model updating in the context of nonlinear mechanical behaviors using a unified formulation of the modified constitutive relation error concept[END_REF] et est adapté au cas où la loi de comportement est décrite par des réseaux de neurones.

Le problème inverse consiste à trouver les paramètres optimaux p opt (impliqués dans les potentiels thermodynamiques) tels que :

p opt = argmin p min ŝ∈A d E 2 mCRE (ŝ ; p) (14)
avec

E 2 mCRE (ŝ ; p) = E 2 CRE (ŝ) + α 2 T 0 ||Πu -u obs || 2 (15)
avec A d = (U ad × S ad ) (bien qu'il soit possible de considérer les conditions aux limites comme une information incertaine), α une pondération entre les termes de la fonction coût, et Π un projecteur de u sur les quantités mesurées.

La minimisation de la mCRE est effectuée avec un processus itératif dans lequel, à chaque itération :

-dans une première étape, une solution admissible ŝ est calculée pour les paramètres actuels p (n) de telle sorte que :

ŝ(n) = argmin ŝ∈A d E 2 mCRE (ŝ ; p (n) ) (16)
Cette étape, la plus coûteuse en terme de temps de calcul, est effectuée selon une stratégie similaire à celle utilisée dans [START_REF] Marchand | Parameter identification and model updating in the context of nonlinear mechanical behaviors using a unified formulation of the modified constitutive relation error concept[END_REF]. Cette stratégie s'inspire de la méthode LATIN [START_REF] Ladevèze | Principles of the Method of Large Time Increments[END_REF], qui est non-incrémentale (c'est-à-dire globale en temps) et convient bien à la structure mathématique de la mCRE. Le choix fait ici est de diviser la mCRE en deux parties positives, E 2 ψ (ε e , σ, X, Y) et E 2 φ (ε p , σ, Ẋ, Y), définies par :

E 2 ψ (ε e , σ, X, Y) = T 0 Ω η ψ (ε e , σ, X, Y) + α 2 T 0 ||Πu -uobs|| 2 (17) E 2 φ (ε p , σ, Ẋ, Y) = T 0 t 0 Ω η φ (ε p , σ, Ẋ, Y) (18)
Cette séparation en deux parties positives permet de minimiser chaque terme alternativement. L'équation (17) est un compromis entre le résidu des équations d'état et l'écart par rapport aux mesures. Sa minimisation est un problème global en espace et linéaire. La minimisation de (18), réalisée ici avec un schéma d'Euler, est locale en espace et correspond à l'intégration des lois d'évolution. -dans une deuxième étape, les paramètres de la loi de comportement sont mis à jour en suivant une étape de descente de gradient :

p (n+1) = p (n) -l r E 2 mCRE (ŝ (n) ; p (n) ) dp (19)
Dans cette étape, le gradient est calculé par la méthode de l'état adjoint.

4 Réseaux de neurones thermodynamiquement consistents pour l'apprentissage de loi de comportement

La précédente section a brievement détaillé le principe de l'erreur en relation de comportement modifiée ainsi que sa minimisation pour l'identification de paramètres de lois de comportement. Ce cadre est adapté aussi bien à l'identification de paramètres d'un modèle donné qu'à l'entrainement d'un réseau de neurones (ce qui revient à identifier les poids et le biais). L'utilisation d'un réseau de neurones pour représenter le comportement permet de relacher la forme du modèle (car un réseau de neurones est un approximant universel) et donc le biais de modèle. Néanmoins, le passage de l'identification de paramètres à l'entrainement de réseau de neurones, bien que conceptuellement proche, soulève des questions importantes :

-la question de la cohérence lors de la phase d'inférence par rapport aux exigences physiques telles que la thermodynamique ; -En comparaison avec l'identification de paramètres d'un modèle donné, un réseau de neurones comporte bien plus de paramètres à identifier, ce qui rend la tâche d'optimisation plus difficile, notamment à cause du nombre important d'hyperparamètres à régler tels que le taux d'apprentissage, le nombre d'epochs, la taille de batch, etc. Ces questions sont traitées dans les deux prochaines sous-sections.

Imposer les contraintes thermodynamiques dans l'architecture

Dans la Section 2, un cadre a été défini pour la formulation de lois de comportement compatibles avec la thermodynamique. Les propriétés permettant de satisfaire automatiquement la thermodynamique peuvent s'imposer dans le choix d'architecture du réseau. Ceci permet, contrairement à une stratégie de pénalisation en phase d'entrainement, de garantir le respect de ces propriétés en phase d'inférence. La positivité est imposée grâce à l'utilisation de fonctions d'activation positives. En ce qui concerne la convexité, l'architecture "input convex neural network " (ICNN) proposée dans [START_REF] Amos | Input convex neural networks[END_REF] est utilisée. Cette architecture utilise des fonctions d'activation convexes et non décroissantes, ainsi que des contraintes de positivité sur les poids intermédiaires. Comme la composition d'une fonction convexe et d'une fonction convexe non décroissante est convexe, et que la somme de fonctions convexes est également convexe, cette architecture garantit la convexité. Enfin, pour imposer une valeur nulle à l'origine, pour n'importe quelle fonction g, on définit :

g(x) = g N N (x) -g N N (x = 0) (20)
où g N N (x) est la sortie du réseau de neurones évaluée en x.

Reglage automatique des hyperparamètres

Dans le contexte actuel de jumeau numérique, dans lequel le jumeau hybride interagit avec la structure physique sans intervention humaine, il n'est pas possible d'imaginer une stratégie basée sur le réglage manuel des hyperparamètres. Ainsi, une stratégie de sélection automatique et adaptative des hyperparamètres est proposée.

Initialisation du réseau

Dans le cas de l'apprentissage de lois de comportement et en particulier lorsque les variables internes sont choisies a priori, une connaissance partielle du modèle peut être disponible. Par conséquent, avant de commencer l'entraînement visant à minimiser la mCRE, le réseau est initialisé de manière à représenter cette connaissance partielle. Ceci permet de réduire la sensibilité d'une initialisation aléatoire ainsi que d'initialiser le réseau plus proche de la solution à apprendre.

Réglage de la pondération entre les fonctions coût

La fonctionnelle mCRE offre une interprétation physique intéressante qui peut être utilisée pour régler le paramètre α. L'idée derrière ce réglage est que le modèle doit être mis à jour de manière jusqu'à ce que l'écart entre les prédictions et les observations atteigne le niveau de bruit présent dans les données, conformément au principe de Morozov. Cette hypothèse repose sur la connaissance a priori du niveau de bruit. Le paramètre α est réglé de manière adaptative au cours de l'apprentissage de manière à vérifier ce critère.

Réglage du taux d'apprentissage

Un autre paramètre essentiel dans l'entraînement des réseaux de neurones est la sélection du taux d'apprentissage, noté l r . Une valeur excessivement faible pour l r peut entraîner un processus d'entraînement lent, tandis qu'une valeur trop grande peut introduire de l'instabilité pendant l'entraînement, potentiellement empêchant le réseau d'atteindre la convergence. Le choix du taux d'apprentissage peut être motivé une analyse des progrès réalisés au cours d'une seule étape d'entraînement. Dans ce contexte, un indicateur de la vitesse de mise à jour de l'entraînement peut être obtenu à partir de la procédure de minimisation en deux étapes. Le nombre d'itérations effectuées à l'étape 1 est directement lié à la valeur de la mise à jour effectuée à l'étape 2 de l'epoch précédente. Par conséquent, à partir de la deuxième epoch, l'apprentissage est mis à jour de manière à ce que le nombre d'itérations effectuées soit proche d'un nombre d'itérations cible défini de manière empirique.

Critère d'arrêt de l'entrainement

L'entraînement des réseaux de neurones nécessite l'établissement d'un critère d'arrêt. Le concept de mCRE offre une base physique intéressante qui peut être utilisée pour établir une condition d'arrêt. En effet, le terme de CRE peut-être interprété comme une erreur de modèle, et étant homogène à une énergie il peut donc être comparé à l'énergie dans la structure. Le crière d'arrêt développé consiste à choisir une valeur d'erreur de modèle normalisée cible et d'arrêter l'entrainement quand cette valeur est atteinte.

Résultats

Cette section vise à illustrer la méthode sur un cas test mettant en oeuvre l'apprentissage d'une loi d'écrouissage non linéaire.

Pour simplifier, une partie du modèle est supposée connue, comme indiqué dans le Tableau 1. Il s'agit d'un choix de modélisation en fonction de la fiabilité des connaissances sur le matériau : d'autres stratégies peuvent être envisagées, telles que l'identification de certains paramètres du modèle ou même la représentation de chaque partie du modèle avec des réseaux de neurones. La stratégie mCRE s'adapte facilement à ces cas. Le cas étudié concerne une poutre en 2D soumise à un chargement de traction variable. Les mesures synthétiques proviennent de deux fibres optiques orientées le long des deux axes principaux et positionnées au milieu de la poutre.

La Figure 1 représente l'idée initiale, la vraie solution et le résultat après l'apprentissage. Initialement, une erreur de modèle est présente car la loi d'écrouissage supposée est linéaire alors que la vraie solution a une forme exponentielle. Après l'apprentissage, la forme de modèle est bien apprise. Ceci n'aurait pas été possible si la forme de modèle n'avait pas été relachée comme c'est la cas avec l'utilisation de réseau de neurones.

La Figure 2 présente la convergence des différents termes composant la mCRE. À la fin de l'entraînement, l'écart aux données normalisée est d'environ 1, ce qui signifie que le critère de

ψ p (ε p ) = Ap + A B (e -Bp -1)
Partie plastique de la loi d'état : représentée par un réseau de neurones initialisé avec

ψ p (ε p ) = 1 2 hp 2
Limite du domaine d'élasticité :

f = σ eq -(R + R 0 )
Limite du domaine d'élasticité :

supposée connue avec les bons paramètres

I 0 W E x v l y 4 R k h C o J t C 3 N l 3 5 G A Z K I = " > A A A I c X i c r V X b b t N A E J 2 2 Q E 2 4 p f C E e F k 1 Q i p C i u I K W h 4 r E B K P Q S F p p a a K v O 4 6 s e o b v g Q q K 3 + B x N f w C t / A H 8 A z P 8 D M 2 C F N 4 6 R b C U d J Z s + c c z w z 6 5 V l 5 L l J 2 m r 9 X F v f u H H z 1 q Z x u 3 b n 7 r 3 7 D + p b D 3 t J m M W 2 6 t q h F 8 Z H 0 k q U 5 w a q m 7 q p p 4 6 i W F m + 9 N S h P H t D + c O x i h M 3 D D 6 k 5 5 E 6 8 a 1 h 4 D q u b a U I D e p m X 8 r E H f r W I D f N i e j 7 M v y c C z c Q T h b Y R B G h I x C X U k U J c w b 1 R q v Z 4 k s s B m Y Z N K C 8 2 u H W x h f o w y m E Y E M G P i g I I M X Y A w s S / B y D C S 2 I E D u B H L E Y I 5 f z C i Z Q Q 2 2 G L I U M C 9 E z / B 3 i 6 r h E A 1 y T Z 8 J q G + / i 4 T d G p Y C n q A m R F 2 N M d x O c z 9 i Z 0 G X e O X t S b e f 4 L 0 s v H 9 E U R o h e p Z s y d X X U S w o O v O I e X O w p Y o S 6 s 0 u X j K d C l Y s L X a X o E C F G 8 S n m Y 4 x t V k 7 n L F i T c O 8 0 W 4 v z v 5 h J K K 3 t k p v B 7 5 V V S v w M W a k Q 9 b j a a X e U o 9 k v 3 7 E c G f 4 V e f L + i L W s Z k 1 9 i B M g 9 o n 3 1 e d a A p x D j n j K y o R 3 g d Y j r p M i i Z N + i + h k p Z 7 u U j g U + o t 9 0 r r I u f y 0 B V p u x B 4 y u l j T Z f c Z d 1 J 6 r 3 Y v O p p 3 c n C t 0 6 c D n U p t R 0 v b r t S 2 t b R R p T b S 0 m a V 2 k x L O 6 7 U j r W 0 3 U p t V

Conclusion

Une méthode a été proposée pour l'apprentissage non-supervisé de lois de comportement décrites par des réseaux de neurones. Cette méthode s'appuie sur l'utilisation de réseaux de neurones qui garantissent par construction le respect de la thermodynamique, grâce à la formulation de lois sous la forme de potentiels convexes. Le cadre de la minimisation de l'erreur en relation de comportement modifiée est particulièrement adapté à la correction de bias de modèle de part la structure de sa fonctionnelle. Par ailleurs, le sens riche sens physique de cette fonctionnelle facilite grandement l'automatisation des réglages des hyperparamètres. La méthode a été illustrée sur le cas test d'apprentissage d'une loi d'écrouissage non-linéaire.

Financement du projet

Ce projet a reçu un financement du European Research Council (ERC) dans le cadre du programme de recherche et d'innovation Horizon 2020 de l'Union européenne (convention de subvention n°101002857).

  avec les variables e e = [ε e , X], e p = [ε p , -X] et des forces thermodynamiques globales s = [σ, Y].

  Modèle utilisé pour générer les données Modèle à apprendre Partie élastique : ψ e (ε e ) = 1 2 (λ(tr ε e ) 2 + 2µε e : ε e ) Partie élastique : supposée connue avec les bons paramètres Partie plastique de la loi d'état :

0 v b q 9 T 2 N

 92 H e 4 s / K 5 m e 2 Y 7 n M 9 n e A M 2 a m c 5 7 N r O I 4 q PZ f V W + y W v v 9 Y o + L x t R x n 5 / f / n P F C 4 f F 7 Q C J u / a u Y T u g O R z 4 j 8 5 y i 4 h q + W 8 3 L b 9 L F o L f b N P e a L 9 + / a B y 8 L t + y B j y B b f Q 3 Y R 8 O 4 B 2 e v y 5 W + B W + w X f 4 s f n H e G w I Y 7 u g r q + V m k c w d x n P / w K s J L o w < / l a te x i t > 11 in function of ✏ 11

Figure 1 -

 1 Figure 1 -Comportement avant et après l'entrainement

Figure 2 -

 2 Figure2-Évolution de la CRE normalisée au cours de l'entrainement (gauche), évolution de l'écart aux mesures normalisé par le bruit de mesure (milieu), évolution de la mCRE durant l'entrainement (droite).

  Pour définir les notations du problème, considérons un corps dans une configuration initiale Ω ⊂ R d (d = 1, 2, 3) avec une frontière ∂Ω et un environnement isotherme, observé sur une période [0, T ] et soumis à dess petites perturbations. Des conditions aux limites de Dirichlet sont imposées sur ∂Ω 1 ⊂ ∂Ω au moyen d'un champ de déplacement dépendant du temps u d .Des conditions aux limites de Neumann sont prescrites sur ∂Ω 2 ⊂ ∂Ω au moyen d'un champ de forces surfaciques dépendant du temps f d s . Un champ de forces volumiques dépendant du temps f d v peut également être prescrit dans Ω. De plus, des mesures de déformation ε obs (dans le cas d'observations à partir de fibres optiques) ou des mesures de déplacement u obs (dans le cas de la corrélation d'images numériques) sont disponibles.La solution du problème mécanique direct est le couple (u, σ) de champs de déplacement et de contrainte de Cauchy qui satisfont les trois groupes d'équations suivants pour chaque instant t ∈ [0, T ] :

	3.1 Définition du problème d'apprentissage

Table 1 -

 1 Résumé du comportement à apprendre, avec λ = et σ eq la contrainte équivalente de Von Mises.Morozov est satisfait. Pour la CRE normalisée, la valeur est inférieure à la valeur cible définie par l'utilisateur. Comme ces deux critères sont satisfaits l'entraînement s'arrête, conformément à la procédure décrite en Section 4. Après l'epoch 500, la CRE normalisée présente un comportement oscillatoire qui est expliqué par la mise à jour du paramètre α.

	Eν (1+ν)(1-2ν) , µ = E 2(1+ν)

< l a t e x i t s h a 1 _ b a s e 6 4 = "


Harmonic control of three-phase AC/DC converter with time domain guarantees

Maxime Grosso, Pierre Riedinger, Jamal Daafouz, Serge Pierfederici, Hicham Janati Idrissi, Blaise Lapôtre

To cite this version:

Maxime Grosso, Pierre Riedinger, Jamal Daafouz, Serge Pierfederici, Hicham Janati Idrissi, et al..

  Harmonic control of three-phase AC/DC converter with time domain guarantees. 2023. hal