The amplification of thermohaline circulation (THC) anomalies resulting from heat and freshwater forcing at the ocean surface is investigated in a zonally averaged coupled ocean–atmosphere model. Optimal initial conditions of surface temperature and salinity leading to the largest THC growth are computed, and so are the structures of stochastic surface temperature and salinity forcing that excite maximum THC variance (stochastic optimals). When the THC amplitude is defined as its sum of squares (equivalent to using the standard L2 norm), the nonnormal linearized dynamics lead to an amplification with a time scale on the order of 100 yr. The optimal initial conditions have a vanishing THC anomaly, and the complex amplification mechanism involves the advection of both temperature and salinity anomalies by the mean flow and of the mean temperature and salinity by the anomaly flow. The L2 characterization of THC anomalies leads to physically interesting results, yet to a mathematically singular problem. A novel alternative characterizing the THC amplitude by its maximum value, as often done in general circulation model studies, is therefore introduced. This complementary method is shown to be equivalent to using the L-infinity norm, and the needed mathematical approach is developed and applied to the THC problem. Under this norm, an amplification occurs within 10 yr explained by the classic salinity advective feedback mechanism. The analysis of the stochastic optimals shows that the character of the THC variability may be very sensitive to the spatial pattern of the surface forcing. In particular, a maximum THC variance and long-time-scale variability are excited by a basin-scale surface forcing pattern, while a significantly higher frequency and to some extent a weaker variability are induced by a smooth and large-scale, yet mostly concentrated in polar areas, surface forcing pattern. Overall, the results suggest that a large THC variability can be efficiently excited by atmospheric surface forcing, and the simple model used here makes several predictions that would be interesting to test using more complex models.
The North Atlantic Ocean sea surface temperature (SST) and salinity (SSS) exhibit variability on different time scales from interannual (e.g., Levitus 1989) to interdecadal and decadal (Kushnir 1994) that are often attributed to the variability of the thermohaline circulation (THC). Numerical models also show THC oscillations on decadal time scales (e.g., Bryan et al. 2006; Delworth and Greatbatch 2000; Jungclaus et al. 2005) that can possibly be explained as a small-amplitude damped oscillatory linear mode excited by atmospheric stochastic forcing (e.g., Griffies and Tziperman 1995).
While a large number of studies have dealt with the stability and variability of the THC, the nonnormal dynamics of the THC have received only very little attention. In stable linear systems, small perturbations may undergo a strong amplification (before eventually decaying) due to the interaction of several nonorthogonal eigenmodes if the dynamical system is nonnormal (e.g., Farrell and Ioannou 1996). The initial conditions leading to such transient growth are referred to as optimal initial conditions, and the spatial structure of the stochastic forcing leading to maximum variance of the model solution is referred to as stochastic optimal (Farrell and Ioannou 1996; Kleeman and Moore 1997). These concepts of optimal initial conditions and stochastic optimals were applied to several weather and climate variability phenomena (e.g., Chhak et al. 2006; Farrell 1988, 1989; Moore 1999; Moore and Kleeman 1997a, b; Penland and Sardeshmukh 1995). The nonnormal dynamics and predictability of the Stommel box model were discussed by Lohmann and Schneider (1999), and in addition, transient amplification and stochastic optimals emerging from the nonnormality of the THC dynamics were investigated by Tziperman and Ioannou (2002) and Zanna and Tziperman (2005, hereinafter ZT).
Optimal initial conditions that extend throughout the ocean depth leading to maximum THC growth may be excited by internal ocean variability due to the mesoscale eddy field. However, the ocean circulation is more typically forced by stochastic atmospheric forcing at its surface (Battisti et al. 1995; Eden and Willebrand 2001; Halliwell 1998; Seager et al. 2000). The role of atmospheric surface forcing in efficiently exciting THC variability through nonnormal dynamics has not been explored yet. Sevellec et al. (2007) have recently calculated in a 2D ocean model the stochastic optimals of surface salinity only. Air–sea heat and freshwater fluxes affect the temperature and salinity at the surface of the ocean and thus result in the surface excitation of the THC variability. The first objective of the present work is therefore to investigate the effect of surface forcing on the nonnormal dynamics of the THC in a simple coupled ocean–atmosphere model.
The calculation of the optimal initial conditions using a standard sum-of-squares measure for the amplitude of the THC perturbations at all spatial locations (i.e., L2 norm) leads to physically meaningful results yet to a singular mathematical problem (Tziperman and Ioannou 2002; Zanna and Tziperman 2005). A second objective of this paper is therefore to introduce an alternative in which the amplitude of the THC is measured by its spatial maximum value. We show that this alternative is equivalent to using the L∞ norm and then introduce the needed mathematical machinery to find the optimal initial conditions under this norm, and finally apply it to the THC problem. We explain that using the two different norms corresponds to asking two different and complementary questions about the dynamics of the THC growth. This new formalism may be useful in many other transient growth physical problems where the maximum norm is of interest. The introduction of the L∞ norm here completes a set of tools for the study of linear dynamics that so far included the standard L2 norm (e.g., Farrell and Ioannou 1996) and L1 norm (sum of absolute values; Farrell and Ioannou 2000).
Section 2 provides a brief description of the model. We then derive and discuss the optimal surface initial conditions of temperature and salinity leading to an amplification of the THC in section 3 using the L2 and the L∞ norms. The THC stochastic optimals in our model are presented in section 4, followed by an analysis of the results obtained when forcing the model with the first two stochastic optimals of the THC. We conclude in section 5.
2. Model description
The two-dimensional zonally averaged (3° resolution and 61 meridional grid points per level) coupled ocean–atmosphere box model, described in detail by ZT and developed by Sayag et al. (2004), was used for this study. The atmosphere, composed of a single vertically averaged layer, is an energy balance model for the potential temperature θ and humidity q. The ocean model has two vertical layers, and the temperature T and salinity S are governed by advection–diffusion equations.
The linearized model may be written as
where 𝗠 is the linearized model operator and P is the perturbation vector specifying the time-dependent model anomalies with respect to the steady state: atmospheric temperature θ and moisture q, as well as temperature T and salinity S for the surface and deep ocean layers at all grid points (P has dimensions Ny × 1, where Ny = 366).
A linearized equation of state was used so that the THC intensity varies linearly with temperature and salinity and can therefore be rewritten as
3. Surface optimal initial conditions
a. L2 norm
The explicit solution to the linearized model in (1) is
where 𝗕(t, 0) = e𝗠t is the propagator matrix and P(0) ≡ P0 is the perturbation state vector at t = 0. The system is linearly stable, such that all eigenvalues of 𝗠 have negative real parts.
We wish to evaluate the initial temperature and salinity anomalies leading to maximum amplification of the spatial sum of the squares of the THC anomaly at time τ,
while requiring the initial temperature and salinity anomalies to be restricted to the surface ocean. The surface restriction of the initial anomalies is obtained by demanding that ‖P0‖2Y = PT0𝗬P0 = 1, where ‖ ‖ is the norm of a vector and the norm kernel 𝗬 is given by 𝗬 = 𝗫 + 𝗙S, with 𝗙S being a diagonal matrix such that the elements corresponding to the surface temperature and salinity are zero and all the other elements have large values. The matrix 𝗙S penalizes the initial atmospheric and deep ocean anomalies and forces the optimal initial conditions to be confined to the surface of the ocean.
We find that the maximum amplification of the THC occurs for τ = 121 yr. Because of the singularity of the norm kernel 𝗫, the initial contributions of the temperature and salinity to the THC cancel each other, resulting in a vanishing initial THC anomaly and therefore in an infinite amplification (ZT; Tziperman and Ioannou 2002). The optimal initial conditions P0 leading to this amplification are characterized by a nearly constant surface salinity anomaly as a function of latitude with a slightly stronger amplitude at higher latitudes (Fig. 1, dashed line) and a surface temperature anomaly that switches sign around 30° (Fig. 1, solid line). Moreover, we notice that the gradients of surface temperature and salinity are of the same sign (positive between −90° and −45° and negative elsewhere). The initial contributions of atmospheric temperature, moisture, and deep temperature and salinity vanish as required by the problem formulation.
The meridional gradients of the optimal initial temperature and salinity driving the THC are different from those found in ZT even though their meridional structures may seem similar. Specifically, we found that the gradients of density created by the surface optimals are smaller than those found in ZT. The contributions of surface temperature and salinity to the density gradients initially cancel each other, leading to the zero initial THC anomaly. This is unlike the case in ZT where the initial THC vanishes because of the contribution of both the surface and the deep anomalies that were equally important in the initial cancellation of the temperature and salinity effects on the THC.
Atmospheric and surface ocean temperature anomalies appear to have similar time evolutions (not shown) because of their interaction via air–sea heat fluxes. The salinity contribution to the surface density dominates that of the temperature at all times. Even though the initial deep temperature and salinity anomalies vanish, they are found to be later amplified by the nonnormal growth mechanism. Last, the THC anomaly grows and oscillates in the Southern Hemisphere (SH) while exhibiting a monotonic large-scale growth in the Northern Hemisphere (NH).
In addition to the amplification of the THC, the root-mean-square (rms) of the surface temperature is amplified by a factor of 6 while the rms of the surface salinity is reduced by a factor of 2 at the time of maximum amplification (left panels of Fig. 3, described below). These results are different from those of ZT in which the initial conditions were not restricted to the ocean surface. The amplification time scale of about 120 yr for the present work is longer than the 40-yr growth time scale found in ZT. The optimal initial temperature and salinity in the deep layer contributed to the excitation of relatively faster decaying modes and had an active role in setting the growth time scale and mechanism. The growth time scale of 121 yr found here, which is ⅓ as fast as that found when deep perturbations were allowed as well, is due to the interaction of relatively fast decaying eigenmodes with decay times up to 87 yr with less damped modes with decay times of 173 yr and longer (up to 780 yr). In addition, the rms amplification of the other model variables is not as dramatic as the corresponding 650 and 20 amplification factors found for the SST and SSS when deep initial anomalies are allowed. However, the nonnormal dynamics in the two-layer ocean model used here consistently lead to far stronger amplification of temperature and salinity than obtained with a three-box model (Tziperman and Ioannou 2002).
To identify the physical mechanism resulting in the transient amplification of the THC and other model variables, we ran a series of model experiments (using both the linearized and the full nonlinear model) eliminating one physical process at a time. The growth mechanism appears to be more complex than the one found in ZT, with distinct mechanisms for the SH and NH. We found that the advection of the mean temperature by the anomalous velocity, υ′∇T, plays a dominant role in the decay mechanism. The advection of the mean salinity by the anomalous velocity (υ′∇S) has two different impacts on the growth. The first impact is to reduce the amplitude of the growth by changing the surface salinity in the NH and the deep salinity in the SH. The second impact of υ′∇S is on the growth time scale acting with opposite effects in the SH and NH. The advection of the temperature anomaly by the mean flow, ∇T ′, speeds up and reduces the NH amplification. Eliminating the advection of the salinity anomaly by the mean velocity, ∇S′, results in a reduction of the magnitude of the growth in both hemisphere. As in ZT, the advection of the mean temperature and salinity by the anomalous velocity seems to be a major player in the transient amplification. However, the advection of the temperature and salinity anomalies by the mean flow seems here to participate in the transient growth and not only in the decay mechanism as in ZT.
In addition to the norm kernel 𝗫 (in 4) representing the sum of squares of the THC anomaly, the problem of transient amplification was explored under two additional norm kernels. The first norm is an energy norm kernel reflecting the normalized sum of the square of all the model variables (equivalent to the ocean available potential energy), and the second reflecting the linearized meridional ocean heat transport over both hemispheres. The optimal initial conditions and associated transient amplification mechanisms were found to be similar to the results obtained using the norm kernel representing the sum of squares of the THC. Thus our results do not appear to be sensitive to the norm kernel used.
b. L∞ norm
In the previous section, the infinite amplification found solving the generalized eigenproblem was a result of the singularity of the norm kernel (ZT; Tziperman and Ioannou 2002). This motivated us to propose and develop an alternative approach to transient growth based on the L∞ norm (as opposed to using different norm kernels within the L2 norm as discussed at the end of section 3a). The proposed new formalism is then used to obtain surface optimals leading to the growth of the maximum of the THC amplitude.
The infinity norm of a vector is simply the value of its largest absolute value component. The infinity norm of a given n × n matrix 𝗔 is defined as the maximum infinity norm of the product of 𝗔 with any vector x whose infinity norm is unity (Burden and Faires 2005, p. 425). The maximum norm of 𝗔 can then be shown to be the largest row sum of 𝗔,
where aij are the elements of 𝗔.
Let 𝗥THC be a matrix of dimensions Ny × Ny whose product with the state vector is equal to the THC at all latitudes [and whose columns are the vectors R(y) from (2)], and let 𝗟 be an Ny × Ny diagonal matrix with 0 everywhere along the diagonal except for elements corresponding to the surface temperature and salinity, which are equal to 1. The matrix 𝗟 is a projection matrix used to restrict the optimals to the surface ocean only. Consider a vector U0 defined as the initial THC anomaly at all latitudes induced by surface temperature and salinity only; it may be written as
where we define 𝗤 = 𝗥THC𝗟.
By successively multiplying (6) by 𝗤−1, 𝗕(τ), and finally 𝗥THC, we obtain
Both matrices 𝗟 and 𝗥THC are singular, leading to 𝗤 being singular as well, and therefore not invertible. As a result, the matrix 𝗤 is not full rank and possesses a nonzero null-space, leading, as we shall see below, to a possibly infinite number of solutions for the optimal initial conditions. To address this issue, we choose to leave out the null space and replace 𝗤−1 in all our calculation by the singular value decomposition–based pseudoinverse (Noble and Daniel 1988). We will return to this point later. The right-hand side of (7) is clearly identifiable as a vector containing the THC at all latitudes at a time τ, and we define the matrix 𝗚(τ) ≡ 𝗥THC𝗕(τ)𝗤−1 such that 𝗚(τ)U0 is again the THC at all latitudes at time τ. The L∞ norm of 𝗚(τ) is, by definition (Burden and Faires 2005, p. 425),
Therefore, calculating the infinity norm of 𝗚(τ) by looking for the vector U0 that maximizes the above expression is equivalent to looking for the optimal initial conditions U0 of the THC leading to the maximum THC amplitude at a time τ. Based on the expression for the L∞ norm of 𝗚(τ) (Burden and Faires 2005), the optimal initial conditions U0 of unit norm are given by
where the index i denotes the index of the row of maximum sum of the matrix 𝗚(τ), and Gij are the elements of the matrix 𝗚(τ). As we mentioned previously, we choose to leave out the null space of (8) by replacing the inverse of the matrix 𝗤 by its pseudoinverse. In physical terms, this corresponds to finding the simplest initial temperature and salinity anomalies that are consistent with the optimal initial THC anomaly. These anomalies are the simplest in the sense that adding initial temperature and salinity anomalies (which are in the null space of the above matrix) to the optimal initial conditions may increase the magnitude of the temperature and salinity anomalies without affecting the initial THC anomaly.
To recover the initial perturbations of surface temperature and salinity represented by P0 from the optimal initial conditions for the THC, U0, we use the inverse of the previous transformation (6):
The maximum amplification of the maximum THC anomaly is found to occur for τ = 9 yr, and the amplification factor is 1.6 (Fig. 2a). The optimal initial conditions of surface temperature and salinity shown in Fig. 2b are linear functions of latitude, and their gradient drives initial THC anomalies whose values are 1 everywhere [Fig. 3; (9)], consistent with our problem formulation requiring the maximum initial THC anomaly amplitude to be 1. The initial surface salinity creates a positive density gradient that dominates the initial THC anomaly.
When the model is started from the optimal initial conditions, the relatively small initial surface temperature anomaly is rapidly damped out by air–sea interactions. However, the initial salinity perturbations at high latitudes in the NH creates a positive THC anomaly, enhancing the mean THC there. During the first several years, the NH therefore becomes saltier and warmer because of the advection of the mean temperature and salinity gradients (i.e., of warm and high-salinity water from midlatitudes) by the anomalous velocity. This is precisely the mechanism of the salinity advective THC instability (Marotzke et al. 1988). In the SH, the positive THC anomaly driven by the salinity anomalies reduces the mean THC, leading to the freshening due to reduced salt advection, and the same process as in the NH occurs there with the opposite sign. Eventually the perturbations decay in the NH because of the warmer water advected from the equator overcoming the atmospheric cooling there. The warm water reduces the density signal at high latitudes and causes the positive perturbation of THC to decay. A similar process occurs in the SH, again with an opposite sign.
c. A discussion of the difference between the two norms
Using the two norms to calculate the optimal initial conditions corresponds to asking two different but complementary physical questions. With the L2 norm, we look for the initial conditions leading to maximum amplification of the sum of squares of the THC at all locations. Because this leads to a singular mathematical problem, the solution is characterized by a vanishing initial THC. We therefore find initial conditions for the temperature and salinity corresponding to a vanishing initial THC anomaly and leading to the efficient development of such a THC anomaly. The vanishing initial THC implies that the THC amplification is infinite, yet all model variables remain physical and finite.
With the maximum norm (L∞), we are looking for an initial THC anomaly with a prescribed maximum initial value that then leads to an optimal growth at later times. In this case, we still encounter singular matrices yet, by using their generalized inverse, the initial THC anomaly has the required unit amplitude and the amplification is therefore finite. The advantage of this norm is that it is closely related to THC indices typically used in GCMs and based on the maximum THC value obtained, for example, anywhere within the North Atlantic (e.g., Delworth et al. 1993).
The physical amplification mechanisms obtained using the two norms are also different. The maximum norm leads to a nonvanishing initial THC so that the classic advective feedback, relying on the advection of the mean salinity by the perturbation flow, is therefore active from the initial time of the perturbation. Because the L2 norm leads to a vanishing initial THC anomaly, the growth mechanism needs to rely on other physical factors as analyzed in detail in section 3a.
4. Stochastic optimals of the THC
We find here that the spatial structure of the stochastic forcing plays an important role in exciting the THC variability on different time scales. Consider the linearized model in (1) when stochastically forced:
where ξ(t) is the noise vector and the diagonal matrix 𝗟 is defined such that the noise appears only in the equations for the ocean surface as defined in section 3b. The stochastic forcing can be viewed as being due to freshwater and air–sea heat fluxes. Following Hasselmann (1976) and Frankignoul and Hasselmann (1977), we assume the forcing to be a white-noise process with zero mean and uncorrelated in time, and the statistics of the noise to be separable in time and space:
where angle brackets denote an ensemble average, δmn is the Kronecker delta function, and 𝗦 is the spatial covariance matrix of the noise.
Using a norm kernel 𝗫 reflecting the sum of squares of the THC anomaly over the entire domain, the variance of the THC is then given by (Farrell and Ioannou 1996; Kleeman and Moore 1997; Tziperman and Ioannou 2002)
where Tr denotes the trace of the matrix, the propagator matrix from time s to time τ is defined by 𝗕(τ, s) = e𝗠(τ−s) and the forcing correlation matrix 𝗭 is
The eigenvectors zi of the matrix 𝗭 corresponding to the largest eigenvalues of 𝗭 are the stochastic optimals of the stochastically forced system (Farrell and Ioannou 1996; Kleeman and Moore 1997). These are the spatial patterns of the stochastic forcing that are the most efficient in maximizing and sustaining the variance of the THC anomalies.
The stochastic optimals are found by integrating the forced dynamical system over all excitation times [hence the integral in (14)], whereas the optimal initial conditions are found for a given single response time [(4)]. The eigenvalues of 𝗭, αi, represent the fraction of the variance excited by the corresponding stochastic optimal zi.
The fraction of the variance, αi/Tr(𝗭𝗖), excited by each of the corresponding stochastic optimal zi is shown in Fig. 4a, sorted by increasing order, where 𝗖 is the correlation matrix of the response whose elements are Cij = 〈PiP*j〉 (e.g., Farrell and Ioannou 1993). The structure of the 61 stochastic optimals contributing to most of the variance is associated with the surface salinity, indicating that freshwater forcing is more effective than heat fluxes in driving THC variability. The first and second stochastic optimals explaining, respectively, 22% and 4% of the stochastically induced variance are shown in Figs. 4b and 4c. Both of the first two stochastic optimals exhibit a large-scale spatial structure, although the first one is nonzero throughout the basin, whereas the second is concentrated mostly at high latitudes.
The spatial structure of the first stochastic optimal (for both the surface temperature and surface salinity) shows a large resemblance with the least damped mode of the adjoint operator. In addition, the surface salinity structure seems to be similar to the optimal initial conditions found in the previous section (shown in Fig. 1). Figures 5a,c show time series of the THC at 72°N from two model runs forced with the structure of the first and second stochastic optimals, correspondingly. The associated power spectra based on both the model output and from the theoretical expression (Farrell and Ioannou 1996) are shown in Figs. 5b,d. The power spectrum is mostly red, and no significant peaks are found in either case even though oscillatory modes are present in the linearized model.
The variability excited by the first stochastic optimal has significantly longer time scales, because it excites modes with longer decay times than the second stochastic optimal. The structure of the response (i.e., the EOFs) of the linearized model when forced with the first stochastic optimal is found to have the same spatial structure as the least damped mode of the model (decay time of 734 yr and period of oscillation of 1008 yr). The model response when forced with the second stochastic optimal resembles an eigenmode with a much faster decay time of 87 yr and an oscillation period of 65 yr. In either case, the spatial structure of the response of the model (EOF) is radically different from the stochastic optimal used to force the model, a sign of efficient nonnormal dynamics (Farrell and Ioannou 1999).
Interestingly, even though the two stochastic optimals are found to be large-scale, they excite very different ocean responses. The simple model thus predicts that high-latitude stochastic forcing may excite higher-frequency variability than stochastic forcing with nonvanishing signal throughout the domain. The more general, and perhaps more robust, lesson from this calculation of the stochastic optimals is that somewhat different large-scale stochastic forcing patterns may result in very different temporal response of the THC. Hemispheric large-scale patterns similar to the two first stochastic optimals may possibly be related to actual modes of atmospheric or climate variability, such as annular modes, that have been shown to excite THC variability on interannual to interdecadal time scales (e.g., Marsh 2000; Selten et al. 1999; Timmermann et al.1998). It would be interesting to study this prediction in a more realistic GCM.
We analyzed the optimal excitation of thermohaline circulation anomalies and variability by surface temperature and salinity forcing in a coupled ocean–atmosphere box model. Because of the nonnormality of the model, a large transient amplification was possible even though the model is linearly stable and any initial perturbation eventually decays.
In the case of a fully deterministic linearized dynamics, we analyzed the transient growth using two different norms. First, we used the standard sum-of-squares (L2) norm, and then we introduced a novel approach based on characterizing the THC amplitude by its maximum value throughout the basin as is often done in GCM studies, and showed that it is equivalent to using the maximum (L∞) norm. These two approaches provide complementary views of the problem, and allowed us to gain additional insight into the dynamics of nonnormal THC excitation by surface forcing.
Using the standard L2 norm, we found optimal initial conditions of surface temperature and salinity leading to a transient growth of the THC anomalies with a time scale of about 120 yr. The initial THC anomaly vanishes, then grows because of the effects of the surface salinity anomalies, and finally anomalies rapidly develop in the deep ocean as well, providing a further driving of the THC anomalies. The growth mechanism involves the advection of the mean temperature and salinity by the anomalous velocity as well as the advection of the anomaly fields by the mean velocity. This mechanism is different from that involved when the forcing is not limited to the surface and is also distinct from the advective instability mechanism due to salinity anomalies (e.g., Marotzke et al. 1988; Tziperman et al. 1994). The maximum norm, on the other hand, resulted in a nonvanishing initial THC anomaly and in a different amplification mechanism similar to the classic advective THC instability feedback.
We find that limiting the initial anomalies to the surface leads to a less efficient amplification than the excitation of the THC by both surface and deep ocean temperature and salinity anomalies, as studied by ZT. However, surface excitation by atmospheric forcing is more likely than deep excitation (e.g., by mesoscale eddies), making the mechanism studied here more relevant to the excitation of THC anomalies in the ocean.
In the case of a stochastically forced THC, we found that salinity anomalies were the most efficient structures at sustaining THC variability. Our perhaps most interesting result regarding the stochastic surface forcing is that two different spatial patterns of stochastic forcing that are both large scale and smooth can lead to very different temporal THC response. More specifically, we find that the first stochastic optimal, which has a basin-scale nonzero structure, excites very slow variability, while the second optimal, also having a large-scale smooth structure but concentrated mostly at high latitudes, excites a much faster THC variability. The spectrum was mostly red in both cases without any spectral peak even though some of the eigenmodes of the linearized model are oscillatory. This differs from the results of Griffies and Tziperman (1995) where the atmospheric stochastic forcing excites an oscillatory mode of the linearized model and results in a corresponding spectral peak.
Our results suggest that nonnormal THC dynamics may possibly have a nonnegligible impact on present-day climate variability. It would be interesting to examine the predictions of the simple model used here in a more realistic GCM where the optimal forcing patterns can be compared to observed patterns of atmospheric variability.
The authors are grateful to Andrew M. Moore and Joseph Bernstein for helpful discussions and to two anonymous reviewers for their useful comments on the manuscript. This work was supported by the McDonnell Foundation and by the NSF Paleoclimate Program, Grant ATM-0502482. Author ET thanks the Weizmann Institute, at which some of this work was carried out during a sabbatical.
Corresponding author address: Laure Zanna, Dept. of Earth and Planetary Sciences, and School of Engineering and Applied Sciences, Harvard University, Cambridge, MA 02138. Email: email@example.com