1 INTRODUCTION

A thermocline is defined as a transition layer between the warmer mixed water of the upper ocean and the cooler subsurface water below. Its spatiotemporal variations can greatly influence climate change, marine fishery and underwater communication. As a semi-closed basin, the South China Sea (SCS) lies in the East Asian monsoon region (Fig. 1) and the study area for this paper is deeper than 200 m. Previous investigators have carried out studies on the spatio-temporal variations of the thermocline, mainly focusing on seasonal variation of the thermocline depth. Liu et al. (2001) found that the thermocline in the central SCS becomes deeper and thinner in winter. Based on the global observational dataset MOODS (Master Observational Oceanographic Data Set) and the Generalized Digital Environmental Model (GDEM), Lan et al. (2006) investigated the thermocline depth in the SCS and found that, in January (winter), the thermocline in the northwestern SCS is deeper than that in the southeast, with two shallower centers located in the northwest of the Luzon Island and the Kalimantan Island, and in April (spring), the thermocline depth is about 30 m almost in the entire SCS. In July and October (summer and autumn), however, the thermocline depth in the northwestern SCS is shallower than that in the southeast. According to the Simple Ocean Data Assimilation (SODA) data in the SCS, Fang et al.(2013) obtained the similar seasonal variation results of the thermocline depth as Lan et al. (2006), indicating that the SODA data, though a global product, are suitable for studying the thermocline in the SCS.

Regarding the mechanism for the seasonal variation of thermocline depth, Liu et al. (2000), and Zhou and Gao (2001) suggested that local wind stress can affect the thermocline through Ekman pumping. Liu et al.(2001) pointed out that the variation of the thermocline is mainly induced by the seasonal cycle of heat flux and wind stress. Besides, Lan et al. (2006) put forward that the ocean circulation and multi-eddy structure in the SCS have significant effects on the thermocline variability. Hao et al. (2012) regarded the surface buoyancy flux (caused by net heat flux and fresh water flux) and the wind stress as the main mechanisms for the seasonal variability of the thermocline.

Although the thermocline depth in the SCS in some seasonal months (January, April, July and October) have been investigated, the thermocline depth in the other months remains unclear, which is our motivation to investigate the monthly evolution of the thermocline depth in the SCS. The lower boundary depth (Z_low), the thickness (∆Z) and the intensity (T_z) of the thermocline in the SCS and the monthly variability of them have rarely been studied. In this study, we use 51-years (1960–2010) monthly seawater temperature data of SODA (Carton et al., 2000a, b, 2005) to calculate the upper boundary depth of thermocline (Z_up), Z_low, ∆Z and T_z in the SCS, focusing on their monthly variability and the mechanisms of monthly variability of Z_up. It is clear that this monthly variability study of thermocline in the SCS can be worthwhile, for example, for regional modellers who would like to benchmark their own results against those from a global product. However, in order to avoid too long a text, we only describe the odd-monthly means of the thermocline parameters and the mechanisms.

The paper is organized as follows: Data and analysis methods are shown in Section 2. Monthly variability of the thermocline in the SCS is presented in Section 3. Mechanisms of the monthly variability of Z_up are discussed in Section 4, and the conclusions and discussion are listed in Section 5.

2 DATA AND METHOD 2.1 Data

The monthly seawater temperature, sea surface salinity and wind stress products of SODA (Carton et al., 2000a, b, 2005) are used for calculations of the thermocline, the sea surface buoyancy flux and the wind stress curl, respectively. The ocean data assimilation model is based on observations including virtually all available hydrographic profile data, as well as ocean station data, mooring temperature and salinity time series, surface temperature and salinity observations of various types. The basic temperature and salinity observation sets consist of approximately 7×10⁶ profiles over the global ocean, of which two thirds have been obtained from the World Ocean Database 2009 (Giese and Ray, 2011). Output variables are averaged every five days, and then mapped onto a uniform global 0.5°×0.5° horizontal grid (a total of 720×330×40-level grid points), spanning the latitude range from 79.25°S to 89.25°N and the longitude range from 179.75°W to 179.75°E, using the horizontal grid spherical coordinate remapping and interpolation package (Carton and Giese, 2008). We download 50×50 profiles that bounded by 0.25°–24.75°N and 100.25°–124.75°E and use 470 profiles in which the depth in the SCS is deeper than 200 m (Fig. 1), and we download the upper 25 vertical layers of the seawater temperature, i.e., 5, 15, 25, 35, 46, 57, 70, 82, 97, 112, 129, 149, 171, 198, 229, 268, 317, 381, 465, 579, 729, 918, 1 140, 1 378, and 1 625 m. The version of SODA used in our study is SODA v2.2.4 (Giese and Ray, 2011) and the variables include ocean temperature (including sea surface temperature (SST)) (℃), salinity, horizontal and vertical ocean velocity (m/s), sea level (m) and wind stress (N/m²).

The monthly net long wave radiation, net short wave radiation, sensible heat flux and latent heat flux data (1960-2010, 1.875°×1.875°) used for the sea surface net heat flux calculation are downloaded from the National Centers for Environmental Prediction (NCEP) reanalysis (Kalnay et al., 1996) (ftp://ftp.cdc.noaa.gov/Datasets/ncep.reanalysis.derived/surface_gauss/). Monthly precipitation (1960-2010, 2.5°×2.5°) and evaporation data (1960-2010, 1°×1°) used for the sea surface buoyancy flux calculation are downloaded from the National Oceanic and Atmospheric Administration (NOAA)'s Precipitation Reconstruction (PREC) Dataset (Chen et al., 2002)(http://www.esrl.noaa.gov/psd/data/gridded/data.prec.html) and the Objectively Analyzed air-sea Fluxes (OAFlux) project (Yu et al., 2008) of Woods Hole Oceanographic Institution (WHOI) (ftp://ftp.whoi.edu/pub/science/oaflux/data_v3/monthly/evaporation/), respectively. The oceanic precipitation analysis (PREC/O) is produced by EOF (Empirical Orthogonal Function) reconstruction of historical gauge observations over islands and land areas and the OAFlux project uses the objective analysis to obtain optimal estimates of flux-related surface meteorology and then computes the global fluxes by using the state-of-the-art bulk flux parameterizations.

2.2 Method

There are several methods to determine the thermocline, such as using the quasi-step-function approximation (QFA) (Ge et al., 2003), the 20℃ isotherm and the vertical temperature gradient criterion (GC) methods. However, the QFA method cannot be applied to the region off the shelf (Hao et al., 2008). The 20℃ isotherm method takes the depth of the 20℃ isotherm as the thermocline depth. On the other hand, the GC method defines the thermocline to be the layer with the vertical gradient of temperature continuously larger than a given value (Fig. 2), which can not only determine the upper boundary depth but also the lower boundary depth of thermocline and is properly used off the shelf (Hao et al., 2012). Therefore, we use the GC method to determine thermocline in this study. The gradient criterion of thermocline is 0.05℃/m (CSBTS, 2008).

We define thermocline intensity T_z as

where T_up and T_low are the temperatures of Z_up and Z_low, respectively.

The surface net heat flux is calculated by

where Q is net heat flux and a positive value means absorption of heat by seawater, Q₁, Q₂, Q₃ and Q₄ are the net short wave radiation, net long wave radiation, latent heat flux and sensible heat flux, respectively.

The buoyancy flux B (kg/(m∙s³)) is calculated by

where B_q is the thermal buoyancy due to the net heat flux, B_p is the haline buoyancy due to the net fresh water flux, g is the gravitational acceleration (9.8 m/s²) and ∆ρ is the variation of seawater density, Q is the net heat flux (downward positive; W/m²), ρ is the reference water density (1 024 kg/m³) and C_p is the specific heat of water (3 986.3 J/(kg∙℃)) (Fofonoff and Millard, 1983), α is the thermal coefficient of expansion (1.7×10^-4/℃) and β is the coefficient of unit haline contraction (7.5×10^-4) (McDougall, 1987), S₀, E, and P are surface salinity, evaporation and precipitation, respectively.

The wind stress curl τ_c (N/m³) is defined as

where τ_x and τ_y (N/m²) are zonal and meridional wind stress, respectively.

3 MONTHLY VARIABILITY OF THE THERMOCLINE 3.1 Spatial distribution of thermocline parameters 3.1.1 Upper boundary depth

The odd-monthly distribution of upper boundary depth Z_up in the SCS is shown in Fig. 3a–f. One can see that in January, Z_up in the northwestern SCS is 50–120 m, which is deeper than that in southeast (30– 50 m), with two shallower centers (< 40 m) located in the northwest of the Luzon Island and the Kalimantan Island, respectively. In March (Fig. 3b), Z_up becomes shallower than that in January and continues to shoal in most regions in May (Fig. 3c). On the contrary to January, Z_up starts to deepen in July (Fig. 3d). It becomes deeper in the southeastern SCS (30–50 m) than that in the northwestern SCS (10–30 m), and continues to deepen in September (Fig. 3e). The distribution patterns are the same as that in July with a shallower center (< 20 m) located in the east of the Indo-China Peninsula. In November (Fig. 3f), Z_up becomes deeper than that in September in the northwestern SCS, but shallower in the southwestern SCS.

3.1.2 Lower boundary depth

The odd-monthly distribution of lower boundary depth Z_low in the SCS is shown in Fig. 3g–l. In January (Fig. 3g), Z_low is about 150–190 m in the northwestern SCS and about 150–180 m in the southeastern SCS. In March (Fig. 3h), Z_low in western SCS is 160–180 m, which is deeper than that in east (140–160 m). In May (Fig. 3i), Z_low is between 150 and 170 m in most areas, showing a homogeneous distribution. In July (Fig. 3j), Z_low in the northwestern SCS is 130–160 m, which is shallower than that in the southeastern SCS (160– 190 m). In September (Fig. 3k), Z_lowin the northwestern SCS with a shallower center (120–150 m) east of the Indo-China Peninsula is shallower than that in the southeastern SCS. In November (Fig. 3l), Z_low in northwestern SCS becomes deeper than that in September.

3.1.3 Thickness

The odd-monthly distribution of thickness ∆Z in the SCS is shown in Fig. 3m–r. In January (Fig. 3m), ∆Z is from 120 m to 150 m in the east of the IndoChina Peninsula, while ∆Z is smaller than 110 m in most regions. In March (Fig. 3n), ∆Z increases about 20 m in the areas between the Indo-China Peninsula and the Luzon Island, and about 10 m in the southern SCS in comparison with that in January. In May (Fig. 3o), ∆Z is between 140 m and 160 m in most regions, exhibiting a homogeneous distribution feature. In July (Fig. 3p), ∆Z in the northeastern SCS is 140–160 m, which is a little thicker than that in the southwestern SCS (120–140 m) with a small thinner center (< 120 m) west of the Palawan Island. In September (Fig. 3q) and November (Fig. 3r), ∆Z is generally thinner than that in July.

3.1.4 Intensity

The odd-monthly distribution of intensity T_z in the SCS is shown in Fig. 3s–x. In January (Fig. 3s), from the continental shelf break of the northern SCS to the central SCS, T_z gradually increases from 0.06 to 0.10℃/m. In comparison to January, T_z generally becomes weaker than that in January. In May (Fig. 3u), T_z in the northwestern SCS is between 0.08 and 0.09℃/m, which is weaker than that in the southeastern SCS (0.09–0.11℃/m). In July (Fig. 3v), T_z is between 0.08 and 0.10℃/m in most regions. In September (Fig. 3w), T_z is generally stronger than that in July, especially east of the Indo-China Peninsula. In November (Fig. 3x), the values of T_z are widely smaller than that in September.

3.2 Spatial means of thermocline parameters

In order to analyze the monthly variability of Z_up, Z_low, ∆Z and T_z, we calculate their spatial means and the standard deviations as shown in Fig. 4. One can see that Z_up is the deepest in January (about 54 m) and the shallowest in May (about 17 m) (Fig. 4a, blue bars). It increases gradually by approximately 5.4 m from May to the January of the following year, and decreases by about 12 m from February to May. On the other hand, Z_low remains 162 m (±2 m) throughout the whole year (Fig. 4a, green bars). Monthly variability of ∆Z, shown as light purple bars in Fig. 4b, is in antiphase with Z_up. It is thickest in May (about 148 m) and thinnest in January (about 110 m).

The spatial mean and the standard deviation of monthly T_z values are shown as orange bars in Fig. 4b. One can see that T_z is weakest in March (0.07–0.08℃/m) but strongest in September (0.09– 0.10℃/m). T_z gradually strengthens from March to September, but weakens from September to the following March.

4 MECHANISMS OF MONTHLY VARIABILITY OF Z_up 4.1 Factors influencing monthly variability of Z_up

In this section, we investigate the physical processes responsible for the monthly variability of Z_up. Previous studies have revealed that surface net heat flux and wind stress play important roles in seasonal variation of the thermocline depth. Heat and fresh water flux at the sea surface may change the stability of upper layer stratification, thus affecting Z_up. The sea surface buoyancy consists of the thermal buoyancy (caused by the net heat flux) and the haline buoyancy flux (caused by the fresh water flux). The buoyancy flux, a comprehensive index of thermodynamics, reflects the coherent role of air-sea heat and haline exchanges (Gill, 1982; Schmitt et al., 1989). Wind stress curl also affects the thermocline depth through Ekman pumping. Besides, Z_up is inseparable from the lower boundary of mixed layer. Lozovatsky et al. (2005) suggested that the sea surface buoyancy and the wind stress are the main driving mechanisms leading to the mixed layer variation in the North Atlantic Ocean. Based on this analysis we examine connections of Z_up to the buoyancy flux and the wind stress curl over the SCS.

4.2 Relationships of Z_up with buoyancy flux and wind stress curl

According to Eq.3, the variation of the sea surface buoyancy flux (B) depends on the variability of seawater density (∆ρ). If B is positive, i.e., the loss of buoyancy for surface seawater, then the seawater density increases and the seawater will decline, so that Z_up will deepen without regarding to the other outside forces. If B is negative, the seawater density decreases and Z_up shoals. In the case of B=0, the buoyancy flux does not affect the thermocline. On the other hand, the wind stress also drives the movement of seawater. The positive wind stress curl (τ_c, Eq.4) can generate the cyclonic vorticity on the sea surface and then the upwelling, which can lift the thermocline. On the contrary, the negative wind stress curl can generate the anticyclonic vorticity on the sea surface and then the downwelling, which can sink the thermocline. To sum up, Z_up deepens while B is positive or τ_c is negative, and shoals while B is negative or τ_c is positive.

4.3 Interpretation of monthly variability of Z_up 4.3.1 Qualitative analysis

Based on the relationships of Z_up with B and τ_c, we examine the correlation of Z_up to B and τ_c. The results indicate that Z_up has positive correlation with B and negative correlation with τ_c. Z_up presents the best correlation with one-month-advanced B and twomonth-advanced τ_c, with the correlation coefficient being 0.90 and -0.73 (Table 1), respectively. Accordingly, we conduct a comparative analysis of odd-monthly Z_up (January to November, Fig. 3a–f) with even-monthly B (December to the following October, Fig. 5a–f) and odd-monthly τ_c (November to the following September, Fig. 5g–l).

On the continental shelf break of the northern SCS, B reaches a positive maximum in December (Fig. 5a) and τ_c demonstrates a larger negative value in November (Fig. 5g), causing the deepest Z_up in January (Fig. 3a). As for the two shallower centers northwest of the Luzon Island and the Kalimantan Island, they are mainly resulted from the largest positive τ_c and both the negative B and positive τ_c, respectively.

In February (Fig. 5b), B is smaller than that in December, with a negative value in most areas, so that Z_upin March (Fig. 3b) is shallower than that in January. Meanwhile, the two shallower centers of Z_up become shallower due to the positive τ_c (Fig. 5h).

In May, Z_up (Fig. 3c) is the shallowest throughout the whole year, mainly owing to the larger negative B in April (Fig. 5c), which plays a dominant role in comparison with τ_c because of the small |τ_c| (Fig. 5i).

In July and September, Z_up (Fig. 3d, e) in the northwestern SCS is shallower than that in the southeastern SCS. This is because B in the northeastern SCS is smaller than that in the southeastern SCS in June and August (Fig. 5d, e), while τ_c is positive in the northwestern SCS and larger negative in the southeastern SCS in May and July (Fig. 5j, k).

In October, positive B (Fig. 5f) in the northernmost study region results in the deeper Z_up in November (Fig. 3f). In the region between the northwest of the Kalimantan Island and the west of the Luzon Island, τ_c is positive and a little less than zero (Fig. 5l) on either side of that region. All of these values lead to the fact that Z_up in the northern and southeastern SCS is deeper than that in the other areas (Fig. 3f).

4.3.2 Quantitative interpretation

From the above qualitative analyses, the buoyancy flux and the wind stress curl can well account for the monthly variability of Z_up. In order to analyze the influence level of the buoyancy flux and the wind stress curl on Z_up, we explore the relationships of Z_up to the buoyancy flux or the wind stress curl without considering the wind stress curl (|τ_c| < 1×10^-9 N/m³, close to zero) or the buoyancy flux (|B| < 1×10^-7kg/(m·s³), close to zero), respectively. From the linear regression of their results (Fig. 6), it is evident that Z_up shows direct proportion to the increasing B, but inverse proportion to the increasing τ_c. Their empirical relations are given by Eqs.5 and 6, where Z_{up_B} and Z_{up_τc} denote the variations of Z_up only resulted from B and τ_c, respectively. From Eq.5, Z_up deepens from 4.2 m to 5.0 m (mean 4.6 m), if B increases by 1×10^-5 kg/(m·s³) without regarding τ_c. From Eq.6, however, Z_up shoals from 1.9 m to 3.1 m (mean 2.5 m), if τ_c increases by 1×10^-7 N/m³ without taking B into account. In addition, Z_up is 45 m without considering either B or τ_c according to Eqs.5 and 6, resulting in the linear relationship of Z_up to B and τ_c as expressed in Eq.7. Figure 7 describes the odd-monthly mean of Z_up resulted from Eq.7. The patterns and values (Fig. 7a–f) are very similar and close to the Z_up as shown in Fig. 3a–f, expect in May when Z_up resulted from Eq.7 is about 25–30 m (Fig. 7c), which is larger than Z_up resulted from SODA (10–15 m) (Fig. 3c) in the northwestern SCS.

As mentioned above, the factors influencing Z_up are the buoyancy flux (B) and the wind stress curl (τ_c). To further explore the relative importance of B and τ_c on Z_up, we calculate the contribution of B and τ_c to Z_up using Eq.8, where p means the contribution from the buoyancy flux and 100% minus p represents the contribution from the wind stress curl.

As shown in Fig. 8, most of the Z_up values are influenced by both the buoyancy flux and the wind stress curl from December to the following February when the wind stress curl dominates the SCS deep basin. The buoyancy flux controls Z_up in most regions from April to November (Fig. 8d–k), among which Z_up is mainly affected by the buoyancy flux in May and June, while from July to September, Z_up in the north of 12°N is mainly influenced by the buoyancy flux with the southern part being dominated by both the buoyancy flux and the wind stress curl.

5 CONCLUSION AND DISCUSSION

This study examines variability of the upper and lower boundaries of thermocline (Z_up and Z_low), its thickness (∆Z) and intensity (T_z) in the SCS, and explores the mechanisms responsible for the monthly variability of Z_up. The major results are summarized as follows.

The climatological monthly mean Z_up, Z_low, ∆Z and T_z of the SCS are directly or indirectly determined by the vertical temperature GC method from monthly seawater temperature data (1960–2010). The results show that the spatial mean of Z_up gradually shoals from February to May and deepens from May to the following January, ∆Z is out of phase with Z_up because Z_low remains unchanged all year round.

This study reveals that both the surface buoyancy flux and the wind stress curl play dominant roles in monthly variability of Z_up. The variability of surface buoyancy flux (B) resulting from net heat and fresh water flux variations can lead to mix the surface water and change the thickness of the mixed layer, so as to affect the depth of thermocline. The wind stress curl (τ_c) has the ability to influence Z_up through generating upwelling or downwelling. In order to explore how B and τ_c quantitatively affect Z_up, we investigate the monthly B based on the heat flux and precipitation data from NCEP, the evaporation data provided by the Woods Hole Oceanographic Institution, and the τ_c according to wind stress data from SODA. The results reveal that Z_up shows the best correlation to onemonth-advanced B and two-month-advanced τ_c.

We have analyzed the effects of B and τ_c on Z_up and obtained the regression equation:Z_up=4.6(±0.4)×10⁵B– 2.5(±0.6)×10⁷τ_c+45, which means that Z_up has positive correlation with B and negative correlation with τ_c. Z_up deepens from 4.2 m to 5.0 m (mean 4.6 m) when B increases by 1×10^-5 kg/(m·s³) without considering τ_c, and shoals from 1.9 m to 3.1 m (mean 2.5 m) when τ_c increases by 1×10^-7 N/m³ regardless of B. The relative importance of the buoyancy flux and the wind stress curl to Z_up is also examined, indicating that Z_up is mainly controlled by the buoyancy flux throughout the whole year.

Although the Eq.7, resulting from a linear method, can quantitatively account for the influence of B and τ_c on Z_up, there are still some points remaining for discussion. Before using a linear method, the grid points of Z_up and B, Z_up and τ_c in the SCS must be the same. Z_up and τ_c which are calculated from SODA products with a 0.5°×0.5° horizontal resolution have the same grid points, but B which is calculated from NCEP (1.875°×1.875°), NOAA (2.5°×2.5°) and WHOI (1°×1°) products has a much fewer grid points than Z_up. Therefore, we use the 2-D interpolation method to make B a 0.5°×0.5° horizontal resolution before analyzing the relationship between B and Z_up. This may have distorted the Eq.7 somewhat because some grid points of B are from interpolation which could not represent the observations. Moreover, the wind stress (not the wind stress curl) on the sea surface can generate and maintain turbulence to the upper ocean, which would also influence the depth of mixed layer and thermocline. We will take into account the influence of wind stress on Z_up in our next step.

6 ACKNOWLEDGEMENT

We thank professor John Hodgkiss of the University of Hong Kong for his help with English.