Centennial Changes in the Indonesian Throughflow Connected to the Atlantic Meridional Overturning Circulation: The Ocean's Transient Conveyor Belt

Climate models consistently project a robust weakening of the Indonesian Throughflow (ITF) and the Atlantic meridional overturning circulation (AMOC) in response to greenhouse gas forcing. Previous studies of ITF variability have largely focused on local processes in the Indo‐Pacific Basin. Here, we propose that much of the centennial‐scale ITF weakening is dynamically linked to changes in the Atlantic Basin and communicated between basins via wave processes. In response to an AMOC slowdown, the Indian Ocean develops a northward surface transport anomaly that converges mass and modifies sea surface height in the Indian Ocean, which weakens the ITF. We illustrate these dynamic interbasin connections using a 1.5‐layer reduced gravity model and then validate the responses in a comprehensive general circulation model. Our results highlight the importance of transient volume exchanges between the Atlantic and Indo‐Pacific basins in regulating the global ocean circulation in a changing climate.


Introduction
As the only low-latitude oceanic pathway for freshwater and heat exchange between major ocean basins today, the Indonesian Throughflow (ITF) is an important component of the climate system (e.g., Godfrey, 1996;Gordon, 2005;Lee et al., 2015;Schneider, 1998). Climate models consistently project a weakening of the ITF in response to enhanced greenhouse gas forcing (e.g., Feng et al., 2017;Hu et al., 2015;Sen Gupta et al., 2016). In this study, we show that the Atlantic meridional overturning circulation (AMOC), which climate models also consistently project to decline in a warming climate (Cheng et al., 2013;Weijer et al., 2020), may be a primary cause of the ITF weakening over centennial timescales. This link between the AMOC and the ITF highlights that dynamical processes governing interbasin transport and exchange are critical for representing the transient behavior of the global ocean overturning circulation (Figure 1a).
The ITF transport, a westward flow sustained by a lateral pressure gradient between the western Pacific (high sea surface hight) and the eastern Indian Ocean (low sea surface height) (Wyrtki, 1987), varies over a range of timescales (e.g., Feng et al., 2018;Gordon, 2005;Sprintall et al., 2019). On short timescales, from subseasonal to decadal, low-latitude surface forcing, including surface wind stress forcing and precipitation, dominates ITF transport variability (e.g., Hu & Sprintall, 2017;Lee et al., 2019;Meyers, 1996;Sprintall et al., 2009). On longer timescales, decadal to centennial, the basin-scale wind stress curl determines the sea surface height distribution, which provides a strong constraint on the ITF transport through the Island Rule (Godfrey, 1989), described in section 2. This decadal timescale arises from the transit time for first-mode baroclinic Rossby Waves to cross the Pacific Ocean (Godfrey, 1996). Critically, changes to the surface wind stress in response to greenhouse gas forcing are too small to account for the projected centennial changes in ITF transport in climate models (Feng et al., 2017;Hu et al., 2015;Sen Gupta et al., 2016). Instead, the centennial ITF weakening has been attributed to a reduction in diapycnal upwelling below the thermocline in the Pacific Ocean (Feng et al., 2017). Here, we argue that this interpretation is inconsistent with the processes modifying the Pacific stratification, which are better described by an adiabatic downward displacement of isopycnals. The deepening of Pacific isopycnals is a result of a weakened ITF responding to a variable AMOC. Gordon (1986) first highlighted the ITF as a critical pathway for upwelled Pacific Deep Water (PDW) to return to the Atlantic Ocean, closing the global ocean overturning circulation. The classical "conveyor belt" analogy of the global ocean overturning circulation (Broecker et al., 1991) highlights connections between the ITF and the AMOC in the mean state. However, later observational studies suggested that deep waters, including PDW, mainly return to the surface via along-isopycnal pathways in the Southern Ocean (e.g., Marshall & Speer, 2012). In this paradigm, the ITF is relegated to a component of the circum-Australia circulation with a small role in the global overturning circulation (Rousselet et al., 2020;Sloyan & Rintoul, 2001). Thus, most studies of ITF transport variability have focused on local processes in the Indo-Pacific Basin (e.g., Feng et al., 2018;Godfrey, 1996;Sprintall et al., 2019). Although the "conveyor belt" is not an accurate representation of the mean-state global ocean overturning circulation, here we argue that it is a key component of the overturning's transient response to surface forcing perturbations.
The dynamics of transient, interbasin exchange between Atlantic and Indo-Pacific basins were recently discussed by Sun et al. (2020): in response to a weakened AMOC and an associated southward surface transport anomaly ( T ATL < 0), the Indo-Pacific develops a northward surface transport anomaly ( T IP > 0) that opposes changes in the Atlantic. The Indo-Pacific almost fully compensates AMOC changes on decadal to centennial timescales, quantified as a time-dependent interbasin compensation, − T IP / T ATL , that peaks at around 0.8. Modifications to Southern Ocean upwelling that result from and compensate AMOC changes only become important on longer timescales. Here, we extend these results by resolving separate Indian and Pacific basins and show that the Indo-Pacific northward surface transport response occurs almost exclusively in the Indian Ocean. This northward transport anomaly raises Indian Ocean sea level and weakens the ITF. We illustrate the key dynamics using a 1.5-layer reduced gravity model in section 2. Since the reduced gravity model makes a number of simplifications and omits important physical components of the global ocean Geophysical Research Letters 10.1029/2020GL090615 overturning circulation, we also explore whether the proposed dynamics are a robust feature of the ocean circulation in a more comprehensive general circulation model (GCM) simulation in section 3. We diagnose how much of the ITF weakening in the 21st century can be explained by AMOC changes and discuss the intermodel spread in the Coupled Model Intercomparison Project, phase 6 (CMIP6) (Eyring et al., 2016) in section 4. A brief summary is provided in section 5.

Model and Experiment Descriptions
The 1.5-layer reduced gravity model is an idealized representation of the upper branch of the global ocean overturning circulation, defined as the layer above the isopycnal that separates Intermediate Water from Deep Water (see schematic in Figure 1a). Reduced gravity models have proven to be useful tools in guiding theoretical understanding of the controls on the large-scale ocean circulation (e.g., Allison et al., 2011;Johnson & Marshall, 2004;Sun et al., 2020).
The model domain includes three idealized ocean basins representative of the Atlantic, Indian, and Pacific. The total longitudinal extent is 220 • wide, and it extends from 72 • S to 72 • N in latitude. The Southern Ocean is represented by a zonally re-entrant channel between 45 • S and the southern boundary. A 5 • (∼550 km) opening near the equator represents the low-latitude passages that connect the Indo-Pacific basins. The results discussed in this paper are not sensitive to the width of the ITF, remaining essentially unchanged in a simulation with a 3 • -wide channel. The model is forced at the surface by a zonally uniform wind stress ( Figure 1b). The model evolves the upper layer thickness, h(x, y, t), and is discretized on a B-Grid with a horizontal resolution of 1 • × 1 • . Lateral mixing by mesoscale eddies is parameterized as a layer thickness diffusion with diffusivity K GM = 1,000m 2 s −1 (Gent & Mcwilliams, 1990). Interior diapycnal mixing is parameterized as a diapycnal upwelling velocity, w diap = ∕h, with = 2.0 × 10 −5 m 2 s −1 . Surface water mass transformation in the Southern Ocean is represented as a relaxation of the upper layer thickness to 10 m in the hatched area close to the southern boundary ( Figure 1a). The relaxation timescale increases from 10 days at the southern boundary to 100 days at 62 • S. The formation of North Atlantic Deep Water (NADW) is represented as a prescribed constant downwelling velocity, w NADW , in the hatched area close to the northern boundary ( Figure 1a). Details of the model, including the evolution equations and definitions of the transport components, are provided in Text S1, as well as in Sun et al. (2020).
As a control simulation, we prescribe a 12 Sv NADW formation rate, T NADW (Equation S5), and evolve h for 3,000 years to achieve an approximately steady state, defined as global-mean upper-layer thickness changes less than 1 m over 100 years. In this equilibrium state, meridional transports across 30 • S (Equation S6) have the following magnitudes in each basin: Atlantic, T ATL = 11.1 Sv; Indian, T IND = −14.2 Sv; and Pacific T PAC = 12.4 Sv, where positive values are northward. A majority of the 12-Sv NADW formation is balanced by Southern Ocean water mass transformation, which is approximately equal to T ATL + T IND + T PAC = 9.3 Sv, with the remaining due to interior diapycnal upwelling. This partitioning is consistent with the current understanding of the global ocean overturning circulation (e.g., Cessi, 2019;Marshall & Speer, 2012). The equilibrium state also supports an ITF transport (T ITF ) of 13.7 Sv (Equation S7; positive westward). This value is consistent with the estimated zonal ITF transport that would arise from the Island Rule considerations based on the model's wind stress (Text S2).
Using the control run as initial conditions, we conduct two types of simulations in which NADW formation is modified to represent changes in the North Atlantic surface forcing. In the first experiment, we reduce the NADW formation rate, T NADW , from 12 Sv in the control run to 8 Sv and hold it constant, to explore the dynamical processes involved in the overturning circulation's adjustment to this perturbation. In the second set of experiments, we impose time-dependent perturbations to the NADW formation rate. We will focus on the first simulation and briefly discuss the others in section 2.2.

Response of the ITF to NADW Perturbations
An abrupt, step-change, reduction in NADW formation rate leads to a local deepening of the layer interface h in the high-latitude North Atlantic. This interface deepening signal propagates equatorward along the western boundary, eastward along the equator, and poleward (both north and south) along the eastern boundary. In the Southern Hemisphere, the deepening signal moves eastward around the southern tip of the continent into the Indo-Pacific basins (red lines in Figure S1) (e.g., Huang   Here, we focus on the response of the meridional transport across 30 • S in each basin. This geostrophic transport is supported by a change in interface depth between the western and eastern boundaries (e.g., Ferrari et al., 2017;Jones & Cessi, 2016;Thompson et al., 2016).
Kelvin waves propagating from the North Atlantic deepen the interface on the eastern boundary of the South Atlantic, which produces a southward transport anomaly across 30 • S in the Atlantic Ocean, T ATL < 0 (Figures 2a and 2b). In contrast, the deepening of the interface on the western boundary of the Indian Ocean produces a northward transport anomaly in the Indian Ocean, T IND > 0, which is mainly confined to the western boundary current (Figures 2a and 2b). Due to the low-latitude Indo-Pacific passage, the Indian Ocean eastern boundary and the Pacific eastern boundary at 30 • S are connected by Kelvin waves that propagate along the eastern boundaries and the equator, where viscous dissipation can be neglected. The Kelvin waves allow the interface depth on the eastern boundaries of both the Indian Ocean and Pacific Ocean to evolve similarly, such that the transport across 30 • S in the Pacific Ocean remains approximately constant in response to NADW perturbations, i.e., T PAC ≈ 0. This invariant response of the Pacific transport is also consistent with the Island Rule, in which T PAC is constrained by the constant wind stress forcing in the reduced gravity model (Text S2). Therefore, there is a convergence of volume transport into the upper layer of the Indo-Pacific, which is balanced by a deepening of the interface (e.g., Figure 2a; Equation S12).
At timescales longer than the Rossby wave propagation across the Indo-Pacific basin, the interface deepens at roughly the same rate in the Indian and Pacific Oceans, which implies that the changes in ITF transport, T ITF , are linearly proportional to the Indian Ocean transport response, T IND . This can be expressed as

10.1029/2020GL090615
with the ratio r determined by the basin areas ( Figure S2a), Here, S IND and S PAC denote the horizontal area of the Indian and Pacific basins (see derivation in Text S2).
Values of this ratio are r ≈ 0.76 and r ≈ 0.70 for the reduced gravity model and the real ocean, respectively. On decadal to centennial timescales, Sun et al. (2020) showed that the Indo-Pacific transport compensates around 80% the Atlantic changes, i.e., −( T IND + T PAC )/ T ATL ≈ 0.8. Therefore, with T PAC ≈ 0 from above, Equation 1 predicts that the ITF response should be 0.61 times the AMOC changes in the reduced gravity model and 0.56 of the AMOC changes in the real ocean.
The overturning circulation response intensifies over the first two decades, associated with the spin-up of a gyre circulation in the North Atlantic linked to the NADW perturbation ( Figure 2b) (Sun et al., 2020). During this fast response, the Indian Ocean transport anomaly compensates much of the Atlantic changes ( Figure S2b), with the ITF transport response relative to the Atlantic changes T IND / T ATL close to 0.6 ( Figure S2c). The fast response on decadal timescales is followed by a slower adjustment, over millennial timescales, during which the Atlantic southward transport anomaly continues to increase (AMOC continues to weaken), but the Indian northward transport anomaly and the ITF transport anomaly weaken (ITF strengthens) (Figures 2b and S2). As a result, both the transient interbasin compensation, − T IND / T ATL , and T ITF / T ATL decay ( Figure S2). This occurs because the continuous deepening of the upper layer interface on centennial to millennial timescales steepens the interface slope across the Southern Ocean, producing an anomalous southward eddy transport in this region (Sun et al., 2020). The slow Southern Ocean response eventually relieves the burden on the Indo-Pacific to accommodate Atlantic overturning changes, causing the Indo-Pacific to yield some of the transport changes from its peak response. Approaching equilibrium, both the ITF transport and transport into the Indian Ocean eventually recover their values from the control run (Figure 2b; Text S2).
In Text S3, which draws on a number of previous studies (Feng et al., 2018;Godfrey, 1996; Otto-Bliesner & Brady, 2010), we describe a separate set of simulations with time-dependent NADW formation rates. Consistent with the above discussion, we show that both the transient interbasin compensation level and the relative amplitude of the ITF transport anomaly decreases as the NADW forcing period increases ( Figures S3, S4, and S5). Given that the ITF is strongly affected by processes local to the Indo-Pacific Basin on decadal and shorter timescales (e.g., Feng et al., 2018;Godfrey, 1996), we suggest that the imprint of AMOC variability on ITF transport is most pronounced on centennial timescales.

Reduced Gravity Model-GCM Comparison
Despite the idealized nature of the reduced gravity model, the physical processes linking changes in the AMOC to the ITF discussed above appears to be also relevant in more realistic ocean simulations. Specifically, we compare the reduced gravity model to output from an abrupt CO 2 quadrupling (4xCO 2 ) experiment by the NCAR Community Climate System Model, version 4 (CCSM4; Gent et al., 2011), as part of CMIP5 (CMIP5; Taylor et al., 2012). The CCSM4 4xCO 2 experiment is initialized from an approximately equilibrated CCSM4 preindustrial run at year 1850, but with the atmospheric CO 2 instantaneously quadrupled. Both the preindustrial and the 4xCO 2 simulations are continued from 1850 for another 250 years. Throughout this section, we show the difference between the 4xCO 2 experiment and the preindustrial run.
Following Sun et al. (2020), we quantify the AMOC strength at 30 • S as the maximum value of the residual-mean overturning circulation streamfunction at that latitude in the Atlantic Ocean. In order to highlight interbasin exchange, we focus on the AMOC strength at 30 • S rather than in the North Atlantic, although the latter is more commonly used in the literature (e.g., Cheng et al., 2013). We calculate the ITF transport as the difference in the barotropic streamfunction between the coast of Southeast Asia and northwest Australia, which is equal to the total volume transport through all the passages that connect the Indian Ocean and Pacific Ocean. We also quantify the surface meridional transport across 30 • S in each of the three basins as Figure 3. (a) Changes in the AMOC strength at 30 • S (red; Sv) and ITF volume transport (blue; Sv) in the CCSM4 abrupt 4xCO 2 experiments. The gray dash-dotted line represents the ITF volume transport anomaly in the reduced gravity model ("ITF-RGM") in response to a perturbation to the NADW formation rate that is prescribed to follow the maximum value of the AMOC streamfunction in the North Atlantic from the CCSM4 abrupt 4xCO 2 experiment (discussed in section 3). A scatter plot of the ITF volume transport vs. the AMOC strength anomaly from the CCSM4 abrupt 4xCO 2 experiment is provided in the inset to (a), with the linear slope represented as a green straight line. Note that the AMOC strength is slightly different from the surface transport in Figure 2d (blue line). (b) Changes in the upwelling rates in the Pacific basin at 800-m depth (purple) and its contribution due to isopycnal movement (aqua; see definition in Text S4). The difference between the purple and aqua lines is due to changes in diapycnal upwelling.
where the subscript i indicates the basin, z d = 800 m is approximately the maximum depth that connects the Indo-Pacific basins through the ITF passage in this model, x w and x e represent the western and eastern boundary in each basin, and v represents the residual velocity that includes both the Eulerian-mean velocity and the parameterized eddy bolus velocity.
In response to the abrupt 4xCO 2 , the AMOC weakens roughly from 20 to 12 Sv during the first 100 years, followed by a partial recovery to around 16 Sv during the next 150 years (red line in Figure 3a). Consistent with the AMOC influencing the ITF, the ITF transport covaries with the AMOC on centennial timescales (Figure 3a). The ITF also undergoes strong interannual and decadal fluctuations, likely forced by local processes (e.g., surface forcing) within the Indo-Pacific Basin (e.g., Feng et al., 2018;Godfrey, 1996). The ratio of ITF to AMOC transport changes is 0.51 (inset of Figure 3a), slightly lower than the estimated 0.56 from the reduced gravity model in section 2.2. This overestimate in the reduced gravity model is likely related to the deeper depth of the AMOC maximum streamfunction as compared to the ITF, which is not resolved by the 1.5-layer model. We carry out an additional reduced gravity model simulation, in which we prescribe the NADW formation rate (Equation S12) using the maximum value of the AMOC streamfunction in the North Atlantic below 500 m in the CCSM4 abrupt 4xCO 2 experiment. The ITF volume transport anomaly from this reduced gravity model simulation (gray dash-dotted line in Figure 3a) largely reproduces the centennial ITF changes in the CCSM4 abrupt 4xCO 2 experiment (blue line in Figure 3a).
The similar isopycnal structure between the CCSM4 and the reduced gravity model simulations provides confidence that the ITF is connected to AMOC changes via the same dynamical processes discussed in section 2.2 (Figures 2 and S6). Associated with the weakened AMOC following the 4xCO 2 , the density anomaly along 30 • S in the Atlantic Ocean between 1 and 3 km depths has a zonal gradient consistent with a deepening of isopycnals on the eastern boundary ( Figure 2c) and an anomalous southward surface transport (blue line in Figure 2d). This isopycnal deepening signal on the eastern boundary radiates into the interior via Rossby Waves on decadal timescales and causes isopycnal deepening in the Atlantic interior that weakens westward (Figures 2a and 2c). Consistent with the propagation of coastal Kelvin waves into the Indian Ocean, the density anomaly along 30 • S in the Indian Ocean has a strong zonal gradient that is largely confined to the western boundary (inset of Figure 2c), associated with an anomalous northward surface transport (orange line in Figure 2d). In contrast, there are only weak gradients in density anomaly in the Pacific, as constrained by the basin-scale wind stress forcing through the Island Rule ( Figure S7).
With an approximately invariant surface meridional transport across the southern boundary in the Pacific (green line in Figure 3d) and a response in the Bering Strait transport by less than 0.5 Sv in CCSM4, the centennial changes in the ITF transport can only be balanced by a change in the Pacific upwelling based on volume conservation. Previous studies have attributed this change in the Pacific upwelling to diapycnal processes (Feng et al., 2017;Sen Gupta et al., 2016). However, in the reduced gravity model, the reduced ITF transport is mainly balanced by an adiabatic deepening of the interface in the Pacific (Equation S12). In the CCSM4 experiment, we perform a similar assessment by partitioning the changes in the Pacific upwelling (w) into an adiabatic component due to isopycnal movement (w isop ) and a diabatic component due to changes in diapycnal mixing (w diap ): where and with b the buoyancy of seawater and diapycnal diffusivity. This is discussed in details in Text S4, which draws on (Munk, 1966;Sun et al., 2018). Consistent with the reduced gravity model, we find that the changes in the Pacific upwelling needed to balance the ITF changes are approximately equal to the adiabatic components at 800-m depth (Figure 3b). The adiabatic component exceeds the diabatic component throughout the water column, and the latter makes a significant contributions only at depths greater than 1,500 m. This suggests that much of the subsurface changes in tracer (e.g., temperature and salinity) distributions on centennial timescales can be attributed to horizontal interbasin exchanges rather than mixing between different water masses in the vertical (e.g., Huang, 2015).

ITF and AMOC Changes in CMIP6
Both the reduced gravity model and CCSM4 4xCO 2 experiment highlight the connections between the AMOC and ITF on centennial timescales. In response to a weakened AMOC, our results suggest a decline in the ITF transport ( T ITF ) that is around half of the AMOC strength changes ( T AMOC ), i.e., T ITF ≈ 0.5 T AMOC (Figure 3a). Now, we use this relationship to estimate how much of the ITF weakening during the 21st century can be explained by the AMOC changes by examining the CMIP6 simulations under the high-end emission scenario ("SSP585").
In response to the increasing greenhouse gas forcing, both the AMOC and ITF transports weaken between 2015 and 2100 ( Figures S8 and S9) in all the CMIP6 models. Yet, these changes in AMOC strength and ITF transport, diagnosed as a difference between 2015-2024 and 2091-2100 (Figure 4), have significant intermodel spread among the 15 CMIP6 models analyzed in this study. The AMOC changes explain around 100% of the ITF transport weakening during the 21st century in five models, but this percentage is only around 40% in the other models (Figure 4). This intermodel spread implies differences in their simulated surface forcing changes, which account for additional changes to the ITF transport. For example, Sen Gupta et al. (2016) show that there is a strong intermodel spread in wind changes in CMIP5, such that a portion of the centennial changes in the ITF transport may be explained by the wind. While application of the Island Rule should help with the attribution of changes in the ITF transport to wind stress variations, in practice, we found that due to the complicated continental geometry and bathymetry in the GCMs, the results are highly sensitive to small changes in the integral path. We leave a detailed analysis of this sensitivity to future work.

Summary and Discussion
Climate models consistently project a robust weakening in the AMOC and the ITF during the 21st century in response to greenhouse gas forcing. Here, we propose that the ITF is dynamically linked to the AMOC, and the latter is a primary driver of changes in ITF transport on centennial timescales. In a previous study, Sun et al. (2020) showed that there is a transient overturning compensation between the Atlantic and Indo-Pacific basins. In a warming climate, the AMOC weakens, but the Indo-Pacific develops an opposing overturning circulation anomaly, characterized by an anomalous northward surface transport. This earlier study neglected the potential for low-latitude exchange between the Indian and Pacific basins and the key role of the ITF. By resolving this additional transport pathway, we show that the Indo-Pacific northward surface transport anomaly is almost exclusively confined to the Indian Ocean. The Pacific surface transport is instead constrained by the basin-scale wind stress through the Island Rule. This asymmetry in surface transport between the Indian and Pacific basins has a direct impact on the zonal sea surface height differences between basins, a mechanism that explains the weakening ITF. This asymmetry is also important for constraining oceanic heat and dissolved gas budgets as transport into the northern basins across 30 • S is fed by Antarctic Intermediate Waters responsible for significant uptake of heat (Armour et al., 2016) and carbon dioxide (Gruber et al., 2019).
This transient version of the ocean's conveyor belt circulation provides an oceanic pathway for changes in the high-latitude North Atlantic to affect the low-latitude Indo-Pacific, a teleconnection that could play an important role in regulating the climate system. In response to reduced NADW formation, a reduced ITF transport converges more heat into the Pacific Ocean (e.g., Garuba & Klinger, 2016). Effectively, a weakened ITF and the associated deepening of the isopycnals in the Pacific Basin provide more heat below the mixed layer that could modify the tropical atmosphere-ocean interactions and boost occurrence of extreme El Niño/La Niña events (e.g., Cai et al., 2015).
Our results suggest an ITF transport response that is around half of the AMOC changes. An intriguing corollary of this relationship is the potential to use the ITF (e.g., Susanto & Song, 2015) to monitor or interpret long-term trends in the overturning circulation. However, a diagnosis of the CMIP6 simulations finds an intermodel spread with regard to the relative magnitude of the ITF's weakening in response to AMOC changes during the 21st century. While this study has emphasized the dynamics that enable remote forcing to influence the ITF, the relative importance of local (e.g., wind and surface buoyancy forcing) and remote processes, and why they might differ between models, requires further study.