Atlantic Warming Since the Little Ice Age

s. Technical Report 14-04, Ministry of Climate and Energy, Copenhagen, Denmark, 90 pp. Deacon, G.E.R. 1937. The Hydrology of the Southern Ocean. Discovery Reports Series 15, Cambridge University Press, 124 pp. DeVries, T., and F. Primeau. 2011. Dynamically and observationally constrained estimates of water-mass distributions and ages in the global ocean. Journal of Physical Oceanography 41(12):2,381–2,401, https://doi.org/ 10.1175/JPO-D-10-05011.1. Dickson, R.R., H.H. Lamb, S.-A. Malmberg, and J.M. Colebrook. 1975. Climatic reversal in northern North Atlantic. Nature 256:479–482, https://doi.org/ 10.1038/256479a0. Dickson, R.R., J. Meincke, S.A. Malmberg, and A.J. Lee. 1988. The “great salinity anomaly” in the northern North Atlantic 1968–1982. Progress in Oceanography 20:103–151, https://doi.org/ 10.1016/ 0079-6611(88)90049-3. Durack, P.J., P.J. Gleckler, F.W. Landerer, and K.E. Taylor. 2014. Quantifying underestimates of long-term upper-ocean warming. Nature Climate Change 4(11):999–1,005, https://doi.org/10.1038/ nclimate2389. Fukumori, I. 2002. A partitioned Kalman filter and smoother. Monthly Weather Review 130(5):1,370–1,383, https://doi.org/ 10.1175/ 1520-0493 (2002)130 <1370:APKFAS>2.0.CO;2. Gebbie, G., and P. Huybers. 2010. Total matrix intercomparison: A method for resolving the geometry of water-mass pathways. Journal of Physical Oceanography 40(8):1,710–1,728, https://doi.org/ 10.1175/2010JPO4272.1. Gebbie, G., and P. Huybers. 2011. How is the ocean filled? Geophysical Research Letters 38(6), https://doi.org/ 10.1029/2011GL046769. Gebbie, G., and P. Huybers. 2012. The mean age of ocean waters inferred from radiocarbon observations: Sensitivity to surface sources and accounting for mixing histories. Journal of Physical Oceanography 42(2):291–305, https://doi.org/ 10.1175/ JPO-D-11-043.1. Gebbie, G., and P. Huybers. 2019. The Little Ice Age and 20th-century deep Pacific cooling. Science 363(6422):70–74, https://doi.org/10.1126/ science.aar8413. Giese, B.S., and S. Ray. 2011. El Niño variability in simple ocean data assimilation (SODA), 1871– 2008. Journal of Geophysical Research 116(C2), https://doi.org/ 10.1029/ 2010JC006695. Giese, B.S., H.F. Seidel, G.P. Compo, and P.D. Sardeshmukh. 2016. An ensemble of ocean reanalyses for 1815–2013 with sparse observational input. Journal of Geophysical Research 121(9):6,891–6,910, https://doi.org/ 10.1002/ 2016JC012079. Gleckler, P.J., P.J. Durack, R.J. Stouffer, G.C. Johnson, and C.E. Forest. 2016. Industrial-era global ocean heat uptake doubles in recent decades. Nature Climate Change 6:394–398, https://doi.org/ 10.1038/nclimate2915. Gouretski, V., and K. Koltermann. 2004. WOCE Global Hydrographic Climatology. Technical Report 35, Berichte des Bundesamtes für Seeschifffahrt und


INTRODUCTION
The Little Ice Age (LIA) was a period of cold surface temperatures over roughly the years 1400−1800 CE (e.g., Paasche and Bakke, 2010), the time of cold European winters immortalized in the works of European master painters.A period of modern warming has followed the Little Ice Age, but research is still being conducted to determine what fraction of modern warming is due to natural versus anthropogenic causes (e.g., Abram et al., 2016).To determine what physical processes are at play, key information includes the timing of the end of the Little Ice Age, the magnitude of modern warming, and whether these quantities are globally coherent or have regional variations (Bradley and Jones, 1993).Measurements of ocean temperature using modern instruments only go back as far as the 1870s HMS Challenger expedition (Murray, 1895), however, and probably do not capture the coldest sea surface conditions.Paleoceanographic records recovered from sediment cores provide information about sea surface temperature over the last millennium (e.g., Oppo et al., 2009), but their utility is sometimes limited by dating uncertainty and lack of temporal resolution.
Modern-day radiocarbon observations indicate that waters in the deep ocean were last in contact with the atmosphere several hundred years ago (Key et al., 2004).Atlantic deep waters are not as isolated from the surface as Pacific waters due to vigorous deep water formation, but their radiocarbon concentrations suggest radioactive decay over 200−500 years since they had atmospheric values.This timescale is long enough for some of the deep waters of the Atlantic to have last been in contact with the sea surface during the Little Ice Age.That the Atlantic could have such long timescales may be surprising, as the response of the Deep Western Boundary Current (DWBC) to surface perturbations is a few decades or less (e.g., Jackson et al., 2016), but this DWBC response is just the arrival of an initial signal, and the full equilibrium response takes much longer (Khatiwala et al., 2001;Gebbie and Huybers, 2012).In order to simulate the post-LIA deep ocean evolution accurately, a numerical model may require knowledge of the long-term surface history, or the initial conditions may need to contain the disequilibrium effects that result from this surface history.
The goal of this manuscript is to review what is known about centennialtime scale surface temperature changes and to use that information to set expectations about warming in the subsurface Atlantic since the end of the Little Ice Age.The relevant surface temperature history includes times that precede the instrumental record, and thus we reconstruct surface temperature by taking into account both instrumental and paleoceanographic proxy observations.We explore how well interior ocean variations on multidecadal and longer timescales can be inferred by the passive advection and diffusion of surface temperature anomalies.Due to the long memory of the deep ocean, the signal of the Little Ice Age occurs later in the subsurface Atlantic than at the surface.The remnant effect of the Little Ice Age has consequences for both the Atlantic heat uptake and air-sea fluxes in the modern warming period.

COMMON ERA SEA SURFACE TEMPERATURE EVOLUTION Global-Mean Variations
Here, we reconstruct a plausible surface temperature history of the Common Era (i.e., the last 2,000 years) from the Ocean2k paleoceanographic proxy data compilation (McGregor et al., 2015) and the HadISST 1.1 product derived from modern thermometers (Rayner et al., 2003).For the years 1870−2015, the average sea surface temperature (SST) of the three coldest months in the HadISST product is selected because these waters are preferentially communicated to the ocean interior by the mixedlayer "demon" (Stommel, 1979;Williams et al., 1995).This quantity should be thought of as the "subduction temperature, " or the temperature of subducted water, rather than a strict reconstruction of SST.Temperature reconstructions from late nineteenth century instrumental data are especially uncertain, as the ABSTRACT.Radiocarbon observations suggest that the deep Atlantic Ocean takes up to several centuries to fully respond to changes at the sea surface.Thus, the ocean's memory is longer than the modern instrumental period of oceanography, and the determination of modern warming of the subsurface Atlantic requires information from paleoceanographic data sets.In particular, paleoceanographic proxy data compiled by the Ocean2k project indicate that there was a global cooling from the Medieval Warm Period to the Little Ice Age over the years 900−1800, followed by modern warming that began around 1850.An ocean simulation that is forced by a combined instrumentalproxy reconstruction of surface temperatures over the last 2,000 years shows that the deep Atlantic continues to cool even after the surface starts warming.As a consequence of the multicentury surface climate history, the ocean simulation suggests that the deep Atlantic doesn't take up as much heat during the modern warming era as the case where the ocean was in equilibrium at 1750.Both historical hydrographic observations and proxy records of the subsurface Atlantic are needed to determine whether the effects of the Little Ice Age did indeed persist well after the surface climate had already shifted to warmer conditions.
HadISST data set shows increased divergence from other instrumentally based products during this period (e.g., Huang et al., 2015).The Ocean2k project reconstructed global-mean sea surface temperature for the years 100−1900 with 200-year resolution, and its statistical tests show a coherent global signal (McGregor et al., 2015).This proxy reconstruction must be interpreted carefully because the proxy data may best reflect spring and summer temperatures, the temperatures are anomalies that are not placed on an absolute temperature scale, and no estimate of twentieth century warming is provided.To produce a seamless simulation of the Common Era, these three issues must be addressed.
To address potential seasonal biases, we average the 57 individual marine paleoceanographic records in the Ocean2k SST network to produce a temperature anomaly time series up to 1950.By computing anomalies, any warm signal due to a spring or summer bias in the proxies is removed so long as there are no changes in seasonality.Furthermore, it is assumed that this time series well approximates the variations in winter temperature that subduct into the subsurface.
An improved estimate of modern warming can be produced by blending the instrumental and proxy data, but first the two data sets must be placed on the same absolute temperature scale.The missing quantity is the Common Era-average surface temperature that needs to be added to the proxy-data anomaly time series.We solve for the Common Era-average temperature that leads to a match with the global-mean surface temperature difference between HMS Challenger expedition measurements of the 1870s and the World Ocean Circulation Experiment (WOCE) data of the 1990s, where the instrumental and proxy data sets are blended with a linearly varying weight from 1870 and 1950.The weighting is such that the proxy record has 100% weight at 1870 and the instrumental product is given 100% weight at 1950.This process yields a single blended surface temperature history in which signals from the instrumental and proxy data sets are compared in a consistent manner.
The global-average subduction temperature is reconstructed to have cooled by ~0.5°C between the Medieval Warm Period (~900 CE) and the end of the Little Ice Age, followed by almost 0.9°C of modern warming since 1850 (Figure 1).The coldest subduction temperatures of the Common Era occurred in 1750 and again around 1850, suggesting that the Little Ice Age was not a stable cold climate, but instead was characterized by variability around a cold state.It is difficult to define the Medieval Warm Period because of decadal and centennial variability in the years 600−900.This variability arises in the average of 57 paleoceanographic proxy data sets but is not emphasized in the reconstruction provided by the Ocean2k project, which instead provides a smoother time series by averaging over 200-year intervals.Our reconstruction is probably best thought of as one possible realization of this climate interval, where the timing and magnitude of the multidecadal and centennial events are uncertain.
The global temperature reconstruction shows a number of strong warmings, such as those beginning in 1850, 1920, and 1970 (Figure 2).There is conflicting evidence about the physical mechanisms responsible for these events.The event that begins in the 1920s is detected in several surface meteorological measurements in Greenland as a rapid 4°C warming (Cappelen, 2014).A climate model study suggests that stochastic atmospheric forcing can cause Greenland climate transitions (Kleppin et al., 2015) along with other climate signals around the globe, including the tropical Pacific (Giese and Ray, 2011).This proposed mechanism does not require anthropogenic forcing.Furthermore, the magnitude of the earliest warming that starts around 1850 is uncertain.Paleoceanographic data have been interpreted such that this warming is a consequence of the early onset of human activities (Abram et al., 2016), but its magnitude is usually FIGURE 1. Reconstructed surface temperature, θ b , from blended paleoceanographic and instrumental data products.Surface temperature time series are relative to a baseline value at the year 15 CE (i.e., EQ-0015), including the global area-weighted average (black line) and seven major surface regions of the Atlantic sector (colored lines).The surface regions are defined in the text.Prior to 1870, the regional anomalies collapse onto the global mean.A general definition of the Little Ice Age is included (LIA, horizontal bar; Paasche and Bakke, 2010).ARC = Arctic.MED = Mediterranean.WED = Weddell Sea.LAB = Labrador Sea.GIN = Greenland-Iceland-Norwegian Sea.SUBANT = The subantarctic region of the Atlantic.TROP = The remaining subtropical and tropical regions of the Atlantic.reconstructed to be smaller than shown in Figures 1 and 2. Many of the Coupled Model Intercomparison Project version 5 (CMIP5) models are initialized around this time period and do not simulate any significant late nineteenth century warming (Gleckler et al., 2016).Here, the strong warming is corroborated by the surface HMS Challenger data, and additional historical hydrography from the Norwegian Sea could provide additional evidence (Helland-Hansen and Nansen, 1909).

Regional Temperature Variations
We reconstruct regional temperature variations after 1870 based upon the HadISST product.To mitigate the uncertainty of the small-scale structures in the HadISST product in the late nineteenth and early twentieth centuries, the surface temperature time series are regionally averaged in 14 oceanographically defined surface patches following Gebbie and Huybers (2010).Of the 14 regions, we consider seven to be the focus of this work: (1) Arctic (ARC), ( 2) Mediterranean (MED), (3) Weddell Sea (WED), (4) Labrador Sea (LAB), (5) Greenland-Iceland-Norwegian Sea (GIN), ( 6) the subantarctic region of the Atlantic (SUBANT), and (7) the remaining subtropical and tropical regions of the Atlantic (TROP).The Weddell Sea region is defined to be south of the southern extent of the Antarctic Circumpolar Current at a boundary given by the σ o = 27.55 line (i.e., surface density of 1,027.55 kg m -3 ; Sievers and Nowlin, 1984;Orsi et al., 1995), and bounded zonally by 70°W and 30°E.The subantarctic and the subtropics are divided by the 34.8 isohaline on the practical salinity scale that approximates the subtropical front (Deacon, 1937).The Labrador Sea region is defined to be north of the polar front marked by the 35.4 isohaline (Tomczak and Godfrey, 1994) and south of the Greenland-Scotland Ridge, and thus includes much of the North Atlantic subpolar gyre.The GIN Sea region extends northward from the LAB region as far as 82°N and 30°E, and any waters beyond those boundaries are considered Arctic.The Mediterranean is defined as all waters interior to the Strait of Gibraltar.Over 82% of the original variability in the HadISST product is retained after averaging, and we suggest that the remaining signals are more reliable due to averaging over large regions, although some regions such as the Arctic are severely data limited.
Temperature in the Labrador Sea region has not monotonically increased since the end of the Little Ice Age, but has been influenced by strong multidecadal variability ("LAB, " Figure 2).This inference is consistent with the annual time series of temperature, salinity, and density of the central Labrador Sea that is available back to 1938 (Yashayaev and Clarke, 2008).The suggested cause of the multidecadal variability is intermittent deep convection that creates large quantities of Labrador Sea Water (e.g., Lazier, 1980;McCartney and Talley, 1982).Temperature change is often coincident with changes in salinity and other North Atlantic climatic variables (e.g., Dickson et al., 1975Dickson et al., , 1988)).When considering Labrador Sea surface temperature over the entire 1850−2015 interval, however, the warming has a similar magnitude as the global-mean signal.Our reconstruction suggests significant Labrador Sea warming prior to 1940, consistent with recent inferences from paleoceanographic data (Thornalley et al., 2018).An outstanding question is the extent to which the North Atlantic temperature variability is related to changes in the meridional overturning circulation, global warming, or both (e.g., Caesar et al., 2018).
Paleoceanographic proxies of the last 2,000 years show coherent global changes, but a recent analysis did not detect significant differences in sea surface temperature when compiled into regional bins (McGregor et al., 2015).For this reason, our reconstruction of surface temperature prior to 1870 makes the first guess that all surface regions covary with globalmean surface temperature (collapsed lines, Figure 1).Such an assumption does not imply that the sea surface was homogeneous prior to 1870, but instead that the spatial pattern is fixed.The reconstructed regional anomalies are small relative to the background temperature differences between regions, so the present-day spatial pattern is only slightly modified.Additional information is available in subsurface data from both paleoceanographic proxies (e.g., Mjell et al., 2016) and historical instrumental observations (e.g., Murray, 1895).There are thought to be shifts in the North Atlantic Oscillation, El Niño-Southern Oscillation, and the Pacific Decadal Oscillation that covary with the Medieval Warm Period and the Little Ice Age (e.g., IPCC, 2005).Efforts to fit a model to these types of observations have been successful over the time period 1815−2013 (Giese et al., 2016).It would be especially useful to produce an estimate that is also consistent with the pre-1815 surface history, and to make an estimate that conserves heat.To do so, an optimal control problem needs to be solved (e.g., Fukumori, 2002;Wunsch and Heimbach, 2007), but computational barriers to such a method would have to be overcome.Even in cases where surface temperature varies on monthly timescales, for example, the stable numerical integration of an ocean model typically requires a sub-hourly time step.At this level of time resolution, the equations for either the adjoint method or the Kalman smoother, two methods for solving the control problem, become computationally intractable.

SIMULATION OF THE COMMON ERA Circulation Model
The interior ocean response to surface temperature variability includes the passive advection and diffusion of anomalies along the existing large-scale circulation, and a dynamic response that alters the circulation through the temperature effect on density (e.g., Banks and Gregory, 2006).Surface temperature variability may reflect a change in water mass formation rates, and this dynamical signal is communicated rapidly throughout the deep ocean (Kawase, 1987).In a realistic simulation with an ocean general circulation model, however, Marshall et al. (2015) find that anthropogenic warming is communicated to the interior on the centennial timescale primarily through the passive response to surface anomalies.The large-scale patterns of ocean heat uptake and storage are primarily controlled by the advection and diffusion of a background circulation, and changes in ocean circulation play a secondary role so long as warming is sufficiently small.Our reconstructed surface history indicates temperature perturbations smaller than 1°C and even smaller subsurface anomalies due to ocean mixing.In addition, comparison of hydrographic observations taken decades apart shows only minor changes in density perturbations (Roemmich and Wunsch, 1984).As the large-scale ocean circulation is in balance with the density structure through geostrophy, these observations suggest that the large-scale circulation did not undergo major reorganizations in the recent past (see also, Zanna et al., 2019).This assumption should be reassessed with additional data that are now available, and the potential effect of changes in circulation is later tested in the section on Penetration of Surface Anomalies.
Our goal here is to set expectations for the multidecadal to centennial temperature changes that characterize Atlantic warming after the Little Ice Age.Given the suggestion of stable circulation patterns, we represent interior ocean circulation by the advective and diffusive fluxes empirically derived from observations collected during the WOCE field campaign of the 1990s (Gebbie and Huybers, 2012).Distributions of temperature, salinity, nutrients, oxygen, and radiocarbon were inverted for the circulation and mixing rates while respecting conservation equations for mass and all other properties.Circulation is represented at 2° × 2° horizontal resolution across 33 vertical layers, leading to 291,556 grid cells.At this resolution, diffusive fluxes include the net effect of subgridscale processes such as mesoscale eddies.
The ocean circulation model was expressly derived from WOCE data to accurately represent connections between the surface and the deep.A demonstration of the model's skill is that it accurately reproduces independent carbon isotope values (Gebbie and Huybers, 2011) and deep-ocean radiocarbon concentrations while producing realistic estimates of the fraction of Antarctic Bottom Water and North Atlantic Deep Water in the deep ocean (Gebbie and Huybers, 2010;DeVries and Primeau, 2011).We suggest that issues inherent to models that parameterize small-scale processes are mitigated here by empirically training the model with data.
By using this simplified model of ocean physics, we can simulate back far enough into the past to capture the disequilibrium processes due to the Little Ice Age.While this manuscript focuses on the Atlantic Ocean, a global model domain is used because of the remote influences that can affect the Atlantic on these long timescales.We refer to the simulation as EQ-0015 because it is a 2,000-year simulation that assumes equilibrium in year 15 of the Common Era.Even starting at the year 15 CE does not completely remove sensitivity to initial conditions.In the Atlantic, there is little memory from earlier than 15 CE, as the percentage of water with age greater than 2,000 years is nowhere greater than 6%.Due to uncertainties in small-scale temperature patterns during the Common Era, we force the simulation with the average temperature in 14 surface regions described earlier.The surface boundary conditions are also averaged with a five-year timescale.In order to accurately resolve tracer transport, we use an adaptive time step that is much shorter than the five-year timescale of the surface boundary conditions and is adjusted depending upon spatial gradients.
The evolution of subsurface temperature is also affected by potential changes in the winds, tides, and freshwater forcing that are not considered in the model.Variations in wind forcing, in particular, are known to be transmitted to the deep ocean more rapidly than the advective-diffusive processes in the model.Unfortunately, the magnitude of centennial-scale variability in winds is not known, nor is the sensitivity of internal redistributions of heat to these potential wind shifts.As a major effect of winds is to set off planetary waves that lead to isopycnal heave and a redistribution of ocean waters, large-scale averages can be taken to mitigate the net effect of wind forcing (Roemmich et al., 2015).Here we emphasize that the dynamics of our simulation are expected to be most accurate on the largest spatial and temporal scales.

Penetration of Surface Anomalies
Our 2,000-year simulation permits centennial-scale deep Atlantic temperature anomalies to be coherently tracked back to surface climate variability (Figure 3).The signals of the warm anomaly of the Medieval Warm Period, the cold anomaly of the Little Ice Age, and modern warming all propagate from the surface to the seafloor, as can be seen when an Atlanticwide average temperature profile is diagnosed.On the scale of the Common Era, the penetration of these anomalies appears rapid, but the deep Atlantic generally lags the surface by a century or more.The range of deep Atlantic temperature over the last 2,000 years is expected to be greater than 0.2°C, consistent with the relatively rapid propagation to depth being only slightly damped by ocean mixing processes.
The three-dimensional, global temperature distribution at the beginning of the Common Era is not known, but this field is not needed to initialize the simulation under the assumptions given above.To initialize the simulation in the year 15, the temperature anomaly vanishes under the assumption of equilibrium.Then, the surface boundary conditions are translated to anomalies relative to the year 15.The evolution of the full temperature field can be backed out after the anomaly simulation is completed by calibrating with the WOCE-era (1990s) temperature distribution (Gouretski and Koltermann, 2004).
According to the simulation, the amount of Atlantic-average warming since the Little Ice Age varies strongly as a function of depth.The coldest conditions of the Common Era occur later as depth increases, and the deep Atlantic doesn't begin warming until about 1950.Thus, the simulation suggests that there has been a long-term warming trend at all depths in the Atlantic since at least the mid-twentieth century.The amount of warming is larger than 0.5°C at the surface, and diminishes to less than 0.05°C at the seafloor.
The simulation indicates that cold anomalies persist in the deep Atlantic below 2,200 m depth until the end of the simulation.In this case, cold anomalies do not indicate a present-day cooling trend, but rather that the Atlantic is cooler today than in the year 15 at those depths.If the advective-diffusive processes of the model are representative, the cumulative effect of surface cooling from the Medieval Warm Period to the end of the Little Ice Age is expected to have some signal today.The signal of warm anomalies, however, is penetrating rapidly into the deep Atlantic, especially after 1980.
The meridional overturning circulation may have changed by up to ±25% during the Common Era (Lund et al., 2006;Rahmstorf et al., 2015), and the potential influence of these changes is not captured in the previous simulation of this work.To address this issue, we note that the simulation of temperature anomalies is accomplished through a set of linear equations.Under this assumption, the potential temperature at an interior point is a linear function of the boundary temperatures at previous times.Only the surface-to-interior lag and water-mass decomposition are important for determining the response.These linear functions are often called the multiple-source boundary propagator (Haine and Hall,FIGURE 3. Simulation of Atlantic Ocean subsurface temperature anomaly.Time evolution of the vertical profile of Atlantic-average potential temperature anomaly relative to a baseline value at the year 15 CE.The units of the contour labels are hundredths of a degree Celsius.2002) or the boundary Green's function (Wunsch, 2002).A uniform speedup of the circulation acts to contract the timeaxis of the boundary Green's functions, leading primarily to decreased lag times in the interior Atlantic.We explore the effect of a time-varying circulation by performing simulations of the Common Era where the globalmean surface temperature is assumed to linearly covary with ocean circulation strength (not shown).To approximate how this forcing might affect ocean circulation, we performed one simulation in which the Little Ice Age circulation was 25% slower than found for the 1990s, and another that was 25% faster.Despite the large opposing effects of the time-varying circulation in these two additional simulations, Common Era temperature evolution at 3,500 m depth shows the telltale signs of the Medieval Warm Period, the Little Ice Age, and modern warming, in both cases much like Figure 3.The time-varying circulation modifies the timing of the coldest conditions from the Little Ice Age and the rate of vertical propagation.In all simulations presented here, the present-day Atlantic is warming at all depths, although there are some slight differences in the rate of warming in the deep Atlantic.

Spatial Pattern of Deep Atlantic Warming
The expected spatial pattern of deep Atlantic warming is highlighted by diagnosing the simulated temperature difference between two decades.The 1870s are chosen due to the maximum depth-integrated effect of the Little Ice Age at that time, and the 1990s are chosen because they represent the most recent top-to-bottom global assessment of ocean temperatures (Gouretski and Koltermann, 2004).Greater warming is indicated in the North Atlantic relative to the South Atlantic, as well as the western basin relative to the eastern basin (Figure 4).Although the model is only forced by surface data, the simulation suggests a warming of North Atlantic Deep Water and an enhanced propagation of the signal along the Atlantic deep western boundary current.Ocean reanalyses show a similar pattern of warming over the last few decades (Palmer et al., 2017), indicating that this ocean response is similar in more complex general circulation models.In this study, the model produces a deep western boundary current with a realistic boundary tracer plume (Figure 4).The core of strong currents in the deep western boundary is smaller than the 2° × 2° horizontal resolution of our model, but the net effect of tracer transport along the boundary is enforced through fitting several million modern-day ocean observations in the derivation of the empirical circulation.
While the simulation shows Atlantic warming nearly everywhere after 1870, sediment core data from an intermediate depth on the West African margin shows cooling and strong interdecadal variability (Morley et al., 2011).This cooling trend may be a continuation of the cooling after the Holocene climatic optimum as observed in both Atlantic (Morley et al., 2014) and Pacific (Rosenthal et al., 2013) cores.This cooling could be explained by an intermediate water mass with source waters that do not warm in step with the rest of the ocean, and suggests that the model assumption of globally coherent surface temperature anomalies should be relaxed as in the study of Gebbie and Huybers (2019).

ATLANTIC ENERGY BALANCE Atlantic Heat Uptake
Ocean heat uptake is diagnosed by taking the change in ocean heat content over time.Ocean heat content is most accurately computed as a function of the thermodynamic property of potential enthalpy (McDougall, 2003).
Here, we calculate ocean heat content by taking into account the density, specific heat capacity, potential temperature, and volume of the 291,556 grid cells of the model.Our diagnostic is likely to have errors that approach 1% due to the assumption of a fixed density field through time, and the replacement of Conservative Temperature with the model variable, potential temperature.Other assumptions made in deriving the passive temperature anomaly model likely cause greater errors.Globally, about 160 ZJ (1 ZJ 10 21 J) of heat are taken up in the upper 700 m of the model between 1970 and 2005, consistent with observational estimates (Levitus et al., 2012) and the CMIP5 simulations (Durack et al., 2014).As our model is not as complex as those of the CMIP5 project, reconstructing the exact amount of heat uptake is not the focus here.
After the Little Ice Age, most of the Atlantic heat uptake occurred above 700 m depth, followed in importance by mid-depth (700−2,000 m) and deep (2,000−4,000 m) layers (Figure 5).The order of importance is related to the response time of each layer, as the upper layer was the first to increase heat content.The deep layer did not uptake heat until 1870.Despite the late start and the relatively small size of the Atlantic basin, the deep Atlantic layer uptakes a significant 70 ZJ of heat after 1870.When integrating below 700 m depth, the model simulates 180 ZJ of heat uptake over 130 years in the Atlantic.This corresponds to a nearly 0.1 W m -2 energy flux over the entire surface of the planet and is a significant term in the global energy budget.
Changes in heat content in the deep Atlantic were larger than in the surface ocean prior to 1750.This importance of the deep layer is due to the long timescale of the Little Ice Age cooling trend, where the deep Atlantic had enough time to fully respond to the surface.Our findings highlight how the memory of the ocean preserves the signal of temperature variability in intermediate and deep waters, as previously inferred from paleoceanographic observations over the Holocene (Rosenthal et al., 2017).
The EQ-1750 simulation describes a counterfactual case where the ocean is prescribed to be in equilibrium in the year 1750.After 1750, the surface conditions follow those in the EQ-0015 simulation.To accomplish this, the 1750 initial conditions are set to zero anomaly everywhere, and the post-1750 surface conditions are input to the model as anomalies relative to 1750.Differences in EQ-1750 and EQ-0015 heat content evolution reflect the memory of surface conditions before 1750 (dashed lines, Figure 5).In the upper layer, only small differences exist because the ocean memory is relatively short and the surface history before 1750 is quickly forgotten.In the middepth and deep layers, however, the heat uptake in EQ-1750 is greater than in EQ-0015.We hypothesize that such a bias in heat uptake may exist in any model that is initialized at the end of the Little Ice Age, including the state-of-the-art CMIP5 models (Gleckler et al., 2016).The decrease in heat uptake due to the pre-1750 surface history can be understood in the context of Figure 3.In 1750, the EQ-0015 simulation shows cold anomalies penetrating into the deep Atlantic.The EQ-1750 simulation instead starts with a period lacking any deep ocean changes due to being initialized at equilibrium and the slow response to surface conditions.Consequently, the EQ-1750 heat content doesn't change much for almost a century.During this period, the EQ-0015 and EQ-1750 simulations diverge due to the subsurface response to the LIA-related surface cooling.
This analysis focuses on heat uptake since 1750 by independently setting the temperature anomalies in each simulation to a 1750 baseline, but it doesn't provide any information about the temperature difference between EQ-0015 and EQ-1750.The year 1750 has nearly the coldest surface temperature of the Common Era, and thus the initial conditions of EQ-1750 reflect these cold conditions throughout the globally equilibrated ocean.Consequently, EQ-1750 is colder than EQ-0015 at all times (1750−2015).0−700 m 700−2,000 m 2,000−6,000 m FIGURE 5. Simulated Atlantic heat content changes, 15−2015 CE.Time series of the Atlantic heat content anomaly relative to 1750 CE as decomposed into the upper (magenta, 0-700 m), mid-depth (green, 700-2,000 m) and deep (orange, 2,000-6,000 m) ocean from simulation EQ-0015.Atlantic heat content from an equilibrium simulation initialized at 1750 CE (EQ-1750, dashed lines) diverges from the EQ-0015 simulation.Heat content anomaly is in units of zettajoules (1 ZJ 10 21 J).
As these cold, deep waters upwell to the surface in EQ-1750, a larger heat flux is required to produce the observed surface temperature warming, in accordance with the heat content evolution.

Inferred Heat Fluxes
The anomalous global heat flux is calculated by distributing the rate of change in ocean heat content over Earth's surface area.We take the entire planetary area, rather than the oceanic area, in order to facilitate comparison with other contributions to the planetary energy balance.The modern warming period is marked by significant decadal variability, although there is clearly a consistent net heat flux into the ocean.The implied heat flux is between 0.2 W m -2 and 0.8 W m -2 , with an average of 0.47 W m -2 between 1965 and 2015 (top panel, Figure 6).The simulated global heat flux is variable but generally consistent with recent estimates (e.g., Lyman and Johnson, 2014).
The global heat flux evolution makes clear that there was no time before the modern warming period when the ocean could be considered in equilibrium.The cooling period of the Little Ice Age is described as a time when the number of years with heat flux out of the ocean exceeded those years when the ocean absorbed heat.Whether the decadal variability of global heat flux during these years reflects actual climatic variation or noise in paleoceanographic records deserves further study.
To what extent does the North Atlantic Ocean contribute to the global energy imbalance?To compute the impact of the North Atlantic, we first diagnose the heat content in all subsurface waters that originated from either the LAB or GIN surface regions, here denoted the NATL region.The source of water is pre-computed using the assumption of an equilibrium circulation (Gebbie and Huybers, 2011), with the result being a water-mass distribution for the world ocean.Heat content can increase due to air-sea fluxes through the surface of the NATL region, or through advective-diffusive fluxes through the subsurface water-mass boundaries.This calculation is similar to a water-mass transformation budget (Walin, 1982), but we use a water-mass distribution rather than an isopycnal layer to define the control volume.Like the global calculation above, the rate of heat content change is translated to heat flux by dividing by an appropriate area.Here, we use the area of the NATL surface region as we suggest that air-sea fluxes will dominate the water-mass transformation.
During Atlantic warming after the Little Ice Age, 1−5 W m -2 of anomalous heat flux must enter the North Atlantic to explain the simulated increase in heat content (bottom panel, Figure 6).The implied heat flux reverses from almost 1 W m -2 out of the ocean in 1825 to 5 W m -2 into the ocean in 1995.These values are much greater than those found in the global calculation because these fluxes must replace and counteract the heat transported to depth by deep water formation.While the decadal-average North Atlantic fluxes are large relative  to other oceanic regions, these fluxes are small compared to those that are extracted from the ocean during deep water formation (e.g., Våge et al., 2008).Therein lies one major challenge in this problem: deep ocean variations result from the small residual of episodic high-intensity surface events.

CONCLUSION
The modern warming period is occurring at the end of a long cooling period extending from the Medieval Warm Period through the Little Ice Age.The long memory of the ocean dictates that the modern warming is influenced to some degree by the preceding surface history.Thus, some error is incurred if it is assumed that the ocean was in equilibrium prior to modern warming.The size of the error depends upon the response time of the ocean.Even though the deep Atlantic communicates with the atmosphere more rapidly than does the deep Pacific, as it hosts deep water formation in both its northern and southern extremes, the full response of the deep Atlantic can take several hundred years.Influences on temperature from the ocean's long memory are small, but the expected temperature trends are long-lasting and thus can potentially be detected from data that are far enough apart in time (Gebbie and Huybers, 2019).In addition, these small temperature trends represent large amounts of energy when integrated over the Atlantic basin for these long periods of time.
To quantify modern Atlantic warming, several important issues must be addressed.For example, we must account for the surface temperature of the ocean during times that pre-date our instrumental record as well as potential ocean circulation change during times in which the wind field was uncertain.Here, we proceed by using paleoceanographic observations to fill in gaps in the surface temperature history and to make the plausible assumption that the large-scale circulation was not dramatically altered given the relatively small temperature perturbations over the Common Era.A simulation of the last 2,000 years has basic features that are expected to be robust: (1) the cumulative effect of the Little Ice Age occurring in the deep Atlantic well after the surface climate forcing occurred, and (2) the reduction of heat uptake in the deep Atlantic due to the effect of ocean disequilibrium at the onset of modern warming.

FIGURE 4 .
FIGURE 4. Simulated Atlantic temperature change from the 1870s to 1990s at 2,500 m depth.Temperature differences are expressed in hundredths of a degree Celsius.90°W 60°W 30°W 0°30°E

FIGURE 6 .
FIGURE 6.Heat fluxes from the EQ-0015 simulation.(top) Global air-sea fluxes with a 10-year (dashed) and 50-year (solid) averaging timescale, where the reference area is the global surface area in order to facilitate comparison with planetary energy balance calculations.(bottom) Similar, but the heat flux into all waters originating from the NATL (LAB + GIN) region of the subpolar North Atlantic, where a regional value of the reference area is used.The NATL region includes both the Labrador Sea and GIN Sea subregions.