Episodic Extrema of Surface Stress Energy Input to the Western Arctic Ocean Contributed to Step Changes of Freshwater Content in the Beaufort Gyre

The recent dramatic decline of sea ice in the western Arctic Ocean changes the transfer of momentum across the ice‐ocean boundary layer. The surface stress energy input through the surface geostrophic current in the Beaufort Gyre (BG) based on a numerical model is 0.03 mW/m2 in 1992–2004 versus 0.23 mW/m2 in 2005–2017. This energy input is primarily concentrated over the southern Canada Basin and the Chukchi Sea. It is 1.38 × 1016 J in observations versus 4.90 × 1016 J in the model in the BG during 2003–2014. We find that some well‐known freshwater changes in the BG over 1992–2017 resulted from episodic extrema of energy input in 2007, 2012, and 2016. In particular, most of the energy input in 2007 was transformed into potential energy (57%) which resulted in a new state of freshwater budget. Our study suggests that as of 2016, the BG had not yet reached a saturated freshwater state.


Introduction
The Beaufort Gyre (BG) is a wind-driven large-scale ocean circulation in the western Arctic Ocean. Currently, it contains about 23,300 km 3 of freshwater in the upper ocean (Proshutinsky & Krishfield, 2019) and is the largest oceanic freshwater reservoir in the northern hemisphere (Carmack et al., 2008). The liquid freshwater budget in the BG showed a dramatic increase in the 2000s (e.g., Giles et al., 2012;Krishfield et al., 2014;Morison et al., 2012;Proshutinsky et al., 2009;Rabe et al., 2014). An interesting phenomenon is that since 2010, the increasing freshwater content in the BG region seemed to be stabilized (e.g., Krishfield et al., 2014;Zhang et al., 2016), while at the same time, the surface forcing shows a weakening trend (Petty et al., 2016). Theoretical studies suggested that lateral eddy flux is the key factor to balance the wind-driven Ekman pumping and thus the stabilization of BG Yang et al., 2016), while other studies have considered the important role of geostrophic current in mediating the ice-ocean stress which also contributes to the equilibration of BG (Dewey et al., 2018;Meneghello, Marshall, Timmermans, & Scott, 2018;Zhong et al., 2018). Further, Meneghello, Marshall, Campin, et al., 2018) and Doddridge et al. (2019) noted that the ice-ocean stress feedback is more important than the lateral eddy flux for the equilibration of BG on interannual timescales.
As the BG is a wind-driven system, here we try to address the freshwater issues from a different point of view by investigating its mechanical energy partition. In the ice-free regions, wind energy input through the surface geostrophic current is formulated as: W ¼ τ ⇀ air−ocean ⋅u ⇀ g (e.g., Wunsch, 1998;Ferrari & Wunsch, 2009;Zhai et al., 2012;Zhai, 2013), where the air-ocean stress τ air density ρ air = 1.25kg m −3 , the air-ocean drag coefficientC d = 0.00125, and U ⇀ 10 is the wind speed at 10 m above the sea surface. However, in the BG region where it is partially covered by sea ice, the total surface stress is considered, which is the combination of air-ocean stress and ice-ocean stress scaled by sea ice concentration, that is, τ (here the ocean surface velocity of~5 cm/s is negligible comparing with the wind velocity of~5 m/s or larger), the ice-ocean stress τ , the ocean density ρ ocean = 1023kg m −3 , the ice-ocean drag coefficient C iw = 0.0055, and R o is rotation matrices for ocean (McPhee, 1975(McPhee, , 1980. Thus, it is the surface stress energy (SSE) input through the surface geostrophic current should be considered in the BG as: W ¼ τ ⇀ total ⋅u ⇀ g , which we denoted as SSE hereafter. Previous studies focusing on the global ocean instead of the Arctic revealed that when the surface current is non-negligible comparing with the wind speed, neglecting the ocean surface current would result in a positive bias in calculations of wind energy input (about 20%-30%) (e.g., Duhaut & Straub, 2006;Hughes & Wilson, 2008;Scott & Xu, 2009). The inclusion of ocean surface current is also important when considering the Ekman dynamics in the BG (Dewey et al., 2018;Meneghello, Marshall, Timmermans, & Scott, 2018;Zhong et al., 2018). However, so far as we know, there are still no studies considering the surface stress energy input in the BG region including the effect of increasing geostrophic current. The basic objective of this study is to investigate the interannual variability of surface stress energy input in the BG region, as well as its role in the interannual variability of freshwater content.
The majority of the wind energy input to surface waves and Ekman layer will be lost in the surface layer and will not enter the deeper ocean (e.g., Huang et al., 2006;von Storch et al., 2007;Zhai et al., 2009). How the surface stress energy input is transferred to the deeper BG interior has yet to be addressed. (By saying deeper BG interior, we refer to the depth that deeper than the Ekman layer.) And its linkage with the accumulation of freshwater in the BG is still unclear. Here, we will address these questions by studying the surface stress energy input and mechanical energy partition in the BG. We used a numerical model and an observational method to quantify the surface stress energy input in the BG. The next section describes the model and observational data that we used. Section 3.1 gives the spatial distribution and temporal variability of SSE. We then evaluate the model's performance in simulating the freshwater content variability of the BG in section 3.2. In section 3.3, the SSE is linked with the mechanical energy partition in the BG and the implication for the freshwater content variability. Our major findings and discussion are presented in section 4.

Model and Observations
The monthly mean-derived surface stress τ ⇀ total and surface geostrophic current u ⇀ g are used to calculate SSE as: W ¼ τ ⇀ total ⋅u ⇀ g . Note that although the SSE (low-frequency energy) considered here is only a portion of the total surface stress energy input (including both low-and high-frequency energy) to the ocean, it plays a major role in regulating the large-scale ocean circulation like the BG (e.g., Huang et al., 2006;Zhai et al., 2012). We focus here on monthly mean SSE (this is also the timescale of the observed geostrophic currents) to study the interannual variability of BG, rather than the more high-frequency energy input (e.g., energy input to surface waves, Ekman layer, and near-inertial current) which will mostly dissipate in the surface Ekman layer. We also used the observed geostrophic currents via dynamic ocean topography data (Armitage et al., 2016(Armitage et al., , 2017 to calculate the observation-based SSE, which is denoted as "OBSm" here. We refer readers to Zhong et al. (2018) for the details of surface stress calculation using observations.
In order to extend the study period to 1992-2017, model results from the Marginal Ice Zone Modeling and Assimilation System (MIZMAS) are used. MIZMAS is a coupled ice-ocean model that assimilates satellite observations of sea ice concentration and sea surface temperature. The ocean model was originally based on the Parallel Ocean Program (Smith et al., 1992) and modified by Zhang and Steele (2007). More detailed information about MIZMAS model components, domain, and grid configuration can be found in Zhang et al. (2016), where extensive validation against in situ observations in the BG was also presented. A caveat is that the model we used here has an average horizontal resolution in the Chukchi and Beaufort Seas of~10 km and is thus not fully eddy-resolving. However, the model is sufficient to study the large-scale circulations such as the BG (~10 3 km) with a timescale longer than months. (This had been validated in Zhang et al. (2016).) A caveat is that our model is not an eddy-resolving model, so that the eddy kinetic energy (KE) or eddy dissipation might be underestimated. In the model, the surface geostrophic current is derived from the model sea surface height. The model has a vertical resolution of 5 m over the upper 80 m. Both the observed method and model are using the same National Center for Environmental Prediction/National Center for Atmospheric Research reanalysis.

The Surface Stress Energy Input Through the Surface Geostrophic Current in the BG Region
Figures 1a and 1b show that the surface stress energy input through the surface geostrophic current estimated from OBSm and the model are generally similar both in the spatial distribution and magnitude. Significant and strong SSE appears in the northern Chukchi Sea and the southern Canada Basin which indicates strong mechanical energy input into the ocean. The strong SSE in the southern Canada Basin acts like an engine to speed up the flywheel of BG (Proshutinsky et al., 2002), while the strong SSE along the Chukchi Slope and the northern Chukchi Sea might contribute to a stronger Chukchi slope current (Corlett & Pickart, 2017;Li et al., 2019;Spall et al., 2018). In contrast, negative SSE appears in the southern Chukchi Sea and the northern Canada Basin, which indicates areas where the surface stress represents a brake on oceanic KE. How can this happen? There are two possibilities: (i) relatively faster ocean currents rubbing against slower sea ice (Dewey et al., 2018;Doddridge et al., 2019;Meneghello, Marshall, Campin, et al., 2018) or (ii) wind stress forcing that opposes the northward inflow of Pacific Water in the southern Chukchi Sea. The different sign of SSE in the southern versus the northern Canada Basin indicates that mechanical energy goes into the southern limb of BG, while it is lost in the northern limb of BG. The net effect of this inhomogeneous distribution of SSE might contribute to the asymmetric shape of the BG (Regan et al., 2019; Zhong et al., 2018). With regard to variability, there is high SSE standard deviation in the Chukchi Sea and the southern Canada Basin (Figure S1), where the surface geostrophic current shows a significant increasing trend (Armitage et al., 2017;Zhong et al., 2018) and is dominated by strong seasonal variability of SSE ( Figure S2). The north-south SSE difference in the Canada Basin is larger in the later period (Figures 1c  and 1d), relative to the earlier period. Significant increases of SSE also appear in the northern Chukchi Sea.
The MIZMAS model has a high correlation coefficient with the OBSm regarding the surface stress magnitude and surface geostrophic current both in the monthly and annual variability (Figures 2a and 2c). Both the surface stress and surface geostrophic current reveal a strong seasonal variability, with the maximum value usually appearing in fall ( Figure S3, also see Meneghello, Marshall, Timmermans, & Scott, 2018;Zhong et al., 2018). The surface geostrophic current shows a dramatic increase of~2 cm/s since 2007. The model surface geostrophic current is generally a bit larger than the observational results, while the model's surface stress magnitude is generally consistent with observations. No clear increasing trend is revealed in surface stress variance, although there is an exceptional peak in 2007 (Figures 2a and 2c). The SSE shows contrasting variabilities before 2002 and after 2002 (Figures 2b and 2d), with weak fluctuations (generally

Geophysical Research Letters
from −0.5 to 1 mW/m 2 ) before 2002 and strong fluctuations after (generally from −2 to 2 mW/m 2 ). Stronger fluctuations indicate the combined effects of faster and retreating sea ice and increasing geostrophic current (e.g., Petty et al., 2016;Zhang et al., 2016). The 2007 SSE peak (with~9 mW/m 2 for monthly mean and 1.5 mW/m 2 for annual mean) is consistent with the peak of surface stress magnitude and surface geostrophic current (Figures 2a and 2b). The SSE from model is generally larger than that from the OBSm. This may be either due to an overestimation of the surface geostrophic current by the model or due to an underestimation from observations with coarse resolution in space and time. However, the interannual variability of the SSE from the MIZMAS model is similar to that of the observational derived SSE over 1992-2017. This indicates that we can evaluate the interannual variability of SSE from the model.
The most significant increase of SSE appears in fall (November) followed by winter (February) during 2005-2017 compared with 1992-2004 ( Figures S2 and S3). This indicates that the retreat of sea ice in the southern Canada Basin during these seasons contributes to a larger surface stress energy input into the ocean. There is also a moderate drop of SSE during spring (May) in the second period, which suggests that the surface stress takes the energy out of the ocean when the ocean drives the ice (Dewey et al., 2018).

The Step Change of Freshwater Content in the BG
In this section, we evaluate the model's robustness in simulating the freshwater accumulation in the BG. Figure 3 shows that the model-simulated freshwater content has a very high correlation coefficient with the observations on interannual timescales, despite the model's freshwater content bias. Although the model did not capture the exact magnitude of freshwater content volume, it did capture the interannual variability of this wind-driven freshwater build up and release process. Both the model and observations indicate that the freshwater content dramatically increases in 2007 and 2008, followed by a slightly varying plateau of 4 years. The observations show a steep drop of freshwater content in 2013 which then gradually rebounds. In contrast, the model freshwater content seems to be more or less stabilized during 2008-2015 (with a moderate drop in 2013). Both the model and observation show that the freshwater content reaches maxima in 2016 and 2017. Overall, the model is able to represent the interannual variability of freshwater content in the BG, when compared to the observation. This provides us a basic condition to study the mechanical energy partition in the BG system in our following discussions. The well-known increase of freshwater content in 2007/2008 and the abrupt decrease of freshwater content in 2012/2013 were reported to connect with some extreme events in the BG region. Here we defined the step changes of freshwater content as the change of annual mean freshwater content volume that is larger than 10 3 km 3 (which is about half of the standard deviation of freshwater content volume). The first event was the historic 2007 sea ice minimum with a large freshwater availability (Wang et al., 2018), enhanced by strong Ekman convergence (Meneghello, Marshall, Timmermans, & Scott, 2018;Zhong et al., 2018) which led to a sharp increase of freshwater content in 2008. The second event was strong cyclonic winds (storms) that appeared in the Canada Basin in summer of 2012, which limited Ekman convergence in the surface freshwater layer, resulting in a drop of freshwater content in the following year 2013 (Wang et al., 2018). These results are confirmed in Figure 3. In addition, our results reveal another step increase of freshwater content that appears in 2016-2017. These three events are each preconditioned by SSE extrema (Figures S4 and S5

The Interior Adjustment of Upper Ocean to the Episodic Extrema of Surface Stress Energy Input and Mechanical Energy Partition in the BG
In this section, we will address the mechanical energy partitions in the upper 80 m of the BG. Note that the typical Ekman layer depth in the Arctic Ocean is about 20 m (Cole et al., 2017). We focused on the upper 80 m to reduce the effects of the lateral volumetric injection of Pacific Water (e.g., Shimada et al., 2005;Zhong et al., 2019) since the Pacific Water also provides momentum flux to the BG. Note that the 80 m would eliminate most of the Pacific Winter Water but not the Pacific Summer Water. We employed the basic momentum equation and the rate of change of potential energy Φ to study the surface stress energy input and its mechanical energy partitions in the upper BG (the MIZMAS model results are used for this analysis), that is, , A h and A v are the horizontal and vertical viscosity, respectively. ∇p is the pressure gradient. After some manipulations by combining the above equations we have: , the potential energy (PE) ρΦ = ρgz, the area that is used for the integration on the right side of the equation is shown in Figure 1a. Because we are focusing on the upper 80 m, τ ⇀ bot is the interfacial shear stress (O'Brien & Hurlburt, 1972) between two vertical adjacent model grids at 80 m and 86 m: the divergence terms is related with the elevation or depression of ocean surface, −A h ∇ 2 E h is the lateral energy input from the Pacific Water injection (or other lateral water injection), and ρεis the dissipation of KE (the eddies' effect). To consider the large-scale BG circulation with a timescale longer than month, u ⇀ top is substituted by the surface geostrophic current. Equation (3) indicates that the surface stress energy input τ ⇀ top ⋅u ⇀ top is used for the following processes: (i) acceleration of the ocean current KE, (ii) enhancement of PE via collection of freshwater which changes the density distribution, (iii) dissipation by the interfacial shear stress τ ⇀ bot ⋅u ⇀ bot , and (iv) dissipation of KE ρεand the divergence terms in the residual term R * . We do not attempt to evaluate the residual terms in equation (3), because there are many uncertainties in these terms due to the model we used (especially the eddy-related terms). Instead, our study constitutes a first attempt to evaluate the changes of large-scale surface energy input and its fate in the ocean interior which is linked with the freshwater content variability. Figures 4a and 4b show that the KE from OBSm and from the model is generally consistent in the spatial distribution, although the model shows large KE and more spatial variability. The variability of the model surface KE is consistent with the observed derived KE, despite that the model values are generally larger (Figure 4c). The KE of BG is increased by~1 J/m 3 since 2007.
Because we are considering the mechanical energy input and partitions since 1992, the SSE and the interfacial shear stress work W int are integrated over 1992-2017 (Figure 4d). The change of potential energy ΔPE and kinetic energy ΔKE since 1992 are then analyzed, which both shows significant interannual variability (Figure 4d). The ΔPE is rather stable during 1992-1999, gradually increases during 2000-2006, increases steeply during 2007-2008, and then becomes stable again since 2010. The ΔPE and ΔKE also have a strong seasonal variability, with maxima in fall and minima in midyear. The maxima of ΔPE and ΔKE in fall are associated with the maxima of surface stress and surface geostrophic current during this season. Figure 4d shows that the change of potential energy ΔPE in BG is of order 10 2 -10 3 times the KE (although both of them increase or decrease together); this is consistent with the findings of Gill et al. (1974) who used a simple two-layer model. There was a large step change of the integrated SSE minus W int in 2007 (an increase of 3.38 × 10 16 J), which corresponds to a dramatic increase of ΔPE (1.94 × 10 16 J) during 2007-2008; thus 57% of SSE is transformed into ΔPE (Figure 4d). This indicates that most of the surface stress energy input to the BG is transformed into PE (of order 10 16 J), with the interfacial shear stress work W int (of order 10 12 J) being negligible. The ΔPE has a high correlation coefficient with the freshwater content volume, which indicates that the increase of freshwater content in BG is equivalent to the increase of PE (also see Polyakov et al., 2018). This is reasonable as the accumulation of freshwater steepens the isopycnal slope and thus increases the PE. In this way, the change of potential energy is a reflection of the change of freshwater content volume. The interannual variability of ΔPE could be roughly divided into four periods: 1992-1999 as the stable period, 2000-2006 as the period of gradually increase, 2007-2008 as the period of dramatic increase, and 2009-2017 as the new stable period but with higher ΔPE, which is consistent with periods defined by some previous studies (e.g., Zhang et al., 2016;Zhong et al., 2018). Because our model is not an eddy-resolving model, the energy partition to the eddy terms might change if we were to use a fully eddy-resolving model. This remains to be investigated. In addition, for equation (3), integrating over a thicker depth range results in a higher percentage of mechanical energy (PE + KE) partition (because the surface energy input penetrates deeper than 80 m). Thus, the percentage of 57% is a lower limit. But note that integrating over larger depth could induce a stronger lateral water injection which we try to avoid.

Discussion and Conclusions
For the first time, to our knowledge, we give a robust estimation of the surface stress energy input in the BG from a view of the large-scale gyre circulation. The interannual variability of freshwater budget in the BG depends on low-frequency energy input, which could be quantified by the monthly derived SSE. Before 2003, the SSE is rather small, with larger fluctuations thereafter. There are two significant peaks of SSE in 2007 and 2016 which correspond to step increases of freshwater content in the following years (i.e., 2008 and 2017). During these events, surface stress forcing is generally in the same direction as surface geostrophic currents in the BG, so that the ocean gains mechanical energy ( Figure S4). As a first attempt, our results show how much of the surface energy input had transformed into the PE which is equivalent to the freshwater content changes. The extreme of energy input in 2007 is the key contribution for the new state of high freshwater content since 2007 and persisted through 2017. The resulting high level of freshwater content may remain for many years, which Johnson et al. (2018) described as the decadal memory. Thus, the increase and maintenance of BG freshwater content are associated with episodic events of strong SSE. In contrast, the year 2012 is dominated by cyclonic wind and a negative minimum of SSE which results in a step decrease of freshwater in 2013. At this time, the direction of surface stress forcing is generally opposite to that of surface geostrophic current in the BG, so that the ocean loses mechanical energy. These results suggest that the change of freshwater content is a 1 year or less time-lagged response to episodic SSE change (with the understanding that the changes of freshwater content since 2008 are smaller positive or negative freshwater content anomalies (e.g., 2012 and 2016) on the new freshwater state). To some extent, this provides a useful way to predict future change of BG freshwater if we know the SSE of the previous year. The multiyear mean of SSE in the BG was 0.03 mW/m 2 in 1992-2004 versus 0.23 mW/m 2 in 2005-2017 ( Figure 2). This is a dramatic increase of mechanical energy input in the BG and has thus forced the BG into a "new normal" state with stronger geostrophic current and higher freshwater content.
The total integrated mechanical energy input was 1.38 × 10 16 J in observation versus 4.90 × 10 16 J in model in the BG during 2003-2014, while the MIZMAS model captures the interannual variability of the mechanical energy input very well (i.e., the change of SSE). MIZMAS results showed that most (57%) of the BG surface stress energy input (SSE) in 2007 is transformed into PE, with the rest dissipated through eddies and contributing to the divergence terms in the residual term of equation (3) which is related with the elevation or depression of ocean surface. In this case, the model can capture the interannual variability of freshwater content very well. Because of the strong baroclinic effects in the BG system , the topographic roughness for the energy dissipation of the upper layer is relatively small. On the other hand, some results  suggested that the topographic Rossby waves might be important for the energy dissipation of the BG. This remains to be further investigated. An interesting and important question is "How much freshwater does it take to saturate the BG?" Doddridge et al. (2019) proposed a three-way balance between the wind stress, the ice-ocean governor (the ice-ocean stress feedback), and the eddy fluxes. Ideally, the freshwater content volume should be determined by this three-way balance. The modelintegrated SSE reached a maximum in the early months of 2016, with no further net energy input into the BG after that time (Figure 4d). However, the PE (or freshwater content) continued to increase. This indicates that the peak of SSE in early 2016 is gradually transformed into the PE of BG which is also reflected in the historically high freshwater content in 2017. This suggests that by the year 2016, the BG has not yet reached a state of freshwater content saturation, although 2016 saw a strong SSE. The lagged response of PE to the surface energy input is a result of the inertia of the ocean. The interior ocean would take time to adjust to the changing surface forcing. Our results suggest that with 1-year surface forcing anomaly, the response of the interior ocean would manifest 1 year later. In addition, with an extreme surface energy input like that of 2007, the freshwater content of BG subjects to a decadal memory (Johnson et al., 2018).
For the present study, we investigated on a defined domain of BG and try to evaluate its freshwater state in this domain. However, for the real BG, it could expend its domain that allow more freshwater to be accumulated instead of developing into deeper ocean to reach the saturated state. It is very likely that the BG expands more easily on horizontal direction than to expand to deeper ocean to reach the saturated state. In fact, the observed dynamic ocean topography revealed that the BG is expanding (Regan et al., 2019). This cannot be simulated with those theoretical idealized model simulations (e.g., Doddridge et al., 2019;Yang et al., 2016) which the study region is restricted to a basin domain.
Results have showed that the Ekman pumping is relatively stable in recent years (e.g., Zhang et al., 2016;Zhong et al., 2018). And the area of Ekman downwelling is expanding (Figures 7i and 7l in Zhang et al., 2016). While it would take much more energy for the BG to pump the water down (stronger Ekman downwelling), the expanding Ekman downwelling area would favor the BG to converge the surface freshwater from a larger area without a solid boundary. All these lead to an expansion of BG in horizontal direction.
The increases of cyclone activities tend to decrease the surface stress energy input into the BG. That is because the cyclonic wind forcing opposes the anticyclonic BG circulation, thus taking energy out of the ocean and releasing freshwater that resides in the BG. There is an increasing trend of cyclone activity during winter in the Pacific side of the Arctic Ocean (Moore et al., 2018;Zahn et al., 2018), which might contribute to future releases of BG freshwater.