The ocean absorbs >90% of anthropogenic heat in the Earth system, moderating global atmospheric warming. However, it remains unclear how this heat uptake is distributed by basin and across water masses. Here we analyze historical and recent observations to show that ocean heat uptake has accelerated dramatically since the 1990s, nearly doubling during 2010–2020 relative to 1990–2000. Of the total ocean heat uptake over the Argo era 2005–2020, about 89% can be found in global mode and intermediate water layers, spanning both hemispheres and both subtropical and subpolar mode waters. Due to anthropogenic warming, there are significant changes in the volume of these water-mass layers as they warm and freshen. After factoring out volumetric changes, the combined warming of these layers accounts for ~76% of global ocean warming. We further decompose these water-mass layers into regional water masses over the subtropical Pacific and Atlantic Oceans and in the Southern Ocean. This shows that regional mode and intermediate waters are responsible for a disproportionate fraction of total heat uptake compared to their volume, with important implications for understanding ongoing ocean warming, sea-level rise, and climate impacts.

Mixing by mesoscale eddies profoundly impacts climate and ecosystems by redistributing and storing dissolved tracers such as heat and carbon. Eddy mixing is parameterized in most numerical models of the ocean and climate. To reduce known sensitivity to such parameterizations, observational estimates of mixing are needed. However, logistical and technological limitations obstruct our ability to measure global time‐varying mixing rates. Here, we extend mixing length theory with mean‐flow suppression theory, and first surface modes, to estimate mixing from readily available observational‐based climatological data, of salinity, temperature, pressure, and eddy kinetic energy at the sea surface. The resulting full‐depth global maps of eddy mixing can reproduce the few available direct estimates and confirm the importance of mean‐flow suppression of mixing. The results also emphasize the significant effect of eddy surface intensification and its relation to the vertical density stratification. These new insights in mixing dynamics will improve future mesoscale eddy mixing parameterizations.

This study provides observation-based estimates, determined by inverse methods, of horizontal and isopycnal eddy diffusion coefficients K_H and K_I, respectively, the small-scale mixing coefficient D, and the diathermohaline streamfunction Ψ. The inverse solution of Ψ represents the ocean circulation in Absolute Salinity S_A and Conservative Temperature Θ coordinates. The authors suggest that the observation-based estimate of Ψ will be useful for comparison with equivalent diagnostics from numerical climate models. The estimates of K_H and K_I represent horizontal eddy mixing in the mixed layer and isopycnal eddy mixing in the ocean interior, respectively. This study finds that the solution for D and K_H are comparable to existing estimates. The solution for K_I is one of the first observation-based global and full-depth constrained estimates of isopycnal mixing and indicates that K_I is an order of magnitude smaller than K_H. This suggests that there is a large vertical variation in the eddy mixing coefficient, which is generally not included in ocean models. With ocean models being very sensitive to the choice of isopycnal mixing, this result suggests that further investigation of the spatial structure of isopycnal eddy mixing from observations is required.

Mesoscale eddies stir along the neutral plane, and the resulting neutral diffusion is a fundamental aspect of subgrid‐scale tracer transport in ocean models. Calculating neutral diffusion traditionally involves calculating neutral slopes and three‐dimensional tracer gradients. The calculation of the neutral slope traditionally occurs by computing the ratio of the horizontal to vertical locally referenced potential density derivative. However, this approach is problematic in regions of weak vertical stratification, prompting the use of a variety of ad hoc regularization methods that can lead to rather nonphysical dependencies for the resulting neutral tracer gradients. Here we use a VErtical Non‐local Method “VENM,” a search algorithm that requires no ad hoc regularization and significantly improves the numerical accuracy of calculating neutral slopes, neutral tracer gradients, and associated neutral diffusive fluxes. We compare and contrast VENM against a more traditional method, using an independent objective neutrality condition combined with estimates of spurious diffusion, heat transport, and water mass transformation rates. VENM is more accurate, both physically and numerically, and should form the basis for future efforts involving neutral diffusion calculations from observations and possibly numerical model simulations.

Subantarctic Mode Water (SAMW) forms in deep mixed layers just north of the Antarctic Circumpolar Current in winter, playing a fundamental role in the ocean uptake of heat and carbon. Using a gridded Argo product and the ERA-Interim reanalysis for years 2004-2018, the seasonal evolution of the SAMW volume is analyzed using both a kinematic estimate of the subduction rate and a thermodynamic estimate of the air-sea formation rate. The seasonal SAMW volume changes are separately estimated within the monthly mixed layer and in the interior below it. We find that the variability of SAMW volume is dominated by changes in SAMW volume in the mixed layer. The seasonal variability of SAMW volume in the mixed layer is governed by formation due to air-sea buoyancy fluxes (45%, lasting from July to August), entrainment (35%), and northward Ekman transport across the Subantarctic Front (10%). The interior SAMW formation is entirely controlled by exchanges between the mixed layer and the interior (i.e. instantaneous subduction), which occurs mainly during August-October. The annual mean subduction estimate from a Lagrangian approach shows strong regional variability with hotspots of large SAMW subduction. The SAMW subduction hotspots are consistent with the distribution and export pathways of SAMW over the central and eastern parts of the south Indian and Pacific Oceans. Hotspots in the south Indian Ocean produce strong subduction of 8 and 9 Sv for the light and southeast Indian SAMW, respectively, while SAMW subduction of 6 and 4 Sv occurs for the central and southeast Pacific SAMW, respectively.