# Unraveling the Impacts of Submesoscale Thermal and Current Feedbacks on the Low-Level Winds and Oceanic Submesoscale Currents

Lionel Renault, Marcela Contreras, Patrick Marchesiello, Carlos Conejero,
Igor Uchoa, Jacob Wenegrat

► To cite this version:

Lionel Renault, Marcela Contreras, Patrick Marchesiello, Carlos Conejero, Igor Uchoa, et al. Unraveling the Impacts of Submesoscale Thermal and Current Feedbacks on the Low-Level Winds and Oceanic Submesoscale Currents. Journal of Physical Oceanography, 2024, 54 (12), pp.2463-2486. 10.1175/JPO-D-24-0097.1 . ird-04901096

HAL Id: ird-04901096

https://ird.hal.science/ird-04901096v1

Submitted on 21 Jan 2025

**HAL** is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

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

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
Unraveling the Impacts of Submesoscale Thermal and Current Feedbacks on the Low-Level Winds and Oceanic Submesoscale Currents

LIONEL RENAULT,ª MARCELA CONTRERAS,ª PATRICK MARCHESIELLO,ª CARLOS CONEJERO,ª IGOR UCHOA,ᵇ AND JACOB WENEGRATᵇ

(Manuscript received 4 June 2024, in final form 12 September 2024, accepted 21 October 2024)

**ABSTRACT:** In this study, high-resolution coupled ocean–atmosphere simulations are performed over the Gulf Stream to investigate the influence of submesoscale [O(10) km] thermal feedback (TFB) and current feedback (CFB) on the low-level atmosphere and the oceanic submesoscale kinetic energy (SKE). At the submesoscale, TFB and CFB exhibit constructive and destructive effects on wind and surface stress, making this a more complex problem than for the mesoscale [O(100) km]. This coinfluence alters classical coupling coefficients, posing a challenge to isolate individual coupling mechanisms. Here, the feedbacks are isolated separately by removing their imprint on the air–sea coupling fields in dedicated simulations. Both submesoscale TFB and CFB lead to a damping of the SKE. CFB causes eddy killing by drag friction between currents and the atmosphere. However, while eddy killing should be more efficient than its mesoscale counterpart due to a weaker wind response (less re-energization), its effect is hampered by an energization from TFB and by the highly transient nature of submesoscale flow, resulting in a modest 10% reduction in SKE. TFB also contributes to a reduction of SKE, mainly by causing a potential energy sink, associated with turbulent heat fluxes, especially at scales < 10 km. The potential energy sink affects SKE through a decrease of baroclinic energy conversion although this is slightly modulated by an increase in Ekman pumping by submesoscale CFB. Our results emphasize the importance of considering both TFB and CFB at the submesoscale and highlight the limitations of mesoscale CFB parameterizations for submesoscale applications. Future parameterizations should be scale-aware and account for both TFB and CFB effects on momentum and heat fluxes.

**KEYWORDS:** Ocean; North Atlantic Ocean; Ocean dynamics; Small scale processes; Atmosphere-ocean interaction

**1. Introduction**

The ocean supports motions from a few meters to thousands of kilometers. Mesoscale eddies—the ocean equivalent of the synoptic weather systems—with horizontal scales of a few hundred kilometers (e.g., Wunsch and Ferrari 2004) are the most energetic flows. In recent decades, numerous studies have demonstrated their importance not only in regulating regional marine ecosystems (McGillicuddy et al. 2007; Mahadevan 2016; Lévy et al. 2018; Renault et al. 2016a; Kessouri et al. 2020), in modulating the physical carbon pump (Omand et al. 2015; Harrison et al. 2018), but also in determining the characteristics of the western boundary currents (McWilliams 2008; Chassignet and Marshall 2008; Renault et al. 2019b). Submesoscale ocean currents, on the other hand, are a relatively recently studied phenomenon, but their popularity as a research topic is growing rapidly (McWilliams 2016). Submesoscale currents arise from various mechanisms such as frontogenesis, frontal instabilities, topographic wakes, or interactions between internal waves and currents (Mahadevan and Tandon 2006; Capet et al. 2008c,a,b; Thomas et al. 2008; McWilliams 2016; Srinivasan et al. 2017; Contreras et al. 2023b,a). They occur on an intermediate scale on the order of meters to tens

of kilometers horizontally, 10–100 m vertically, and hours to days temporally, i.e., about an order of magnitude smaller than the mesoscale currents, but still large enough for rotation and density stratification to matter. Submesoscale currents are characterized by large Rossby numbers [O(1)], drive vertical velocities one to two orders of magnitude greater than those associated with mesoscale currents (Capet et al. 2008c; Su et al. 2020; Siegelman et al. 2020; McWilliams 2021), and can influence momentum, buoyancy, biogeochemical material, and ocean–atmosphere gas exchange (Su et al. 2018; Lévy et al. 2009; Kessouri et al. 2020, 2022). They can also affect the oceanic forward and inverse cascades, i.e., the interior energy pathway (McWilliams et al. 2001; Boccaletti et al. 2007; Capet et al. 2008b; Thomas et al. 2013; Schubert et al. 2020; Contreras et al. 2023b). While some of the submesoscale activity occurs at depth (Wenegrat et al. 2018a; Siegelman et al. 2020), most of it occurs within the oceanic surface layer and thus near the ocean–atmosphere interface, where intense exchanges of heat, freshwater, and momentum take place.

In recent decades, satellite observations and numerical models have shown that mesoscale ocean–atmosphere interactions affect both the atmosphere and the ocean, modifying winds, clouds, humidity, and ocean properties (Chelton et al. 2001; Renault et al. 2017a, 2016d; Seo et al. 2016; Ma et al. 2016; Frenger et al. 2013; Wenegrat and Arthur 2018; Seo et al. 2023). So far, the focus has been on two main ocean–atmosphere interaction processes: the thermal feedback (TFB) and the current feedback (CFB) of the ocean to the atmosphere. TFB is related to the influence of sea surface temperature (SST) on

<footnote>ª Université de Toulouse, LEGOS (CNES/CNRS/IRD/UPS), Toulouse, France</footnote>

<footnote>ᵇ Department of Atmospheric and Oceanic Science, University of Maryland, College Park, College Park, Maryland</footnote>

<footnote>Denotes content that is immediately available upon publication as open access.</footnote>

<footnote>Corresponding author: Lionel Renault, lionel.renault@ird.fr</footnote>
the atmosphere [see Small et al. (2008) and Seo et al. (2023) for a review]. TFB modifies directly the turbulent heat fluxes (Moreton et al. 2021) and has a “top-down” effect on the low-level wind and the surface stress (Chelton et al. 2004; Skyllingstad et al. 2007; Renault et al. 2019c) by altering the atmospheric boundary layer turbulence through downward momentum mixing and pressure gradient adjustment (Seo et al. 2023; Desbiolles et al. 2023). CFB is related to the oceanic surface currents and their influence on surface stress and heat fluxes. Unlike TFB, CFB directly modifies the surface stress (Bye 1985; Chelton et al. 2001) and has a “bottom-up” effect on the wind: a positive current anomaly produces a negative stress anomaly, which in turn generates a positive wind anomaly (Renault et al. 2016d). At the mesoscale, the wind response is important to consider because it partially counteracts the direct surface stress response to the currents (Renault et al. 2016d). To quantify the efficiency of the mesoscale coupling between SST, surface currents, surface stress, and low-level wind, coupling coefficients were defined and then estimated from both satellite observations and coupled simulations. TFB mesoscale coupling coefficients are estimated from the linear relationships between the wind stress curl and crosswind component of SST gradients and between the wind stress divergence and downwind component of SST gradients (Chelton et al. 2001, 2004; Renault et al. 2019c). CFB mesoscale coupling coefficients are only related to the curl of mesoscale currents, surface stress, and low-level wind because mesoscale currents are mostly nondivergent (Renault et al. 2016d, 2019c).

Several studies (Ma et al. 2016; Bishop et al. 2020; Renault et al. 2023b; Holmes et al. 2024) suggest that mesoscale TFB, by causing heat flux anomalies, reduces the potential energy available in the ocean in favor of the atmosphere, affecting energy conversion by baroclinic instability, and hence mesoscale activity. Other studies show that lower eddy activity slows the rate of restratification by eddies, thus cooling the regional SST (Shan et al. 2020). Alternatively, CFB slows down the mean oceanic currents (e.g., Luo et al. 2005; Renault et al. 2016a) and induces damping of the mesoscale activity by ≈30% by acting as an *eddy killer*, that is, by causing a sink of kinetic energy from the eddies to the atmosphere (e.g., Dewar and Flierl 1987; Duhaut and Straub 2006; Eden and Dietze 2009; Renault et al. 2016d,c; Xu and Scott 2008). Mesoscale CFB also exerts significant control over the dynamics of western boundary currents, as lower mesoscale activity translates into weaker eddy–mean flow interaction, i.e., a reduced inverse cascade of energy (Renault et al. 2016c, 2017b, 2019b, 2023a; Larrañaga et al. 2022). Finally, both TFB and CFB cause additional Ekman pumping (e.g., Gaube et al. 2015). However, at the mesoscale, these additional vertical velocities do not significantly modify baroclinic energy conversion (Renault et al. 2023b).

Very little is known about air-sea interactions at the sub-mesoscale. Renault et al. (2018), for the U.S. West Coast, suggest that the eddy killing effect is still active at the submesoscale but at a lower level than at the mesoscale (17% vs more than 30% at the mesoscale for the U.S. West Coast during summer). At these scales, the coupling-induced surface

stress anomaly causes significant Ekman pumping, which can alter baroclinic energy conversion (Renault et al. 2018). CFB-induced surface stress anomaly can also alter the Ekman transport of buoyancy at fronts and thus the potential vorticity budget of the surface oceanic mixed layer (Wenegrat 2023). In addition, Uchoa et al. (2024, manuscript submitted to *J. Phys. Oceanogr.*) show that TFB causes a 10% damping of submesoscale currents by inducing a potential energy sink, an effect hardly present at the mesoscale.

Nevertheless, the response of the atmosphere to both sub-mesoscale TFB and CFB is poorly understood, and it is still unclear whether the coupling coefficients mentioned above behave similarly at the submesoscale as they do at the mesoscale. Recently, using a submesoscale permitting ($dx \approx 3$ km) coupled simulation in the Southern Ocean, Bai et al. (2023) suggest that submesoscale TFB and CFB affect wind and surface stress anomalies, either reinforcing or opposing each other. This could modify the coupling coefficients and the eddy killing process. Conejero et al. (2024) find a similar result for the Caribbean Sea using a submesoscale permitting (1 km) coupled configuration. Not surprisingly, but in contrast to the mesoscale, the authors found a linear relationship between submesoscale surface current divergence and surface stress divergence. They also show that the wind response to CFB (re-energization, secondary to the direct effect on stress) is weaker at the submesoscale than at the mesoscale, which tends to increase the stress response. This would imply more efficient eddy killing but would then contradict the results of Renault et al. (2018) for the U.S. West Coast.

This study focuses on the Gulf Stream (GS), a western boundary current that has a significant impact on our climate. At Cape Hatteras, the GS separates from the coastal boundary and extends horizontally into the offshore ocean. This region is known for its strong air-sea interaction (e.g., Czaja et al. 2019; Masunaga et al. 2020) and intense mesoscale and submesoscale activities (Contreras et al. 2023b; Callies et al. 2015), making it a hotspot for submesoscale interaction with the atmosphere. These interactions are also important for the seasonal cycle of subtropical mode water subduction, linking surface processes to interior variability (Wenegrat et al. 2018b). In this study, submesoscale ocean–atmosphere simulations are used with a spatial resolution of 700 m for the ocean and 2 km for the atmosphere. They are run over the GS for a period of 1 year. The simulations account for different degrees of coupling by smoothing or not smoothing the submesoscale SST or surface current fields sent to the atmosphere. The aim is to properly unravel the effects of the TFB and CFB on both the low-level atmosphere and the ocean.

Based on this methodology, our objectives are twofold and are centered around two main questions: Does the low-level wind and surface stress respond to oceanic submesoscale features, and if so, how? What are the mechanisms by which TFB and CFB influence oceanic submesoscale activity, and do these mechanisms work together or against each other? The paper is organized as follows: section 2 sums up our main methodology. Section 3 evaluates the realism of the simulations, characterizes the prevalent submesoscale activity within the study region, and assesses the degree to which it is
<figure><img src="image_2.png" /><figcaption>FIG. 1. Domain configuration and example of online filtering for (a), (b), (c) SST and (d) surface zonal current. (left) Snapshot of the SST from NOTFB_SUB simulated by (a) CROCO and sent to (b) WRF. (c), (d) As in (a) and (b), but for the zonal flow from the NOCFB_SUB simulation. (e) Spatial spectrum of SST from NOTFB_SUB simulated by CROCO (black) and seen by WRF (red). Using an online spatial filter, the NOTFB_SUB and NOCFB_SUB simulations ignore submesoscale TFB and CFB, respectively.</figcaption></figure>

damped by submesoscale air-sea interactions. Section 4 analyzes the response of low-level winds and surface stress to both CFB and TFB and assesses how these combined effects influence classical coupling coefficients. Section 5 explores the mechanisms driving the response of submesoscale currents to submesoscale air-sea interactions. In particular, we will reconcile the contradiction between a theoretically more effective eddy killing and a weaker damping at the submesoscale than at the mesoscale by evaluating the submesoscale damping time scale associated with CFB. Finally, section 6 synthesizes the key findings and discusses their limitations and broader implications.

## 2. Methods

### *a. Coupled ocean-atmosphere model*

#### 1) MAIN OCEANIC CONFIGURATION

This study uses the Coastal and Regional Ocean Community (CROCO), version 1.3.1, model (Shchepetkin and McWilliams 2005; Debreu et al. 2012) for performing oceanic simulations. CROCO is a free-surface, terrain-following, hydrostatic and Boussinesq model utilizing a split-explicit time-stepping scheme. Its accuracy is achieved using a third-

order predictor-corrector time-stepping algorithm and high-order numerical discretization for pressure gradients and momentum advection (fifth-order upstream scheme). For tracer advection, a third-order upstream-biased advection scheme is implemented in a split mode, the diffusive part being rotated along isopycnal surfaces to minimize excessive diapycnal mixing (Marchesiello et al. 2009; Lemarié et al. 2012). The K-profile parameterization (KPP) scheme by Large et al. (1994) is used to represent unresolved turbulent diffusion at the surface, bottom, and interior of the ocean (Large et al. 1994).

CROCO is implemented over the GS region, focusing on its zonal extent, as depicted in Fig. 1. The computational domain spans from 78.6° to 55.27°W and 34° to 42.9°N, discretized on a horizontal grid with 2942 × 1430 points, resulting in a nominal grid spacing of approximately 700 m. The model bathymetry is derived from the SRTM global bathymetry and elevation data at 15-arc-s resolution with data voids filled (SRTM15_PLUS) data set which is based on the 15-arc-s global dataset by Sandwell and Smith (1997). To minimize aliasing and ensure smoothness at the grid scale, a Gaussian smoothing kernel with a width 4 times the grid spacing is applied to the bathymetry. Mitigating pressure gradient errors is achieved through a dedicated numerical scheme by Shchepetkin and McWilliams (2003). Additionally, local smoothing of the
bottom topography is applied where the slope exceeds an r-factor value of 0.2. In the vertical dimension, the model employs 80 sigma levels with stretching parameters of $\theta_s = 7$, $\theta_b = 2$, and $h_{\text{cline}} = 200$ m. The bottom drag parameterization assumes a classic vertical logarithmic profile for flow within the bottom boundary layer, defined by a friction velocity $u^*$ and a bottom roughness length $Z_{0b} = 0.5$ cm (Renault et al. 2023a).

The model uses a combination of boundary condition schemes with mixed active-passive properties (Marchesiello et al. 2001). It includes a modified Flather-type scheme for the barotropic mode and an Orlanski-type radiation scheme for the baroclinic mode and tracers, which involves boundary forcing from external data. Boundary forcing data are interpolated from a high-resolution (2 km) simulation with tides from Contreras et al. (2023a). They have a temporal resolution of 3 h, allowing the propagation of external and internal tides into the computational domain.

The numerical solution is first run for a 3-month spinup period in uncoupled mode starting from 1 January 2005 of the 2-km simulation with tides from Contreras et al. (2023a). Three months is sufficient for the ocean model to adjust from a resolution of 2 km to a resolution of 700 m. In this case, atmospheric forcing is derived from hourly fifth major global re-analysis produced by ECMWF (ERA5) data (Heiderich and Todd 2020), and turbulent freshwater, heat, and momentum fluxes are estimated with the bulk parameterization scheme of Fairall et al. (2003). CFB is parameterized following the stress correction approach described in Renault et al. (2020). After spinup, the ocean model is still forced at its lateral boundaries by the 2-km simulation of Contreras et al. (2023a), but it is then coupled at the surface with the Weather Research and Forecasting (WRF) Model through the OASIS3-MCTV3 coupler (Craig et al. 2017) from April 2005 for a period of 15 months and with a coupling frequency of 1 h. The first 3 months are considered a spinup period for the coupled model and are excluded from further analysis. During the coupled simulation, WRF provides hourly momentum, heat, and freshwater fluxes to drive the ocean model, which are estimated using the same bulk formulas mentioned above.

## 2) MAIN ATMOSPHERIC CONFIGURATION

The WRF Model (version 4.2; Skamarock et al. 2008) is implemented in a nested grid configuration with two domains. The larger domain encompasses the northwest Atlantic Ocean, spanning from 22° to 49.1°N and 82.5° to 35.4°W, with a horizontal resolution of approximately 6.2 km (Contreras 2023). This parent domain simulation runs from January 2005 to December 2007, forced by hourly lateral boundary conditions from the ERA5 dataset. SST data are obtained from the Operational Sea Surface Temperature and Sea Ice Analysis (OSTIA) satellite product (Donlon et al. 2012). A finer child domain nested within the parent domain focuses on the GS region, with a horizontal resolution of 2 km. Both WRF simulations utilize 50 hybrid eta levels with a model top pressure set at 1000 Pa and a first level at 10 m over the ocean.

The physics suite employed in the parent domain simulation consists of the WRF single-moment 6-class (WSM6)

microphysics scheme with droplet concentration inclusion (Hong and Lim 2006; Jousse et al. 2016), the Korea Institute of Atmospheric Prediction Systems simplified Arakawa-Schubert (KIAPS SAS) (KSAS) convective scheme (Han and Pan 2011; Kwon and Hong 2017), the Rapid Radiative Transfer Model longwave radiation scheme (Mlawer et al. 1997), the Dudhia shortwave radiation scheme (Dudhia 1989), the Yonsei University (YSU) planetary boundary layer scheme (Hong et al. 2006), the revised MM5 surface layer scheme (Jiménez et al. 2012), and the Noah land surface model (Tewari et al. 2004). CFB effect is implemented in both the surface layer and planetary boundary layer schemes, following the approach of Renault et al. (2019a). The finer 2-km domain simulation is driven by 3-hourly lateral boundary conditions from the coarser simulation and coupled to the CROCO ocean model. The parameterization schemes are nearly identical to the parent configuration, except for cumulus parameterization. In the 2-km simulation, only a shallow convection parameterization is retained, in this case from the Tiedtke scheme (Tiedtke 1989; Zhang et al. 2011).

### b. Online filtering procedure and experiments

Three 15-month coupled simulations are conducted starting from April 2005. All experiments use surface stress estimated with the wind relative to the oceanic current ($U_{10rel}$). They differ only in the extent of submesoscale coupling between the ocean and atmosphere:

* In the control (CTRL) simulation, CROCO provides WRF with hourly averages of surface current and SST data. Conversely, WRF provides CROCO with hourly momentum (surface stress), heat, and freshwater flux data.

* In the NOTFB_SUB experiment, the atmosphere receives an SST field smoothed by a two-dimensional spatial Gaussian filter with a cutoff of around 30 km, which effectively eliminates the submesoscale signal. Figures 1a and 1b show an example of the SST simulated by CROCO in the NOTFB_SUB experiments as perceived by WRF after smoothing. The corresponding spatial spectrum is shown in Fig. 1e. This simulation removes the submesoscale TFB effect, emphasizing the CFB effect.

* NOCFB_SUB is similar to NOTFB_SUB, but instead of SST, the zonal and meridional surface currents are smoothed before being sent to WRF (Figs. 1c,d show an example of surface zonal current in NOCFB_SUB). This simulation removes the submesoscale CFB effect, emphasizing the TFB effect.

Note that NOCFB_SUB and the difference between CTRL and NOTFB_SUB both highlight the effect of submesoscale TFB, but in a way that is only statistically consistent since the simulations are run independently and therefore show slight differences in eddy activity due to nonlinear behavior (the same remark holds for the effect of CFB). We tested several filter sizes and chose the 30-km threshold to effectively remove the submesoscale signal without affecting the mesoscale signal too much.

### c. Submesoscale anomalies

To calculate some of the offline diagnostics, the submesoscale anomalies are isolated from other signals, including
<table><caption>TABLE 1. Coupling coefficients.</caption><thead><tr><td>Coefficient</td><td>Description</td></tr></thead><tbody><tr><td>s<sub>τ</sub></td><td>Submesoscale surface current vorticity with surface stress curl</td></tr><tr><td>s<sub>w</sub></td><td>Submesoscale surface current vorticity with 10-m wind curl</td></tr><tr><td>s<sub>τ<sub>δ</sub></sub></td><td>Submesoscale surface current divergence with surface stress divergence</td></tr><tr><td>s<sub>w<sub>δ</sub></sub></td><td>Submesoscale surface current divergence with 10-m wind divergence</td></tr><tr><td>sC<sub>str</sub></td><td>Submesoscale cross-stress SST gradient with surface stress curl</td></tr><tr><td>sC<sub>u</sub></td><td>Submesoscale crosswind SST gradient with 10-m wind curl</td></tr><tr><td>s<sub>thf</sub></td><td>Submesoscale SST and turbulent heat fluxes</td></tr></tbody></table>

mesoscale, large-scale, and inertial waves, by using both a temporal filter and a spatial filter. The temporal filter detects deviations from a 2-day running mean, while the spatial filter is the same as the online Gaussian filter (cutoff of around 30 km) used for SST or surface currents. Submesoscale anomalies are then designated by the prime sign ('). Note that the use of a 1-day running mean does not change qualitatively the results.

### d. Coupling coefficients

To evaluate the impact of TFB and CFB on low-level wind, surface stress, and turbulent heat fluxes, we calculated several coupling coefficients as outlined in Table 1, following the methodology of, e.g., Conejero et al. (2024). The calculation details are provided below.

As previous studies have pointed out (Renault et al. 2018; Conejero et al. 2024), to understand submesoscale coupling processes, it is essential to isolate the atmospheric response to oceanic submesoscale features from other sources of atmospheric variability. This includes sorting out fine-scale atmospheric features that may have similar spatial scales to oceanic submesoscale structures but propagate significantly faster. To estimate the coupling coefficients between the oceanic and atmospheric submesoscale variations, we use a space-time filtering approach similar to the one described in section 2b. However, for temporal smoothing, a 6-h running window is applied to the data. Conejero et al. (2024) show that similar results are obtained using 12- or 24-h running windows. The coupling coefficients are subsequently evaluated for each grid point across the entire simulation period (excluding spinup).

#### 1) SUBMESOSCALE CFB COUPLING COEFFICIENTS

Following established practices in previous studies (Renault et al. 2016d, 2019c; Conejero et al. 2024), we calculate two key CFB coupling coefficients to quantify the interaction between oceanic submesoscale dynamics and atmospheric processes:

* $s_{\tau}$: This coefficient represents the slope of the linear relationship between surface stress curl ($\tau$) and surface current vorticity ($\xi$) anomalies. A more negative value for $s_{\tau}$ indicates a stronger damping effect on submesoscale eddies due to increased eddy killing by the atmosphere. Conversely, a value closer to zero indicates a weaker damping effect, and a positive $s_{\tau}$ indicates a strengthening of the submesoscale flows.

* $s_w$: This coefficient represents the slope of the linear relationship between 10-m wind curl and surface current vorticity anomalies. A more positive value for $s_w$ indicates a stronger wind response to submesoscale ocean currents, implying a more pronounced partial re-energization of the ocean by the atmosphere (Renault et al. 2016d).

Building on the work of Conejero et al. (2024), we additionally compute submesoscale coupling coefficients associated with surface current divergence:

* $s_{\tau_{\delta}}$: This coefficient represents the slope of the linear robust regression (Li 1985) between surface stress divergence and surface current divergence anomalies.

* $s_{w_i}$: This coefficient represents the slope of the linear robust regression between 10-m wind divergence and surface current divergence anomalies.

Those coupling coefficients can be interpreted as the efficiency of how surface current divergence impacts the surface stress and low-level wind divergence. Note that these divergence-related coupling coefficients are not calculated at the mesoscale because mesoscale ocean currents are essentially geostrophic (horizontally nondivergent).

#### 2) TFB COUPLING COEFFICIENTS

To assess the influence of submesoscale SST on the surface stress and on the surface turbulent heat fluxes, following O'Neill et al. (2012), Chelton et al. (2007), and Conejero et al. (2024), we estimate the TFB coupling coefficients as follows:

* $s_{C_{str}}(s_{C_u})$: The slope of the linear robust regression between surface stress (wind) curl and cross-stress (crosswind) SST gradient.

* $s_{d_{str}}(s_{d_u})$: The slope of the linear robust regression between surface stress (wind) divergence and down-stress (downwind) SST gradient.

* $s_{thf}$: The slope of the linear robust regression between turbulent heat flux (sensible plus latent, defined as positive upward) and SST anomalies.

### e. Spectral analysis

In the following, the cospectra and spectra of the data are focused on the GS region, highlighted by a black box in Fig. 4a. To prepare the data for fast Fourier transform (FFT) spectral analysis, we use a two-step preprocessing approach. The first step is to remove the spatial average of the data over the GS (black box). Next, windowing is applied to deal with spectral leakage, due to the abrupt start and end of the data segment creating unwanted spectral components polluting the frequency domain. A Hanning window is applied to the data to smoothly taper it to zero at both ends of the region of interest.
<figure><img src="image_3.png" /><figcaption>FIG. 2. CTRL realistically represents the mean state of the GS region. (top) Mean sea surface height from (a) CTRL and (b) AVISO, (middle) mean SST from (c) CTRL and (d) OSTIA, and (bottom) mean surface stress curl (color) as well as surface stress intensity and direction (vectors) from (e) CTRL and (f) SCOW.</figcaption></figure>

**3. Submesoscale activity and damping by submesoscale air-sea interactions**

a. *Evaluation of the CTRL experiment*

1) GENERAL EVALUATION

An evaluation of the 2-km parent solution to the CTRL
ocean-atmosphere simulation is presented against available
observational data in Contreras (2023) and Contreras et al.
(2023a). Here, we briefly evaluate the higher-resolution
CTRL simulation. As shown in Figs. 2a and 2b, CTRL effec-
tively reproduces the mean sea surface height obtained from
the AVISO/ Developing Use of Altimetry for Climate Studies
(DUACS) climatology dataset (CNES-CLS18) (Taburet et al.
2019) and captures the stable path of the GS as it turns away
from Cape Hatteras. Some discrepancies between CTRL and
observations are visible, however, but may be due to the lim-
ited period considered (1 year) as much as to model biases.
Note that the simulation adequately captures mesoscale activity
(not shown), consistent with the results presented in Contreras
et al. (2023b). In particular, intense mesoscale activity, reaching
up to 3000 m² s⁻², is produced over the GS path after separa-
tion (Contreras 2023). Figures 2c and 2d show the mean SST as
simulated in CTRL alongside data from the OSTIA product
(corresponding to the same time period as the simulation)
(Donlon et al. 2012). Again, a general agreement between
CTRL and observational data is observed, with warm water

in the southern domain and cooler water in the northern
domain. The GS shows the characteristic advection of warm
water from the Florida coast toward the northwest of the do-
main. However, there is a warm bias in the northern domain,
likely due to inherited characteristics from the 2-km solution
used to drive the model.

Figures 2e and 2f show the mean surface stress curl (colors) as well as the surface stress intensity and direction (vectors) as simulated by CTRL compared to climatological data from the Scatterometer Climatology of Winds (SCOW) product (Risien and Chelton 2008). CTRL effectively captures the primary atmospheric low-level circulation characterized by northwesterly winds blowing from the northwest to the southeast before veering eastward at ≈63°W. Consistent with findings by Chelton et al. (2001) and Renault et al. (2016c), the negative surface stress curl along the GS reflects the influence of the mean surface currents. It is worth noting that SCOW appears smoother than CTRL due to its coarser resolution and also because SCOW is a climatology based on a longer time averaging (1999–2009 compared to 1 year for CTRL).

Figure 3 shows the different components of the heat fluxes
simulated by CTRL compared to data from the objectively
analyzed air-sea fluxes (OAFlux; spatial resolution of 1°)
product for the same time period (Jin and Weller 2008).
CTRL realistically simulates turbulent heat fluxes (latent and
sensible) as well as longwave and shortwave fluxes, suggesting
<figure><img src="image_4.png" /><figcaption>FIG. 3. The heat fluxes simulated by CTRL are in good agreement with OAFlux. The positive values indicate a heat input to the ocean. (left) The simulated heat fluxes from CTRL and (right) the observed heat fluxes from OAFlux. (a),(b) Sensible heat flux, (c),(d) latent heat flux, (e),(f) shortwave heat flux, (g),(h) longwave heat flux, and (i),(j) total net heat flux.</figcaption></figure>

a reasonable representation of cloud cover. However, differences between CTRL and OAFlux are visible, such as an excessive latent heat flux over the GS, possibly due to model biases or the coarse spatial resolution of OAFlux. In addition, the shortwave radiation over the GS appears to be slightly overestimated (by about 20 W m⁻²). Overall, the realistic representation of the mean state of the GS system provides confidence in the model's ability to represent air-sea coupling processes. Similar conclusions are drawn for the other experiments.

## 2) SUBMESOSCALE ACTIVITY

As a measure of the intensity of the submesoscale activity, we computed the surface submesoscale kinetic energy (SKE) as

$$ \text{SKE} = \frac{u_o'^2 + v_o'^2}{2}, \qquad (1) $$

where $u_o$ and $v_o$ are the zonal and meridional surface currents, respectively, and the ' represents anomalies computed from a 30-km Gaussian spatial filter and a 2-day running mean to remove large-scale features (see section 2). Note that the spatial filter is the same as that used in the numerical experiments NOTFB_SUB and NOCFB_SUB.

Figure 4a shows the mean SKE estimated from CTRL. As expected, CTRL is characterized by intense submesoscale activity, mainly located along the GS path, which is consistent with the previous results of Callies et al. (2015) and Contreras et al. (2023b). The kinetic energy (KE) spectra is furthermore estimated as

$$ \text{CS}_{\text{KE}} = \frac{\int_{-100}^{0} 0.5(\mathbb{R}[u_o \hat{u}_o^* + v_o \hat{v}_o^*])}{\int_{-100}^{0} 1}, \qquad (2) $$

where the caret (^) denotes the two-dimensional Fourier transform, the symbol R is the real component of the spectra, and the asterisk (*) represents the complex conjugate operator. The CROCO simulations have a spatial resolution of 700 m and use a fifth-order momentum advection scheme, which allows for resolving a large part of the submesoscale range. Consequently, the KE spectrum integrated over 100-m depth (Fig. 4d), where the submesoscale activity is most pronounced (Contreras et al. 2023b), and exhibits a spectral slope between $k^{-3}$ and $k^{-2}$ as expected (Capet et al. 2008b). A discontinuity in the spectrum around 3-4 km (dashed line in Fig. 4d) roughly indicates the effective resolution of CROCO simulations, consistent with the expected capabilities of the
<figure><img src="image_5.png" /><figcaption>FIG. 4. SKE in the different experiments: (a) CTRL, (b) NOCFB_SUB, and (c) NOTFB_SUB. (d) Kinetic energy cospectrum from CTRL. The slopes $k^{-2}$ and $k^{-3}$ are shown in dot-dashed and dashed lines, respectively. (e) Kinetic energy cospectrum difference between NOCFB_SUB (NOTFB_SUB) and CTRL. The dashed line represents approximately the effective resolution of the oceanic simulations. All the experiments are characterized by a large submesoscale activity, mainly situated along the GS path. The simulation that neglects the submesoscale CFB effect (NOCFB_SUB) has a higher SKE. The simulation that neglects the submesoscale TFB effect (NOTFB_SUB) has a higher SKE at scales lower than 10 km.</figcaption></figure>

fifth-order advection scheme (Ménesguen et al. 2018; Soufflet et al. 2016).

Figure 5a shows a snapshot of surface relative vorticity (Rossby number, $\zeta/f$), and Fig. 5b shows a snapshot of surface divergence ($\delta/f$) during winter, the most active season for

submesoscale currents (Callies et al. 2015; Contreras et al. 2023b). The GS north wall is clearly visible in both the vorticity and divergence fields. Large anticyclonic and cyclonic vortices separate from the GS and propagate westward, surrounded by submesoscale eddies, fronts, and filaments.

<figure><img src="image_6.png" /><figcaption>FIG. 5. Characteristics of curl and divergence of submesoscale currents. Snapshots of (a) $\zeta/f$ and (b) $\delta/f$ from CTRL. The PDFs of (c) $\zeta/f$ and (d) $\delta/f$ for CTRL, NOCFB_SUB, and NOTFB_SUB. Both submesoscale TFB and CFB reduce the intensity of strong (large Rossby radius) submesoscale currents.</figcaption></figure>
**Figures 5c and 5d** show the probability density function (PDF) of the surface $\zeta/f$ and $\delta/f$ over the analyzed 1-yr period (without any filtering). Consistent with the literature (e.g., Shcherbina et al. 2013; Srinivasan et al. 2017), $\zeta/f$ shows a strong positive skewness, indicating a predominance of cyclical vorticity, especially in filaments. The $\zeta/f$ can reach values greater than 5, indicating a very strong unbalanced flow. The $\delta/f$ is also characterized by negative skewness although less pronounced than for $\zeta/f$. In summary, the CTRL simulation successfully reproduces intense submesoscale activity around the GS extension.

### b. Damping of submesoscale activity by both TFB and CFB

**Figures 4a-c** show the mean surface SKE from CTRL, NOCFB_SUB, and NOTFB_SUB simulations. Compared to NOCFB_SUB, CTRL shows a damping of surface SKE of about 10%, which is significantly lower than that of mesoscale eddies (about 30%, Renault et al. 2016c). The damping extends below the surface (down to 100-m depth, not shown), as illustrated by the difference ratio of the 100-m depth-integrated KE spectra between NOCFB_SUB and CTRL (Fig. 4e). The weaker damping of KE by CFB in CTRL compared to its effect on the mesoscale implies a reduced eddy killing activity at submesoscale.

The surface SKE in NOTFB_SUB appears similar to CTRL. However, a closer look at the 100-m depth-integrated KE cospectra in Fig. 4e reveals a different story. While the 10-30-km range in $CS_{KE}$ of CTRL and NOTFB_SUB shows similar energy levels, NOTFB_SUB shows a progressive increase in energy at scales finer than 10 km. This increase reaches 20% at 4 km and up to 50% at even finer scales. Note that the model's numerical dissipation range starts at these scales (<4 km). In conclusion, TFB appears to damp submesoscale currents but only at scales finer than approximately 10 km. This damping effect is even stronger at these scales (but not on the whole) than that induced by CFB.

**Figures 5c and 5d** show the PDFs of surface $\zeta/f$ and $\delta/f$ for the CTRL, NOCFB_SUB, and NOTFB_SUB simulations (note that the mesoscale signal is not included). As expected, all three simulations have similar overall characteristics, with a positive skewness in the $\zeta/f$ distribution (around 2) and a negative skewness in the $\delta/f$ distribution. Consistent with the damping effect of CFB on submesoscale flow, the distribution of $\zeta/f$ in NOCFB_SUB is slightly more skewed toward higher Rossby numbers (values above 2 and up to 10). The effect on the $\delta/f$ distribution is less pronounced. Neglecting submesoscale TFB effect in NOTFB_SUB increases the occurrence of very strong positive Rossby numbers (exceeding 5), which is consistent with the increased KE at scales finer than 10 km observed in Fig. 4e. The absence of the submesoscale TFB effect also has a small influence on $\delta/f$, promoting a slightly weaker current divergence. This could be linked to wind or heat flux anomalies associated with TFB. As a result, the skewness of $\zeta/f$ in both NOCFB_SUB and NOTFB_SUB is greater than that in CTRL. However, NOTFB_SUB is slightly more skewed, which can be explained by the higher frequency

of large Rossby numbers (>6). For $\delta/f$, because CFB tends to increase the occurrence of negative values, while TFB tends to increase the occurrence of positive values, the skewness of NOCFB_SUB and NOTFB_SUB are similar.

In the following sections, the mechanisms governing the reduction of SKE by CFB and TFB are discussed, with particular emphasis on the atmospheric response to these submesoscale forcing.

## 4. Atmospheric response to oceanic submesoscale features

### a. Coaction of TFB and CFB on surface stress and low-level wind

Our set of simulations enables us to isolate and analyze the wind and surface stress responses to CFB and TFB. In the following, we examine how these processes interact, similar to modeling studies of the Antarctic Circumpolar Current (Bai et al. 2023) and the Caribbean Sea (Conejero et al. 2024). Following Bai et al. (2023), we evaluate these responses by bin-averaging hourly submesoscale surface stress (wind) curl data based on both ocean vorticity and cross-stress (cross-wind) SST gradients. To isolate the impact of submesoscale TFB, we first analyze Figs. 6a and 6b, which depict data from the NOCFB_SUB simulation. This simulation neglects submesoscale CFB, allowing us to assess TFB's influence on wind and surface stress in isolation. As expected, a positive correlation exists between SST gradients and surface stress/wind curl. Negative cross-stress (wind) SST gradients induce negative stress (wind) curl, and vice versa. This is because TFB has a top-down effect, initially impacting the atmospheric boundary layer and wind, which then influences surface stress. Notably, the absence of submesoscale CFB in the NOCFB_SUB simulation results in no dependence of the results on surface vorticity (vertical white line). While downward momentum mixing and pressure adjustments are likely mechanisms, their detailed contributions are outside the scope of this study.

We now shift our focus to the NOTFB_SUB experiment, which deliberately excludes the submesoscale TFB effect, allowing us to isolate the wind and surface stress response to submesoscale surface currents (Figs. 6c,d). As expected, the wind and stress curl anomalies in this scenario show no dependence on SST gradients, confirming the absence of TFB. In contrast to the top-down influence of TFB, CFB operates in a bottom-up manner at the submesoscale. Here, as for the mesoscale (Renault et al. 2016d), negative (positive) surface current vorticity induces positive (negative) stress curl, resulting in negative (positive) 10-m wind curl. There is therefore a negative correlation between submesoscale surface currents and surface stress curl, which further extends to a positive correlation with wind curl. Notably, the zero response, in this case, is horizontal, reflecting the absence of submesoscale TFB influence (unlike the vertical white line in Figs. 6a,b).

The CTRL simulation shown in Figs. 6e and 6f illustrates the surface stress and low-level wind responses to both CFB and TFB. Consistent with the previous results, their effects on
<figure><img src="image_7.png" /><figcaption>FIG. 6. Submesoscale TFB and CFB relationship to surface stress and wind curl. (top) The plots show the submesoscale stress curl anomalies as a function of cross-stress SST gradients (x axis) and surface ocean vorticity (y axis) from (a) NOCFB_SUB, (c) NOTFB_SUB, and (e) CTRL. (b), (d), and (f) As in (a), (c), and (e), respectively, but for the 10-m wind curl anomalies as a function of crosswind SST gradients (x axis) and surface ocean vorticity (y axis).</figcaption></figure>

the wind response tend to be out of phase. Similar to the NOCFB_SUB experiment, TFB in CTRL induces a positive relationship between surface stress (wind) curl and cross-stress (crosswind) SST gradient for a given surface vorticity. This translates into larger stress and wind curl for larger SST gradients, and as expected, the wind response remains positively correlated with the stress response due to the top-down influence of TFB. However, CFB introduces a contrasting effect. For a fixed SST gradient, vorticity is anticorrelated to stress curl, but the wind response to CFB (anticorrelated to the stress response) is opposite to that of TFB. Consequently, the zero response line in both stress and wind curl (white line) tilts between the dominant regimes of TFB and CFB. This finding is consistent with modeling studies by Bai et al. (2023) and Conejero et al. (2024): TFB and CFB can have a constructive or a destructive combined effect on the surface stress, and their combination tends to phase out the wind and stress responses. For example, when acting constructively on surface stress, TFB and

CFB synergistically produce the most anomalous stress values but with a damped wind response.

Similar results are found for the divergence of the surface stress (wind) curl and down-stress (downwind) SST gradients (not shown). Given the combined effect of submesoscale TFB and CFB on the surface stress and wind curl, we expect contamination of the CFB coupling coefficients by TFB and contamination of the TFB coupling coefficients by CFB.

### b. CFB coupling coefficients

Figure 7 shows the CFB submesoscale coupling coefficients: $s_τ$ for surface stress (left column) and $s_w$ for wind (right column)—estimated from the different experiments. In CTRL, $s_τ$ has negative offshore values (Fig. 7a) ranging from -0.005 to -0.02 N s m⁻³, in contrast to the positive $s_w$ values that range from 0.1 to 0.2 (Fig. 7b). This pattern is consistent with previous mesoscale results (Renault et al. 2019c) and recent submesoscale results in the Caribbean Sea by Conejero et al. (2024).
<figure><img src="image_8.png" /><figcaption>FIG. 7. Coupling coefficients reveal that, similar to the mesoscale, submesoscale currents vorticity induces anticorrelated surface stress curl anomalies and correlated 10-m wind curl anomalies. This is illustrated by (left) $s_τ$ and (right) $s_w$ coupling coefficients from (a),(b) CTRL, (c),(d) NOCFB_SUB, and (e),(f) NOTFB_SUB. Both $s_τ$ and $s_w$ are modified by submesoscale TFB, decreasing $s_τ$ (less negative), and increasing $s_w$ (more positive). Note the different colorscale in (c).</figcaption></figure>

The positive values along the coast in both $s_τ$ and $s_w$ are not related to CFB but rather likely to small-scale orographic wind features that drive coastal currents (Renault et al. 2016b; Kessouri et al. 2022). More importantly here, the submesoscale $s_τ$ values are more negative than the mesoscale values, with an average value over the GS of ≈-0.019 N s m⁻³ vs ≈-0.015 N s m⁻³ (Renault et al. 2019c), corresponding to a difference of ≈27%, which is similar to the finding of Conejero et al. (2024) for the Caribbean Sea (30% difference). This indicates a weaker wind response (re-energization), reflected in $s_w$ hovering around 0.15 (compared to 0.3 at the mesoscale). As a result, we expect that a more negative $s_τ$ would imply more efficient eddy damping at submesoscale. However, this assertion needs to be moderated, and we will discuss eddy damping and, in particular, a submesoscale damping time scale in detail in the next section.

The nonuniform spatial distribution of the two coefficients in CTRL is probably due to the influence of TFB. To confirm this hypothesis, Figs. 7c and 7d show $s_τ$ and $s_w$ from the NOCFB_SUB experiment, where submesoscale CFB is excluded. Here, positive values of both coefficients dominate (energization of submesoscale currents), highlighting the influence of SST gradients on the wind and subsequently on the surface stress. This is consistent with Fig. 6, where TFB and CFB are mainly constructive for wind but destructive for surface stress. However, some negative values are present in the southern edge of the GS, reflecting the dominance of negative SST gradients in this region. In the NOTFB_SUB simulation (Figs. 7e,f) where submesoscale SST gradients are not felt by the atmosphere (CTRL to

NOTFB_SUB), both $s_τ$ and $s_w$ are more uniform, with negative $s_τ$ and positive $s_w$.

In contrast to their predominantly nondivergent behavior at the mesoscale, submesoscale currents often exhibit a divergent component (see Fig. 5). This property requires additional coupling coefficients, as proposed by Conejero et al. (2024). These coefficients, $s_{δτ}$ and $s_{δw}$, quantify the relationship between surface current divergence and surface stress/low-level wind divergence, respectively. In particular, 10-m wind divergence affects atmospheric vertical velocities, which in turn can affect cloud cover and precipitation patterns.

Figure 8 shows the spatial distribution of $s_{δτ}$ and $s_{δw}$ derived from our simulations. In the CTRL run, negative $s_{δτ}$ and positive $s_{δw}$ values are not surprising (Figs. 8a,b). This indicates, e.g., that positive surface flow divergence leads to negative surface stress divergence, which subsequently induces positive 10-m wind divergence. Interestingly, the wind divergence response ($s_{δw}$) is stronger than the curl response ($s_w$), with average values of 0.3 and 0.15, respectively. This is consistent with Chelton et al. (2001), who observed that the wind divergence response to SST gradients is often 30%–50% larger than the curl response. A more detailed examination of the mechanisms behind the CFB and TFB wind responses is beyond the scope of this study.

In CTRL, both $s_{δτ}$ and $s_{δw}$ show spatial variations. In particular, while larger values are observed along the GS for the 10-m wind, weaker values are present for the surface stress. Again, this can be attributed to the influence of SST gradients on these coupling coefficients. Figures 8c and 8d show $s_{δτ}$ and
<figure><img src="image_9.png" /><figcaption>FIG. 8. As in Fig. 7, but for divergence instead of curl relationship. Submesoscale currents cause divergence anomalies in both surface stress (negative correlation) and wind (positive correlation).</figcaption></figure>

$s_{\delta W}$, respectively, from the NOCFB_SUB experiment, where the direct influence of submesoscale surface flows on the atmosphere is excluded. Consistent with a thermal response, both coefficients become positive in NOCFB_SUB, with larger values along the GS path. This explains the spatial heterogeneity in CTRL and is further confirmed by NOTFB_SUB (Figs. 8e,f). Here, the influence of submesoscale SST on the atmosphere contributes to the heterogeneity in $s_{\delta T}$ and $s_{\delta W}$.

### c. TFB coupling coefficients

#### 1) STRESS AND WIND COUPLING COEFFICIENTS

This section explores the TFB coupling coefficients, $s_{C_{\text{str}}}$ and $s_{C_{\text{w}}}$, representing the surface stress (wind) curl response to submesoscale SST cross-stress (crosswind) gradients (see section 2). Similar to mesoscale dynamics, TFB acts as a top-down process, influencing atmospheric boundary layer turbulence

<figure><img src="image_10.png" /><figcaption>FIG. 9. Surface stress and 10-m wind curl response to submesoscale SST anomalies. Coupling coefficient (left) $s_{C_{\text{str}}}$ and (right) $s_{C_{\text{w}}}$ from (a),(d) CTRL, (b),(e) NOCFB_SUB, and (c),(f) NOTFB_SUB.</figcaption></figure>
<figure><img src="image_11.png" /><figcaption>FIG. 10. Influence of oceanic submesoscale features on the heat fluxes. (a) Snapshot of submesoscale turbulent heat fluxes during winter. (b) RMS of the submesoscale turbulent heat fluxes. (c)–(e) Coupling coefficients between submesoscale SST anomalies and submesoscale turbulent heat fluxes from CTRL, NOCFB_SUB, and NOTFB_SUB [positive values indicate an input for the atmosphere as in Conejero et al. (2024)]. Note the different scale range in (e).</figcaption></figure>

and wind, which then affect surface stress. As defined by Chelton et al. (2004), in CTRL, positive values for $sC_{str}$ and $sC_w$ reflect the positive correlation between wind and surface stress curl anomalies and SST gradients (Figs. 9a,b). Specifically, in CTRL, $sC_{str}$ ranges from 0.005 to 0.015 N m$^{-2}$ s$^{-1}$ while $sC_w$ ranges from -0.05 to 0.15 m s$^{-1}$ $^{\circ}$C$^{-1}$ (Figs. 9a,d). Notably, these submesoscale wind and surface stress curl responses to SST gradients are weaker than those found for the mesoscale in Renault et al. (2019c) ($sC_w$ ≈ 0.2 m s$^{-1}$ $^{\circ}$C$^{-1}$) but align with results on submesoscales from Conejero et al. (2024). Similar patterns hold for divergence (not shown here).

The spatial distribution of $sC_{str}$ in CTRL also shows heterogeneity, particularly pronounced in the southern domain where SKE is highest (Fig. 4). Although $sC_w$ also shows some spatial variation, it is less pronounced than $sC_{str}$. Examination of the other simulations confirms the co influence of CFB on TFB coupling coefficients. In the NOCFB_SUB experiment, where submesoscale surface current effects are excluded, $sC_{str}$ becomes more homogeneous, especially where SKE is highest, where it decreases by about 30% (Figs. 9c,d). This decrease highlights the significant contribution of CFB to the $sC_{str}$ values. The NOTFB_SUB simulation (Figs. 9e,f) further confirms this influence. Since NOTFB_SUB removes submesoscale TFB, $sC_{str}$ should be approximately zero. However, positive values persist in the southern domain, indicating the influence of submesoscale surface currents. Interestingly, the wind curl response ($sC_w$) appears to be less sensitive to surface currents. Values in CTRL and NOCFB_SUB show remarkable similarity, suggesting that weak submesoscale wind responses to CFB may not be fully aligned with SST gradients. Similar results are observed for the divergence of stress and

wind with respect to down-stress/downwind SST gradients (not shown).

## 2) SUBMESOSCALE TURBULENT HEAT FLUXES AND ASSOCIATED COUPLING COEFFICIENT

Beyond its influence on surface stress, wind curl, and divergence anomalies, TFB also acts as a significant driver of mesoscale turbulent heat and moisture fluxes across the air-sea interface. This has been well documented in previous work (e.g., Small et al. 2008; Bishop et al. 2017; Moreton et al. 2021; Renault et al. 2023b; Desbiolles et al. 2023). At the oceanic mesoscale, these turbulent fluxes play a critical role in damping SST anomalies, essentially acting as a sink for heat and oceanic potential energy that are transferred to the atmosphere (Bishop et al. 2020; Renault et al. 2023b; Uchoa et al. 2024, manuscript submitted to J. Phys. Oceanogr.) (discussed further in the next section). Submesoscale features within the ocean are characterized in particular by large SST anomalies (as shown in Figs. 5a and 6). The associated turbulent heat flux anomalies are particularly pronounced. This is illustrated in Fig. 10a, which shows a snapshot of winter submesoscale turbulent heat fluxes, where turbulent heat flux values can exceed -100 W m$^{-2}$, i.e., a transfer of heat from the ocean to the atmosphere. To comprehensively assess the overall importance of submesoscale turbulent heat fluxes, we present the RMS of this quantity in Fig. 10b. Examination of the RMS distribution reveals a range from about 0 to 30 W m$^{-2}$, highlighting the spatial heterogeneity inherent in these fluxes. In particular, the spatial distribution of the RMS is more pronounced over the GS and its northern vicinity,
which reveals regions characterized by intense submesoscale SST variability.

To assess the efficiency of turbulent heat fluxes in damping submesoscale SST anomalies, in Figs. 10c-e, we derived the coupling coefficient $s_{thr}$ following the methodology suggested by Conejero et al. (2024) (see methods, Fig. 10). In the CTRL experiment (Fig. 10a), $s_{thr}$ exhibits strong spatial heterogeneity, with pronounced values concentrated over the GS and south of the GS. The coefficient ranges from 7 to 62 W m$^{-2}$ °C$^{-1}$, with a spatial mean of ≈32 W m$^{-2}$ °C$^{-1}$. Two regions can be distinguished: north of the GS where the values are close to those documented by Uchoa et al. (2024, manuscript submitted to J. Phys. Oceanogr.) for the California region and south of the GS where the values are larger and mirror those found by Conejero et al. (2024) for the Caribbean Sea (around 53 m$^{-2}$ °C$^{-1}$). In the absence of submesoscale CFB, as simulated in NOCFB_SUB and presented in Fig. 10d, the values of $s_{thr}$ show a marginal decrease of up to 10 W m$^{-2}$ °C$^{-1}$. This decrease is particularly noticeable over the GS, where SKE is highest. This finding is further corroborated by estimates from NOTFB_SUB simulations (Fig. 10e). Despite the absence of submesoscale SST anomalies in the calculation of heat fluxes by the atmospheric model, $s_{thr}$ remains slightly positive, with values reaching up to 10 W m$^{-2}$ °C$^{-1}$ over regions characterized by the strongest submesoscale flows (see Fig. 4). However, these values are weaker than those induced by submesoscale TFB, suggesting that for this region, submesoscale CFB is a second-order mechanism in driving the transfer of heat from the ocean to the atmosphere.

## 5. Main oceanic mechanisms

SKE variation can be derived from the total kinetic energy equation, considering the submesoscale range:

$$ \frac{\partial \mathrm{SKE}}{\partial t} = \underbrace{-\nabla \cdot [\overline{\mathbf{u} \mathrm{SKE}}]}_{AK_s} - \underbrace{\frac{1}{\rho_0} \nabla \cdot (\overline{\mathbf{u}' p'})}_{p_s} - \underbrace{\frac{1}{\rho_0} g \overline{p' w'}}_{p_s K_s} - \underbrace{\overline{(\mathbf{u}' \mathbf{u}' \cdot \overline{\nabla \overline{\mathbf{u}}}} + \overline{\mathbf{v}' \mathbf{u}' \cdot \overline{\nabla \overline{\mathbf{v}}}})}_{B_s K} + DK_s + F_s K_s, \quad (3) $$

where prime denotes the submesoscale part of the signal estimated as 2-day running-mean anomalies, also using a spatial Gaussian high-pass filter (see methods). The first and second terms on the right-hand side represent the transport of SKE due to advection ($AK_s$, including the horizontal and vertical terms) and pressure work ($P_s$). Note that these terms can be distributed energy from one region to another, but they cannot generate or dissipate it (Kundu et al. 2015). The fourth term, $B_s K$, is the cross-scale transfer of energy from the background kinetic energy reservoir, resulting from Reynolds stress work, including vertical and horizontal shear terms. The terms $DK_s$ and $F_s K_s$ are the dissipation and forcing terms for SKE, respectively. Finally, $P_s K_s$ represents the conversion of energy from submesoscale potential energy (SPE) to SKE through buoyancy production. SPE is defined as $g^2 \rho'^2/(2\rho_0 N^2)$, where $\rho(T, S, z)$ is a perturbation density relative to a mean density

profile $\rho_r$. The SPE budget is derived by multiplying the budget equation for $\rho$ by $g^2 \rho'/(\rho_0 N^2)$:

$$ \frac{\partial \mathrm{SPE}}{\partial t} = \underbrace{-\nabla \cdot [\overline{\mathbf{u} \mathrm{SPE}}]}_{AP_s} + \underbrace{\frac{1}{\rho_0} g \overline{\rho' w'}}_{-P_s K_s} - \underbrace{\frac{g^2}{\rho_0 N^2} (\overline{\rho' \mathbf{u}'} \cdot \overline{\nabla \rho})}_{B_s P} + DP_s + GP_s, \quad (4) $$

The first term, $AP_s$, is the transport of SPE by advection. The term $B_s P$ is the cross-scale flux of potential energy from the background energy reservoir. The term $P_s K_s$ is the conversion between potential and kinetic energy, already defined in Eq. (3). The terms $DP_s$ and $GP_s$ refer to the dissipation and forcing of potential energy.

As shown in the previous section, at the submesoscale, both TFB and CFB influence the surface stress and turbulent heat fluxes. Consequently, differences in SKE between simulations can be attributed to two main mechanisms: (i) the damping of eddies by $F_s K_s$, i.e., a transfer of kinetic energy from the ocean to the atmosphere, commonly referred to as eddy killing (Renault et al. 2016d), and (ii) changes in baroclinic energy conversion ($P_s K_s$) triggered by a potential energy sink ($GP_s$), as emphasized by Uchoa et al. (2024, manuscript submitted to J. Phys. Oceanogr.), modulated by Ekman pumping resulting from CFB-induced surface stress curl anomalies (Renault et al. 2018). Notably, while the forward cascade of kinetic energy can play a significant role (although second order with respect to the potential energy pathway, Contreras 2023) in generating submesoscale flows, its contribution remains relatively consistent across simulations (not shown). The primary objective of this section is to assess the extent to which these mechanisms explain the coupling effect on SKE outlined in section 3.

### a. Submesoscale kinetic energy flux

As shown in the literature (e.g., Renault et al. 2018), one direct effect of CFB is the transfer of kinetic energy from submesoscale currents to the atmosphere, denoted by $F_s K_s$ in Eq. (3). The term $F_s K_s$ ($m^3 s^{-3}$) is estimated as

$$ F_s K_s = \frac{\overline{(\tau_x u_o' + \tau_y v_o')}}{\rho_0}, \quad (5) $$

where $\rho_0$ is the ocean surface density; $\tau_x$ and $\tau_y$ are the zonal and meridional surface stress, respectively; and $u_o$ and $v_o$ are the zonal and meridional currents, respectively. Figure 11a illustrates $F_s K_s$ estimated from the CTRL experiment. Consistent with the results of Renault et al. (2018) for the California region, this analysis reveals robust pathways of energy transfer from oceanic submesoscale currents to the atmosphere, particularly concentrated over the GS, coinciding with regions of highest SKE (see Fig. 4). When submesoscale CFB is neglected in NOCFB_SUB (Fig. 11b), the energy sinks mostly disappear and are generally replaced by energy inputs from the atmosphere to the submesoscale flow. In NOCFB_SUB, as indicated by the positive coupling coefficient $s_r$, the surface
<figure><img src="image_12.png" /><figcaption>FIG. 11. Submesoscale CFB induces a sink of kinetic energy from submesoscale flows to the atmosphere (negative $F_sK_s$), which is slightly attenuated by the submesoscale TFB influence on surface stress. (a)–(c) $F_sK_s$ from CTRL, NOCFB_SUB, and NOTFB_SUB. (d) Cospectrum of windwork from CTRL, NOCFB_SUB, and NOTFB_SUB. The dashed line represents approximately the effective resolution of the oceanic simulations.</figcaption></figure>

stress anomalies induced by TFB show spatial coherence with the surface flow that drives a positive $F_sK_s$. This has a tapering effect on CFB energy sinks, which is absent (or barely present) at the mesoscale. As confirmation, in NOTFB_SUB, when only the influence of submesoscale currents on the atmosphere is considered, $F_sK_s$ generally exceeds CTRL values due to the absence of TFB energization (Fig. 11c). To further confirm this finding, following, e.g., Jullien et al. (2020) and Uchoa et al. (2024, manuscript submitted to J. Phys. Oceanogr.), we estimated in Fig. 11d the spatial cospectrum of wind work as

$$CS_{WW} = \frac{1}{\rho_0} \left( \mathbb{R}[\hat{\tau}_x \hat{u}_o^* + \hat{\tau}_y \hat{v}_o^*] \right). \quad (6)$$

Consistent with the spatial estimates of $F_sK_s$, the analysis of $CS_{WW}$ in the CTRL experiment reveals a clear transfer of energy from submesoscale currents to the atmosphere, i.e., the eddy killing process. However, in NOCFB_SUB, where the influence of submesoscale surface currents is neglected, the cospectral analysis confirms a slight energization of submesoscale currents by TFB-induced surface stress anomalies. In NOTFB_SUB, which neglects the submesoscale TFB and, thus, this energization effect, there is a more pronounced negative transfer of energy from submesoscale currents to the atmosphere.

At the oceanic mesoscale, the eddy killing process induces a damping of the mesoscale activity by about 30% over the GS (Renault et al. 2016c). However, at the submesoscale

level, despite a more negative $s_τ$, the SKE is only damped by about 10%. While the discrepancy between the two indicators may be due to other processes at play (see next section), it also reflects the highly transitory nature of submesoscale currents. A useful and relevant assessment of the eddy killing effect may be given as in Renault et al. (2016d) by an eddy damping time scale derived from energy considerations:

$$t_{\text{decay}} = \frac{\rho_0 D}{2s_\tau} / 2, \quad (7)$$

where $t_{\text{decay}}$ is primarily determined by D, representing the depth of the eddy, and $s_τ$, the current coupling coefficient. Assuming a typical depth of submesoscale currents of 100 m, a standard seawater density $\rho_0 = 1025 \text{ kg m}^{-3}$, and a typical offshore coupling coefficient $s_τ = 0.017 \text{ N s m}^{-3}$, the eddy decay time scale is about 17.5 days, contrasting with a mesoscale decay time scale of about 125 days (Renault et al. 2016d). At submesoscale, the decay time scale is thus 8.75 times longer than the typical lifetime of submesoscale flow (about 2 days). Note that some submesoscale eddies that last longer will be more affected by CFB. For example, considering a typical, e.g., lifetime of 10 days results in a $t_{\text{decay}}$ that is 1.75 times longer. In contrast, the decay time scale of long-lived mesoscale eddies is only 2.2 times longer than their lifetime exceeding 16 weeks (Chelton et al. 2011). This difference suggests that, despite a more efficient eddy killing mechanism, there is not enough time for a submesoscale eddy to be significantly affected by interaction with the atmosphere before it is
<figure><img src="image_13.png" /><figcaption>FIG. 12. Submesoscale TFB causes a transfer of PE ($GP_s$) between ocean and atmosphere, which is negative over the GS and positive around the GS. Negative $GP_s$ is driven by turbulent heat fluxes, while positive $GP_s$ is driven by salt fluxes. (a)–(c) $GP_s$ from CTRL, NOCFB_SUB, and NOTFB_SUB, respectively. (d) $GP_s$ estimated from CTRL but considering only surface salt fluxes and thus ignoring surface heat fluxes. (e) CSPE from CTRL, NOCFB_SUB, and NOTFB_SUB and, in addition, from CTRL and NOTFB_SUB when only the contribution of salt flux is considered. The dashed line represents approximately the effective resolution of the oceanic simulations.</figcaption></figure>

dissipated by other mechanisms. This likely accounts for the relatively modest 10% reduction of SKE due to CFB over the GS region and reconciles the contradiction between a more efficient eddy killing and a weaker eddy damping.

b. Submesoscale potential energy flux

TFB exerts a significant influence on oceanic buoyancy variability through heat and salt fluxes, manifesting as evaporation and cooling of relatively lighter water masses, or precipitation and heating of denser water masses. This process directly affects the energy budget of mesoscale and submesoscale currents by altering the available potential energy and, consequently, the baroclinic conversion of energy near the surface (von Storch et al. 2012; Bishop et al. 2020; Renault et al. 2023b; Uchoa et al. 2024, manuscript submitted to J. Phys. Oceanogr.).

Following Cronin and Sprintall (2009) and Uchoa et al. (2024, manuscript submitted to J. Phys. Oceanogr.), we first approximate the submesoscale buoyancy flux, denoted $B_0$, considering contributions from both heat flux ($B_{0T}$) and salt flux ($B_{0S}$) components:

$$B'_{0} = B'_{0T} + B'_{0S}, \qquad (8)$$

$$B'_{0} = \frac{\alpha_{\theta} g}{\rho_{0} C_{p}} Q'_{net} - \beta_{s} g [SSS(E - P)]', \qquad (9)$$

where g is the gravity, SSS is the sea surface salinity, $C_p$ is the specific heat of water, $Q_{net}$ is the net surface heat flux

(positive values represent a flux into the ocean), and E and P are the evaporation and precipitation, respectively. Furthermore, $\alpha_\theta$ and $\beta_s$ represent the thermal expansion and salinity contraction coefficients calculated at each point. The submesoscale flux of potential energy $GP_s$ is estimated as

$$GP_s = \frac{\overline{b'_{0}B'_{0}}}{\rho_{0} N_{r}^{2}}, \qquad (10)$$

where $b'_{0}$ is the surface submesoscale buoyancy ($b' = g\rho'/\rho_0$) and $N_r^2 = \partial\langle b_r \rangle / \partial z$ is the reference squared Brunt-Väisälä frequency, where $\langle b_r \rangle = g\rho_r/\rho_0$ is the horizontally and temporally averaged buoyancy over the domain and 1 year. A negative $GP_s$ indicates a loss of potential energy (PE) to the atmosphere, while a positive $GP_s$ indicates a gain by the ocean.

Figure 12a depicts $GP_s$ estimated from the CTRL experiment. It reveals two distinct regions: the GS path, where $GP_s$ is negative (PE loss for the ocean) with a mean value of $-9.1 \times 10^{-9} \text{ m}^4 \text{ s}^{-3}$, and regions beyond the GS, where $GP_s$ is positive (PE gain) with a mean value of $1.5 \times 10^{-10} \text{ m}^4 \text{ s}^{-3}$ in the northern part. This spatial variability differs from that estimated in the California region (Uchoa et al. 2024, manuscript submitted to J. Phys. Oceanogr.), where $GP_s$ is generally positive. Uchoa et al. (2024, manuscript submitted to J. Phys. Oceanogr.) suggest that negative $GP_s$ values primarily arise from heat flux anomalies responding to SST perturbations. These heat flux anomalies reduce $b'$ values and consequently the potential energy available at the ocean surface. To confirm this hypothesis, we computed $GP_s$ in two
additional experiments: NOCFB_SUB and NOTFB_SUB, depicted in Figs. 12b and 12c, respectively. In the NOCFB_SUB experiment, negative $GP_s$ values still prevail along the GS path, indicating a persistent sink of potential energy. However, in the NOTFB_SUB experiment, neglecting submesoscale TFB removes the heat flux response to submesoscale SST, suppressing the potential energy sink. The remaining negative values in NOTFB_SUB are due to the effect of relative wind on heat fluxes, which also contributes to the positive values of $s_{THF}$ in this case (Fig. 10). In addition, we estimated $GP_s$ from the CTRL experiment, but considering only the salt flux contribution $B_{0S}$ to the buoyancy flux (Fig. 12d). In line with the potential energy sink driven by heat fluxes, negative values disappear in $GP_s$ (EMP).

Beyond the GS, a regime change is observed: all simulations show positive $GP_s$, pointing to a transfer of potential energy from the atmosphere to the ocean. The large positive values in the northern part of the domain, evident in both CTRL and NOCFB_SUB, disappear when $GP_s$ is estimated from NOTFB_SUB and from CTRL using only $B_{0S}$ (Figs. 12c,d). Here, the SSS adjusts to the evaporation caused by latent heat flux and the resulting larger $b'$ anomalies lead to positive $GP_s$. The remaining positive values in NOTFB_SUB (and CTRL using only $B_0$) are primarily due to the salt flux contribution [via the freshwater flux, i.e., evaporation minus precipitation (EMP)], which tends to increase the oceanic potential energy (Uchoa et al. 2024, manuscript submitted to J. Phys. Oceanogr.). It is worth noting that some of the large positive values observed in both NOTFB_SUB and CTRL $GP_s$ (EMP) can also be attributed to the interaction between mesoscale (convective) atmospheric vortices and submesoscale oceanic SSS anomalies. These interactions can lead to localized variations in the potential energy, influencing the overall energy budget of the system.

To confirm our results, following Uchoa et al. (2024, manuscript submitted to J. Phys. Oceanogr.), we now evaluate the cospectrum of the potential energy flux by estimating

$$CS_{PE} = \frac{1}{\rho_0} \left[ \frac{\hat{b}_0^* \hat{B}_0}{N_r^2} \right]. \quad (11)$$

This analysis has the advantage of displaying the flux at each wavelength, but the disadvantage of hiding the spatial heterogeneity of the potential energy flux is revealed in Figs. 12a-c. As a result, CTRL and NOCFB_SUB spectra show a similar loss of available potential energy in the submesoscale range. This is consistent with the average negative values in CTRL over the entire range shown in Fig. 12a. In contrast, the NOTFB_SUB experiment shows a gain in potential energy, consistent with the positive values of $CS_{PE}$ in this case (Fig. 12c).

### c. Modulation of baroclinic energy conversion

We expect that modulation of potential energy by submesoscale air-sea coupling may modify SKE through baroclinic energy conversion. Therefore, this subsection investigates the effect of TFB and CFB on submesoscale baroclinic energy conversion $P_s K_s$, in order to assess its role in reducing SKE.

$P_s K_s$ is estimated over the first 100 m:

$$P_s K_s = \int_{-100}^{0} w' b', \quad (12)$$

where $w$ is the vertical velocity. Note that $w'b'$ reaches a maximum around 50-m depth and becomes negligible around 100 m (not shown). Figures 13a-c show $P_s K_s$ estimated from all the experiments, revealing the intensity and spatial distribution. Notably, in all experiments, energy conversion is prominently concentrated along the GS path, consistent with the spatial distribution of SKE. Smoothing out the submesoscale current feedback to the atmosphere in NOCFB_SUB shows minimal impact on $P_s K_s$, suggesting only a weak CFB effect on submesoscale baroclinic energy conversion. This is consistent with the low impact of CFB on submesoscale heat fluxes and the potential energy sink shown in the previous section (see Fig. 12). In contrast, neglecting the submesoscale TFB in NOTFB_SUB results in a significant increase in $w'b'$, averaging 15% over the whole domain, which supports the previous assessment of a decrease in SKE by TFB. This is consistent with the results of Uchoa et al. (2024, manuscript submitted to J. Phys. Oceanogr.) for the California upwelling system.

We have already noted (Fig. 4e) a strong dependence of SKE reduction by TFB on the wavelength considered, with the greatest reduction below 10 km. To confirm this spectral behavior on $P_s K_s$, following Uchoa et al. (2024, manuscript submitted to J. Phys. Oceanogr.), the cospectrum of the baroclinic energy conversion is estimated as

$$CSP_s K_s = \overline{\mathbb{R}[\hat{w}^* \hat{b}^*]} \quad (13)$$

Figures 13g and 13h show the results from all the experiments, as well as the percentage difference between NOCFB_SUB and CTRL (green) and NOTFB_SUB and CTRL (red). The striking similarity between CTRL and NOCFB_SUB cospectra confirms a minimal alteration of $P_s K_s$ by CFB, whatever the wavelength. In contrast, the difference between CTRL and NOTFB_SUB reveals a $CSP_s K_s$ signal reduced by TFB, particularly pronounced below 10 km, where the difference can reach 75%. This is in agreement with a potential energy sink attributed to TFB.

While the observed difference in $P_s K_s$ can be attributed to the change in $b$ (through heat fluxes), they may be modulated by a change in oceanic vertical velocities induced by the surface stress curl response to both TFB and CFB (Renault et al. 2018), i.e., a change in near-surface Ekman pumping. Following Renault et al. (2023b), to estimate an upper bound on the importance of submesoscale Ekman pumping, constant Ekman pumping is considered over the first 100 m of depth. This assumption is simplistic and has many caveats because it considers a constant Ekman layer depth (constant turbulent viscosity), where the vertical Ekman velocity is taken equal to the maximum value proposed by Stern (1965) and Wenegrat and Thomas (2017):

$$w'_{\text{Ekman}} = \frac{1}{\rho_0} \nabla \times \frac{\tau}{f + \zeta_b}, \quad (14)$$
<figure><img src="image_14.png" /><figcaption>FIG. 13. (a)–(c) Baroclinic energy conversion ($\int_{-100}^{0} w'b'$) from CTRL, NOCFB_SUB, and NOTFB_SUB. (d)–(f) As in (a)–(c), but for that induced by submesoscale Ekman pumping. Cospectrum of (g) $wb$ and (h) difference between CTRL and other experiments. The dashed line represents approximately the effective resolution of the oceanic simulations. Note that the nearshore values naturally are zero.</figcaption></figure>

where $f$ is the Coriolis parameter and $\zeta_b$ is the vorticity of the surface balanced flow. The conversion of energy by baroclinic instability $P_s K_{s\_OA\_Ekman}$ associated with $w_{Ekman}$ is then estimated as

$$P_s K_{s\_OA\_Ekman} = \int_{z=100}^{0} \overline{w'_{Ekman} b'}. \quad (15)$$

Figure 13d shows the result from CTRL. The term $P_s K_{s\_OA\_Ekman}$ is generally positive and thus contributes to an increase of SKE by $P_s K_s$, as opposed to the overall TFB effect. It is significant primarily over the GS path, i.e., where the surface stress curl response to submesoscale TFB and CFB is strongest. However, $P_s K_{s\_OA\_Ekman}$ is significantly weaker than $P_s K_s$, representing less than 3% of the total $P_s K_s$. Besides, this value is probably overestimated, given the abovementioned caveats. To unravel

the contributions of TFB and CFB in $P_s K_{s\_OA\_Ekman}$, it is further estimated from NOCFB_SUB and NOTFB_SUB (Figs. 13e,f). Neglecting submesoscale CFB in NOCFB_SUB leads to a drastic reduction in $P_s K_{s\_OA\_Ekman}$, suggesting that most of the additional baroclinic energy conversion induced by Ekman pumping is caused by CFB. This is confirmed by $P_s K_{s\_OA\_Ekman}$ estimated in NOTFB_SUB, which is very similar to CTRL. This suggests that the TFB-induced surface stress anomaly and associated Ekman pumping are only weakly spatially coherent with buoyancy anomaly. It does not cause any significant additional baroclinic energy conversion. This is consistent with the mesoscale behavior shown by Gaube et al. (2015) and Oerder et al. (2018). They show that Ekman pumping induced by TFB over an eddy has a dipole pattern that cancels its effect once averaged over the eddy.
<figure><img src="image_15.png" /><figcaption>FIG. 14. (a)–(b) Schematic representation of the (a) CFB effects and of the (b) TFB effects. The yellow area represents the GS region, where the submesoscale activity is the most intense. The black and green arrows are related to CFB and TFB, respectively.</figcaption></figure>

**6. Conclusions and discussion**

This study uses coupled submesoscale ocean–atmosphere simulations to explore the intricate influence of submesoscale thermal feedback (TFB) and current feedback (CFB) on various atmospheric and oceanic processes. Specifically, we investigate their impacts on low-level winds and surface stress, the exchange of heat, kinetic and potential energy between the ocean and atmosphere, and the submesoscale kinetic energy (SKE). Our results show that both submesoscale TFB and CFB lead to a damping of the SKE by about 10%, which is weaker than eddy damping at the mesoscale (about 30%). The difference can be explained by time-scale considerations and also by the interplay between submesoscale TFB and CFB, resulting in constructive or destructive effects on surface drag, depending on the specific characteristics of SST gradients and vorticity.

The coupling processes dominant at the mesoscale remain dominant at the submesoscale. CFB appears to be the main process for depleting kinetic energy from the ocean, resulting in negative wind work (Fig. 14a, black arrows). Surprisingly, our study tends to show not a reduction but an increase in the coupling relationship between currents and surface stress at submesoscale, due to a weaker wind response—the secondary process acting as negative feedback from the atmosphere, responsible for a partial re-energization of surface currents. However, the short lifetime of submesoscale eddies, compared to the frictional time scale associated with surface drag, limits the efficiency of the momentum coupling process, resulting in a relatively modest reduction in SKE, of the order of 10%. This reconciles the results from Renault et al. (2018) and Conejero et al. (2024); the weaker damping is explained by the short lifetime of submesoscale currents compared with mesoscale currents. Another difference for CFB at the submesoscale is its effect on flow divergence that was not prevalent at the mesoscale. The resulting wind divergence and associated vertical airflow can affect cloud cover.

On the other hand, TFB, by altering the low wind and the surface stress, has an energizing effect on ocean currents,

favoring positive wind work (Fig. 14b, black arrows), but this remains a secondary process that cannot counteract the dominant role of CFB in the exchange of kinetic energy between the ocean and the atmosphere. TFB is more significant for the potential energy budget, where it emerges as a key player in submesoscale air–sea exchanges. Specifically, the response of submesoscale SST to heat flux anomalies results in a reduction of surface ocean buoyancy and a sink of potential energy from the ocean to the atmosphere (Fig. 14b, green arrows). This is in agreement with Uchoa et al. (2024, manuscript submitted to J. Phys. Oceanogr.) that shows the presence of a similar potential energy sink for the California region. Interestingly, in the vicinity of the GS, we show that freshwater fluxes resulting from evaporation can cause a reverse transfer of potential energy to the ocean (Fig. 14b). However, the magnitude of this process seems second order compared to the loss of potential energy from the ocean to the atmosphere caused by heat fluxes over the GS. Overall, while CFB primarily affects kinetic energy exchange, TFB plays a key role in modulating potential energy through its effect on surface heat fluxes.

The loss of potential energy due to TFB has an effect on submesoscale kinetic energy through baroclinic energy conversion. The latter is weakened, particularly at scales below 10 km, where the decrease in SKE can be locally substantial (in spectral space), reaching 40% below 5 km. We note that CFB has an inverse effect on baroclinic energy conversion since it tends to increase it through Ekman pumping, but this effect remains secondary (Fig. 14). This leaves us with the central roles of submesoscale CFB and TFB in the kinetic and potential energy sinks, respectively, both leading to SKE damping.

While our study provides valuable insights into submesoscale processes, several limitations must be considered. The first concern arises from spatial resolution. The resolution used in our simulations for the ocean ($dx = 700$ m, corresponding to an effective resolution of about 3 km), while sufficient to capture many submesoscale features (fronts and mixed layer instabilities), may still miss or underrepresent
others (e.g., ageostrophic instabilities, internal waves). In addition, TFB appears to be most active at scales below 10 km, particularly at scales below 3 km, i.e., in a gray zone of our simulations where submesoscale features are only partially (or not at all) resolved. Future coupled simulations should further increase the horizontal resolution to verify the role of TFB at the kilometer scale. Another potential concern is the mismatch between oceanic and atmospheric spatial resolutions. Jullien et al. (2020), from simulations at the oceanic mesoscale resolution, show that an ocean-atmosphere resolution ratio of 3 does not affect the exchange of kinetic energy between the ocean and the atmosphere. This is mainly due to the effective resolution of ocean models, which is generally greater than 5dx (see, e.g., Soufflet et al. 2016). Consequently, when interpolating ocean fields onto the atmospheric grid, there is relatively little loss of structure if the resolution ratio is between 1 and 4. In our submesoscale study, we expect, for the same argument of effective resolution, that using a resolution ratio of 3 will have little impact on the results. Increasing the spatial resolution of the atmosphere could also improve the representation of coastal atmospheric features, but this should not affect the results presented here as we focus on offshore regions. Atmospheric fronts must also be well resolved, but we assume that a mesoscale resolution of the atmosphere is sufficient for our purposes. Second, our model does not account for surface gravity waves, which can directly affect submesoscale currents and indirectly the air-sea heat and momentum exchanges. Although our simulations provide valuable information on submesoscale dynamics, the short duration of the simulation limits our ability to fully assess how these processes affect the large-scale, long-term behavior of oceanic and atmospheric dynamics. Further research with long-term simulations is needed to assess these effects.

Our results highlight the shortcomings of current parameterizations for capturing CFB at submesoscale resolutions (Renault et al. 2020). These shortcomings are due to the magnitude of the wind response, the comparative time scales of the flow and coupling processes, and, to a lesser extent, the combined effects of CFB and TFB on momentum and buoyancy fluxes. New parameterizations should be inherently scale-aware, capable of adapting to the specific spatial scales under consideration to ensure accurate representation of air-sea interactions across a range of scales.

Finally, there is a pressing need for progress in the simultaneous measurement of surface currents, SST and wind, to better understand and quantify mesoscale and submesoscale oceanic processes and their interaction with the atmosphere. While existing satellite missions provide valuable data, their coverage and resolution are limited. The recent Surface Water and Ocean Topography (SWOT; Morrow et al. 2019) mission is promising in this regard, providing measurements of fine-scale geostrophic currents that provide a more detailed description of ocean circulation. However, it is essential to complement these surface current measurements with wind data to assess air-sea interactions. In addition to SWOT, other satellite missions such as ODYSEA (Rodríguez et al. 2019) can provide simultaneous wind and surface current measurements. The Harmony satellite, with its very high spatial resolution of

surface currents, SST and wind, also offers an opportunity, albeit with nonglobal coverage. Nevertheless, these advances in satellite technology hold great promise for improving our understanding of fine-scale processes and their impact on air-sea interactions.

*Acknowledgments.* L. R. appreciates support from CNES TOSCA I-CASCADE, M-ODYSEA, SWOT-POSEIDON, LEFE VENUS, and SHOM project 22CP06. I. U. and J. W. acknowledge support from NASA Grants 80NSSC21K0554 and 80NSSC24K0412. This work used the GENCI (project 13051) computing resources.

*Data availability statement.* The CROCO model and its toolbox are available at https://www.croco-ocean.org/. The WRF version we used can be downloaded from https://github.com/wrf-croco/. AVISO and OSTIA data are available on the Copernicus website at https://marine.copernicus.eu/, SCOW at https://chapman.ceoas.oregonstate.edu/scow/, and OAFlux at https://oaflux.whoi.edu.

## REFERENCES

Bai, Y., A. F. Thompson, A. B. Villas Bôas, P. Klein, H. S. Torres, and D. Menemenlis, 2023: Sub-mesoscale wind-front interactions: The combined impact of thermal and current feedback. *Geophys. Res. Lett.*, **50**, e2023GL104807, https://doi.org/10.1029/2023GL104807.

Bishop, S. P., R. J. Small, F. O. Bryan, and R. A. Tomas, 2017: Scale dependence of midlatitude air-sea interaction. *J. Climate*, **30**, 8207–8221, https://doi.org/10.1175/JCLI-D-17-0159.1.
—, —, and —, 2020: The global sink of available potential energy by mesoscale air-sea interaction. *J. Adv. Model. Earth Syst.*, **12**, e2020MS002118, https://doi.org/10.1029/2020MS002118.

Boccaletti, G., R. Ferrari, and B. Fox-Kemper, 2007: Mixed layer instabilities and restratification. *J. Phys. Oceanogr.*, **37**, 2228–2250, https://doi.org/10.1175/JPO3101.1.

Bye, J. A. T., 1985: Large-scale momentum exchange in the coupled atmosphere-ocean. *Coupled Ocean-Atmosphere Models*, J. C. J. Nihoul, Ed., Elsevier Oceanography Series, Vol. 40, Elsevier, 51–61.

Callies, J., R. Ferrari, J. M. Klymak, and J. Gula, 2015: Seasonality in submesoscale turbulence. *Nat. Commun.*, **6**, 6862, https://doi.org/10.1038/ncomms7862.

Capet, X., J. C. McWilliams, M. J. Molemaker, and A. F. Shchepetkin, 2008a: Mesoscale to submesoscale transition in the California current system. Part II: Frontal processes. *J. Phys. Oceanogr.*, **38**, 44–64, https://doi.org/10.1175/2007JPO3672.1.
—, —, —, and —, 2008b: Mesoscale to submesoscale transition in the California current system. Part III: Energy balance and flux. *J. Phys. Oceanogr.*, **38**, 2256–2269, https://doi.org/10.1175/2008JPO3810.1.

—, —, —, and —, 2008c: Mesoscale to submesoscale transition in the California current system. Part I: Flow structure, eddy flux, and observational tests. *J. Phys. Oceanogr.*, **38**, 29–43, https://doi.org/10.1175/2007JPO3671.1.

Chassignet, E. P., and D. P. Marshall, 2008: Gulf Stream separation in numerical ocean models. *Ocean Modeling in an*
Eddying Regime, Geophys. Monogr., Vol. 177, Amer. Geo-phys. Union, 39–61, https://doi.org/10.1029/177GM05.

Chelton, D. B., and Coauthors, 2001: Observations of coupling between surface wind stress and sea surface temperature in the eastern tropical Pacific. *J. Climate*, **14**, 1479–1498, https://doi.org/10.1175/1520-0442(2001)014<1479:OOCBSW>2.0.CO;2.

—, M. G. Schlax, M. H. Freilich, and R. F. Milliff, 2004: Satellite measurements reveal persistent small-scale features in ocean winds. *Science*, **303**, 978–983, https://doi.org/10.1126/science.1091901.

—, —, and R. M. Samelson, 2007: Summertime coupling between sea surface temperature and wind stress in the California Current System. *J. Phys. Oceanogr.*, **37**, 495–517, https://doi.org/10.1175/JPO3025.1.

—, —, and —, 2011: Global observations of nonlinear mesoscale eddies. *Prog. Oceanogr.*, **91**, 167–216, https://doi.org/10.1016/j.pocean.2011.01.002.

Conejero, C., L. Renault, F. Desbiolles, J. C. McWilliams, and H. Giordani, 2024: Near-surface atmospheric response to meso- and submesoscale current and thermal feedbacks. *J. Phys. Oceanogr.*, **54**, 823–848, https://doi.org/10.1175/JPO-D-23-0211.1.

Contreras, M., 2023: Study of the kinetic energy cycle and dissipation pathways in the Gulf Stream. Ph.D. thesis, Université Paul Sabatier, 193 pp.

—, L. Renault, and P. Marchesiello, 2023a: Tidal modulation of energy dissipation routes in the Gulf Stream. *Geophys. Res. Lett.*, **50**, e2023GL104946, https://doi.org/10.1029/2023GL104946.

—, —, and —, 2023b: Understanding energy pathways in the Gulf Stream. *J. Phys. Oceanogr.*, **53**, 719–736, https://doi.org/10.1175/JPO-D-22-0146.1.

Craig, A., S. Valcke, and L. Coquart, 2017: Development and performance of a new version of the OASIS coupler, OASIS3-MCT_3.0. *Geosci. Model Dev.*, **10**, 3297–3308, https://doi.org/10.5194/gmd-10-3297-2017.

Cronin, M., and J. Sprintall, 2009: Wind-and buoyancy-forced upper ocean. *Elements of Physical Oceanography: A Derivative of the Encyclopedia of Ocean Sciences*, J. H. Steel, Ed., Academic Press, 3217–3224.

Czaja, A., C. Frankignoul, S. Minobe, and B. Vannière, 2019: Simulating the midlatitude atmospheric circulation: What might we gain from high-resolution modeling of air-sea interactions? *Curr. Climate Change Rep.*, **5**, 390–406, https://doi.org/10.1007/s40641-019-00148-5.

Debreu, L., P. Marchesiello, P. Penven, and G. Cambon, 2012: Two-way nesting in split-explicit ocean models: Algorithms, implementation and validation. *Ocean Modell.*, **49–50**, 1–21, https://doi.org/10.1016/j.ocemod.2012.03.003.

Desbiolles, F., A. N. Meroni, L. Renault, and C. Pasquero, 2023: Environmental control of wind response to sea surface temperature patterns in reanalysis dataset. *J. Climate*, **36**, 3881–3893, https://doi.org/10.1175/JCLI-D-22-0373.1.

Dewar, W. K., and G. R. Flierl, 1987: Some effects of the wind on rings. *J. Phys. Oceanogr.*, **17**, 1653–1667, https://doi.org/10.1175/1520-0485(1987)017<1653:SEOTWO>2.0.CO;2.

Donlon, C. J., M. Martin, J. Stark, J. Roberts-Jones, E. Fiedler, and W. Wimmer, 2012: The Operational Sea Surface Temperature and Sea Ice Analysis (Ostia) system. *Remote Sens. Environ.*, **116**, 140–158, https://doi.org/10.1016/j.rse.2010.10.017.

Dudhia, J., 1989: Numerical study of convection observed during the winter monsoon experiment using a mesoscale two-

dimensional model. *J. Atmos. Sci.*, **46**, 3077–3107, https://doi.org/10.1175/1520-0469(1989)046<3077:NSOCOD>2.0.CO;2.

Duhaut, T. H. A., and D. N. Straub, 2006: Wind stress dependence on ocean surface velocity: Implications for mechanical energy input to ocean circulation. *J. Phys. Oceanogr.*, **36**, 202–211, https://doi.org/10.1175/JPO2842.1.

Eden, C., and H. Dietze, 2009: Effects of mesoscale eddy/wind interactions on biological new production and eddy kinetic energy. *J. Geophys. Res.*, **114**, C05023, https://doi.org/10.1029/2008JC005129.

Fairall, C. W., E. F. Bradley, J. Hare, A. A. Grachev, and J. B. Edson, 2003: Bulk parameterization of air-sea fluxes: Updates and verification for the COARE algorithm. *J. Climate*, **16**, 571–591, https://doi.org/10.1175/1520-0442(2003)016<0571: BPOASF>2.0.CO;2.

Frenger, I., N. Gruber, R. Knutti, and M. Münnich, 2013: Imprint of Southern Ocean eddies on winds, clouds and rainfall. *Nat. Geosci.*, **6**, 608–612, https://doi.org/10.1038/ngeo1863.

Gaube, P., D. B. Chelton, R. M. Samelson, M. G. Schlax, and L. W. O'Neill, 2015: Satellite observations of mesoscale eddy-induced Ekman pumping. *J. Phys. Oceanogr.*, **45**, 104–132, https://doi.org/10.1175/JPO-D-14-0032.1.

Han, J., and H.-L. Pan, 2011: Revision of convection and vertical diffusion schemes in the NCEP Global Forecast System. *Wea. Forecasting*, **26**, 520–533, https://doi.org/10.1175/WAF-D-10-05038.1.

Harrison, C. S., M. C. Long, N. S. Lovenduski, and J. K. Moore, 2018: Mesoscale effects on carbon export: A global perspective. *Global Biogeochem. Cycles*, **32**, 680–703, https://doi.org/10.1002/2017GB005751.

Heiderich, J., and R. E. Todd, 2020: Along-stream evolution of Gulf Stream volume transport. *J. Phys. Oceanogr.*, **50**, 2251–2270, https://doi.org/10.1175/JPO-D-19-0303.1.

Holmes, R. M., L. Renault, L. Maillard, and J. Boucharel, 2024: Air-sea coupling feedbacks over tropical instability waves. *J. Phys. Oceanogr.*, **54**, 2165–2184, https://doi.org/10.1175/JPO-D-24-0010.1.

Hong, S.-Y., and J.-O. J. Lim, 2006: The WRF Single-Moment 6-class Microphysics scheme (WSM6). *Asia-Pac. J. Atmos. Sci.*, **42**, 129–151.

—, Y. Noh, and J. Dudhia, 2006: A new vertical diffusion package with an explicit treatment of entrainment processes. *Mon. Wea. Rev.*, **134**, 2318–2341, https://doi.org/10.1175/MWR3199.1.

Jiménez, P. A., J. Dudhia, J. F. González-Rouco, J. Navarro, J. P. Montávez, and E. García-Bustamante, 2012: A revised scheme for the WRF surface layer formulation. *Mon. Wea. Rev.*, **140**, 898–918, https://doi.org/10.1175/MWR-D-11-00056.1.

Jin, X., and R. A. Weller, 2008: Multidecade global flux datasets from the objectively analyzed air-sea fluxes (OAFlux) project: Latent and sensible heat fluxes, ocean evaporation, and related surface meteorological variables. *Lisan Yu. OAFlux Project Tech. Rep. OA-2008-01*, 64 pp., https://apdrc.soest.hawaii.edu/doc/OAFlux_TechReport_3rd_release.pdf.

Jousse, A., A. Hall, F. Sun, and J. Teixeira, 2016: Causes of WRF surface energy fluxes biases in a stratocumulus region. *Climate Dyn.*, **46**, 571–584, https://doi.org/10.1007/s00382-015-2599-9.

Jullien, S., S. Masson, V. Oerder, G. Samson, F. Colas, and L. Renault, 2020: Impact of ocean-atmosphere current feedback on ocean mesoscale activity: Regional variations and sensitivity to model resolution. *J. Climate*, **33**, 2585–2602, https://doi.org/10.1175/JCLI-D-19-0484.1.
Kessouri, F., D. Bianchi, L. Renault, J. C. McWilliams, H. Frenzel, and C. A. Deutsch, 2020: Submesoscale currents modulate the seasonal cycle of nutrients and productivity in the California Current System. *Global Biogeochem. Cycles*, **34**, e2020GB006578, https://doi.org/10.1029/2020GB006578.

—, L. Renault, J. C. McWilliams, P. Damien, and D. Bianchi, 2022: Enhancement of oceanic eddy activity by fine-scale orographic winds drives high productivity, low oxygen, and low pH conditions in the Santa Barbara channel. *J. Geophys. Res. Oceans*, **127**, e2022JC018947, https://doi.org/10.1029/2022JC018947.

Kundu, P. K., I. M. Cohen, and D. R. Dowling, 2015: *Fluid Mechanics*. Academic Press, 928 pp.

Kwon, Y. C., and S.-Y. Hong, 2017: A mass-flux cumulus parameterization scheme across gray-zone resolutions. *Mon. Wea. Rev.*, **145**, 583–598, https://doi.org/10.1175/MWR-D-16-0034.1.

Large, W. G., J. C. McWilliams, and S. C. Doney, 1994: Oceanic vertical mixing: A review and a model with a nonlocal boundary layer parameterization. *Rev. Geophys.*, **32**, 363–403, https://doi.org/10.1029/94RG01872.

Larrañaga, M., L. Renault, and J. Jouanno, 2022: Partial control of the Gulf of Mexico dynamics by the current feedback to the atmosphere. *J. Phys. Oceanogr.*, **52**, 2515–2530, https://doi.org/10.1175/JPO-D-21-0271.1.

Lemarié, F., J. Kurian, A. F. Shchepetkin, M. J. Molemaker, F. Colas, and J. C. McWilliams, 2012: Are there inescapable issues prohibiting the use of terrain-following coordinates in climate models? *Ocean Modell.*, **42**, 57–79, https://doi.org/10.1016/j.ocemod.2011.11.007.

Lévy, M., P. Klein, and M. Ben Jelloul, 2009: New production stimulated by high-frequency winds in a turbulent mesoscale eddy field. *Geophys. Res. Lett.*, **36**, L16603, https://doi.org/10.1029/2009GL039490.

—, P. J. Franks, and K. S. Smith, 2018: The role of submesoscale currents in structuring marine ecosystems. *Nat. Commun.*, **9**, 4758, https://doi.org/10.1038/s41467-018-07059-3.

Li, G., 1985: Robust regression. *Exploring Data Tables, Trends, and Shapes*, D. C. Hoalin, F. M. Mosteller, and J. W. Tukey, Eds., John Wiley and Sons, 281–340.

Luo, J.-J., S. Masson, E. Roeckner, G. Madec, and T. Yamagata, 2005: Reducing climatology bias in an ocean-atmosphere CGCM with improved coupling physics. *J. Climate*, **18**, 2344–2360, https://doi.org/10.1175/JCLI3404.1.

Ma, X., and Coauthors, 2016: Western boundary currents regulated by interaction between ocean eddies and the atmosphere. *Nature*, **535**, 533–537, https://doi.org/10.1038/nature18640.

Mahadevan, A., 2016: The impact of submesoscale physics on primary productivity of Plankton. *Annu. Rev. Mar. Sci.*, **8**, 161–184, https://doi.org/10.1146/annurev-marine-010814-015912.

—, and A. Tandon, 2006: An analysis of mechanisms for submesoscale vertical motion at ocean fronts. *Ocean Modell.*, **14**, 241–256, https://doi.org/10.1016/j.ocemod.2006.05.006.

Marchesiello, P., J. C. McWilliams, and A. Shchepetkin, 2001: Open boundary conditions for long-term integration of regional oceanic models. *Ocean Modell.*, **3** (1–2), 1–20, https://doi.org/10.1016/S1463-5003(00)00013-5.

—, L. Debreu, and X. Couvelard, 2009: Spurious diapycnal mixing in terrain-following coordinate models: The problem and a solution. *Ocean Modell.*, **26**, 156–169, https://doi.org/10.1016/j.ocemod.2008.09.004.

Masunaga, R., H. Nakamura, B. Taguchi, and T. Miyasaka, 2020: Processes shaping the frontal-scale time-mean surface wind convergence patterns around the Kuroshio extension in

winter. *J. Climate*, **33**, 3–25, https://doi.org/10.1175/JCLI-D-19-0097.1.

McGillicuddy, D. J., and Coauthors, 2007: Eddy/wind interactions stimulate extraordinary mid-ocean Plankton blooms. *Science*, **316**, 1021–1026, https://doi.org/10.1126/science.1136256.

McWilliams, J. C., 2008: The nature and consequences of oceanic eddies. *Ocean Modeling in an Eddying Regime, Geophys. Monogr.*, Vol. 177, Amer. Geophys. Union, 5–15, https://doi.org/10.1029/177GM03.

—, 2016: Submesoscale currents in the ocean. *Proc. Roy. Soc., 472A*, 20160117, https://doi.org/10.1098/rspa.2016.0117.

—, 2021: Oceanic frontogenesis. *Annu. Rev. Mar. Sci.*, **13**, 227–253, https://doi.org/10.1146/annurev-marine-032320-120725.

—, M. Molemaker, and I. Yavneh, 2001: From stirring to mixing of momentum: Cascades from balanced flows to dissipation in the oceanic interior. *Aha Huliko'a Proc. 2001*, Honolulu, HI, University of Hawaii, 59–66, https://apps.dtic.mil/sti/pdfs/ADA412459.pdf#page=62.

Ménesguen, C., S. Le Gentil, P. Marchesiello, and N. Ducousso, 2018: Destabilization of an oceanic meddy-like vortex: Energy transfers and significance of numerical settings. *J. Phys. Oceanogr.*, **48**, 1151–1168, https://doi.org/10.1175/jpo-d-17-0126.1.

Mlawer, E. J., S. J. Taubman, P. D. Brown, M. J. Iacono, and S. A. Clough, 1997: Radiative transfer for inhomogeneous atmospheres: RRTM, a validated correlated-k model for the longwave. *J. Geophys. Res.*, **102**, 16663–16682, https://doi.org/10.1029/97JD00237.

Moreton, S., D. Ferreira, M. Roberts, and H. Hewitt, 2021: Air-sea turbulent heat flux feedback over mesoscale eddies. *Geophys. Res. Lett.*, **48**, e2021GL095407, https://doi.org/10.1029/2021GL095407.

Morrow, R., and Coauthors, 2019: Global observations of fine-scale ocean surface topography with the Surface Water and Ocean Topography (SWOT) mission. *Front. Mar. Sci.*, **6**, 232, https://doi.org/10.3389/fmars.2019.00232.

Oerder, V., F. Colas, V. Echevin, S. Masson, and F. Lemarié, 2018: Impacts of the mesoscale ocean-atmosphere coupling on the Peru-Chile ocean dynamics: The current-induced wind stress modulation. *J. Geophys. Res. Oceans*, **123**, 812–833, https://doi.org/10.1002/2017JC013294.

Omand, M. M., E. A. D'Asaro, C. M. Lee, M. J. Perry, N. Briggs, I. Cetinić, and A. Mahadevan, 2015: Eddy-driven subduction exports particulate organic carbon from the spring bloom. *Science*, **348**, 222–225, https://doi.org/10.1126/science.1260062.

O'Neill, L. W., D. B. Chelton, and S. K. Esbensen, 2012: Covariance of surface wind and stress responses to sea surface temperature fronts. *J. Climate*, **25**, 5916–5942, https://doi.org/10.1175/JCLI-D-11-00230.1.

Renault, L., C. Deutsch, J. C. McWilliams, H. Frenzel, J.-H. Liang, and F. Colas, 2016a: Partial decoupling of primary productivity from upwelling in the California Current system. *Nat. Geosci.*, **9**, 505–508, https://doi.org/10.1038/ngeo2722.

—, A. Hall, and J. C. McWilliams, 2016b: Orographic shaping of US West Coast wind profiles during the upwelling season. *Climate Dyn.*, **46**, 273–289, https://doi.org/10.1007/s00382-015-2583-4.

—, M. J. Molemaker, J. Gula, S. Masson, and J. C. McWilliams, 2016c: Control and stabilization of the Gulf Stream by oceanic current interaction with the atmosphere. *J. Phys. Oceanogr.*, **46**, 3439–3453, https://doi.org/10.1175/JPO-D-16-0115.1.

—, —, J. C. McWilliams, A. F. Shchepetkin, F. Lemarié, D. Chelton, S. Illig, and A. Hall, 2016d: Modulation of wind
work by oceanic current interaction with the atmosphere. *J. Phys. Oceanogr.*, **46**, 1685-1704, https://doi.org/10.1175/JPO-D-15-0232.1.

—, J. C. McWilliams, and S. Masson, 2017a: Satellite observations of imprint of oceanic current on wind stress by air-sea coupling. *Sci. Rep.*, **7**, 17747, https://doi.org/10.1038/s41598-017-17939-1.

—, —, and P. Penven, 2017b: Modulation of the Agulhas Current retroflection and leakage by oceanic current interaction with the atmosphere in coupled simulations. *J. Phys. Oceanogr.*, **47**, 2077-2100, https://doi.org/10.1175/JPO-D-16-0168.1.

—, —, and J. Gula, 2018: Dampening of submesoscale currents by air-sea stress coupling in the Californian upwelling system. *Sci. Rep.*, **8**, 13388, https://doi.org/10.1038/s41598-018-31602-3.

—, F. Lemarié, and T. Arsouze, 2019a: On the implementation and consequences of the oceanic currents feedback in ocean-atmosphere coupled models. *Ocean Modell.*, **141**, 101423, https://doi.org/10.1016/j.ocemod.2019.101423.

—, P. Marchesiello, S. Masson, and J. C. McWilliams, 2019b: Remarkable control of western boundary currents by eddy killing, a mechanical air-sea coupling process. *Geophys. Res. Lett.*, **46**, 2743-2751, https://doi.org/10.1029/2018GL081211.

—, S. Masson, V. Oerder, S. Jullien, and F. Colas, 2019c: Disentangling the mesoscale ocean-atmosphere interactions. *J. Geophys. Res. Oceans*, **124**, 2164-2178, https://doi.org/10.1029/2018JC014628.

—, —, T. Arsouze, G. Madec, and J. C. McWilliams, 2020: Recipes for how to force oceanic model dynamics. *J. Adv. Model. Earth Syst.*, **12**, e2019MS001715, https://doi.org/10.1029/2019MS001715.

—, P. Marchesiello, and M. Contreras, 2023a: Coaction of top and bottom drags in Gulf Stream dynamics. *J. Geophys. Res. Oceans*, **128**, e2022JC018939, https://doi.org/10.1029/2022JC018939.

—, S. Masson, V. Oerder, F. Colas, and J. McWilliams, 2023b: Modulation of the oceanic mesoscale activity by the mesoscale thermal feedback to the atmosphere. *J. Phys. Oceanogr.*, **53**, 1651-1667, https://doi.org/10.1175/JPO-D-22-0256.1.

Risien, C. M., and D. B. Chelton, 2008: A global climatology of surface wind and wind stress fields from eight years of QuikSCAT scatterometer data. *J. Phys. Oceanogr.*, **38**, 2379-2413, https://doi.org/10.1175/2008JPO3881.1.

Rodríguez, E., M. Bourassa, D. Chelton, J. T. Farrar, D. Long, D. Perkovic-Martin, and R. Samelson, 2019: The winds and currents mission concept. *Front. Mar. Sci.*, **6**, 438, https://doi.org/10.3389/fmars.2019.00438.

Sandwell, D. T., and W. H. Smith, 1997: Marine gravity anomaly from Geosat and ERS 1 satellite altimetry. *J. Geophys. Res.*, **102**, 10 039-10 054, https://doi.org/10.1029/96JB03223.

Schubert, R., J. Gula, R. J. Greatbatch, B. Baschek, and A. Biastoch, 2020: The submesoscale kinetic energy cascade: Mesoscale absorption of submesoscale mixed layer eddies and frontal downscale fluxes. *J. Phys. Oceanogr.*, **50**, 2573-2589, https://doi.org/10.1175/JPO-D-19-0311.1.

Seo, H., A. J. Miller, and J. R. Norris, 2016: Eddy-wind interaction in the California Current System: Dynamics and impacts. *J. Phys. Oceanogr.*, **46**, 439-459, https://doi.org/10.1175/JPO-D-15-0086.1.

—, and Coauthors, 2023: Ocean mesoscale and frontal-scale ocean-atmosphere interactions and influence on large-scale climate: A review. *J. Climate*, **36**, 1981-2013, https://doi.org/10.1175/JCLI-D-21-0982.1.

Shan, X., and Coauthors, 2020: Surface heat flux induced by mesoscale eddies cools the Kuroshio-Oyashio extension region. *Geophys. Res. Lett.*, **47**, e2019GL086050, https://doi.org/10.1029/2019GL086050.

Shechepetkin, A., and J. McWilliams, 2003: A method for computing horizontal pressure-gradient force in an oceanic model with a nonaligned vertical coordinate. *J. Geophys. Res.*, **108**, 3090, https://doi.org/10.1029/2001JC001047.

—, —, 2005: The regional oceanic modeling system (ROMS): A split-explicit, free-surface, topography-following-coordinate oceanic model. *Ocean Modell.*, **9**, 347-404, https://doi.org/10.1016/j.ocemod.2004.08.002.

Shcherbina, A. Y., E. A. D'Asaro, C. M. Lee, J. M. Klymak, M. J. Molemaker, and J. C. McWilliams, 2013: Statistics of vertical vorticity, divergence, and strain in a developed submesoscale turbulence field. *Geophys. Res. Lett.*, **40**, 4706-4711, https://doi.org/10.1002/grl.50919.

Siegelman, L., P. Klein, P. Rivière, A. F. Thompson, H. S. Torres, M. Flexas, and D. Menemenlis, 2020: Enhanced upward heat transport at deep submesoscale ocean fronts. *Nat. Geosci.*, **13**, 50-55, https://doi.org/10.1038/s41561-019-0489-1.

Skamarock, W. C., and Coauthors, 2008: A description of the Advanced Research WRF version 3. NCAR Tech. Note NCAR/TN-475+STR, 113 pp., https://doi.org/10.5065/D68S4MVH.

Skyllingstad, E. D., D. Vickers, L. Mahrt, and R. Samelson, 2007: Effects of mesoscale sea-surface temperature fronts on the marine atmospheric boundary layer. *Bound.-Layer Meteor.*, **123**, 219-237, https://doi.org/10.1007/s10546-006-9127-8.

Small, R. J., and Coauthors, 2008: Air-sea interaction over ocean fronts and eddies. *Dyn. Atmos. Oceans*, **45**, 274-319, https://doi.org/10.1016/j.dynatmoce.2008.01.001.

Soufflet, Y., P. Marchesiello, F. Lemarié, J. Jouanno, X. Capet, L. Brebu, and R. Benshila, 2016: On effective resolution in ocean models. *Ocean Modell.*, **9**, 36-50, https://doi.org/10.1016/j.ocemod.2015.12.004.

Srinivasan, K., J. C. McWilliams, L. Renault, H. G. Hristova, J. Molemaker, and W. S. Kessler, 2017: Topographic and mixed layer submesoscale currents in the near-surface southwestern tropical Pacific. *J. Phys. Oceanogr.*, **47**, 1221-1242, https://doi.org/10.1175/JPO-D-16-0216.1.

Stern, M. E., 1965: Interaction of a uniform wind stress with a geostrophic vortex. *Deep-Sea Res. Oceanogr. Abstr.*, **12**, 355-367, https://doi.org/10.1016/0011-7471(65)90007-0.

Su, Z., J. Wang, P. Klein, A. F. Thompson, and D. Menemenlis, 2018: Ocean submesoscales as a key component of the global heat budget. *Nat. Commun.*, **9**, 775, https://doi.org/10.1038/s41467-018-02983-w.

—, H. Torres, P. Klein, A. F. Thompson, L. Siegelman, J. Wang, D. Menemenlis, and C. Hill, 2020: High-frequency submesoscale motions enhance the upward vertical heat transport in the global ocean. *J. Geophys. Res. Oceans*, **125**, e2020JC016544, https://doi.org/10.1029/2020JC016544.

Taburet, G., A. Sanchez-Roman, M. Ballarotta, M.-I. Pujol, J.-F. Legeais, F. Fournier, Y. Faugere, and G. Dibarboure, 2019: DUACS DT2018: 25 years of reprocessed sea level altimetry products. *Ocean Sci.*, **15**, 1207-1224, https://doi.org/10.5194/os-15-1207-2019.

Tewari, M., and Coauthors, 2004: Implementation and verification of the unified NOAH Land Surface Model in the WRF model (formerly paper number 17.5). *Proc. 20th Conf. on Weather Analysis and Forecasting/16th Conf. on Numerical Weather Prediction*, Vol. 14, Seattle, WA, Amer. Meteor. Soc., 14.2a, https://ams.confex.com/ams/84Annual/techprogram/paper_69061.htm.
Thomas, L. N., A. Tandon, and A. Mahadevan, 2008: Submeso-scale processes and dynamics. *Ocean Modeling Eddying Regime, Geophys. Monogr.*, Vol. 177, Amer. Geophys. Union, 17–38, https://doi.org/10.1029/177GM04.

—, J. R. Taylor, R. Ferrari, and T. M. Joyce, 2013: Symmetric instability in the Gulf Stream. *Deep-Sea Res. II*, **91**, 96–110, https://doi.org/10.1016/j.dsr2.2013.02.025.

Tiedtke, M., 1989: A comprehensive mass flux scheme for cumulus parameterization in large-scale models. *Mon. Wea. Rev.*, **117**, 1779–1800, https://doi.org/10.1175/1520-0493(1989)117<1779:ACMFSF>2.0.CO;2.

von Storch, J.-S., C. Eden, I. Fast, H. Haak, D. Hernández-Deckers, E. Maier-Reimer, J. Marotzke, and D. Stammer, 2012: An estimate of the Lorenz energy cycle for the World Ocean based on the STORM/NCEP simulation. *J. Phys. Oceanogr.*, **42**, 2185–2205, https://doi.org/10.1175/JPO-D-12-079.1.

Wenegrat, J. O., 2023: The current feedback on stress modifies the Ekman buoyancy flux at fronts. *J. Phys. Oceanogr.*, **53**, 2737–2749, https://doi.org/10.1175/JPO-D-23-0005.1.

—, and L. N. Thomas, 2017: Ekman transport in balanced currents with curvature. *J. Phys. Oceanogr.*, **47**, 1189–1203, https://doi.org/10.1175/JPO-D-16-0239.1.

—, and R. S. Arthur. 2018: Response of the atmospheric boundary layer to submesoscale sea surface temperature fronts. *Geophys. Res. Lett.*, **45**, 13505–13512, https://doi.org/10.1029/2018GL081034.

—, J. Callies, and L. N. Thomas, 2018a: Submesoscale baroclinic instability in the bottom boundary layer. *J. Phys. Oceanogr.*, **48**, 2571–2592, https://doi.org/10.1175/JPO-D-17-0264.1.

—, L. N. Thomas, J. Gula, and J. C. McWilliams, 2018b: Effects of the submesoscale on the potential vorticity budget of ocean mode waters. *J. Phys. Oceanogr.*, **48**, 2141–2165, https://doi.org/10.1175/JPO-D-17-0219.1.

Wunsch, C., and R. Ferrari, 2004: Vertical mixing, energy, and the general circulation of the oceans. *Annu. Rev. Fluid Mech.*, **36**, 281–314, https://doi.org/10.1146/annurev.fluid.36.050802.122121.

Xu, Y., and R. B. Scott, 2008: Subtleties in forcing eddy resolving ocean models with satellite wind data. *Ocean Modell.*, **20**, 240–251, https://doi.org/10.1016/j.ocemod.2007.09.003.

Zhang, C., Y. Wang, and K. Hamilton, 2011: Improved representation of boundary layer clouds over the Southeast Pacific in ARW-WRF using a modified Tiedtke cumulus parameterization scheme. *Mon. Wea. Rev.*, **139**, 3489–3513, https://doi.org/10.1175/MWR-D-10-05091.1.