In this study, the MLD was defined as the depth at which σθ (potential density) is 0.125 kg m-3 denser than at the sea surface, and the maximum MLD in a given year was used to calculate the global subduction rate [Liu and Huang, 2011]. In the North Pacific, the MLD front is located to the east of Japan, where it tilts northeastward, indicating an eastward advection effect of the Kuroshio Extension.
In the 20C3M simulation results (Fig. 1a), the MLD is about 300 m deep to the east of Japan, with a maximum of about 350 m at (40°N, 150°E). To the north of the subtropical gyre of the North Atlantic, the equatorward shoaling MLD extends northeastward across the central North Atlantic, with a maximum of 400 m. The deep mixed layer extends into the subpolar region and Labrador Sea, with the depth varying from roughly 600 m in the east to more than 800 m in the west. In the east of the subtropical gyre of the South Pacific, an MLD maximum of about 200 m is centered at (30°S, 100°W). Two zonal northwest slant bands of the deep MLD are found in the Southern Ocean—one is located in the east of the South Indian Ocean [Tsubouchi et al., 2010] and south of New Zealand (30°-60°S), with a maximum of 600 m and extending from 50°-180°E; while the other is located in the Southeast Pacific (45°-60°S), with a maximum of 700 m and extending from about 130°-50°W.
In the SRES A1B runs (Fig. 1b), the spatial pattern of the MLD is similar to that in 20C3M; however, the MLD becomes shallower in most parts of the global ocean, particularly in the Southern Ocean. This significant shoaling of the MLD results from a more stratified upper layer due to greater ocean warming near the surface. In comparison with 20C3M, the MLD is about 100 m shallower to the southeast of Japan and about 50 m shallower along the path of the Gulf Stream. In the Southern Ocean, the deep MLD areas are also seen to shrink significantly under the warming climate, being about 100 m shallower to the south of the New Zealand and in the east of the South Indian Ocean. In the southeast of the South Pacific, the MLD is 250 m shallower than in the 20C3M results.
Figure 2 shows the distribution of maximum MLD in the World Ocean Atlas 2009 (WOA09) and, as can be seen, the pattern is similar to that derived from the 20C3M experiment. However, the MLD in the subtropical gyre of the Northern Hemisphere is roughly 50 m deeper in the results of the IPCC AR4 models than in WOA09, and in the former the mixed layer maximum occupies a larger region. In the Antarctic Circumpolar Current (ACC) region, the MLD in the IPCC model results is relatively shallower and shifts poleward in comparison to WOA09 (Fig. 2), and this may result from the different MLD distributions in the ACC region in the IPCC AR4 models.
The annual subduction rate in the global ocean was calculated for both the 20C3M and the SRES A1B scenarios. From the kinematic perspective, by following a particle released at the base of the mixed layer, the Lagrangian subduction rate is calculated at the time when the effective detrainment starts throughout the entire seasonal cycle (Huang and Qiu, 1994). The annual subduction rate (Sann) is defined as follows:
where hm denotes the mixed layer depth and wnb is the vertical velocity of a water parcel at the base of the mixed layer. By integrating from the end of the first winter (t1) to that of the second winter (t2) over one year (T), the Sann could be obtained. On the right-hand side of Eq. (1), the first term denotes the contribution of vertical pumping at the base of the mixed layer. The second term, which is called lateral induction, represents the contribution due to the slope of the mixed layer base. Compared to previous studies (Huang and Qiu, 1994), the Ekman pumping velocity and a reduction term due to the meridional velocity in the surface layer is replaced by the vertical velocity.
As mentioned above, regions of high subduction rate are generally characterized by thick mixed layers. In the Southern Hemisphere in the 20C3M runs (Fig. 3a), strong subduction occurs around the northern region of the Southern Ocean in the Pacific and eastern Indian sectors, with maximums of more than 300 m yr-1 and 250 m yr-1, respectively. In the subtropical region of both the Pacific and Atlantic, the subduction rates are roughly 100 m yr-1, corresponding to local MLD fronts. In the subpolar region of the North Atlantic, there is also a maximum subduction rate of about 250 m yr-1 at around (55°N, 40°W) corresponding to the local deep MLD [Marshall et al., 1993].
Similar to the response of the MLD to the SRES A1B global warming scenario, the spatial pattern of global subduction rate does not change (Fig. 3b), but there are reductions of subduction rates in most parts of the global ocean. In the subpolar region of the North Atlantic, the subduction rate is around 75 m yr-1, which is 25 m yr-1 less than in 20C3M. A decrease in subduction can also be identified in the subtropical region of both hemispheres, especially in the North Pacific, as found by Luo et al., (2009). Compared to the subtropical region, there is noticeable weakening of subduction under the SRES A1B scenario in the northern region of the Southern Ocean in the Pacific sector, with a maximum of about 250 m yr-1, and in the eastern Indian sector, with a decrease of about 50 m yr-1 (Fig. 3b).
In order to evaluate the uncertainties of the model results, we integrated the subduction volume in the North Pacific (20°-50°N, 120°-240°E), North Atlantic (20°-50°N, 80°-0°W) and South Pacific (10°-60°S, 140°-290°E) for each model, respectively (Fig. 4)—regions for which there are observational estimates. The model mean subduction volume is 44.2 Sv, 25.5 Sv and 65.3 Sv in the North Pacific, North Atlantic and South Pacific, respectively. Although these values are slightly different from the 35.2 Sv and 27.0 Sv reported by [Qiu and Huang, 1995] for the North Pacific and North Atlantic, respectively, and the 48.8 Sv reported by [Qu et al., 2008] for the South Pacific, all of which were derived from observational data, the multi-model mean results still capture the major subduction process in the main subduction regions. Figure 5 presents the subduction volume of the global ocean in the 20C3M and SRES A1B experiments separately for each model. As can be seen, all of the models present a substantial decrease of the subduction rate under a warming climate (Fig. 5a), and most of them also depict a decrease in subduction in the Southern Ocean, except for Cnrm_cm3 and Mri_cgcm2_3_2a (Fig. 5b). The standard deviation for the change in subduction volume in the global ocean and Southern Ocean is 15.2 Sv and 18.2 Sv respectively, reflecting the fact that this decrease in subduction is significant. We also checked the spatial pattern of the response of the subduction process in SRES A1B for each model, and the results revealed that all the models produce a similar decreasing spatial pattern, except for CNRM_CM3 and MRI_CGCM2_3_2a in the Southern Ocean, with an increase of subduction (not shown).
Figure 6 shows that the contribution of lateral induction, which is closely associated with mixed layer fronts (Xie et al., 2000), also has a decreasing over the global ocean, especially at around 45°-60°S in the Southern Ocean. The lateral induction terms for the 20C3M runs (Fig. 6a) present a similar spatial pattern as the subduction rate, suggesting that lateral induction may play a dominant role in subduction in most parts of the global ocean, except in the tropical region. In contrast to lateral induction, no significant differences can be seen in the vertical pumping terms between the 20C3M (Fig. 7a) and SRES A1B (Fig. 7b) results. In the Northern Hemisphere, the subduction associated with vertical pumping presents a northwest slant in the subtropical area, while in the Southern Hemisphere it presents a southeast slant in the subtropical area. As suggested by (Huang and Qiu, 1994), the vertical pumping term is predominantly controlled by the wind stress curl. A negative value in the tropical region of the Pacific and Indian oceans is apparently induced by wind-driven upwelling, which may imply that subduction would be inhibited by the vertical pumping term in the tropical area.
The total volumetric subduction rate integrated over the global ocean and averaged from 1951-2000 is estimated at 308.6 Sv from the 20C3M runs, which is about 40 Sv more than that (268.8 Sv) averaged from 2051-2100 from the SRES A1B runs (Fig. 8a), with the standard deviation of the subduction response being 15.1 Sv. The integrated subduction associated with the lateral induction term over the global ocean and averaged from 1951-2000 is estimated at 294.1 Sv from the 20C3M runs, which is about 50 Sv more than that (244.4 Sv) averaged from 2051-2100 from the SRES A1B runs (Fig. 8b). Note that, according to the calculation suggested by Liu and Huang(2011), the negative value was set to zero, so that the differences in lateral induction between the 20C3M and SRES A1B runs may be larger than those of net subduction volume due to the absence of cancellation between the lateral induction and vertical pumping. The integrated subduction associated with the vertical pumping term is estimated at 134.0 Sv from the 20C3M runs, which is only about 5 Sv more than that (128.8 Sv) averaged from 2051-2100 from the SRES A1B runs (Fig. 8c). This indicates that the significant decrease in the global subduction volume is dominated by the lateral induction process, and the vertical pumping term plays a secondary role in this reduction. Figure 8d depicts the subduction volume of the Southern Ocean in different runs. The volumetric subduction rate integrated over the Southern Ocean averaged from 1951-2000 is estimated at 163.0 Sv from the 20C3M runs, which is about 24 Sv more than that (139.0 Sv) averaged from 2051-2100 from the SRES A1B runs (Fig. 8d), with the standard deviation of the subduction response being 18.2 Sv. These results show that subduction of the Southern Ocean plays the most important role in the subduction of the global ocean, and thus the reduction of subduction in the Southern Ocean plays a key role in the weakening of global subduction.
As noted previously, subduction consists of a vertical pumping term and a lateral induction term. However, lateral induction, through which mixed layer waters are swept into the permanent thermocline, may play the dominant role in reduction of subduction under global warming. Since lateral induction of subduction is typically associated with MLD horizontal gradients, the change of subduction under a warming climate may be closely related to the weakened MLD front. Furthermore, this may indicate that the winter MLD front could make a significant contribution to the evolution of subduction. Under a global warming scenario, wintertime convection can be inhibited by a more stratified upper ocean, resulting in a shallower MLD and thus a reduction of lateral induction.
The stratification of the upper ocean is influenced by the buoyancy flux (Bnet) (Downes et al., 2009; Downes et al.,2010), which is defined as follows:
In Eq. (2), a negative Bnet indicates buoyancy loss (i.e., making the density of surface waters weighted) and a positive value indicates a buoyancy gain. The variableρis sea water density at the surface, and the gravitational force is given by g. The first term in Eq. (2) is the buoyancy flux induced by surface heat flux (H), with Cw being the heat capacity of water andαbeing the thermal expansion coefficient. The second term is the buoyancy flux induced by the surface freshwater flux (W), with S being the mixed layer salinity and β being the haline contraction coefficient. The third term represents the buoyancy gain or loss amassed by the Ekman drift, where k is the unit vertical vector, τis the wind stress, and f is the Coriolis parameter.
Since MLD reaches a maximum in late winter, we diagnosed each term in Eq. (2) at that time. Figure 9 depicts the three terms of the buoyancy flux in late winter from the results of the 20C3M experiment. In the subtropical region of both hemispheres, the term induced by surface heat flux shows a minimum <-10×10-5 N m-2 s-1. In the ACC region, it also reveals a buoyancy loss >-5×10-5 N m-2 s-1 through heat loss from ocean to atmosphere. The term induced by freshwater flux (Fig. 9b) tends to make the surface water lighter in the equatorial region, with a maximum >1.5×10-5 N m-2 s-1, while in the subtropical region more evaporation than precipitation tends to make the surface water denser, with a minimum <-1×10 N m-2 s-1. The term associated with Ekman drift (Fig. 9c) contributes less outside of the equatorial region, while the term induced by surface heat flux may play a dominant role in the global ocean in later winter.
Figure 10 displays the response of each term in Eq. (2) under a warming climate. In the subtropical region of both hemispheres, the term induced by surface heat flux tends to make the surface water lighter, with a maximum >1× 10-5 N m-2 s-1 (Fig. 10a). In the ACC region and the subpolar region of the North Atlantic, it also reveals a buoyancy gain, with a maximum >1.5×10-5 N m-2 s-1. The term induced by freshwater flux (Fig. 10b) and Ekman drift (Fig. 10c) all contribute a positive buoyancy flux in the equatorial region, with maximums of roughly 5× 10-6 N m-2 s-1 for the freshwater flux and 2× 10-6 N m-2 s-1 for the Ekman drift. The term induced by freshwater flux also makes a positive contribution in the Southern Ocean through more precipitation than evaporation. However, both the term induced by freshwater flux and that by Ekman drift contribute much less to the buoyancy flux than the term induced by surface heat flux (Fig. 10a). This implies that surface heat flux may play a dominant role in the intensification of the stratification of the surface layer (Fig. 10a).
The stratification of the upper ocean is mainly influenced by changes in atmosphere-ocean heat flux under global warming, which can be seen to have notable differences in the simulation results of 20C3M and SRES A1B. The climatological heat flux in 20C3M (Fig. 11a) has a similar spatial pattern to that in SRES A1B (not shown). However, the areas of strong ocean-to-atmosphere heat loss shrink significantly as a response to global warming. Large ocean-to-atmosphere heat loss is seen in the west of the subtropical gyre of the North Pacific and the North Atlantic due to the influence of local winter monsoon. In response to the global warming scenario, the heat loss in the SRES A1B runs is about 50 W m-2 less than that in the 20C3M runs in those regions (Fig. 11b). In the subpolar region of the North Atlantic, the heat loss in the SRES A1B runs is about 100 W m-2 less than that in the 20C3M runs. In the Southern Ocean, the heat loss areas also shrink significantly under the warming climate. In the Southern Indian Ocean, the heat loss in the SRES A1B runs is about 40 W m-2 less than in 20C3M, and about 50 W m-2 less in the southern Pacific Ocean around 60°S.