Multidecadal variability of the Atlantic meridional overturning circulation (MOC) is investigated diagnostically in the NCAR Community Climate System Model version 3 (CCSM3) present-day simulations, using the highest (T85 × 1) resolution version. This variability has a 21-yr period and is present in many other ocean fields in the North Atlantic. In MOC, the oscillation amplitude is about 4.5 Sv (1 Sv ≡ 106 m3 s−1), corresponding to 20% of the mean maximum MOC transport. The northward heat transport (NHT) variability has an amplitude of about 0.12 PW, representing 10% of the mean maximum NHT. In sea surface temperature (SST) and sea surface salinity (SSS), the peak-to-peak changes can be as large as 6°–7°C and 3 psu, respectively. The Labrador Sea region is identified as the deep-water formation (DWF) site associated with the MOC oscillations. In contrast with some previous studies, temperature and salinity contributions to the total density in this DWF region are almost equal and in phase. The heat and freshwater budget analyses performed for the DWF site indicate a complex relationship between the DWF, MOC, North Atlantic Oscillation (NAO), and subpolar gyre circulation anomalies. Their complicated interactions appear to be responsible for the maintenance of this multidecadal oscillation. In these interactions, the atmospheric variability associated with the model’s NAO plays a prominent role. In particular, the NAO modulates the subpolar gyre strength and contributes to the formation of the temperature and salinity anomalies that lead to positive/negative density anomalies at the DWF site. In addition, the wind stress curl anomalies occurring during the transition phase between the positive and negative NAO states produce fluctuations of the subtropical–subpolar gyre boundary, thus creating midlatitude SST and SSS anomalies. Comparisons with observations show that neither the pattern nor the magnitude of this dominant SST variability is realistic.
Many coupled general circulation models (CGCMs) used in climate studies exhibit multidecadal oscillations in their meridional overturning circulations (MOCs) in the Atlantic Ocean (e.g., Delworth et al. 1993; Timmermann et al. 1998; Cheng et al. 2004; Dong and Sutton 2005; Dai et al. 2005; Jungclaus et al. 2005). These oscillations are mostly irregular and their periods change considerably among models. For example, while ECHAM5/Max Planck Institute Ocean Model (MPI-OM; Jungclaus et al. 2005) has one of the longest periods with 70–80 yr, the Third Hadley Centre Coupled Ocean–Atmosphere General Circulation Model (HadCM3; Dong and Sutton 2005) and Parallel Climate Model (PCM; Dai et al. 2005) show the shortest periods with about 25 yr. The middle range includes the ECHAM3/large-scale geostrophic ocean GCM (LSG; Timmermann et al. 1998) and an earlier version of the Geophysical Fluid Dynamics Laboratory (GFDL) model (Delworth et al. 1993) with periods of about 35 and 50 yr, respectively. This oscillation is also present in other ocean fields, particularly in the North Atlantic, including the northward heat transport (NHT) and sea surface temperatures (SSTs), which are two of the climatically most important ocean fields.
Recent observational studies based on instrumental and proxy data also show distinct multidecadal variability in SSTs with periods of about 40–70 yr (e.g., Kushnir 1994; Delworth and Mann 2000). The associated spatial pattern is particularly dominant in the North Atlantic, but it has suggestions of a broader hemispheric–global pattern. In the North Atlantic, it is largely basin scale, indicating broad warming and cooling, and its maximum local amplitude is about 0.5°C. This multidecadal variability is sometimes referred to as the Atlantic Multidecadal Oscillation (AMO) and it has been associated with multidecadal variations of the North American and western European summertime climate (Sutton and Hodson 2005).
A broad resemblance between the CGCM simulated and observed multidecadal SST variability patterns in the North Atlantic is shown in Delworth et al. (1993) and Timmermann et al. (1998). This SST variability is usually associated with the MOC variability in CGCM studies. Griffies and Bryan (1997) indicate that the MOC variability may be predictable on either decadal or longer time scales, implying the potential predictability of the associated climate changes in North America and western Europe. Such predictability may also have some implications for hurricane forecasts in the North Atlantic, because there are small SST changes (0.1°–0.2°C) in the tropical Atlantic associated with the MOC oscillations. However, the influence of such multidecadal natural variability on hurricane activity remains controversial (see, e.g., Trenberth and Shea 2006; Vimont and Kossin 2007). Clearly the presence of such a multidecadal intrinsic, that is, unforced, variability also complicates climate studies investigating anthropogenic effects. Therefore, it is important to understand the details of this oscillation.
Similar (midlatitude) multidecadal or longer time-scale MOC oscillations have been identified and analyzed in numerous other studies, using simple, idealized models (e.g., Weaver and Sarachik 1991; Greatbatch and Zhang 1995; Saravanan and Mc Williams 1997; Capotondi and Holland 1997; Neelin and Weng 1999; Colin de Verdiere and Huck 1999; Te Raa and Dijkstra 2002; Dijkstra et al. 2006). Saravanan et al. (2000) briefly summarize various mechanisms that have been proposed to explain this variability. While we acknowledge that such simpler models are useful to test hypotheses and isolate mechanisms, we mostly discuss results from the more complex CGCMs in this study.
Delworth et al. (1993) show that the density anomalies in the sinking region of the overturning circulation drive these multidecadal oscillations. Reduced heat transport associated with a weak MOC leads to a cold, dense pool in the middle of the North Atlantic. This cold pool has an associated cyclonic circulation that transports salt into this sinking region, thus increasing further the density there. As a result, MOC strengthens, leading to the transport of warmer, less dense waters into the sinking area. In turn, MOC weakens again, accompanied by reduced heat transport. In Delworth et al. (1993), the existence of the oscillation crucially depends on the phase lag between the temperature and salinity contributions to the total density in the deep-water formation (DWF) regions. In this early, flux-corrected GFDL model, this variability is interpreted as a damped ocean-only mode excited by atmospheric noise (Delworth et al. 1993; Griffies and Tziperman 1995). Further analyses of the same model by Weaver and Valcke (1998) and Delworth and Greatbatch (2000, hereafter DG00) produce two rather different conclusions. While the former argues that this variability is a coupled mode, the latter reinforces the view of a damped ocean-only mode, continuously excited by low-frequency atmospheric forcing. This major discrepancy in these conclusions can be partly attributed to differences in surface forcings used in the respective sensitivity experiments as well as to the differences in the definition of a coupled mode; our view of a coupled mode follows that of DG00, namely, a mode as represented by the El Niño–Southern Oscillation phenomenon. We note that DG00 also show that surface heat flux variations with patterns resembling the ones associated with the North Atlantic Oscillation (NAO) play a dominant role in driving these oscillations. Both Dai et al. (2005) and Dong and Sutton (2005) find a very similar mechanism in their CGCMs to that of Delworth et al. (1993). In particular, the lagged-phase relationship between temperature and salinity contributions to the total density plays a prominent role. However, both of these studies suggest stronger ties with the NAO than those in the Delworth et al. (1993) and DG00 studies.
In another study using the same GFDL model output, Delworth et al. (1997) suggest a role for the enhanced transports of relatively fresh water and sea ice from the Arctic into the sinking regions to weaken MOC in the North Atlantic. Here, the Greenland Sea oscillations are likely implicated in generating the MOC oscillations. Another possibility is a large-scale atmospheric response to the MOC oscillations, creating the Greenland Sea variability. This latter possibility is in contrast with their previous interpretation of this oscillation. Similar to Delworth et al. (1997), Jungclaus et al. (2005) implicate the storage and release of freshwater from the central Arctic and circulation changes in the Nordic Seas in the MOC variability in their model. Nevertheless, they concur with Delworth et al.’s (1993) conclusion that this is a damped ocean-only mode excited by the atmosphere.
In contrast to above studies, Timmermann et al. (1998) indicate that the Atlantic and Pacific Oceans are coupled via an atmospheric teleconnection pattern and interpret the multidecadal oscillations as an inherently coupled atmosphere–ocean mode. However, the mechanism that provides the oscillation remains the same as in the above studies, that is, the phase delay between temperature and salinity contributions to the total density in the DWF regions.
We note that the role of the flux corrections used in Delworth et al. (1993) and Timmermann et al. (1998) in affecting the period and suggested mechanisms of this variability as well as the patterns and amplitude of the associated SST oscillations remains unclear. In particular, as indicated earlier, these two studies present SST anomaly patterns and amplitudes that broadly resemble those of the observed SST multidecadal variability associated with the AMO. In contrast, the CGCM studies that use no flux corrections tend to produce SST anomaly patterns and amplitudes that differ noticeably from those associated with the AMO (e.g., Dai et al. 2005; Jungclaus et al. 2005). We speculate that such flux corrections constrain the path of a model’s North Atlantic Current (NAC), thus eliminating a persistent model bias that appears to be highly relevant in an analysis of this variability (see section 3b).
Multidecadal MOC oscillations in the Atlantic also exist in all present-day simulations of the National Center for Atmospheric Research (NCAR) Community Climate System Model version 3 (CCSM3), as discussed in Bryan et al. (2006). This study shows that both the mean value of the MOC and amplitude of its variability have significant resolution dependency, with the highest resolution (T85 × 1) simulation having the largest values for both. In particular, the standard deviations of the maximum MOC time series differ by about a factor of 2 [1.24 versus 0.72 Sv (1 Sv ≡ 106 m3 s−1)] between the T85 × 1 and T42 × 1 resolution versions where the ocean model resolution is identical, but the atmospheric model resolutions are 1.4° and 2.8°, respectively. Another preliminary spectral analysis of the annual-mean North Atlantic MOC maximum time series shows dominant and significant (at 99% level) peaks at 21.4 and 20 yr in T85 × 1 and T42 × 1, respectively. In contrast, the coarsest-resolution simulation reveals a dominant and significant peak at a period of about 100 yr. Given these sensitivities to the model resolution and possible role of the atmospheric resolution in determining its response to midlatitude SST anomalies, we have chosen to analyze the T85 × 1 simulation in the current study, because it is the most widely used version of the model, including the Intergovernmental Panel on Climate Change Fourth Assessment Report. A detailed comparison of resolution dependencies of the MOC variability will be the subject of a future work. Therefore, the primary foci of the present diagnostic study are to document this multidecadal variability and investigate its mechanism in this CCSM3 version. Because of its role in the air–sea interactions, we also focus on the SST anomalies and investigate how they are created. Furthermore, we assess the contributions of temperature and salinity to the total density in the related DWF region in comparison with the previous CGCM studies. The role of the atmospheric variability is also discussed. The paper is organized as follows: In section 2, we give a brief description of the model as well as the experiment. The results are presented in section 3. A summary and discussion are given in section 4.
2. Model and experiment descriptions
The NCAR CCSM3 is a coupled climate model that uses no flux corrections (Collins et al. 2006a). The model components are the Community Atmosphere Model version 3 (CAM3; Collins et al. 2006b), the Community Land Model version 3 (CLM3; Oleson et al. 2004; Dickinson et al. 2006), the Parallel Ocean Program version 1.4 (POP1.4; Smith and Gent 2004; Danabasoglu et al. 2006), and the Community Sea Ice Model (CSIM; Briegleb et al. 2004; Holland et al. 2006).
In CCSM3, the land model is on the same horizontal grid as CAM3 and the sea ice model shares the same horizontal grid as the ocean model. In the T85 × 1 resolution version, the atmospheric model uses T85 spectral truncation in the horizontal (about 1.4° resolution) with 26 vertical levels. The ocean model has a nominal 1° horizontal resolution (constant at 1.125° in longitude and varying from 0.27° at the equator to about 0.64° in the far northwest Pacific in latitude) with 40 vertical levels.
The present-day control integration was integrated for about 700 yr, using a global- and annual-mean CO2 mixing ratio of 355 ppmv that corresponds to the 1990 observed value. The ocean model was initialized with the January-mean climatological potential temperature and salinity [a blending of the Levitus et al. (1998) and Steele et al. (2001) datasets] and zero velocities. The remaining components were initialized with January conditions obtained from stand-alone integrations.
As discussed below in section 3a (see Fig. 1a), the MOC multidecadal variability is very regular, particularly between model years 151 and 450. Therefore, we use this 300-yr segment in the present analysis. This choice also avoids the initial adjustment period because of coupling. However, during the analysis period (indeed throughout the entire 700-yr integration), the ocean model potential temperature and salinity fields show continued drift, particularly below 1000-m depth. Consequently, we mostly restrict our analysis to the upper ocean where such drifts are less pronounced.
We use annual-mean fields in the present study except for the heat and salt budget analysis where monthly mean fields are utilized. The time-mean distributions show the 300-yr mean (years 151–450) fields. We employ standard correlation, regression, spectral analysis, and empirical orthogonal function (EOF) methods. All of the time series are detrended using a linear least squares fit prior to analysis. Unless otherwise noted, no time filtering is applied. The EOF time series is normalized to have unit variance, so that the EOF spatial pattern magnitudes can be directly multiplied by the corresponding time series values to obtain the magnitude of an anomaly. The power spectra use the Hanning Window. The reference red noise spectra with the same total variances, computed from the lag-1 autocorrelations, and the associated 95% and 99% confidence levels are shown in all related plots. The significance of the correlations is tested using a two-sided Student’s t test. Here, in situations where two time series have significant differences in their number of degrees of freedom (DOF), based on lag-1 autocorrelations, the time series with the higher DOF is smoothed to match the DOF of the other time series.
a. MOC and northward heat transport variability
The time series of the maximum MOC in the North Atlantic is shown in Fig. 1a. The MOC rapidly increases from about 18 to 28 Sv during the first 40 yr. Thereafter, it starts to develop a decadal time-scale oscillation. Particularly between years 151 and 450, this oscillation becomes regular, varying between 17.5 and 26.5 Sv. This range corresponds to an amplitude of about 4.5 Sv. After year 450, the oscillation is somewhat less regular and its amplitude is lower in comparison to the year 151–450 segment. Given this regularity and larger amplitude between years 151 and 450, we restrict our analysis to this 300-yr segment. Figure 1b shows the corresponding time-mean MOC. The circulation associated with the North Atlantic Deep Water (NADW) has a maximum transport of >20 Sv. Thus, the amplitude of the MOC variability between years 151 and 450 represents a rather significant fraction of the mean maximum circulation, that is, 20%. The circulation associated with the Antarctic Bottom Water (AABW) has a maximum of >2 Sv.
Figures 2a,b present the MOC EOF1 spatial pattern and its time series, respectively, accounting for 71.2% of its total variance. The EOF1 has a single cell pattern, covering the entire Atlantic domain south of 60°N. This (positive) pattern indicates an overall strengthening and deeper penetration of the NADW cell. There is also a corresponding weakening in the AABW circulation. The associated time series is significantly (>99%) correlated with the maximum MOC time series of Fig. 1a, with a simultaneous correlation coefficient value of 0.94. Therefore, we choose to use the MOC EOF1 time series as our reference time series in our analysis below. The power spectrum of the EOF1 time series (Fig. 2c) shows three periods that are significant at the 99% level—16.6, 21.4, and 25 yr. Among these, the 21.4-yr period (hereafter referred to as the 21-yr period) has the largest power. The autocorrelation function for the EOF1 time series (Fig. 2d) reveals that the biggest negative correlation occurs at a lag of 10–11 yr, indicating again that the dominant period is 21 yr. We note that although the power spectrum of the EOF1 time series that includes the latter part of the integration also produces the same three peaks as significant periods, with the 21-yr peak showing the highest power, the power spectrum only for years after 450 does not show a multidecadal peak, but instead an approximately 50-yr period at about 95% significance level appears.
The time-mean NHT in the Atlantic (Fig. 3a) has a maximum transport of 1.14 PW located at 23°N. The direct estimates for this maximum transport are 1.2–1.3 PW between 14° and 24°N, with estimated errors of 0.3 PW (Bryden and Imawaki 2001). Therefore, the model maximum transport is in agreement with these estimates. The local minimum at 43°N is due to a positive heat flux bias, that is, into the ocean, that tries to counteract a cold SST bias discussed in section 3b (see also Figs. 5a,c). The regressions of the NHT time series with those of the MOC EOF1 (Fig. 3b) show enhanced (reduced) NHT with larger (smaller) MOC magnitudes. Using a mean amplitude of 2 for the MOC EOF1 time series (Fig. 2b), the amplitude of the NHT variability is about 0.12 PW at 23°N, corresponding to about 10% of the mean maximum transport. The peak amplitude, however, is collocated at 43°N with the latitude of the mean NHT local minimum.
b. SST and SSS variability
We show the spatial pattern of the SST EOF1 and its time series in Figs. 4a,b, respectively. This EOF accounts for 24.3% of the total SST variance, and the power spectrum of its time series is very similar to the MOC EOF1 spectrum, with the largest power occurring at the 21-yr period (not shown). The correlation function between the MOC EOF1 and SST EOF1 time series (Fig. 4c) shows that the largest correlations occur when the MOC maximum leads the SST time series by 2–3 yr with a significance level of >99%. In Fig. 4d, we present the composite-mean SST difference distribution between the high and low MOC states, shifted by 2 yr to capture the maximum correlation indicated in Fig. 4c. Each composite mean contains six episodes for a total of 24 yr of data. The resulting pattern closely matches that of the SST EOF1 given in Fig. 4a. The regressions of the SST anomalies with the MOC EOF1 time series (not shown) also produce similar patterns as those in Figs. 4a,d, with the most significant correlations (>99%) again occurring when the MOC EOF1 leads the SST anomaly time series by 2–3 yr.
The spatial patterns presented in Figs. 4a,d have a dipole structure with similar positive and negative anomaly magnitudes. However, the negative anomaly pattern is much smaller in its spatial extent than the positive one. Similar potential temperature anomaly patterns with reduced magnitudes exist down to about 1000-m depth, suggesting that they have a barotropic structure (not shown). The peak-to-peak changes in the EOF1 time series (Fig. 4b) exceed 4, thus indicating that the peak-to-peak SST changes can be as large as 6°–7°C, consistent with the magnitudes of the composite-mean difference shown in Fig. 4d. These anomalies represent rather large SST variability on multidecadal time scales. If such anomalies exist in nature, then they should be readily present in the observations. Therefore, for comparison, we follow Kushnir (1994) and plot the observed SST differences between a warm (1950–64) and a cold (1970–84) period in Fig. 4e, using the Hadley Centre Sea Ice and Sea Surface Temperature, version 1 (HadISST1) dataset (Rayner et al. 2003). The observations largely show a basin-scale pattern with differences rarely exceeding 0.8°C. Indeed, in regions where the SST EOF1 anomalies and composite-mean differences have large magnitudes, the observational differences are only 0.4°–0.6°C, which is an order of magnitude smaller then the model anomalies. Another comparison is made with the dipole SST mode of Deser and Blackmon (1993) with a 9–12-yr period. Their dipole mode has more of a north–south orientation as opposed to the east–west structure of Figs. 4a,d. Moreover, Deser and Blackmon (1993) show about 40% smaller amplitude compared to the present values. Therefore, we conclude that neither the spatial pattern nor the magnitude of this dominant SST variability in CCSM3 is realistic as is also the case in some other CGCM studies (e.g., Jungclaus et al. 2005). Alexander et al. (2006) suggest a similar conclusion based on their analysis of wintertime SST anomalies in CCSM3 (but they do not show any related figures). We note that the SST EOF2, accounting for 13.6% of the total SST variance, has a spatial pattern (not shown) more closely resembling that of Deser and Blackmon (1993).
The time-mean SST and sea surface salinity (SSS) distributions in the North Atlantic, along with their differences from observations, are presented in Fig. 5. As documented in Large and Danabasoglu (2006), in CCSM3 the largest time-mean SST bias occurs in the North Atlantic. This cold bias is about 9°C (Fig. 5c), roughly collocated with the regions of the largest SST anomalies (about 40°–50°N, 50°–20°W) shown in Figs. 4a,d. There is a density-compensating large, fresh bias in SSS (about 4 psu) in the same region (Fig. 5d). Large and Danabasoglu (2006) attribute these mean biases to the incorrect path of the model NAC as shown in Fig. 5e. Although the separation point of the model Gulf Stream is realistic, its extension (i.e., NAC) remains too zonal as it crosses the North Atlantic, and it does not become northeastward until after 40°W. Consequently, the subpolar gyre penetrates farther south, thus producing these big mean biases in the presence of large SST and SSS meridional gradients (Figs. 5a,b). As in earlier versions of the CCSM (Danabasoglu 1998), the interior gyre transports in CCSM3 are generally set by the wind stress curl (WSC), following Sverdrup dynamics. We show the time-mean WSC distribution in Fig. 5f, revealing that WSC magnitudes remain rather small (indeed near zero) in regions where these large SST and SSS biases and anomalies exist. Given the sharp meridional gradients in SST and SSS in this region, small north–south excursions of the subtropical–subpolar gyre boundary can lead to large SST and SSS anomalies. These sorts of subtropical–subpolar gyre boundary fluctuations can result from small fluctuations in WSC because a positive (negative) WSC anomaly imparts (extracts) vorticity to (from) the ocean, thus allowing water parcels to move northward (southward) crossing planetary vorticity contours.
To test the above hypothesis, we first create a time series to represent these north–south excursions of the subtropical–subpolar gyre boundary. Here, we simply use the latitudinal changes of the zero contour line of the barotropic streamfunction (BSF) along 30°W (shown in Fig. 5e). The resulting time series (not shown) has a maximum range of 4° with typical peak-to-peak oscillations of 2°–3°. The associated power spectrum (Fig. 6a) shows that the highest power occurs at a period of 21 yr. We next present the simultaneous regressions of this time series with WSC in Fig. 6b, where a 5-yr running mean is applied to WSC prior to analysis. For reference purposes, the figure includes the most northward and southward mean positions of the gyre boundary as represented by the zero BSF contour line. These mean positions are based on six episodes, each 4 yr long, when the gyre boundary is at its most northern (NORTH) and southern (SOUTH) positions, and thus represent 24-yr means each. The figure clearly indicates that the NAC path shifts southward and northward in response to the negative and positive WSC anomalies centered at about 43°N, 45°W and 45°N, 35°W, respectively, in agreement with our hypothesis. As stated earlier, because these shifts occur in regions with large meridional SST and SSS gradients, the resulting SST and SSS changes are rather large as demonstrated by the NORTH–SOUTH difference distributions of Figs. 6c,d. During a NORTH phase, the NAC takes a more northeasterly path east of 40°W, bringing warmer waters from south and creating a positive SST anomaly of 6°–7°C (Fig. 6c). The negative SST anomaly reaches −5°C and is a direct consequence of the negative WSC anomaly pushing colder waters southward. The presence of a narrow patch of relatively fresh water (e.g., 34 psu in Fig. 5b) in the mean SSS distributions results in the dipole pattern of the SSS difference distributions (Fig. 6d). These SSS differences can be as large as 3 psu. We note that both the patterns and magnitudes of these difference distributions are very similar to their respective EOF1 patterns and magnitudes (shown for SST in Figs. 4a,b). The negative WSC anomaly centered at about 41°N, 35°W also contributes to local transport changes, implying a larger southward transport, and hence affects SST and SSS distributions. However, the meridional gradients, particularly in SST, are not as large here.
The WSC anomalies implicated with the gyre boundary fluctuations (Fig. 6b) are quite local and exhibit rather small-scale features. The associated correlation function coefficients (not shown) are the largest at lag 0, that is, when they are simultaneous, and when the WSC anomalies lead the gyre boundary fluctuations by 1 yr, suggesting a fast barotropic ocean response. We construct a representative time series for these WSC fluctuations by spatially averaging the time-filtered WSC in 44°–46°N, 40°–30°W, designed to include the positive WSC pattern centered at 45°N, 35°W in Fig. 6b. Its power spectrum (Fig. 7) reveals several significant peaks at the 99% level, including the one at 21 yr. Based on these analyses, we believe that the large midlatitude SST and SSS anomalies are not directly associated with the MOC oscillations. Instead, these anomalies are actually created by the oscillations of the subtropical–subpolar gyre boundary driven by the WSC anomalies. Thus, the statistical relationship between the SST EOF1 and MOC EOF1 time series (Fig. 4c) does not indicate causality.
Recently, Zhang and Vallis (2007) suggested a more direct mechanism that ties the strength of the MOC to the fluctuations of the subtropical–subpolar gyre boundary. This mechanism depends on bottom vortex stretching induced by a downsloping deep western boundary current (DWBC) and subsequent formation of a northern recirculation gyre (NRG). In our simulations, neither such a NRG exists nor are there significant shifts of the gyre boundary west of 50°W (Fig. 5e).
An EOF analysis of the WSC time series shows that these rather local, small-scale patterns do not appear in the first two EOFs, but are revealed as EOF3 and EOF4, accounting for about 10% of the total variance (not shown). The associated time series show only weak and insignificant correlations with the model NAO time series, when the latter leads by 5–7 yr. In this study, we define the NAO as the wintertime [December–March (DJFM)] sea level pressure EOF1. In an analysis of the North Pacific in CCSM2, Kwon and Deser (2007) also report similar meridional shifts of the Kuroshio Extension resulting from some WSC anomalies, leading to decadal SST variability. However, their WSC anomalies remain basin scale. We defer further discussion to section 4.
c. Exploration for a mechanism and role of the atmosphere
To investigate the cause of the oscillations in MOC, we now focus on the northern North Atlantic where DWF occurs. We use the boundary layer depth (BLD) to determine both the formation sites and the variability of the DWF in the model. Because BLD attains its maximum depths in March, our analysis is based on the March-mean values. BLD is determined by the K-profile parameterization (KPP) vertical mixing scheme (Large et al. 1994) as the shallowest depth at which a bulk Richardson number exceeds a specified critical Richardson number for the first time. The March-mean BLD climatology (Fig. 8a) shows essentially the following three formation sites: in the Greenland–Iceland–Norway (GIN) Sea between Iceland and Spitsbergen, mostly along the ice edge (Holland et al. 2006); south of Iceland between Greenland and Scotland; and south of Greenland, mostly in the Labrador Sea basin. These time-mean formation sites are the same as those shown in Bryan et al. (2006) where the maximum mixed layer depths are used instead. While the GIN Sea and south of Iceland sites show maximum depths of only about 500–600 m, BLD exceeds 1300 m in the Labrador Sea site. The BLD EOF1 distributions (Fig. 8b) clearly indicate that the Labrador Sea region is the only DWF site with substantial variability that penetrates the abyssal ocean. This EOF accounts for 44.4% of the total variance. The associated EOF1 time series (not shown) has a range between −1.5 and 3. Therefore, in the Labrador Sea site, the March-mean BLD can be about 200 m at its shallowest and reach the ocean bottom (about 3500 m) at its deepest. In contrast, the range of variability is only from about −30 to 60 m about their respective means at the other two DWF sites. The power spectrum of the BLD EOF1 time series has its largest significant peak at a period of 21 yr (Fig. 8c). Figure 8d shows the correlation function between the MOC EOF1 and BLD EOF1 time series, indicating that the maximum and minimum BLDs occur about 4–5 yr prior to a MOC maximum and minimum, respectively.
Figure 9a presents the density (ρ) regressions with the MOC EOF1 time series for the Labrador Sea DWF site identified by the boxed region in Fig. 8b. The regression coefficients are simply volume averaged in this box for the upper 0–212-m depth range. The figure also includes the individual temperature (ρT) and salinity (ρS) contributions to ρ. The densest upper-ocean waters form about 3–5 yr prior to a maximum in MOC. As expected, this indicates that the deepest BLDs (Fig. 8) occur in response to these increased surface densities. Similarly, the BLD and ρ minimums are approximately simultaneous, that is, both lag MOC maximum by about 5 yr. In contrast to some previous studies (e.g., Delworth et al. 1993; Dong and Sutton 2005; Dai et al. 2005), we find that ρT and ρS contributions to ρ are almost equal and lag each other by only 2–3 yr, particularly during the minimum and maximum phases of ρ.
While this phase lag represents about 1/8 of the 21-yr period, it is substantially smaller then in the previous studies in which the mechanism for the MOC oscillation is largely attributed to a phase lag of at least 1/4 of the MOC period, or sometimes an entirely out-of-phase relationship, between ρT and ρS in their contributions to ρ. To explain this apparent discrepancy between the present study and previous ones, we next evaluate the same regression coefficients over a much larger horizontal area (50°–80°N, 70°W–20°E) and the entire ocean depth. The resulting plots are given in Fig. 9b, clearly showing a largely out-of-phase relationship between ρT and ρS. Another average based on the same larger horizontal domain, but only for the upper 0–212-m depth range (not shown), also reveals similar distributions as in Fig. 9b, suggesting that the horizontal domain used in averaging plays a primary role here. In our view, such large horizontal and vertical averaging domains do not properly reflect the phase relationships and balance of terms in a model’s actual DWF sites, particularly in higher-resolution models as demonstrated in Fig. 9a. Moreover, nonlinear drifts at depth in potential temperature and salinity may further adversely affect analysis when regressions are performed over the full depth of an ocean model. Therefore, we advocate performing this sort of analysis for the actual DWF sites that feed the MOC, and that are found to be associated with its decadal variability. Here, we conclude that a Delworth et al. (1993) type of ocean mode does not exist in CCSM3.
We note that the delay time (about 5 yr) between the maximum of ρ and BLD and the MOC maximum corresponds to about a quarter of the 21-yr period in CCSM3 and may play a role in setting the time scale of this variability. Similar lag times (3–6 yr) are also reported in previous studies, regardless of the period of their dominant oscillation. For example, the lag time and the oscillation period are 4 and 50 yr, respectively, in Delworth et al. (1993). Eden and Willebrand (2001) give a detailed discussion of some possible mechanisms for this delay. These include boundary wave propagation and advection of density anomalies by the DWBC. They also suggest that the resolution of a model can influence the delay time.
Undoubtedly, the upper-ocean heat and freshwater budget analyses for the DWF site identified in Fig. 8 provide the most comprehensive picture of how the ρT and ρS anomalies that are associated with the BLD and MOC oscillations are created. However, before performing such analyses, it is useful to examine the variability and associated forcing of the northern North Atlantic circulation. Figure 10 shows the first two BSF EOFs, accounting for about 29% and 18% of the total variance, respectively. Figure 10a shows that, north of 45°N, the counterclockwise and clockwise circulation patterns imply strengthening and weakening of the mean subpolar gyre in the Labrador Sea basin and to the east of 40°W, respectively, in a seesaw pattern. When multiplied by the corresponding time series (not shown), these positive and negative anomalies correspond to peak-to-peak changes of about 9 and 21 Sv, respectively. The EOF2 distribution (Fig. 10b), in contrast to EOF1, shows a single counterclockwise circulation pattern north of 45°N, indicating an overall weakening and strengthening of the subpolar gyre. The peak-to-peak circulation changes reach about 18 and 9 Sv for the local maxima at about 50° and 25°W, respectively. The power spectrum of the BSF EOF1 time series (Fig. 10c) shows a distinct peak at the 21-yr period, while the BSF EOF2 time series (Fig. 10d) has an additional peak at a period of 25 yr. The correlation functions between the MOC EOF1 and BSF EOF1 and EOF2 time series are presented in Fig. 10e. The figure indicates that the BSF EOF1 and EOF2 are in quadrature. While the BSF EOF1 time series lag the MOC EOF1 time series by 1–2 yr, BSF EOF2 time series lead the MOC EOF1 time series by 3–4 yr. Consequently, the subpolar gyre circulation is enhanced between lag −6 and lag +5 in the Labrador Sea. At 35°W, where the eastern boundary of the DWF budget region resides (see Fig. 8b), however, the weakening of the subpolar gyre transport starts at about lag +3.
Figure 11 presents the WSC EOF1 distribution along with the power spectrum of its time series. A 5-yr running mean is applied to WSC prior to analysis. Under ice-covered regions, for example, east of Greenland, this field represents the curl associated with the stress between the ocean and the sea ice. The WSC EOF1 accounts for 25.7% of the total variance and it has a significant peak again at the 21-yr period. The WSC EOF1 pattern shows a broad dipole structure with positive anomalies north of 60°N and negative anomalies between 45° and 60°N. The time series of the WSC EOF1 is significantly (>99% level) correlated with the time series of the NAO at lag 0 (not shown) and the anomaly pattern of Fig. 11a corresponds to the positive phase of the NAO. The positive phase of the NAO shows a stronger-than-normal subtropical high pressure center.
The regressions of the WSC with the BSF EOF1 time series (not shown) produce very similar anomaly patterns and magnitudes north of about 50°N, as in Fig. 11a, when the WSC leads the BSF time series by 4–5 yr. This apparent lagged response in ocean circulation to the WSC anomalies is verified by the correlation function between the WSC EOF1 and the BSF EOF1 time series given in Fig. 12a. We note that the significance is only at about a 90% level. To match the degrees of freedom between the BSF and WSC time series, a 9-yr running mean is applied to the latter for the regression and correlation analysis. Thus, the analysis suggests a delayed relationship between WSC EOF1 (or NAO) and BSF EOF1. Such a delayed response of the ocean circulation to the NAO is reported by Eden and Willebrand (2001) in their ocean-only experiments forced with realistic surface fluxes. They speculate that the mechanism is likely related to baroclinic wave adjustment that leads to a slow barotropic mode adjustment via an interaction between baroclinic and barotropic modes in the presence of topography. A similar delayed response has also recently been found in the HadCM3 coupled model by Dong and Sutton (2005). However, while Eden and Willebrand (2001) report a 6–8-yr delay time, it is only 1–2 yr in Dong and Sutton (2005). Unfortunately, the observational evidence for this type of delayed response is inconclusive. While Taylor and Stephens (1998) support a delay of 2 yr in the Gulf Stream response, Joyce et al. (2000) favor a synchronous response.
A similar analysis between the WSC and BSF EOF2 time series shows neither a recognizable WSC pattern nor a statistically significant relationship. Nevertheless, the strengthening of the subpolar gyre circulation indicated by the BSF EOF2 starts at about lag −6 and occurs during the positive phase of the NAO. Another possibility is that BSF EOF2 is related to the density anomalies in the vicinity of the DWF site. As discussed in Delworth et al. (1993) and Dong and Sutton (2005), an enhanced subpolar gyre can be due to the positive density anomalies in this region because the resulting surface height (dynamic topography) anomalies are related to the (surface) currents through geostrophy. Although this relationship is very appealing because it provides a simple and direct link between the MOC and the subpolar gyre circulation variability, we believe that the WSC anomalies (or NAO) play a role in creating the density anomalies in the first place (as discussed below) as well as directly modulating the subpolar gyre strength.
For completeness, the correlation function between the WSC EOF1 (or NAO) and MOC EOF1 time series is included in Fig. 12b, showing that the WSC leads MOC time series by 2 yr. Thus, Figs. 10e and 12 summarize the following phase relationships between NAO, MOC, and BSF EOF1 time series: a positive NAO occurs 2 yr prior to a MOC maximum, and the delayed BSF ocean response to this NAO phase occurs about 2 yr after a MOC maximum.
We now return to the computation of the upper-ocean heat and freshwater budgets for the DWF site identified by the boxed region in Fig. 8b to examine the creation of the density anomalies (see the appendix for details). The results are presented in Fig. 13 as regressions with the MOC EOF1 time series. The positive (negative) fluxes indicate heat and freshwater input (loss) to (from) the budget region. In addition to the budget terms, the temperature (T) and salinity (Sa) regressions are included in Figs. 13b,d, respectively. Although the salinity anomaly amplitude (about 0.02 psu) is much smaller than the temperature anomaly amplitude (about 0.1°C), their contributions to density remain about equal, because the salinity contraction coefficient is much larger than the thermal expansion coefficient at these low temperatures.
Figures 13a,c clearly show a largely out-of-phase relationship between the two largest advective components of the heat and freshwater budgets in our model. Namely, the advective transport through the southern boundary (S) is largely compensated for by the advective transport through the eastern and northern (E + N) boundaries of the budget region. For example, as the MOC spins up, S brings warmer and saltier waters northward, starting at lag −7. These S transports peak 1–2 yr prior to a MOC maximum. In contrast, E + N fluxes show increased cooling and freshening of the upper ocean during the same period, that is, from lag −7 to +2, in which the subpolar gyre transport is larger, substantially canceling S. Later, reversal of this cooling and freshening trend in E + N starts at about lag +3, coinciding with the weakening of the subpolar gyre transport at the eastern boundary of the budget region, as depicted in Fig. 10a for BSF EOF1, and continuing through the negative phase of the NAO. Again, this warming and increased salt trend in E + N is opposing the cooling and freshening trend in S during this later period. Further analysis of all advective fluxes indicates that the advection of the mean temperature and salinity by the anomalous circulation is the primary contributor to the advective terms (not shown). We note that the net contributions of the BSF fluctuations to the creation of density anomalies in the DWF region are included in the advective budget terms. An increase (decrease) in the subpolar gyre circulation indicates increased (reduced) transport both into and out of the budget region. Therefore, it is crucial to consider the net effect of these transport changes in a budget analysis.
As a result of this significant cancellation between S and E + N, the remaining budget terms and how their phase relationships are determined become important. For example, both the sign and phase of the vertical advective fluxes at 212-m depth (B) reflect those of the vertical velocities associated with the WSC anomalies. The negative (positive) WSC anomalies produce anomalous Ekman downwelling (upwelling), taking (bringing) heat and salt away from (to) the budget volume. These negative (positive) WSC anomalies are associated with the positive (negative) phase of the NAO, and they peak 1–2 yr before (7–8 yr after or 10–11 yr before) a MOC maximum. As shown in Fig. 13c, B is particularly responsible for the initiation of positive (negative) salinity anomalies during a MOC minimum (maximum). The total advective heat flux has a minimum at lag −3 (Fig. 13b), in contrast with the total advective freshwater flux showing a maximum at lag −1 (Fig. 13d). Therefore, while these total advective fluxes remain largely out of phase in their respective contributions to density, surface and diffusive (DIFF) fluxes become important. We believe that the warming and reduced freshening resulting from DIFF before a MOC maximum are attributable to increased convection during this period.
The sea ice plays a role in determining the phases of the surface flux components through its response to the NAO and MOC variability. The regressions of the sea ice extent time series with those of the NAO and MOC EOF1 (not shown) indicate opposing sea ice responses. While the positive NAO phase results in increased sea ice cover in the Labrador Sea basin (see also Holland 2003), including our budget region, increased NHT resulting from enhanced MOC reduces sea ice extent. The time lag between these opposing effects is about 2 yr (the NAO leading). The freshwater fluxes associated with melting sea ice are by far the dominant component in surface freshwater fluxes (not shown). Nevertheless, the surface freshwater flux contributes little to its budget in the DWF region (Fig. 13d) in contrast with the surface heat flux (Fig. 13b). These surface heat flux anomalies, which also include the effects of the opposing responses of the sea ice extent to the NAO and MOC, first reinforce the negative SST anomalies prior to a maximum in BLD and density. Then, they have a damping effect after about lag −4 (Fig. 13b). This heat gain by the ocean along with the positive diffusive fluxes stop the cooling trend at about lag −2. In CCSM3, the sea ice cover in the Labrador Sea region is somewhat too extensive compared to observations (Holland et al. 2006). Consequently, the anomalies associated with the NAO are not confined to the shelf/ice-edge regions, but extend over the Labrador Sea basin (Holland 2003). It is unclear how these model biases influence the budget terms, particularly surface fluxes, discussed here.
4. Summary and discussion
The multidecadal variability in the Atlantic Ocean that exists in the CCSM3 present-day simulations is investigated, using the highest (T85 × 1) resolution version. Many ocean fields in the North Atlantic show statistically significant oscillations at a 21-yr period. In MOC, the oscillation amplitude is about 4.5 Sv, corresponding to 20% of the mean maximum MOC transport. The associated EOF has a single-cell pattern, and in its positive phase it indicates an overall strengthening and deeper abyssal penetration of the NADW circulation. The variability in the NHT has an amplitude of about 0.12 PW, representing 10% of the mean maximum NHT. In SST and SSS, the peak-to-peak changes can be as large as 6°–7°C and 3 psu, respectively. Similar to some other CGCM studies, comparisons with observations show that neither the pattern nor the magnitude of this dominant SST variability is realistic. The regions with the highest SST and SSS variability are roughly collocated with the regions of large mean SST and SSS model biases. Because these biases are attributed to the incorrect path of the model NAC, we believe that the elimination of this model bias can significantly improve the simulation of SST and SSS variability in the North Atlantic. In the tropical Atlantic, the peak-to-peak SST changes remain modest (0.1°–0.2°C), but may be large enough to prompt an atmospheric response (Dong and Sutton 2005). These oscillations become more regular and, perhaps, predictable only after the initial transient period that can last 150–200 yr in CCSM3.
We identify the Labrador Sea region as the DWF site associated with the MOC oscillations. The positive density anomalies and the resulting deep BLD reach their maxima 5 yr before a MOC maximum. The heat and freshwater budget analyses performed for this DWF region indicate a complex relationship between the DWF, MOC, NAO, and subpolar gyre circulation anomalies. We believe that their complicated interactions are responsible for the maintenance of this multidecadal oscillation. Because the beginning point of such an oscillation is arbitrary, we can start with a minimum in MOC at about lag −10. This state occurs about 2 yr after a negative NAO phase and is associated with rather weak subpolar gyre circulation. During this phase, the oceanic heat loss increases, accompanying negative surface heat flux anomalies, and vertical upwelling produces negative freshwater flux anomalies. In turn, they initiate lower temperature and higher salinity anomalies. Here, the negative phase of the NAO appears to be particularly responsible for the upwelled freshwater flux anomalies. These cold and saline waters start to form the positive density anomalies in the DWF site. At lag −7, the WSC anomalies become rather weak, indicating a transition from a negative to positive phase of the NAO. Partly associated with this transition, the subpolar gyre starts to gain strength. At about lag −5, the density anomalies reach their maximum, producing deepest BLD anomalies. At lag −2 the positive NAO phase peaks. During this phase, the positive surface heat flux anomaly, together with the diffusive fluxes, initiates the reversal of the cooling trend. At the same time, the vertical advective fluxes contribute to the freshening of the upper-ocean waters. As the MOC attains its maximum transport at lag 0, both temperature and salinity anomalies are already reversing their signs to produce negative density anomalies, thus completing a half-cycle within about 10–11 yr. We note that the phasing of the surface heat flux anomalies also include the opposing responses of the sea ice extent to the NAO and MOC. The resulting heat fluxes first reinforce the negative SST anomalies prior to a maximum in BLD and density. Then, they have a damping effect after about lag −4.
We find that the large, midlatitude SST and SSS anomalies are not directly associated with the MOC oscillations. Instead, they are created by the fluctuations of the subtropical–subpolar gyre boundary driven by small-scale WSC anomalies that are quite different in their spatial patterns than the ones related to the NAO. The SST and SSS anomalies are simultaneously correlated with these WSC anomalies that exist about 5 yr after (before) a positive (negative) phase of the NAO. Therefore, we deduce that these WSC anomalies occur during the transition phase between the positive and negative NAO states. We conclude that the high correlations between the SST and SSS anomalies and the MOC time series do not indicate causality, in contrast to the conclusion in Dai et al. (2005). With the caveat that the model SST and SSS anomaly patterns and magnitudes do not resemble observations (e.g., those associated with the AMO), they also appear to contradict Bjerknes (1964), where decadal SST changes are related to basin-scale circulation changes.
In contrast with some previous CGCM studies (e.g., Delworth et al. 1993; Dong and Sutton 2005; Dai et al. 2005), our study does not support the existence of a Delworth et al. (1993) type of ocean mode that relies on the phase-lagged relationship between the temperature and salinity contributions to the total density in the DWF regions. Instead, the atmospheric variability associated with the model’s NAO appears to play a prominent role in maintaining the multidecadal oscillations. In particular, the NAO modulates the subpolar gyre strength and contributes to the formation of the temperature and salinity anomalies that lead to positive/negative density anomalies at the DWF site. In addition, the WSC anomalies occurring during the transition phase between the positive and negative NAO states are responsible for the creation of midlatitude SST and SSS anomalies. The power spectrum of the observed NAO index is slightly red, with somewhat enhanced, but statistically insignificant, variances between 2–3- and 8–10-yr bands (Hurrell et al. 2003). We present the power spectrum of the model-filtered NAO time series in Fig. 14. It shows that there are several statistically significant peaks, including the one at the 21-yr period. Interestingly, the model spectrum also shows a significant peak at an 8-yr period. Although the NAO is considered to be an internal mode of the atmosphere, as discussed in Hurrell et al. (2003), interactions between the atmosphere and ocean can be argued to modulate the NAO variability, particularly on long time scales. It is likely that the processes setting the 21-yr time scale have oceanic origins and the corresponding peak in the model NAO is merely a reflection of these processes. For example, an advective mechanism that is partially excited by atmospheric forcing, as suggested by Saravanan et al. (2000), can represent such an oceanic process.
The level of atmospheric and oceanic interactions is partially dictated by the response of an atmospheric general circulation model (AGCM) to midlatitude SST anomalies. This response appears to vary considerably among AGCMs, possibly depending on the model resolution and the region, magnitude, and duration of the imposed SST anomalies [see Saravanan and McWilliams (1997) and Neelin and Weng (1999) for related discussions]. In the North Atlantic, Magnusdottir et al. (2004) document that an earlier version of the CCSM atmosphere component is on the weak-response side of this spectrum. Therefore, we believe that the large midlatitude SST anomalies associated with the multidecadal oscillation do not significantly influence the atmospheric circulation. However, using a simple model, Neelin and Weng (1999) argue that even such a weak atmospheric response to decadal SST anomalies can in turn reexcite these SST anomalies, thus maintaining a coupled mode. In addition, there may be the possibility of a remotely forced coupled mode, because the Northern Hemisphere SST regressions with the MOC EOF1 time series (not shown) reveal broad SST anomalies in the midlatitude Pacific Ocean with order 1°C peak-to-peak changes. In the Pacific Ocean, the CCSM3 atmosphere component shows a somewhat stronger response to imposed midlatitude anomalies (Kwon and Deser 2007). Consequently, an atmospheric teleconnection between the Pacific and Atlantic may exist. This is broadly similar to the proposed coupled mode of Timmermann et al. (1998). The pattern of these Pacific SST anomalies resembles that of the Pacific decadal oscillation. Although the implication of such a Pacific–Atlantic coupled relationship is interesting, we believe that it is unlikely to be a mode in CCSM3 for the following reasons. First, the simultaneous SST regression coefficients show insignificant correlations with the MOC EOF1 time series. Second, these correlations get much larger, but still remain insignificant (at about an 80% level), when the Pacific SST anomalies lead those of the Atlantic by about 3 yr. Finally, the Pacific SST anomalies show variability at a period of 23 yr, rather than the 21-yr period seen in the North Atlantic.
We believe that the mean model biases in the North Atlantic in CCSM3 need to be addressed prior to the investigations of the nature of the oscillatory mode, that is, either coupled or ocean only, and how the 21-yr period is set. A final remark concerns the associations between various processes in the presence of a rather regular oscillation discussed here. In such an oscillation, processes may have indirect relationships, for example, between MOC and SST anomalies; MOC oscillations may be responsible for the presence of the small-scale WSC anomalies that create the midlatitude SST and SSS variability through an atmospheric response.
The author would like to thank Drs. F. O. Bryan, C. Deser, P. R. Gent, M. M. Holland, Y.-O. Kwon, W. G. Large, R. Saravanan, and J. J. Tribbia for helpful discussions and three anonymous reviewers for their constructive comments and suggestions.
Heat and Freshwater Budget Details
The budget volume is defined by the boxed area in Fig. 8b and extends from the surface to a 212-m depth. The computations are based on the model tracer equations, which can be written as
where TEN is tendency, S is the horizontal advective flux from the southern face, E + N is the horizontal advective flux from the eastern and northern (segment east of Greenland) faces, LAB is the horizontal advective flux from the Labrador Seaside, B is the vertical advective flux from the bottom, that is, 212-m depth, and SFLX is the total surface flux. The positive (negative) fluxes indicate heat and freshwater input (loss) to (from) the budget region. Our computation of individual budget terms is based on the monthly mean model output. We then create annual-mean, detrended time series for each component and regress these with the MOC EOF1 time series. TEN is evaluated based on the differences of two subsequent January-mean values and the diffusive fluxes (DIFF) are computed as residuals. Therefore, DIFF contains both isopycnal and vertical (including convection) diffusion as well as the eddy-induced advection contributions. The advective fluxes are obtained based on the conservative form of the model equations, and hence involve tracer times plane-normal velocity component multiplications. Also, the model temperature is in degrees Celsius. All fluxes (except SFLX) are normalized by the surface area of the analysis region so that all terms have the same units (i.e., W m−2 and Kg m−2 s−1, respectively) for heat and freshwater fluxes. TEN is not shown in Fig. 13.
* The National Center for Atmospheric Research is sponsored by the National Science Foundation.
Corresponding author address: Dr. Gokhan Danabasoglu, National Center for Atmospheric Research, P.O. Box 3000, Boulder, CO 80307. Email: email@example.com