# Changes in Global Ocean Circulation due to Isopycnal
Diffusion

Ashwita Chouksey, Alexa Griesel, Manitoba Chouksey, Carsten Eden

► To cite this version:

Ashwita Chouksey, Alexa Griesel, Manitoba Chouksey, Carsten Eden. Changes in Global Ocean Circulation due to Isopycnal Diffusion. Journal of Physical Oceanography, 2022, 52 (9), pp.2219-2235. 10.1175/JPO-D-21-0205.1 . hal-05039819

HAL Id: hal-05039819

https://hal.univ-brest.fr/hal-05039819v1

Submitted on 19 Apr 2025

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

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

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

Distributed under a Creative Commons Attribution 4.0 International License
# Changes in Global Ocean Circulation due to Isopycnal Diffusion

ASHWITA CHOUKSEY,ª,ᵇ ALEXA GRIESEL,ª MANITA CHOUKSEY,ª,ᶜ AND CARSTEN EDENª

(Manuscript received 23 September 2021, in final form 3 June 2022)

**ABSTRACT:** We investigate changes in the ocean circulation due to the variation of isopycnal diffusivity ($\kappa_{iso}$) in a global non-eddy-resolving model. Although isopycnal diffusion is thought to have minor effects on interior density gradients, the model circulation shows a surprisingly large sensitivity to the changes: with increasing $\kappa_{iso}$, the strength of the Atlantic residual overturning circulation (AMOC) and the Antarctic Circumpolar Current (ACC) transport weaken. At high latitudes, the isopycnal diffusion diffuses temperature and salinity upward and poleward, and at low latitudes downward close to the surface. Increasing isopycnal diffusivity increases the meridional isopycnal fluxes whose meridional gradient is equatorward, hence leading to a negative contribution to the flux divergence in the tracer equations and predominant cooling and freshening equatorward of 40°. The effect on temperature overcompensates the countering effect of salinity diffusion, such that the meridional density differences decrease, along with which ACC and AMOC decrease. We diagnose the adjustment process to the new equilibrium with increased isopycnal diffusion to assess how the other terms in the tracer equations react to the increased $\kappa_{iso}$. It reveals that around $\pm40°$ latitude, the cooling induced by the increased isopycnal flux is only partly compensated by warming by advection, explaining the net cooling. Overall, the results emphasize the importance of isopycnal diffusion on ocean circulation and dynamics, and hence the necessity of its careful representation in models.

**SIGNIFICANCE STATEMENT:** The effect of mixing by mesoscale eddies, represented as diffusion along surfaces of constant density in models, on the ocean circulation is not well understood. Here, we show that an increase in the eddy diffusivity in different setups of a global ocean model leads to a surprisingly large change of the ocean circulation. The strength of the Atlantic overturning circulation and the Antarctic Circumpolar Current decrease. We find that the interior ocean becomes cooler and fresher and that the temperature effect on density dominates over salinity, resulting in a decrease in the density gradients. Our results point out the importance of eddy diffusion on ocean circulation, and hence the necessity of its correct representation in ocean and climate models.

**KEYWORDS:** Meridional overturning circulation; Ocean circulation; Isopycnal mixing

## 1. Introduction

The large-scale global ocean circulation is driven by winds, tidal forces, and density mixing induced to a large extent by the breaking of internal gravity waves. In addition to the omn-present internal wave field, geostrophically balanced energetic mesoscale eddies are present in the ocean and they mix properties along rather than across isopycnals [see, e.g., McDougall et al. (2017) for a recent discussion of this issue]. The eddy-induced isopycnal mixing becomes most pro-nounced at high latitudes, western boundary currents, and in frontal regions such as the Antarctic Circumpolar Current (ACC). The mesoscale eddies, with typical horizontal scales of less than 100 km, cannot be resolved by climate models and coarse ocean models which have too coarse a resolution to explicitly represent any but the largest eddies. Conse-quently, many present-day climate models still rely on diffu-sive parameterizations of the effects of mesoscale eddies, usually involving a diffusion tensor consisting of an advective (antisymmetric) and a diffusive (symmetric) part. The advective

part seeks to parameterize the advection of properties by an additional eddy driven velocity (sometimes called skew diffusion) and is often parameterized with the Gent-McWilliams (GM) closure (Gent et al. 1995) involving the so-called GM coefficient $\kappa_{gm}$. The diffusive part of the mixing tensor can be interpreted as turbulent mixing of properties along isopycnal surfaces, and is often (as we do here) parameterized with the isotropic, so-called Redi scheme (Redi 1982) involving a single scalar isopycnal diffusivity $\kappa_{iso}$.

The effect of isopycnal diffusion on the oceanic tracer uptake is well studied (e.g., Gregory 2000; Kuhlbrodt et al. 2015; Abernathey and Ferreira 2015; Gnanadesikan et al. 2007, 2015a,b; Jones and Abernathey 2019), but its role for the ocean's dynamics, however, remains unclear. As the mesoscale eddies, predominantly generated by baroclinic instability, are of leading-order importance in setting up the ocean stratification, tracer transport, and uptake (e.g., Karsten et al. 2002; Booth and Kamenkovich 2008; McWilliams 2013; Kamenkovich et al. 2017), a better understanding of the influence of isopycnal diffusivity on ocean circulation becomes important. In this paper, we investigate the effect of changing isopycnal diffusivity in a model with global ocean configuration.

<footnote>ª Institut für Meereskunde, Universität Hamburg, Hamburg, Germany</footnote>

<footnote>ᵇ Univ. Brest, CNRS, Laboratoire d'Océanographie Physique et Spatiale, IUEM, Brest, France</footnote>

<footnote>ᶜ Institut für Umweltphysik and MARUM, Universität Bremen, Bremen, Germany</footnote>

<footnote>*Corresponding author:* Ashwita Chouksey, ashwita.chouksey@univ-brest.fr</footnote>
The sensitivity of the ocean circulation to the GM coefficient $\kappa_{gm}$ has extensively been studied (Danabasoglu and McWilliams 1995; Kamenkovich and Sarachik 2004; Griesel and Morales Maqueda 2006; Viebahn and Eden 2010; Farneti and Gent 2011; Kuhlbrodt et al. 2012; Jochum and Eden 2015). However, since along-isopycnal diffusion of temperature and salinity in the ocean interior has to first order no effect on density gradients, apart from supposedly minor nonlinear effects, its potential effect on ocean circulation has received less attention and studies have focused mainly on the sensitivity of passive tracer and heat uptake: Pradal and Gnanadesikan (2014) discuss the effect of an increase in isopycnal diffusivity in a global climate model and find that the poleward transport of salt increases and suggest that this would enhance convection, with further discussion on feedback loops involving sea ice. Abernathey and Ferreira (2015) find that stronger Southern Ocean wind stress increases eddy kinetic energy resulting in an increase in isopycnal diffusivity, which is a crucial factor controlling Southern Ocean ventilation. Using a coupled Earth system model with a constant isopycnal diffusivity, Gnanadesikan et al. (2015a) show that a higher isopycnal diffusivity $\kappa_{iso}$ leads to a significant increase in oceanic carbon uptake.

Studies have also looked at the connection between isopycnal diffusivity and circulation and water mass formation. Sijp et al. (2006) and Sijp and England (2009) discuss that isopycnal diffusion affects the stability of North Atlantic Deep Water (NADW) formation making it more stable to freshwater perturbations and attribute this to the increased diffusion of salinity, preventing a halocline, that otherwise suppresses sinking. Sévellec and Fedorov (2011) also study the stability of the Atlantic meridional overturning circulation (AMOC) and describe a weakening of the circulation with increasing isopycnal diffusivity in a zonally averaged model, but they did not further explain this decrease. A recent study by Jones and Abernathey (2019) shows a reduction in the fraction of NADW in the deep Pacific with an increase in isopycnal diffusivity in an idealized-geometry ocean model while keeping the circulation fixed. Ragen et al. (2020) explicitly address circulation changes in the Southern Ocean due to increased isopycnal diffusion and find a decrease in ACC transport in a coupled climate model for increased isopycnal diffusivity. They suggest a role for feedback between changed SST and increased wind stress curl contributing to the decrease. However, it is difficult to understand the primary driver of the coupled model results, whether it is the atmosphere-ocean feedbacks related to wind stress and sea ice [as suggested by Ragen et al. (2020)], or the transmission of the surface condition changes to the interior that changes the pressure gradients, and thus the circulation. Here, in an uncoupled model where the wind stress is fixed, we find similar changes in the ACC transport, thus narrowing the possible mechanisms, and show explicitly how the diffusive fluxes along isopycnals are affected, and also discuss a decrease in the meridional overturning circulation.

In this study, we employ a global ocean model with realistic geometry and allow the circulation to change with the changing isopycnal diffusivity $\kappa_{iso}$. We focus on the role of isopycnal

diffusivity on ocean circulation by changing isopycnal diffusivities $\kappa_{iso}$ in the range of observed values in the ocean and assess its effects on the two active tracers' temperature and salinity. Estimates of isopycnal diffusivities from observations vary widely and are still uncertain. Groeskamp et al. (2017) estimate an isotropic diffusivity $\kappa_{iso}$ with an average value of $350 \text{ m}^2 \text{ s}^{-1}$ from observation-based inverse methods in the ocean interior, while larger values up to $10,000 \text{ m}^2 \text{ s}^{-1}$ are found closer to the ocean surface from drifter observations (Zhurbas et al. 2014), consistent with a general decrease of eddy diffusivities with depth as inferred from Argo float data (Cole et al. 2015; Roach et al. 2016). Regional estimates of isotropic isopycnal eddy diffusivities around $1000 \text{ m}^2 \text{ s}^{-1}$ have been estimated from floats and tracer releases in the ACC and North Atlantic (Ledwell et al. 1998; Tulloch et al. 2014). However, not much is known in general about the anisotropy of isopycnal mixing so far, though there is evidence that in the presence of strong jets the mixing is reduced in the cross-stream direction (e.g., Ferrari and Nikurashin 2010; Griesel et al. 2010; Klocker and Abernathey 2014; Griesel et al. 2014). Therefore (and also in this study), ocean models use isotropic GM and Redi schemes, with values close to $1000 \text{ m}^2 \text{ s}^{-1}$ for both the Redi and GM coefficients, $\kappa_{iso}$ and $\kappa_{gm}$, but this equality has no further justification. In fact, linearized quasi-geostrophic theory predicts the PV diffusivity (which can be a surrogate of isopycnal diffusivity) to not be equal to the buoyancy (GM) diffusivity (Smith and Marshall 2009). We will vary only the Redi coefficient, i.e., $\kappa_{iso}$ and keep the GM coefficient $\kappa_{gm}$ constant in our model experiments. To compare the circulation changes attributed to changes in $\kappa_{iso}$, we also assess the impact of changing the formulation of the surface fluxes of heat and freshwater on the ocean circulation in a suite of model experiments. Our results suggest that, in addition to the tracer uptake, the along-isopycnal diffusion controlled by the isopycnal diffusivity plays a more dominant role than the surface flux restoring time scales in controlling the circulation, at least with a fivefold increase of isopycnal diffusivity.

The paper is arranged as follows: section 2 describes the numerical model used and different model experiments; section 3 compiles the results on changes in circulation, isopycnal tracer gradients, and isopycnal diffusive fluxes, including the effect on meridional heat and freshwater transports as well as the effect of fixed fluxes; and a summary and conclusions follow in section 4.

## 2. Model and experiments

We use the realistic global non-eddy-resolving ocean model by Eden (2016) with a horizontal resolution of ~2° × 2° and 45 vertical levels. The across isopycnal mixing by internal wave breaking is provided by the IDEMIX closure (Internal Gravity Wave Dissipation, Energy and Mixing) (Olbers and Eden 2013) and the prognostic closure by Gaspar et al. (1990) is used to parameterize small-scale turbulent mixing. Mesoscale eddy mixing is implemented with the mixing operator following Griffies (1998), where the additional eddy-driven advection velocity is implemented in form of an
<table><caption>TABLE 1. Description of the model setups along with the initial conditions, the prescribed along isopycnal diffusivity ($\kappa_{iso}$) values and the length of each run. The value of $\kappa_{gm}$, representing effect of eddy advection, is fixed to 1000 m$^2$ s$^{-1}$ for all the model setups.</caption><thead><tr><th>Acronym</th><th>Surface flux formulation</th><th>Initial conditions</th><th>Length of the run (yr)</th><th>κ<sub>iso</sub> (m<sup>2</sup> s<sup>-1</sup>)</th></tr></thead><tbody><tr><td>K05_REST</td><td>Restoring to SST, SSS</td><td>State of rest</td><td>8218</td><td>500</td></tr><tr><td>K07_REST</td><td>Restoring to SST, SSS</td><td>State of rest</td><td>8205</td><td>700</td></tr><tr><td>K1_REST</td><td>Restoring to SST, SSS</td><td>State of rest</td><td>6789</td><td>1000</td></tr><tr><td>K2_REST</td><td>Restoring to SST, SSS</td><td>State of rest</td><td>5726</td><td>2000</td></tr><tr><td>K5_REST</td><td>Restoring to SST, SSS</td><td>State of rest</td><td>5221</td><td>5000</td></tr><tr><td>K5_REST_K1TOK5</td><td>Restoring to SST, SSS</td><td>K1_REST after 6500 yr</td><td>7500</td><td>5000</td></tr><tr><td>K1_FORC1</td><td>FK1, fixed flux diag. from K1_REST</td><td>K1_REST after 5500 yr</td><td>6500</td><td>1000</td></tr><tr><td>K5_FORC1_RESTλ0</td><td>Eqs. (1), (2) with λ = 0</td><td>K1_REST after 5500 yr</td><td>6500</td><td>5000</td></tr><tr><td>K5_FORC1_RESTλ0.01</td><td>Eqs. (1), (2) with λ = 0.01</td><td>K1_REST after 5500 yr</td><td>3623</td><td>5000</td></tr><tr><td>K5_FORC1_RESTλ0.5</td><td>Eqs. (1), (2) with λ = 0.1</td><td>K1_REST after 5500 yr</td><td>3660</td><td>5000</td></tr><tr><td>K5_FORC1_RESTλ0.5</td><td>Eqs. (1), (2) with λ = 0.5</td><td>K1_REST after 5500 yr</td><td>3450</td><td>5000</td></tr><tr><td>K5_FORC1_RESTλ1</td><td>Eqs. (1), (2) with λ = 1</td><td>K1_REST after 5500 yr</td><td>3880</td><td>5000</td></tr></tbody></table>

antisymmetric (skew) component. In the case of isopycnal slopes being too steep, the mixing scheme by Griffies (1998) becomes unstable. Therefore, the isopycnal diffusivity is multiplied by the factor $d_{taper} = \{1 + \tanh[(s_c - |s|)/s_d]\}/2$, where $s_c = 0.001$ and $s_d = 0.001$ are parameters controlling the maximal allowed isopycnal slopes, and the isopycnal slope $s = -\nabla_h\rho/\partial_z\rho$ with $\rho$ being the locally referenced potential density. In regions where the isopycnal slopes are too steep, isopycnal mixing is replaced by lateral mixing with a prescribed constant diffusivity of $\kappa_{steep} = 500 \text{ m}^2 \text{ s}^{-1}$ multiplied with the factor $(1 - d_{taper})$. We set $\kappa_{gm} = 1000 \text{ m}^2 \text{ s}^{-1}$ for all the model experiments and isopycnal diffusivity $\kappa_{iso}$ is varied from 1000 to 5000 $\text{ m}^2 \text{ s}^{-1}$ (see Table 1). We choose $\kappa_{iso} = 1000 \text{ m}^2 \text{ s}^{-1}$ as an average value commonly used in ocean models and $\kappa_{iso} = 5000 \text{ m}^2 \text{ s}^{-1}$ which is within the range of isopycnal diffusivities as derived from observations or eddying ocean models (e.g., Vollmer and Eden 2013; Tulloch et al. 2014; Griesel et al. 2015; Cole et al. 2015).

We also found the surface boundary condition formulation to be important in the model experiments. The model by Eden (2016) uses a monthly climatology of realistic forcing datasets for momentum and heat fluxes and a restoring boundary condition for surface salinity. The surface heat flux boundary condition follows Barnier et al. (1995) and is also called a Haney-type surface boundary condition (Haney 1971). The forcing as in Barnier et al. (1995) is based on the linearized bulk formulas and assumes an infinite heat reservoir of the atmosphere, but allows for an evolving ocean circulation. The restoring time scale for surface salinity is 30 days for the 10-m-thick surface grid box. There is no explicit sea ice model. In the case of surface temperatures below the freezing point, surface heat fluxes out of the ocean and salinity restoring are set to zero (but not the momentum fluxes). In the remainder of this paper, we will refer to the combination of Haney-type boundary conditions for temperature and restoring boundary conditions for salinity as restoring boundary conditions.

Since the restoring formulation for the surface fluxes by Barnier et al. (1995) allows the heat and freshwater fluxes to change with increasing isopycnal diffusivity, we test whether similar changes in ocean circulation are obtained when the surface fluxes are not allowed to change as much by

increasing the restoring time scales. Therefore, as alternative surface boundary condition formulation, we diagnosed the resulting heat and freshwater fluxes in the simulations using the restoring conditions in form of monthly means over a 10-yr time period at the end of the model experiments. These diagnosed surface fluxes are then applied directly in a subsequent simulation as monthly means to the model, which we call fixed surface fluxes in the following. To allow for varying degrees of restoring time scales we use the following formulation:

$$F(T) = F_1(T) + \lambda \times \frac{Q}{\rho_0 c_p} (T_s - T^*), \quad (1)$$

$$F(S) = F_1(S) + \lambda \times \frac{\Delta z_s}{t_s} (S_s - S^*), \quad (2)$$

where $F_1(T)$ denotes the fixed monthly mean heat flux diagnosed from the simulation with $\kappa_{iso} = 1000 \text{ m}^2 \text{ s}^{-1}$ using restoring conditions, and $F_1(S)$ denotes the corresponding fixed monthly mean surface salt flux. The second terms on the RHS of Eqs. (1) and (2) are the restoring terms for temperature and salinity, where $T^*$ and $S^*$ are the surface climatologies from Barnier et al. (1995); $T_s$ and $S_s$ are the model's sea surface temperature and salinity; $t_s = 30$ days; and $\Delta z_s$ the uppermost level thickness. The term $Q/(\rho_0 c_p)$ can be interpreted as a piston velocity that restores the surface temperatures using the heat flux climatology $Q$ leading to restoring time scales between 40 and 70 days (Barnier et al. 1995). The factor $0 \le \lambda \le 1$ allows to vary the restoring time scales. For decreasing $\lambda$, the restoring time scale is increased allowing the surface flux to adjust less to a changing isopycnal diffusivity, $\lambda = 0$ means the fluxes are fixed to $F_{1,S}$, while $\lambda = 1$ allows the surface fluxes to adjust toward the climatologies.

All model setups are described below in Table 1. Naming conventions for the experiments (Table 1 and hereafter in the paper) are such that the first part refers to the chosen value of the isopycnal diffusivity and the second part refers to the surface flux formulation. The first set of experiments use the restoring formulation for the heat flux and freshwater flux. Initial conditions are climatologies for temperature and
salinity and a state of rest. These sets of experiments are called K1_REST and K5_REST with isopycnal diffusivities $\kappa_{iso} = 1000 \text{ m}^2 \text{ s}^{-1}$ and $\kappa_{iso} = 5000 \text{ m}^2 \text{ s}^{-1}$, respectively. We also diagnose the adjustment process to an ocean state with increased isopycnal diffusivity, by performing an additional simulation in which the K1_REST ocean state after 6500 years of integration is used as the initial condition and the diffusivity is abruptly increased to $\kappa_{iso} = 5000 \text{ m}^2 \text{ s}^{-1}$. This simulation is then further integrated for another 1000 years to reach a new equilibrium state with $\kappa_{iso} = 5000 \text{ m}^2 \text{ s}^{-1}$ using the restoring formulation for the surface heat and freshwater fluxes. This setup is called K5_REST_K1TOK5. Five model runs using fixed surface fluxes in addition to varying degrees of restoring time scales will be discussed: K1_FORC1 and K5_FORC1 both use a fixed forcing derived from K1_REST (corresponding to $\lambda = 0$), but in K1_FORC1, $\kappa_{iso}$ is left unchanged at $1000 \text{ m}^2 \text{ s}^{-1}$, while for K5_FORC1 $\kappa_{iso}$ is increased to $5000 \text{ m}^2 \text{ s}^{-1}$. The K5_FORC1_REST$\lambda$ experiments use the surface flux formulation given by Eqs. (1) and (2) with the restoring time scale from the REST experiments ($\lambda = 1$, K5_FORC1_REST$\lambda$1), the restoring time scales increased by a factor 2 ($\lambda = 0.5$, K5_FORC1_REST$\lambda$0.5), by a factor 10 ($\lambda = 0.1$, K5_FORC1_REST$\lambda$0.1) and by a factor 100 ($\lambda = 0.01$, K5_FORC1_REST$\lambda$0.01). All the analyses presented in this study are conducted after all the model experiments have reached a steady state. The next sections describe the results from the different model experiments.

## 3. Results

We first consider the changes in ocean circulation and density in experiments K1_REST and K5_REST in section 3a, discuss the changing isopycnal tracer gradients in section 3b, and the direction of the isopycnal diffusive fluxes in section 3c. We show the effect on heat and freshwater transports in section 3d, consider the adjustment process to the increased isopycnal diffusion in section 3e, and discuss the role of the surface forcing in section 3f explaining different sensitivity experiments enlisted in the Table 1. We focus on large-scale temperature and salinity changes, zonal transports, and the residual meridional overturning circulation. The latter is here taken as the sum of the Eulerian mean meridional overturning and the parameterized eddy-induced one (Olbers et al. 2012).

### a. The effect of increased isopycnal diffusion on circulation and density

We begin by comparing the changes in circulation in the experiments with restoring flux formulation, i.e., the experiments K1_REST and K5_REST. By increasing the isopycnal diffusivity $\kappa_{iso}$ from $1000 \text{ m}^2 \text{ s}^{-1}$, the residual meridional overturning circulation in the North Atlantic decreases from 20.79 to 14.33 Sv (1 Sv $\equiv 10^6 \text{ m}^3 \text{ s}^{-1}$; Fig. 1 and Table 2). In the Southern Ocean the volume transport of the residual circulation is not affected very much and decreases only slightly from 28.62 to 27.12 Sv which is within the limits of the standard deviation, but the circulation cell extends in depth. The positive cell in the Southern Ocean associated with

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

FIG. 1. Residual meridional overturning circulation for the model setups (a) K1_REST and (b) K5_REST. The color shading indicates the volume transport (Sv) with overlaid black contours indicating the same. The negative contours indicate the volume transport in the clockwise direction, and the positive contours indicate the volume transport in anticlockwise direction. The minimum in the Northern Hemisphere streamfunction between 40° and 80°N corresponds to the maximum AMOC, the minimum of the Southern Hemisphere clockwise circulation between 40° and 60°S corresponds to the maximum Southern Ocean residual circulation (cf. Table 2). The counterclockwise cell south of 60°S is related to AABW formation.

Antarctic Bottom Water (AABW) formation increases slightly. The most significant change for the Southern Ocean is that the ACC volume transport decreases from 149.5 to 107.60 Sv (Table 2). Understanding these surprisingly large changes is the task of this study. We also note that the sensitivity of the circulation to isopycnal diffusivity is reduced for values < $1000 \text{ m}^2 \text{ s}^{-1}$. We find that in those experiments there still is a slow drift of increasing ACC and AMOC and larger values will be reached after the circulation has reached a steady state. The ACC and AMOC in the experiment with $K_{iso} = 2000 \text{ m}^2 \text{ s}^{-1}$ are reduced by about 14 and 2 Sv, respectively.

Temperature and salinity for K1_REST and the difference between the experiments K5_REST and K1_REST are shown for the zonally averaged temperature and salinity in Fig. 2. The strongest signals are the cooling and freshening of the deep ocean below 500 m in the NADW and the cooling and freshening in the upper ocean between latitudes 20° and 40°. This means that around the 40° latitudes in the Northern and
<table><caption>TABLE 2. Volume transports (Sv; 10-yr average and associated standard deviation) for the residual Atlantic meridional overturning circulation (AMOC), the Southern Ocean meridional residual circulation, and the Antarctic Circumpolar Current (ACC) (through Drake Passage) for the different model setups (first column) described in Table 1. The ↑ indicates still increasing trend.</caption><thead><tr><th colspan="2">Model setups</th><th colspan="3">Volume transport (Sv)</th></tr><tr><th>Acronym</th><th>κ<sub>iso</sub> (m<sup>2</sup> s<sup>-1</sup>)</th><th>AMOC</th><th>Southern Ocean</th><th>ACC</th></tr></thead><tbody><tr><td>K05_REST</td><td>500</td><td>20.08 ± 1.91 ↑</td><td>28.14 ± 3.17</td><td>148.2 ± 3.45 ↑</td></tr><tr><td>K07_REST</td><td>700</td><td>20.83 ± 1.96 ↑</td><td>28.57 ± 3.13</td><td>149.7 ± 3.47 ↑</td></tr><tr><td>K1_REST</td><td>1000</td><td>20.79 ± 1.95</td><td>28.62 ± 3.17</td><td>149.5 ± 3.47</td></tr><tr><td>K2_REST</td><td>2000</td><td>18.64 ± 1.73</td><td>28.09 ± 3.29</td><td>135.7 ± 3.57</td></tr><tr><td>K5_REST</td><td>5000</td><td>14.33 ± 1.38</td><td>27.12 ± 3.22</td><td>107.6 ± 3.39</td></tr><tr><td>K5_REST_K1TOK5</td><td>5000</td><td>14.30 ± 0.49</td><td>27.61 ± 3.09</td><td>111.60 ± 3.56</td></tr><tr><td>K1_FORC1</td><td>1000</td><td>22.35 ± 0.52</td><td>28.23 ± 2.83</td><td>142.20 ± 3.09</td></tr><tr><td>K5_FORC1_RESTλ0</td><td>5000</td><td>17.89 ± 0.60</td><td>27.10 ± 2.16</td><td>86.54 ± 3.18</td></tr><tr><td>K5_FORC1_RESTλ0.01</td><td>5000</td><td>15.44 ± 1.36</td><td>27.59 ± 3.23</td><td>105.3 ± 3.02</td></tr><tr><td>K5_FORC1_RESTλ0.1</td><td>5000</td><td>13.92 ± 1.24</td><td>26.93 ± 3.32</td><td>107.7 ± 3.04</td></tr><tr><td>K5_FORC1_RESTλ0.5</td><td>5000</td><td>11.67 ± 1.01</td><td>27.67 ± 3.55</td><td>109.2 ± 2.88</td></tr><tr><td>K5_FORC1_RESTλ1</td><td>5000</td><td>11.22 ± 0.93</td><td>28.25 ± 3.77</td><td>111.5 ± 2.59</td></tr></tbody></table>

Southern hemispheres the meridional temperature gradients decrease, whose effect is a decrease in meridional density gradient, while the decrease in salinity would lead to an increase in meridional density gradients.

To analyze the degree of compensation between both factors, we show the effects of the temperature and salinity changes on the (dynamically relevant) in situ density in Figs. 3a and 3b. Figure 3a shows the density $ρ(T, S, z)$, using temperature $T$ from K5_REST but salinity $S$ from K1_REST, where $ρ(T, S, z)$ using $T$ and $S$ from K1_REST is removed, and then zonally averaged. That is, it shows the density change is only due to the temperature change.

Between 40°S and 40°N for depths below 200 m, water masses become all denser. At depths between 250 and 1000 m around latitudes of ±40°, the maximum density increases by up to ≈0.24 kg m⁻³ in the Southern Hemisphere and by more than ≈0.26 kg m⁻³ in the Northern Hemisphere. The

corresponding density changes due to salinity while keeping temperature at its value from K1_REST (Fig. 3b), show a general lightening of all water masses between 40°S and 40°N of similar structure, but note that the minima associated with the salinity decrease are not exactly collocated with the maxima due to the temperature decrease. The water becomes lighter in the upper ocean by only about 0.18 kg m⁻³. Hence the effect of the temperature change on density dominates, and the net effect on the in situ density—$ρ(T, S, z)$ from K5_REST minus $ρ(T, S, z)$ from K1_REST shown in Fig. 3c—is a decrease in the meridional density gradient around the ±40° in the upper 200–1500 m. In particular between ±20° and ±50°, density increases by about 0.16 kg m⁻³ in the Southern Ocean between depths of 500 and 1200 m and by about 0.18 kg m⁻³ in the Northern Hemisphere between depths 300 and 700 m. The density poleward of those latitudes decreases slightly. The effect of temperature and salinity

<figure><img src="image_3.png" /><figcaption>FIG. 2. Zonally averaged latitude–depth sections for (a),(b) temperature and (c),(d) salinity. (left) The results from the setup K1_REST and (right) the difference $K5_REST - K1_REST$. The color shading and overlaid black contours indicate temperature (°C) in (a) and (b) and salinity (g kg⁻¹) in (c) and (d).</figcaption></figure>
<figure><img src="image_4.png" /><figcaption>FIG. 3. Zonally and yearly averaged in situ density differences (kg m⁻³) for K5_REST – K1_REST as a function of latitude and depth. (a) In situ density difference due to temperature change with model nonlinear equation of state, (b) in situ density difference due to salinity change with model nonlinear equation of state, (c) actual in situ density change with model nonlinear equation of state, and (d) density change using a linear equation of state (with constant thermal expansion coefficient, α = 1.65 × 10⁻⁴ K⁻¹ and contraction coefficient, β = 0.78 × 10⁻³ kg g⁻¹). (e) Contribution of cabbeling to density difference and (f) contribution of thermobaricity to density difference. The overbar denotes zonal and time average.</figcaption></figure>

cancel each other to a larger part in the deep ocean below 1500-m depth, where the density decreases only by up to 0.1 kg m⁻³.

Other than the GM diffusivity, the isopycnal diffusivity does not directly change density in the ocean interior, unless effects of the nonlinear dependency of ρ on T and S play a role, or where isopycnals outcrop and there is an interaction with the surface fluxes. To explore the effect of the nonlinear equation of state, we show in Fig. 3d the corresponding density change due to temperature and salinity using a linear equation of state

$$ \Delta\rho_{\text{lin}} = \rho_0(\beta_0\Delta S - \alpha_0\Delta T), \quad (3) $$

with β₀ = 0.78 × 10⁻³ kg g⁻¹, α₀ = 1.65 × 10⁻⁴ K⁻¹, and ρ₀ = 1027 kg m⁻³ (Vallis 2006), instead of the nonlinear equation of state used in the model integration (IOC et al. 2010). The Δρₗₗₙ is similar to the density difference using the model's

nonlinear equation of state (cf. Figs. 3c,d), but there are also differences, with the largest deviation occurring in the deep ocean below 1000-m depth. These differences can be partly but not completely decreased by different choices of β₀ and α₀ (not shown). The effect of thermobaricity Δρₗₙ (increase of the thermal expansion coefficient α₀ with depth) and cabbeling Δρₗₙ (increase of α₀ with temperature) can be quantified using the formulations of Vallis (2006) and Roquet et al. (2015) as

$$ \begin{aligned} \Delta\rho_{\text{tb}} &= \rho_0^2 \alpha_0 \gamma_s g z \Delta T, \\ \Delta\rho_{\text{cb}} &= -\frac{1}{2} \rho_0 \beta_{\text{TS}} \left[ (T_5 - T_0)^2 - (T_1 - T_0)^2 \right] \end{aligned} \quad (4) $$

with γs = 1.11 × 10⁻⁸ s² m kg⁻¹, βₗₙ = 10⁻⁵ K⁻², T₀ = 9.85 K, S₀ = 35. The Δρₗₙ and Δρₗₙ are the most important contributions by nonlinearities in the equation of state, other contributions are of much smaller magnitude than Δρₗₙ and Δρₗₙ (not shown).
The negative effect of cabbeling on the density difference shown in Fig. 3e is compensated to a large extent by the positive effect of thermobaricity shown in Fig. 3f in the deep ocean below 1000 m. In the upper ocean, as expected, there is no large effect of the thermobaricity, but the cabbeling contributes with a positive density correction in the shallow subtropical thermocline, where the deviation from the basic state temperature used for the linearization is largest. We conclude that even though the exact contribution of the nonlinearities is difficult to quantify due to the dependence on the choice of coefficients and background state, the contribution from cabbeling to the density difference could be significant in the shallow subtropical thermocline. For the dominant density difference occurring between depths of 250 and 1000 m around the 40° latitudes (Fig. 3c) the contribution from nonlinearities is less significant.

### b. The effect of increased isopycnal diffusion on isopycnal tracer gradients and horizontal distributions of density changes

We now look at the contributions of the different ocean basins to the zonally averaged tracer gradient changes. Figures 4a–d show temperature and salinity averaged over 200–1000 m depth. It can be seen that the isopycnal gradients of temperature (Figs. 4a,b) and salinity (Figs. 4c,d) decrease from K1_REST (left panels) to K5_REST (right panels) as expected for stronger isopycnal diffusion. It can also be seen that the strongest isopycnal temperature and salinity gradients occur in the Southern Ocean around 40°S and in the North Atlantic between 40° and 60°N, most notably in the Western boundary currents. Associated with these changes in isopycnal tracer gradients are cooling and freshening around the 20°–40° latitudes in the upper 600 m (Figs. 2b,d) in the western boundary currents of all basins (Figs. 5a,b). The strongest cooling and freshening are seen in the Brazil Current, the Gulf Stream, and the Agulhas Current. The cooling is clearly associated with increased densities in the western boundary currents (Fig. 5d).

### c. Direction of isopycnal diffusive fluxes

The direction and magnitude of the isopycnal diffusive fluxes of temperature and salinity in the upper ocean can be seen from the positioning of isotherms and isohalines with respect to the isopycnals, shown for the setup K1_REST in Fig. 6. The slopes of the isotherms are more aligned with isopycnal slopes than for the isohalines, indicating that the temperature is to first order proportional to the density. The slopes of isohalines show much more misalignment with the isopycnal slopes, hence large isopycnal gradients of salinity. The diffusive flux is from higher concentration to lower concentration along the isopycnals and Fig. 6 indicates that for both temperature and salinity the diffusive flux is from the oceanic interior toward the surface for latitudes poleward of about 20°–30°. Equatorward of 30°, where the thermocline depth is deepest, the fluxes change direction and are equatorward and downward. The diffusive flux hence tends to transport heat and salt along the isopycnals toward the surface in

higher latitudes, consistent with the cooling/freshening and warming/salinification seen equatorward and poleward of 40°, respectively (Fig. 2).

Within the small-slope approximation ($|s| \ll 1$), the isothermal diffusive flux of temperature is given by

$$ 
\begin{aligned}
F_{\text{iso}}^{y} &= -\kappa_{\text{iso}}(\partial_x T + s_x \partial_z T), & F_{\text{iso}}^{y} &= -\kappa_{\text{iso}}(\partial_y T + s_y \partial_z T), \\
F_{\text{iso}}^{z} &= -\kappa_{\text{iso}}(\mathbf{s} \cdot \nabla_h T + |\mathbf{s}|^2 \partial_z T), & & (5)
\end{aligned}
$$

where $\mathbf{s} = (s_x, s_y) = -(\beta\nabla_h S - \alpha\nabla_h T)/(\beta\partial_z S - \alpha\partial_z T)$ is the isoneutral slope vector. The divergence of the flux $F_{\text{iso}} = (F_{\text{iso}}^x, F_{\text{iso}}^y, F_{\text{iso}}^z)$ changes temperature, and equivalently for salinity. The zonally averaged meridional isopycnal diffusive flux $F_{\text{iso}}^y$ diagnosed from K1_REST for the upper 1500 m is shown in Figs. 7a and 7b. We see a predominantly poleward flux poleward of 30°–40°, consistent with Fig. 6. Equatorward of about 40°, the flux increases with decreasing latitude and there is a sign change poleward of 30°, which is due to the thermocline “dip” in isolines of temperature, salinity, and density, changing the sign of $\partial_y T$. The northward increase in isopycnal diffusive fluxes leads to cooling ($-\partial_y F_{\text{iso}}^y$), which explains the cooling equatorward from 40° to about 20° in both hemispheres. The cooling is largest where the flux changes sign. Equivalently, the meridional gradients in the isopycnal salt fluxes (Fig. 7b) explain the freshening between 40° and 20°.

### d. Surface fluxes and meridional heat and freshwater transports

The along isopycnal fluxes of temperature and salinity are related to the isopycnal gradients and connect to the surface fluxes when the isopycnals outcrop. The zonally averaged surface heat and salt fluxes are shown in Figs. 7c and 7d, where negative salt flux indicates freshwater gain by the ocean surface and negative heat flux indicates heat loss of the ocean surface. Poleward of 40° (and both heat and salt fluxes predominantly decrease with increasing isopycnal diffusivity, whereas equatorward of about 40° heat and salt fluxes increase. This is consistent with Fig. 6, since the increased isopycnal diffusive fluxes of temperature and salinity poleward of 40° make the ocean surface warmer and saltier, hence leading to a decrease in the surface restoring fluxes. The opposite occurs equatorward of 40° where the increased downward diffusive fluxes tend to make the ocean surface colder and fresher leading to an increase in the predominantly positive surface restoring fluxes of temperature and salinity.

In equilibrium, the zonally averaged surface heat and freshwater fluxes are directly related to the zonal and vertical integral of the meridional heat and freshwater transports. Figures 7e and 7f show those meridional transports for K1_REST and K5_REST split up into the advective (including the GM eddy advection) and the isopycnal diffusive contributions. The increase of the predominantly poleward isopycnal diffusive transport with increasing isopycnal diffusivity is apparent for both hemispheres (purple lines). The northward advective heat and salt fluxes in the Northern Hemisphere decrease, related to the decrease in overturning circulation. In the Southern Hemisphere, the
<figure><img src="image_5.png" /><figcaption>FIG. 4. (a),(b) Temperature (°C) and (c),(d) salinity (g kg⁻¹) averaged over 200–1000 m depth for (left) K1_REST and (right) K5_REST.</figcaption></figure>

southward advective heat and salt transports increase (more precisely south of 40°S for the salt transport), consistent with the net surface heat and salt loss (Figs. 7c,d). As a result of the increased isopycnal diffusion, mixed layers become deeper in the Southern Ocean, consistent with enhanced heat loss south of 60°S. The advective transports are dominant, but the isopycnal transports contribute significantly south of 20°S.

### e. Diagnosing the adjustment process to the ocean state with increased isopycnal diffusion

The comparison of the two equilibrium states in the previous sections gave insight into the changes after the adjustment to the new equilibrium. In particular, the density changes shown in Figs. 2 and 5 are also the result of the adjusted circulation. Those changes are the consequence of an interplay of the advective, isopycnal and diapycnal diffusive fluxes, but they are initiated by the increase in isopycnal diffusion. To isolate the effect of those individual processes we therefore now discuss the adjustment to the equilibrium with increased isopycnal diffusion to $κ_{iso} = 5000 \text{ m}^2 \text{ s}^{-1}$ with the experiment K5_REST_K1TOK5, in which isopycnal diffusivity is abruptly increased starting from the equilibrium of K1_REST after 6500 years and then integrated for another 1000 years with $κ_{iso} = 5000 \text{ m}^2 \text{ s}^{-1}$. We note that this experiment has the same residual overturning circulation in the new equilibrium

for the North Atlantic and Southern Ocean as K5_REST but a slightly larger ACC of 111 Sv compared to 107 Sv for K5_REST. This difference occurs because high-latitude surface fluxes in the Southern Ocean adapt to this sudden increase in isopycnal diffusivity and less cold water is mixed downward leading to a less cold AABW on the poleward side of the ACC below 1500 m. This reduces the meridional density gradient in the deep ocean contributing to a smaller ACC decrease than in K5_REST. We note that the salinity difference distributions across the full depth and temperature difference distributions in the upper 1500 m are not affected and we focus in the following discussion on the upper 1500 m.

The local time rate of change of tracer concentration $T$ can be written as

$$ \partial_t T = -\nabla \cdot (\mathbf{u}T) - \nabla \cdot [(\kappa_{iso} + \kappa_{gm}) \cdot \nabla T] + \partial_z(\kappa \partial_z T) \quad (6) $$

$$ = -\nabla \cdot \mathbf{F}_{adv} - \nabla \cdot (\mathbf{F}_{iso} + \mathbf{F}_{gm}) + \partial_z F_{vdiff}, \quad (7) $$

where (the tensor) $κ_{iso} + κ_{gm}$ is the Redi-GM isoneutral mixing tensor, e.g., given by Olbers et al. (2012), that includes the contributions both from the diffusive fluxes $\mathbf{F}_{iso} = (F_{iso}^x, F_{iso}^y, F_{iso}^z)$ related to $κ_{iso}$, and the skew diffusive fluxes $\mathbf{F}_{gm}$ related to $κ_{gm}$. The last term in the two rows is the cross-isopycnal diffusive flux with vertical diffusivity $κ$. The change in temperature
<figure><img src="image_6.png" /><figcaption>FIG. 5. (a) Temperature difference K5_REST-K1_REST averaged over 200–1000 m depth, (b) salinity difference K5_REST-K1_REST averaged over 200–1000 m depth, (c) density of K1_REST averaged over 200–1000 m depth, (d) density difference K5_REST-K1_REST averaged over 200–1000 m depth.</figcaption></figure>

and salinity due to the individual flux terms **F** in Eq. (7) is given by the integral in time $T_F(t_0) - T_F(t_1) = -\int_{t_1}^{t_0} \nabla \cdot \mathbf{F} dt$, where the divergence of the flux **F** is integrated in time from the equilibrium with $\kappa_{iso} = 1000$ at time $t_0$ to the equilibrium with $\kappa_{iso} = 5000$ at time $t_1$.

Figures 8a–c and 8e show the integral of the temperature change due to advection [first term on the RHS of Eq. (7)], isopycnal advection and diffusion [second term on the RHS of Eq. (7)] and vertical diffusion [third term on the RHS of Eq. (7)]. We first note that there is a large cancellation between the advective and isopycnal flux contributions, and the individual flux terms integrated over time (Figs. 8a–c) are three orders of magnitude larger than the sum of the three terms (Fig. 8e), which gives the actual temperature change. The isopycnal term is the only one that leads to cooling between 20° and 40° south and north, consistent with the equatorward increase and sign change of the isopycnal-diffusive fluxes (Fig. 7a), and it leads to warming equatorward of 20°. This cooling/warming pattern is consistent with the meridional isopycnal diffusive flux changes shown in Fig. 7e and is balanced by the opposite warming/cooling pattern due to the advective term (Fig. 8a). Overall, the combined effect of the advective (adv) and isopycnal (iso) flux terms is a cooling equatorward of 40° and some warming poleward (Fig. 8c). Vertical mixing (vmix) contributes strongest in the upper 200 m and otherwise leads to

cooling poleward of 40° and minor warming equatorward (not shown).

We now focus on the prominent cooling between 20° and 40°S and between depths of 400–1000 m (also evident in Fig. 2b) and show the evolution with time of the three terms (the divergence of the advective (adv), isopycnal (iso), and vertical diffusive (vmix) fluxes) averaged zonally and over these latitude and depth ranges for both temperature and salinity (Figs. 8e,f). At time $t = t_0$, the local time rate of change of temperature and salinity is zero in the annual mean and so the sum of the three terms adds up to zero (blue line). Then the change to increased isopycnal diffusivity occurs leading to an abrupt increase in the cooling rate due to the isoneutral mixing term (black line). The advective term (red line) reacts with an increase in the rate of warming, while the time rate of change due to vertical mixing (green line) is unaffected by the change in isopycnal diffusivity. The total time rate of change (i.e., the sum of all three terms, blue line) means cooling due to the excess cooling caused by the increased isopycnal diffusivity. After about 400 years, the sum of the three terms is zero and the new equilibrium with higher isopycnal diffusion is reached (not shown). It is one in which the sum of the advective and isopycnal terms lead to continuing (small) cooling, balanced by (small) warming from vertical mixing. We see a similar scenario for the salinity time tendencies with freshening by the isopycnal fluxes,
<figure><img src="image_7.png" /><figcaption>FIG. 6. Zonally averaged isosurfaces of (a) temperature (isotherms) and (b) salinity (isohalines) shown along with the isosurfaces of density (isopycnals). Isopycnals are shown in black contours, and isotherms and isohalines in red contours, shown here for the setup K1_REST for the upper 280 m.</figcaption></figure>

balanced by salinification from the advective tendency term (Fig. 2f).

To summarize: increasing isopycnal diffusivity increases the meridional isopycnal fluxes whose gradient is equatorward, leading to predominantly cooling and freshening equatorward of 40°, that is only partly compensated by advection. Since the effect of temperature on density dominates over the effect of salinity on density, meridional density gradients are reduced. Along with the density changes, both ACC and North Atlantic residual overturning transports decrease. Nonlinear effects of the equation of state may play a significant role in explaining those density changes particularly in the shallow subtropical thermocline. To shed light on the influence of using restoring surface fluxes, we next consider experiments where the fluxes have been fixed.

### f. Toward fixed flux experiments

In this section, we discuss the model experiments where the surface fluxes of freshwater and heat are fixed to varying degrees as stated in section 2. While a restoring to sea surface temperature may be justified (it is realistic that the atmosphere-ocean heat flux responds to a warming ocean surface), a restoring surface flux for salinity is not realistic, since the

freshwater flux should not depend on the sea surface salinity. Hence the restoring surface fluxes respond to the changing isopycnal diffusivity and may introduce feedbacks that over or underestimate the sensitivity to isopycnal diffusivity. Before we compare the sensitivity to the isopycnal diffusivity for varying restoring time scales, we test whether the circulation changes when we change from restoring to fixed surface fluxes. For example, the sea ice edge position in the prescribed surface fluxes could be inconsistent with the changed model dynamics and lead to model drift. We see that the AMOC and ACC volume transports change by about 2 and 7 Sv, respectively (K1_REST compared to K1_FORC1). This ACC change is still not as large as the changes we see when the isopycnal diffusivity is increased under varying surface forcing (K5_FORC1_RESTλ experiments compared to K1_FORC1). We first discuss the simulation where a fixed surface flux diagnosed from K1_REST ($λ = 0$) is used with $κ_{iso} = 5000 \text{ m}^2 \text{ s}^{-1}$. The ACC in this simulation is drastically decreased to about 86 Sv, while the sensitivity of the AMOC is reduced and changes only by about 4 Sv (Table 2). As in the experiments with restoring surface boundary condition, we see an increase in the density for latitudes equatorward of 40°, but the
<figure><img src="image_8.png" /><figcaption>FIG. 7. (left) Zonally averaged latitude–depth section of meridional isopycnal diffusive flux of (a) temperature (m s⁻¹ K) and (b) salinity (m s⁻¹) for the upper 1500 m in K1_REST. The overlaid black contours indicate potential density referenced to the surface (kg m⁻³) (center) Zonally averaged surface (c) heat flux (°C m s⁻¹) and (d) salt flux (g kg⁻¹ m s⁻¹) for the setups K1_REST and K5_REST shown with solid and dashed lines, respectively. Positive freshwater flux is equivalent to negative salt flux. (right) Zonally and vertically integrated meridional (e) heat transport (PW) and (f) salt transport (Sv). Solid and dashed lines are for K1_REST and K5_REST, respectively. The purple curve is the contribution from isopycnal diffusion and the green curve from advection.</figcaption></figure>

density difference is now much larger, exceeding 0.24 kg m⁻³ in the upper ocean (Fig. 10a, note the different color scale as compared to Fig. 10b). This experiment is allowed to drift far away from realistic surface temperature and salinity distributions and hence unrealistic upper-ocean density.

With a short restoring time scale ($\lambda = 1$), the surface flux adjusts to the change in isopycnal diffusivity and is closer to the restoring flux from K5_REST (the difference to the K5_REST runs is 3 Sv for AMOC and 4 Sv for ACC). Around $\pm 40°$ latitude, the increased diffusion transports warm waters toward the surface and the restoring surface heat flux reacts by pumping cold water into the surface (cf. section 3d) that subsequently mixes thereby increasing the density of the isopycnal layer that outcropped there. For $\lambda = 0.01$, the restoring time scale is 100 times larger and the surface fluxes are closer to the fixed fluxes diagnosed from K1_REST so they are allowed to vary less in response to the increased isopycnal diffusion. However, even with this very long restoring time scale, the surface fluxes are still significantly different from K1_FORC1 (not shown). The upper-ocean density differences for $\lambda = 0.01$ (Fig. 9b) are much smaller than the ones with the fixed flux and are more similar to K5_REST-K1_REST (cf. Fig. 3c to Fig. 9b). We find that the AMOC is reduced by about 7 Sv for $\lambda = 0.01$ and it is reduced by about 11 Sv when the flux is allowed to adjust ($\lambda = 1$). The ACC is reduced by 30 Sv for $\lambda = 1$ and by 34 Sv for $\lambda = 0.01$. We see a consistent decrease in the sensitivity to isopycnal diffusivity for the AMOC when going toward a fixed flux formulation, while the sensitivity of the ACC increases with decreasing $\lambda$ with a sudden decrease from $\lambda = 0.01$ to $\lambda = 0$.

Overall, the impact of the changed surface flux formulation on the ACC is small as compared to the impact of increasing the isopycnal diffusivity, and the surface flux formulation does not seem to play a significant role for the sensitivity to the increase in isopycnal diffusivity. In particular, the sensitivity of the ACC to increased isopycnal diffusion is not reduced when the restoring time scale is increased as one might expect when surface fluxes are held fixed. This leaves nonlinear effects in the equation of state to explain the density changes. Figure 9 shows that the contribution from nonlinear terms is very similar for the varying length of restoring time scales, except for an effect at high latitudes. The largest density difference occurs again close to the surface where the effect of cabbeling is strong, but between 250- and 1000-m depth, and equatorward of $40°$, the contribution from the nonlinear terms is again less significant. There is no dependence of nonlinear effects on the surface flux formulation equatorward of $60°$, and in particular we do not see nonlinear terms to become more significant when $\lambda$ is decreased.

Overall, the analyses with the different restoring time scales show that the largest changes occur in the ACC, whose response is highly sensitive to the increase in isopycnal diffusivity. The AMOC on the other hand is less sensitive. Given that it is the isopycnal diffusion of temperature that is crucial in determining the density gradient changes, allowing the heat fluxes to change as in the restoring surface heat flux might actually more realistically represent the sensitivity to isopycnal diffusion.

## 4. Discussion and conclusions

The analysis presented in this paper shows a previously not much discussed sensitivity of ocean circulation to isopycnal
<figure><img src="image_9.png" /><figcaption>FIG. 8. Zonally averaged temperature change $\int_{j_0}^{f_1} \nabla \cdot Fdt$, where $t_1 - T_0 = 1000$ years, due to (a) the divergence of the advective flux (10³ K; adv), (b) the divergence of the isopycnal fluxes $F_{iso} + F_{gm}$ (10³ K; iso), and (c) the sum of the divergence of the advective, isopycnal, and vertical diffusive fluxes (K; adv + iso + vmix), as a function of depth and latitude. (d) As in (c), but for salinity (g kg⁻¹). (e) Rate of change of temperature (K yr⁻¹) zonally averaged and averaged over the black box in (c) due to advection (red), isopycnal fluxes (black), vertical diffusion (pink), the sum of the advective and isopycnal fluxes (green) and total time rate of change (blue). The switch to $\kappa_{iso} = 5000 \ m^2 \ s^{-1}$ occurs at year 10. The time rate of change for only the first 100 years are shown. (f) As in (e), but for salinity (g kg⁻¹ yr⁻¹).</figcaption></figure>

diffusivity. Using a ~2° × 2° non-eddy-resolving global ocean configuration and changing isopycnal diffusivity $\kappa_{iso}$ along with different formulations of surface forcing, we assess the changes in ocean circulation resulting from the changes in horizontal tracer gradients. A suite of model experiments is set up with varying initial conditions, restoring surface fluxes and fixed surface fluxes, and with different values of isopycnal diffusivity; most prominent choices are $\kappa_{iso} = 1000 \ m^2 \ s^{-1}$ and $\kappa_{iso} = 5000 \ m^2 \ s^{-1}$. Changes in isopycnal diffusivity affect the tracer gradients, density, and the resulting circulation. The ACC transport consistently decreases with increasing isopycnal diffusivity using surface boundary fluxes with different restoring

time scales. The AMOC also decreases with increasing isopycnal diffusivity but becomes less sensitive once the surface fluxes become more fixed. The effect of the isopycnal diffusivity increase is to decrease the isopycnal gradients of temperature and salinity, leading to a shallowing of isopycnals and a decrease in the in situ density gradients around 40° latitude in the Northern and Southern Hemisphere. The density in the experiments with restoring heat and salt fluxes as well as in the experiments where the surface fluxes are varied using different restoring time scales has decreased poleward of 40° latitude (upward/poleward salinity and temperature diffusion) and increased equatorward (downward/equatorward salinity and temperature
<figure><img src="image_10.png" /><figcaption>FIG. 9. In situ density difference $\Delta\rho$ (kg m$^{-3}$) for experiments (a) K5_FORC1$\lambda$0 – K1_FORC1 and (b) K5_FORC1$\lambda$0.01 – K1_FORC1.</figcaption></figure>

diffusion). Even though there is a compensating effect on density gradients between the diffusive fluxes of temperature and salinity inducing cooling and freshening or warming and salinification, the net effect is dominated by the temperature: the along-isopycnal poleward of around 40° latitude and upward diffusion of heat leads to cooling equatorward of 40° and warming poleward of 40°. This leads to a significant decrease in meridional density gradients around 40° latitude in both hemispheres. Note that the way how such gradients effect the circulation (AMOC) is beyond the scope of the current discussion, but it is clear that density changes will affect it. Since the ACC is in geostrophic balance with the meridional density gradients, the related decrease of the ACC is not surprising. The relation of meridional density gradients and the overturning circulation is less clear as discussed by, e.g., Straub (1996), Greatbatch and Lu (2003), and Brüggemann et al. (2011). References to a parallel

recent discussion can be found in Johnson et al. (2019), but the issue is not the scope of this study.

Figure 11 summarizes the processes that influence density in an isopycnal surface around ±40° latitude. The isoneutral diffusion scheme does not introduce any spurious diapycnal mixing, except in regions with steep slopes—essentially the surface mixed layer—where weak horizontal mixing is introduced for stability with a horizontal diffusivity that is the same in all experiments. We find that the difference in the number of grid points along any latitude between K1_REST and K5_REST where this slope threshold is exceeded, is smaller than 2% (not shown), and occurs mainly in the Southern Ocean at high latitudes. We thus conclude that the tapering does not contribute to a significant change in diapycnal flux. This leaves the surface fluxes (and their formulation) and the nonlinear equation of state as the two candidates responsible for the density changes described above and illustrated in

<figure><img src="image_11.png" /><figcaption>FIG. 10. In situ density difference $\Delta\rho$ – $\Delta\rho_{\text{lin}}$ (kg m$^{-3}$) for experiments (a) K5_FORC1$\lambda$1 – K1_FORC1, (b) K5_FORC1$\lambda$0.5 – K1_FORC1, (c) K5_FORC1$\lambda$0.1 – K1_FORC1, and (d) K5_FORC1$\lambda$0.01 – K1_FORC1.</figcaption></figure>
<figure><img src="image_12.png" /><figcaption>FIG. 11. Summary of the processes that influence density in the isopycnal layers around ±40° latitude under increasing isopycnal diffusivity.</figcaption></figure>

Fig. 11. We find that nonlinear effects in the equation of state likely play a role for explaining the changes in density. Studies comparing the ocean circulation when running ocean models using linear and nonlinear equations of state find significant differences (Nycander et al. 2015; Roquet et al. 2015). In these studies, nonlinear effects are of importance owing to their influence on the static stability of the water column thereby influencing the atmosphere-interior ocean exchanges. In fact, stratification, water mass distributions, and circulation are significantly different from realistic conditions using a linear equation of state in those simulations. Nycander et al. (2015), for example, show that cabbeling is crucial for the formation of AAIW. Groeskamp et al. (2016) show that enhanced isopycnal diffusivity on the flanks of the ACC contributes to enhanced AAIW transformation by cabbeling. Likely, the sensitivity to isopycnal diffusivity will be different in an ocean simulation run with a linear equation of state, but the scope of this study was to investigate the sensitivity in a realistic model simulation with nonlinear equation of state. We quantified the contribution from the most significant nonlinear terms with a formulation as in Roquet et al. (2015) for given temperature and salinity distributions. With this diagnostic we find that the contribution from cabbeling to the density difference could be significant in the shallow subtropical thermocline, but that for the dominant density difference occurring between depths of 250–1000 m around ±40° latitude the contribution from nonlinearities is less significant. Below 1000 m, the effect of thermobaricity and cabbeling on the density change tend to cancel each other to some extent.

The restoring surface flux is the other factor that influences density; the increased isopycnal diffusivity diffuses heat and salt poleward toward the surface and the restoring heat and freshwater fluxes react by pumping cold and freshwater into the surface waters that subsequently mix downward and

increase the density of the outcropping isopycnal layer around the ±40° latitude since temperature dominates the effect on density. We found that different restoring time scales for the surface flux formulation do not significantly affect the sensitivity of the circulation to increased isopycnal diffusivity. The experiment with fixed surface flux exhibits the largest sensitivity of the ACC, pointing at the role of nonlinearities for the density difference, but it also produces unrealistic upper-ocean densities.

While we find the dominance of temperature diffusion for the effect on circulation, Pradal and Gnanadesikan (2014) also discuss an increased poleward diffusion of heat and salt with a sixfold increase in isopycnal diffusivity employing a global climate model. They stress the importance of the salt diffusion toward the surface, suggesting that this then exerts a strong control over Antarctic Sea ice. Without an explicit sea ice component, we can exclude in our experiments such feedbacks on the circulation response to increased $\kappa_{iso}$. In the same coupled model Ragen et al. (2020) report a 20% decrease in Drake Passage transport with a sixfold increase in isopycnal diffusivity, which is comparable to the decrease seen here of 25%–30% with a fivefold increase in isopycnal diffusivity. While Ragen et al. (2020) attribute this decrease to both an increase in the density of the intermediate and deep waters on the northern edge of the ACC as well as an increase in the wind stress curl, we find similar reduction in ACC transport in our uncoupled model, suggesting that the feedback on the wind stress must be of minor importance for the ACC transport changes.

We note that once isopycnal diffusivity is increased, there is a possible feedback of the eddy-induced overturning where even though $\kappa_{gm}$ is kept constant, skew fluxes decrease when density gradients decrease. Using a smaller GM coefficient $\kappa_{gm}$ could potentially increase the sensitivity to the isopycnal diffusivity. Ocean models sometimes set the GM coefficient equal to the isopycnal diffusivity primarily for computational convenience. Estimates of both coefficients based on quasi-geostrophic theory, however, suggest that this may be unrealistic (Smith and Marshall 2009; Vollmer and Eden 2013). Kamenkovich and Sarachik (2004) discuss density and circulation changes similar to what we find here by varying the GM and isopycnal diffusivity, but setting them both equal. Our study suggests that part of the circulation sensitivity can be attributed to the increase in isopycnal diffusivity.

Our study shows that the strength of isopycnal diffusion not only controls the uptake of carbon and other passive tracers as shown in previous studies, but is also as important as the surface fluxes in shaping the circulation in the ocean. The results emphasize that constraining both $\kappa_{gm}$ and $\kappa_{iso}$ is crucial to improve climate models. Estimates from observations and eddying ocean models suggest that the eddy diffusivities are large where eddy kinetic energy is large and that mixing is suppressed across strong currents (Eden and Greatbatch 2008; Ferrari and Nikurashin 2010; Bates et al. 2014; Griesel et al. 2015). Flow-dependent parameterizations can induce additional feedbacks not only on the passive tracer concentrations but also on the circulation.
Groeskamp, S., R. P. Abernathey, and A. Kletterer, 2016: Water mass transformation by cabbeling and thermobaricity. *Geo-phys. Res. Lett.*, **43**, 10 835–10 845, https://doi.org/10.1002/2016GL070860.

—, B. M. Sloyan, J. D. Zika, and T. J. McDougall, 2017: Mixing inferred from an ocean climatology and surface fluxes. *J. Phys. Oceanogr.*, **47**, 667–687, https://doi.org/10.1175/JPO-D-16-0125.1.

Hallberg, R. W., and A. Gnanadesikan, 2006: The role of eddies in determining the structure and response of the wind-driven Southern hemisphere overturning: Results from the Modeling Eddies in the Southern Ocean (MESO) project. *J. Phys. Oceanogr.*, **36**, 2232–2252, https://doi.org/10.1175/JPO2980.1.

Haney, R. L., 1971: Surface thermal boundary condition for ocean circulation models. *J. Phys. Oceanogr.*, **1**, 241–248, https://doi.org/10.1175/1520-0485(1971)001<0241:STBCFO>2.0.CO;2.

Hofmann, M., and M. A. Morales Maqueda, 2011: The response of Southern Ocean eddies to increased midlatitude westerlies: A non-eddy resolving model study. *Geophys. Res. Lett.*, **38**, L03605, https://doi.org/10.1029/2010GL045972.

IOC, SCOR, and IAPSO, 2010: The international thermodynamic equation of seawater – 2010: Calculation and use of thermodynamic properties. Intergovernmental Oceanographic Commission, Manuals and Guides 56, UNESCO, 196 pp., http://www.teos-10.org/pubs/TEOS-10_Manual.pdf.

Jochum, M., and C. Eden, 2015: The connection between Southern Ocean winds, the Atlantic Meridional Overturning Circulation, and Indo-Pacific upwelling. *J. Climate*, **28**, 9250–9257, https://doi.org/10.1175/JCLI-D-15-0263.1.

Johnson, H. L., P. Cessi, D. P. Marshall, F. Schloesser, and M. A. Spall, 2019: Recent contributions of theory to our understanding of the Atlantic meridional overturning circulation. *J. Geophys. Res. Oceans*, **124**, 5376–5399, https://doi.org/10.1029/2019JC015330.

Jones, C. S., and R. P. Abernathey, 2019: Isopycnal mixing controls deep ocean ventilation. *Geophys. Res. Lett.*, **46**, 13 144–13 151, https://doi.org/10.1029/2019GL085208.

Kamenkovich, I., and E. Sarachik, 2004: Mechanisms controlling the sensitivity of the Atlantic thermohaline circulation to the parameterization of eddy transports in ocean GCMs. *J. Phys. Oceanogr.*, **34**, 1628–1647, https://doi.org/10.1175/1520-0485(2004)034<1628:MCTSOT>2.0.CO;2.

—, Z. Garraffo, R. Pennel, and R. A. Fine, 2017: Importance of mesoscale eddies and mean circulation in ventilation of the Southern Ocean. *J. Geophys. Res. Oceans*, **122**, 2724–2741, https://doi.org/10.1002/2016JC012292.

Karsten, R., H. Jones, and J. Marshall, 2002: The role of eddy transfer in setting the stratification and transport of a circum-polar current. *J. Phys. Oceanogr.*, **32**, 39–54, https://doi.org/10.1175/1520-0485(2002)032<0039:TROETI>2.0.CO;2.

Kletterer, A., and R. Abernathey, 2014: Global patterns of mesoscale eddy properties and diffusivities. *J. Phys. Oceanogr.*, **44**, 1030–1046, https://doi.org/10.1175/JPO-D-13-0159.1.

Kuhlbrodt, T., R. S. Smith, Z. Wang, and J. M. Gregory, 2012: The influence of eddy parameterizations on the transport of the Antarctic Circumpolar Current in coupled climate models. *Ocean Modell.*, **52–53**, 1–8, https://doi.org/10.1016/j.ocemod.2012.04.006.

—, J. M. Gregory, and L. C. Shaffrey, 2015: A process-based analysis of ocean heat uptake in an AOGCM with an eddy-permitting ocean component. *Climate Dyn.*, **45**, 3205–3226, https://doi.org/10.1007/s00382-015-2534-0.

Ledwell, J. R., A. J. Watson, and C. S. Law, 1998: Mixing of a tracer in the pycnocline. *J. Geophys. Res.*, **103**, 21 499–21 529, https://doi.org/10.1029/98JC01738.

McDougall, T. J., S. Groeskamp, and S. M. Griffies, 2017: Comment on Tailleux, R. Neutrality versus Materiality: A thermodynamic theory of neutral surfaces. *Fluids* 2016, **1**, 32. *Fluids*, **2**, 19, https://doi.org/10.3390/fluids2020019.

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

Nycander, J., M. Hieronymus, and F. Roquet, 2015: The nonlinear equation of state of sea water and the global water mass distribution. *Geophys. Res. Lett.*, **42**, 7714–7721, https://doi.org/10.1002/2015GL065525.

Olbers, D., and C. Eden, 2013: A global model for the diapycnal diffusivity induced by internal gravity waves. *J. Phys. Oceanogr.*, **43**, 1759–1779, https://doi.org/10.1175/JPO-D-12-0207.1.

—, J. Willebrand, and C. Eden, 2012: *Ocean Dynamics*. Springer, 703 pp.

Pradal, M.-A., and A. Gnanadesikan, 2014: How does the Redi parameter for mesoscale mixing impact global climate in an Earth system model? *J. Adv. Model. Earth Syst.*, **6**, 586–601, https://doi.org/10.1002/2013MS000273.

Ragen, S., M.-A. Pradal, and A. Gnanadesikan, 2020: The impact of parameterized lateral mixing on the Antarctic Circumpolar Current in a coupled climate model. *J. Phys. Oceanogr.*, **50**, 965–982, https://doi.org/10.1175/JPO-D-19-0249.1.

Redi, M. H., 1982: Oceanic isopycnal mixing by coordinate rotation. *J. Phys. Oceanogr.*, **12**, 1154–1158, https://doi.org/10.1175/1520-0485(1982)012<1154:OIMBCR>2.0.CO;2.

Roach, C. J., D. Balwada, and K. Speer, 2016: Horizontal mixing in the Southern Ocean from Argo float trajectories. *J. Geophys. Res. Oceans*, **121**, 5570–5586, https://doi.org/10.1002/2015JC011440.

Roquet, F., G. Madec, L. Brodeau, and J. Nycander, 2015: Defining a simplified yet “realistic” equation of state for seawater. *J. Phys. Oceanogr.*, **45**, 2564–2579, https://doi.org/10.1175/JPO-D-15-0080.1.

Sévellec, F., and A. V. Fedorov, 2011: Stability of the Atlantic meridional overturning circulation and stratification in a zonally averaged ocean model: Effects of freshwater flux, Southern Ocean winds, and diapycnal diffusion. *Deep-Sea Res. II*, **58**, 1927–1943, https://doi.org/10.1016/j.dsr2.2010.10.070.

Sijp, W. P., and M. H. England, 2009: The control of polar haloclines by along-isopycnal diffusion in climate models. *J. Climate*, **22**, 486–498, https://doi.org/10.1175/2008JCLI2513.1.

—, M. Bates, and M. H. England, 2006: Can isopycnal mixing control the stability of the thermohaline circulation in ocean climate models? *J. Climate*, **19**, 5637–5651, https://doi.org/10.1175/JCLI3890.1.

Smith, K., and J. Marshall, 2009: Evidence for deep eddy mixing in the Southern Ocean. *J. Phys. Oceanogr.*, **39**, 50–69, https://doi.org/10.1175/2008JPO3880.1.

Spence, P., O. A. Saenko, M. Eby, and A. J. Weaver, 2009: The Southern Ocean Overturning: Parameterized versus permitted eddies. *J. Phys. Oceanogr.*, **39**, 1634–1651, https://doi.org/10.1175/2009JPO4120.1.

Straub, D. N., 1996: An inconsistency between two classical models of the ocean buoyancy driven circulation. *Tellus*, **48A**, 477–481, https://doi.org/10.3402/tellusa.v48i3.12073.
Tulloch, R., and Coauthors, 2014: Direct estimate of lateral eddy
diffusivity upstream of Drake Passage. J. Phys. Oceanogr.,
**44**, 2593–2616, https://doi.org/10.1175/JPO-D-13-0120.1.

Vallis, G. K., 2006: *Atmospheric and Oceanic Fluid Dynamics*. Cambridge University Press, 964 pp.

Viebahn, J., and C. Eden, 2010: Towards the impact of eddies on the response of the Southern Ocean to climate change. *Ocean Modell.*, **34**, 150–165, https://doi.org/10.1016/j.ocemod.2010.05.005.

Vollmer, L., and C. Eden, 2013: A global map of meso-scale eddy diffusivities based on linear stability analysis. *Ocean Modell.*, **72**, 198–209, https://doi.org/10.1016/j.ocemod.2013.09.006.

Zhurbas, V., D. Lyzhkov, and N. Kuzmina, 2014: Drifter-derived estimates of lateral eddy diffusivity in the world ocean with emphasis on the Indian Ocean and problems of parameterisation. *Deep-Sea Res. I*, **83**, 1–11, https://doi.org/10.1016/j.dsr.2013.09.001.