Siegelman Lia

Klein Patrice

Rivière Pascal

Andrew F Thompson

Hector S Torres

Flexas Mar

Menemenlis Dimitris

Lia Siegelman
email: lsiegelman@caltech.edu

Patrice Klein

Pascal Rivi Ère

Dimitris Menemenlis

Enhanced upward heat transport at deep submesoscale ocean fronts

The ocean is the largest solar energy collector on Earth. The amount of heat it can store is modulated by its complex circulation, which spans a broad range of spatial scales, from metres to thousands of kilometres. In the classical paradigm, fine oceanic scales, less than 20 km in size, are thought to drive a significant downward heat transport from the surface to the ocean interior, which increases oceanic heat uptake. Here we use a combination of satellite and in situ observations in the Antarctic Circumpolar Current to diagnose oceanic vertical heat transport. The results explicitly demonstrate how deep-reaching submesoscale fronts, with a size smaller than 20 km, are generated by mesoscale eddies of size 50-300 km. In contrast to the classical paradigm, these submesoscale fronts are shown to drive an anomalous upward heat transport from the ocean interior back to the surface that is larger than other contributions to vertical heat transport and of comparable magnitude to air-sea fluxes. This effect can remarkably alter the oceanic heat uptake and will be strongest in eddy-rich regions, such as the Antarctic Circumpolar Current, the Kuroshio Extension and the Gulf Stream, all of which are key players in the climate system.

Vertical heat transport is one of the key mechanisms by which the ocean regulates Earth's climate. Munk and Wunsch 1 first postulated that vertical heat transport is balanced by an upward component due to the large-scale mean flow (>300 km) and a downward component explained by fine-scale diffusive processes (<20 km). However, recent studies 2,3 highlighted the importance of mesoscale eddies (50-300 km) for vertical heat transport. Thus, global vertical heat transport by the mean flow and diffusion are now thought to be both downwards and balanced by an upward eddy heat flux, with the mean and eddy components generally being the largest contributors 2 .

Mesoscale eddies are intensified in energetic areas such as the Antarctic Circumpolar Current (ACC), Kuroshio Extension and Gulf Stream 4 . They are known to drive the production of submesoscale fronts (<20 km), i.e. regions that separate waters of different densities [5][START_REF] Ferrari | A frontal challenge for climate models[END_REF][7][START_REF] Mcwilliams | Submesoscale currents in the ocean[END_REF][START_REF] Lévy | The role of submesoscale currents in structuring marine ecosystems[END_REF] , ubiquitous on satellite images of ocean color (Fig. 1). To date, submesoscale fronts are mainly thought to be confined to the ocean surface mixed-layer (⇠100 m deep) [START_REF] Mcwilliams | Submesoscale currents in the ocean[END_REF][START_REF] Lévy | The role of submesoscale currents in structuring marine ecosystems[END_REF] . This is because, in the classical paradigm, motions below the mixed-layer are broadly assumed to be in geostrophic balance (a balance between Coriolis and horizontal pressure forces), preventing the formation of strong density gradients [START_REF] Vallis | Atmospheric and oceanic fluid dynamics[END_REF] . As a consequence, fine scales below the mixed-layer are still principally associated with diffusive processes that drive downward heat transport. However, a handful of studies suggest that submesoscale fronts associated with upward heat transport may also penetrate below the mixed-layer [START_REF] Taylor | Ocean fronts trigger high latitude phytoplankton blooms[END_REF][START_REF] Thomas | Symmetric instability in the gulfstream[END_REF][START_REF] Bachman | Modelling of partially-resolved oceanic symmetric instability[END_REF][START_REF] Yu | An annual cycle of submesoscale vertical flow and restratification in the upper ocean[END_REF] . Yet, these studies do not explain the formation mechanism of these submesoscale fronts because they are confined to small domains (<50 km). Observational evidence over domains sufficiently large to capture multiple mesoscale eddies and their associated submesoscale fronts are needed for the closure of oceanic vertical heat transport. Such measurements are still lacking due to the fact that capturing fine spatial and temporal scales over extended domains lie at the edge of our observational capabilities.

Here, we meet this observational challenge by analyzing a unique dataset collected by an instrumented elephant seal in the ACC.

The seal data reveal the presence of numerous deep-reaching submesoscale fronts that extend well below the mixed-layer and are principally found on eddies' edges. We then develop a pioneering methodology, combining satellite and seal observations, to retrieve a 3-D synoptic view of ocean dynamics from meso-to submesoscale. Results explicitly demonstrate that deep-reaching submesocale fronts are generated by the strain field associated with co-interacting mesoscale eddies. By invoking properties of mesoscale turbulence, we explain why deep-reaching fronts, counter-intuitively, drive enhanced vertical heat transport below the surface mixed-layer that is directed upward, which is also supported by a numerical model analysis (see Supplementary Information). Associated vertical heat transport is larger than the mean flow and eddy contributions 2 , and of similar magnitude to air-sea fluxes [START_REF] Large | The global climatology of an interannually varying air-sea flux data set[END_REF] . This effect is argued to crucially limit oceanic heat uptake and therefore to alter the ocean heat storage capacity.

Frontal region sampled by elephant seal and satellite

Measurements were collected by a southern elephant seal (Mirounga leonina) in the Kerguelen area (Indian sector of the Southern Ocean, Fig. 1) during the austral spring and summer. The dataset has a horizontal resolution of 1 km, a vertical resolution of 1 m and extends from the surface down to 400 m in the ocean interior. Over a period of more than 3 months and a distance exceeding 5000 km, the seal continuously recorded temperature, conductivity and pressure, from which buoyancy b (of opposite sign to density) is derived (see Methods). An example of the buoyancy field is shown in Fig. 2a-c. In addition, satellite altimeter observations supply daily horizontal fields of Sea Surface Height (SSH), a proxy of surface pressure, with an effective resolution sufficient to resolve mesoscale eddies with a size of 50 to 100 km.

The area sampled by the seal from km 1000 to 2300 (in red in Fig. 1), i.e., from 11/02-24/2014, is rich in eddies and lies within the energetic ACC [START_REF] Hogg | Recent trends in the southern ocean eddy field[END_REF] as illustrated by the SSH contours in Fig. 1. The seal crossed numerous co-interacting cyclonic (low SSH) and anti-cyclonic (high SSH) mesoscale eddies. Geostrophic currents derived from SSH reach values of 1 m/s around the eddy edges (Fig. 2a), similar to those found in the Gulf Stream and the Kuroshio Extension, the two other most energetic ocean currents 4 . This area is thus representative of the flow field encountered within eddy-rich regions of the world oceans, and will receive particular emphasis in this study. The other areas crossed by the seal (Fig. 1) are much less energetic (Extended Data Fig. 1).

The comparison between satellite SSH and buoyancy anomalies sampled by the seal in the turbulent area (Figs. 2a-c and S1,a-c) highlights the expected vertical structure of the mesoscale eddies 10 : buoyancy anomalies are positive and bowl-shaped in anticyclonic eddies (high SSH) and negative and reversed bowl-shaped in cyclonic eddies (low SSH, Figs. 2a-c and Extended Data Fig. 1a-c). The combination of satellite observations of SSH and seal measurements of buoyancy provides a synoptic 3-D view of the flow field encountered by the seal, and in particular of the eddy field, that can extend down to depths of at least 400 m. We refer the reader to the Supplementary Information for a more detailed analysis of the consistency at mesoscale between these two independent datasets.

Frontal structures, or buoyancy fronts, are identified by the along-track derivative of buoyancy and are defined as b s = ∂ b/∂ s, with s the curvilinear abscissa, i.e. the along-track direction. The fronts, shown in Fig. 2d, have a width between 5 and 20 km and are thus submesoscale features. Indeed, the Rossby radius of deformation in the Kerguelen area is of ⇠15 km [START_REF] Chelton | Geographical variability of the first baroclinic rossby radius of deformation[END_REF] , which corresponds to an eddy radius of ⇠50 km [START_REF] Chelton | Geographical variability of the first baroclinic rossby radius of deformation[END_REF][START_REF] Chelton | Satellite measurements reveal persistent small-scale features in ocean winds[END_REF] , consistent with the SSH observations (Fig. 2a). These fronts are well resolved by the seal's 1-km horizontal resolution measurements. They are more numerous below the mixed-layer and extend down to at least 400 m. Large buoyancy gradients are preferentially found at the edges of the mesoscale eddies and in-between them.

Buoyancy gradient magnitudes reach values larger than 5 ⇥ 10 7 s 2 . The root mean square of the lateral gradient of buoyancy, used as an indicator of the gradient magnitude, is larger than 0.5 ⇥ 10 7 s 2 , regardless of depth (Extended Data Fig. 2a).

Compared to the few other existing submesoscale-resolving datasets [START_REF] Thompson | Open-ocean submesoscale motions: A full seasonal cycle of mixed layer instabilities from gliders[END_REF][START_REF] Du Plessis | Submesoscale processes promote seasonal restratification in the subantarctic ocean[END_REF][START_REF] Viglione | Abrupt transitions in submesoscale structure in southern drake passage: glider observations and model results[END_REF] , these values are of the same order of magnitude, highlighting the important frontal activity of the area.

Frontal dynamics in the ocean interior

Fig. 3 illustrates how the production of horizontal gradients of buoyancy is driven by a pure strain field, which corresponds to the hyperbolic regions in-between co-interacting eddies 5 , as can be seen in Fig. 1. In the schematic Fig. 3a, the strain field (black arrows) stretches a tracer patch in the y-direction and compresses it in the x-direction. This leads to the formation of strong horizontal gradients of buoyancy, or fronts, at submesoscale that are associated with a growth rate related to the strain rate u x = ∂ u ∂ x (Fig. 3b).

To assess the relevance of this mechanism in the Kerguelen area, we analyze the background strain field in relation to the observed buoyancy gradients. We use daily Finite Size Lyapunov Exponents (FSLE) [START_REF] Ovidio | Mixing structures in the mediterranean sea from finite-size lyapunov exponents[END_REF] computed from satellite-derived geostrophic velocities to characterize the strain field properties (see Methods). FSLE indicate both the orientation and time scale (colorbar in Fig. 4a) of the stretching and compression induced by the strain field (respectively the red and blue curves in Figs. 3a and4a). As illustrated in Fig. 3a, we expect tracer patch, or equivalently buoyancy anomalies, to be aligned with stretching FSLE [START_REF] Ovidio | Mixing structures in the mediterranean sea from finite-size lyapunov exponents[END_REF] . The spatial distribution of FSLE (Fig. 4a and Extended Data Fig. 3) confirms that regions in-between eddies and on eddy edges are prone to the formation of strong horizontal buoyancy gradients. However, the seal's trajectory is often oblique rather than perpendicular to the buoyancy fronts, which may lead to an underestimation of the front's magnitude. Therefore, to correct for the seal's orientation with respect to the FSLE it encounters, buoyancy gradients are first normalized by the angle between the seal's trajectory and the FSLE direction (see Methods). Note that normalized buoyancy gradients are now referred to as b x , where the subscript x refers to the direction normal to the stretching FSLE (Fig. 3a). We then compare the along-track time series of stretching FSLE estimated along the seal's path with the normalized buoyancy gradients at 300 m (comparison with buoyancy gradients at different depths below the mixed-layer produces similar results). There is good agreement between FSLE and normalized buoyancy gradients at 300 m (Fig. 4 and Extended Data Fig. 4a), as well as a consistent relationship linking the two (Fig. 4c). The two time-series are strongly intermittent. Peaks of FSLE and buoyancy gradients are often co-located (Fig. 4a-b and Extended Data Fig. 4a). Note that similar results are obtained with non-normalized gradients, which demonstrates that no bias is added via the normalization. Combined together, these two datasets provide observational evidence of how an eddy strain field, diagnosed at the sea surface, drives strong buoyancy gradients at submesoscale in the ocean interior.

Comparison between panels c and d in Fig. 2 further reveals that, whereas mesoscale buoyancy anomalies are slanted on the vertical, as evidenced by the bowl-shaped structures of buoyancy (panel c), buoyancy fronts at submesoscale are almost vertical (panel d). This density slope property is related to the dynamical regime that drives these structures and can be quantified by the non-dimensional Richardson number Ri (see Methods). The small Ri-values observed in Fig. 2d suggest an energetic ageostrophic regime associated with intense vertical currents. The emergence of these vertical currents, w, counter-balances the formation of sharp submesoscale fronts generated by the ambient strain field. This mechanism is referred to as frontogenesis [START_REF] Hoskins | A new look at the w-equation[END_REF] and is sketched in Fig. 3b. It leads to an equilibrium that is captured by the classical omega equation, which involves the vertical current field w, the strain field u x , and the lateral gradient of buoyancy field b x . Here, we diagnose the vertical currents w by solving the 2-D (x,z) quasi-geostrophic version of the omega equation [START_REF] Hoskins | A new look at the w-equation[END_REF] (see Methods), with b x obtained from seal observations and u x estimated from stretching FSLE (Fig. 4). However, compared to the full omega equation, there is a tendency for the vertical currents diagnosed here to be underestimated by ⇠ 1.4 (see Methods).

The vertical section of w shown in Fig. 5b reveals positive and negative w-values with large magnitudes of up to 100 m/day, i.e.

almost an order of magnitude greater than what is attributed to mesoscale eddies alone [START_REF] Mcgillicuddy | Mechanisms of physical-biological-biogeochemical interaction at the oceanic mesoscale[END_REF] . Vertical currents have a width of 5 to 10 km. They are intensified in the ocean interior, below the mixed-layer down to at least 400 m, and do not necessarily penetrate into the mixed-layer. The continuous vertical extent of these vertical currents highlights the crucial role played by deep-reaching submesoscale fronts in connecting the ocean interior to the surface. The w-field is characterized by the same intermittency present in the FSLE and buoyancy gradient fields (Figs. 2d, 4c , 5b and Extended Data Fig. 4b). Large w-magnitudes are co-located with strong buoyancy gradients and intense FSLE, and are therefore mostly found on the edges of, and in-between, eddies.

The analysis of the vertical currents presented here is supported by a comparison with a realistic regional numerical model with a 1.5-km horizontal resolution (see Supplementary Information for more details on the model). The values of w derived from the observations are comparable to, although smaller than, what is obtained in the model (Extended Data Fig. 5a). Overall, these results suggest that the vertical pathway provided by deep ocean fronts in spring and summer is likely a generic mechanism throughout the ocean. This deep-reaching vertical pathway has important consequences for the vertical transfer of heat between the ocean interior and the surface mixed-layer, as explored in the next section.

Vertical heat transport at deep submesoscale fronts

Oceanic heat transport is estimated from temperature and vertical velocity anomalies, from the surface down to 400 m depth (Fig. 5a-b). A vertical section of this transport (see Methods for the calculation) is shown in Fig. 5c, where a positive (negative)

value indicates an upward (downward) heat transport. Positive values result from frontogenesis processes, i.e the production of fronts, as illustrated in Fig. 3b. On the other hand, negative values arise from frontolysis processes, i.e. the destruction of fronts (see Supplementary Information). Heat transport is strongly enhanced at the location of submesoscale fronts generated by the background strain field and has a local amplitude that reaches 2000 W/m 2 below the surface mixed-layer, and extending down to at least 400 m, consistent with the high-resolution regional model (Extended Data Figs. 5b and6) and, surprisingly, of the same order of magnitude as instantaneous air-sea heat fluxes [START_REF] Thompson | Open-ocean submesoscale motions: A full seasonal cycle of mixed layer instabilities from gliders[END_REF] .

We now explore the contribution of fine oceanic scales to the domain-averaged vertical heat transport. Averaged vertical heat transport within the eddy-rich area of the ACC sampled by the seal (km 1000 to 2000, in red in Fig. 1) is directed upward, i.e. from the ocean interior back to the surface (Fig. 5f). This direction is strikingly opposite to the one induced by the diffusive processes used in the classical paradigm. However, this result can be understood in terms of the direct cascade

3/23
of potential energy, a well-known property of mesoscale eddy turbulence [START_REF] Salmon | Baroclinic instability and geostrophic turbulence[END_REF] . The direct cascade implies that frontogenesis statistically dominates frontolysis, and thus that the net vertical heat transport associated to submesoscale fronts is positive (see Supplementary Information), which is also striking in the model outputs (Extended Data Fig. 6). Furthermore, the heat transport magnitude reaches an averaged value of ⇠100 W/m 2 at 200 m (Fig. 5f). Remarkably, this value is an order of magnitude larger than that associated with mesoscale eddies alone [START_REF] Su | Ocean submesoscales as a key component of the global heat budget[END_REF] , which is, along with the mean flow, traditionally thought to be the main contributor to vertical heat transport 2 . It is also of comparable magnitude to domain-averaged air-sea heat fluxes [START_REF] Large | The global climatology of an interannually varying air-sea flux data set[END_REF] , illustrating a potential coherent pathway linking the ocean interior to the atmosphere. Note that the percentage of seal dives decreases with depth (gray line in Fig. 5f), from 100 to 40% between 200 and 400 m. This is likely to account for a part of the RMS decrease (stagnation) of temperature (velocity) anomalies below 200 m, as one would have expected 5 (Fig. 5d,e). As a consequence, the already unexpected strong heat transport (Fig. 5c,f) is likely to be underestimated because of it, in addition to the conservative estimate of the vertical currents discussed in the previous section. Indeed, the domain-integrated heat transport from the high-resolution numerical simulation yields a similar -although stronger -positive vertical heat transport (Extended Data Fig. 7) with, for instance, a value of ⇠140 W/m 2 at 200 m. As such, these results strongly contrast with the classical paradigm based on diffusive heat transport, as they emphasize the existence of intense and upward heat transport in the ocean interior well below the ocean surface mixed-layer that are preferentially localized in strain dominated areas (Extended Data Fig. 6).

The data presented here provide the first observational evidence of large and anomalous, i.e. upward, heat transport associated with deep-reaching submesoscale fronts in an eddy-rich area of the world ocean. The observations, supported by a high-resolution regional model (see Supplementary Information), highlight the crucial role played by submesoscale frontal dynamics, in the ocean interior well below the ocean surface mixed-layer, for oceanic heat transport.

In summary, The deep-reaching frontal dynamics and its associated large positive vertical heat transport, observed here in numerous sharp fronts of the ACC, are driven by mesoscale eddies and are likely to occur widely in the ocean, such as in the eddy-rich Gulf Stream and Kuroshio Extension, all of which are key players in the climate system. It potentially plays an important role by, for instance, exacerbating restratification processes as warm (cold) waters become warmer (colder). Furthermore, a first order estimation indicates these deep ocean fronts lead to an additional increase of the sea surface temperature ranging between ⇠0.2 and ⇠0.6 C within a month (see Methods for the calculation and caveats), illustrating the potential impact of deep-reaching ocean fronts on air-sea fluxes and how they can substantially limit oceanic heat uptake from the atmosphere.

These observational results suggest revisiting current estimates of the Earth's heat budget and stress the need to account for small-scale physics, not only within, but also below the ocean surface mixed-layer, in the prediction of future climate states.

Inaccurate representation of these physics could considerably underestimate the amount of heat transferred from the ocean interior back to the surface and, as a consequence, potentially overestimate the amount of heat the ocean can absorb. Finally, these results may have a far greater scope as the evidence for intense vertical currents associated to the deep-reaching ocean fronts presented here also provide an efficient pathway for the transport of chemical and biological tracers, with potential major implications for biogeochemical systems.

Methods

Sea Surface Height and geostrophic currents

The Sea Surface Height (SSH) used to retrieve the mesoscale dynamics of the Kerguelen area (67-91 E, 48-58 S) is composed of the Mean Dynamic Topography CNES/CLS 2015 [START_REF] Rio | Improving the altimeter-derived surface currents using high-resolution sea surface temperature data: A feasability study based on model outputs[END_REF] and the Sea Level Anomaly maps produced by the DUACS processing chain merging the multi-altimeter along-track data. This specific regional dataset was produced in the context of the DUACS-HR project aiming at increasing the resolution of altimetry maps by applying recently developed methods such as dynamic interpolation [START_REF] Ubelmann | Dynamic interpolation of sea surface height and potential applications for future high-resolution altimetry mapping[END_REF] . The period from October 2014 to January 2015 benefits from an extremely favorable satellite coverage of four altimeters: AltiKa, Jason-2, Cryosat-2 and HY-2, allowing a very good spatial and temporal sampling. The temporal resolution is daily and, as estimated in Ballarotta et al. (personal communication), the spatial effective resolution of the gridded anomalies approaches 40 km in wavelength. Using the geostrophic approximation, geostrophic surface currents (u, v) are derived from SSH following 10 :

u = g f ∂ SSH ∂ x , v = g f ∂ SSH ∂ y , ( 1
)
where g is gravity and f the Coriolis parameter.

Instead of using a direct estimate of the strain rate (i.e., du/dx calculated directly from SSH, Fig. 3) that only describes the buoyancy gradient's growth rate, we use Finite Size Lyapunov Exponents (FSLE). Indeed, FSLE have the advantage of exploiting both the spatial and temporal variability of the velocity field deduced from SSH, and consequently provide information about both the growth rate and orientation of elongated buoyancy gradients [START_REF] Waugh | Spatial variations of stirring in the surface ocean: A case study of the tasman sea[END_REF][START_REF] Ovidio | Comparison between eulerian diagnostics and finite-size lyapunov exponents computed from altimetry in the algerian basin[END_REF] . FSLE is a Lagrangian diagnostic that measures the separation of close initial particles embedded in a given flow field. The separation's growth rate is defined as An along-track time series of stretching (i.e. positive) FSLE was extracted along the seal's path in order to capture the stretching direction and intensity of the strain field. A subsequent 5-km window moving average was applied to remove the strong intermittency present in the raw data (Fig. 4 and Extended Data Fig. 4a). Stretching FSLEs are then compared to lateral gradients of buoyancy sampled by the seal.

l (d 0 , d f ) = 1 t log( d f d 0 ), (2)

Southern Elephant Seal measurements

A newly released in situ dataset collected in the austral summer by a female southern elephant seal during her three-month post-breeding trip (20 October 2014 to 16 January 2015) east of the Kerguelen Islands (Indian sector of the Southern Ocean, see Fig. 1) is analyzed. Compared to previous tags mounted on elephant seals, this one records every single dive realised by the animal (>80 dives/day) at high-resolution as opposed to 2 to 5 dives/day for previous tags. The seal is localized through the Argos satellite system and is equipped with sensors recording conductivity, temperature and pressure (CTD-Satellite Relay

Data Logger) at a continuous frequency of 0.5 Hz. Only the ascending phase of a dive is used because it is more uniform in speed and direction compared to the descent when the seal dives sinuously to forage [START_REF] Siegelman | Correction and accuracy of high-and low-resolution ctd data from animal-borne instruments[END_REF] . The dataset is comprised of 6333 dives, which corresponds to a cumulative length of 5270 km with a median spacing between two dives of 700 m (Extended Data Fig. 8a). Dives, that can be as deep as 500 to 1000 m, last less than 25 minutes and are separated by intervals of a few minutes where the seal breathes but does not transit. More than 80% of the dives reach a depth of at least 200 m, 50% reach 300 m and 35% 400 m or deeper.

To ensure a better accuracy of the temperature and conductivity-derived salinity data, two additional steps are applied: first a thermal cell effect correction is applied to the temperature and conductivity fields. A density inversion removal algorithm, which seeks the minimum adjustment of the salinity profile to achieve neutral stability is subsequently applied to the salinity field. The accuracy of the final post-processed data is of ±0.02 C for temperature and ±0.03 g/kg for salinity. More details on the post-processing method and final data accuracy can be found in Siegelman et al. (2019) [START_REF] Siegelman | Correction and accuracy of high-and low-resolution ctd data from animal-borne instruments[END_REF] . Potential density is calculated from corrected conservative temperature and absolute salinity with the TEOS-10 equation [START_REF]The international thermodynamic equation of seawater-2010: calculation and use of thermodynamic properties[END_REF] .

The

Buoyancy

Along-track time series of buoyancy, b = g(1 r/r 0 ), where g is gravity, r is potential density, and r 0 = 1025 kg m 3 is a reference density reveals variability covering both meso-and submesoscales. For the analysis, in particular the calculation of lateral buoyancy gradients b s = ∂ b/∂ s with s the curvilinear abscissa, buoyancy was first linearly interpolated along the seal's path onto a regular grid of 100 m resolution, corresponding to the shortest along-track distance between two dives (Extended Data Fig. 6). A moving average with a 1 km window was then applied such that the final dataset has a horizontal resolution of 1 km and a vertical resolution of 1 m. Buoyancy anomalies are resolved by multiple vertical profiles, such that the structures are

5/23
not related to aliasing of the along-track data. Surface buoyancy anomalies in Fig. 2a,b is defined by the time series of the mean buoyancy from 15 to 50 m minus its mean value calculated from 15 to 50 m along the entire trajectory.

Mixed layer depth

The mixed layer depth (MLD) is defined as the level of a 0.03 kg m 3 density increase from 15 m depth.

Normalization of buoyancy gradients

Since the seal's trajectory is more often oblique to the stretching FSLE it crosses rather than perpendicular to them (Extended Data Fig. 8b) and buoyancy fronts are assumed to be aligned with stretching FSLE, buoyancy gradients sampled by the seal need to be corrected in order to account for the seal's orientation with respect to the FSLE it encounters. To do so, buoyancy gradients are divided by the sine of the angle q between the seal's trajectory and the FSLE direction. More precisely, q is the angle between the FSLE eigenvector and the seal's path. To focus on the regions prone to the formation of submesoscales, only b s associated with large FSLE (>0.15 day 1 ) have been normalized and are then referred to as b x . However, a sensitivity analysis (ranging from 0.1 to 0.3 day 1 , not shown here) to this threshold value leads to similar results.

Link between strain, frontogenesis and vertical heat transport (Fig. 3)

The ambient strain field in Fig. 3a acts to elongate the tracer patch in the stretching direction (y-direction) and to compress it in the x-direction since the area of the tracer patch is conserved to leading order [START_REF] Vallis | Atmospheric and oceanic fluid dynamics[END_REF] . This increases the tracer gradient in the x-direction. The orientation and time scale of the strain field can be captured by the Finite Size Lyapunov Exponents (FSLE, see main text and Methods). Red (blue) FSLE lines identify the horizontal stretching (compressing) direction. This mechanism is particularly relevant for the formation of submesoscale fronts, that are thus aligned with the stretching FSLE (red line).Fig. 3b shows a schematic illustration of the frontogenesis process. The As a consequence, the vertical heat transport associated to frontogenesis is directed upward because of the positive correlation between temperature and vertical velocity anomalies.

Richardson number

The non-dimensional Richardson number Ri, estimated from the seal measurements, is defined as Ri ⌘ f

= f 2 N 2 ⇣ ∂ b/∂ z ∂ b/∂ x ⌘ 2 = f 2 N 2 ⇣ ∂ x ∂ z ⌘ 2 .
Thus, for Ri >> 1, which corresponds to the quasi-geostrophic regime [START_REF] Molemaker | Baroclinic instability and loss of balance[END_REF][START_REF] Thomas | Submesoscale processes and dynamics[END_REF] , the steepness of the isopycnal slope, ∂ b/∂ x ∂ b/∂ z , is small.

For Ri close to one (i.e. Ri  4), which corresponds to the ageostrophic regime [START_REF] Molemaker | Baroclinic instability and loss of balance[END_REF][START_REF] Thomas | Submesoscale processes and dynamics[END_REF] , the steepness of the isopycnal slope,

∂ b/∂ x ∂ b/∂ z ,
is large. Seal observations suggest a strong ageostrophic regime since instances of strong b x , |b x | 2.5 ⇥ 10 7 s 2 , coincide with Ri < 2 (Extended Data Fig. 2b). Small Ri-values (Ri < 2) indicate an ageostrophic regime where vertical currents are large [START_REF] Molemaker | Baroclinic instability and loss of balance[END_REF][START_REF] Thomas | Submesoscale processes and dynamics[END_REF] . Instances of strong b x , |b x | 2.5 ⇥ 10 7 s 2 , coincide with Ri < 2 (Fig. 2 and Extended Data Fig. 2b).

Vertical velocities

Classical methods to diagnose vertical velocities are based on the Q-vector version of the omega equation [START_REF] Hoskins | A new look at the w-equation[END_REF][START_REF] Giordani | Advanced insights into sources of vertical velocity in the ocean[END_REF] . In this study, we use the quasi-geostrophic (QG) version of the omega equation (see the Supplementary Information section for a discussion about this choice). Buoyancy fronts are assumed to be elongated along a stretching direction (red FSLE in Fig. 3a) such that the along-front gradient of buoyancy in the y-direction is negligible with respect to the cross-front one (x-direction). Thus, we consider the 2-D version (x,z) of the QG omega equation. This equation assumes that the trajectory is normal to the front, i.e.

that the front is embedded in a pure strain field, which is achieved through the buoyancy gradients normalization mentioned above. This equation reads :

N 2 w xx + f 2 w zz = 2 (u x b x ) x , (3)
where subscripts indicate derivatives. u x is estimated from stretching FSLE derived from satellite altimetry and b x , N 2 and f 2 from the seal's measurements.

Equation ( 3) is solved using the flexible framework for spectrally solving differential equations provided by dedalus [START_REF] Burns | Flexible framework for spectrally solving differential equations[END_REF] .

6/23

Note that when the Richardson number is close to one, equation (3) underestimates vertical velocities. When Ri ⇠ 2, the magnitude of this underestimation is of ⇠ 0.7, implying that the vertical velocities w diagnosed at sharp fronts should be closer to 1.4 ⇥ w (see Supplementary Information for more details).

Note also that this study focuses on strain dominated regions. However, in other regions, for instance, within mesoscale eddies, the strain rate is weaker and its impact on the formation of horizontal gradients is counterbalanced by the impact of the relative vorticity, which leads to the formation of weaker gradients of buoyancy. In such regions, even though the 2-D version of the omega equation is no longer appropriate and a 3-D version needs to be used, the resulting vertical velocities are 7.5 times weaker than the ones associated to submesoscale fronts in strain dominated regions 5 , like the ones considered in the present study.

Vertical heat transport

The vertical heat transport is defined as rC p w 0 T 0 , where C p = 3985 J kg 1 K 1 is the specific heat capacity of seawater, r = 1025 kg m 3 is the density of the fluid, w 0 and T 0 the vertical velocity and temperature anomalies, respectively.

Impact of deep-reaching ocean fronts on the mixed-layer temperature

Impact of deep-reaching ocean fronts on the mixed-layer temperature, and therefore the sea surface temperature, is derived from a thermodynamic equation that captures the evolution of the mixed-layer temperature T . This approach has already been used and validated at leading order [START_REF] Su | Ocean submesoscales as a key component of the global heat budget[END_REF] . The equation, integrated over the mixed-layer depth, reads:

C dT dt = S l T, (4)
where C = rC p H is the total heat capacity of the mixed-layer. H is the mixed-layer depth, S is the heat transport at the base of the mixed-layer due to deep-reaching fronts, and -l T is the negative feedback from air-sea heat fluxes. From equation ( 4), positive S causes a higher temperature T , and hence a larger upward air-sea heat exchange l T .

A sensitivity analysis is carried out to assess the range of the mean mixed-layer temperature change due to submesoscales over a period of one month, which corresponds to the time spent by the seal in the turbulent region and the time span of the high-resolution model. We set S to vary between 50 and 100 W/m 2 , as obtained from both the observational (Fig. 5f) and model (Extended Data Fig. 7) results, l to vary between 15 and 25 W m 2 K 1 , according to Vallis (2012) [START_REF] Vallis | Climate and the Oceans[END_REF] , and the mixed-layer H to vary between 100 and 200 m, implying that C varies between ⇠ 4 ⇥ 10 8 and ⇠ 8 ⇥ 10 8 J m 2 K 1 . As a consequence, equation (4) integrated over one month indicates that the mean mixed-layer temperature increase ranges between ⇠ 0.2 and

⇠ 0.6 C.
Note this is a first order estimation that includes several caveats. In particular, this estimation does not take into account any subsequent atmospheric feedback on the ocean that may arise in response to these surface heat fluxes (l T ). For instance, such feedback may include interactions between the windstress and SST anomalies at submesoscales [START_REF] Foussard | Response of surface wind divergence to mesoscale sst anomalies under different wind conditions[END_REF] , which may limit the mean temperature increase. FSLE are greatly enhanced in the strongly turbulent region (black rectangle and in red in Fig. 1) compared to the rest of the domain.

Figure captions

Extended Data Fig. 4 Finite size Lyapunov exponents and horizontal gradient of buoyancy, vertical velocities and vertical heat transport at 300 m. Times series of a) Horizontal gradients of buoyancy at 300 m sampled by the seal (in black)

and FSLE derived from satellite altimetry along the seal's track (in blue). b) Vertical velocities at 300 m derived from the seal and satellite data by solving the omega equation (see main text and Methods). c) Vertical heat transport (see Methods). The areas described in the main text and in Fig. 1 are highlighted by the colored rectangles. VHT is enhanced at depth and follows the isotherms on the eddy edges, and its averaged value is directed upward (positive value), all of which is consistent with the observational results presented in the main text.

22/23

Extended Data Fig. 7. Domain averaged vertical heat transport from the high-resolution numerical simulation.

Monthly averaged vertical heat transport (<VHT>) as a function of depth over the entire domain from the high-resolution numerical simulation. VHT is directed upwards (positive values) and its magnitude is similar -although even higher -than what is derived from the observational data presented in the main text.

Link between the strain field, Finite Size Lyapunov Exponents and lateral gradients of buoyancy

As illustrated in Fig. 3a, buoyancy anomalies embedded in a strain field are elongated in one direction and compressed in the perpendicular direction. This mechanism can be understood in terms of particle dispersion, and thus in terms of Finite Size Lyapunov Exponents (FSLE); two particles, initially close and embedded in a strain field, become quickly separated in one particular direction (the stretching direction, in red in Figs 3a and4a). As such, the time scale of their separation can be characterized by λ -1 (see equation 2 in Methods), with large λ being indicative of an intense strain field. More precisely, FSLE characterize both the time scale (via the largest FLSE eigenvalue λ ) and the direction (via the FSLE eigenvector associated to λ ) of the elongation of buoyancy anomalies.

Thus, the time scale λ (shown in the colorbar of Fig. 4a for instance) can be used to characterize the production of buoyancy gradient expressed in the RHS of the omega equation (equation 3 in Methods). In particular, the intersection of intense compressing and stretching FSLE lines, respectively in blue and red in Figs 3a and4a, identify Lagrangian hyperbolic points, where particles (or tracers) are simultaneously being stretched along one direction and compressed along the other one. Areas surrounding hyperbolic points are particularly prone to frontogenesis.

Vertical velocities underestimation

At steeply titled isopycnals' fronts, the vertical velocity field is more accurately captured by the full omega equation, known as the Sawyer-Eliassen (SE) equation, than by its QG version 2, 3 . The 2-D (x,z) SE equation, assuming the front is embedded in a pure strain field and diffusive processes are negligible, reads 2 :

N 2 ∂ 2 ψ ∂ x 2 + F 2 ∂ 2 ψ ∂ z 2 + 2S 2 ∂ 2 ψ ∂ z∂ x = 2u x b x , (2)
with ψ a meridional streamfunction related to the ageostrophic circulation (i.e. defined as

v a = -∂ ψ ∂ z , w = ∂ ψ ∂ x , with v a the ageostrophic component of the meridional component of the horizontal velocity field. N 2 = N 2 (x, z), S 2 = -b x and F 2 = f ( f + ∂ v ∂ x )
, where ∂ v ∂ x is the relative vorticity associated with the front 3 .

Compared with the QG version (equation 3 in Methods) rewritten in the form of the ageostrophic stream function, equation

2) involves F 2 = f ( f + ∂ v ∂ x (
) instead of f 2 as well as N 2 (x, z) instead of N 2 (z). In addition, equation (2) involves the extra term 2S 2 ∂ 2 ψ ∂ z∂ x . A dimensional analysis, using L and D as respectively the horizontal and depth scales of the front, indicates that this extra term is of the order of Ri -1 , with Ri the Richardson number, compared to the other two terms on the LHS 2, 3 . Since Ri -1 is not small but close to one at locations of sharp fronts, this term should not be ignored. Nevertheless, a real solution for equation (2) exists if and only if F 2 N 2 -S 4 > 0. This condition can be written as f 2 N 2 (1+Ro-Ri -1 )>0, with Ro the Rossby number. As such, it places a strong constraint on Ri -1 , which is postive and has to be smaller than |1+Ro|. This leads to a limitation of the steepness of the front and also to a limitation of the importance of the S 2 term relatively to others. However, we do not have access to F 2 and in particular to ∂ v ∂ x , the along-front relative vorticity, that can be of order f . SSH observations have a spatial resolution too low to resolve ∂ v ∂ x at these scales. This is why we use the QG version (equation 3 in Methods). Nevertheless, many studies, starting with Hoskins (1982) 4 and more recently confirmed by Hakim and Keyser (2001) 2 , emphasize that the SE and QG solutions are qualitatively similar when the condition F 2 N 2 -S 4 > 0 is met. In particular, the SE and QG ageostrophic circulation, and therefore the w-field, have the same sign. What differs is the amplitude of the w-field, which is larger in the SE solution, as well as the shape of the ageostrophic circulation, which is more tilted and oriented parallel to the isopycnals in the SE solution. Thus, the QG omega solution gives the correct sign for w but the magnitude may be underestimated at low Ri.

This underestimation can however be quantified using the analytical solutions for the QG and SE versions of the omega equation derived in Hakim and Keyser (2001) PV PE ] 1/2 = [ 1 1-Ri -1 ] 1/2 . When Ri = 2, which corresponds to the strongest fronts in the seal measurements (see main text), we get w SA ∼ 1.4 × w QG . In other words, the maximum w-values found in our study should be closer to 140 m/day than to 100 m/day.

2/5

Dominance of frontogenesis versus frontolysis understood via the direct cascade of potential energy

The classical frontogenesis process 5 emphasizes that the increase tendency of a buoyancy gradient embedded in a strong strain field is balanced by a vertical velocity field that tends to decrease this buoyancy gradient (black arrows in Fig. 3b). The reason is that a buoyancy gradient increase destroys the thermal wind balance and this balance is restored by the emergence of the vertical velocity field. These mechanisms lead to the omega equation (see Methods section "Vertical velocities"). As illustrated in Fig. 3b, the frontogenesis process leads to an upward vertical heat flux (red arrows in Fig. 3b) because of the positive correlation between temperature and vertical velocities anomalies.

On the other hand, the frontolysis process, or the destruction of front, occurs when the strain rate decreases and/or changes its direction. As a consequence, the buoyancy gradient relaxes, which is then balanced by a vertical velocity field (and therefore a vertical heat flux) in the opposite direction for the same reason as before (thermal wind balance). The frontolysis process is also captured by the omega equation and is associated with a downward vertical heat flux because of the negative correlation between temperature and vertical velocities anomalies (Fig. 5c). Frontogenesis (frontolysis) processes explain occurrences of positive (negative) vertical heat fluxes displayed in Fig. 5c and Extended Data Fig. 5b. However, a well-known property of mesoscale eddy turbulence is the direct cascade of potential energy driven by the background strain field. The direct cascade refers to the creation of buoyancy anomalies at smaller and smaller scales and thus to the continuous production of submesoscale fronts [START_REF] Ferrari | A frontal challenge for climate models[END_REF] . This is why frontogenesis processes statistically dominate frontolysis ones, as can be seen in Fig. S2c where positive buoyancy fluxes are more numerous and of greater magnitude than negative ones. As a consequence, we propose that the dominance of positive vertical heat fluxes at submesoscale as demonstrated in the main text is consistent with the strong background eddy field (Fig. 5f).

Numerical simulation comparison

The vertical velocity field and vertical heat transport (VHT) calculated from the observations are compared to the same fields obtained from a high-resolution numerical simulation in the Kerguelen area in late spring-early summer (November 2011) performed with the Massachusetts Institute of Technology general circulation model (MITgcm). This high-resolution simulation has a horizontal resolution of 1/48 • , 90 vertical levels and 10 minute outputs and is described in many papers (see for example Su et al. (2018) 7 ). The model domain is 45 • S-55 • S/ 85 • E-95 • E (or ∼ 1000 × 1000 km). This domain is embedded in the MITgcm LLC 4320 global numerical simulation performed at the same resolution but with hourly outputs 7 . Boundary conditions and forcings are supplied by the global simulation.

Vertical sections of daily averaged vertical velocities and vertical heat transport obtained from this simulation are presented in Extended Data Fig. 5. At depths of 100-600 m, intense vertical velocities with a width of 5-10 km are present below the mixed layer and reach 300 m/day. Similar features are observed for vertical heat transport with values locally reaching 4000 W/m 2 . Two-dimensional slices (x,y) of vertical heat transport averaged over 10 days are shown in Extended Data Fig. 6 at 50 and 200 m. This Fig. clearly highlights that the organization of VHT is driven by the background mesoscale strain field, intensified on the eddy edges. This Fig. also shows the dominance of positive VHT in the domain as well as its intensification at depth. Finally, Extended Data Fig. 7 presents the domain-averaged VHT over one month, which reveals similar -although largervalues than the ones derived from the observational data. Overall, the similar shape, distribution within the water column, and magnitude of both quantities in the model and in the observations strengthen the results presented in this study, and further highlight the impact of deep reaching ocean fronts on oceanic heat transport.

  where d 0 (d f ) is the initial (final) separation distance and t the first time at which a separation d f is reached. FSLE has the dimension of time1 . FSLE is particularly suited to diagnose the properties of a strain field. Positive (negative) FSLE indicate that patches of particles are being stretched (compressed) in a given direction determined by the background strain field. As such, large positive (negative) FSLE values indicate regions of strong stretching (compression) (see Fig.3aand Supplementary Information for more details). Here, positive (negative) FSLE are computed backward (forward) in time[START_REF] Ovidio | Mixing structures in the mediterranean sea from finite-size lyapunov exponents[END_REF] .88 daily maps of altimetry-based FSLEs were computed for the Kerguelen area following d'Ovidio et al. (2004) [START_REF] Ovidio | Mixing structures in the mediterranean sea from finite-size lyapunov exponents[END_REF] and using the geostrophic velocities derived from SSH. Parameters were chosen close to d'Ovidio et al. (2004) 22 with d 0 = 0.04 and d f = 1 , i.e. a final separation distances of about 110 km.

  animal in this study was handled in accordance with the Institut polaire francais Paul-Emile Victor (IPEV) ethical and Polar Environment Committees guidelines as part of the SNO-MEMO and IPEV program 109 (PI. H. Weimerskirch). The experimental protocols were approved by the Ethics Committee of IPEV and Polar Environment Committees.

  Fig. corresponds to a 3-D slice through a submesoscale front generated by a background mesoscale strain field (curved black arrows) like the one in a). Thin black lines are isotherms. Straight black arrows correspond to the vertical velocities that develops in response to the front intensification. Warmer (or equivalently lighter) fluid (light red) is on the right of the Fig. and colder (or equivalently heavier) fluid (light blue) is on the left.

2 N 2

 22 /b 2 x , with N 2 = ∂ b/∂ z the Brunt-Väissälä frequency and b x = ∂ b/∂ x the along-track lateral gradient of buoyancy. Ri characterizes the dynamical regime and can be interpreted as the steepness of the isopycnal slopes relatively to N/ f , since Ri

Fig. 1 / 23 Fig. 2

 1232 Fig.1 Study region. Three-month seal's trajectory (black line) superimposed with a snapshot of chlorophyll a (colorscale) and SSH (white contours). Three regions are highlighted: the highly turbulent area (red), the weakly turbulent area (grey) and the southern eddy edge (orange). Cyclones (C) and anti-cyclones (A) are identified from elliptic SSH contours. Hyperbolic SSH contours located in-between eddies identify the strain field (see main text and Fig.3a). Red arrows indicate the direction of the seal. Bathymetry contours of 0.5, 2 and 3 km (thin black lines) from NOAA ETOPO5 outline the Kerguelen plateau. The inset shows the Kerguelen region (red polygon).

Fig. 3

 3 Fig.3 Strain field, frontogenesis and vertical heat transport. a) Horizontal slice (x,y) of tracer patch (light blue) in a strain field (black arrows). The strain elongates (compresses) the tracer in the y(x)-direction[START_REF] Vallis | Atmospheric and oceanic fluid dynamics[END_REF] . Red (blue) Finite Size Lyapunov Exponents (FSLE, proxy of strain) identify horizontal stretching (compressing) directions. Fronts are aligned with stretching FSLE. b) 3-D slice of strain-induced submesoscale front. Thin black lines are isotherms. Vertical velocities (w, straight black arrows) develop in response to the front intensification. Since temperature and w anomalies are positively correlated, frontogenesis-induced vertical heat transport is upward (see Methods for a full description).

Fig. 4 FSLE

 4 Fig. 4 FSLE and horizontal buoyancy gradients (|b x |). a) Snapshot of FSLE (colorscale) and SSH contours (dashed gray) for 11/13/2014 in the strongly turbulent region (red in Fig. 1), with surface |b x | from 11/02/2014 to 11/24/2014. b) Scatterplot of along-seal-track stretching FSLE and |b x | at 300 m, over the entire trajectory, after a 30km-moving average and averaging the y-axis over 100 bins (black dots), unaveraged data in gray. Relationship obtained between FSLE and |b x | in red. c) Along-seal-track time series of stretching FSLE (blue) and |b x | at 300 m (black) after a 5km-moving average, from km 1600 to 2100.

Fig. 5 Extended Data Fig. 1 Extended Data Fig. 3

 513 Fig. 5 Temperature anomaly, vertical velocity and vertical heat transport (VHT). Vertical section of a) seal-sampled temperature anomaly, b) vertical velocity derived from seal and satellite data with equation 3, c) VHT, as defined in the Methods. In a)-c), the mixed layer depth is in bold black. RMS(z) from km 1000 to 2000 of d) temperature anomaly, e) vertical velocity. f) averaged VHT (<VHT>). <VHT> reaches ⇠ +110 W/m 2 at 200 m. The percentage of profiles (in gray) goes from 100 % at 200 m to 40% at 400 m.

10 / 23 Extended Data Fig. 5 Extended Data Fig. 6 Extended Data Fig. 7 Figure 1 . 23 Figure 2 .Figure 3 .Figure 4 . 23 Figure 5 . 23 Extended Data Fig. 1 .. 2 .

 10235671232342352312 Figure1. Study region. Three-month seal's trajectory (black line) superimposed with a snapshot of chlorophyll a (colorscale) and SSH (white contours). Three regions are highlighted: the highly turbulent area (red), the weakly turbulent area (grey) and the southern eddy edge (orange). Cyclones (C) and anti-cyclones (A) are identified from elliptic SSH contours. Hyperbolic SSH contours located in-between eddies identify the strain field (see main text and Fig.3a). Red arrows indicate the direction of the seal. Bathymetry contours of 0.5, 2 and 3 km (thin black lines) from NOAA ETOPO5 outline the Kerguelen plateau. The inset shows the Kerguelen region (red polygon).

  Figure S5

  2 . These authors show that w SA ∼ w QG [PV QG PV PE ]1/2 with w SA and w QG respectively the solution of the QG and SA equations and PV QG and PV PE respectively the QG potential vorticity and the Ertel potential vorticity. Using the approximations PV QG ∼ f N 2 and PV PE ∼ f N 2 [1 -Ri -1 ], as detailed by Thomas et al. 2008 3 , leads to [

	PV QG

Acknowledgements

Thanks to Kelvin Richards for his insightful comments, Francesco d'Ovidio for providing the code to compute FSLE and two anonymous reviewers whose comments improved the manuscript. This work was supported by the CNES-TOSCA project Elephant seals as Oceanographic Samplers of submesoscale features led by C. Guinet with support of the French Polar Institute (Progams 109 and 1201). This research was carried out, in part, at the Jet Propulsion Laboratory, California Institute of Technology, under a contract with the National Aeronautics and Space Administration (NASA). High-end computing resources for the numerical simulation were provided by the NASA Advanced Supercomputing Division at the Ames Research Center. This work was partly funded by the CNES (OSTST-OSIW) and the Laboratoire d'Excellence LabexMER (ANR-10-LABX-19). L.S. is a NASA-JVSRP affiliate and is supported by a joint CNES-Région Bretagne doctoral grant. P.K. is supported by the NASA-CNES SWOT mission and a NASA Senior NPP Fellowship. A.F.T. is supported by the David and Lucille Packard Foundation and NASA grant NNX16AG42G. M.F. is supported by NASA grant NNX15AG42G.

Data availability

The marine mammal data were collected and made freely available by the International MEOP Consortium and the national programs that contribute to it, and is available at http://www.meop.net/database/meop-databases/ meop-sms-database-submesosc.html. The Ssalto/Duacs altimeter products were produced and distributed by the Copernicus Marine and Environment Monitoring Service with support from CNES, and is available at http://marine. copernicus.eu/services-portfolio/access-to-products/.

Author contributions statement L.S. and P.K. conceived the experiments, analysed the results and wrote the manuscript. D.M. and H.S.T. ran the numerical simulations. H.S.T helped with analyzing the regional simulation. L.S., P.K., P.R., A.F.T., H.S.T. and M.F. reviewed the manuscript.

Competing interests

The authors declare no competing interests. Extended Data Fig. 3. Map of finite size Lyapunov exponents Map of FSLE over the entire domain on 13 November 2014. FSLE are greatly enhanced in the strongly turbulent region (black rectangle and in red in Fig. 1) compared to the rest of the domain.

19/23

Consistency between seal and satellite observations at meso-and submesoscale

The contribution of the vertical stratification to the SSH can be estimated by integrating the hydrostatic equation dp/dz = -ρg. The SSH (η) referenced from z = 0 is then given by 1

where p is the pressure, ρ the potential density, ρ 0 the reference potential density (1025 kg m -3 ), ρ the potential density anomaly (ρ = ρ 0 +ρ ), g the gravity of Earth, p a the atmospheric surface pressure, H the depth of the ocean and p b = p b -ρ 0 gH the bottom pressure anomaly with p b the ocean bottom pressure. The term η 0 ρ /ρ 0 dz has been neglected because η H. The three terms on the right hand side of equation ( 1) represent respectively the contributions from the bottom pressure, the atmospheric pressure loading, and the steric height. The steric height (ζ ) is computed with the seal dataset from the surface down to H = 400 m for dives reaching at least 400 m (3100 dives) as ζ 400 = 0 -400 ρ /ρ 0 dz. ζ 400 is then compared to the SSH observed by altimetry.

SSH obtained from satellite data and mixed layer depth (MLD) and buoyancy anomalies, sampled by the seal, have a larger magnitude in the strongly turbulent region (Fig. 2) than in the weakly turbulent one (Extended Data Fig. 1i). Anticyclonic structures, or positive SSH anomalies, are associated with an increase in MLD and bowl-shaped positive buoyancy anomalies from the surface down to 400 m depth (Fig. 2 and Extended Data Fig. 1). Shallow MLD and negative buoyancy anomalies are associated with cyclonic structures, or negative SSH anomalies, confirming the fact that buoyancy anomalies compensate SSH anomalies. To further quantify this compensation over the water column, we have compared SSH with the depth-integrated buoyancy also called steric height. Since seal observations mostly sample the upper ocean, the steric height contribution has been estimated from surface down to 400 m (which concerns at least 35% of the dives) as ζ 400 = 0 -400 bdz with b the buoyancy (see Methods). As emphasized by the time-series in Fig. S1a (SSH in black and ζ 400 in blue), SSH and steric height fluctuations correlate well over the entire trajectory, but only for scales smaller than 100 km. SSH reveals a large-scale signal (∼1000 km) not present for ζ 400 . This indicates that SSH at larger scales is likely compensated by the buoyancy field below 400 m and/or by bottom pressure, The contribution of ζ 400 to SSH is no more than 30% as revealed by the linear regression coefficient of 3.58 linking the two-time series (not shown).

To remove this large-scale contribution and focus only on the meso/sub-mesoscale band, time-series of the horizontal gradients of SSH and ζ 400 are computed (a gradient operator overemphasizes contribution of smaller scales). This is the equivalent of comparing horizontal currents from SSH using the geostrophic approximation, with those explained by the buoyancy field contribution in the water column using the thermal wind balance 1 . The two time-series are now closer in terms of amplitude differences (Fig. S1b). ζ 400 gradients statistically explain almost 75% of the SSH gradients as revealed by the linear regression coefficient of 1.35 linking the two time-series (not shown). This further emphasizes the consistency between the two independent datasets at meso/sub-mesoscale despite their different spatial resolution. The main differences concern the extrema of SSH gradients (mostly located on the eddy edges) that overcome those of ζ 400 gradients by a factor of two to three. This suggests the steric height is not integrated deep enough to capture SSH in the upper mesoscale band and therefore that dynamics of this band affects at least the first 400 m below the surface.


Y Wu

Z Sun

F Brisset

T Baudin

A L Helbert

D Retraint

In-situ EBSD investigation of thermal stability of a 316L stainless steel nanocrystallized by surface mechanical attrition treatment

Keywords: Nanocrystalline materials, SMAT, In-situ heating, EBSD, Thermal analysis

In-situ Electron BackScatter Diffraction (EBSD) was used to study the thermal stability of a stainless steel nanocrystallized by Surface Mechanical Attrition Treatment (SMAT). A grain size gradient was generated after SMAT from the surface to the interior of the specimen.

Observations at different temperatures were performed to study the thermal stability of the gradient microstructure. During the investigation, oxidation was detected which greatly affected the indexation quality especially for high temperatures. After an ionic polishing, the grains could be properly revealed. No obvious microstructure changes were highlighted up to 720°C, which indicated a good thermal stability of the nanocrystalline grains.

Introduction

Surface Mechanical Attrition Treatment (SMAT) is one promising surface modification technique. It can transform the surface layer of coarse grains to ultrafine grains by means of severe plastic deformation [1][2][3]. The nanocrystallization process is determined by high strain rate multi-directional impacts between shot and the part. During SMAT, a large quantity of crystallographic defects such as plastic slips can be generated. The multiplication of plastic slips can progressively lead to ultrafine grains in the near surface region if the treatment intensity is high enough. However, in the bulk region, the material state is not modified, and the microstructure remains unchanged. Between the ultrafine grained layer and the bulk, a transition region is present and characterized by a high strain hardening gradient. Some metallurgical parameters generated by SMAT including ultrafine grains can be beneficial to enhance mechanical properties such as wear [4] or fatigue [5][6][7][8].

For components exposed to high temperature in service, the thermal stability of gradient microstructure generated by SMAT is an important factor as it determines the contribution of SMAT for such conditions. During the last decade, the thermal stability has been studied for various SMATed materials [9][10][11]. The results show that the thermal stability is dependent on materials. For example, nanometer grains of alloys with precipitates are more stable when subjected to temperature [9], whereas nanometer grains of pure metals are less stable [10,12].

In these studies, different techniques such as Transmission Electron Microscopy (TEM), Scanning Electron Microscopy (SEM) were used, and all the investigations were carried out ex-situ. The heating process must thus be interrupted in order to observe the microstructure.

Compared to in-situ observations, real time microstructure evolution cannot be monitored with ex-situ techniques.

In this work, in-situ Electron BackScatter Diffraction (EBSD) was used to follow eventual microstructure changes of a SMATed steel during thermal exposure. The temperature was increased up to 720°C to study the thermal stability, especially of the ultrafine grains in the surface layer.

Material and experiments

The investigated material is a 316L stainless steel. Its chemical composition is as follows (wt.%): 17.37Cr, 2.80Mo, 14.52Ni, 1.70Mn, 0.26Si, 0.17P, 0.08Cu, 0.07V and Fe balance. For SMAT, shot of 3 mm diameter in 100Cr6 steel was boosted by an ultrasonic generator working at 20 kHz. The specimen was first treated during 15 minutes with a generator power of 30%, then treated during 5 minutes with a power of 50%.

In-situ EBSD technique was used to characterize the microstructure under a temperature increase up to 720°C. The observations were performed using a FEG-SEM SUPRA 55-VP equipped with a platinum based heating system. The temperature was controlled in order to follow a programmed increase path. It was increased with a rate of about 40°C/minute. For different values, the temperature was kept constant so as to perform EBSD measurements.

EBSD scanning was performed in an area of 200 µm×50 µm with a step size of 0.5 µm. About 4 minutes were needed to achieve a cartography. For the specimen preparation, a sample was successively SMATed, cut along cross-section, molded, mechanically ground, polished to a mirror-like finish and then polished with an OPS solution.

Results and discussion

EBSD observations on the cross-sectional surface reveal a microstructure gradient in the SMAT affected region, similar to that obtained in a previous study [3]. A grain refinement down to the nanometer scale occurred in the top surface region. The thickness of the nanostructured layer is about 10 µm, followed by a mixed grain size region, where both ultrafine grains and relatively large grains are present. A typical microstructure obtained at 400°C is illustrated in Fig. 1a. This microstructure can also represent the as-treated state of the material, since at low temperature the microstructure is rather stable. In addition, SMAT doesn't induce any austenitic phase transformation into martensite in this 316L alloy, as found in previous studies [11,13].

Fig. 1 illustrates the obtained maps for temperatures from 400°C to 720°C. It can be seen that with an increase of temperature, some areas become dark (Figs. 1b to 1e), which means that these areas are no longer well indexed due to the degradation of the Kikuchi diagram quality. This phenomenon is more significant in the areas with small grains, especially in the nanostructured layer. At 720°C, the indexation is totally lost (Fig. 1f). Fig. 2 presents some characteristic areas after a thermal exposure at 720°C and cooling down to room temperature. As it can be seen, all areas give rise to Kikuchi diagrams with more or less acceptable quality, which indicates that over a certain temperature, electron diffraction is deteriorated. From 560°C to 650°C, a loss of Kikuchi pattern quality occurs at grain boundaries (Fig. 1) and especially at the SMATed surface. According to Fig. 2d, this loss of indexation quality is due to an oxidation effect, as indicated by the Energy Dispersive X-ray Spectrometry (EDS) analysis. This phenomenon occurs more easily on small grains. It seems strange that the sample is oxidized even in a secondary vacuum (around 10 -10 atm), i.e. the working condition of a conventional SEM. In fact, according to the Ellingham diagram [13], chromium oxide can be formed at 500°C even when the pressure of oxygen is as low as 10 -42 atm. However, the oxidation kinetics is quite low at 500°C and the temperature holding is not long enough to notice any significant changes in the indexation quality (see Fig. 1b). That is why a noticeable indexation quality degradation can be observed only beyond 560°C (Fig. 1c). Fig. 2a shows an oxidation gradient from the surface to the centre of the sample. In the transition area where the alloy was strongly plastically deformed, oxidation takes place and the indexation quality is degraded (Fig. 2c). Finally, in the area far from the surface, the material is subjected to thermal grooving while the oxidation level remains low: the indexation quality remains thus good (Fig. 2b). Note that the thermal grooving of grain boundaries is observed from 600°C. The observed heterogeneous oxidation can be related to the poor oxidation resistance of fine grained materials. As a matter of fact, in the region where the grain size is small, the fraction of grain boundaries is relatively high. It is well known that at grain boundaries, the atomic arrangement is less regular. Grain boundaries can thus provide short circuit diffusion path for oxygen to penetrate more easily into the material. Oxidation preferably begins at the SMATed surface where grain boundary density is high and extends to the bulk where the diffusion kinetics is relatively low. In the case of austenitic stainless steels, another mechanism can be involved in the oxidation process. When these steels are heated to 425 ̵ 815°C, chromium carbides can form at grain boundaries. As a result, the narrow areas adjacent to the grain boundaries are highly depleted in chromium and oxidation may easily occur there [START_REF] Liang | [END_REF]15]. It is worth indicating that the formation of chromium oxide is easier than that of iron oxide for these pressure and temperature conditions, according to the Ellingham diagram [13]. Contrary to what was initially expected, the microstructure evolution is negligible despite the presence of areas with high stored energy favorable for recrystallization and grain growth. As thermal grooving and oxidation can both hamper the grain boundary migration at the surface, it was decided to verify the microstructure under the oxidized surface. Ionic polishing was therefore carried out a posteriori on the sample heated at 720°C to remove the oxidized layer.

After ionic polishing, the indexation quality became much better according to the Kikuchi diffraction pattern shown in Fig. 3. The microstructure looks very similar to that obtained after SMAT, and no notable change in grain size (Fig. 3) and misorientation is observed according to the analyses of the obtained cartographies. This means that the microstructure is not obviously changed by the temperature. The ultrafine grain size remains constant at 0.7µm +/-0.2µm. This investigation is consistent with the results presented in the literature concerning the thermal stability of nanostructured layer generated by SMAT. For example, in the studies on 316L steels, grain size change seems to begin at about 600°C [16] or beyond 0.5Tm (Tm is the melting temperature) [17].

Conclusion

In this work, the feasibility of the in-situ EBSD was illustrated by investigating the thermal stability of a 316L stainless steel gradient microstructure generated by SMAT. The initial microstructure obtained after SMAT shows a gradient microstructure from the top surface to the interior of the specimen. The thermal stability investigation shows that the microstructure does not significantly evolve up to 720°C, and especially the small grains are maintained. An oxidation phenomenon takes place and becomes more pronounced when the temperature exceeds 560°C.

Figure 1 -

 1 Figure 1 -EBSD maps from (a) to (f) show microstructure at different temperatures. The appearance of dark areas signifies a progressive degradation of indexation quality, especially for the treated surface (indicated by an arrow).

Figure 2 -

 2 Figure 2 -Micrographs obtained after heating up to 720°C and cooling, EDS analyses and Kikuchi diffraction patterns: (a) global view of a large region (the arrow indicates the surface), (b) bulk region, (c) strongly work-hardened region, and (d) very near surface region.

Figure 3 -

 3 Figure 3 -EBSD observation performed after the removal of the oxidized layer (720°C) using ionic polishing. Ultrafine grain size distribution at the SMATed surface for both 400 and 720°C.


Afef Mabrouki

Patricia Le Nahenec-Martel

Abdelkader Kriaa

Ahmed Hedhli

Pierre-Yves Renard

Cyrille Sabot
email: cyrille.sabot@univ-rouen.fr

Dearomatization of a 3-Hydroxypyridine through an unexpected Oxidative Deformylation Process: an Entry to Aza-cyclohexadienones

Keywords: hydroxypyridine, oxidation, hypervalent iodine reagant, PIDA, azacyclohexadienone Scheme 2 Oxidative hydroxypyridine dearomatization

Phenols are well-known precursors of cyclohexadienones, which have found useful widespread applications in organic synthesis. In contrast, their hydroxypyridine counterparts have not been explored yet. An unprecedented oxidative dearomatization of 3-hydroxypyridine involving an unexpected deformylation step, is reported. The chemical reactivity of the subsequent unreported aza-cylohexadiene-type compound, was also explored.

The oxidative dearomatization of phenols offers a straightforward entry to a variety of 4-alkoxy-2,5-cyclohexadien-1-ones (also named p-quinol ethers) mostly from commercially or readily available starting materials. 1,[START_REF] Quideau | Phenol Dearomatization with Hypervalent Iodine Reagents[END_REF] This oxidation step notably increases the derivatization ability of phenols, as their cyclohexadienone derivatives can undergo subsequent nucleophilic 1,2 or 1,4 additions, 35 α-alkylations or -arylations, 6 DielsAlder reactions, 7,8 photochemical ring contractions, 9 or dienone-phenol rearrangements, 6,9 to cite a few (Scheme 1, X = CH).

Scheme 1 Oxidative dearomatization of phenols and hydroxypyridines.

Whilst this strategy has extensively been studied and exploited through the expedient synthesis of natural products, 1014 on the other hand the oxidative dearomatization of their hydroxypyridine counterparts has not been reported yet, despite promising synthesis applications (Scheme 1, X = N). 1519 Although the oxidative homocoupling of 3-hydroxypyridines was recently described, 20 the competitive oxidation of the pyridine nitrogen and hydroxyl group is a plausible issue 21 which may have prevented further dearomatization studies and reports on hydroxypyridine systems.

Herein, as a-proof-of-concept study, we report the unexpected formation of a 5-aza-2,5-cyclohexadien-1-one via an oxidative deformylation/1,3-alkoxy rearrangement process of commercially available 2,6-lutidine-α2,3-diol. The chemical reactivity of the corresponding unreported dearomatized product was investigated.

Hypervalent iodine reagents, such as phenyliodine(III) diacetate (PIDA) PhI(OAc)2, are reagents of choice to promote phenols oxidation under mild and environmentally friendly conditions. [22][23][24][25][26][27] Accordingly, the oxidation of commercially available 2,6-lutidine-α2,3-diol 1 was investigated under standard alcoholic conditions with the aim to provide the corresponding aza-p-quinol ether 3. To this end, the oxidation was carried out in the presence of 1 equivalent of PIDA in MeOH at room temperature for 15 minutes (Scheme 2). mixture using flash column chromatography provided compound 2a in 27% yield, whose structure was confirmed by spectroscopic and spectrometric data. Unexpectedly, the hydroxy-methyl group in 1 was replaced by a methoxy substitutent in 2a. Compound 2a being two oxidation states higher than the parent compound 1, 2.14 equivalents of PIDA were used in order to drive the reaction to completion and improve the overall efficiency of this transformation. As expected, under these optimized conditions, compound 2a was isolated in 60% yield (Scheme 2b).

Following these preliminary results, the reaction was conducted in other low-boiling-point alcoholic solvents, and obtained results are presented in Table 1. [START_REF]Aza-p-quinol ether (2a) Typical procedure To a stirred solution of 2,6-lutidine-α2,3-diol (100 mg, 0.719 mmol) in 1.5 ml of solvent (MeOH, EtOH or TFE), phenyliodine(III) diacetate (PIDA) PhI(OAc)2 (500 mg, 1.538 mmol, 2.14 equiv) in 1,5 ml of solvent (MeOH, EtOH or TFE) was added dropwise (using pipette) under air atmosphere at room temperature. After stirring for 15 min, the reaction which came yellow and was completed by TLC analysis, was concentrated under reduced pressure. The crude product was purified by flash-column chromatography on silica gel using cyclohexane and ethyl acetate (7:3) as eluent to give a yellow oil product[END_REF] The yields of dearomatized 3-hydroxypyridines decreased as the bulkiness of the solvent increased. In fact, a moderate yield of 43% of 2b was obtained in ethanol (entry 2), whereas only traces of product 2c could be identified in isopropanol (entry Besides, the use of electron-deficient alcohols such as 2,2,2trifluoroethanol (TFE) markedly decreased the efficiency of the transformation with the isolation of 11% of dearomatized compound 2d, while product 2b was obtained in 43% in nonfluorinated analog, ethanol (entry 2versus 4). Not surprisingly, no product was isolated when the reaction was conducted in 1,1,1,3,3,3-hexa-fluoroisopropanol (HFIP), a bulky and nonnucleophilic fluorinated alcohol. The reaction was also explored in the non-protic and chemically inert solvent, dichloromethane, in order to investigate whether the acetic acid released in the course of the reaction could trap the oxidized intermediates as alcoholic solvents did (Scheme 3).

Scheme 3 Oxidation of compound 1 in CH2Cl2

Under these conditions, only the deformylated mono-oxidized compound 4 was obtained, albeit in a low 14% yield. In contrast to alcoholic solvents, the formation of this compound would involve a mono-oxidation process.

No example of such unfused aza-p-quinol ether systems having being reported yet, their reactivity was investigated using the representative compound 2a. First, the direct reductive aromatization of 2a to provide an electron-rich pyridine, was carried out by the action of Pd/C mediated hydrogenation, which delivered compound 5 in a low 37% yield. This volatile material possessed a pleasant fragrance (Scheme 4).

Scheme 4 Direct aromatization of aza-p-quinol ether 2a

Then, compound 2a was subject to 1,2-nucleophilic addition conditions. In view of the operational simplicity associated with the use of non toxic allyl indium reagents, compound 2a was treated with indium (In) in the presence of allyl bromide to produce the corresponding allylic alcohol 6 in 46% yield as a mixture of inseparable diastereoisomers in a 1:1 ratio (Scheme 5).

Scheme 5 Allylindium-mediated 1,2-alkylation of 2a

It is worth mentioning that, in contrast to what was reported with p-quinol ether systems, 29 no further anionic vinylogous semipinacol type rearrangement of 6 into 7 was observed in the present case.

The conjugate addition of in situ generated organocuprate reagent to the α,β-unsaturated system of 2a was not observed (Scheme 6a). In fact, the reaction of methyl magnesium bromide in the presence of copper (I) bromide favored the formation of carbonyl addition product 8, which was isolated in 33% yield as a ~3:1 mixture of inseparable isomers.

Scheme 6 Investigation of conjugate additions with dimethyl malonate

Thus, the reaction was next investigated with the Michael donor dimethyl malonate, which enabled this time the formation of the 1,4-addition product 9, isolated in 70% yield as a single diastereoisomer (Scheme 6b).

Next, the reactivity of aza-p-quinol ether as a diene in a Diels-Alder reaction was investigated. To this end, compound 2a was reacted with freshly distilled cyclopentadiene in the presence of phosphoric acid as the catalyst (Scheme 7). 30 The corresponding cycloadduct 10 was formed in 23% isolated yield. This low yield may partly be ascribed to the formation of fused bicyclic pyridine 11, which can be formed cleanly in 77% yield from isolated compound 10, in the presence of methanesulfonic acid in methanol.

Scheme 7 Diels-Alder reaction of aza-p-quinol ether 2a

Finally, a plausible reaction mechanism for the oxidative deformylation process of 3-hydroxypyridine 1 can be proposed in Scheme 8. The oxidation of the aromatic hydroxyl function promoted by the first equivalent of hypervalent iodine reagent would lead to the formation of the phenoxonium intermediate 12, which would subsequently be trapped by methanol to furnish the transient intermediate 13.

Scheme 8 Proposed reaction mechanism

Then, compound 13 would undergo a retroaldolizationdriven aromatization process to provide the 2-methoxy-3hydroxy-pyridine 14. This reaction would occur with the concomitant release of formaldehyde. It is worth noting that the putative intermediate 14 has not been isolated, even in the presence of sub-stoichiometric amounts of PIDA, presumably because this electron-rich pyridine is more readily oxidized than the starting 3-hydroxypyridine 1. Then, the second equivalent PIDA would oxidize 14 to generate the aza-o-quinone monoketal 16 through the regioselective addition of methanol to the methoxy-stabilized carbocation at the α-position of the carbonyl of the intermediate 15. 31 This compound would be converted into enone 17 through a 1,3-methoxy rearrangement process, similarly to what was reported with o-quinone monoketals under Lewis acid catalysis. 32 Accordingly, the presence of the temporary hydroxy-methyl group at the position 2 of the pyridine ring system of 1 would play an important role in favoring the methanol addition to occur at this position due to the enhanced stability of this tertiary carbocation. To investigate this point, the 3-hydroxypyridine 18, unsubstituted at position 2, was oxidized under identical optimized PIDA conditions (Scheme 9). This reaction provided 2a in only 20% yield, instead of 60% when the position 2 was functionalized by the hydroxy-methyl substituent. This result shows that the presence of a removable hydroxy-methyl group, was required to ensure effective formation of aza-p-quinol ethers 2.

In conclusion, the first oxidative dearomatization of a hydroxypyridine system leading to the formation of aza-p-quinol ethers, is reported. This reaction proceeded under mild conditions, and involved a deformylation/1,3-migration process favoring the formation of aza-p-quinol ethers. Subsequent modifications of this particular cyclohexadienone through 1,2-, 1,4-addition, and Diels-Alder reaction were investigated.

Scheme 9

 9 Scheme 9 Oxidation of 2-unsubstituted-3-hydroxypyridine 18

Table 1

 1 Solvent-scope of the reaction

	Entry	Solvent	R	Cmpds	Yield (%)
	1	MeOH	Me	2a	60
	2	EtOH	Et	2b	43
	3	iPrOH	iPr	2c	traces
	4	TFE	CF2CH2	2d	11
	5	HFIP	(CF3)2CH	2e	-

a Reaction conditions: PIDA (2.14 equiv.), solvent, rt, 15 min.

Acknowledgment

We thank Albert Marcual (CNRS) for HRMS analyses.

Funding Information

This work was funded by the Centre National de la Recherche Scientifique (CNRS), INSA Rouen, Normandie Rouen University, Région Normandie and the Labex SynOrg (ANR-11-LABX-0029).

Supporting Information

YES (this text will be updated with links prior to publication)

Primary Data

NO (this text will be deleted prior to publication)