The authors characterize impacts on heat in the ocean climate system from transient ocean mesoscale eddies. Their tool is a suite of centennial-scale 1990 radiatively forced numerical climate simulations from three GFDL coupled models comprising the Climate Model, version 2.0–Ocean (CM2-O), model suite. CM2-O models differ in their ocean resolution: CM2.6 uses a 0.1° ocean grid, CM2.5 uses an intermediate grid with 0.25° spacing, and CM2-1deg uses a nominal 1.0° grid.
Analysis of the ocean heat budget reveals that mesoscale eddies act to transport heat upward in a manner that partially compensates (or offsets) for the downward heat transport from the time-mean currents. Stronger vertical eddy heat transport in CM2.6 relative to CM2.5 accounts for the significantly smaller temperature drift in CM2.6. The mesoscale eddy parameterization used in CM2-1deg also imparts an upward heat transport, yet it differs systematically from that found in CM2.6. This analysis points to the fundamental role that ocean mesoscale features play in transient ocean heat uptake. In general, the more accurate simulation found in CM2.6 provides an argument for either including a rich representation of the ocean mesoscale in model simulations of the mean and transient climate or for employing parameterizations that faithfully reflect the role of eddies in both lateral and vertical heat transport.
The ocean plays a fundamental role in establishing the planetary heat budget, particularly on decadal and longer time scales in which net heating or cooling of the planet is largely reflected in a net heating or cooling of the depth-integrated ocean. In particular, the ocean has absorbed an estimated 90% of the additional radiative heating because of industrialization (Otto et al. 2013; Church et al. 2013). Transient climate change, generally defined in terms of changes in the global mean surface air temperature, is thus linked fundamentally to how the ocean transports heat, both vertically and laterally. This role for the ocean in the earth’s heat budget is the key reason that global climate models have incorporated an ocean component since the early days of climate modeling (Bryan and Cox 1967; Bryan 1969).
a. Lateral and vertical heat transport in the ocean
Sea surface temperature (SST) establishes the lower boundary condition for most of the atmosphere, with SST in turn impacting the terrestrial climate state. SST evolution is determined by heat exchange across the ocean surface, as well as lateral and vertical heat transport within the ocean. The importance of lateral, more specifically poleward, heat transport has long been highlighted in ocean climate studies (e.g., Bryan 1996; Jayne and Marotzke 2002). The scale and structure of the poleward ocean heat transport is largely tied to that implied by the atmospheric poleward heat transport. Additionally, there is a growing appreciation for the role that both transient eddies and time-mean currents play on the vertical transport of ocean heat (Gregory 2000; Gnanadesikan et al. 2005; Wolfe et al. 2008; Gregory and Tailleux 2011; Delworth et al. 2012; Morrison et al. 2013; Hieronymus and Nycander 2013; Zika et al. 2013).
Wind-driven vertical circulation tends to transport cold water upward and warm water downward (e.g., Wunsch and Ferrari 2004; Gnanadesikan et al. 2005). In particular, in an ocean that is stably stratified according to temperature, upwelling from Ekman suction brings cold water toward the surface, whereas downwelling from Ekman pumping brings warm water to depths. Consequently, when averaging horizontally over a region, the time-mean circulation acts on the time-mean temperature to impart a downward ocean heat transport. Lateral variations in the vertical heat transport steepen isopycnals and thus increase baroclinicity, particularly in the western boundary current regions and the Southern Ocean. Increases in baroclinicity raise the available potential energy, which in turn enhances mesoscale eddy activity. As emphasized in the analysis of eddying simulations by Wolfe et al. (2008), Delworth et al. (2012), Morrison et al. (2013), and Zika et al. (2013), transient mesoscale eddies respond to the downward heat transport from the mean flow by imparting an offsetting (or compensating) upward eddy-induced heat transport in the process of releasing the available potential energy imparted by the mean circulation. Mesoscale eddies therefore act as a regulator or “gatekeeper” for wind-driven ocean heat uptake, so that interior ocean heating by mean flow does not continue unabated. We, in turn, expect that how eddies are either parameterized or represented in a simulation will greatly impact how much heat enters the numerical ocean. This point serves as the theme of our study.
The fundamental role that mesoscale processes play in climate has motivated a suite of model studies characterizing and quantifying their role, with the mesoscale either parameterized, as in Gregory (2000) and Hieronymus and Nycander (2013), or explicitly represented, as in Yokohata et al. (2007), Shaffrey et al. (2009), Roberts et al. (2009), Bryan et al. (2010), McClean et al. (2011), Delworth et al. (2012), Kirtman et al. (2012), Palter et al. (2014), Bryan et al. (2014), and Winton et al. (2014, manuscript submitted to Geophys. Res. Lett.). Our study contributes to this literature by identifying mechanisms for how ocean mesoscale eddies impact the evolution of ocean heat within a hierarchy of climate models, with an emphasis on the vertical heat transport.
b. Aim and content of this paper
The scientific aim of this paper is to characterize how ocean mesoscale eddies impact the evolution of heat within the ocean climate system. In particular, Delworth et al. (2012) hypothesized that transient mesoscale eddy heat transport plays a primary role in determining climate model drift. Here, we explicitly evaluate the ocean eddy heat transport and ocean heat budget. Our analysis of this budget supports their hypothesis and provides details for the mechanism.
We start the paper in section 2 by introducing the Geophysical Fluid Dynamics Laboratory (GFDL) Climate Model, version 2.0–Ocean (CM2-O) climate model suite used for our study. Besides summarizing model details, we introduce elements of the simulated transient eddy activity by comparing the dynamic sea level variance in the CM2-O suite to that from satellite estimates. In section 3, we characterize features of 140-yr simulations using the CM2-O models. Notably, this simulation length is too short to reach thermodynamic equilibrium. However, it is long enough to capture decadal to multidecadal adjustments, with these transient features the focus of our study. We present a detailed ocean heat budget in section 4 and follow with a discussion of vertical heat transport in section 5. The analysis in sections 4 and 5 quantifies the role of mean and transient features in transporting ocean heat both laterally and vertically. We close the main portion of the paper with discussion and conclusions in section 6. Appendix A summarizes details of the ocean model configurations, and appendix B presents the method used to partition ocean heat transport into transient and time-mean contributions.
2. The CM2-O climate model suite
The present-day ocean thermal state results from accumulating the previous millennia of boundary fluxes, most of which occurred in preindustrial times. The ocean is therefore not in thermal equilibrium with fluxes based on present-day radiative conditions. Hence, we expect an ocean model initialized from modern observations [e.g., Conkright et al. (2002) used here] to accumulate heat, although the spatial patterns and rate of heat uptake are not known a priori. Furthermore, there are errors in numerical methods, physical parameterizations, and initialization that modify the evolution of ocean heat during such transient climate simulations, with these errors generally leading to “model drift.” It is difficult to disentangle model drift from disequilibrium, with both leading to biases (differences) relative to present-day ocean measures. Nonetheless, we are aided in understanding the drift versus disequilibrium issue by considering climate model simulations where only the ocean component differs, such as those considered here.
Circulation features occurring at the ocean mesoscale form the focus of this study, with the mesoscale corresponding to currents and transient eddy motions having spatial scales proportional to the first baroclinic Rossby radius (see, e.g., McWilliams 1996; Smith et al. 2000; Hecht and Hasumi 2008). In particular, we examine how the representation and/or parameterization of the ocean mesoscale impacts on heat transport in the global climate system. For this purpose, we analyze centennial-scale (140-yr) simulations forced with “present-day” 1990 atmospheric radiation that is based on a constant globally averaged CO2 mixing ratio of 355 parts per million by volume (ppmv). In this section, we introduce the CM2-O climate model suite and discuss some questions regarding the resolution of ocean mesoscale features.
a. A model suite to characterize impacts from the ocean mesoscale
We consider results from three coupled climate models forming the GFDL CM2-O suite. All CM2-O models use the same atmosphere of roughly 50-km resolution, as well as identical sea ice and land configurations. Each of these components is described in Delworth et al. (2012), as is the initialization procedure [identical for the CM2-O models, as detailed in section 2f of Delworth et al. (2012)]. The sole difference between models in the CM2-O suite concerns how the ocean is represented.
Among the CM2-O model suite, CM2-1deg is configured with the coarsest ocean, using a nominal 1° horizontal grid that is identical to that used in earlier GFDL climate models, such as the CM2.1 of Delworth et al. (2006), CM3 of Griffies et al. (2011), and Earth System Model with MOM, version 4 component (ESM2M) of Dunne et al. (2012). We make use of a mesoscale eddy parameterization in the tracer budgets, closely configured according to Dunne et al. (2012) (see appendix A for summary). The intermediate- and finest resolution ocean models are the CM2.5 and CM2.6 climate models described by Delworth et al. (2012). CM2.5 uses an ocean grid with a nominal 0.25° spacing, whereas CM2.6 uses a nominal 0.10° grid spacing. The horizontal grid in each of the ocean models uses the tripolar configuration from Murray (1996). The meridional equatorial waveguide region is refined to ⅓° in CM2-1deg, but there is no extra refinement for the CM2.5 and CM2.6 grids. The grid cell thicknesses are a function of time, through fluctuations of the free surface, as per the use of the z* generalized vertical coordinate proposed by Stacey et al. (1995) and Adcroft and Campin (2004). The ocean models in the CM2-O suite use the same vertical grid configuration, with 50 cells distributed throughout the ocean column down to 5500 m. When the ocean is at rest, the upper ocean in the CM2-O models has 10-m thick grid spacing, whereas the deepest cell is 210-m thick.
Neither CM2.5 nor CM2.6 have ocean grids sufficiently refined to justify the absence of a mesoscale eddy parameterization. Indeed, as shown in this paper, CM2.5 is far from converged in regards to its heat budget. However, the question of how to parameterize mesoscale eddies in a model that also permits the eddies remains a topic of ongoing research [see Fox-Kemper and Menemenlis (2008) for a review]. Given the uncertainty in how to best formulate a mesoscale closure for these models, we did not include a mesoscale eddy parameterization for the tracer equations in CM2.5 and CM2.6. Instead, we only implemented a biharmonic friction operator (Griffies and Hallberg 2000) in the momentum equation, aimed at dissipating the downscale enstrophy cascade near the grid scale. Additionally, and central to our study, by not using a mesoscale eddy parameterization in the tracer equation, we can focus on comparing and contrasting the mesoscale features that are represented in CM2.5 and CM2.6.
b. Sea surface variability
We show in Fig. 1 the temporal standard deviation of the dynamic sea level (DSL) computed from daily mean DSL values from the CM2-O models and DSL from Archiving, Validation, and Interpretation of Satellite Oceanographic Data (AVISO) over years 1993–2010, each referenced to their respective climatological annual cycle (see appendix B). This figure provides a climatological measure of the transient behavior of the upper ocean and serves to introduce the CM2-O model suite. Large variability occurs in regions of enhanced mesoscale eddy activity, which also correspond to regions where the surface eddy kinetic energy is larger (see Fig. 14 of Delworth et al. 2012). The zonal-mean standard deviation (Fig. 1e) reveals more quantitatively the comparison between the CM2-O models and satellite measurements.
CM2-1deg (Fig. 1b) generally has very small transient variability, with an exception being the deep tropics of the central and west Pacific, some of which may reflect tropical instability wave and/or interannual variability. Additionally, there is high northern latitude shelf variability in all models. Animations indicate that the shelf variability from the Arctic and into the Baltic arises from fast barotropic wave fluctuations, with the barotropic flow dominating the depth-integrated kinetic energy in these regions. This variability is largely outside the latitude range for satellite measurements.
In general, CM2.6 (Fig. 1d) captures more of the broad regions of modest variance in subtropical gyres of the Pacific and Indian Oceans, whereas this variability is muted in CM2.5 (Fig. 1c) and virtually nonexistent in CM2-1deg (Fig. 1b). We also notice that CM2.6 captures some of the variance in the Gulf of Mexico associated with the loop current, whereas this variance is absent in CM2-1deg and CM2.5. Variability in both CM2.5 and CM2.6 is enhanced in the Kuroshio region of the west Pacific (Sasaki et al. 2014). However, the models show a region of large variance upstream of the separation latitude, with this enhanced variance absent from AVISO. A similar “coastal trapping” of variance is seen in the Gulf Stream region of the North Atlantic in CM2.5 and CM2.6, whereas the AVISO analysis shows variance aligned as a jetlike structure extending away from the coast and turning northwest when it reaches the Grand Banks. This comparison suggests that the Gulf Stream becomes unstable farther away from the coast in the satellite measurements than the simulations. These simulation features represent well-known biases (e.g., Barnier et al. 2006; Bryan et al. 2007; Chassignet and Marshall 2008) that contribute to the large biases in the Atlantic sector (see SST biases in Fig. 5 discussed later in section 3c).
The Agulhas region of the South Atlantic and south Indian Ocean has extensive eddy activity associated with the Agulhas Current and the Antarctic Circumpolar Current (ACC). The eddies that enter the Atlantic in CM2.5 penetrate westward into the central South Atlantic, whereas they dissipate closer to South Africa in CM2.6 as well as in AVISO. This feature of the overly coherent eddies in CM2.5 is reflected in other models of similar resolution (Barnier et al. 2006; Biastoch et al. 2008b). There is particular importance ascribed to accurately simulating this eddy variability because of impacts on Atlantic overturning variability associated with Agulhas mediated interbasin exchange (Biastoch et al. 2008a). Another striking feature of the South Atlantic is the C-shaped variability surrounding the Zapiola region of the southwest Atlantic. CM2.6 well captures this variability as compared to AVISO, whereas the feature is somewhat muted in CM2.5, particularly in its zonal extension away from the coast. Note that CM2.6 shows variability larger than AVISO between 40° and 60°S. We are unsure if this difference represents a limitation of the model or the coarseness of the satellite measurements.
c. Resolution compared to Rossby radius
This correspondence to satellite measures seen for CM2.6 in Fig. 1 suggests that the 0.1° class of ocean models need not employ a mesoscale eddy closure, as such models appear to be “eddy resolving” by this measure, at least in the middle latitudes (Smith et al. 2000). Hallberg (2013) further examines this resolution question via the grid spacing needed to represent the first baroclinic Rossby radius with two grid points.1 Referring to Fig. 1 in Hallberg (2013), the required grid spacing extends from roughly 100 km in the deep tropics to finer than 2 km in the high-latitude shelf regions. Correspondingly, for the CM2.5 and CM2.6 grids used here, a 0.25° Mercator spacing resolves the Rossby radius in most regions equatorward of 30°, whereas a 0.10° Mercator spacing is sufficient for most regions equatorward of 50°.2
3. Ocean temperature evolution and bias patterns
In this section, we characterize features of the temperature drift and associated bias patterns as realized in 140-yr simulations with the CM2-O model suite.
a. Time series
We exhibit in Fig. 2a the time series for annual mean global mean SST. Although we expect the volume mean ocean to warm during these transient simulations, the SST is seen here to cool during the initial 10 years in each of the three models. After the first 10 years, the SST trajectories differ. SST in CM2.5 undergoes a cooling phase for roughly 40 years, after which it displays a gradual rise for the remainder of the 140-yr simulation. SST in CM2-1deg continues to slowly cool until around 60–70 years, after which it reaches a rough steady state that is ≈0.5°C cooler than the initial state. SST in the CM2.6 simulation halts cooling after the first decade and then warms slightly. Starting around year 40, SST in CM2.6 fluctuates slightly below the initial SST value for the remainder of the 140-yr simulation.
Our expectation of ocean warming is reflected in the time series of the global mean heat flux crossing the ocean surface boundary (Fig. 2b), with a positive flux indicating that heat is accumulating in the ocean. Estimates for the net heat uptake were presented by Otto et al. (2013). Observed ocean heat uptake during the decade of the 2000s is estimated as 0.9 ± 0.4 W m−2, as normalized by the ocean area.3 This ocean heat uptake estimate compares to the roughly 1.0 W m−2 in CM2-1deg and CM2.6, whereas the heat uptake in CM2.5 is 1.6 W m−2 (Fig. 2b). CM2.5 experiences a correspondingly higher rate of rising global volume mean potential temperature and steric sea level (Figs. 2c,d). Differences in heat uptake are solely due to different representations and/or parameterizations of ocean mesoscale processes.
The climatological annual cycle for SST based on years 101–140 is shown in Fig. 3. The annual cycle is an important reflection of the interactions between the ocean and atmosphere over the course of a year, and it represents the primary coupled air–sea mode [see Bates et al. (2012) for a recent analysis of the annual cycle in a climate model]. It is therefore of interest to diagnose this cycle after a century spinup to see how the model drift may have impacted upper-ocean seasonal climate. For the Northern Hemisphere (Fig. 3a), CM2.6 is consistently closer to the observation-based estimates. However, all three models are cooler than observations. For the Southern Hemisphere (Fig. 3b), which model is closer to observations depends on the month. Finally, when averaged over the World Ocean (Fig. 3c), CM2.6 agrees more to observations during most of the year, though it is too warm during the Southern Hemisphere summer.
b. Temperature drift as a function of depth and basin
We link the SST and volume mean time series by exhibiting the horizontal mean drift in ocean potential temperature as a function of depth (Fig. 4) and ocean basin. For all basins, the models typically exhibit an upper-ocean cooling accompanied by a deeper warming. This result indicates that heat is transported from the ocean surface vertically downward into the interior of all basins. For CM2-1deg, there is a sizable compensation between upper-ocean cooling and deeper-ocean warming, especially in the Atlantic and Pacific (Figs. 4g,m), that accounts for the modest global volume mean ocean warming (see Fig. 2c). The deep warming in the Southern Ocean (south of the latitude of South Africa) reflects the deep penetration of heat in this region of strong deep currents and generally weak vertical stratification (see Figs. 4p–r).
A striking feature in Fig. 4 is the generally smaller drift in CM2.6 relative to CM2.5 and CM2-1deg, with the Southern Ocean (final row of Fig. 4) a modest exception. Mean downward flow tends to be warmer than the mean upward flow, so that the mean circulation renders a downward heat transport when area averaging. In contrast, eddy heat transport is generally upward and so counteracts the mean (we explore this point in sections 4 and 5). Enhanced ocean heat uptake may therefore be caused either by an overly dominant mean flow or too weak mesoscale eddy field. Delworth et al. (2012) hypothesized that the reduced drift in CM2.6 is due to the improved representation of the ocean mesoscale. We provide support for this hypothesis in section 4.
Delworth et al. (2012) also noted that mesoscale eddy parameterizations impact the interior ocean warming. The global volume mean warming found in CM2-1deg is indeed smaller than CM2.5 (Fig. 2c), suggesting that the eddy parameterization is acting to ameliorate the drift in a manner reflective of the eddies represented in CM2.6. We support this perspective in section 4b, where we explore the physical mechanisms for eddy heat transport. However, the eddy parameterization implemented in CM2-1deg is unable to reflect the small drift over basin scales found in CM2.6. Indeed, the drift is largest for the CM2-1deg simulation in some basins, such as the Atlantic (cf. Figs. 4g–i).
c. Surface and interior temperature bias maps
In Fig. 5 we provide geographical context for the time series in Figs. 2–4. Here, we exhibit maps of the SST and 730-m potential temperature biases averaged over years 101–140 relative to an observation-based analysis; 730 m is the model depth roughly sitting in the middle of the maximum warming found in the global ocean (see Figs. 4a–c). CM2-1deg settles into a cool SST state in the Northern Hemisphere, with widespread sea ice extent (not shown). Furthermore, SST in the subtropical regions of both hemispheres is generally cool. In contrast, CM2-1deg shows widespread warming at 730 m, particularly in the subtropical gyres of the Atlantic and Pacific basins. That is, where the upper portion of the gyres is cool, the deeper portion is warm. This pattern again reflects the transport of heat from the upper ocean and its deposition into the interior.
CM2.5 also shows relatively cool SST in the subtropical gyres, though with a generally warmer Northern Hemisphere than CM2-1deg. The Gulf Stream extension of the North Atlantic also shows a strong cold spot characteristic of biases simulating the Gulf Stream path (e.g., Bryan et al. 2007). The work of Zhang et al. (2011) suggests that much of this bias is related to the relatively poor representation of the overflows in the models. Namely, the poor overflows are associated with weak deep stratification in the Labrador Sea (see also Danabasoglu et al. 2010) and a weak deep western boundary current, each of which impact the North Atlantic current system. In general, these drifts have the potential to adversely impact North Atlantic climate related to air–sea interactions and atmospheric blocking (Weese and Bryan 2006; Joyce et al. 2009; Minobe et al. 2008; Scaife et al. 2011). Upon moving into the ocean interior, CM2.5 shows a warm drift analogous to CM2-1deg. The CM2.5 interior drift is more widespread than CM2-1deg, yet with a notably smaller drift in the Atlantic.
CM2.6 generally shows the least bias among the three models. The cool SST drift in the subtropical gyres found in CM2-1deg and CM2.5 is nearly absent. There is, however, a notable warm drift in the Southern Ocean SST. This drift, seen in both CM2.5 and CM2.6, is associated with a polynya-like event in the Weddell Sea initiated around year 125, in which there is a large release of interior heat that causes a transient surface warming that is advected downstream in the ACC. Such Weddell Sea polynya-like events have been seen in other coarse-resolution climate model simulations, such as those discussed recently by Martin et al. (2013) and de Lavergne et al. (2014).
The CM2-O models exhibit a modest warm SST drift off the west coast of North and South America. As noted by Delworth et al. (2012), the magnitude of this bias pattern is reduced relative to other GFDL climate model simulations [e.g., see Fig. 2 in the CM2.1 simulation from Delworth et al. (2006), as well as Figs. 2 and 5 of Wittenberg et al. (2006)]. The reduced bias is associated with the use of a relatively fine (50 km) atmospheric grid, thus allowing the CM2-O models to more accurately represent the coastal winds off South America giving rise to coastal upwelling [e.g., see Fig. 3 in Gent et al. (2010) or Fig. 6 in Griffies et al. (2011)]. However, unlike the Atlantic improvements found by Gent et al. (2010), who also used a 50-km atmospheric model, all three of the CM2-O models show a warm SST drift of roughly 4°C off the west coast of Southern Africa, which are comparable to the GFDL models with coarser atmospheric components (e.g., see Fig. 2 in Delworth et al. 2006).
d. Zonal-mean bias patterns
In Fig. 6, we show the zonal-mean deviation of temperature averaged over years 101–140 relative to the initial decade (years 1–10) in the World Ocean, Indian–Pacific, and Atlantic–Arctic. In general, as for the other patterns of temperature drift, CM2.6 shows a muted bias for each basin relative to CM2.5 and CM2-1deg.
As seen in Fig. 4, warming in the various basins predominantly occurs within the middepth regions beneath roughly 100–200 m down to 1000–1500 m, whereas there is generally a smaller-amplitude cooling in the near-surface ocean. The Atlantic patterns are rather complex, with CM2-1deg exhibiting a broad warming reaching down to around 1500 m under a relatively cold upper ocean. CM2.5 and CM2.6 have a generally smaller-amplitude Atlantic drift than CM2-1deg. Furthermore, CM2.5 and CM2.6 show cooling that reaches from the mode water regions near the surface around 40°N and extends southward as it descends toward 1000 m. In the high-latitude North Atlantic and Arctic region, there is some warming drift reaching below 1500 m, with CM2-1deg showing the most warming in this region. Near Antarctica, CM2.6 shows a slight warming drift that penetrates to the bottom and toward 40°S, with this warming less coherent in CM2.5 and CM2-1deg.
4. Ocean heat budget and the role of transient eddies
To help understand mechanisms for the different temperature drifts exhibited by the CM2-O models, we make use of an online diagnosis of the ocean heat budget. In particular, we decompose heat advection into transient eddy and mean flow contributions (see appendix B).
a. Depth-integrated ocean heat budget
The depth-integrated ocean heat budget is given by4
This expression states that the time tendency of heat vertically integrated over an ocean column (left-hand side) equals the convergence of heat transport from lateral advective fluxes (ρ dzuΘ) and subgrid-scale processes (ρ dzF), plus contributions from surface boundary heat fluxes.5 Transport processes redistribute heat crossing the ocean surface; they do not modify the net heat in the global ocean. In a steady state [zero left-hand side of Eq. (1)], boundary heat fluxes locally balance the horizontal divergence of depth-integrated heat transport. The CM2-O simulations are not in steady state, with mechanisms for the drift the focus of this section.
b. Poleward heat transport
In Fig. 7, we display the poleward ocean heat transport from the CM2-O model suite for the global ocean, Atlantic–Arctic region, and Pacific–Indian region. For the Pacific–Indian region (Fig. 7c), the CM2-O simulations generally have less southward export of heat from the tropics in the Southern Hemisphere as compared to the reanalysis products, and they are at the lower end of the in situ measurements. Export of heat from the tropics to the north is at the upper margins of the Ganachaud and Wunsch (2003) observation-based analysis.
The simulations generally underestimate the Atlantic heat transport relative to in situ measurements (Fig. 7b). In particular, the cross-equatorial Atlantic heat transport in all three models (about 0.6 PW) is much lower than the in situ measurements (about 1 PW), and this asymmetric meridional heat transport in the tropical Atlantic is mainly due to the asymmetric Atlantic meridional overturning circulation (AMOC). The weaker Atlantic equatorial heat transport extends to the latitude of the Rapid Climate Change (RAPID) measurements at 26.5°N (Johns et al. 2011). In general, the weak heat transport is consistent with the weaker and shallower AMOC in all three models (not shown), with the slightly larger value in CM2-1deg probably due to its slightly stronger AMOC.
In the middle- and high-latitude North Atlantic, there is an agreement of the ocean poleward heat transport in CM2.5 and CM2.6 with Ganachaud and Wunsch (2003) at 45°N, whereas CM2-1deg is consistent with the reanalysis products at this latitude and northward. Specifically, the heat transport at 45°N from Ganachaud and Wunsch (2003) is 0.6 PW, which compares well to the 0.57 PW in both CM2.5 and CM2.6, whereas CM2-1deg and the reanalysis products show 0.45 PW. As noted in sections 14 and D2 of Griffies et al. (2009), there are many limitations of using implied transports from reanalysis products to estimate ocean heat transports, thus motivating us to place more emphasis on the in situ measurements in those regions where the two estimates differ. We also give more credence to the larger transports in CM2.5 and CM2.6 given that CM2-1deg has drifted into a relatively cold North Atlantic surface climate that is outside the range of observational estimates. Decomposing the heat transport into overturning and gyre components (not shown) reveals that the enhanced poleward heat transport in CM2.5 and CM2.6 north of 40°N is largely due to a larger gyre transport rather than overturning.
c. Eddy-mean decomposition of depth-integrated heat advection
We here characterize the mean and eddy advective contributions to the depth-integrated heating in CM2.5 and CM2.6. We also consider how the mesoscale eddy parameterization in CM2-1deg transports heat. As discussed in appendix A, the CM2-1deg mesoscale eddy parameterization contains a contribution from neutral diffusion (Redi 1982; Griffies et al. 1998) and eddy-induced advection (Gent and McWilliams 1990; Gent et al. 1995) [implemented in CM2-1deg as a skew diffusion according to Griffies (1998) and Ferrari et al. (2010)]. Neutral diffusion and eddy-induced advection in CM2-1deg aim to parameterize missing transient mesoscale eddy effects. Both neutral diffusion and eddy advection are explicitly represented, to differing degrees, by the mesoscale eddy contribution to heat advection in CM2.5 and CM2.6 [see Lee et al. (2007)].
Our analysis focuses on the 20-yr period 101–120, with similar results holding for other periods. In Fig. 8, we present the convergence of depth-integrated advective heat transport from CM2.5 and CM2.6 (Figs. 8a,b), and for CM2-1deg (Fig. 8c) we also include the contribution from the mesoscale eddy parameterization (neutral diffusion plus skew diffusion). We decompose the time-mean advection operator in all three models into a contribution from advection due to the time-mean fields (Figs. 8d–f) and advection due to eddy correlations (Figs. 8g–i). Details for the decomposition are provided in appendix B. This analysis is conducted on constant depth surfaces, so that adiabatic processes appear as a compensation between mean and transient. This comment applies in particular to CM2-1deg, where sizable transient eddy contributions appear in the tropics because of both interannual variability and tropical instability waves, with a modest level of mixing associated with wave breaking. Further transient eddy contributions in CM2-1deg appear in western boundaries and the Southern Ocean.
There is a striking resemblance across the three models for contributions from the mean advection, with the exception of the North Atlantic, where CM2-1deg shows distinct behavior associated with its relatively cold surface climate state (see SST bias map in Fig. 5a). Likewise, the sum of the eddy parameterization plus transient eddy term in CM2-1deg shows a close correspondence to patterns of the transient eddies in CM2.5 and CM2.6. We see contributions from eddy transport especially in the Southern Ocean, boundary currents of all basins, and throughout the North Atlantic and North Pacific. The amplitude of the eddy transport generally increases when refining the grid spacing. Also, the eddy variability extends further away from boundary regions as the resolution is refined.
In many regions, there is a compensation between eddy- and time-mean contributions, whereby the eddy and mean heat transports are nearly equal and oppositely directed, thus leading to a relatively small net transport. This compensation is clearly reflected in the poleward heat transport shown in Fig. 9 for the World Ocean, Atlantic–Arctic region, and Indian–Pacific region. Mesoscale eddy compensation of the time-mean heat transport was discussed in idealized studies by Cox (1985), Bryan (1991), Böning and Budich (1992), Drijfhout (1994), Bryan (1996), and Jayne and Marotzke (2002), where they noted that the degree to which eddy heat transport compensates mean heat transport is a function of surface diabatic forcing and interior diapycnal mixing.
The equatorial region in Figs. 8g–i exhibits broad zonally oriented eddy contributions associated with tropical instability waves (TIWs) (Legeckis 1977; Philander et al. 1985), as well as interannual variability. Note the slight north equatorial preference for heating from TIWs arises from their dynamical association with lateral shear instabilities between the South Equatorial Current and the North Equatorial Countercurrent (between 2° and 7°N). Heating from the TIWs mutes the cooling from the Ekman divergence acting on the mean (Figs. 8d–f). Consequently, the equatorward heat transport from TIWs partially compensates the larger poleward transport of heat by the mean flow, as revealed in the poleward heat transport in Figs. 9g–i.
The eddy parameterization plus transient eddies in CM2-1deg, and the transient eddies in CM2.5 and CM2.6, provide significant depth-integrated heating/cooling in the western boundary current regions in the North Atlantic, as well as in the Southern Ocean (Figs. 8g–i). The transient eddies in CM2.5 and CM2.6 also play a role in interbasin transport, particularly between the Indian and Atlantic basins around South Africa. In regions where the convergence patterns are zonally aligned, such as the Southern Ocean, cooling occurs predominantly on the equatorial side of the pattern and warming to the poleward side. This feature reflects a property of the eddy parameterization to transport heat poleward, and it is also represented by the eddies in CM2.5 and CM2.6. A summary of this effect of eddy heat transport, particularly in the Southern Ocean, is seen in the poleward heat transport of Figs. 9a–c. The central and equatorial Pacific in CM2-1deg are notable regions with a negligible contribution from the eddy parameterization, but where transient eddies associated with interannual variability and tropical instability waves are explicitly represented by the grid (roughly 30 km meridionally; e.g., Roberts et al. 2004, 2009).
For the region around 40°N in the North Atlantic, the eddy parameterization convergence patterns in CM2-1deg are largely aligned meridionally (Fig. 8g), so that contributions to poleward transport are small. Nonetheless, it is notable that the eddy parameterization is very active in the North Atlantic, as expected given the baroclinicity associated with the Gulf Stream and North Atlantic current. For CM2.5 and CM2.6, the Gulf Stream is more zonally aligned (Figs. 8h,i), so these models express far more poleward eddy heat transport around 40°N than the eddy parameterization in CM2-1deg (see Figs. 9d–f).
d. Global horizontally integrated heat budget
Differing ocean temperature drifts for the CM2-O model suite are well illustrated by the time–depth drift shown in Fig. 4. This drift is directly related to the evolution of horizontally integrated ocean heat, determined according to the budget
In this equation, is a horizontal integral over the global ocean, with dA the area of a grid cell; −Δk[ρwΘ + ρFz] is the thickness weighted convergence of vertical heat transport from advection plus subgrid scale processes; and δk,1 Qsurf is the surface heat flux directly impacting just the k = 1 grid cell. As noted in section 3b, the surface boundary heat fluxes Qsurf are impacted by the rate that ocean processes move heat from the surface to the interior. Notably, all lateral transport terms are removed from the budget (2) because of the periodic and/or solid wall boundary conditions. Hence, the horizontally integrated heat budget is determined solely by vertical transport processes and boundary fluxes.
The vertical heat flux ρFz is associated with vertical subgrid-scale parameterizations, penetrative shortwave radiation, and nonlocal redistribution of surface boundary heat fluxes (generally cooling fluxes) according to the K-profile parameterization (KPP) scheme of Large et al. (1994). Convective gravitational instability is parameterized by enhancing the vertical diffusivity (Klinger et al. 1996). Convection leads to an upward heat flux, which is dominated in our budgets by the downward heat flux from vertical diffusion in stably stratified regions [see Fig. 2 in Hieronymus and Nycander (2013) for a similar result]. The shortwave term is generally smaller in our budgets than the KPP redistribution of surface boundary fluxes, although more detailed analyses than considered here have identified the importance of shortwave radiation for the upper-ocean heat budget (see Sweeney et al. 2005; Iudicone et al. 2008; Hieronymus and Nycander 2013). The vertical component of parameterized mixed-layer submesoscale eddy processes (Fox-Kemper et al. 2011) contributes a restratification term within the mixed layer. Additionally, CM2-1deg contains a contribution from parameterized mesoscale eddy fluxes via neutral diffusion and skew diffusion as summarized in appendix A.
We decompose the climatological time-mean advective heat flux in CM2.5 and CM2.6 according to appendix B and ask whether and how vertical eddy and time-mean heat flux convergences compare. We perform the same decomposition in CM2-1deg, with the transient eddy contribution smaller than the eddy parameterization in all regions except the tropics. These terms, and those from other processes, are shown in Fig. 10, which displays the time-averaged horizontally integrated heat budget given by Eq. (2) for the CM2-O suite. In the following, we detail the various contributions to this budget, with an aim to identify the dominant terms contributing to the differing drifts shown in Fig. 4.
1) Heat budget above 200 m
We start by considering the heat budget over the upper 200 m, which is dominated by vertical diffusion and the redistribution of surface boundary fluxes from the nonlocal KPP term. In a column that is stably stratified in temperature, vertical diffusion cools the surface grid cell as it weakens vertical stratification. The KPP nonlocal term acts only under negative (i.e., destabilizing) surface buoyancy forcing, which occurs predominantly under cooling surface fluxes. It acts to transport a portion of the negative surface boundary heat flux into the deeper cells within the ocean boundary layer, which in turn warms the surface grid cell and cools the deeper boundary-layer cells.
The parameterized submesoscale eddy restratification acts solely in the mixed layer, where it warms the upper portion and cools the lower, thus increasing mixed-layer stratification. The submesoscale parameterization becomes systematically smaller in magnitude as the model resolution is refined (cf. contributions in Figs. 10a–c), with this behavior arising from the grid-scale dependence built into the scheme (see section 2.1 of Fox-Kemper et al. 2011).
Vertical advection due to the mean generally acts to cool the upper 100 m, moving heat from the upper ocean into the interior below. This interior warming by the mean advection reflects the generally warmer waters that move downward and the generally cooler waters that move upward.
Effects in the upper ocean from both parameterized and represented mesoscale eddies are generally smaller than other processes. For the eddies in CM2.5 and CM2.6 (Figs. 10b,c), they provide a warming in the very near surface cells, slight cooling below, then warming until near 200 m. The parameterized eddy term in CM2-1deg is very small for the surface cells, but then generally warms the upper 200 m (Fig. 10a). One expects the mesoscale eddies and parameterized eddies to warm in the upper ocean and cool below, in order to release available potential energy and in turn to increase stratification. In general, the warming from transient eddies in CM2.6 is larger than CM2.5, as well as the transient plus parameterized mesoscale eddy effects in CM2-1deg. This characteristic helps to greatly reduce the cold near-surface drift in CM2.6 relative to CM2.5 and CM2-1deg.
2) Heat budget below 200 m
Below 200 m, the magnitude of the heat flux convergences are much smaller than in the upper ocean, largely because of the absence of strong boundary-layer mixing processes. Recall from Fig. 4 that it is in the region between 200 and 1600 m that differences in temperature drift among the CM2-O model suite are most visible. Figure 10d in particular shows that the net heating (drift) for CM2-1deg and CM2.5 is indeed larger than in CM2.6, with CM2.5 showing nearly twice the heating. We now detail how the budget terms in Figs. 10a–c account for these differences.
The KPP nonlocal redistribution and submesoscale restratification both cool the interior. Impacts from this cooling get systematically smaller and shallower as the ocean model grid is refined. Warming from vertical diffusion in CM2.5 contributes the smallest level of heating among the CM2-O suite, perhaps due to this model losing the most vertical stratification over the global ocean (see Fig. 4).
Mesoscale eddies (CM2.5 and CM2.6) or parameterized eddies (CM2-1deg) generally cool the ocean interior, whereas mean advection warms. Warming in CM2.6 from mean advection is the largest among the CM2-O suite, perhaps due to this model maintaining a stronger mean vertical stratification on which mean advection acts. The magnitude of interior cooling from the advection of transient eddies in CM2.6 increases when moving upward from the interior into the region shallower than 200 m. In contrast, eddy-induced cooling is relatively constant in CM2.5 between 600 and 300 m and reaches a maximum cooling around 300 m. Interior cooling from the mesoscale eddy parameterization in CM2-1deg reaches a maximum at the deeper level of roughly 600 m and switches to a warming near 300 m and shallower.
We conclude that CM2.5 exhibits a larger rate of warming below 200 m than CM2.6 because of the differing contributions from transient eddies. Namely, the relatively strong cooling by mesoscale eddies in CM2.6 is far more effective at counteracting the warming by mean advection. For example, at 500-m depth, which is near the maximum of the warming drift for the global ocean in Fig. 4, the cooling tendency of −0.057°C yr−1 from eddies in CM2.6 nearly balances the warming from the mean, whereas in CM2.5, the eddy cooling is roughly 25% weaker than the corresponding mean warming.
In general, there is a very close balance between the advection of the mean and eddies in CM2.6 that extends from 200 m into the deep ocean (Fig. 10c), implying that the net advective transport (sum of mean and eddies) is close to adiabatic (i.e., zero net heating). In contrast, heat advection in CM2.5 (Fig. 10b) has a notable imbalance, particularly for regions shallower than 800 m. We further reemphasize that warming in this depth range from vertical diffusion and advection acting on the time-mean fields is largest in CM2.6, whereas the magnitude of cooling from KPP and submesoscale is smallest. Together, these effects suggest that CM2.6 should warm the most. However, cooling by advection of transient eddies is sufficient to bring the net heating down in CM2.6 to the smallest in the CM2-O suite, in particular to a rate a bit more than half that of CM2.5 (Fig. 10d). The comparison to CM2-1deg is largely analogous to CM2.5, though the magnitude of heating/cooling from both the mean and eddy parameterization in CM2-1deg are somewhat less than the eddying simulations.
This analysis provides direct evidence that the enhanced cooling effect from transient eddies in CM2.6 is the reason this model experiences a smaller warming drift in the region below 200 m. Stronger eddies more ably counteract the heating from the mean flow, which is precisely the hypothesis proposed by Delworth et al. (2012). In so doing, eddies regulate how much heat is pumped into the interior ocean and retain more heat in the upper ocean. This process in turn feeds back to the atmosphere by reducing the heat entering the ocean through the air–sea interface.
e. Regional heat budgets
Heat budgets for the Pacific, Atlantic, and Southern Oceans are shown in Fig. 11, with results plotted to the ocean bottom to illustrate the different depth penetration of the mean and eddies across the ocean basins. Note that lateral transport terms contribute to the regional budgets along the region boundaries. The global horizontal averaged budget largely reflects the behavior in the different basins, with no single basin as an outlier, at least insofar as the contributions from mean and eddy advection. However, impacts from transient eddies reach only into the upper few hundred meters in the Indian (not shown) and Pacific (Figs. 11a–c), slightly deeper in the Arctic (not shown), much deeper in the Atlantic (Figs. 11d–f), and to the bottom in the Southern Ocean (Figs. 11g–i).
Note the deep penetration of mean flow in the Southern Ocean (Figs. 11g–i), associated with deep penetration of the wind-driven circulation (e.g., Rintoul and Naveira-Garabato 2013). The temperature tendency below roughly 800 m is driven almost entirely by a large warming through mean advection, with the cooling effect from mesoscale eddies only partially offsetting the warming. Below the sills around 4000 m, the eddy and mean contribution to ocean heating switch signs. It is for these depths that Antarctic Bottom Water (AABW) flows equatorward, acting to cool the deep Southern Ocean. Wolfe et al. (2008) also noted this characteristic in their eddying simulations, in which their mean and eddy vertical heat fluxes changed sign below roughly 3500 m. A similar sign switch in tendencies is seen in the Atlantic basin (Figs. 11d–f), presumably because of the impacts from AABWs cooling the deep Atlantic basin.
The Atlantic basin (Figs. 11d–f) shows some striking distinctions for the advective and mesoscale eddy parameterization contributions in the region below 200 m. For both CM2-1deg and CM2.6, mean advection exhibits a strong warming, which is somewhat weaker in CM2.5. For both CM2.5 and CM2.6, the eddy cooling largely compensates the mean warming below roughly 600 m, whereas the mean warming in CM2-1deg dominates the transient eddies plus mesoscale parameterization until roughly 1200 m. This feature of the parameterization in CM2-1deg largely explains the large warming seen in the time–depth drifts in Fig. 4, as well as the zonal-mean drifts in Fig. 6.
5. Vertical heat transport
We now consider the vertical ocean heat transport in the simulations and provide a conceptual picture for why we should expect the transport from mesoscale eddies to be upward.
a. Vertical heat transport in the simulations
Vertical heat transport is computed by performing a vertical running sum of the heat budget (2) from the bottom upward.6 Results shown in Fig. 12 provide a vertical analog to the poleward heat transports shown in Figs. 7 and 9.
For all three of the CM2-O models, vertical heat transport in Fig. 12 from the submesoscale scheme is positive, indicating that the scheme moves heat upward. It has a maximum in the near-surface ocean, with cooling below the maximum and warming above, again reflecting the restratifying effects from this scheme. The vertical heat transport from the KPP nonlocal redistribution is also positive, which is associated with the transport of negative heat flux downward throughout the boundary layer. In contrast, vertical diffusion transports heat downward into the ocean interior. Likewise, there is a downward heat transport from mean advection, with this transport largest in CM2.6.
A maximum upward heat transport from the mesoscale eddy parameterization plus transient eddies in CM2-1deg is ≈ 1.3 PW, with this maximum occurring around 300 m (Fig. 12a). The transient eddy heat transport in CM2.5 is ≈ 1.3 PW near 200 m (Fig. 12b). The heat transport is larger in CM2.6, with a peak of ≈ 2.1 PW (50% larger than CM2.5) occurring slightly shallower than 200 m (Fig. 12c). The extra 0.7 PW of upward heat transport in CM2.6 is entirely due to the stronger eddy heat fluxes. The eddy parameterization in CM2-1deg produces a parameterized vertical heat transport that is both significantly weaker than CM2.6 and has a maximum that is more than twice as deep. The net downward heat transport associated with the mean plus eddy advection (plus eddy parameterization in CM2-1deg) occurs around 100 m in all the models and has values of ≈ 2.3 PW in CM2-1deg, whereas CM2.5 has ≈ 2.1 PW and CM2.6 has ≈ 1.4 PW. These differences largely arise from variations in the eddy parameterization (in CM2-1deg) or transient eddy contributions.
b. Why we expect mesoscale eddy heat fluxes to be upward
Consider regions where potential density stratification is dominated by temperature (e.g., middle and low latitudes). An eddy-mean flow decomposition of the temperature budget indicates that in the steady state, correlations between vertical velocity and temperature fluctuations lead to an upward eddy heat flux in order to balance small-scale dissipation. The associated slumping of isopycnals/isotherms by the eddy fluxes reflects the adiabatic extraction of available potential energy by baroclinic eddies, in which heavy (cold) water is rearranged to sit under lighter (warm) water. This vertically upward heat transport is a basic property of baroclinic eddies that is built into the Gent et al. (1995) parameterizations of mesoscale eddies (e.g., see Fig. 1 of Griffies 1998). In particular, parameterized vertical heat transport is enhanced in strongly baroclinic regions where eddies are vigorous, such as western boundary currents and the Southern Ocean.
Neutral directions are generally misaligned with conservative/potential temperature isolines. In such regions, there is a component of the eddy heat flux that is generally parameterized as downgradient diffusion oriented along neutral directions (Solomon 1971; Redi 1982; McDougall and Church 1986; Griffies et al. 1998). In a study based on observational data, Osborn (1998) noted the generally upward heat flux associated with neutral diffusion. Earlier, Davis (1994, his section 3c) pointed to the general cooling effects from neutral diffusive fluxes over broad regions of the ocean interior. Our particular region of interest concerns the Southern Ocean, where isotherms are typically sloped more steeply than neutral directions in the upper 1000 m. An example of this orientation is shown in Fig. 13, in which case downgradient neutral diffusion of temperature fluxes heat poleward and upward. This property of Southern Ocean eddy heat fluxes was emphasized by Gregory (2000) and Hieronymus and Nycander (2013) in coarse-resolution climate models, whereas Morrison et al. (2013) discussed it in an eddy-permitting process study. Each of these studies highlighted the primary role played by the Southern Ocean in effecting the global ocean heat budget, with neutral diffusive fluxes the key contributor in that region.
6. Discussion and conclusions
Global climate models are transitioning toward refined ocean grid resolutions sufficient to admit time-mean and transient ocean mesoscale features. During this transition, the community must develop robust and accurate numerical methods, such as suggested by Adcroft and Hallberg (2006) for vertical coordinates, and transport methods that minimize spurious mixing, as emphasized by Ilicak et al. (2012). Additional work is needed for methods that parameterize the unresolved transient mesoscale without overdamping the resolved mesoscale, such as through ideas proposed by Hallberg (2013) and Jansen and Held (2014). As the numerics and parameterizations are improved, so too must the characterization of how the ocean mesoscale impacts the climate system. It is toward the latter goal that the present paper is focused. Here, we identified roles for the mesoscale in establishing the distribution of heat within the ocean and by implication within the climate system.
Before summarizing our analysis that focused on the ocean heat budget and associated temperature distributions, we acknowledge that accurate salinity distributions are also critical for maintaining ocean water mass integrity. Analysis of the salinity biases (not shown) indicates that the finest-resolution ocean in CM2.6 retains the least biases among the CM2-O climate model suite, consistent with its behavior for temperature.
a. The importance of vertical heat transport by transient mesoscale eddies
The structure and scale of poleward ocean heat transport is largely tied to the poleward transport in the atmosphere. We may thus expect similar general features of oceanic poleward heat transport among the CM2-O models, given that they each use the same atmospheric model component. This expectation is supported by the poleward heat transports shown in Figs. 7 and 9, which show only subtle differences between the models, particularly between CM2.5 and CM2.6. However, as shown in this study, differences in oceanic processes have a leading order impact on vertical heat transport.
The predominant impact from vertical circulation acting on the time-mean fields is to transport warm water from the upper ocean into the interior, and cool water from the interior to the upper ocean, thus rendering a net downward heat flux when averaging over a horizontal region. Wunsch and Ferrari (2004) and Gnanadesikan et al. (2005) emphasized the role of winds in imparting the necessary mechanical energy needed to support this heat transport. Correspondingly, as noted by Delworth et al. (2006), biases in the wind patterns can lead to significant biases in ocean heating.
Mechanical forcing from winds enhances ocean baroclinicity, which increases available potential energy (APE). Transient mesoscale eddies, particularly in the western boundaries and Southern Ocean, feed off the APE. As identified by Osborn (1998), Gregory (2000), Wolfe et al. (2008), Delworth et al. (2012), Morrison et al. (2013), Hieronymus and Nycander (2013), and Zika et al. (2013), mesoscale eddies transport heat upward, thus partially compensating (or offsetting) the downward heat transport from the mean. Figure 12 indicates this role for the eddies (and eddy parameterization) in the CM2-O model suite, with differences in the strength and vertical structure of the eddy transport affecting the rates of heat drift. In particular, we found that the interior (200–1500 m) cooling from transient mesoscale eddies is stronger in CM2.6 than CM2.5, thus contributing to a smaller interior warming drift in CM2.6. Conversely, stronger warming in CM2.6 from eddies in the region shallower than 200 m contributed to its relatively small upper-ocean temperature drift. This heat budget analysis thus provides a mechanistic explanation for the very different heat drifts found in CM2.5 versus CM2.6, and in so doing, it supports the eddy hypothesis of Delworth et al. (2012).
Although we make use of a relatively high-order advection scheme in CM2-O, there is no perfect discretization. In particular, dispersion errors are smoothed with an upwind bias that mixes water analogously to diffusion. Ilicak et al. (2012), building on earlier work from Griffies et al. (2000), emphasized the role of the grid Reynolds number and the representation of features near the grid scale. In general, poorly resolved features lead to increases in spurious diapycnal mixing, with such spurious mixing increasing interior advective heating whether impacting on the mean or eddy contributions. Correspondingly, one may conjecture that CM2.6 better resolves the dominant energy scales than CM2.5, thus reducing the role of spurious interior ocean heating in CM2.6. Unfortunately, we do not know how to determine the spatial patterns of spurious mixing in realistic simulations, with the energetic methods of Ilicak et al. (2012) providing only a global diagnostic indicative of spurious mixing levels. It therefore remains a topic of ongoing research to determine accurate spatial information about spurious mixing, thus allowing for investigations into its relevance for ocean tracer budgets.
b. Mesoscale eddy parameterizations
Disentangling the many issues associated with climate model drift and model initialization remains a difficult task for the modeling community. In particular, the results here point to the importance of vertical mesoscale eddy heat transport (either parameterized or represented) in effecting model drift associated with ocean heating in 1990 radiatively forced simulations. We conjecture that these results are relevant for more realistic, historically forced simulations where the radiative forcing changes in time. In general, the more accurate simulation found in CM2.6 provides an argument for either including a rich representation of the ocean mesoscale in model simulations of the mean and transient climate, or for employing parameterizations that faithfully reflect the role of eddies in vertical heat transport.
As noted by Bryan et al. (2014), mesoscale eddy parameterizations based on Redi (1982) and Gent and McWilliams (1990) capture a good part of the impacts from mesoscale eddies through their upward heat transport that compensates, at least partially, for the downward mean advective transport. Nonetheless, our implementation of the eddy parameterization in CM2-1deg failed to fully capture the regional heat uptake found in CM2.6. In particular, it underestimated the strength and structure of the upward heat transport by transient eddies represented in CM2.6.
Part of the shortcoming in CM2-1deg may be an inadequate implementation of its mesoscale eddy parameterization. One issue may be related to the maximum neutral slope used by the neutral diffusion scheme (Griffies et al. 1998), with Gnanadesikan et al. (2007) indicating an impact from this parameter especially on high-latitude ventilation. More testing is needed to determine a suitable value for this parameter. Another issue may be related to the use of a depth-independent eddy diffusivity (Griffies et al. 2005; Gnanadesikan et al. 2006) in CM2-1deg, where recent studies suggest that depth structure is important to better capture impacts from wind perturbations in the Southern Ocean (Gent and Danabasoglu 2011; Hofmann and Morales Maqueda 2011; Farneti and Gent 2011). A depth-dependent eddy diffusivity was also used by Bryan et al. (2014), who focused on climate change in the Southern Ocean, though they did not analyze the heat budget. We conjecture that depth dependence to the eddy diffusivity may offer a better means for CM2-1deg to emulate the vertical structure of ocean heating found in CM2.6. This question is the topic of ongoing research.
It remains a research question how to formulate a suitable mesoscale eddy parameterization that captures effects from unresolved transient eddies in a model that also allows for a portion of the eddy spectrum to be represented. For such models, which characterize CM2.5 and CM2.6 particularly for the higher latitudes (see Fig. 1 from Hallberg 2013), the use of a traditional parameterization often comes at the price of overdamping the finescale currents, also important for the mean climate. Roberts and Marshall (1998) and Smith and Gent (2004) propose closures applied to the tracer equation aimed at overcoming this limitation, whereas Hallberg (2013) proposes a resolution function tailoring application of the traditional closures only to those regions where needed. Guided by the inverse cascade of geostrophic turbulence, Jansen and Held (2014) consider a “negative Laplacian” operator in the momentum equation as a means to reintroduce energy to the large scales that is otherwise dissipated by the biharmonic friction operator. We suggest that research into the fidelity of these, and other, methods should be gauged in part by their impact on vertical heat transport. This is an important practical issue, particularly since many modeling groups for the next few years are unlikely to be able to afford ocean resolutions finer than the ≈ 25 km used in CM2.5.
c. What resolution is sufficient?
The ocean is absorbing heat because of increasing anthropogenic greenhouse gases. A warming ocean reduces the amount of heat in the atmosphere while it increases steric sea level (Church et al. 2013). Correspondingly, a key goal of climate science is to produce models that accurately and faithfully simulate ocean heat transport and its convergence. Unfortunately, inaccurate numerical methods and/or uncertain physical parameterizations generally lead to climate model drift. Model drift, and differences across the suite of climate models, compromise our ability to provide deductive statements regarding physical mechanisms and reduce our ability to produce skillful dynamically based climate predictions and long-term projections. Understanding mechanisms for climate model drift therefore remains one of the key issues for climate model development and analysis. The results of this paper shed some light on this issue by emphasizing roles for the ocean mesoscale in affecting how heat is distributed laterally and vertically in the ocean. By extension, we presume that this understanding also aids in uncovering mechanisms for the wide differences in climate model–projected thermosteric sea level rise (Kuhlbrodt and Gregory 2012).
An assumption underlying our study is that an explicit representation of ocean features, the mesoscale in particular, provides a more accurate simulation than parameterizing these features. This assumption leads us to consider CM2.6 as the most accurate simulation in the CM2-O suite. In turn, we place more trust in the CM2.6 simulation of the heat budget than the eddy parameterized simulation from CM2-1deg. It also allowed us to deduce why CM2.5 underrepresented the impacts of mesoscale eddies on ocean heat transport, in particular the vertical transport. Even so, there is no study that determines what resolution is sufficient to accurately simulate ocean tracer budgets for purposes of global climate modeling. That is, we do not know what resolution is required for the budgets to remain unchanged when resolution is refined. We are therefore unable to address just how accurate CM2.6 is compared to a climate model with an even finer ocean resolution. This question remains a topic for future studies.
We thank the GFDL community for their support of this work, which required a tremendous amount of computational hardware and software resources. We thank Alistair Adcroft, Ivy Frenger, Robert Hallberg, Andy Hogg, Malte Jansen, Rym Msadek, Jorge Sarmiento, Ron Stouffer, and three anonymous reviewers for helpful and encouraging comments on drafts of this paper. Many figures were produced by Ferret developed at NOAA/PMEL. Carolina Dufour was supported by the U.S. Department of Energy under Contract DE-SC0006848. Adele Morrison was supported by the Carbon Mitigation Initiative (CMI) project at Princeton University, sponsored by BP.
Summary of the Ocean Model Formulations
The Modular Ocean Model, version 5 (MOM5), code of Griffies (2012) provides the ocean component of the CM2-O climate model suite. Notably, MOM5 conserves scalar fields (e.g., heat, salt, and biogeochemical tracers) to machine precision, with details presented in sections 12.5 and 12.6 of Griffies (2004). This conservation property facilitates our analysis of the ocean heat budget using online diagnostics in section 4.
Delworth et al. (2012) detail physical and numerical choices for the CM2.5 and CM2.6 ocean configurations, as well as the atmosphere, land, and sea ice components. Ocean lateral tracer transport in CM2.5 and CM2.6 is impacted by resolved advection, which is based on a piecewise parabolic method, and the Fox-Kemper et al. (2011) parameterization of mixed layer submesoscale eddies. CM2.5 and CM2.6 do not use a mesoscale eddy parameterization in their tracer equation. Parameterized vertical mixing processes include the internal gravity wave breaking scheme of Simmons et al. (2004), the coastal mixing scheme of Lee et al. (2006), and the internal shear mixing and KPP surface boundary-layer scheme of Large et al. (1994). Regions of gravitational instability are stabilized by enhancing both vertical diffusivity and viscosity (Klinger et al. 1996). None of the CM2-O ocean configurations make use of an overflow parameterization.
CM2-1deg makes use of the same atmosphere, land, and sea ice configurations as CM2.5 and CM2.6, with the ocean also employing the same vertical grid arrangement and vertical physical parameterizations. For the lateral physical parameterizations, CM2-1deg deviates from CM2.5 and CM2.6 in the following ways.
Lateral friction is based on the Laplacian plus biharmonic form chosen for the 1° ocean component of the earth system model ESM2M detailed by Dunne et al. (2012).
The mesoscale eddy parameterization also follows the ESM2M setup in Dunne et al. (2012). In particular, we employ the skew flux formulation of Ferrari et al. (2010), which is a modified version of the Gent et al. (1995) parameterization. We use here a larger maximum diffusivity of 1200 m2 s−1 rather than 800 m2 s−1 used by Dunne et al. (2012). The larger maximum diffusivity is motivated by studies of Farneti and Gent (2011) and Gent and Danabasoglu (2011), who investigated the Southern Ocean response to increasing winds. For neutral diffusion, we retain the ESM2M choice with a constant neutral diffusivity of 600 m2 s−1.
Marginal seas are connected to the open ocean through the cross-land mixing scheme used in ESM2M, with details given by Griffies et al. (2005).
Table 1 exhibits time steps for the CM2-O ocean models and coupling time steps for interactions with the atmosphere and sea ice. We chose time steps small enough to ensure model stability even in the presence of a relatively energetic 50-km atmosphere model [far more energetic than the 200-km atmosphere used in CM2.1/ESM2M (Delworth et al. 2006) or in CM3 (Donner et al. 2011)], along with diurnal air–sea fluxes standard in GFDL climate models since Delworth et al. (2006). In particular, the synoptic atmospheric storms can sporadically input a tremendous amount of inertial energy to the ocean, much more than that provided by ocean-ice forcing such as that from Large and Yeager (2009).
Vecchi et al. (2014) employ a prediction model similar to CM2-1deg, but with differing ocean vertical grid resolution [corresponding to CM2.1 of Delworth et al. (2006)], differing mesoscale eddy parameterizations, and differing lateral friction. These differences lead to the surface ocean climate of this “CM2.5-FLOR” model of Vecchi et al. (2014) to be warmer, and more realistic, in the North Atlantic and North Pacific than CM2-1deg. However, the interior ocean warm drift in CM2.5-FLOR is far larger than in CM2-1deg. Preliminary studies (not shown) suggest that a key difference arises from the reduced lateral viscosity in CM2.5-FLOR relative to CM2-1deg.
Tracer Advection Decomposed into Eddy and Mean
We introduce a semidiscrete expression for the convergence of advective tracer fluxes,
which imparts a tendency for tracer mass per horizontal area in a discrete model grid cell. In this equation, C is either the tracer concentration, potential temperature θ (as used in the CM2-O suite), or conservative temperature Θ (IOC et al. 2010). The thickness-weighted convergence of the vertical advective flux is
with k a discrete vertical label that increases downward. Likewise, −∇ · [ρ dzuC] is the convergence of the horizontal advective fluxes. The ocean configurations as part of the CM2-O climate models assume a Boussinesq ocean, so that the in situ density ρ reduces in these expressions to the constant reference density ρo. However, the generalized level coordinates used in CM2-O oceans means that the grid cell thickness dz is a function of (x, y, z, t) for all grid cells.
Introduce the thickness weighted horizontal advective mass transport, and the vertical advective mass transport, each crossing the respective faces of a tracer cell:
in which case the advection operator (B1) takes the form
Next, decompose the mass transport and tracer concentration into time-mean and transient eddy components:
where angle brackets denote a climatological time-mean and transient fluctuations are denoted by the prime. Following Bryan et al. (2014), we define a transient fluctuation relative to its climatological monthly value. In this way, time “mean” contributions incorporate a time-mean plus climatological seasonal cycle. Eddies, in turn, represent fluctuations associated with mesoscale eddies, tropical instability waves, and interannual variability. Interannual variability is large in the tropics (e.g., associated with El Niño), whereas transient fluctuations are dominated by mesoscale eddies in the middle and high latitudes.
As so defined, a fluctuating eddy quantity is defined relative to a mean field, so that the mean of the fluctuating quantity vanishes. Hence, decomposing the time-mean advection operator yields
is the advection operator based on fluxes constructed with the climatological time-mean transport and time-mean tracer concentration (the mean advection), whereas
is the advection operator based on the convergence of the climatological transient eddy fluxes (the eddy advection). By focusing on the convergence of advective fluxes, rather than components to the advective fluxes, we can deduce effects on the tracer tendency from eddies without concern for the arbitrary rotational portion of the eddy flux field (Fox-Kemper et al. 2003). This approach was also taken by Bryan et al. (1999).
The Rossby radius R is related to the baroclinic Rossby wavelength λ via R = λ/2π. If the grid spacing Δ satisfies Δ < R, then 2πΔ < λ, which is the traditional criteria for resolving a wave on a discrete grid. Hallberg’s criteria is even more conservative, in which 2Δ < R so that 4πΔ < λ.
The 0.25° Mercator spacing is ~ (28 cosϕ) km, where ϕ is the latitude, whereas 0.10° Mercator spacing is ~ (11 cosϕ) km.
We consider the decade of the 2000s as more representative than the 1990s since our simulations do not consider volcanic forcing, which plays an important role in the observed climate system during the 1990s.
The variable Θ is conservative temperature (IOC et al. 2010), which we approximate as potential temperature θ in our simulations, and ρ the in situ density, which we set to the constant ρo = 1035 kg m−3 according to the Boussinesq approximation. The seawater heat capacity is taken as the constant Cp = 3992 J kg−1 K−1. The vertical sum in Eq. (1) extends over an ocean column with space–time–dependent grid cells of thickness dz.
There are no bottom geothermal heat fluxes in the CM2-O suite of models.
The vertical heat flux vanishes at the bottom because of the absence of geothermal heating in the CM2-O models.