The authors show that a simple three-dimensional ocean model linearized about a resting basic state can accurately simulate the dynamical ocean response to wind forcing by the Madden–Julian oscillation (MJO). This includes the propagation of equatorial waves in the Indian Ocean, from the generation of oceanic equatorial Kelvin waves to the arrival of downwelling oceanic equatorial Rossby waves in the western Indian Ocean, where they have been shown to trigger MJO convective activity. Simulations with idealized wind forcing suggest that the latitudinal width of this forcing plays a crucial role in determining the potential for such feedbacks. Forcing the model with composite MJO winds accurately captures the global ocean response, demonstrating that the observed ocean dynamical response to the MJO can be interpreted as a linear response to surface wind forcing. The model is then applied to study “primary” Madden–Julian events, which are not immediately preceded by any MJO activity or by any apparent atmospheric triggers, but have been shown to coincide with the arrival of downwelling oceanic equatorial Rossby waves. Case study simulations show how this oceanic equatorial Rossby wave activity is partly forced by reflection of an oceanic equatorial Kelvin wave triggered by a westerly wind burst 140 days previously, and partly directly forced by easterly wind stress anomalies around 40 days prior to the event. This suggests predictability for primary Madden–Julian events on times scales of up to five months, following the reemergence of oceanic anomalies forced by winds almost half a year earlier.
The dynamics of the equatorial oceans are dominated by equatorial Kelvin and Rossby waves on time scales from intraseasonal (Kessler et al. 1995; Hendon et al. 1998) to interannual (Battisti 1988; McPhaden 1999). These are planetary-scale internal waves with wavelengths of thousands of kilometers that affect sea surface height (SSH) and pycnocline depth. Variations in pycnocline depth are several orders of magnitude larger than the SSH anomalies and have the opposite sign; thus positive SSH anomalies correspond to a deeper mixed layer depth (MLD) and downwelling anomalies. This variability modifies the upper-ocean heat content and the mixing of cold subsurface waters, leading to changes in sea surface temperature (SST) (McCreary 1983; Battisti 1988). These SST anomalies subsequently modulate atmospheric convection; such processes have been shown to be important for the El Niño–Southern Oscillation (ENSO) (Battisti 1988; Kessler and McPhaden 1995; McPhaden 1999) and the Madden–Julian oscillation (MJO) (Han et al. 2001; Webber et al. 2010).
Equatorial Kelvin waves exhibit a single equatorial SSH maximum, while equatorial Rossby waves have more complex meridional structures, the simplest and most commonly observed being the first meridional mode with symmetric off-equatorial SSH maxima (Chelton et al. 2003). A westerly wind burst on the equator will generate Ekman convergence, thus forcing a positive SSH anomaly on the equator, along with downwelling anomalies and an eastward current. These anomalies propagate along the equator as a downwelling Kelvin wave (Giese and Harrison 1990). Westerly wind bursts will typically also trigger upwelling Rossby waves due to off-equatorial Ekman divergence, although this is somewhat dependent on the meridional structure of the wind stress (Chelton et al. 2003).
For a continuously stratified ocean there are an infinite number of vertical (baroclinic) modes for both equatorial waves, in addition to the barotropic mode (sometimes referred to as the zeroth baroclinic mode). However, typically only the barotropic and first few baroclinic modes are observed to be important (Giese and Harrison 1990; Chelton and Schlax 1996; Wunsch 1997). The propagation speed of such waves depends upon both the baroclinic mode and the stratification of the fluid through which they propagate. Equatorial Kelvin waves are nondispersive with the first baroclinic mode (n = 1) phase speed ce between around 2.4 and 3.0 m s−1, depending on the local stratification (Chelton et al. 1998). The first baroclinic mode equatorial Rossby wave phase speed, c, also depends on the meridional mode (Chelton et al. 2003), and follows the dispersion relation:
where β ≈ 2.3 × 10−11 m−1 s−1 is the meridional gradient of planetary vorticity, k is the zonal wavenumber, and m is the meridional mode number. Thus, long first baroclinic, first meridional mode Rossby waves propagate westward at approximately or around 0.8–1.0 m s−1 in the Indian Ocean.
Kelvin waves can also propagate along coastal waveguides, with the coastline to the right (left) of the wave in the Northern (Southern) Hemisphere. Coastal Kelvin waves can be generated by local wind stress forcing or by equatorial waves incident on a meridionally oriented coastal boundary. However, the effect of sloping coastal bathymetry is to produce a coastally trapped wave that combines the properties of Kelvin waves and barotropic shelf waves, with modified phase speed and increased dispersion (Huthnance 1975; Brink 1982, 1991).
The Madden–Julian oscillation is the dominant atmospheric mode of intraseasonal variability in the tropics (Madden and Julian 1971, 1972). It is associated with convective rainfall and wind variability with a periodicity of around 30–60 days, although it is only quasiperiodic and thus has a broadband spectral signal extending as far as 100 days (Salby and Hendon 1994). See Zhang (2005) and Lau and Waliser (2005) for a full review of the MJO and its impacts.
The MJO is an important source of surface wind stress variability and thus generates substantial oceanic equatorial wave activity (Hendon et al. 1998; Han et al. 2001; Webber et al. 2010). Observations have linked the MJO-generated oceanic equatorial Kelvin waves to the triggering of ENSO events (McPhaden 1999). There is also growing evidence of the potential for oceanic equatorial Rossby waves in the Indian Ocean to force variability in the MJO (Webber et al. 2010, 2012) and for feedbacks between the two (Han et al. 2001; Han 2005; Fu 2007).
The MJO is sporadic, with periods of cyclical activity (“successive” events) interspersed with intervals of relative quiescence. “Primary” Madden–Julian (MJ) events are those that are not preceded by any coherent MJO activity. Here, we follow Matthews (2008) by examining the subset of primary events that begin in the western Indian Ocean, which accounts for the largest portion of such events. Matthews found no triggers or precursors from the atmospheric dynamics or thermodynamic interaction with the ocean for such events. However, Webber et al. (2012) showed, through a combination of case studies and composite analysis, that these events do coincide with the arrival of a downwelling oceanic equatorial Rossby wave in the western Indian Ocean, implying that such waves could act as a trigger. They suggested that the associated warm SST anomalies act to destabilize the atmospheric boundary layer, thus priming the atmosphere for the large-scale convection associated with a MJ event. The western Indian Ocean appears to be a favorable region for such forcing to occur, possibly because of the relatively shallow thermocline (Xie et al. 2002), higher rates of entrainment cooling (McCreary et al. 1993), and the relatively weak intraseasonal surface heat fluxes here.
This paper examines the dynamical response of the Indian Ocean to the MJO and its potential to couple to the MJO through simulations in a linearized ocean circulation model. A combination of model runs forced by both idealized and realistic surface wind stress forcing are used. The pathways and time scales of the wave propagation are examined using simple idealized wind patches and the effects of stratification, bathymetry, and forcing scale are analyzed. Model runs with composite forcing are compared with the observational results of Webber et al. (2010) to validate the model performance and further investigate the global dynamical ocean response. The model is then run with observed winds from the period prior to a case study primary event on 24 September 2004. Comparison with observations from that case study allows the key atmospheric forcing events to be identified and their relative importance to be tested. Thus, we identify the time scales that are important for this primary MJ event triggering mechanism.
2. Observational data and model forcing
The surface wind stress data used in this study are calculated using daily global 10-m wind velocities v from the European Centre for Medium-Range Weather Forecasts (ECMWF) Era-Interim dataset. For compositing relative to the MJO, the period of 1 January 1990 to 31 December 2008 is used. The wind stress τ is then approximated using the standard bulk formula τ = ρCD|v|v, with the air density ρ = 1.23 kg m−3 and the drag coefficient CD = 1.5 × 10−3. This is calculated separately for the zonal and meridional components of the wind field and then spatially interpolated using cubic splines onto the ⅓° grid used in the ocean model.
The SSH anomalies used in this study are from the merged Ocean Topography Experiment (TOPEX)/Poseidon–Earth Remote Sensing (T/P–ERS) satellite altimetry product (Fu et al. 1994; Ducet et al. 2000; Le Traon et al. 2001). Weekly data on a 0.25° grid were obtained for the period from 14 October 1992 to 23 January 2008 and then interpolated to daily values using cubic splines. SST data from the Tropical Rainfall Measuring Mission Microwave Imaging (TMI) satellite (Kummerow et al. 2000; Gentemann et al. 2004) were extracted from 1 January 1998 to 17 December 2008. The data are available as daily fields of the 3-day running-mean SST at 0.25° resolution. Outgoing longwave radiation (OLR) is used as a proxy for tropical deep convective precipitation. The OLR data for this study were obtained as 2.5° resolution gridded daily output from the optimally interpolated Liebmann and Smith (1996) dataset. For comparison with the subsurface variability in the model, we use density data from the Estimating the Circulation and Climate of the Ocean-Global Ocean Data Assimilation Experiment (ECCO-GODAE) ocean state estimate (Wunsch and Heimbach 2007), which we treat as “observations” in the same sense as an atmospheric reanalysis product.
To construct composites of observational fields for surface forcing or comparison with model output, we use the eight phases of the MJO as defined by the Wheeler–Hendon index (Wheeler and Hendon 2004). Phase 1 corresponds to minimum convection over the Maritime Continent and the initiation of active convective anomalies in the western Indian Ocean. These anomalies move eastward in each successive phase; by phase 4 the positive convective anomalies overlie the Maritime Continent while the suppressed convection has propagated into the western Pacific. A new region of suppressed convection originates in the western Indian Ocean at phase 5; phases 5–8 are thus equivalent to phases 1–4 with anomalies of the opposite sign.
The linear trend, mean, and first three harmonics of the annual cycle were removed pointwise from each dataset to produce detrended anomaly fields, which were then bandpass filtered using a 20–200-day Lanczos filter. The use of such a broadband filter is motivated by its ability to exclude low-frequency climate signals as well as high-frequency noise, while better retaining the MJO signal and the distinction between individual MJ events than a narrower filter window (Matthews 2000). The filter uses 241 symmetric weights, meaning that 120 days of data are lost at each end of the dataset. All datasets were further truncated so that an integer number of calendar years are retained, thus avoiding seasonal bias. To investigate the oceanic response to the MJO, composites are created separately for each MJO phase, using all days within that phase when the amplitude of the Wheeler–Hendon index exceeds an arbitrary threshold value of 1.5, along with counterclockwise rotation in the Wheeler–Hendon phase space to ensure eastward propagation. The threshold value was chosen to select the stronger MJO events (whose impact on the ocean should be stronger and more coherent), as opposed to the threshold value of 1 chosen by Wheeler and Hendon (2004).
3. Model description and setup
The ocean model used for this study is a linearized version of the general circulation model described by Webb (1996). A similar approach was used by Anderson et al. (1979) to study the transient response in the North Atlantic to a change in wind forcing. There are no mean currents and the mean vertical potential density stratification, , is independent of horizontal position. The equations of motion are
The variables φ, λ, z, t, u, υ, w, p, ρ, represent latitude, longitude, depth (negative of), time, zonal velocity, meridional velocity, vertical velocity, pressure, and density anomaly, respectively. The radius of the earth is a, g is the acceleration due to gravity, and ρ0 is a reference potential density. The Coriolis parameter f = 2Ω sinφ where Ω is the speed of angular rotation of the earth. The horizontal eddy viscosity is Ah = 103 m2 s−1 and the vertical eddy viscosity, Km, and diffusivity, Kh, are both set to 10−4 m2 s−1. At the ocean surface a wind stress is applied and related to the ocean velocities by
and no flux of density is prescribed.
The model is configured quasi-globally (the Arctic and the Nordic Seas north of 67°N are omitted), with a horizontal resolution of ⅓° × ⅓° and 40 uneven levels in the vertical that vary smoothly from 10-m thickness at the sea surface to 300-m thickness at depth. The mean vertical potential density profile, , is taken from the Indian Ocean, calculated from the World Ocean Database 2001 (Levitus et al. 2002) over 10°N–10°S, 60°–90°E (Fig. 1). The model can also be configured with only two density layers, where the upper 12 model layers are given the surface potential density (1022.1 kg m−3) and the remaining 28 layers the deepest layer potential density (1027.8 kg m−3). The bathymetry of the model is shown in Fig. 2 along with key geographical regions for this study and the pathway for the Indonesian Throughflow (ITF).
4. Idealized westerly wind burst forcing experiments
a. Control run
Several studies have shown the importance of westerly wind bursts (WWBs) in forcing equatorial oceanic Kelvin waves (Giese and Harrison 1990; McPhaden et al. 1992; Kessler et al. 1995). The MJO is a major source of relatively long-lived WWBs at the same spatial scale as the resulting Kelvin waves, and is thus an effective forcing mechanism (Hendon et al. 1998). Here we idealize this forcing mechanism by applying a Gaussian WWB centered on the equator (φ0 = 0) at λ0 = 75°E, with a decay scale of 10° in both longitude (Wλ) and latitude (Wφ). The wind forcing is applied constantly with a peak westerly magnitude of τ0 = 0.1 N m−2 (corresponding to a wind speed of 7.4 m s−1) while t ≤ 5 days and is zero thereafter. The surface zonal wind stress is thus given by
The spatial and temporal scale of this forcing is broadly characteristic of the MJO, as is the peak magnitude (Hendon et al. 1998). This experiment is qualitatively similar to that of McCreary et al. (1993) and Valsala (2008) who used a 2½-layer model to examine the response of the Indian Ocean to seasonal to interannual wind anomalies.
Figure 3a shows the model SSH anomalies for this control simulation after the initial five days of forcing. The positive SSH anomalies indicate a downwelling equatorial Kelvin wave being forced at 75°E, while the negative off-equatorial SSH anomalies indicate the upwelling Rossby wave response to the wind forcing. By day 35 the Kelvin wave has propagated eastward, impinging on the coast of Sumatra where coastal waves are triggered along with reflected downwelling Rossby waves (Fig. 3b). The coastally trapped waves propagate both northward and southward, with some of the energy from the latter portion propagating into the Maritime Continent. This energy spreads in all directions, with coastal wave signals propagating around the island coastlines (e.g., the positive SSH anomalies along the east coast of Java and Sumatra and then along the west coast of Borneo) and through the various channels. The surface wave signal is not apparently dependent on channel depth and is capable of propagating through channels with a maximum depth of less than 20 m, such as between Borneo and Sumatra. Note that the vertical resolution is around 10 m in the surface layers of the model, so such waves will be imperfectly resolved in such shallow depths.
The subsequent westward propagation of the reflected Rossby wave along the off-equatorial waveguide can be seen in the transition from day 35 to 65 (Figs. 3b,c) and is similar to that observed (Oliver and Thompson 2010; Webber et al. 2010). At the same time, the upwelling Rossby waves (negative SSH anomalies) have reached the coast of Africa and have begun to reflect into an upwelling equatorial Kelvin wave. The coastally trapped waves propagating along the northern coastal waveguide have circumnavigated the Bay of Bengal and the Arabian Sea, with the westernmost positive SSH anomalies reaching the southern coast of the Arabian Peninsula, as seen in models (Valsala 2008) and observations (Oliver and Thompson 2010; Webber et al. 2010). The downwelling coastal wave that propagated southward along the west coast of Australia has generated westward-propagating Rossby waves in the southern Indian Ocean, while the anomalies that propagated through the Maritime Continent have formed a downwelling equatorial Kelvin wave in the Pacific (positive SSH anomalies at 140°–160°E on the equator).
Figure 3d shows the eventual SSH anomalies at 95 days along with the propagation paths of the various wave responses. The thick black arrow shows the equatorial Kelvin waveguide and subsequently the northern coastal waveguide, emphasized because this waveguide retains energy within the tropical Indian Ocean. The thin black lines show the southern coastal waveguide including the path through the Lombok Straits into the Maritime Continent. The gray arrows and dashed boxes show the off-equatorial Rossby waveguide. The arrival of the downwelling equatorial Rossby wave in the western Indian Ocean can be seen in this figure (positive off-equatorial SSH anomalies at 45°–60°E). Coastal Kelvin waves have propagated all the way around the northern coastal waveguide (shown by the thick black arrow) to arrive in the western Indian Ocean at the same time as the equatorial Rossby wave. Once they reach the equatorial region they may contribute to and reinforce the positive SSH anomalies there at around 100 days after the initial wind forcing; we return to this question later. The propagation and timing of both the equatorial Rossby and coastally trapped waves agrees well with observational composites of the free wave response to MJO forcing (Webber et al. 2010).
b. Sensitivity to latitudinal width of westerly wind forcing
The sensitivity of the ocean response to the latitudinal Gaussian decay scale (Wφ) of the wind patch is examined. Four experiments were carried out with latitudinal decay scales of 20°, 10°, (the control run), 5°, and 1° (Fig. 4). The peak amplitude τ0 is the same for all four simulations (0.1 N m−2); hence the wider wind bursts will input more total energy to the ocean, as well as projecting differently onto the equatorial wave meridional modes. As would be expected, there are larger SSH signals away from the equator in the wider wind burst simulations, leading to stronger variability in the Bay of Bengal and the Arabian Sea. However, in the equatorial waveguide the sign of the SSH anomalies changes between the simulations, with positive anomalies in the western Indian Ocean in the 20° experiment but predominantly negative SSH anomalies in the same region for the 1° and 5° simulations. The key result is therefore that the latitudinal decay scale of the westerly wind burst can alter not just the magnitude, but also the sign of the resulting equatorial anomalies 95 days later.
The relative strength of the various baroclinic modes (with their different vertical structures and zonal propagation speeds) is also different between the simulations. Within the equatorial waveguide, there are strong positive anomalies at 40°–70° and 80°–90°E clearly separated zonally by negative anomalies at 70°–80°E in the 20° experiment (Fig. 4a), while narrower wind bursts lead to much less distinction and separation between the modes (zonally uniform negative near-equatorial anomalies in the 5° and 1° experiments). It therefore appears that varying the latitudinal decay scale of the wind forcing leads to a different combination of baroclinic and meridional modes in the model and subsequently a different equatorial SSH signal 95 days after the initial forcing. The implication is that a westerly wind burst at the equator needs to have a latitudinal decay scale greater than 5° latitude in order for the ocean forcing mechanism postulated by Webber et al. (2010) to be important.
c. Sensitivity to stratification and bathymetry
We now construct simulations with the same wind forcing as the control run but with varying bathymetry and stratification to examine the interaction between Kelvin and shelf wave modes along the coastal waveguide. The bathymetry is either the default model bathymetry used in the control run (see Fig. 1) or a “flat” bathymetry, set to zero where it is shallower than 200 m and 4000 m elsewhere. The stratification is either the same as the control, or a two-layer stratification comprising the surface potential density (ρ1 = 1022.1 kg m−3) in the top 12 model levels (~200 m) and the bottom layer potential density (ρ2 = 1027.8 kg m−3) in the lower 28 model levels (Fig. 1). For the two-layer configuration, we expect the ocean dynamics to consist solely of the barotropic and first baroclinic mode waves, with coastal waves resembling pure Kelvin waves. The propagation speed of the first baroclinic mode Kelvin wave for this configuration is
where He is the equivalent depth, that is, the depth of the surface density layer. This phase speed is in good agreement with the calculated first baroclinic mode equatorial Kelvin wave phase speed for the Indian Ocean, which is typically between 2.6 and 2.8 m s−1 (Chelton et al. 1998). The estimated Kelvin wave phase speed for the model with full 40-layer stratification and the realistic (control run) bathymetry (Fig. 2) is around 2.5 m s−1, so we expect the waves to propagate slightly faster in the two-layer configuration. To isolate the Kelvin wave response we construct Hovmöller diagrams along the equatorial Kelvin and coastal waveguides, as shown by the thick black arrow in Fig. 3d. Note that some of the reflected Rossby wave propagation can be seen in the equatorial waveguide, but this Rossby wave signal is relatively weak.
The first experiment is the most idealized, with just the two-layer stratification and flat bathymetry (Fig. 5a). The downwelling first baroclinic mode Kelvin wave (positive SSH anomalies) propagates eastward along the equatorial waveguide to the coast of Sumatra (point B) then anticlockwise around the Bay of Bengal to the southern tip of India (point C) and onward around the Arabian Sea to the Somali coast at the equator (point A). This propagation is very clear and coherent, with minimal energy loss. The arrival of the coastal Kelvin wave in the western Indian Ocean (point A) coincides almost perfectly with that of the reflected Rossby wave which can be seen to a certain extent in the equatorial waveguide (left-hand side of Fig. 5a). Interestingly, some of the energy initially associated with the upwelling equatorial Rossby waves escapes into the coastal waveguide near the tip of India (point C) where the coastal waveguide extends into the equatorial domain. There may also be direct forcing of coastal waves here, although the alongshore component of the westerly wind forcing will be small.
After the initial downwelling coastal Kelvin wave there are several wave signals of both signs that can be seen propagating along the coastal waveguide. Some of these can be traced back to reflection of equatorial Rossby waves at the western boundary (point A); for example, the negative SSH anomalies here at around day 40 reflect into an upwelling equatorial Kelvin wave, which reaches the eastern boundary (point B) around day 75. This signal then propagates coherently around the coastal waveguide as an upwelling equatorial Kelvin wave (negative SSH anomalies) to arrive just past point C at day 120. Similar behavior can be seen in the downwelling wave signal (positive SSH anomalies) that reaches the western boundary around day 100, subsequently reflecting into an equatorial Kelvin wave that arrives at point B at day 120.
The effect of bathymetry is investigated by replacing the flat bathymetry with realistic sloping bathymetry (cf. Fig. 2), while keeping the idealized two-layer stratification (Fig. 5b). This experiment shows less coherent propagation along the coastal waveguide, whereas the equatorial propagation pattern in the open ocean (left-hand portion of the diagrams) is almost identical. The sloping bathymetry will lead to coastal trapped waves with dispersive properties (Huthnance 1975; Brink 1982) and, thus, will attenuate the signal as it propagates around the coastal waveguide. Thus it is unsurprising that the downwelling Kelvin wave signal is relatively spread out and that the signal is weaker by the time it reaches the western Indian Ocean.
The effect of using a realistic 40-layer stratification, but reinstating the idealized flat bathymetry, is shown in Fig. 5c. The propagation of energy along the equatorial waveguide is much more spread out than in the idealized two-layer stratification experiments (Figs. 5a,b), probably due to the partitioning of wave energy between the multiple baroclinic modes that are now allowed by the model. The coastal wave propagation is clearly split into two baroclinic modes with different propagation speeds, as shown by their relative slopes on the Hovmöller diagrams. As expected, the wave propagation is slightly slower for this realistic stratification than for the two-layer configuration (Figs. 5a,b) with the peak in positive SSH at the western boundary (point A) occurring 10–20 days later.
Finally, Fig. 5d shows the corresponding Hovmöller diagram for the control simulation (as in Fig. 3) with both realistic stratification and bathymetry. In this figure, the two factors of energy partitioning between baroclinic modes and increased energy dispersion along the coastal waveguide lead to quite a different picture from Figs. 5a–c. The two baroclinic modes evident in Fig. 5c are no longer apparent. The propagation of the positive SSH anomalies along the equatorial waveguide is less coherent, compared with Fig. 5b, with weaker positive SSH anomalies at the western boundary (point A) at days 100–120. In addition, the coastal wave signal spreads farther in time and space than in the other simulations and becomes insubstantial before arriving in the western Indian Ocean; the weak positive SSH signal on the equator at point A is primarily due to the Rossby wave propagation. Therefore, in the real ocean, it seems likely that the coastal waveguide will be less important than the equatorial Rossby waveguide for forcing MJO variability in the western Indian Ocean.
The simulations described above summarize the dynamical ocean response to idealized WWBs associated with the MJO in the Indian Ocean. This consists of a combination of equatorial Kelvin and Rossby waves along with coastally trapped waves. Some of the equatorial wave energy escapes along the southern coastal waveguide and into the Maritime Continent, which subsequently leads to weak equatorial Kelvin wave activity in the Pacific. The remainder returns across the Indian Ocean in the form of both equatorial Rossby waves and coastally trapped waves propagating along the northern coastal waveguide. For realistic simulations, the interaction of wave modes along the sloping bathymetry of this coastal waveguide leads to dissipation of this coastal signal before it arrives in the western Indian Ocean. In contrast, equatorial Rossby waves do propagate coherently across the width of the Indian Ocean and arrive in the western Indian Ocean around 100 days after the initial downwelling Kelvin wave was forced. There are also upwelling wave signals (negative SSH anomalies) that are approximately half a cycle out of phase with the downwelling waves. Multiple meridional and baroclinic wave signals are seen, the distribution of which partly depends upon the latitudinal decay scale of the initial wind forcing. The combination of these wave modes has the potential to critically alter the dynamic ocean state along the equatorial waveguide, and is thus crucial for determining the potential for triggering atmospheric convection and the MJO.
5. Realistic composite MJO wind forcing
The key processes in the dynamical ocean response to MJO-like wind forcing have been determined in the idealized experiments of section 4. In this section, the dynamical ocean response to realistic MJO wind forcing is examined. We force the model with composite MJO winds globally (see section 2). These composites consist of eight global fields of the zonal and meridional components of the surface wind (one for each of the eight Wheeler–Hendon MJO phases), which are then interpolated onto a daily time scale by assuming a constant 6-day time interval between phases. This gives a 48-day cycle for the MJO, which is consistent with previous definitions of the MJO life cycle (Madden and Julian 1972, 1994). However, this emphasizes the cyclical component of the MJO, which may not always be appropriate during periods of sporadic behavior (Matthews 2008), and will tend to smear out some of the freely propagating oceanic wave signals generated (Webber et al. 2010). Nevertheless, it is useful as a tool to investigate the key components of the dynamical ocean response to the cyclical MJO. The simulation is run for two full cycles of the MJO (i.e., 96 days), to allow processes occurring over the time scale of multiple MJO cycles to be seen.
Global SSH anomalies
Figure 6 shows the SSH anomalies from the second MJO cycle of the simulation (i.e., day 54 to day 96) in the left column, with the equivalent composites of observed SSH anomalies in the right-hand column. In general the agreement between the model and observations is strikingly good. There is more small-scale noise in the observations, which is to be expected as the model is forced by the relatively low resolution 1° ECMWF winds (cf. the 0.25° SSH data) and does not contain currents and their associated eddies. The propagation of waves around the Indian Ocean in the model shows all of the same features as the observations, as does the propagation of SSH anomalies into the Maritime Continent. Several of the features evident in the idealized westerly wind burst experiment (Fig. 3) are also evident in this composite simulation; the reflection of equatorial waves at the eastern boundary is similar, as is the propagation of waves around the coastal waveguide.
The modeled propagation of Kelvin waves across the equatorial Pacific appears consistent with observations, although the subsequent generation of coastal waves and reflected Rossby waves at the eastern boundary is stronger and more coherent in the model than in observations. There is also evidence in the model of equatorial Kelvin wave activity in the Atlantic, which is somewhat evident, but less clear in the observations. A role for the MJO in forcing dynamical variability in the Atlantic has been postulated previously, but not proven (Foltz and McPhaden 2004; Han et al. 2008; Webber et al. 2010). Our results indicate that at least a portion of this observed variability arises through wind forcing associated with the MJO.
The conclusion of this simulation is that the majority of the observed SSH variability composited relative to the MJO can be explained by the linear dynamical response to surface wind forcing. This implies that other, nonlinear processes are less important, although they may account for the decay in the equatorial Kelvin wave signals as they propagate across the Pacific. There are also large coherent regions of SSH anomalies outside of the equatorial waveguide in the Pacific in both hemispheres. Additionally, the model simulates the relatively large (>1.5 cm) variability of SSH in the Red Sea noted by Webber et al. (2010). The positive SSH anomalies in phase 7 coincide with strong southeasterly wind anomalies over the preceding two MJO phases (not shown), which will lead to anomalous inflow of water into the Red Sea and thus raise the sea surface; the reverse process accounts for the negative SSH anomalies in phase 3. The variability in the Gulf of Carpentaria found by Oliver and Thompson (2011) and in the observational composites here is not reproduced in quite the same manner; the observed strong positive (negative) SSH anomalies in phases 6–7 (2–3) correspond to a dipole pattern of SSH anomalies in the model. It is not clear what causes this discrepancy, but it may be due to inaccurate representation of the shallow bathymetry of the area in the model or due to other surface processes not included in this model.
6. Realistic primary Madden–Julian event wind forcing experiment
We now move to the question of how ocean dynamics might force primary Madden–Julian (MJ) events. Webber et al. (2012) investigated this question using a combination of composite and case study analysis applied to observational data. Here, we will examine their case study of the primary event starting on 24 September 2004 in further detail by attempting to simulate the oceanic equatorial Rossby wave activity that preceded the triggering of this event. Four experiments are run. The first uses the observed winds for the 150 days previous to the start of the primary MJ event to drive the ocean circulation, that is, wind forcing from 27 April 2004. The later experiments use the winds for the 120, 50, and 30 days previous to the start of the primary MJ event, respectively. This will allow insight into the key components of the atmospheric forcing necessary to generate this oceanic equatorial Rossby wave.
a. 150-day experiment
Figure 7a shows the observed SSH and SST anomalies in the off-equatorial Rossby waveguide for the 100 days prior to the initiation of the event of 24 September 2004, until 10 days after the event. Note that day zero of the primary event (i.e., 24 September 2004) is defined to be the maximum principal component amplitude during “phase A” of the MJO, following the definition of Matthews (2008). Therefore, the actual initiation of the convection associated with this MJ event occurs around 10–15 days earlier, as shown by the OLR anomalies (dashed contours in Fig. 7b).
The arrival of the downwelling Rossby wave in the western Indian Ocean can be clearly seen in Fig. 7a as can the resultant positive SST anomalies, as shown by the diagonally oriented solid contours in the western Indian Ocean between days −60 and +10. The location and timing of the warm SST anomalies agree well with the initiation of the convective anomalies associated with the primary MJ event (dashed contours in Fig. 7b). There are likely to be other factors influencing the triggering of this event. For example, there are cold SST anomalies (dashed contours in Fig. 7a) in the eastern Indian Ocean that are associated with suppressed convection in this region, which may induce circulation favorable for the initiation of convection in the western Indian Ocean. In addition, there is some evidence for a weak convective signal that circumnavigates the globe following the convective activity in the Pacific at lag −50 (not shown) and may also be favorable for the triggering of the event. Nevertheless, given the strength of the dynamically induced SST anomalies, it seems reasonable to conclude that the arrival of the downwelling Rossby wave has a role in triggering this primary MJ event, following the arguments of Webber et al. (2010, 2012).
When the model is forced using observed winds over the 150 days prior to the primary event of 24 September 2004, it is able to simulate the observed wave propagation with reasonable accuracy (cf. SSH anomalies in Figs. 7a,b). The magnitude of the SSH anomalies is somewhat smaller in the simulation than in the observations but the pattern of anomalies agrees well, albeit with less small-scale noise than the observations. These differences may be due to erroneous noise in the observations, to genuine variability not captured in the temporally filtered 1° wind forcing, or possibly to variations in the propagation speed of the Rossby waves between the model and observations. There appear to be two distinct downwelling Rossby wave signals: one (wave A) that is triggered by reflection of the downwelling Kelvin wave arriving at the eastern boundary (Fig. 7d) and a later one (wave B) that begins near the eastern boundary at day −60 and appears to be augmented by local wind forcing around day −40. The combination of these waves appears to generate the broad region of positive SSH anomalies in the western Indian Ocean around days −20 to −10. The latter wave (wave B) continues to propagate westward after the MJ event is triggered and is further amplified following the subsequent easterly winds (days +10 to +30, not shown). However, wave B may be somewhat less important than wave A for triggering the primary event itself, judging by the location at which convection is triggered (west of 50°E at day −10).
b. Reemergence of dynamical ocean anomalies
Analysis of the zonal wind averaged over 10°N–10°S in the Indian Ocean (Fig. 7c) shows that there was a long interval of relatively weak equatorial winds during boreal summer (approximately day −120 to −50), with only a slight increase in variability in the 50 days preceding the primary event. However, around day −140 there was a large westerly wind burst, especially intense in the region of 70°–80°E (black line in Fig. 7c). This appears to be linked to the triggering of the downwelling equatorial Kelvin wave shown in Fig. 7d. The winds remain westerly for a period of 20 days after this wave is triggered, further intensifying the signal.
There is also an easterly wind burst around day −40 that may be partly responsible for the intensification and rapid westward propagation of wave B around day −40 in Fig. 7b. It is likely that the period of relative quiescence between day −120 and −40 is partly due to the northward displacement of the monsoonal circulation and MJO during boreal summer (Wang and Rui 1990). The westerly wind burst that triggers the initial equatorial Kelvin wave is thus probably associated with the last major MJO cycle of the boreal spring season (140 days prior to 24 September is 7 May).
A key result is, therefore, that a westerly wind burst in the central Indian Ocean associated with a strong MJ event can trigger equatorial wave activity in the Indian Ocean that can affect the MJO 140 days (approximately five months) later. This is associated with the “reemergence” of SST anomalies generated by ocean dynamics forced by wind variability during the preceding boreal spring and subsequently leading to convective anomalies. One hundred, forty days is an extraordinary lead time for the MJO, for which predictability is normally limited to around three weeks (Love and Matthews 2009; Kang and Kim 2010; Vitart and Molteni 2010; Rashid et al. 2011). Thus, even the relatively short lead time associated with the easterly wind forcing at day −40 suggests the potential for unprecedented predictability.
c. Sensitivity to length of wind forcing history
Although we have clearly shown that 150 days of wind forcing simulates the observed wave propagation preceding the primary MJ event of 24 September 2004, it remains an open question whether it is possible to simulate similar dynamics with a shorter period of wind forcing. It is possible that a similar SSH structure could be generated without the westerly wind burst at −140 days. We now turn to this question by examining Hovmöller diagrams of Rossby wave propagation from such model simulations integrated over a shorter period of wind forcing (Fig. 8). The observed SSH anomalies are reproduced in each panel to facilitate comparison. The run with 150 days of wind forcing (Fig. 8a) shows that the model SSH anomalies agree broadly with the location of the observed anomalies, although the fit is not exact, as discussed above.
When the model is forced only with observed winds over the preceding 120 days (Fig. 8b), there is little, if any, evidence of the Rossby wave generated by reflection at the western boundary in Fig. 8a. This is unsurprising given that this wave can be traced back to wind forcing that occurred before the 120 day cutoff. There are some positive SSH anomalies along the propagation path of wave A, which may be locally forced, but these do not propagate coherently as a Rossby wave. However, there is evidence of an upwelling Rossby wave that propagates from the eastern Indian Ocean before being largely overridden by the broad region of positive SSH anomalies between 60° and 90°E, associated with wave B and forced by the easterly wind anomalies around day −40. The positive SSH anomalies associated with this downwelling wave (wave B) do propagate relatively coherently, but do not create a broad region of positive SSH anomalies between days −20 and −10 in time for triggering the primary event.
Figures 8c,d show the results from shorter simulations. For the simulation with 50 days (Fig. 8c), wave B is strongly evident; more so than for the 120-day simulation since the preceding upwelling Rossby wave is not simulated. The SSH anomalies do not agree well with the location and timing of the observed SST anomalies but could still have a role in enhancing the convective anomalies of the primary MJ event. However, from Fig. 8d, it is clear that 30 days of wind forcing is not sufficient to generate any positive SSH anomalies locally in the western Indian Ocean. We therefore conclude that the ocean dynamics are not merely a response to wind forcing in the few weeks before the primary event, but rather a response to wind forcing from more than a month before. Furthermore, the role of the westerly wind burst five months previously is confirmed as being important to the triggering of the primary event.
d. Rossby wave vertical structure
We now analyze the vertical structure of the model Rossby waves in order to investigate the magnitude of the thermocline displacement associated with the SSH perturbations. Figure 9 shows depth–longitude sections for the 450-m surface layer over 45°–95°E in both our model and “observations” from the ECCO ocean state estimate (Wunsch and Heimbach 2007). In both model and observations, the downwelling wave can be traced in the form of negative density anomalies back to the eastern boundary at day −100. These anomalies are maximum around 100-m depth but extend coherently into the deep ocean with a slight westward tilt with increasing depth. This tilt implies upward phase propagation and downward energy propagation and has been observed in equatorial Kelvin waves for the Pacific Ocean forced by the MJO (Matthews et al. 2007, 2010). The coherent propagation across the width of the Indian Ocean reinforces the hypothesis that the wave activity can be traced back to reflection of the equatorial Kelvin wave incident on the eastern boundary prior to day −100.
One difference between the model and the observations is in the near-surface density anomalies, especially from day −60 to 0. The model shows a persistent surface layer of positive density anomalies overlying the negative density anomalies associated with the downwelling Rossby wave, while this layer is much more broken in the observations. This is likely to be associated with surface thermodynamical and precipitation processes not included in our simulations, which will act to modify the density of the surface layers. In addition, nonlinear processes and mixing by near-surface currents are not present in the model and may also play a role.
a. Ocean dynamics
Our results show that the observed wave propagation can be simulated by a simple ocean model linearized about a resting basic state. This realistic behavior is evident even from model simulations forced by an idealized Gaussian westerly wind burst on the equator in the Indian Ocean. The subsequent wave propagation shows many of the features identified in previous studies of the dynamics of this ocean basin in response to wind forcing by the MJO (Oliver and Thompson 2010; Webber et al. 2010). Crucially, the arrival of downwelling oceanic equatorial Rossby waves in the western Indian Ocean around 90 days after the initial wind burst agrees with the findings of Webber et al. and with hypotheses of coupling between the atmosphere and the ocean dynamics on such time scales (Han et al. 2001; Han 2005; Fu 2007). The western Indian Ocean is an important region in which forcing by oceanic equatorial Rossby waves can influence atmospheric convection due to the relatively shallow thermocline, higher rates of entrainment cooling (McCreary et al. 1993), and thus increased correlation between SSH and SST anomalies (Xie et al. 2002). It is worth noting that the model does not contain the complex equatorial current systems of the Indian Ocean, which would be expected to alter the phase speed of Rossby waves through Doppler shifting and by altering the meridional vorticity gradient through which they propagate (McPhaden and Ripa 1990). Therefore, the location of the waves in time and space may vary between model and observations.
The model reproduces the observed propagation of coastally trapped waves, both southward along the coasts of Sumatra and Java into the Maritime Continent and northward around the Bay of Bengal. This latter coastal waveguide continues around the tip of India, where wave energy may escape from the equatorial waveguide or be triggered by near-equatorial winds. Coastally trapped waves continue to propagate around the Arabian Sea and arrive in the western Indian Ocean at around the same time as the equatorial Rossby waves. This coincidence in timing led Webber et al. (2010) to hypothesize that the propagation of anomalies along the coastal waveguide could contribute to feedbacks onto the MJO. This hypothesis is consistent with the results of a two-layer model with flat bathymetry in which the coastal Kelvin wave arrives at the same time as the equatorial Rossby wave. However, simulations with more realistic stratification and bathymetry suggest that interactions between baroclinic modes and barotropic shelf wave modes lead to higher levels of energy dispersion such that the coastally trapped wave signal never reaches the western Indian Ocean. Therefore, we can conclude that it is the oceanic equatorial Rossby wave dynamics that are important for understanding coupling between the MJO and the dynamics of the Indian Ocean.
Another constraint on the potential for feedbacks from the ocean dynamics onto the MJO appears to be the latitudinal decay scale of the westerly wind burst that generates the initial wave activity. In order for the eventual SSH anomalies to be coherent, it appears that the Gaussian westerly wind burst must have a latitudinal decay scale greater than 5° (10° is certainly sufficient). Narrower wind bursts change not just the magnitude but also the sign of the eventual equatorial anomalies, with substantially altered spatial patterns for the SSH anomalies. It is hypothesized that this difference between the simulations is due to different projections of the initial wind forcing onto meridional and baroclinic wave modes. This variation in the dynamic response suggests that atmospheric variability needs to be coherent over a relatively broad latitudinal range, centered on the equator, in order to effectively force dynamical variability in the equatorial oceans. In addition, the boreal summer MJO, which exhibits less coherent variability around the equator (Wang and Rui 1990) will probably be less effective in forcing such dynamics.
When the model is forced by realistic composite MJO winds globally, it is able to replicate the observed composite SSH anomalies. This is an important result as it reinforces the notion that the global SSH response to the MJO is a simple product of the linear dynamic response to surface momentum flux; nonlinear and thermodynamic effects are not required and neither is a realistic ocean circulation. We might, therefore, infer that the interannual variability of the ocean will have a relatively minor effect on the dynamical ocean response to the MJO, which will be predominantly a linear addition to this low-frequency variability.
In addition, this model run highlights some intriguing characteristics of the global dynamic ocean response to the MJO. For example, there is evidence of equatorial Kelvin wave propagation in the Atlantic, which subsequently can be seen to propagate in both directions along the coast of Africa. Such behavior has been seen by Polo et al. (2008), but has not previously been attributable to the MJO; our results suggest that at least a portion of this can be attributed to a local linear response to wind variability associated with the MJO. The model also simulates substantial SSH variability in the Red Sea, in agreement with observations (Webber et al. 2010), but fails to simulate the observed SSH variability in the Gulf of Carpentaria (Oliver and Thompson 2011). Such forms of local variability are likely to have a substantial impact on biological productivity due to variations in the upwelling of nutrients. This would build upon existing research that has established a role for the MJO in generating chlorophyll variability through such mechanisms (Waliser et al. 2005; Isoguchi and Kawamura 2006).
b. Primary events
A model simulation with observed winds over the 150 days preceding the primary MJ event of 24 September 2004 closely reproduces the observed SSH variability previously implicated in triggering this primary event (Webber et al. 2012). A westerly wind burst 140 days before the primary event is implicated in generating a downwelling Kelvin wave that subsequently reflects into the first of two Rossby wave signals seen (wave A). This wave propagates westward to arrive in the western Indian Ocean shortly before the initiation of convection, although there is a slight difference in the propagation speed of this wave between model and observations. Simulations with shorter periods of wind forcing are less successful at replicating the observed dynamics, although 50 days is sufficient to generate the latter portion of the wave activity observed (wave B). From the observations of SST and OLR, it appears that wave A is at least as important as wave B (indeed, probably more so) in generating the warm SST anomalies, which prime the atmosphere for convection. The reduced convection in the eastern Indian Ocean that precedes the initiation of the primary event may also be influential through modifications to the atmospheric circulation; these reduced convective anomalies coincide with cold SST anomalies that may be partly driven by the preceding upwelling Rossby wave.
To put these forcing events into context, Fig. 10 shows the intraseasonal OLR anomalies averaged over 15°N–15°S for the period from 150 days prior to 30 days after the primary event of 24 September 2004. The westerly wind burst around day −140 can be clearly linked to intense convective activity over the eastern Indian Ocean region (negative OLR values below −25 W m−2 between 65° and 100°E in early May). This is one of the last strong cycles of MJO activity of the boreal spring season, which is followed by relatively weak and incoherent anomalies during the summer season when convection shifts north of the equator (Wang and Rui 1990). The MJO becomes organized again in a cyclical sense in late September and early October. This reinvigoration of the MJO begins with the primary event on 24 September, which coincides with the reemergence of ocean anomalies forced by a westerly wind burst before the summer period. There may also be a role for atmospheric variability, including the component of the MJO that exists to the north of the equator over the summer months, but it seems likely that the ocean dynamics at least play a role in triggering the primary event in question.
The key conclusion of this paper is therefore that the ocean dynamics act to integrate intraseasonal wind forcing over as much as five months and can subsequently trigger a primary MJ event. In the case study of 24 September 2004, the wave activity is triggered by the last major episode of MJO variability before the summer season. The equatorial waves propagate almost freely back and forth across the Indian Ocean before being reinforced by a second Rossby wave, triggered by easterly wind anomalies around 40 days prior to the primary event. Thus, an MJ event is able to influence another one four cycles later. Although we have only presented results from one case study, composite analysis suggests that downwelling Rossby waves are consistently important for triggering primary MJ events (Webber et al. 2012). Therefore, although the specifics will vary, it is likely that interaction between intraseasonal atmospheric variability and the ocean dynamics, with their differing time scales, is at least partly responsible for the sporadic nature of the MJO.
This potential for coupling between the intraseasonal atmospheric variability and dynamic ocean processes occurring at a longer time scale also brings the potential for forecasting MJO activity at long lead times. Given that current state-of-the art forecasts of the MJO only exhibit skill up to 20–25 days (Love and Matthews 2009; Rashid et al. 2011)—a process predictable over 140 days is exceptional. The degree to which the ocean dynamics can be used to forecast the MJO remains to be tested, but it is likely to be strongest during periods of reduced MJO activity when forecast skill is low (Jones et al. 2004). These results emphasize the importance of monitoring the equatorial Indian Ocean using combinations of moored buoys and ARGO floats, in combination with satellite data, to potentially facilitate such predictability.
The altimeter products were produced by SSALTO/DUACS and distributed by AVISO with support from CNES and are available online at www.aviso.oceanobs.com/. The TMI SST data were obtained from the SSM/I website at www.ssmi.com/tmi/. The EMCWF surface wind data were obtained through the British Atmospheric Data Centre (BADC) and from the ECMWF website at data-portal.ecmwf.int/data/d/interim_daily. The interpolated OLR data were obtained from the NOAA/OAR/ESRL PSD website at www.cdc.noaa.gov/. Data from the ECCO-GODAE project were obtained through the ftp site: www.ecco-group.org/. The Wheeler–Hendon index was obtained from http://www.bom.gov.au/. BGMW was supported by a NERC Ph.D. studentship. The research presented in this paper was carried out on the High Performance Computing Cluster supported by the Research Computing Service at the University of East Anglia. We thank three anonymous reviewers for comments that helped to improve the manuscript.