# Deep Ocean Ventilation: A Comparison Between a General Circulation Model and Data-Constrained Inverse Models

Bruno Millet, Casimir de Lavergne, William Robert Gray, Christian Éthé, Gurvan Madec, Mark Holzer, Tim Devries, Geoffrey Gebbie, Didier M Roche

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

► To cite this version:

Bruno Millet, Casimir de Lavergne, William Robert Gray, Christian Éthé, Gurvan Madec, et al., Deep Ocean Ventilation: A Comparison Between a General Circulation Model and Data-Constrained Inverse Models. Journal of Advances in Modeling Earth Systems, 2025, 17 (7), 10.1029/2024ms004914. hal-05173206

HAL Id: hal-05173206

https://hal.science/hal-05173206v1

Submitted on 21 Jul 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_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
RESEARCH ARTICLE

10.1029/2024MS004914

**Special Collection:**

Advances in large-scale
hydrological modeling and
prediction under global change

**Key Points:**

* Deep ocean ventilation pathways and timescales are assessed in three state-of-the-art models of ocean transport

* Comparisons reveal broad agreement in ventilation pathways but also important differences in ventilation timescales and bulk diffusiveness

* The mid-depth North Pacific takes the longest to renew and is too isolated in nucleus for European modeling of the ocean, possibly due to biased subarctic Pacific dynamics

**Supporting Information:**

Supporting Information may be found in
the online version of this article.

**Correspondence to:**

B. Millet,
bruno-millet@hotmail.com

**Citation:**

Millet, B., de Lavergne, C., Gray, W. R., Éthé, C., Madec, G., Holzer, M., et al. (2025). Deep ocean ventilation: A comparison between a general circulation model and data-constrained inverse models. *Journal of Advances in Modeling Earth Systems, 17*, e2024MS004914.
https://doi.org/10.1029/2024MS004914

Received 18 DEC 2024

Accepted 19 JUN 2025

© 2025 The Author(s). Journal of Advances in Modeling Earth Systems published by Wiley Periodicals LLC on behalf of American Geophysical Union. This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited.

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

# Deep Ocean Ventilation: A Comparison Between a General Circulation Model and Data-Constrained Inverse Models

Bruno Millet¹,²,³ ID, Casimir de Lavergne⁴, William R. Gray³ ID, Christian Éthé⁴, Gurvan Madec⁴, Mark Holzer⁵ ID, Tim DeVries⁶ ID, Geoffrey Gebbie⁷ ID, and Didier M. Roche³,⁸

¹Climate and Environmental Physics, Physics Institute, University of Bern, Bern, Switzerland, ²Oeschger Centre for Climate Change Research, University of Bern, Bern, Switzerland, ³Laboratoire des Sciences du Climat et de l'Environnement, LSCE/IPSL, CEA-CNRS-UVSQ, Université Paris-Saclay, Gif-sur-Yvette, France, ⁴LOCEAN Laboratory, Sorbonne Université-CNRS-IRD-MNHN, Paris, France, ⁵School of Mathematics and Statistics, University of New South Wales, Sydney, NSW, Australia, ⁶Earth Research Institute, University of California, Santa Barbara, CA, USA, ⁷Department of Physical Oceanography, Woods Hole Oceanographic Institution, Woods Hole, MA, USA, ⁸Faculty of Science, Department of Earth Sciences, Vrije Universiteit Amsterdam, Amsterdam, The Netherlands

**Abstract** Ocean ventilation, or the transfer of tracers from the surface boundary layer into the ocean interior, is a critical process in biogeochemical cycles and the climate system. Here, we assess steady-state ventilation patterns and timescales in three models of ocean transport: a 1° global configuration of the Nucleus for European Modeling of the Ocean (NEMO), a recent 2° solution of the Ocean Circulation Inverse Model (OCIM), and a 2° solution of the Total Matrix Intercomparison (TMI). We release artificial dyes in six surface regions of each model and compare equilibrium dye distributions as well as ideal age distributions. We find good qualitative agreement in large-scale dye distributions across the three models. However, the distributions indicate that TMI and OCIM are more diffusive than NEMO. A shallow bias of North Atlantic ventilation in NEMO contributes to a stronger presence of the North Atlantic dye in the mid-depth Southern Ocean and Pacific. This isopycnal communication between the North Atlantic surface and the mid-depth Pacific is very slow, however, and NEMO simulates a maximum age in the North Pacific (NP) about 900 years higher than the data-constrained models. Overly slow NP ventilation persists across NEMO sensitivity experiments encompassing our current best knowledge of diapycnal and isopycnal mixing, pointing to biases in subarctic Pacific dynamics. This study provides a synoptic picture of deep ocean ventilation and a framework for assessing its representation in general circulation models.

**Plain Language Summary** The ocean can be viewed as composed of two compartments: a rapidly mixing surface layer overlying a slowly renewed interior. The transfer of fluid parcels from the surface layer to the ocean interior is often called ocean ventilation. Ocean ventilation is crucial to global climate and marine biology because it controls the injection and storage of heat, carbon, oxygen and nutrients in the vast ocean interior. In this work, we evaluate ventilation simulated by a state-of-the-art computer model of the global ocean named Nucleus for European Modeling of the Ocean (NEMO). We inject virtual dyes at the surface of the model ocean and track their propagation across the deep ocean. We repeat the same experiment with two different models of ocean transport which are constrained by observations. The comparisons reveal broad agreement in ventilation pathways, but also systematic differences in ventilation times. In particular, observationally constrained models predict faster renewal of the deep Pacific Ocean. This result is robust to various representations of turbulent mixing in NEMO and points to important biases in simulated currents within the North Pacific. In general, this study provides a refined picture of ocean ventilation and a framework for its assessment in computer models of the Earth system.

## 1. Introduction

Global ocean circulation is often pictured as a gigantic conveyor belt moving tracers across basins, latitudes and layers (Broecker, 1991; Gordon, 1986). Accordingly, large-scale tracer distributions are commonly interpreted to result from transport and sources or sinks along the conveyor belt. In fact, the circuit of the conveyor belt is itself largely deduced from observed tracer distributions (Richardson, 2008; Talley, 2013). This simple paradigm has proven powerful to describe some of the ocean's roles in climate and biogeochemical cycles (Broecker, 1991). However it has been challenged by the observation that the ocean is fundamentally chaotic (Gregg, 1987;
Lozier, 2010) and the ensuing demonstration that tracer distributions are shaped by turbulent mixing as much as they are by large-scale currents (Holzer & Primeau, 2006; Naveira Garabato et al., 2017).

The chaotic nature of the ocean represents a formidable challenge for models of ocean tracer transport: such models need to represent mixing along and across density surfaces in addition to capturing the large-scale flows of mass. For example, prognostic Ocean General Circulation Models (OGCMs) employ uncertain parameterizations of isopycnal and diapycnal mixing, which have a prominent influence on simulated tracer distributions and transports (Jones & Abernathey, 2019; Oka & Niwa, 2013) and ultimately on the climate simulated by coupled ocean-atmosphere GCMs (Hieronymus et al., 2019; Pradal & Gnanadesikan, 2014). Similarly, inverse models that infer ocean transports from observed tracer distributions often resort to a priori constraints and adjustable mixing parameters to reconstruct plausible three-dimensional transport pathways (Hernández-Guerra & Talley, 2016; Holzer et al., 2021). The chosen constraints are partly subjective and can have a substantial impact on the obtained solutions (Groeskamp et al., 2017).

These difficulties in mapping and modeling tracer transports limit our comprehension of global ocean ventilation. Ocean ventilation broadly signifies the transport of tracers from the surface boundary layer to deeper layers. Here we distinguish ‘ventilation’ from ‘subduction’ and ‘circulation’, the latter two terms traditionally referring to mass transports. Tracer transports are the result of both circulation (or advection) and turbulent mixing (or diffusion). Ocean ventilation is an essential process in the climate system because it affects the baseline climate state and controls the response of ocean heat and carbon reservoirs to perturbations (Gregory et al., 2024; Kwon et al., 2011; Ödalen et al., 2018). In addition, ventilation shapes the distributions of dissolved oxygen, nutrients and other biochemical substances in the ocean interior that are critical to the functioning of marine ecosystems (Davila et al., 2023; Gupta et al., 2022; Holzer, 2022). These considerations have motivated substantial observational and modeling efforts to assess the rates and locations of tracer injection into the ocean interior (England, 1995b; Fine, 2011; Key et al., 2004; Lester et al., 2024; Marzocchi et al., 2021; Orsi et al., 2001). Comprehensive characterization of ventilation trajectories and timescales throughout the modern ocean has also been sought. Such assessments have typically relied on a single OGCM (England, 1995a; Holzer et al., 2020; Holzer & Primeau, 2008; Primeau, 2005) or on a dedicated inverse model (DeVries & Primeau, 2011; Gebbie & Huybers, 2010).

It is established that the deeper half of the ocean (below about 2 km depth) is ventilated primarily from localized sources in the subpolar North Atlantic and on Antarctic shelves (Orsi et al., 2001). Yet the percentage volume filled from each source remains uncertain: studies estimate that the Antarctic source represents between ~ 70% and 300% of its North Atlantic counterpart (DeVries & Primeau, 2011; Gebbie & Huybers, 2010; Holzer et al., 2021; Johnson, 2008; Rae & Broecker, 2018). The composition of the voluminous mid-depth (1–3 km) Pacific is most uncertain because these waters do not fall on the mixing line (in temperature-salinity space) between waters of Antarctic origin and North Atlantic origin (Worthington, 1981), nor do they lie on the circuit of the global conveyor belt (Holzer et al., 2021; Lavergne et al., 2017).

Relatedly, the pathways that connect the surface North Atlantic to the deep Pacific remain debated. Talley (2013) proposed that the vast majority of dense waters formed in the subpolar North Atlantic collectively referred to as North Atlantic Deep Water (NADW) upwell near Antarctica, where they are converted into colder Antarctic Bottom Water (AABW) which subsequently spreads northward in the abyssal Pacific. A Lagrangian analysis (Rousselet et al., 2021) indicates a more direct route skipping upwelling at high southern (HS) latitudes, whereby NADW advects eastward and gains density along the northern flank of the Antarctic Circumpolar Current (ACC) before moving north in the abyssal Pacific. Recent tracer-based studies highlight two routes and suggest that denser varieties of NADW are incorporated in AABW at southern high latitudes whereas lighter varieties spread diffusively at mid-depths (Holzer et al., 2021; Millet, Gray, et al., 2024). Such contrasting paradigms for the connectivity between the surface North Atlantic and deep Pacific carry large implications for nutrient and carbon inventories and sensitivities in the global deep ocean (Baker et al., 2023; Kwon et al., 2012).

In the upper half of the ocean (above about 2 km depth), ventilation has been traced primarily to the production and subduction of mode and intermediate waters at latitudes between 30° and 60° (Gebbie & Huybers, 2011; Samelson & Vallis, 1997). Yet the spreading processes and pathways of these water masses in the ocean interior are poorly characterized (Sloyan & Rintoul, 2001; Toggweiler et al., 2019; You, 2003). For example, the degree to which Subantarctic Mode Water and Antarctic Intermediate Water spread into Northern Hemisphere oceans or
into the abyssal ocean is largely unknown. The role of mixing in controlling ventilation by mode and intermediate waters in northern and southern basins is also subject to debate (Fu et al., 2020; Urakawa & Hasumi, 2012).

Counterintuitively, ventilation timescales are often better known than ventilation routes and volumes. Indeed, whereas estimating ventilation patterns generally requires combining numerous tracer distributions, ventilation timescales can be constrained more directly using radiocarbon (Gebbie & Huybers, 2012; Holzer et al., 2010) or chlorofluorocarbon (Orsi et al., 2001) measurements. Ventilation timescale is usually quantified by the ideal age of seawater, which is the mean time elapsed since this water was last at the surface. Because a given interior water volume can be fed by different surface-to-interior trajectories of different duration, care must be taken when estimating and interpreting the ideal age of a water volume, which is an effective average across trajectories (England, 1995a; Holzer & Hall, 2000; Wunsch, 2002). Nuclear bomb radiocarbon (Broecker & Olson, 1960), the Suess effect (Friedli et al., 1986), corrections for reservoir age (Gebbie & Huybers, 2012), measurement limitations and temporal variability further complicate ideal age reconstructions. Nevertheless, combining chlorofluorocarbon and radiocarbon measurements provides strong constraints on upper-ocean and deep-ocean ventilation rates (Schlitzer, 2007). Generally, studies indicate rapid ventilation in the deep Atlantic, slower ventilation of the deep NP, reaching a maximum age circa 1,500 years, with the deep Indian Ocean falling in between (DeVries & Primeau, 2011; Wunsch & Heimbach, 2008).

Observational inventories of radiocarbon and chlorofluorocarbon (Cimoli et al., 2023; Key et al., 2004; Lavergne et al., 2017; Olsen et al., 2016) provide an invaluable resource to ground-truth OGCM-simulated ventilation (Dunne et al., 2020; Dutay et al., 2002; England, 1995b; Kawasaki et al., 2021; Orr et al., 2017). Nonetheless, the time dependence of atmospheric concentrations and sparse ocean sampling of these tracers complicate model-data comparisons. In addition, the instrumental record and typical model integration times are too short to allow systematic evaluation of deep ocean ventilation occurring over centuries and longer. When included, biochemical tracers such as oxygen offer additional but indirect views on ventilation rates, because their modeled distributions conflate biases in physics and biogeochemistry (Lévy et al., 2022; Talley et al., 2016). Hence, the deep ocean's physical representation in climate-scale OGCM experiments is usually evaluated through the lenses of temperature, salinity and meridional overturning circulation (e.g., Adcroft et al., 2019; Kiss et al., 2020). Given the climatic implications of—and the influence of mixing on—deep ocean ventilation, model evaluation of deep ocean ventilation that goes beyond overturning streamfunctions is desirable (England, 1995a; Holzer et al., 2020).

Here, we evaluate steady-state ventilation patterns and timescales simulated by a global NEMO configuration whose nominal horizontal resolution is 1° (Lavergne et al., 2024), focusing on the deep ocean. This simulation includes an ideal age tracer and regional dye tracers that track ventilation emanating from six distinct ocean regions. The evaluation is performed by comparing equivalent tracer distributions in two data-constrained inverse models: a recent 2° solution of Ocean Circulation Inverse Model (OCIM) (Holzer et al., 2021), and a 2° solution of Total Matrix Intercomparison (TMI) (Gebbie & Huybers, 2010, 2012). Both inverse models are constrained by multiple tracers including radiocarbon, and their solutions provide a global gridded estimate of advective-diffusive transport in the modern ocean consistent with assimilated observations. The three models are state-of-the-art but structurally different: in essence, NEMO is dynamically constrained, TMI is observationally constrained, while OCIM is both observationally and dynamically constrained. The employed models and tracers are described in more detail in the next section. In Section 3, we examine steady-state ventilation volumes, pathways and timescales predicted by the three models. Section 4 further explores salient questions regarding fundamental aspects of deep ocean ventilation raised by the comparison, namely (a) NADW pathways to the deep Pacific, (b) sensitivity of the maximum ideal age to mixing parameterizations, (c) dynamics of the subarctic Pacific and (d) implications for the glacial ocean. We summarize our findings in Section 5.

## 2. Methods

### 2.1. NEMO General Circulation Model

We employ a global configuration of NEMO (NEMO System Team, 2022b, 2022a) that has a nominal horizontal resolution of 1° and 75 layers in the vertical. We build on the work of Lavergne et al. (2024), who performed and evaluated experiments spanning 1,000 simulated years with this configuration under climatological COREv2 forcing at the surface (Large & Yeager, 2009) and geothermal forcing at the bottom (Lucazeau, 2019). They showed that the simulated stratification and circulation are quasi-steady by model year 1,000. In other words, temperature, salinity and velocity fields display very little drift and little inter-annual variability by year 1,000. In
<figure><img src="image_4.png" /><figcaption>**Figure 1.** Surface regions used to define the six regional dye tracers. From south to north, these regions cover Antarctic polar seas (AA), high southern latitudes (HS), mid southern latitudes (MS), low latitudes (LL), the northern North Pacific (NP), and the northern North Atlantic and Arctic (NA). Region boundaries mostly track longitude circles at 30°N, 30°S, 45°S and 63°S or 68°S.</figcaption></figure>

this study, we focus on experiment 'TRA' of Lavergne et al. (2024), which includes the most comprehensive representation of diapycnal mixing by breaking internal tides. We run this experiment for an additional year with output at 5-day frequency, and then use this output to propagate passive tracers in offline mode. The offline mode implies that velocities, diffusivities and stratification parameters are not computed again but instead read and time-interpolated from 5-day output files (corresponding here to model year 1,001), so as to evolve only the passive tracers. In this way, the seasonal cycle is represented. The numerical scheme for passive tracer advection is the same as that used in the online simulation: a second-order Flux Corrected Transport scheme, which is computationally expensive but generates relatively low levels of spurious numerical mixing (NEMO System Team, 2022b, 2022a).

Seven passive tracers are included in the offline simulation: an ideal age tracer and six regional dye tracers released in different regions that cover the entire ocean surface (Figure 1). The ideal age tracer is restored to zero in the top 10 m, and increases at a rate of 1 year per year deeper than 10 m. This tracer was included in the online simulation and is thus extended offline for another 6,000 years. The regional dye tracers are restored to one in the top 10 m within their source region, and restored to zero elsewhere in the top 10 m. Deeper than 10 m, they are conservative, that is, they are only affected by advection and diffusion. The initial state for these dye tracers is zero at every grid point. Like the ideal age tracer, dye tracers are evolved for 6,000 years in our reference offline simulation, which we call REF hereafter. The 6,000-year duration was found necessary to reach quasi steady annual-mean dye and age distributions throughout the ocean. Even though active tracers (temperature and salinity) initialized from observations (Gouretski, 2019) approached cyclo-stationarity in 1,000 years of online simulation, passive tracers (like ideal age) initialized with uniform concentration can take much longer to equilibrate. Past studies of ocean ventilation resorting to similar passive tracers called them 'dye tracers' (e.g., England, 1995a; Wunsch & Heimbach, 2008) or 'surface-labeled water-mass tracers' (e.g., Galbraith & Lavergne, 2019; Primeau & Holzer, 2006). For simplicity, we will stick to the former terminology, noting that equilibrium dye concentrations discussed here are equivalent to so-called water-mass fractions.

To ascertain that the chosen year (1,001) and 5-day frequency of model output from experiment 'TRA' are adequate to reproduce ventilation as would occur in the online simulation, our passive tracers were also simulated online in an extension of 'TRA' from model year 1,000 to 1,500. For this extension, passive tracers were initialized from their final (year 6,000) distributions in the offline experiment REF. Dye distributions and
<table><caption>**Table 1**
Summary of Nucleus for European Modeling of the Ocean Sensitivity Experiments</caption><thead><tr><th rowspan="2">Simulation name</th><th rowspan="2">K<sub>z</sub> includes subinertial contrib.</th><th colspan="3">K<sub>i</sub> depends on</th><th rowspan="2">Simulation length (years)</th></tr><tr><th>Latitude</th><th>Longitude</th><th>Depth</th></tr></thead><tbody><tr><td>REF</td><td>Yes</td><td>Yes</td><td>Yes</td><td>No</td><td>6,000</td></tr><tr><td>DIALOW</td><td>No</td><td>Yes</td><td>Yes</td><td>No</td><td>3,000</td></tr><tr><td>DIAHIGH</td><td>Yes, tripled in NP</td><td>Yes</td><td>Yes</td><td>No</td><td>3,000</td></tr><tr><td>ISOCST</td><td>Yes</td><td>Yes</td><td>No</td><td>No</td><td>3,000</td></tr><tr><td>ISO3D</td><td>Yes</td><td>Yes</td><td>Yes</td><td>Yes</td><td>3,000</td></tr></tbody></table>

*Note.* Diapycnal mixing does not change across experiments except for the representation of mixing by subinertial internal tides.

inventories showed minimal change over the extra 500 years (Figure S1 in Supporting Information S1), validating the offline methodology and confirming the stability of the online experiment.

In order to study the impact of different mixing scenarios on global ocean ventilation, we perform four additional offline NEMO simulations of 3,000-year duration, as summarized in Table 1. The reference (REF) offline simulation follows the experiment 'TRA' of Lavergne et al. (2024): it has an isopycnal diffusivity that is independent of depth but varies with latitude and longitude according to the simulated stratification and flow (Figure 2b). Two simulations, ISOCST and ISO3D, are used to assess the influence of isopycnal mixing parameterization on ventilation. In the ISOCST experiment, isopycnal diffusivity varies only with latitude (Figure 2a). In the ISO3D experiment, isopycnal diffusivity is set by the static three-dimensional map of Groeskamp et al. (2020), based on observations and theory (Figures 2c and 2d). Note that ISO3D and ISOCST share the same stratification, velocities and diapycnal diffusivities as REF: isopycnal diffusion was not changed

<figure><img src="image_5.png" /><figcaption>**Figure 2.** Maps of the isopycnal diffusivity used in Nucleus for European Modeling of the Ocean (NEMO) sensitivity experiments called ISOCST (a), REF (b) and ISO3D (c, d). The diffusivity is independent of depth in both ISOCST and REF. In ISO3D, it varies with depth and is shown here at depths of 200 m (c) and 1,950 m (d). In ISOCST, isopycnal diffusivity is only a function of grid cell size. In REF, it is proportional to the product of the local Rossby radius squared and the Eady growth rate (NEMO System Team, 2022b, 2022a). In ISO3D, an extended mixing length theory including suppression by mean flows (Ferrari & Nikurashin, 2010) is employed to map diffusivity from observed ocean temperature, salinity, and surface geostrophic velocity anomalies (Groeskamp et al., 2020).</figcaption></figure>
during the online simulations. This choice saves computational cost and simplifies interpretation (Jones & Abernathey, 2019), but implies that dynamics did not adjust to isopycnal diffusion changes in ISO3D and ISOCST.

The last two simulations, DIALOW and DIAHIGH, aim to assess the impact of different diapycnal mixing parameterizations. These simulations were run online for 1,000 years, then offline for 3,000 years; hence the impact of vertical mixing on stratification and circulation is accounted for. The DIALOW experiment is identical to REF except that it excludes the mixing contribution of subinertial internal tides (thus corresponding to experiment 'NEW' of Lavergne et al., 2024). The DIAHIGH experiment is identical to REF except that the power input to subinertial internal tides in the NP is tripled. Subinertial internal tides are bottom-trapped internal waves generated by tide-topography interactions at latitudes where the tidal frequency is lower than the Coriolis frequency (Chapman, 1989). They exist only poleward of 30° and their contribution to mixing is most significant in subpolar and polar oceans (Müller, 2013).

Most of the forthcoming analysis will focus on the REF simulation, which shares the same isopycnal mixing ($K_i$) or the same diapycnal mixing ($K_z$) scheme with other simulations. The different parameterizations will be discussed only in Section 4. Hence, unless otherwise stated, all statements about NEMO refer to experiment REF.

## 2.2. Inverse Models

We use two inverse models, TMI and OCIM, to assess global ocean ventilation and its representation in our NEMO experiments. Total Matrix Intercomparison was developed by Gebbie and Huybers (2010). It connects every point of the ocean interior to every surface point by inverting six tracers: potential temperature, salinity, phosphate, nitrate, oxygen and the oxygen-18 isotopic ratio of seawater. The method relies on the notion that any grid box in the ocean interior exchanges volume with its six neighbors. Assuming steady transports in the ocean, and using volume conservation, volume exchanges (representing net advective-diffusive transport) can then be deduced from gridded fields of at least five independent tracers. This approach reconstructs connectivity across the global ocean but not ventilation rates. Transport rates were later estimated by Gebbie and Huybers (2012) assimilating the gridded natural radiocarbon climatology of Key et al. (2004) into the original solution. We use the code available at https://github.com/ggebbie/TMI.mat to simulate ideal age and regional dye tracers for 6,000 years. The simulated fields are defined on a regular longitude-latitude grid of uniform 2° resolution, with 33 layers in the vertical. An alternative TMI solution at 4° resolution was also explored but deemed less realistic; we thus restrict our analysis to the 2° solution.

Ocean Circulation Inverse Model integrates both dynamical and observational constraints to estimate the mean state of the ocean circulation (DeVries & Primeau, 2011). It uses a linearized dynamical ocean circulation model, with adjustable forcing terms in the horizontal momentum equation determined through minimization of a cost function (DeVries & Primeau, 2011). The density distribution is prescribed following the World Ocean Atlas 2013 climatology (Locarnini et al., 2013; Zweng et al., 2013). In its most recent iteration (Holzer et al., 2021), OCIM deduces isopycnal diffusivities as part of the optimization whereas it prescribes diapycnal diffusivities using the three-dimensional map of Lavergne et al. (2020). The cost function includes quadratic mismatch between modeled and observed tracer concentrations. Six tracers are assimilated: potential temperature, salinity, natural radiocarbon ($\Delta^{14}$C), chlorofluorocarbon-11, chlorofluorocarbon-12 and helium-3 ($\delta^3$He). Chlorofluorocarbon concentrations help to better constrain the injection of tracers into the thermocline (DeVries, 2014) while $\delta^3$He and $\Delta^{14}$C help to constrain the abyssal circulation (DeVries & Holzer, 2019). The version of OCIM employed here is OCIM2-48L (Holzer et al., 2021), whose resolution is 2° in longitude and latitude and 48 levels in the vertical. For brevity, OCIM2-48L will be referred to as OCIM in all the following. As for TMI and NEMO, OCIM-derived transports are used to simulate ideal age and six regional dye tracers for 6,000 years.

The analysis of the three models will focus on steady-state distributions, hence on final model years. For OCIM and TMI, the steady state was also solved directly from transport matrices, without any time stepping. We find that solutions from direct inversion match the final state (year 6,000) of time-stepped simulations. We will only show annual means, noting that the seasonal cycle is absent in TMI and OCIM while present in NEMO experiments.
<figure><img src="image_6.png" /><figcaption>**Figure 3.** Panel (a) shows the percentage of the global ocean volume which has been filled as a function of model time, in the reference (REF) Nucleus for European Modeling of the Ocean experiment (blue), in Ocean Circulation Inverse Model (orange) and in Total Matrix Intercomparison (green). The vertical dashed line marks year 2,500. Panel (b) shows the Pacific zonal mean total dye concentration at year 2,500 of the REF simulation. Even though the maximum zonal mean ideal age is about 2,500 years in REF, the core of the shadow zone is only 60% filled by year 2,500.</figcaption></figure>

## 2.3. Equilibration Times

Dye tracers can be used to estimate the volume of the ocean filled as a function of model time. We define the total dye as the sum of the six regional dyes. This total dye is uniformly equal to zero in the initial state. From the first time step onward, it is restored to 1 everywhere at the ocean surface, and thus gradually fills the ocean until its concentration reaches 1 in every grid cell. The total dye concentration in a grid cell thus represents the fraction of the cell's volume which has been filled or ventilated somewhere at the surface, regardless of where. Integrating the total dye concentration over the global ocean volume yields the global volume filled. Steady state is reached when this filled volume equals the model's global ocean volume, since all water must have been at the surface at some point in the past.

The evolution of the filled volume in TMI, OCIM and the reference (REF) NEMO experiment is shown in Figure 3a. It takes almost 6,000 years to fill over 99% of the global ocean volume in REF, compared to only about 3,000 years in TMI and 4,000 years in OCIM. The differences in filled volume across the three models are relatively small during the first 500 years, then the difference between REF and inverse models increases markedly. This indicates that REF predicts slower ventilation of the weakly ventilated oceanic realms. Indeed, we find that the discrepancy comes largely from the mid-depth Pacific shadow zone, where the mean meridional advection is very weak (Holzer et al., 2021; Lavergne et al., 2017). The phrase ‘shadow zone’ refers here to an oceanic domain largely excluded from the meridional overturning circulation (Lavergne et al., 2017); it is analogous but distinct from upper-ocean shadow zones excluded from the horizontal gyre circulation (Luyten et al., 1983). In REF, the core of this shadow zone, between 1.5 and 3 km in the NP, resists complete ventilation for about 6,000 years (Figure 3b). Strikingly, this filling time (taken as the time when the percentage volume filled surpasses 99%) is much larger than the maximum zonal mean Pacific ideal age simulated by REF of about 2,500 years. Similarly, in TMI, the filling time of the shadow zone is about 3,000 years whereas the maximum zonal mean Pacific ideal age is close to 1,600 years. Ocean Circulation Inverse Model produces a similar shadow zone as in REF (Holzer et al., 2021), albeit with a reduced filling time of about 4,000 years. This filling time again far exceeds the maximum zonal mean ideal age in OCIM's shadow zone, which is close to 1,600 years.

The difference between filling time and ideal age can be understood as follows: ideal age equilibration requires the uniform interior source to be balanced by transport from the surface sink, whereas dye tracer equilibration requires unbroken transport from surface source to surface sink via the ocean interior. Hence, dye tracer equilibration time is more akin to residence time, which is the sum of time since last surface contact (ideal age) and time to next surface contact (resurfacing or re-exposure time; Holzer & Primeau, 2008). The large filling times
found in the three models, exceeding 3,000 years, underscore the need for multi-millennial climate simulations to achieve climate states without artificial deep ocean drift (Gebbie & Huybers, 2019; Jansen et al., 2018).

Despite near-complete filling of the global ocean by year 6,000 in REF, the total dye concentration is not exactly equal to 1 in all grid cells. To estimate steady-state dye concentrations, we divide each dye concentration by the total dye concentration in every grid cell. This correction assumes that dye proportions at year 6,000 are representative of the steady state; an excellent approximation given the attained near-complete equilibration. No such correction is needed for OCIM and TMI, whose tracer distributions are spun up to exactly steady state.

Nucleus for European Modeling of the Ocean sensitivity experiments were not carried for 6,000 years and instead stopped at 3,000 years (Table 1). The different mixing parameterizations alter the rates of ventilation of the ocean: for the same year at the same grid cell, the filled volume will vary across the different simulations. In order to compare the different experiments, we first estimate steady-state dye concentrations for each experiment. By year 3,000, about 90% of the ocean has been filled (Figure 3a) and dye concentrations evolve quasi linearly with the total dye concentration (Figure S2 in Supporting Information S1), so that extrapolation is feasible. We thus estimate the equilibrium concentration $C_i(eq)$ of a dye sourced in region i as a function of local dye concentrations at model years 2,900 and 3,000:

$$C_i(eq) = C_i(3000) + \frac{C_i(3000) - C_i(2900)}{\sum_j [C_j(3000) - C_j(2900)]} \times \left(1 - \sum_j C_j(3000)\right) \quad (1)$$

The equilibrium dye concentration is thus the concentration at year 3,000 plus the product of the centennial variation from year 2,900 to 3,000 and the corresponding number of centuries to equilibrium; the latter number being estimated as the ratio of the volume fraction yet to be filled to the centennial variation of the total dye. This extrapolation ensures that the equilibrium total dye concentration, $\sum_j C_j(eq)$, is equal to 1. We assess the accuracy of this extrapolation by applying it to REF, for which we know (to a very good approximation) the true equilibrium concentrations. We find differences between true and extrapolated equilibrium dye concentrations of less than half a percent in the Pacific zonal mean (Figure S3 in Supporting Information S1). This demonstrates the predictive accuracy of Equation 1.

As with the dye concentrations, the ideal age tracer is not fully equilibrated after 3,000 years of simulation. A similar extrapolation can be performed to estimate the steady-state ideal age $A(eq)$ in NEMO sensitivity experiments:

$$A(eq) = A(3000) + \frac{A(3000) - A(2900)}{\sum_j [C_j(3000) - C_j(2900)]} \times \left(1 - \sum_j C_j(3000)\right) \quad (2)$$

Again, we validate this extrapolation by applying it to the REF experiment which was extended to 6,000 years. We find differences of less than 2.5% between extrapolated and true equilibrium ideal ages (Figure S4 in Supporting Information S1).

## 3. Results

We now compare the steady-state dye and ideal age distributions in REF, OCIM and TMI.

### 3.1. Dye Volumes and Distributions

The contributions of each surface source region to the volume of the global, Atlantic and Pacific oceans are shown in Figure 4. In this and following figures, Atlantic and Pacific oceans are extended to Antarctica: boundaries south of 32°S are 68°W, 25°E for the Atlantic and 150°E, 68°W for the Pacific. There is good agreement overall in the ventilated volumes across the three models. The northern North Atlantic (NA) and Antarctic (AA) regions contribute each about a third of the global ocean volume, while HS latitudes (HS) contribute about a fifth (Figure 4a). In the Atlantic, almost 60% of the volume is ventilated from NA, about two times more than from AA (Figure 4b). Conversely, in the Pacific, about 60% of the volume is ventilated south of 45°S (AA and HS), and only about 25% from the NA region (Figure 4c). The 45°S–30°N latitude band (low latitudes (LL) and mid
<figure><img src="image_7.png" /><figcaption>**Figure 4.** Percentage ocean volume ventilated from each surface region. The contributions of Antarctic (AA), high southern (HS), mid southern (MS), low latitude (LL), North Atlantic (NA) and North Pacific (NP) surface regions—as defined in Figure 1—to the global (a), Atlantic (b) and Pacific (c) ocean volumes are given. The blue circle corresponds to the REF experiment of Nucleus for European Modeling of the Ocean, the orange cross to Ocean Circulation Inverse Model and the green cross to Total Matrix Intercomparison. Atlantic and Pacific basins are extended to Antarctica, thus spanning the 80°S–65°N latitude range. The whole depth range of the ocean (0–6,000 m) is included.</figcaption></figure>

southern (MS)) ventilates about 10% of the Pacific, Atlantic and global ocean volumes. The northern North Pacific (NP) ventilates about 5% of the Pacific volume, and virtually none of the Atlantic.

Differences between REF and either OCIM or TMI are of the same order of magnitude as differences between the two data-constrained models (Figure 4). REF tends to have more water of Antarctic origin: the AA global volume fraction of REF is 0.04 and 0.11 higher than its TMI and OCIM counterparts, respectively. At global scale, this excess is mostly compensated by a smaller volume contribution of the neighboring HS region. However, in the Atlantic, compensation occurs mostly through a deficit of NA water. The most noticeable difference between the two inverse models is the larger contribution of NA in OCIM, which manifests mostly outside of the Atlantic and goes with a weaker influence of Antarctic and subantarctic zones (AA and HS).

Similarity in ventilation volumes does not warrant similarity in ventilation patterns. Steady-state ventilation patterns in the Atlantic are illustrated in Figure 5. In all three models, a clear layering is apparent: the top two hundred meters equatorward of 30° are ventilated by the LL region (LL), thermocline waters immediately beneath (0.2–0.5 km) largely originate from 45°S–30°S (MS), HS latitudes (HS) dominate the ventilation of intermediate layers around 0.5–1.5 km depth, the northern Atlantic (NA) ventilates the bulk of the deep Atlantic between 1.5 and 3.5 km depth, and Antarctic marginal seas (AA) take over in densest layers. The two data-constrained models display strong resemblance of their dye distributions, despite weaker vertical gradients in TMI compared to OCIM. REF shows the strongest layering, including a sharp separation between NA and AA dyes near 3.5 km depth. This stronger layering in REF is associated with overly pure and abundant AA dye in the abyssal Atlantic, and overly pure and abundant NA dye in the overlying mid-depths. This bias likely stems from the lack of dense varieties of NADW overflowing from the Nordic seas (Lavergne et al., 2024), a common deficiency of OGCM simulations (Danabasoglu et al., 2014). It should be noted that REF has little bias in the overall density stratification of the water column (Figure S6 in Supporting Information S1; Lavergne et al., 2024), therefore biases cannot be attributed to overly dense AABW and overestimated deep stratification.

In the Pacific, steady-state zonal mean dye concentrations also show broad qualitative similarity between the three models (Figure 6). The AA dye decreases from the bottom upward and dominates throughout much of the deep Pacific, up to about 2 km depth. The subantarctic (HS) dye is predominant between 0.5 and 2 km depth between 60°S and 30°N; it is also important in deeper layers within TMI. This indicates that ventilation from the subantarctic zone spreads well into the Northern Hemisphere. By contrast, the NP dye is largely confined to the Northern Hemisphere and depths shallower than 1.5 km. The NA dye has a strong but relatively uniform influence throughout depths greater than 1.5 km. Finally, the low-latitude and mid-latitude southern dyes fill mostly the top 500 m of the low-latitude Pacific.
<figure><img src="image_8.png" /><figcaption>**Figure 5.** Steady-state zonal means of the different dye concentrations over the Atlantic basin for Total Matrix Intercomparison (left column), Ocean Circulation Inverse Model (middle column) and the REF Nucleus for European Modeling of the Ocean experiment (right column). Dyes originate from the surface at low latitudes (a–c), southern mid latitudes (mid southern) (d–f), the North Pacific (g–i), high southern (HS) latitudes (HS) (j–l), the North Atlantic and Arctic (NA) (m–o) or the Antarctic margins (AA) (p–r). Semi-enclosed marginal seas are excluded from zonal averages.</figcaption></figure>

Major differences in Pacific ventilation patterns are also apparent. In TMI, the subantarctic dye spreads downward into the abyssal ocean, such that it plays a major role in the deep South Pacific ventilation (Figure 6j). Such vertical penetration is less pronounced in OCIM and even less in REF (Figures 6k and 6l). In the NP, the NP dye
<figure><img src="image_9.png" /><figcaption>**Figure 6.** Steady-state zonal means of the different dye concentrations over the Pacific basin for Total Matrix Intercomparison (left column), Ocean Circulation Inverse Model (middle column) and the REF Nucleus for European Modeling of the Ocean experiment (right column). Dyes originate from the surface at low latitudes/low latitude (a–c), southern mid latitudes (mid southern) (d–f), the North Pacific (g–i), high southern (HS) latitudes (HS) (j–l), the North Atlantic and Arctic (NA) (m–o) or the Antarctic margins (AA) (p–r). Semi-enclosed marginal seas are excluded from zonal averages.</figcaption></figure>

spreads deeper in the inverse models than it does in REF (Figures 6g–6i), suggesting that NEMO simulates
vertical transports that are too weak within this region. REF also stands out for its distinct mid-depth maximum of
the NA dye up to about 20°N (Figure 6o). This mid-depth maximum is visible in the Southern Ocean of OCIM but
<figure><img src="image_10.png" /><figcaption>**Figure 7.** Transient Southern Ocean ventilation in Total Matrix Intercomparison (a, b), Ocean Circulation Inverse Model (c, d) and the REF Nucleus for European Modeling of the Ocean experiment (e, f). Shading shows the total dye concentration minus the NA dye concentration, at model year 100 and at depths of 300 m (a, c, and e) and 1 km (b, d, and f). The red line corresponds to the separation between AA and high southern source regions.</figcaption></figure>

altogether absent in TMI (Figures 6m and 6n). A mid-depth maximum of the NA dye exists in all three models in the Atlantic sector of the Southern Ocean (Figures 5m–5o); however in REF this maximum is higher and less rapidly eroded in the ACC (Figure S5 in Supporting Information S1). The northward penetration of this mid-depth tongue of NA dye in the Pacific is due to isopycnal diffusion rather than advection (Millet, Gray, et al., 2024).

A recurring feature in the analyzed dye distributions is the lower contribution of the subantarctic surface region (HS) to deep ocean ventilation in REF compared to OCIM and especially to TMI. Both TMI and OCIM predict the subantarctic dye to spread into the abyssal ocean, whereas REF does not. This difference may stem from differing levels of deep diapycnal mixing and/or from differing degrees of entrainment of the subantarctic dye in the lower branch of the overturning (see Figures 5j–5l and 6j–6l).

To shed light on these possibilities, we show in Figure 7 the total dye concentration (minus the NA dye concentration) at depths of 300 m and 1 km in the Southern Ocean, after only 100 years of simulation. A ventilation minimum is visible in TMI, OCIM and REF at 300 m depth between the Polar Front and the Antarctic continental slope (Figures 7a, 7c and 7e). We interpret this ventilation desert as the signature of advective upwelling, which is focused in this area of the Southern Ocean (Gray et al., 2023; Styles et al., 2024). The separation between the subantarctic (HS) and Antarctic (AA) surface regions (red line in Figure 7) largely tracks the ventilation minimum, which justifies the choice of boundary. Beyond these shared patterns, the ventilation desert is evidently more pronounced in REF, such that AA and HS dyes barely mix in the subsurface Southern Ocean, each following a distinct physical ventilation route. By contrast, ventilation is more geographically diffuse in OCIM, and still more in TMI, such that AA and HS dyes intertwine at depth and partly share antarctic and subantarctic ventilation routes. Furthermore, the surface-released dye propagates faster and deeper in the ACC in both inverse models (Figures 7b, 7d and 7f), suggesting higher effective mixing there. It is also apparent that TMI lacks deep
<figure><img src="image_11.png" /><figcaption>**Figure 8.** Intercomparison of diffusiveness across the three models. Panel (a) represents the average dye concentration weighted by dye volume, and panel (b) displays the similarity score of the density distribution of water volume and the density distribution of dye volume. Both metrics are shown for each dye (x-axis) and each model (symbols). Blue circles correspond to the REF Nucleus for European Modeling of the Ocean experiment, orange crosses to Ocean Circulation Inverse Model and green crosses to Total Matrix Intercomparison. The distributions have been normalized prior to calculating the similarity score, so that this metric only quantifies differences in the shapes of distributions. Larger similarity scores indicate lower bulk diapycnal mixing. Density is measured with $\sigma_2$, the potential density referenced to 2,000 dbars.</figcaption></figure>

penetration of dye within the Antarctic marginal seas on this centennial timescale (Figure 7b), whereas both NEMO and OCIM capture sinking down the Antarctic continental slope (Figures 7d and 7f).

## 3.2. Quantifying Water-Mass Mixing

The foregoing qualitative analysis of dye distributions suggests that REF predicts stronger layering and weaker mixing of the different dyes than do both inverse models. Here we quantify the degree to which dyes mix with other dyes (Figure 8a) or mix across density surfaces (Figure 8b) in the three models. The first measure of ‘diffusiveness’ does not distinguish between isopycnal and diapycnal mixing: it is simply the average dye concentration weighted by dye volume. At every grid cell $n$, we know the volume of the cell $V(n)$, the concentration $C_i(n)$ of each dye ($i$), and thus the volume occupied by each dye $V_i(n) = C_i(n) \times V(n)$. We can thus compute the average dye concentration weighted by the volume it occupies, $C_{i,avg}$, as follows:

$$C_{i,avg} = \frac{\sum_n C_i(n) \times V_i(n)}{\sum_n V_i(n)} = \frac{\sum_n C_i(n)^2 \times V(n)}{\sum_n C_i(n) \times V(n)} \quad (3)$$

If this weighted average concentration is close to one, it implies that the dye is extremely pure on average; in other words, it remains largely separated from other dyes. If instead the weighted average concentration of a dye is close to the global volume fraction it fills (cf Figure 4a), dilution is maximum. We find that REF is least diffusive in that its dyes tend to have higher weighted average concentrations (Figure 8a). The difference is largest for the low-latitude, mid-latitude southern and Antarctic dyes, hence for the lightest and densest water masses, which are much less diluted in REF than in OCIM or TMI. For the subantarctic, NP and especially the North Atlantic dye, differences across the three models are less pronounced. This is likely because diapycnal mixing plays a smaller role in the dilution of the latter three dyes, which spread through intermediate and mid-depth layers (away from strongly mixing boundaries) largely via horizontal circulation and isopycnal diffusion.

Indeed, the hierarchy between the three models and the different dyes is more obvious when only mixing across density surfaces is considered. To quantify degrees of diapycnal diffusion, we first bin the ocean volume as a function of annual mean $\sigma_2$ (the potential density referenced to 2,000 dbars) in each of the three models. We then repeat this calculation for each dye, that is, we bin dye volumes as a function of $\sigma_2$. We thus obtain the density distribution of total water volume, and the density distribution of dye volume, for each dye and each model. Should a dye be completely mixed across the global ocean, its concentration would be uniform and its distribution as a function of $\sigma_2$ would mirror the distribution of total water volume. To measure the degree of mixing, we thus
normalize all distributions and calculate the similarity score ($\Delta\rho_i$) between the volume distribution for a given dye $V_i(\sigma_2)$ and the distribution for total water volume $V(\sigma_2)$:

$$ \Delta\rho_i = \sqrt{\sum_{\sigma_2} (V(\sigma_2) - V_i(\sigma_2))^2} \quad (4) $$

The larger $\Delta\rho_i$, the lower the bulk diapycnal mixing. We find that this score is systematically higher for REF than for TMI and OCIM (Figure 8b). Differences between TMI and OCIM vary across dyes: low-latitude, mid-latitude southern and NP dyes spread more across isopycnals in OCIM, whereas subantarctic, North Atlantic and Antarctic dyes show stronger diapycnal spreading in TMI. These contrasted differences are consistent with qualitative patterns visible in Figures 5 and 6. Our quantification demonstrates notably that the subantarctic dye features the largest discrepancies in diapycnal spreading across models (Figure 8b), indicating that ventilation from the subantarctic zone remains poorly modeled and constrained. In addition, we find a consistent hierarchy across dyes: the North Atlantic dye is the most well-mixed in the global ocean, whereas low-latitude, mid-latitude southern and NP dyes are the least well-mixed, as might be expected from their more regional realms.

Overall, these diagnostics demonstrate that REF is less diffusive than the two inverse models. In addition, TMI is more diffusive than OCIM in the deep ocean, consistent with previous findings (Davila et al., 2023). Density stratification is very similar in REF and inverse models (Figure S6 in Supporting Information S1). Hence, the diagnosed difference in bulk diffusiveness suggests that in our NEMO configuration (a) dyes are injected within smaller (and less overlapping) density ranges, and/or (b) interior diapycnal mixing is underestimated. The lack of denser varieties of NADW in REF, and the good match of parameterized tidal mixing in REF with available in situ mixing estimates (Lavergne et al., 2020), are in support of hypothesis (a). On the other hand, several sources of deep ocean mixing are not represented in REF, such as wind-induced low-mode near-inertial waves (Alford, 2020) and interactions between geostrophic currents and topography (Evans et al., 2022; Klymak, 2018; Naveira Garabato et al., 2019; Yang et al., 2021). Efforts to faithfully parameterize these processes in NEMO can contribute to a more realistic simulation of ventilation rates and pathways. Indeed, the low bulk diffusiveness across all dyes in NEMO signals that spurious numerical mixing by the model's discrete advection schemes (Burchard & Rennau, 2008; Holmes et al., 2021) is sufficiently weak to make the addition of physical sources of mixing beneficial.

### 3.3. Ideal Age

Steady-state dye distributions illustrate ventilation patterns of the global ocean but do not provide information on ventilation timescales. We now examine steady-state distributions of the ideal age tracer in the three models. The ideal age predicted by TMI, OCIM and REF at 500 m depth, 2 km depth and the seafloor is shown in Figure 9. The models share large-scale patterns. Waters are younger than 300 years at 500 m depth except in upwelling zones such as the equatorial Pacific and polar Southern Ocean (Figures 9a, 9d and 9g). The maximum age at 500 m is actually located in the northwest corner of the Pacific, signaling the upwelling of particularly old deep waters in this region (Millet, Gray, et al., 2024; Sarmiento et al., 2004; Toggweiler et al., 2019). A strong inter-basin contrast at 2 km depth features centenarian waters in the Atlantic and thousand-year-old waters in the Pacific, with the Indian Ocean lying between these extremes (Figures 9b, 9e and 9h). Ideal age increases northward in the deep Indian and Pacific basins, consistent with ventilation from the south (Figures 9b, 9c, 9e, 9f, 9h, and 9i). In the Atlantic, age increases with depth below 2 km, due to slower ventilation by AABW relative to NADW; the reverse is true in the deep Pacific and eastern Indian oceans, where northward-flowing bottom waters are the youngest water mass (Figures 9b, 9c, 9e, 9f, 9h, and 9i).

Major differences between models are also apparent. Most conspicuously, the REF experiment has generally higher ages than the inverse models (Figure 9). The largest difference is in the mid-depth Pacific shadow zone, where ideal age peaks at about 2,500 years in REF, compared to about 1,600 years in TMI and OCIM (Table 2). The location of the maximum age is also different: it lies at the northeast end of the Pacific in REF, but near 35°N in inverse models (Figures 9b, 9e and 9h) like the minimum in measured radiocarbon (Key et al., 2004). The overly sluggish ventilation of the shadow zone in REF contributes to higher ages at shallower depths in upwelling zones, especially in the subarctic Pacific (Figure 9g). In the Atlantic, differences are less pronounced except at
<figure><img src="image_12.png" /><figcaption>**Figure 9.** Steady-state ideal age at 500 m depth (a, d, and g), 2 km depth (b, e, and h) and at the seafloor (c, f, and i) in Total Matrix Intercomparison (a–c), Ocean Circulation Inverse Model (d–f) and the REF Nucleus for European Modeling of the Ocean experiment (g–i).</figcaption></figure>

abyssal depths, where REF has higher ages associated with a deficit of ventilation from the North Atlantic (see Section 3.1).

Figure 10 shows the zonal mean ideal age in the Pacific for the three models. The range of values in TMI and OCIM is similar, but the depth structure is not the same: TMI has a much weaker depth dependence of age in the deep Pacific (Figures 10a and 10b). We attribute this difference to the radiocarbon data used: TMI uses the GLODAP gridded product of bomb-corrected radiocarbon (Key et al., 2004), whereas OCIM assimilates only local radiocarbon measurements from GLODAPv2 (Olsen et al., 2016) that have no anthropogenic imprint (based

on absence of chlorofluorocarbons). The latter method incorporates additional measurements, avoids their interpolation and obviates the correction for bomb-produced ¹⁴C (Rubin & Key, 2002); correction which appears to introduce biases in the deep Pacific Ocean (Burke et al., 2015; Lavergne et al., 2017; Sweeney et al., 2007). Thus, we posit that (a) both inverse models predict realistic age maxima in the Pacific shadow zone, where natural ¹⁴C is directly measured; and (b) the spatial structure of age is more accurate in OCIM.

In the REF NEMO experiment, zonal mean Pacific ideal age has a similar structure to OCIM, but the age exceeds that of OCIM across most of the basin (Figures 10b and 10c). Bottom waters entering and filling the South Pacific

<table><caption>**Table 2**
Maximum Ideal Age (in Years) in the North Pacific Across the Three Models</caption><thead><tr><th></th><th>NEMO</th><th>OCIM</th><th>TMI</th></tr></thead><tbody><tr><td>500 m</td><td>1,492</td><td>902</td><td>878</td></tr><tr><td>2,000 m</td><td>2,519</td><td>1,640</td><td>1,535</td></tr><tr><td>4,000 m</td><td>2,022</td><td>1,386</td><td>1,527</td></tr></tbody></table>

*Note.* The maximum is taken within the Pacific between 0° and 65°N, excluding marginal seas, along three geopotential surfaces: 500 m, 2 and 4 km depth.
<figure><img src="image_13.png" /><figcaption>**Figure 10.** Pacific zonal mean ideal age in Total Matrix Intercomparison (a), Ocean Circulation Inverse Model (b) and the REF Nucleus for European Modeling of the Ocean experiment (c). Semi-enclosed marginal seas are excluded from zonal averages.</figcaption></figure>

abyss display the same age in OCIM and REF, thus the higher mid-depth ages in REF do not stem from waters that are already too old when they enter the Pacific but rather from excess aging within the Pacific. Note that the Pacific density stratification is not overestimated in REF (Figure S6 in Supporting Information S1); however, the Pacific abyssal overturning is slower in REF than in OCIM (Figure S7 in Supporting Information S1; Lavergne et al., 2024). As we discuss in Section 4, a lack of diffusive and/or advective upwelling in the subarctic Pacific of REF is plausibly an important contributor to the age bias. The widespread overestimation of age apparent in Figure 9 further suggests that systematic underestimation of isopycnal and/or diapycnal diffusion in REF is likely. We explore hereafter the sensitivity of simulated age and dye distributions to mixing parameterizations in NEMO.

# 4. Discussion

Comparisons detailed above highlighted three distinctive features of the Pacific Ocean simulated by the REF NEMO experiment: (a) a mid-depth maximum of the North Atlantic dye concentration that pervades the basin; (b) excessive isolation of mid-depth waters; and (c) an age maximum located in the northeast corner of the basin, rather than further south as suggested by radiocarbon observations and inverse models. We now explore ventilation processes underpinning these features and their sensitivity to mixing parameters, resorting to NEMO experiments described in Table 1. We remind the reader that REF, ISO3D and ISOCST experiments share the exact same velocities and stratification: they only differ by the isopycnal diffusivity felt by passive tracers.

## 4.1. From the North Atlantic Surface to the Deep Pacific

Steady-state dye distributions do not always allow a clear identification of ventilation pathways. To aid interpretation, we show in Figure 11 the transient spreading of the NA dye in the Pacific, at year 800 of the REF simulation. Two pathways are clearly apparent: an abyssal pathway that tracks the circulation of AABW, and a mid-depth pathway that transports the NA dye horizontally from its core in the ACC toward the north. The weakness of simulated overturning flows in this mid-depth zone (especially north of the equator; Figure S7 in Supporting Information S1) and the large latitudinal gradient in the NA dye suggest that isopycnal diffusion is the main driver of this mid-depth propagation, consistent with recent proposals (Holzer et al., 2021; Millet, Gray, et al., 2024). Importantly, the northward propagation of the NA dye at mid-depths is at odds with the common view that tracers move toward the south (along the conveyor) in this layer (Chen et al., 2022; Talley, 2013). We argue that tracers can move either northward or southward in the mid-depth Pacific shadow zone, depending on their large-scale isopycnal concentration gradient.

The role of isopycnal diffusion in the mid-depth spreading of the NA dye is further evidenced by analysis of NEMO sensitivity experiments (Figure 12). We focus on steady-state distributions of the zonal mean NA dye concentration in the Pacific. In ISO3D, which has substantially weaker isopycnal diffusivities than REF in the deep ocean (Figures 2b, 2d), a mid-depth maximum of the NA dye is present but only south of 5°N (Figure 12b). By contrast, ISOCST has a slightly stronger mid-depth tongue of NA dye than REF, extending to 40°N, consistent with higher average isopycnal diffusivities across the basin (Figures 2a and 2b). Differences in isopycnal
<figure><img src="image_14.png" /><figcaption>**Figure 11.** Pacific zonal mean NA dye concentration in the REF Nucleus for European Modeling of the Ocean experiment after 800 years of simulation. White arrows are a simplified depiction of the two pathways of NA dye spreading: in the abyss along the abyssal overturning cell, and at mid-depth via isopycnal diffusion. Note the dye concentration scale is different than in other figures, for improved visualization.</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>**Figure 12.** Pacific zonal mean NA dye concentration at steady state, for each Nucleus for European Modeling of the Ocean experiment. Semi-enclosed marginal seas are excluded from zonal averages. White lines mark depths of 1.5 and 2.5 km. Panels have been organized so that isopycnal diffusion in the deep Pacific increases toward the right, from ISO3D (b) via REF (c) to ISOCST (d), while diapycnal diffusion increases toward the top, from DIALOW (e) via REF (c) to DIAHIGH (a). The colorscale saturates at a concentration of 0.6 so that differences across panels are clearer. To aid quantification, in each panel we show the NA dye concentration (in percent) averaged over 1.5–2.5 km depth (dashed black lines) at 50°S and 40°N.</figcaption></figure>
<figure><img src="image_16.png" /><figcaption>**Figure 13.** Pacific zonal mean ideal age at steady state in the different Nucleus for European Modeling of the Ocean experiments: REF (a), DIALOW (b), ISOCST (c), ISO3D (d) and DIAHIGH (e). The equivalent in Ocean Circulation Inverse Model is shown in panel (f). Semi-enclosed marginal seas are excluded from all zonal averages. Gray contour lines show 500 years intervals in the zonal mean ideal age.</figcaption></figure>

diffusivity outside the deep Pacific can also explain some of the differences in dye distributions. In particular, ISO3D has smaller mean isopycnal diffusivities south of 45°S (Figure 2); this leads to weaker resurfacing and removal of the NA dye at the Southern Ocean surface, and thereby to an overall increase in NA dye content within the ACC and across the global deep ocean (Figures 12b, 12c, Figure S5 in Supporting Information S1). Thus, changes in isopycnal mixing alone can modify the balance between northern and southern ventilation of the deep ocean (Jones & Abernathey, 2019).

Diapycnal mixing also influences the northward spreading of the NA dye in the Pacific. Removing the contribution of subinertial internal tides to diapycnal mixing (as in DIALOW) enhances the mid-depth maximum of NA dye concentration, which then spreads to the northern boundary of the Pacific (Figure 12e). This is due to weaker diapycnal diffusion south of 30°S, which reduces removal of NA dye in the polar Southern Ocean, and to weaker diapycnal diffusion in the Pacific north of 30°N, which reduces local erosion of the mid-depth maximum. Tripling mixing by subinertial tides in the NP (as in DIAHIGH) does not significantly alter the NA dye distribution, because of weak vertical gradients already present in REF north of 30°N (Figures 12a and 12c).

Overall, the presence of a marked mid-depth maximum of NA dye in the Pacific sector of the ACC which spreads northward in the basin near 2 km depth is a robust feature of our NEMO sensitivity experiments (Figure 12). This maximum is stronger and shallower than that found in the ACC of OCIM (Figures 6n and 6o). The light and shallow bias of NADW in our NEMO configuration (Lavergne et al., 2024) may largely explain this difference. Specifically, the lack of NADW in the abyssal Atlantic in NEMO may favor the mid-depth pathway over the abyssal pathway connecting the surface North Atlantic to the deep Pacific (Figure 11).

## 4.2. The Maximum Age of the Deep Pacific

Pacific ventilation timescales display stronger sensitivity to mixing parameterizations than do ventilation patterns. Figure 13 compares Pacific zonal mean ideal ages in the five NEMO experiments and in OCIM. The broad structure is similar across all simulations, with a peak zonal mean age located near 2 km depth in the NP. However the magnitude of this peak age varies between 2,300 years (DIAHIGH) and 3,000 years (ISO3D) in NEMO,
<figure><img src="image_17.png" /><figcaption>**Figure 14.** Ideal age at 2 km depth in the North Pacific simulated by REF (a), DIALOW (b), ISOCST (c), ISO3D (d), DIAHIGH (e) and Ocean Circulation Inverse Model (OCIM) (f). Red crosses indicate the location of the maximum value of age in each simulation. The difference in age between DIALOW and REF (g), DIAHIGH and REF (h), and ISO3D and REF (i) is also shown. Note that the ideal age of OCIM is multiplied by two in panel (f) to share the range of other panels. For the same reason the difference between ISO3D and REF is divided by two in panel (i).</figcaption></figure>

compared to only about 1,600 years in OCIM. Hence all NEMO experiments simulate an overly isolated Pacific shadow zone. The degree of isolation decreases as bulk isopycnal diffusivity in the deep Pacific increases (from ISO3D via REF to ISOCST) and as bulk diapycnal diffusivity increases (from DIALOW via REF to DIAHIGH). Puzzlingly, simulation ISO3D, which includes state-of-the-art maps of diapycnal (Lavergne et al., 2024) and isopycnal diffusivity (Groeskamp et al., 2020), possesses the largest age bias. A more realistic three-dimensional structure of isopycnal diffusivity, including a decrease of diffusivity with depth, thus reinforces the regional deficit of ventilation.

This result could imply that diapycnal diffusivity is substantially underestimated in the deep NP. However, tripling the mixing power of subinertial internal tides in the NP reduces the peak age by only about 150 years (compare DIAHIGH and REF, Figures 13a, 13e, and 14h). This reduction is less than that obtained by adding the mixing contribution of subinertial internal tides in the first place (compare REF and DIALOW, Figures 13a, 13b, and 14g). The increase in shadow zone ventilation rate from DIALOW to REF is mostly due to regionally increased mixing in the Pacific north of 30°N, since age differences at the southern entrance of the Pacific are relatively small (Figures 13a, 13b, and 14g). Hence, increasing regional diapycnal mixing further than tested in DIAHIGH may not achieve much further reduction in peak age, unless the increase in diapycnal mixing is
<figure><img src="image_18.png" /><figcaption>**Figure 15.** Pacific zonal mean difference of the North Pacific (a, d), high southern (HS) latitude (HS) (b, e) and Antarctic (AA) (c, f) dye concentrations between DIALOW and REF (a–c) and between DIAHIGH and REF (d–f).</figcaption></figure>

unrealistically strong or widely spread. We hypothesize that our NEMO experiments misrepresent vertical mass transports across the shadow zone in the subarctic Pacific (see Section 4.3), and thereby systematically overestimate deep Pacific age.

Analysis of the spatial pattern of ideal age at 2 km depth corroborates some of these inferences. In REF, DIALOW and ISOCST, the maximum age is found in the northeastern corner of the basin, near 140°W, 55°N (Figures 14a–14c). By contrast, it is situated about 20° further south in ISO3D, DIAHIGH and OCIM (Figures 14d–14f), in accord with radiocarbon observations (Lavergne et al., 2017). The addition of mixing by subinertial internal tides, which is particularly strong along the Aleutian arc and in the Kuril Straits (Müller, 2013), contributes to shift the peak age away from the northern rim of the basin and away from the Bering Sea (Figures 14g and 14h). Reduced isopycnal diffusivity in the deep NP, with lowest values in the eastern half of the subpolar gyre (Figures 2c and 2d) as predicted by Groeskamp et al. (2020), also favors the emergence of peak ages at middle latitudes offshore of North America (Figure 14i). Lower isopycnal diffusivities enhance the footprint of bathymetrically constrained horizontal circulation pathways and diapycnal mixing patterns on tracer distributions (Hautala, 2018), leading to a more complex age distribution. Overall, the magnitude and pattern of age biases across NEMO experiments suggest that REF overestimates large-scale isopycnal homogenization rates but underestimates vertical communication rates, notably in the subarctic Pacific.

### 4.3. Upwelling in the Subarctic Pacific

Maps of ideal age at 500 m depth in TMI, OCIM and REF suggest localized upwelling of old deep waters in the subarctic Pacific (Figures 9a, 9d and 9g). The presence of such upwelling has been inferred previously from in situ observations and regional tracer budgets (Holzer et al., 2021; Millet, Gray, et al., 2024; Sarmiento et al., 2004; Toggweiler et al., 2019; Warren, 1983). However it remains unclear what is the driver of tracer upwelling from abyssal depths to the surface in this understudied region. Candidates include tidal mixing (Sarmiento et al., 2004) and wind-driven upwelling (Toggweiler et al., 2019).

Internal tide-driven mixing in the subarctic Pacific is dominated by subinertial internal tides (Lavergne et al., 2024; Müller, 2013; Tanaka, 2023). Comparison of DIALOW, REF and DIAHIGH experiments may thus
<figure><img src="image_19.png" /><figcaption>**Figure 16.** Annual mean Ekman pumping velocity in REF (a) and in a high-resolution ocean reanalysis (Copernicus Marine Service, 2024) (b) in the North Pacific (NP). The zonally averaged Ekman velocity in both REF and the reanalysis is shown in panel (c). This average was taken over the Pacific basin excluding semi-enclosed marginal seas. Vertical velocity at 300 m depth (d) and 1,000 m depth (e) in REF. In panel (f), we plot the horizontally integrated vertical velocity in the NP (purple), and the same integral for negative-only (red) and positive-only (blue) vertical velocities. The x-axis of panel (f) uses a logarithmic scale and units of Sverdrups (Sv; $1 \text{ Sv} \equiv 10^6 \text{ m}^3 \text{ s}^{-1}$).</figcaption></figure>

shed light on the role of tidal mixing on vertical tracer transports in the northern NP. Figure 15 shows zonal mean differences in dye concentrations across the three experiments. Sensitivity of dye concentrations peaks north of 40°N between 300 m and 2 km depth. Zonal mean concentrations of the Antarctic (AA) dye tracer in this depth range increase by up to 0.03 from DIALOW to REF, and by a similar amount from REF to DIAHIGH (Figures 15c and 15f). The increased upward spreading of the Antarctic dye is accompanied by increased downward spreading of the NP dye; both increases are compensated by a drop of the subantarctic dye (Figures 15a, 15b, 15d, and 15e). This demonstrates that mixing by subinertial internal tides contributes to vertical transports across the mid-depth shadow zone, and in particular to the upwelling of tracers from the abyss to the near-surface in the subarctic Pacific. Nevertheless, the effect on dye distributions is limited: it changes the overall balance of ventilation (Figure 6) by a few percent only.

Radiocarbon observations (Toggweiler et al., 2019) and the simulated age pattern at 500 m depth (Figures 9a, 9d and 9g) point to the area just east of the Kamchatka peninsula as the primary upwelling zone. This location does not coincide with the strongest rates of tidal mixing, suggesting that other—likely advective—processes are at play. To probe the representation of vertical flows in our NEMO experiments, we map the annual mean Ekman pumping velocities implied by the wind stress forcing of REF (Figure 16a) and the equivalent in a high-resolution ocean reanalysis (Figure 16b). Both REF and the reanalysis feature Ekman upwelling across most of the subpolar gyre, including regionally enhanced upwelling in the northeastern corner of the basin and off the Kamchatka peninsula. Averaged zonally over the main Pacific basin, the annual mean Ekman pumping in REF is comparable to that of the reanalysis (Figure 16c). Nonetheless, local maxima in Ekman suction close to the coast—including along the Kamchatka peninsula—are weaker in REF, potentially causing muted vertical transports in these hotspots (Figures 16a and 16b). Hence, wind stress forcing in our NEMO experiments (which is based on the COREv2 normal year data set; Large & Yeager, 2009) does not appear to misrepresent Ekman suction on average in the northern NP, but could contribute to underestimate wind-driven upwelling at some key sites.
The actual vertical velocities in REF do not simply reflect Ekman dynamics: while upwelling dominates at 300 m in the cyclonic subpolar gyre, large vertical velocities of both signs are simulated along the boundary (Figure 16d). Note that vertical velocities include a parameterized eddy-driven component (Gent & McWilliams, 1990) in addition to the resolved component; however we find the parameterized component to be one order of magnitude smaller than the resolved component in this region (not shown). Near-boundary vertical velocities vary substantially and frequently change sign from one grid cell to the next (Figure 16d). The same is true at 1,000 m depth (Figure 16e). Vertical tracer communication in the subarctic Pacific in REF could thus be dominated by these large opposing velocities. Horizontal integration of the vertical velocity across the NP indicates net upwelling of $1 - 3 \times 10^6 \text{ m}^3 \text{ s}^{-1}$ between 0.5 and 5 km depth, which is the residual of much larger upward and downward volume transports (Figure 16f).

Alternating up and down velocities simulated near the bottom boundary (Figures 16d–16f) are not the expected physical response to bottom-intensified diapycnal mixing (Ferrari et al., 2016; Lavergne et al., 2016; McDougall & Ferrari, 2017): they are too large and do not coincide with mixing hotspots. They are also inconsistent with fine-scale dynamics of deep currents that hug and rub continental slopes (Capó et al., 2024; Gula et al., 2016; Naveira Garabato et al., 2019); these dynamics are unresolved here. Rather, they reflect numerical artifacts associated with the model’s staircase-like topography, a common issue in depth-coordinate OGCMs (e.g., Adcroft et al., 1997; Debreu et al., 2020). In sum, our NEMO experiments do not faithfully simulate vertical flows along the northern rim of the Pacific, due to misrepresented dynamics near continental slopes and possibly to local wind forcing biases. As a result, the observationally inferred upwelling in the subarctic Pacific is poorly captured. This likely explains why the deep NP age bias persists regardless of mixing parameterization choices. Our findings imply that improving the fidelity of deep ocean ventilation in OGCMs requires a better representation of currents near sloping topography in addition to improved isopycnal and diapycnal mixing.

## 4.4. Perspectives on Glacial Ocean Ventilation

While focusing on the present-day ocean, this study offers insights of potential relevance to understanding ventilation of the glacial ocean. During the Last Glacial Maximum (LGM) about 21,000 years ago, the salinity difference between bottom and intermediate waters is thought to have been larger than today (Ferrari et al., 2014; Gebbie & Huybers, 2019), leading to a stronger deep ocean stratification (Adkins et al., 2002) and a greater disconnect between NADW and AABW realms (Ferrari et al., 2014; Gebbie, 2014). To draw apart AABW and NADW, models generally require a large increase in salt stratification, typically achieved via enhanced brine rejection in the Southern Ocean (Bouttes et al., 2010; Jansen, 2017). Here, our NEMO experiments simulate a rather strict separation between NADW and AABW (as measured by the NA and AA dye tracers), especially in the Atlantic and Southern Oceans (Figure 5), despite having a bulk density stratification close to present day (Figure S6 in Supporting Information S1). This suggests that sharp separation between northern- and southern-sourced water can potentially be achieved without a major change in density stratification. In the modern climate, formation of slightly lighter NADW than currently observed, as is the case in our NEMO experiments (Lavergne et al., 2024), may indeed suffice to largely separate dense waters of North Atlantic and Antarctic origins.

Millet, Gray, et al. (2024) analyzed the oxygen-18 isotopic ratio of calcite within shells of benthic foraminifera in sediment cores of the NP. They found a stronger mid-depth maximum of this quasi-conservative tracer at LGM relative to the Holocene, suggesting increased layering of the deep (>2 km) NP at LGM. This stronger mid-depth maximum could come from increased influence of isopycnal diffusion relative to vertical transports in the region’s ventilation (Millet, Gray, et al., 2024). The NA dye tracer distributions simulated in present NEMO experiments corroborate these arguments: they demonstrate that a NP mid-depth maximum in conservative tracers can exist via isopycnal diffusion from the south, notably via isopycnal spreading of tracers contained in mid-depth waters of North Atlantic origin. In addition, the NA dye distributions show that the strength of the mid-depth maximum is sensitive to the balance of isopycnal and diapycnal diffusion in the deep Pacific. We hypothesize that reduced density of NADW relative to AABW, together with reduced vertical transports in the deep NP, may underpin the LGM ventilation patterns inferred by Millet, Gray, et al. (2024).

Finally, we have shown that the filling time of the ocean far exceeds the maximum ideal age of the ocean (by a factor of about two; see also Holzer & Primeau, 2008), and varies strongly between basins. In the modern state, assuming steady circulation and mixing, we estimate that near complete (99%) renewal of the mid-depth (1–3 km) Pacific requires between 3,000 and 4,000 years (Figure 3). Under changing boundary conditions, such as during
the last deglaciation (about 20,000 to 10,000 years ago), the mid-depth Pacific would be expected to adjust with multi-millennial lag and to feature persistent disequilibrium with respect to the surface climate (Gebbie & Huybers, 2019; Wunsch & Heimbach, 2008, 2014). Interpretation of proxy records of past climate changes thus needs to account for spatially varying lags between deep sea states and global surface climate (Skinner & Shackleton, 2005). The very notion of climatic equilibrium or steady state is questionable given the presently inferred adjustment time of the mid-depth Pacific.

**5. Conclusion**

In this study, we compared a global configuration of the NEMO general circulation model to two data-constrained inverse models, focusing on steady-state ventilation patterns and timescales. More precisely, we examined the volumes ventilated by different surface regions and the ideal age distribution of the three model oceans. The three models share several important features. First, we find that the time needed to fill (or renew) 99% of the ocean volume is roughly twice the maximum ideal age. The mid-depth Pacific hosts the maximum age and takes the longest to fill. Given a steady modern circulation state as reconstructed by the inverse models, over 3,000 years are needed to entirely renew the mid-depth Pacific. Second, global volumetric contributions of each source region are fairly similar in all three models. The North Atlantic (north of 30°N) and Antarctic marginal seas each ventilate approximately one-third of the global ocean volume (Figure 4a). About one-fifth of the ocean volume originates in the subantarctic zone where Antarctic Intermediate Water and Subantarctic Mode Water are formed (Figure 4a). In the deep (> 1 km) Pacific Ocean (30°S-65°N), 26%–34% of the water volume comes from the surface North Atlantic and 61%–69% comes from the Southern Ocean surface south of 30°S; here ranges bracket the three model estimates. The lion's share of the deep Pacific ventilation comes from Antarctic marginal seas, which contribute 37%–54%. These volumetric censuses are not a mere curiosity: they carry direct implications for the proportion of unused (preformed) nutrients stored in the deep ocean, and thus for the efficiency of the biological carbon pump (Ito & Follows, 2005; Kwon et al., 2012; Rae & Broecker, 2018).

Major differences between the models also reveal important biases and unknowns. First, ideal ages are consistently higher in the NEMO simulations than in the inverse model states. When state-of-the-art maps of isopycnal (Groeskamp et al., 2020) and diapycnal (Lavergne et al., 2024) mixing are implemented in NEMO (experiment ISO3D), the maximum ideal age in the Pacific reaches 3,300 years, twice the inverse estimates of about 1,600 years constrained by radiocarbon data. This difference occurs despite realistic representation in NEMO of the density stratification in the Pacific. A widespread underestimate of both isopycnal and diapycnal mixing rates in the deep ocean of experiment ISO3D is deemed likely. Nonetheless, sensitivity experiments suggest that plausible increases of diapycnal mixing are unlikely to eliminate the age bias. In point of fact, NEMO appears to misrepresent vertical flows along the northern rim of the Pacific, when these flows seem to play a key role in NP tracer budgets and ventilation rates (Millet, Gray, et al., 2024; Sarmiento et al., 2004; Toggweiler et al., 2019). The drivers and pathways of ventilation in the subarctic Pacific deserve further scrutiny. Second, NEMO experiments exhibit stronger layering of the numerical dyes than the inverse models do. Several factors may underpin this difference: (a) overestimation of diapycnal diffusion by the inverse models, especially by TMI in the deep ocean; (b) production of NADW in a too narrow neutral density range (27.8–27.9 kg m⁻³) in NEMO; and (c) underestimated diapycnal mixing in the ocean interior of NEMO. The latter factor suggests that numerical mixing (i.e., spurious mixing due to the discretization of advection) is sufficiently weak in our NEMO experiments that parameterizing additional (non tidal) sources of deep diapycnal mixing would be beneficial.

Overall, present NEMO experiments, benchmarked against the two observationally constrained inverse model solutions, show more skill in the representation of ventilation patterns than ventilation timescales. Indeed, while dye volumes and distributions compare well across the models on basin scales, ideal ages appear to be substantially overestimated in slowly renewed pools of the global ocean in NEMO. Such age biases are of concern when modeling the Earth system's response to anthropogenic change: they affect not only the long-term response of the ocean but also the simulated base state of marine biogeochemistry, which is tightly dependent on deep-ocean ventilation rates (Eggleston & Galbraith, 2018; Ödalen et al., 2018). In other words, large biases in ideal age (and re-exposure time) will affect the calibration of biogeochemical models and thereby the simulated transient response of marine ecosystems and carbon cycling. Our sensitivity analysis suggests that energetically plausible changes in mixing parameterizations may not completely overcome age biases in this NEMO configuration, pointing to biases in resolved-scale dynamics.
More broadly, the present study outlines a methodology to thoroughly assess ocean ventilation simulated by general circulation models. The employed regional dye tracers (Figure 1) are straightforward to implement and permit the analysis of ventilation characteristics by surface origin on all timescales. The comparison with identical tracers in state-of-the-art inverse models allows bias assessment across the global ocean, bypassing complications from time dependence and sparse sampling. In this study, we used the latest gridded solutions of OCIM and TMI, which assimilate multiple observed tracers (DeVries & Primeau, 2011; Gebbie & Huybers, 2010, 2012), providing a well-constrained observational reference. The computational cost of simulating dye and age tracers to equilibrium in an OGCM can be substantial; here it was mitigated by running these passive tracers in offline mode (i.e., without recomputing ocean and sea ice physics). The cost could also be mitigated by restricting the evaluation to transient ventilation over the first few centuries, rather than over millennia. The present framework is thus relatively versatile and applicable to different models and timescales. It is less suited than Lagrangian particles to analysis of advective trajectories (e.g., Rousselet et al., 2021; Sebille et al., 2018; Styles et al., 2024; Tamsitt et al., 2017) but better suited to ground truthing of general circulation models with parameterized diffusive transports. Thus, the approach based on dye and age tracers developed here could complement traditional evaluations of simulated ocean transports and prove helpful for the calibration and evaluation of Earth system models.

## Data Availability Statement

Nucleus for European Modeling of the Ocean code is available at https://forge.nemo-ocean.eu/nemo/nemo/-/releases/ (NEMO System Team, 2022b, 2022a). The employed version is 4.2.0. Total Matrix Intercomparison output was generated from scripts and data sets available at https://github.com/ggebbie/TMI.mat (Gebbie & Huybers, 2010). The OCIM2-48L base output employed here is available in the Figshare database under accession code https://doi.org/10.6084/m9.figshare.14802732 (DeVries, 2021). The code used for generating the figures and performing the analysis is available and can be accessed through the following GitHub repository: https://github.com/brumillet/dyes-evaluation_ventilation/. Steady-state distributions of regional dye tracers and ideal age from each of the three models are provided at: https://doi.org/10.17882/103915 (Millet, Lavergne, et al., 2024), and the transient evolution of dye tracers in the three models at: https://doi.org/10.17882/106965 (Millet et al., 2025). The ocean reanalysis outputs used in Figure 16 is available at https://doi.org/10.48670/moi-00185 (winds, Copernicus Marine Service, 2021a) and at https://doi.org/10.48670/moi-00021 (vertical velocities, Copernicus Marine Service, 2021b).

## Acknowledgments

The authors thank X. Davila, E. McDonagh, D. Ferreira for comments and B. Gulk for help with analysis. We also thank the anonymous reviewers that evaluated the previous version of this manuscript and helped shape it to its final version. This work was granted access to the HPC resources of IDRIS under the allocations 2022-A0120107451 and 2023-A0140107451 made by GENCI. This research was supported by Ocean Cryosphere Exchanges in ANTarctica: Impacts on Climate and the Earth system, OCEAN ICE, which is funded by the European Union, Horizon Europe Funding Programme for research and innovation under grant agreement Nr. 101060452, 10.3030/101060452. Financial support was also provided by the French National LEFE programme through the ROOF project and a thesis grant from the French Alternative Energies and Atomic Energy Commission (CEA).

## References

Adcroft, A., Anderson, W., Balaji, V., Blanton, C., Bushuk, M., Dufour, C. O., et al. (2019). The GFDL global ocean and sea ice model OM4.0: Model description and simulation features. *Journal of Advances in Modeling Earth Systems*, 11(10), 3167–3211. https://doi.org/10.1029/2019MS001726

Adcroft, A., Hill, C., & Marshall, J. (1997). Representation of topography by shaved cells in a height coordinate ocean model. *Monthly Weather Review*, 125(9), 2293–2315. https://doi.org/10.1175/1520-0493(1997)125(2293:ROTBSC)2.0:CO2

Adkins, J. F., McIntyre, K., & Schrag, D. P. (2002). The salinity, temperature, and δ¹⁸O of the glacial deep ocean. *Science*, 298(5599), 1769–1773. https://doi.org/10.1126/science.1076252

Alford, M. H. (2020). Global calculations of local and remote near-inertial-wave dissipation. *Journal of Physical Oceanography*, 50(11), 3157–3164. https://doi.org/10.1175/JPO-D-20-0106.1

Baker, J. A., Bell, M. J., Jackson, L. C., Renshaw, R., Vallis, G. K., Watson, A. J., & Wood, R. A. (2023). Overturning pathways control AMOC weakening in CMIP6 models. *Geophysical Research Letters*, 50(14), e2023GL103381. https://doi.org/10.1029/2023GL103381

Bouttes, N., Paillard, D., & Roche, D. M. (2010). Impact of brine-induced stratification on the glacial carbon cycle. *Climate of the Past*, 6(5), 575–589. https://doi.org/10.5194/cp-6-575-2010

Broecker, W. S. (1991). The great ocean conveyor. *Oceanography*, 4(2), 79–89. https://doi.org/10.5670/oceanog.1991.07

Broecker, W. S., & Olson, E. A. (1960). Radiocarbon from nuclear tests, II. *Science*, 132(3429), 712–721. https://doi.org/10.1126/science.132.3429.712

Burchard, H., & Rennau, H. (2008). Comparative quantification of physically and numerically induced mixing in ocean models. *Ocean Modelling*, 20(3), 293–311. https://doi.org/10.1016/j.ocemod.2007.10.003

Burke, A., Stewart, A. L., Adkins, J. F., Ferrari, R., Jansen, M. F., & Thompson, A. F. (2015). The glacial mid-depth radiocarbon bulge and its implications for the overturning circulation. *Paleoceanography*, 30(7), 1021–1039. https://doi.org/10.1002/2015PA002778

Capó, E., McWilliams, J. C., Gula, J., Molemaker, M. J., Damien, P., & Schubert, R. (2024). Abyssal slope currents. *Journal of Physical Oceanography*, 54(11), 2373–2392. https://doi.org/10.1175/JPO-D-24-0028.1

Chapman, D. C. (1989). Enhanced subinertial diurnal tides over isolated topographic features. *Deep-Sea Research, Part A: Oceanographic Research Papers*, 36(6), 815–824. https://doi.org/10.1016/0198-0149(89)90030-7

Chen, H., Haumann, F. A., Talley, L. D., Johnson, K. S., & Sarmiento, J. L. (2022). The Deep Ocean’s carbon exhaust. *Global Biogeochemical Cycles*, 36(7), e2021GB007156. https://doi.org/10.1029/2021GB007156

Cimoli, L., Gebbie, G., Purkey, S. G., & Smethie, W. M. (2023). Annually resolved propagation of cfc5 and sf6 in the global ocean over eight decades. *Journal of Geophysical Research: Oceans*, 128(3), e2022JC019337. https://doi.org/10.1029/2022JC019337
Copernicus Marine Service. (2021a). Global ocean hourly reprocessed sea surface wind and stress from scatterometer and model. *Copernicus Marine Environment Monitoring Service*. https://doi.org/10.48670/moi-00185

Copernicus Marine Service. (2021b). Global ocean physics reanalysis. *Copernicus Marine Environment Monitoring Service*. https://doi.org/10.48670/moi-00021

Copernicus Marine Service. (2024). Global ocean hourly reprocessed sea surface wind and stress from scatterometer and model [Computer software manual]. *E.U. Copernicus Marine Service Information (CMEMS), Marine Data Store (MDS)*. https://doi.org/10.48670/moi-00185

Danabasoglu, G., Yeager, S. G., Bailey, D., Behrens, E., Bentsen, M., Bi, D., et al. (2014). North Atlantic simulations in coordinated ocean-ice reference experiments phase ii (core-ii). Part I: Mean states. *Ocean Modelling*, **73**, 76–107. https://doi.org/10.1016/j.ocemod.2013.10.005

Davila, X., Olsen, A., Lauvset, S. K., McDonagh, E. L., Brakstad, A., & Gebbie, G. (2023). On the origins of open ocean oxygen minimum zones. *Journal of Geophysical Research: Oceans*, **128**(8), e2023JC019677. https://doi.org/10.1029/2023JC019677

Debreu, L., Kevlahan, N. K. R., & Marchesiello, P. (2020). Brinkman volume penalization for bathymetry in three-dimensional ocean models. *Ocean Modelling*, **145**, 101530. https://doi.org/10.1016/j.ocemod.2019.101530

de Lavergne, C., Madec, G., Le Sommer, J., Nurser, A. J. G., & Naveira Garabato, A. C. (2016). On the consumption of Antarctic bottom water in the abyssal ocean. *Journal of Physical Oceanography*, **46**(2), 635–661. https://doi.org/10.1175/JPO-D-14-0201.1

de Lavergne, C., Madec, G., Roquet, F., Holmes, R. M., & McDougall, T. J. (2017). Abyssal ocean overturning shaped by seafloor distribution. *Nature*, **551**(7679), 181–186. https://doi.org/10.1038/nature24472

de Lavergne, C., Rathore, S., Madec, G., Sallée, J.-B., Ethé, C., Nasser, A., et al. (2024). Effects of improved tidal mixing in NEMO one-degree global ocean model. *Journal of Advances in Modeling Earth Systems*. https://doi.org/10.22541/essoar.173152139.95978362/v1

de Lavergne, C., Vic, C., Madec, G., Roquet, F., Waterhouse, A. F., Whalen, C. B., et al. (2020). A parameterization of local and remote tidal mixing. *Journal of Advances in Modeling Earth Systems*, **12**(5), e2020MS002065. https://doi.org/10.1029/2020MS002065

DeVries, T. (2014). The oceanic anthropogenic co2 sink: Storage, air-sea fluxes, and transports over the industrial era. *Global Biogeochemical Cycles*, **28**(7), 631–647. https://doi.org/10.1002/2013GB004739

DeVries, T. (2021). OCIM2-48L base state model output Dataset. Figshare. Retrieved from https://figshare.com/articles/dataset/OCIM2-48L_base_state_model_output/14802732

DeVries, T., & Holzer, M. (2019). Radiocarbon and helium isotope constraints on deep ocean ventilation and mantle-3he sources. *Journal of Geophysical Research: Oceans*, **124**(5), 3036–3057. https://doi.org/10.1029/2018JC014716

DeVries, T., & Primeau, F. (2011). Dynamically and observationally constrained estimates of water-mass distributions and ages in the global ocean. *Journal of Physical Oceanography*, **41**(12), 2381–2401. https://doi.org/10.1175/JPO-D-10-0501.1

Dunne, J. P., Horowitz, L. W., Adcroft, A. J., Ginoux, P., Held, I. M., John, J. G., et al. (2020). The GFDL earth system model version 4.1 (GFDL-ESM 4.1): Overall coupled model description and simulation characteristics. *Journal of Advances in Modeling Earth Systems*, **12**(11), e2019MS002015. https://doi.org/10.1029/2019MS002015

Dutay, J.-C., Bullister, J., Doney, S., Orr, J., Najjar, R., Caldeira, K., et al. (2002). Evaluation of ocean model ventilation with cfc-11: Comparison of 13 global ocean models. *Ocean Modelling*, **4**(2), 89–120. https://doi.org/10.1016/S1463-5003(01)00013-0

Eggleston, S., & Galbraith, E. D. (2018). The devil’s in the disequilibrium: Multi-component analysis of dissolved carbon and oxygen changes under a broad range of forcings in a general circulation model. *Biogeosciences*, **15**(12), 3761–3777. https://doi.org/10.5194/bg-15-3761-2018

England, M. H. (1995a). The age of water and ventilation timescales in a global ocean model. *Journal of Physical Oceanography*, **25**(11), 2756–2777. https://doi.org/10.1175/1520-0485(1995)025(2756:TAOWAV)2.0.CO;2

England, M. H. (1995b). Using chlorofluorocarbons to assess ocean climate models. *Geophysical Research Letters*, **22**(22), 3051–3054. https://doi.org/10.1029/95GL02670

Evans, D. G., Frajka-Williams, E., & Naveira Garabato, A. C. (2022). Dissipation of mesoscale eddies at a western boundary via a direct energy cascade. *Scientific Reports*, **12**(1), 887. https://doi.org/10.1038/s41598-022-05002-7

Ferrari, R., Jansen, M. F., Adkins, J. F., Burke, A., Stewart, A. L., & Thompson, A. F. (2014). Antarctic sea ice control on ocean circulation in present and glacial climates. *Proceedings of the National Academy of Sciences*, **111**(24), 8753–8758. https://doi.org/10.1073/pnas.1323922211

Ferrari, R., Mashayek, A., McDougall, T. J., Nikurashin, M., & Campin, J.-M. (2016l). Turning Ocean mixing upside down. *Journal of Physical Oceanography*, **46**(7), 2239–2261. https://doi.org/10.1175/JPO-D-15-0244.1

Ferrari, R., & Nikurashin, M. (2010). Suppression of eddy diffusivity across jets in the Southern Ocean. *Journal of Physical Oceanography*, **40**(7), 1501–1519. https://doi.org/10.1175/2010JPO4278.1

Fine, R. A. (2011). Observations of CFCs and SF6 as ocean tracers. *Annual Review of Marine Science*, **3**(3), 173–195. https://doi.org/10.1146/annurev.marine.010908.163933

Friedli, H., Lötscher, H., Oeschger, H., Siegenthaler, U., & Stauffer, B. (1986). Ice core record of the 13c/12c ratio of atmospheric CO₂ in the past two centuries. *Nature*, **324**(6094), 237–238. https://doi.org/10.1038/324237a0

Fu, Y., Li, F., Karstensen, J., & Wang, C. (2020). A stable Atlantic meridional overturning circulation in a changing north Atlantic ocean since the 1990s. *Science Advances*, **6**(48), eabc7836. https://doi.org/10.1126/sciadv.abc7836

Galbraith, E., & de Lavergne, C. (2019). Response of a comprehensive climate model to a broad range of external forcings: Relevance for deep ocean ventilation and the development of late cenozoic ice ages. *Climate Dynamics*, **52**(1), 653–679. https://doi.org/10.1007/s00382-018-4157-8

Gebbie, G. (2014). How much did glacial north Atlantic water shoal? *Paleoceanography*, **29**(3), 190–209. https://doi.org/10.1002/2013PA002557

Gebbie, G., & Huybers, P. (2010). Total matrix intercomparison: A method for determining the geometry of water-mass pathways. *Journal of Physical Oceanography*, **40**(8), 1710–1728. https://doi.org/10.1175/2010JPO4272.1

Gebbie, G., & Huybers, P. (2011). How is the ocean filled? *Geophysical Research Letters*, **38**(6), L06604. https://doi.org/10.1029/2011GL046769

Gebbie, G., & Huybers, P. (2012). The mean age of ocean waters inferred from radiocarbon observations: Sensitivity to surface sources and accounting for mixing histories. *Journal of Physical Oceanography*, **42**(2), 291–305. https://doi.org/10.1175/JPO-D-11-043.1

Gebbie, G., & Huybers, P. (2019). The little ice age and 20th-century deep pacific cooling. *Science*, **363**(6422), 70–74. https://doi.org/10.1126/science.aar8413

Gent, P., & McWilliams, J. (1990). Isopycnal mixing in ocean circulation models. *Journal of Physical Oceanography*, **20**(1), 150–155. https://doi.org/10.1175/1520-0485(1990)020(0150:IMIOCM)2.0.CO;2

Gordon, A. L. (1986). Interocean exchange of thermocline water. *Journal of Geophysical Research*, **91**(C4), 5037–5046. https://doi.org/10.1029/JC091iC04p05037

Gouretski, V. (2019). A new global ocean hydrographic climatology. *Atmospheric and Oceanic Science Letters*, **12**(3), 226–229. https://doi.org/10.1080/16742834.2019.1588066
Gray, W. R., de Lavergne, C., Jnglin Wills, R. C., Menviel, L., Spence, P., Holzer, M., et al. (2023). Poleward shift in the southern hemisphere westerly winds synchronous with the deglacial rise in CO₂. *Paleoceanography and Paleoclimatology, 38*(7), e2023PA004666. https://doi.org/10.1029/2023PA004666

Gregg, M. C. (1987). Diapycnal mixing in the thermocline: A review. *Journal of Geophysical Research*, 92(C5), 5249–5286. https://doi.org/10.1029/JC092iC05p05249

Gregory, J. M., Bloch-Johnson, J., Couldrey, M. P., Exarchou, E., Griffies, S. M., Kuhlbrodt, T., et al. (2024). A new conceptual model of global ocean heat uptake. *Climate Dynamics*, 62(3), 1669–1713. https://doi.org/10.1007/s00382-023-06989-z

Groeskamp, S., LaCasce, J. H., McDougall, T. J., & Rogé, M. (2020). Full-depth global estimates of ocean mesoscale eddy mixing from observations and theory. *Geophysical Research Letters*, 47(18), e2020GL089425. https://doi.org/10.1029/2020GL089425

Groeskamp, S., Sloyan, B. M., Zika, J. D., & McDougall, T. J. (2017). Mixing inferred from an ocean climatology and surface fluxes. *Journal of Physical Oceanography*, 47(3), 667–687. https://doi.org/10.1175/JPO-D-16-0125.1

Gula, J., Molemaker, M. J., & McWilliams, J. C. (2016). Topographic generation of submesoscale centrifugal instability and energy dissipation. *Nature Communications*, 7(1), 12811. https://doi.org/10.1038/ncomms12811

Gupta, M., Williams, R. G., Lauderdale, J. M., Jahn, O., Hill, C., Dutkiewicz, S., & Follows, M. J. (2022). A nutrient relay sustains subtropical ocean productivity. *Proceedings of the National Academy of Sciences*, 119(41), e2206504119. https://doi.org/10.1073/pnas.2206504119

Hautala, S. L. (2018). The abyssal and deep circulation of the Northeast Pacific Basin. *Progress in Oceanography*, 160, 68–82. https://doi.org/10.1016/j.pocean.2017.11.011

Hernández-Guerra, A., & Talley, L. D. (2016). Meridional overturning transports at 30°s in the Indian and pacific oceans in 2002–2003 and 2009. *Progress in Oceanography*, 146, 89–120. https://doi.org/10.1016/j.pocean.2016.06.005

Hieronymus, M., Nycander, J., Nilsson, J., Döös, K., & Hallberg, R. (2019). Oceanic overturning and heat transport: The role of background diffusivity. *Journal of Climate*, 32(3), 701–716. https://doi.org/10.1175/JCLI-D-18-0438.1

Holmes, R. M., Zika, J. D., Griffies, S. M., Hogg, A. M., Kiss, A. E., & England, M. H. (2021). The geography of numerical mixing in a suite of global ocean models. *Journal of Advances in Modeling Earth Systems*, 13(7), e2020MS002333. https://doi.org/10.1029/2020MS002333

Holzer, M. (2022). The fate of oxygen in the ocean and its sensitivity to local changes in biological production. *Journal of Geophysical Research: Oceans*, 127(8), e2022JC018802. https://doi.org/10.1029/2022JC018802

Holzer, M., Chamberlain, M. A., & Matar, R. J. (2020). Climate-driven changes in the Ocean's ventilation pathways and time scales diagnosed from transport matrices. *Journal of Geophysical Research: Oceans*, 125(10), e2020JC016414. https://doi.org/10.1029/2020JC016414

Holzer, M., DeVries, T., & de Lavergne, C. (2021). Diffusion controls the ventilation of a pacific shadow zone above abyssal overturning. *Nature Communications*, 12(1), 4348. https://doi.org/10.1038/s41467-021-24648-x

Holzer, M., & Hall, T. M. (2000). Transit-time and tracer-age distributions in geophysical flows. *Journal of the Atmospheric Sciences*, 57(21), 3539–3558. https://doi.org/10.1175/1520-0469(2000)057(3539:TTATAD)2.0:CO;2

Holzer, M., & Primeau, F. W. (2006). The diffusive ocean conveyor. *Geophysical Research Letters*, 33(14), L14618. https://doi.org/10.1029/2006GL026232

Holzer, M., & Primeau, F. W. (2008). The path-density distribution of oceanic surface-to-surface transport. *Journal of Geophysical Research*, 113(C1), C01018. https://doi.org/10.1029/2006JC003976

Holzer, M., Primeau, F. W., Smethie, W. M., Jr., & Khatiwala, S. (2010). Where and how long ago was water in the western North Atlantic ventilated? Maximum entropy inversions of bottle data from WOCE line A20. *Journal of Geophysical Research*, 115(C7), C07005. https://doi.org/10.1029/2009JC005750

Ito, T., & Follows, M. J. (2005). Preformed phosphate, soft tissue pump and atmospheric CO₂. *Journal of Marine Research*, 63(4), 813–839. https://doi.org/10.1357/0022240054663231

Jansen, M. F. (2017). Glacial ocean circulation and stratification explained by reduced atmospheric temperature. *Proceedings of the National Academy of Sciences*, 114(1), 45–50. https://doi.org/10.1073/pnas.1610438113

Jansen, M. F., Nadeau, L.-P., & Merlis, T. M. (2018). Transient versus equilibrium response of the ocean's overturning circulation to warming. *Journal of Climate*, 31(13), 5147–5163. https://doi.org/10.1175/JCLI-D-17-0797.1

Johnson, G. C. (2008). Quantifying antarctic bottom water and north Atlantic deep water volumes. *Journal of Geophysical Research*, 113(C5), C05027. https://doi.org/10.1029/2007JC004477

Jones, C. S., & Abernathey, R. P. (2019). Isopycnal mixing controls deep ocean ventilation. *Geophysical Research Letters*, 46(22), 13144–13151. https://doi.org/10.1029/2019GL085208

Kawasaki, T., Hasumi, H., & Tanaka, Y. (2021). Role of tide-induced vertical mixing in the deep Pacific Ocean circulation. *Journal of Oceanography*, 77(2), 173–184. https://doi.org/10.1007/s10872-020-00584-0

Key, R. M., Kozyr, A., Sabine, C., Lee, K., Wanninkhof, R., Bullister, J., et al. (2004). A global ocean carbon climatology: Results from global data analysis project (glodap). *Global Biogeochemical Cycles*, 18(4), GB4031. https://doi.org/10.1029/2004GB002247

Kiss, A. E., Hogg, A. M., Hannah, N., Boeira Dias, F., Brassington, G. B., Chamberlain, M. A., et al. (2020). ACCESS-OM2 v1.0: A global ocean–sea ice model at three resolutions. *Geoscientific Model Development*, 13(2), 401–442. https://doi.org/10.5194/gmd-13-401-2020

Klymak, J. M. (2018). Non-propagating form drag and turbulence due to stratified flow over large-scale abyssal hill topography. *Journal of Physical Oceanography*, 48(10), 2383–2395. https://doi.org/10.1175/JPO-D-17-0225.1

Kwon, E. Y., Hain, M. P., Sigman, D. M., Galbraith, E. D., Sarmiento, J. L., & Toggweiler, J. R. (2012). North atlantic ventilation of “southern-sourced” deep water in the glacial ocean. *Paleoceanography*, 27(2). https://doi.org/10.1029/2011PA002211

Kwon, E. Y., Sarmiento, J. L., Toggweiler, J. R., & DeVries, T. (2011). The control of atmospheric pco2 by ocean ventilation change: The effect of the oceanic storage of biogenic carbon. *Global Biogeochemical Cycles*, 25(3). https://doi.org/10.1029/2011GB004059

Large, W. G., & Yeager, S. G. (2009). The global climatology of an interannually varying air–sea flux data set. *Climate Dynamics*, 33(2), 341–364. https://doi.org/10.1007/s00382-008-0441-3

Lester, J. G., Graven, H. D., Khatiwala, S., & McNichol, A. P. (2024). Changes in oceanic radiocarbon and CFCs since the 1990s. *Journal of Geophysical Research: Oceans*, 129(7), e2023JC020387. https://doi.org/10.1029/2023JC020387

Lévy, M., Resplandy, L., Palter, J. B., Couespel, D., & Lachkar, Z. (2022). Chapter 13 - the crucial contribution of mixing to present and future ocean oxygen distribution. In M. Meredith & A. Naveira Garabato (Eds.), *Ocean mixing* (pp. 329–344). Elsevier. https://doi.org/10.1016/B978-0-12-821512-8.00020-7

Locarnini, R. A., Mishonov, A. V., Antonov, J. I., Boyer, T. P., Garcia, H. E., Baranova, O. K., et al. (2013). World ocean atlas 2013. *Temperature*, 1(73), 40.

Lozier, M. S. (2010). Deconstructing the conveyor belt. *Science*, 328(5985), 1507–1511. https://doi.org/10.1126/science.1189250

Lucazeau, F. (2019). Analysis and mapping of an updated terrestrial heat flow data set. *Geochemistry, Geophysics, Geosystems*, 20(8), 4001–4024. https://doi.org/10.1029/2019GC008389
Luyten, J. R., Pedlosky, J., & Stommel, H. (1983). The ventilated thermocline. *Journal of Physical Oceanography*, 13(2), 292–309. https://doi.org/10.1175/1520-0485(1983)013(0292:TVT)2:0.CO;2

Marzocchi, A., Nurser, A. J. G., Clément, L., & McDonagh, E. L. (2021). Surface atmospheric forcing as the driver of long-term pathways and timescales of ocean ventilation. *Ocean Science*, 17(4), 935–952. https://doi.org/10.5194/os-17-935-2021

McDougall, T. J., & Ferrari, R. (2017). Abyssal upwelling and downwelling driven by near-boundary mixing. *Journal of Physical Oceanography*, 47(2), 261–283. https://doi.org/10.1175/JPO-D-16-0082.1

Millet, B., de Lavergne, C., Gray, W. R., Ethé, C., Madec, G., Holzer, M., et al. (2024a). Water-mass fractions and ideal age from three models of ocean transport [Dataset]. *SEANOE*. https://doi.org/10.17882/103915

Millet, B., Divel, H., de Lavergne, C., Gebbie, G., & Holzer, M. (2025). Propagation of regional dye tracers in three models of ocean transport over 6000 years [Dataset]. https://doi.org/10.17882/106965

Millet, B., Gray, W. R., de Lavergne, C., & Roche, D. M. (2024b). Oxygen isotope constraints on the ventilation of the modern and glacial pacific. *Climate Dynamics*, 62(1), 649–664. https://doi.org/10.1007/s00382-02-06910-8

Müller, M. (2013). On the space- and time-dependence of barotropic-to-baroclinic tidal energy conversion. *Ocean Modelling*, 72, 242–252. https://doi.org/10.1016/j.oocemod.2013.09.007

Naveira Garabato, A. C., Frajka-Williams, E. E., Spingys, C. P., Legg, S., Polzin, K. L., Forryan, A., et al. (2019). Rapid mixing and exchange of deep-ocean waters in an abyssal boundary current. *Proceedings of the National Academy of Sciences*, 116(27), 13233–13238. https://doi.org/10.1073/pnas.1904087116

Naveira Garabato, A. C., MacGilchrist, G. A., Brown, P. J., Evans, D. G., Meijers, A. J. S., & Zika, J. D. (2017). High-latitude ocean ventilation and its role in earth's climate transitions. *Philosophical transactions. Series A, Mathematical, physical, and engineering sciences*, 375(2102), 20160324. https://doi.org/10.1098/rsta.2016.0324

NEMO System Team. (2022a). Nemo ocean engine (Computer software manual No. 27). Zenodo. https://doi.org/10.5281/zenodo.6334656

NEMO System Team. (2022b). NEMO ocean engine – Version 4.2. Zenodo. https://doi.org/10.5281/zenodo.7328870

Ödalen, M., Nycander, J., Oliver, K. I. C., Brodeau, L., & Ridgwell, A. (2018). The influence of the ocean circulation state on ocean carbon storage and CO₂ drawdown potential in an earth system model. *Biogeosciences*, 15(5), 1367–1393. https://doi.org/10.5194/bg-15-1367-2018

Oka, A., & Niwa, Y. (2013). Pacific deep circulation and ventilation controlled by tidal mixing away from the sea bottom. *Nature Communications*, 4(1), 2419. https://doi.org/10.1038/ncomms3419

Olsen, A., Key, R. M., van Heuven, S., Lauvset, S. K., Velo, A., Lin, X., et al. (2016). The global ocean data analysis project version 2 (glodapv2) – An internally consistent data product for the World Ocean. *Earth System Science Data*, 8(2), 297–323. https://doi.org/10.5194/essd-8-297-2016

Orr, J. C., Najjar, R. G., Aumont, O., Bopp, L., Bullister, J. L., Danabasoglu, G., et al. (2017). Biogeochemical protocols and diagnostics for the CMIP6 ocean model intercomparison project (OMIP). *Geoscientific Model Development*, 10(6), 2169–2199. https://doi.org/10.5194/gmd-10-2169-2017

Orsi, A. H., Jacobs, S. S., Gordon, A. L., & Visbeck, M. (2001). Cooling and ventilating the abyssal ocean. *Geophysical Research Letters*, 28(15), 2923–2926. https://doi.org/10.1029/2001GL012830

Pradal, M.-A., & Gnanadesikan, A. (2014). How does the redi parameter for mesoscale mixing impact global climate in an earth system model? *Journal of Advances in Modeling Earth Systems*, 6(3), 586–601. https://doi.org/10.1002/2013MS000273

Primeau, F. (2005). Characterizing transport between the surface mixed layer and the ocean interior with a forward and adjoint global ocean transport model. *Journal of Physical Oceanography*, 35(4), 545–564. https://doi.org/10.1175/JPO2699.1

Primeau, F., & Holzer, M. (2006). The Ocean's memory of the atmosphere: Residence-time and ventilation-rate distributions of water masses. *Journal of Physical Oceanography*, 36(7), 1439–1456. https://doi.org/10.1175/JPO2919.1

Rae, J. W. B., & Broecker, W. (2018). What fraction of the pacific and indian oceans' deep water is formed in the southern ocean? *Biogeosciences*, 15(12), 3779–3794. https://doi.org/10.5194/bg-15-3779-2018

Richardson, P. L. (2008). On the history of meridional overturning circulation schematic diagrams. *Progress in Oceanography*, 76(4), 466–486. https://doi.org/10.1016/j.pocean.2008.01.005

Rousselet, L., Cessi, P., & Forget, G. (2021). Coupling of the mid-depth and abyssal components of the global overturning circulation according to a state estimate. *Science Advances*, 7(21), eabf5478. https://doi.org/10.1126/sciadv.abf5478

Rubin, S. I., & Key, R. M. (2002). Separating natural and bomb-produced radiocarbon in the ocean: The potential alkalinity method. *Global Biogeochemical Cycles*, 16(4), 52–1–52–19. https://doi.org/10.1029/2001GB001432

Samelson, R., & Vallis, G. K. (1997). Large-scale circulation with small diapycnal diffusion: The two-thermocline limit. *Journal of Marine Research*, 55(2), 223–275. https://doi.org/10.1357/0022240973224382

Sarmiento, J. L., Gruber, N., Brzezinski, M. A., & Dunne, J. P. (2004). 01). High-latitude controls of thermocline nutrients and low latitude biological productivity. *Nature*, 427(6969), 56–60. https://doi.org/10.1038/nature02127

Schlitzer, R. (2007). Assimilation of radiocarbon and chlorofluorocarbon data to constrain deep and bottom water transports in the world ocean. *Journal of Physical Oceanography*, 37(2), 259–276. https://doi.org/10.1175/JPO3011.1

Skinner, L., & Shackleton, N. (2005). An Atlantic lead over pacific deep-water change across termination i: Implications for the application of the marine isotope stage stratigraphy. *Quaternary Science Reviews*, 24(5), 571–580. https://doi.org/10.1016/j.quascirev.2004.11.008

Sloyan, B. M., & Rintoul, S. R. (2001). The southern ocean limb of the global deep overturning circulation. *Journal of Physical Oceanography*, 31(1), 143–173. https://doi.org/10.1175/1520-0485(2001)031(0143:TSOLOT)2.0.CO;2

Styles, A. F., MacGilchrist, G. A., Bell, M. J., & Marshall, D. P. (2024). Spatial and temporal patterns of southern ocean ventilation. *Geophysical Research Letters*, 51(4), e2023GL106716. https://doi.org/10.1029/2023GL106716

Sweeney, C., Gloor, E., Jacobson, A. R., Key, R. M., McKinley, G., Sarmiento, J. L., & Wanninkhof, R. (2007). Constraining global air-sea gas exchange for CO₂ with recent bomb 14C measurements. *Global Biogeochemical Cycles*, 21(2), https://doi.org/10.1029/2006GB002784

Talley, L. D. (2013). Closure of the global overturning circulation through the Indian, pacific, and southern oceans: Schematics and transports. *Oceanography*, 26(1), 80–97. https://doi.org/10.5670/oceanog.2013.07

Talley, L. D., Feely, R. A., Sloyan, B. M., Wanninkhof, R., Baringer, M. O., Bullister, J. L., et al. (2016). Changes in ocean heat, carbon content, and ventilation: A review of the first decade of GO-SHIP global repeat hydrography. *Annual Review of Marine Science*, 8(1), 185–215. https://doi.org/10.1146/annurev-marine-052915-100829

Tamsitt, V., Drake, H. F., Morrison, A. K., Talley, L. D., Dufour, C. O., Gray, A. R., et al. (2017). Spiraling pathways of global deep waters to the surface of the southern ocean. *Nature Communications*, 8(1), 172. https://doi.org/10.1038/s41467-017-00197-0

Tanaka, Y. (2023). Energy conversion rate from subinertial surface tides to internal tides. *Journal of Physical Oceanography*, 53(5), 1355–1374. https://doi.org/10.1175/JPO-D-22-0201.1
Toggweiler, J. R., Druffel, E. R. M., Key, R. M., & Galbraith, E. D. (2019). Upwelling in the ocean basins north of the acc: 1. On the upwelling exposed by the surface distribution of δ14c. *Journal of Geophysical Research: Oceans*, 124(4), 2591–2608. https://doi.org/10.1029/2018JC014794

Urkawa, L. S., & Hasumi, H. (2012). Eddy-resolving model estimate of the cabbeling effect on the water mass transformation in the southern ocean. *Journal of Physical Oceanography*, 42(8), 1288–1302. https://doi.org/10.1175/JPO-D-11-0173.1

van Sebille, E., Griffies, S. M., Abernathey, R., Adams, T. P., Berloff, P., Biastoch, A., et al. (2018). Lagrangian ocean analysis: Fundamentals and practices. *Ocean Modelling*, 121, 49–75. https://doi.org/10.1016/j.ocemod.2017.11.008

Warren, B. A. (1983). Why is no deep water formed in the north pacific? *Journals of Marine Research*, 41(2), 327–347. https://doi.org/10.1357/002224083788520207

Worthington, L. (1981). The water masses of the world ocean: Some results of a fine-scale census. *Evolution of Physical Oceanography, Scientific Surveys in Honor of Henry Stommel*, 42–60. Retrieved from https://ocw.mit.edu/ans7870/textbooks/Wunsch/Edited/Chapter2.pdf

Wunsch, C. (2002). Oceanic age and transient tracers: Analytical and numerical solutions. *Journal of Geophysical Research*, 107(C6), 1–18. https://doi.org/10.1029/2001JC000797

Wunsch, C., & Heimbach, P. (2008). How long to oceanic tracer and proxy equilibrium? *Quaternary Science Reviews*, 27(7), 637–651. https://doi.org/10.1016/j.quascirev.2008.01.006

Wunsch, C., & Heimbach, P. (2014). Bidecadal thermal changes in the abyssal ocean. *Journal of Physical Oceanography*, 44(8), 2013–2030. https://doi.org/10.1175/JPO-D-13-096.1

Yang, L., Nikurashin, M., Hogg, A. M., & Sloyan, B. M. (2021). The impact of lee waves on the southern ocean circulation. *Journal of Physical Oceanography*, 51(9), 2933–2950. https://doi.org/10.1175/JPO-D-20-0263.1

You, Y. (2003). The pathway and circulation of north pacific intermediate water. *Geophysical Research Letters*, 30(24), 2291. https://doi.org/10.1029/2003GL018561

Zweng, M. M., Reagan, J. R., Antonov, J. I., Locarnini, R. A., Mishonov, A. V., Boyer, T. P., et al. (2013). World ocean atlas 2013. *Salinity*, 2(74), 39.