Ocean stratification and the vertical extent of the mixed layer influence the rate at which the ocean and atmosphere exchange properties. This process has direct impacts for anthropogenic heat and carbon uptake in the Southern Ocean. Submesoscale instabilities that evolve over space (1–10 km) and time (from hours to days) scales directly influence mixed layer variability and are ubiquitous in the Southern Ocean. Mixed layer eddies contribute to mixed layer restratification, while down-front winds, enhanced by strong synoptic storms, can erode stratification by a cross-frontal Ekman buoyancy flux. This study investigates the role of these submesoscale processes on the subseasonal and interannual variability of the mixed layer stratification using four years of high-resolution glider data in the Southern Ocean. An increase of stratification from winter to summer occurs due to a seasonal warming of the mixed layer. However, we observe transient decreases in stratification lasting from days to weeks, which can arrest the seasonal restratification by up to two months after surface heat flux becomes positive. This leads to interannual differences in the timing of seasonal restratification by up to 36 days. Parameterizing the Ekman buoyancy flux in a one-dimensional mixed layer model reduces the magnitude of stratification compared to when the model is run using heat and freshwater fluxes alone. Importantly, the reduced stratification occurs during the spring restratification period, thereby holding important implications for mixed layer dynamics in climate models as well as physical–biological coupling in the Southern Ocean.
In the Southern Ocean, the seasonal cycle dominates the mixed layer depth (MLD) variability (Dong et al. 2008; Sallée et al. 2010). At the ocean surface, buoyancy loss during winter months initiates vertical convection and results in an erosion of the stratification and associated deepening of the mixed layer. Buoyancy gain during spring and summer increases the vertical stratification, shoaling the mixed layer to 50 m all around the Antarctic Circumpolar Current (ACC; Sallée et al. 2010). Subseasonal MLD variability is defined here as vertical variations of the mixed layer occurring within hours to months. Mixed layer deepening at these scales arises from wind-driven processes such as mechanical stirring at the surface and Langmuir turbulence, which has shown to improve biases of the Southern Ocean MLD in climate models (Fan and Griffies 2014; Li et al. 2016). In the Southern Ocean, strong atmospheric storms occurring at the synoptic scale are associated with wind speeds regularly exceeding 20 m s−1 (Yuan et al. 2009). The passage of storms is found to erode the mixed layer stratification and deepen summer mixed layers by as much as 50 m (Swart et al. 2015; Nicholson et al. 2016). These synoptic perturbations of the mixed layer have direct implications for biological processes, where the vertical entrainment of nutrients into the mixed layer from below may sustain phytoplankton production across the summer (Swart et al. 2015; Tagliabue et al. 2014; Carranza and Gille 2015; Nicholson et al. 2016). Furthermore, vertical entrainment of essential climate gasses such as carbon dioxide has direct implications for the global carbon cycle (Sabine et al. 2004). Despite this, global climate models fail to accurately simulate the depth and extent of stratification of the Southern Ocean mixed layer. Current simulations provide mixed layers which are too shallow and stratified compared to observations, which has attributed to excess freshwater at the ocean surface (Sallée et al. 2013). The overstratification leads to mixed layers 3°–4°C warmer than the observations (Belcher et al. 2012). One reason postulated for the overstratification is due to a missing parameterization of surface-wave processes that force Langmuir turbulence, which act to deepen the mixed layer (Belcher et al. 2012).
We owe a significant part of our understanding of mixed layer variations to one-dimensional forcing mechanisms (Niiler and Kraus 1977; Price et al. 1978). However, the ocean is impacted by horizontal processes in response to fronts, eddies, and filaments, which can modify upper-ocean stratification. These potentially important dynamics can occur at small spatial scales, namely submesoscales, O(1–10) km (e.g., Thomas 2005; Mahadevan et al. 2010; D’Asaro et al. 2011; Mahadevan et al. 2012; Thompson et al. 2016). One particular submesoscale process is the formation of baroclinic instabilities within the mixed layer, which grows as a baroclinic wave along a front (Haine and Marshall 1998; Boccaletti et al. 2007). The flow dynamics associated with baroclinic instability approach a regime where the Rossby number Ro = is O(1), where ζ = is the vertical relative vorticity and f is the Coriolis parameter. As Ro ≈ 1, the flow departs from geostrophic balance (Thomas et al. 2008). From this dynamical definition, submesoscale motions will be active in regions of large vorticity. Mixed layer baroclinic instabilities arise at the mixed layer Rossby radius of deformation, L = , where N is the buoyancy frequency in the mixed layer, and H is the MLD. Mixed layer baroclinic instabilities have typical length scales of L ~ 1–10 km and can spin down into submesoscale-sized eddies, referred to as mixed layer eddies (MLEs). MLEs can directly impact the mixed layer stratification by rearranging horizontal buoyancy gradients associated with mixed layer fronts (horizontal change in buoyancy over some spatial distance) to vertical stratification through an ageostrophic secondary circulation with upwelling on the lighter side of the front and downwelling on the denser side (Fox-Kemper et al. 2008). In a previous study occurring west of 0°E in the subantarctic, MLEs are argued to promote the spring mixed layer restratification during periods of weak wind forcing (du Plessis et al. 2017).
Surface winds blowing in the direction of the frontal flow (down-front winds) drive a cross-frontal horizontal Ekman advection from the denser side of the front to the lighter side. The cross-frontal flow can force convective instabilities, enhancing mixing through small-scale turbulence, which can increase dissipation within the mixed layer by up to an order of magnitude compared to wind-driven shear mixing (Thomas 2005; D’Asaro et al. 2011). Conversely, up-front winds (winds directed against the frontal flow) advect the lighter side of the front over the denser side, thus increasing the vertical stratification. The wind-driven Ekman advection at fronts is known as Ekman buoyancy flux (EBF).
Given its remoteness and harsh conditions, multi-month observational studies in the Southern Ocean which sample at the spatial and temporal resolutions necessary to resolve submesoscale dynamics are limited. These lack of observations result in an overreliance on high-resolution numerical modeling (Nikurashin et al. 2013; Rosso et al. 2014; Bachman et al. 2017) and relatively short-duration ship-based measurements (Rocha et al. 2016; Adams et al. 2017) to tease out the role of submesoscale processes impacting mixed layer stratification. Therefore, long-endurance observational platforms, such as profiling gliders, are becoming a common tool to address the data requirements to observe these finescale processes. Gliders have already began to provide quasi-continuous observations in the Southern Ocean at horizontal resolutions of less than 5 km and temporal resolutions of 2–5 h (Schofield et al. 2010; Thompson et al. 2014; Swart et al. 2015; Schofield et al. 2015; Erickson et al. 2016; Miles et al. 2016; du Plessis et al. 2017; Viglione et al. 2018).
In this paper, we use data acquired from Seagliders over four separate years in the Subantarctic Zone region of the Southern Ocean (SAZ). We attempt to elucidate the roles of MLE and EBF impacting the subseasonal variability of the mixed layer stratification. We do this by applying already existing parameterizations which scale MLE and EBF as equivalent heat fluxes. These fluxes are incorporated into a one-dimensional mixed layer model to investigate the potential importance of submesoscale processes impacting the seasonal evolution of stratification. Section 2 describes the field deployments of gliders and supplementary data used. Results from the glider experiments and model simulations are presented in section 3, while section 4 comprises the discussion summarized in section 5.
2. Methods, data, and model simulations
a. Field campaign and regional setting
Seagliders sample the top 1000 m of the ocean in a V-shaped pattern and have been shown to provide an adequate resolution for investigating submesoscale dynamics within the mixed layer (e.g., Ruiz et al. 2012; Baird and Ridgway 2012; Mahadevan et al. 2012; Swart et al. 2015; Todd et al. 2016; Thompson et al. 2016; Erickson et al. 2016; du Plessis et al. 2017; Viglione et al. 2018). The field campaign forms a part of the Southern Ocean Seasonal Cycle Experiment (SOSCEx; Swart et al. 2012), with the aim to understand the seasonal cycle dynamics of mixed layer characteristics in the Southern Ocean. The sampling plan for SOSCEx was to deploy a glider at roughly 43°S and 8°E in the SAZ before the onset of seasonal restratification of the mixed layer for the seasons of 2012, 2013, 2015, and 2016 (Fig. 1). All gliders continually sampled the upper ocean within the SAZ for the duration of each experiment before being retrieved in late summer (February/March). The duration of each mission ranges from 3 to 6 months (Fig. 2). The deployments are labeled incrementally from SOSCEx1 to SOSCEx4. Note that SOSCEx1 sampled approximately 1°N of SOSCEx2–4.
The average time taken for a glider to complete a dive is 5 h, while the horizontal resolution between profiles is 1.4 ± 1.1 km (Fig. 1b). The raw data were initially processed using the University of Washington’s base station processing toolbox, which corrects for thermal lag. We manually remove bad profiles before optimally interpolating to a constant time and depth grid of 2 h and 5 m with a Gaussian correlation function of 1 day. Sensitivity analysis (not shown) indicates that this gridding sufficiently resolves the mesoscale gradients in the mixed layer properties. The geographical position of the glider is mapped onto this grid to produce a monotonically increasing along-track distance. We use the horizontal buoyancy difference between each grid point and along-track distance to calculate the horizontal buoyancy gradient. The definition of the MLD follows the density difference criteria of = 0.03 kg m−3 from a reference depth of 10 m (de Boyer Montégut et al. 2004; Dong et al. 2008). At the location and time of each glider deployment and retrieval, CTD calibration casts were performed and checked for sensor drift and corrected accordingly, as in Swart et al. (2015).
b. Observational bias
Calculating the full magnitude of the horizontal buoyancy gradients for a particular front using a glider is only possible when the glider dives perpendicular to the front sampled. Thompson et al. (2016) perform an analysis where the horizontal buoyancy gradient is calculated for a glider sampling a front at all possible angles. Averaging over all these angles leads to the underestimation of the horizontal buoyancy gradient by a factor of . Applying this method to determine the potential bias in our study requires the direction of the front for each glider dive. To do this, we use the direction of the depth-averaged current to represent the frontal flow direction. The depth-averaged current is obtained by comparing the glider’s dead-reckoning positioning system with the true location. Dead-reckoning relies on estimates of speed and direction to propagate a known position forward in time and thus is subject to integration errors. However, it is ideally suited to profiling gliders where no underwater location reference is available once the glider begins to dive and this method has shown to provide good estimates of glider velocity and heading (Rusello et al. 2012). We note that the depth-averaged current is likely to be dominated by both the mixed layer flow and deeper mesoscale motions. However, the length scale of the fronts evaluated in this study is larger than the submesoscale eddies which dominate the mixed layer flow, and thus we aim to represent the front direction dominated by the geostrophic flow. By determining the difference between the front direction and the glider dive direction, we deduce that across all SOSCEx experiments, the gliders underestimated the true buoyancy gradient of the front by on average 64% (Fig. 3). This value is a similar estimate to the 71% estimated by Thompson et al. (2016) and thus provides confidence that although errors may exist in individual calculations of horizontal buoyancy gradients, we are providing a statistical representation of the magnitude of the fronts which we sample. Furthermore, we indicate that around 46% of the profiles across all experiments capture over 80% of the true buoyancy front. Given the mixed layer Rossby radius of deformation L and typical values within the mixed layer of N ~ 10−3 s−1, MLD of 100 m and f (10−4 s−1), a frontal width of about 1 km exists. Should the glider sample perpendicular to the local front, the horizontal resolution of the gliders should resolve the local submesoscale gradient. However, the by the nature of their sampling, gliders do not do this, and therefore we emphasize that we are providing a statistical representation of the mixed layer fronts for this study.
c. Potential vorticity calculations
where is the vertical component of the absolute vorticity, b = is the buoyancy, and U is the three-dimensional velocity vector (u, υ, w). When , the ocean can become susceptible to a number of flow instabilities occurring through either variations of vertical vorticity, stratification, and/or vertical shear of the velocity (Thomas et al. 2013). In particular, unstable stratification (gravitational instability), horizontally sheared flows (centrifugal or inertial instability) or symmetrical instability may arise when q is positive in the Southern Ocean. Note that baroclinic instability does not require . Unstable stratification may occur as a result of down-front winds, which are known to extract q from the ocean surface (Thomas 2005).
Separating q into vertical and baroclinic components is useful as this helps identify what process may make the flow unstable:
The vertical component of q consists of the vertical component of the absolute vorticity and the vertical stratification:
When horizontal variations in velocity and buoyancy are small, q may be positive (unstable) if the flow is unstably stratified, (N2 = < 0). As the horizontal buoyancy gradients and vorticity become larger, a baroclinic component of the flow must be considered,
Term contains the terms related to the nontraditional component of the Coriolis frequency. Estimating q using gliders has successfully been achieved in numerous studies to date (e.g., Shcherbina et al. 2013; Thompson et al. 2016; Todd et al. 2016; Erickson et al. 2016; Viglione et al. 2018). These studies have shown that although possible, calculating q using glider data requires the following assumptions: (i) we ignore terms in the q calculation that involve the vertical velocity w; (ii) we neglect the term, which only makes a small correction to q; and (iii) we assume the flow to be in thermal wind balance, such that the vertically sheared horizontal velocities can be directly related to the horizontal buoyancy gradients = . When these approximations are applied, appears in a more compact form:
where is the glider-derived horizontal buoyancy gradient . Here x is taken as the along-track distance of the glider’s trajectory.
Following Thompson et al. (2016), and are combined to provide an observational expression for q using glider data,
This expression makes the contribution of horizontal buoyancy gradients to a positive clear, as is a definite quantity. A complicated term in is ζ, which is an approximation of the vertical relative vorticity. To determine the values of ζ, we first linearly interpolate the depth-averaged current to the timestamp of the optimally interpolated data and apply a Gaussian smoothing of 1 day, which provides an estimation of the mesoscale frontal flow. Estimates of ζ obtained from the along-track gradient of the frontal flow have a median of 0.28 × 10−4 s−1 for all experiments combined. We are aware that calculations of PV using and that encompass only two-dimensions of the three-dimensional field are likely to lead to errors in relative to the full PV. Observational PV estimates obtained from a glider experiment in the Drake Passage show that is generally lower than the full PV, with the sign of both PV estimates consistent over 90% of the time (Viglione et al. 2018).
d. Submesoscale buoyancy fluxes
1) Ekman buoyancy flux
By the Coriolis deflection, winds directed along a front drive a cross-frontal flow of water over the Ekman layer depth (Thomas 2005; Thomas and Lee 2005), promoting either mixing (down-front winds) or restratification (up-front winds). The process of EBF can be quantified as an equivalent heat flux comparable directly to surface heat fluxes [; Eq. (7)]. Negative values of represent a negative buoyancy flux, while positive values denote a positive buoyancy flux (W m−2),
where is the alongfront component of wind stress τ. Estimations of require knowledge of the direction for both the wind and the front. The magnitude and direction of τ are determined from the National Centers for Environmental Prediction Reanalysis-2 (NCEP-2) 10 m at 6-hourly time intervals (see section 2e for further details). We assume that the submesoscale eddies are growing off larger mixed layer gradients through mixed layer baroclinic instability (Fox-Kemper et al. 2008). Thus, the fronts are at scales larger than individual submesoscale eddies and are more likely to reflect the flows that dominate the depth-averaged current. Thus, is determined from the angle difference between the direction of the depth-averaged current and the wind. Thompson et al. (2016) perform insightful analysis on the error of misrepresenting the wind-front alignment by taking all possible angles of the glider dive direction and wind orientation with respect to a fixed buoyancy gradient. Their summation is that calculating from glider data is likely to represent almost all or little of the true EBF (Fig. 15 in their paper). They note that the error between the true value of and that calculated from the glider data may exist in both the magnitude and the sign. However, over the period of a season, the mean value of estimated from the glider data is smaller than the mean of the true EBF.
2) Mixed layer eddies
Mixed layer baroclinic instabilities drive a thermally direct vertical buoyancy flux within the mixed layer by upwelling the lighter side of a front over the denser side. Fox-Kemper et al. (2008) provide a parameterization of MLEs to represent the vertical rearrangement of buoyancy [Eq. (8)]. Mahadevan et al. (2012) represent as an equivalent heat flux, which has subsequently been applied to glider observations by Thompson et al. (2016), du Plessis et al. (2017), and Viglione et al. (2018). Term , which is dependent on the horizontal buoyancy gradient and the MLD, is determined as
where H is the MLD. We acknowledge that both and provide a sense of how much heat would be required in the mixed layer to arrive at a similar restratification or mixing and not the diabatic process which occurs due to the response to surface heat flux. 0.06 is an empirically defined coefficient determined by numerical models (Fox-Kemper et al. 2008). We accept this may not be a true representation for the Southern Ocean but is currently our best estimate available. Note that always acts as a restratification flux. For both Eqs. (7) and (8), is averaged over the MLD. The strength of the restratification from MLEs is not uniform throughout the mixed layer, rather it is maximum in the middle of the mixed layer and decreases to zero at the surface and bottom of the mixed layer (Fox-Kemper et al. 2008). For Eq. (8), we estimate the maximum restratification occurring within the mixed layer for each buoyancy flux parameterization. Furthermore, we assume that the depth over which the Ekman advection occurs (Ekman layer) is equivalent to the MLD and thus horizontal advection by Ekman transport generates mixing at the base of the mixed layer.
e. Additional datasets/reanalysis products
Both wind speed and direction are important variables in this study. To obtain collocated wind stress and direction, we use the data from NCEP-2 (https://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanalysis2.html). The NCEP-2 wind stress was compared to in situ observations from Wave Glider deployments at the SOSCEx location, providing the highest correlation to the in situ wind measurements compared with other gridded wind products (Schmidt et al. 2017). Thomson et al. (2018) use Wave Glider measurements of wind direction near the Antarctic Peninsula to show that collocated NCEP-2 wind direction successfully mirrors the in situ observations. The temporal resolution of the NCEP-2 wind product is 6 h. In addition to wind data, we use NCEP-2 for surface heat fluxes (solar, net longwave, latent and sensible) and precipitation.
f. Model description
The Price–Weller–Pinkel (PWP; Price et al. 1986) bulk mixed layer model is used as a diagnostic tool to elucidate the role of one-dimensional mixing and restratification processes. PWP applies a momentum flux induced by winds, which along with cooling and evaporation contribute to the three types of mixing: (i) convective instability, (ii) entrainment from the pycnocline, and (iii) mixing through enhanced vertical current shear. Shortwave radiation is input at the surface and absorbed into the profile with a double exponential depth dependence. The water column restratifies when the buoyancy flux at the surface is positive, for example, through heating and precipitation. The surface net heat flux constitutes shortwave and longwave radiation, as well as latent and sensible heat fluxes. The model is initialized with a mean of the first 20 profiles for each experiment and forced with parameters of and τ that are initially collocated in space and time to the glider surfacing location and then interpolated to the grid of the optimal interpolated glider coordinates. The density difference criteria = 0.03 kg m−3 defines the MLD for each time step (de Boyer Montégut et al. 2004). The simulations have a 1-m depth resolution and 2-h temporal resolution to reflect the resolution of the optimally interpolated data. We do not employ a background diffusivity in our simulations.
PWP simulations using the above criteria are the one-dimensional simulations (PWP1D). We repeat the four simulations described above with the addition of and as equivalent heat fluxes. These fluxes are treated the same way as and are applied at each model time step equally throughout the mixed layer. We refer to these simulations as PWPSM. Direct comparisons of the seasonal evolution of stratification between PWP1D and PWPSM are used to diagnose the role of MLEs and EBF impacting the seasonal development of mixed layer stratification.
a. Mixed layer seasonality
Multi-month (ranging from midwinter to late summer) glider deployments over four separate years show an evident seasonal evolution of mixed layer temperature and salinity (Figs. 4a,c). The most substantial increase in mixed layer temperature from winter to summer exhibits a range of 4.7°C (6.7°–11.4°C). This warming occurs over 147 days between August 2015 and January 2016, equating to a daily mean surface heat flux of 146 W m−2 over a mixed layer of 100 m. In contrast, a seasonal freshening of the mixed layer from December to late summer (most prominent in SOSCEx2 and SOSCEx3) displays a salinity decrease from around 34.3 to 34.1, equating to an approximately 1°C change in temperature. This freshening of the mixed layer is likely to be a signature of an equatorward freshwater flux driven by seasonal ice melt (Haumann et al. 2016). Embedded within the seasonal cycle of thermohaline variability are density-compensating features (1°C, 0.2 psu) occurring over the order of a day (see Figs. 4c,d). We apply the density ratio to quantify the relative effect of horizontal variations of mixed layer temperature and salinity on density,
where α and β are the thermal expansion and haline contraction coefficients, and the horizontal differences of temperature and salinity, and , are taken across the spatial interval between consecutive glider dives. For values of > 1, horizontal gradients of mixed layer temperature contribute to larger changes in density than salinity. Values whereby < 1 represent the opposite. The median value of monthly compositions of are consistently above 1 for all glider experiments (Fig. 5), indicating the tendency for horizontal gradients in temperature to impact density fronts over salinity gradients. Our observations of corroborate with early work by Stommel (1993) who proposed a density ratio of 2 based on regionality of the thermal and haline atmospheric forcing of the upper ocean. However, we observe that exhibits seasonality, with the lowest values that are generally closest to 1 in winter (defined here as July–September) when the surface ocean undergoes atmospheric cooling. Monthly median values of increase from late winter, reaching a maximum of = 3.6 during November–December before decreasing again in late summer. The seasonality of compensation is likely a response to the mixing and distribution of heat within the mixed layer. Rudnick and Martin (2002) show that a stronger density compensation exists during deep mixed layers, where mixing distributes temperature and salinity vertically and horizontally. Meanwhile, the ocean tends not to be compensated during shallow mixed layers. Increasing mixed layer temperatures due to heating during summer ultimately leads to a seasonal shoaling of the mixed layer, and therefore likely allows for less vertical mixing and horizontal temperature gradients to propagate in the mixed layer. Differences in the timing of the maximum median between the four years studied suggest that mixing/restratifying processes that allow for density compensation are variable on interannual time scales.
A composite of the mixed layer horizontal buoyancy gradients from all four glider experiments shows a seasonal signal where the lowest horizontal buoyancy gradients occur in winter and highest in summer. We represent the underestimation in the horizontal buoyancy gradient given that, in the mean, gliders underestimate the true front gradient by 64%. The upper limit of the winter horizontal buoyancy gradients are lower than the spring and summer horizontal buoyancy gradients observed from the glider. Thus, we are confident that the seasonality of mixed layer fronts seen by the glider exists. During strong thermohaline compensation (winter), only 3% of the horizontal buoyancy gradients exceed 10−7 s−2 (Fig. 6). Meanwhile, during spring (October–November) and summer (December–March) when density compensation breaks down and temperature fronts dominate mixed layer density fronts, mixed layer horizontal buoyancy gradients exceed 10−7 s−2 during 12% and 13% of the profiles, respectively.
b. Seasonality of summer restratification
Figure 7 shows the temporal evolution of the upper-ocean buoyancy frequency and MLD for all SOSCEx studies. The MLD reaches a maximum of around 220 m during August, while the shallowest mixed layers are above 100 m during late November and early December. Over subseasonal scales, episodes of mixed layer restratification are signified by a rapid shallowing of the MLD (~50 m day−1), which occur via the formation of new stratification within the top 20 m of the ocean (N2 ~ 0.3 × 10−5 s−2; e.g., SOSCEx3 at the end of August and SOSCEx4 at the end of July). These restratification events occur 2–3 times per month and can remain for periods from one day up to a week.
The physical process by which the mixed layer undergoes restratification is through the emergence of a new pycnocline from the surface, which we refer to as the seasonal pycnocline. The seasonal pycnocline forms about 100 m above the winter pycnocline, thus creating a layering of stratification in the upper ocean. This layering is particularly evident during SOSCEx3 when the formation of the seasonal pycnocline in late November is superseded at the surface by a second seasonal pycnocline in January to generate three separate layers of stratification within the upper 300 m. A common feature of the seasonal pycnocline between the different glider experiments is after the initial formation at the surface; it gradually deepens to around 100 m over the period of a month. We use the metric of the seasonal pycnocline as the seasonal restratification, which allows us to separate the winter/spring (before) and summer (after) periods.
To objectively determine a seasonal restratification date, N2 is averaged from the surface to the depth of the 1026.75 kg m−3 isopycnal (H26.75). The H26.75 occurs within the winter pycnocline, and thus stratification above this depth during winter is low (Fig. 8). The seasonal restratification date is the first day when the mean N2 above H26.75, denoted here , remains continuously above 10−5 s−2 for the duration of the respective experiment. The = 10−5 s−2 represents the maximum stratification observed within the mixed layer during winter. Thus, we argue that increasing , be it due to atmospheric forcing or internal ocean dynamics such baroclinic instability, will define the restratification of the mixed layer (Fig. 8). By this definition, the seasonal restratification dates are only obtainable for the latter three glider experiments (SOSCEx2 to SOSCEx4) as restratification of the mixed layer occurred before the glider deployment for SOSCEx1. The dates of seasonal restratification, shown as vertical bars in Fig. 8a, for SOSCEx2 to SOSCEx4 are 22 October, 26 November, and 14 November, respectively, indicating a 36-day interannual difference in the timing of mixed layer restratification.
The onset date of net positive surface heat flux occurs during late September, corresponding to a delay in mixed layer restratification by up to 1–2 months after (Fig. 8b). Interestingly, from mid-October, the does reveal a response by increasing to above the maximum winter stratification of = 10−5 s−2. However, subsequently erodes to <10−5 s−2, thereby arresting the seasonal restratification.
c. Potential vorticity structure
We now consider sections of glider-derived PV [Eq. (6)] to investigate the seasonality of flow stability in the upper ocean (Fig. 9). Negative PV dominates at the mixed layer pycnocline throughout all seasons, indicating stable flow. However, a seasonality exists whereby weak PV during the winter ( ~ −0.5 × 10−8 s−3) become more negative ( < −1 × 10−8 s−3) in response to the formation of the seasonal pycnocline during summer. Interannual differences in the magnitude of PV at the base of the mixed layer indicate that slightly less negative values during SOSCEx4 ( from ~−0.1 to −0.3 × 10−8 s−3) compared to SOSCEx3 ( ~ −0.5 × 10−8 s−3) are associated with mixed layers that are deeper by about 40 m. Throughout the experiments, the mixed layer itself is prominent with PV that regularly take the opposite sign of f, suggesting the potential for mixed layer instabilities to occur regularly throughout all seasons. We note that baroclinic instabilities do not require , and therefore their presence is not limited to the instances discussed here. We find that between 9% (SOSCEx1) and 26% (SOSCEx4) of the instances of , the vertical stratification was stable. These intermittent periods whereby the baroclinic component of PV generates the condition for instabilities reveals that the mixed layer does experience conditions for other classes of instabilities, such as symmetric and centrifugal.
d. Submesoscale instabilities: Wind–front interactions and mixed layer baroclinic instabilities
The estimation of EBF [Eq. (7)] requires knowledge of the wind-front alignment. For this analysis, the direction of the depth-averaged current acquired from the glider dive cycle represents the direction of the mixed layer front (Fig. 10). The frontal flow direction is predominantly between 45° and 90°, occurring between 33% and 42% over the four experiments. Eastward frontal flow (between 0° and 180°) is observed between 81% (SOSCEx1) and 94% (SOSCEx3), meaning that flow reversals toward the west are most often found during SOSCEx1 (19%) and least often in SOSCEx3 (6%).
The wind direction is strongly dominant toward the east (89% during SOSCEx3 and 4 and 95% during SOSCEx1; Fig. 11). In particular, the wind is predominantly toward the east and southeast (90°–135°), accounting for 42% and 31% of the wind direction during SOSCEx1 and SOSCEx4. SOSCEx3 and SOSCEx4 experienced the most westward wind reversals (11%). The coherent alignment of westerly winds and frontal direction toward the east promotes the occurrence of down-front winds.
Calculating requires both the mixed layer horizontal buoyancy gradient and the component of the wind stress aligned with the front [Eq. (7)]. The wind stress is collocated to each glider profile to provide time series of wind stress values for each experiment. The angle difference between the front and wind direction for each glider profile then provides the alongfront wind stress component (Fig. 12). The wind-front alignment is predominantly down-front, ranging between 81% and 85% of the time. Notably, there is a positive skewness in the distribution for positive (down-front) wind stress. On average across all the experiments, down-front winds stress exceeding 0.2 N m−2 occurs 22% of the time, while only 1% for the up-front winds.
Equivalent heat flux estimates of and are compared directly to the collocated surface heat fluxes (Fig. 13). The datasets which contain the wintertime series (SOSCEx3 and SOSCEx4) reveal the seasonal cycle of , where late winter cooling gradually changes sign, alternating between cooling (~−200 W m−2) and warming (~200 W m−2) throughout August and September. The temporal variability of and appear as intermittent spikes, whereby the magnitude of regularly exceeds −500 W m−2 throughout spring and summer. The is comparatively weaker than and does not surpass 500 W m−2. Furthermore, undergoes sustained periods (>1 week) of negative buoyancy flux (−500 W m−2), which exceeds the amplitude of and can thereby change the sign of the net buoyancy flux. Integrating the three fluxes ( + + ) across each experiment reveals that the contribution of submesoscale equivalent heat fluxes is to reduce the warming provided by alone by 53%, 41%, and 58% for the first three experiments, respectively. For SOSCEx4, the equivalent cooling flux by results in a net cooling of the mixed layer across the study period.
e. Model comparison
The mean stratification above the 1026.75 kg m−3 isopycnal () is computed for both PWP model runs (PWP1D and PWPSM) and directly compared to the glider observations (Fig. 14). The seasonal evolution of from both model simulations and observations increase from winter to late summer. However, at the end of the four PWP1D experiments approaches or exceeds 10−4 s−2, between 4 × 10−5 and 7.5 × 10−5 s−2 larger than the glider . Meanwhile, values of PWPSM exceed the observations by between only 0.2 × 10−5 and 2.7 × 10−5 s−2. Although the evolution of in response to from September is to increase in both simulations as well as the data, it is during this period where the over stratification by PWP1D begins. Overall, in PWPSM appears to remain low relative to PWP1D. The main contribution of EBF is to reduce the mean difference of between PWPSM and the data by between 38% and 88% when compared to the mean difference between PWP1D and the data. What is notably important is that during October and November, when the PWP1D begins to ramp significantly, remains low in the PWPSM simulation and the data.
a. Seasonal cycle of the Subantarctic mixed layer
This study investigates interannual variations of the subseasonal evolution of stratification using four seasonal cycles of upper-ocean glider data from the Subantarctic Zone of the Southern Ocean. Evidence of interannual variability in the timing of seasonal mixed layer restratification exists. Two of the four experiments encompass the austral winter (August and September) when atmospheric cooling promotes convective instabilities and deep mixed layers (Fig. 7). During this time, periodic events of mixed layer restratification occur across the order of a day, synonymous with the time scale of restratification by baroclinic instabilities (Boccaletti et al. 2007).
The magnitude of the mixed layer horizontal buoyancy gradients are generally an order of magnitude weaker than regions of the Southern Ocean preconditioned for strong mesoscale eddy fields (Viglione et al. 2018) and topographical influence (Rosso et al. 2014), but are comparable with the open ocean conditions of the North Atlantic (Thompson et al. 2016). The horizontal buoyancy gradients undergo a seasonal cycle, where weaker gradients occur during the winter months, contrasting the observations from Callies et al. (2015) and Thompson et al. (2016), where horizontal buoyancy gradients are stronger during winter. We associate the strengthening horizontal buoyancy gradients during summer with the seasonal warming and subsequent increase in the contribution of horizontal temperature gradients to density fronts.
Calculations of PV using gliders reveal a seasonality in the role of PV in the upper ocean. We show that PV is weak in the winter pycnocline, strengthening during the summer. The weak PV at the base of winter mixed layer allows for deeper mixing, as indicated by a small reduction in the strength of the stable PV layer during SOSCEx4 relative to the year before resulting in mixed layers deeper by around 40 m. The interannual difference (SOSCEx3 and SOSCEx4) between the PV layer is a result of weakened vertical stratification, which suggests that these differences may translate to variability in the vertical transfer of properties at the base of the mixed layer (Erickson and Thompson 2018). This may be an important consideration for biogeochemical dynamics considering the importance of the upward vertical flux of iron into the mixed layer during winter in the Southern Ocean (Tagliabue et al. 2014). Furthermore, the weakening of the vertical component of PV during winter may allow for deeper mixed layers (as is observed in this study), which can increase the potential energy of the mixed layer when lateral buoyancy gradients are present, and possibly enhance the baroclinic component of PV (e.g., Thomas et al. 2013). We find the mixed layer to be predominantly susceptible to gravitational instabilities, although we do see evidence for the baroclinic component of PV to reverse the sign of PV to the opposite of f. This is a key finding as a number of studies in the Southern Ocean have shown that symmetric instabilities can arise from instances where , however, these studies occur in regions preconditioned for submesoscale instabilities, such as downstream of the Shackleton Fracture Zone (Viglione et al. 2018), downstream of the Kerguelen Plateau (Rosso et al. 2015), or on the edges of mesoscale eddies (Adams et al. 2017). Our observations show that although limited, these instabilities may be present in the open-ocean Southern Ocean. Full seasonal cycle measurements and a focused study on symmetric, centrifugal, and gravitational instabilities are required to understand their occurrence and seasonality further.
The seasonal mixed layer restratification occurs through the emergence of the seasonal pycnocline from the surface. Our observations indicate that a requirement for seasonal restratification is a positive surface heat flux, consistent with previous observations (Dong et al. 2008; Sallée et al. 2010). We define the onset date of seasonal restratification as a continued increase of the mean stratification above H26.75. This increased stratification is due to either mixed layer waters getting lighter by heating or freshening and thereby increasing the vertical density gradient, or horizontal advective processes such as the slumping of isopycnals due to baroclinic instability. We find that the restratification is likely to be a combination of both processes (not shown). In addition to the restratification mechanisms, we observe transient mixing events throughout spring and as a result, the timing of restratification becomes highly variable between different years (up to 36 days). This arrest of seasonal restratification by the mixing may allow for a prolonged vertical exchange of properties between the mixed layer and below, which may directly influence mixed layer heat budget estimations (Dong et al. 2007).
Furthermore, Swart et al. (2015) show that seasonal restratification directly results in a bloom of biological activity when phytoplankton growth is light limited, as is the case in the SAZ. Thomalla et al. (2011) indicate the presence of spatial heterogeneity of phytoplankton bloom initiation dates in the Southern Ocean. We consider that interannual variability of mixed layer restratification observed here may partially be responsible for these discrepancies.
b. Submesoscale impacts on seasonal restratification
Parameterizations of MLE and EBF require information of the mixed layer horizontal buoyancy gradient, the MLD and the alongfront wind stress. The horizontal buoyancy gradient and MLD are calculated directly from glider measurements. The alongfront wind stress is obtained using the frontal direction inferred by the depth-averaged current and the wind direction. The consistent eastward alignment of winds and upper-ocean flow in the SAZ is indicative of a down-front dominant regime. The propagation of cyclonic storms in the Southern Ocean (Yuan et al. 2009) is associated with periods of 4–10 days in the SAZ (Swart et al. 2015). Estimates of a negative buoyancy flux by EBF suggest that the westerly winds drive enhanced down-front Ekman flow, which manifests as enhanced gravitational mixing exceeding the buoyancy input by a positive surface heat flux. In contrast, calculations of MLE for our experiments do not provide a significant contribution to the upper-ocean buoyancy flux compared to surface heat flux and EBF. We recognize this may be the result of (i) the glider not sampling the fronts perpendicularly and thus underestimating the magnitude of the horizontal buoyancy gradient; (ii) baroclinic instabilities spin off the mesoscale horizontal buoyancy gradient, which is unable to be determined as the gliders remain in a localized region; and (iii) relatively shallow winter mixed layers (150–200 m) and weak horizontal buoyancy gradients (order 10−7 s−2) compared to other regions where submesoscales are shown to be active (200–300 m and order 10−6 s−2 in the Drake Passage; Viglione et al. 2018).
By parameterizing EBF and MLE as buoyancy fluxes into the PWP one-dimensional mixed layer model, the seasonal evolution of stratification within the mixed layer dramatically improves compared to the model run with surface heat and freshwater fluxes alone. An important consideration for the model analysis is that the glider underestimates the magnitude of the true mixed layer fronts thus reducing the potential contribution of and to the extent and variability of the mixed layer stratification and its evolution. We are aware there are other processes not included in the model which also impact the stratification of the mixed layer such as ocean wave effects (e.g., Li et al. 2016) and that uncertainties in the surface forcing will likely impact our estimates of stratification as well. Therefore, we emphasize that the evolution of the mixed layer stratification in PWPSM provides a statistical evaluation of the evolution of the mixed layer in response to MLE and EBF, and not all mixing and restratifying processes in the ocean. Additionally, the equivalent fluxes provide a sense of how much heat or cooling is required to arrive at a similar restratification or mixing, and indicate that a one-dimensional perspective on mixed layer dynamics in the Subantarctic is insufficient to explain the seasonal timing of restratification. We provide evidence that parameterizations for MLEs and EBF are required to represent the stratification in the Southern Ocean adequately.
A similar analysis using the Monthly Isopycnal/Mixed Layer Ocean Climatology (MIMOC) provides global maps of and for the month before seasonal restratification of the mixed layer (Johnson et al. 2016). They show that for the month before in the Southern Ocean, provides a mean restratifying flux of around 30 W m−2, while destratifies the upper ocean with a flux around −20 W m−2. They provide the caveat that their large-scale monthly climatologies of wind and density gradients likely misrepresent the impact of EBF at localized fronts. By averaging the submesoscale equivalent heat fluxes from the glider time series over the same period for SOSCEx3 and SOSCEx4, we find similar values (26 and 16 W m−2), but larger (−83 and −105 W m−2), revealing that the lower-resolution climatological maps represent the restratification by submesoscale eddies, but are not able to capture the sizeable destratifying flux by EBF.
The role of down-front wind mixing in other regions of the global ocean (D’Asaro et al. 2011) shows an enhancement in the rate of energy dissipation of the upper ocean by an order of magnitude. Although these estimates represent a region of strong mesoscale frontal activity, we show that submesoscale horizontal buoyancy gradients are ubiquitous in the open ocean Southern Ocean, despite exhibiting weaker horizontal buoyancy gradients compared to regions preconditioned for submesoscale activity. We show that the ubiquity of horizontal buoyancy gradients in the Subantarctic responds to a predominantly down-front Southern Ocean wind field, leading to episodic EBF-induced mixing. The wind-front alignment observed in this study is not persistent for all regions in the Southern Ocean, where topographical features may steer the flow and therefore periodically misalign the frontal flow and the wind field. For these regions, the impact of EBF mixing may be reduced (e.g., Viglione et al. 2018). The favorable wind-front alignment in the Subantarctic may contribute to enhancing turbulence in the mixed layer in addition to other mixing processes such as shear-driven mixing and Langmuir turbulence (Fan and Griffies 2014; Li et al. 2016). We show that EBF may be an important mixing process to the synoptic modulation of the SAZ mixed layer (Nicholson et al. 2016).
Furthermore, the delay of seasonal restratification is likely to result in interannual variability of phytoplankton bloom initiation dates and general bloom heterogeneity observed in the SAZ by Swart et al. (2015) and elsewhere in the Southern Ocean by Thomalla et al. (2011) and Carranza and Gille (2015). Three-dimensional processes, and in particular MLEs, which directly impact the winter to spring restratification of the mixed layer are becoming studied more frequently (Mahadevan et al. 2012; Johnson et al. 2016; du Plessis et al. 2017). Mahadevan et al. (2012) argue that a significant shift in the timing of the spring bloom occurs due to the onset of restratification by MLEs before . Furthermore, du Plessis et al. (2017) propose that MLEs promote an increase of the stratification during spring which cannot be explained by surface heating alone. Our results indicate that in addition to the restratification by MLEs, a destratifying flux by EBF may also provide a shift in the timing of seasonal restratification. We propose that EBF in the Southern Ocean can be dominant in regions where there are alongfront wind alignment and strong wind stress. In these regions, we suggest that EBF may also be an essential process determining the onset of springtime phytoplankton blooms, providing an exciting avenue for future studies to explore.
This work forms a part of a growing body of literature, which continues to show the presence and importance of submesoscale processes in the Southern Ocean (e.g., Rosso et al. 2014; Swart et al. 2015; Rocha et al. 2016; Adams et al. 2017; Bachman et al. 2017; Erickson et al. 2016; du Plessis et al. 2017; Viglione et al. 2018). Our observations have shown that for climate models to correctly simulate the seasonal restratification, the wind direction and fronts in the Southern Ocean need to be adequately represented. A further step in improving this field would be to distinguish the discrepancies in the distributions of critical submesoscale parameters sampled when using various glider sampling patterns. Furthermore, obtaining an understanding of the relative importance of EBF across the SAZ, or even the entire Southern Ocean, would be useful going forward.
Over four separate years, ocean gliders were deployed in the Subantarctic Zone of the Southern Ocean to investigate the subseasonal and interseasonal variability of mixed layer stratification. Observational studies, which elucidate the role of submesoscale motions in the Southern Ocean are rare, while those that cover multiple consecutive seasons have not previously existed. The datasets presented here range between winter and late summer, capturing the transition from deep winter mixed layers to strongly stratified and shallow summer mixed layers. From these valuable datasets the major conclusions are as follows:
Horizontal fronts within the mixed layer exhibit strong seasonality and are driven primarily by changes in temperature, which exhibits the most substantial influence in early summer. Mixed layer buoyancy gradients are lowest in winter and highest in summer.
Winter-to-summer glider time series shows that the restratification of the mixed layer can occur up to 2 months after the onset of seasonal surface heat flux warming. This is an important consideration given that restratification regulates the exchange of properties between the mixed layer and ocean interior as well as the vertical control of tracer properties important for biological production.
The magnitude of the estimated Ekman buoyancy flux is large enough to cause the observed delay in the onset of restratification. The conditions which promote EBF are dominant in the Southern Ocean—strong westerly winds promote down-front conditions, which interact with an ocean substrate of prevalent horizontal buoyancy gradients induced by fronts, meandering jets, and eddies.
The net effect of EBF on mixed layer stratification is to dampen the seasonal evolution of restratification by about half the amplitude of that generated by the surface heat flux alone using simulations of a one-dimensional bulk mixing model.
This study has shown that the combination of the dominantly westerly winds of the Southern Ocean and submesoscale motions may enhance the periodic input of energetic vertical motions and directly impact the production of upper-ocean biomass. Therefore, we propose that the effect of submesoscale processes need to be considered when constraining global climate models. The intermittency of these mixing events suggests that this may be difficult to incorporate. Enhanced mixing by EBF may explain part of the inadequacies to represent the Southern Ocean MLD in GCMs accurately.
MdP acknowledges numerous research visits to the Department of Marine Science, University of Gothenburg, and a visit to Woods Hole Oceanographic Institution, which greatly enhanced this work. We thank SANAP and the captain and crew of the S.A. Agulhas II for their assistance in the deployment and retrieval of the gliders. We acknowledge the work of SAMERC-STS for housing, managing, and piloting the gliders. SS was supported by NRF-SANAP Grant SNA14071475720 and a Wallenberg Academy Fellowship (WAF 2015.0186). Lastly, SS thanks the numerous technical assistance, advice, and IOP hosting provided by Geoff Shilling and Craig Lee of the Applied Physics Laboratory, University of Washington.