Waves are an important form of motion at the air-sea interface and are generated by wind acting directly on the sea surface. Wave breaking, Stokes drift and other phenomena associated with waves have an important impact on the mixed layer in the ocean (Longuet-Higgins et al. 1953; Hasselmann et al. 1970; Li et al. 1995; Lewis et al. 2004; Craik et al. 1997; Huang et al. 2008; Li et al. 2008; Deng et al. 2012; Zhang et al. 2014). Recently, increasing attention has been paid to the impact of Stokes drift on the upper ocean.
Stokes put forward the “Stokes drift” concept in 1847. Due to the existence of surface gravity waves, the nonlinear effect of Stokes drift leads to the non-closed trajectory of sea surface water points. If the periodic fluctuations of irrotational, inviscid and finite amplitudes are considered, the net displacement of waves in their propagation direction, i.e. Stokes drift, can vary with depth (Stokes 1847). When considering the seawater viscosity, the Stokes drift not only causes the net mass transport but also affects the average current at the sea surface (Longuet-Higgins et al. 1953; Kenyon et al. 1970). Weber et al. (1983) found that viscosity is present in the ocean and that no matter how small, the Stokes drift causes nonzero mass transport if viscosity is taken into account, which is consistent with the findings of Longuet-Higgins et al. (1953). The interaction between the Stokes drift and planetary vorticity induces the Coriolis-Stokes force under the Eulerian average, and the Coriolis-Stokes force changes the energy input of the mixed layer in the ocean, the velocity profile structure and the surface current field (Hasselmann et al. 1970). The Stokes-vortex force, which is generated by the interaction between the Stokes drift and the wind-driven current, is the main source of Langmuir circulation (Craik et al. 1976). Langmuir circulation enhances the vertical shear instability of the upper circulation and strengthens the upper turbulence effect. Therefore, the mixing effect in the mixed layer is intensified and the mixing depth is deepened (Li et al. 1995; Polton et al. 2007). McWilliams et al. (1999) compared the global Stokes transport with Ekman transport and concluded that Stokes transport can reach a magnitude comparable to that of wind-driven circulation transport in high wind speed areas. The maximum percentage of Stokes transport in relation to Ekman transport can reach 50%, and Bi (2013) arrived at the same conclusion. The Coriolis-Stokes force was introduced into the classical Ekman model, and the analytical solution obtained is closer to the actual situation than the solution without waves, and the Coriolis-Stokes force changes the dynamic process as well as the velocity profile in the mixed layer (Lewis et al. 2004; Polton et al. 2005; Huang et al. 2008). Restrepo et al. (2007) concluded that in the study of wave-current interactions, the Stokes drift affects the advection transport together with the average current. Li et al. (2008) used the parameterization method to introduce the Stokes drift into the Mellor-Yamada 2.5 turbulent closure model and studied the impact of Stokes drift on the upper ocean mixing. The results showed that the Stokes drift can penetrate the whole mixed layer, while wave breaking affects only the surface of the ocean. The Stokes drift can enhance turbulent mixing and increase energy dissipation during mixing. Tamura et al. (2012), Zhang et al. (2015) and Wang et al. (2015) calculated the sea surface Stokes drift velocity using simulated results from ocean wave models. In these calculations, the Stokes drift velocity at high latitudes was larger than that at low latitudes. Wang (2012) introduced the Stokes drift into the Princeton Ocean Model (POM). The results showed that Stokes drift enhances the shear instability of the upper ocean and strengthens the mixing effect. Zhang et al. (2014) introduced the Stokes drift into the temperature advection equation of the HYCOM (Hybrid Coordinate Ocean Model) model. The results showed that the contribution of the Stokes drift and the mean current to the mixed-layer temperature changes is comparable.
Many studies have dealt with the effects of the Stokes drift on the mixed layer in the ocean, but few studies have dealt with its impact on sea surface temperature (SST). The main factors affecting the SST are the air-sea heat flux and ocean dynamic processes, while ocean dynamic processes include mainly advection, convection and mixing. Deng et al. (2012) introduced the Coriolis-Stokes force as a boundary condition into the momentum equation of the POM and investigated whether both SST and mixed-layer thickness were affected by the Coriolis-Stokes force. Based on the NEMO (Nucleus for European Modelling of the Ocean) ocean circulation model, Breivik et al. (2015) discussed the impact of sea surface stress, the turbulent energy flux caused by wave breaking and the Coriolis-Stokes force on the mixed layer. The three factors played an important role in reducing the SST simulation biases.
In our experimental work, based on the existing research, the wave-induced Stokes drift is calculated using the simulated results from the WAVEWATCH III wave model. The Stokes drift velocity is taken into account in the effect of the wave-induced Stokes drift on the SST simulation in the North Pacific Ocean. The numerical simulation results of SST using the Stokes drift spectrum parameterization scheme and three approximate parameterization schemes are compared and analyzed. The paper consists of four sections. The study area and datasets, parameterization schemes for Stokes drift calculations and description of models and settings are described in the Materials and methods. The wave simulation results, distributions of the Stokes drift and its effect as well as SST simulation results are presented in the Results. Simulated SST validation and a comparison of approximate parameterization schemes are presented in the Discussion. Finally, conclusions are given in the Conclusions section.
The study area is located in the Pacific between 0° and 66°N latitude and between 100°E and 80°W longitude (Fig. 1). In the middle of the study area, the sea area is wide and the sea waves are large.
Geographical location of the North Pacific
Topographic data used in the model were ETOPO5 data with a resolution of 5’ × 5’. Daily averaged wind data and heat flux data (shortwave radiation, longwave radiation, sensible heat flux and latent heat flux) with a resolution of 1.875° × 1.905° came from the National Centers for Environment Prediction (NCEP). The initial temperature field and the salinity field, as well as the lateral boundary were derived from HYCOM data with a resolution of 0.08° × 0.08°. The Optimum Interpolation Sea Surface Temperature (OISST) is a converged product with a spatial resolution of 25 km and a temporal resolution of 1 day, came from the Advanced Very High Resolution Radiometer (AVHRR). This product is derived from the AVHRR satellite-measured data and includes a large-scale adjustment of satellite biases with respect to the in situ data from ships and buoys and uses the optimal interpolation (OI) algorithm for the calculation (Sun et al. 2007; Reynolds et al. 2007; Yan et al. 2016). The European Center for Medium-Range Weather Forecast (ECMWF) reanalysis product ERA-40 with a resolution of 1.5° distinguishes between wind waves and swells, covers the globe and performs well in the North Pacific (Zheng et al. 2015). Buoy data were obtained from the National Data Buoy Center (NDBC) and included wind direction, wind velocity, significant wave height (SWH), wave period and sea temperature. TAO buoys are observational arrays located to the north and south of the equator in the Pacific Ocean, and their datasets are updated daily and are freely available to researchers (Mcphaden 1995; Kara 2006).
The spectral profile parameterization scheme is considered to be a relatively accurate Stokes drift parameterization scheme (Rascle et al. 2008; Webb & Foxkemper 2011; Rascle & Ardhuin 2013; Webb & Fox-Kemper 2015) and was therefore mainly used in this work, whereas the approximate parameterization schemes were introduced only to compare them with the former. The Stokes drift profile based on a one-dimensional wave spectrum is as follows (Breivik et al. 2016):
where
where E(
The sea surface Stokes drift velocity is as follows:
The Stokes transport is as follows:
For deep water,
In Equations (1), (3) and (5), the spectral formula (EL spectrum) proposed by Elfouhaily et al. (1997) is used:
where
where α
where
where
The term
where α
where
The relationship between the frequency spectrum
where
The spectral profile is inconvenient as it is numerically challenging and expensive to integrate the wave spectrum at every vertical level. Therefore, an approximate parameterization scheme is more suitable for the calculation of the Stokes drift. Three approximate parameterization schemes can be alternatives to the spectral profile.
(1) The monochromatic profile is given by Breivik et al. (2014)
where
(2) The exponential profile is given by Breivik et al. (2014):
where
(3) The Phillips spectrum profile is given by Breivik et al. (2016):
where β = 1 and
When
For a vast number of calculations of the spectral profile, it has been customary to use the monochromatic profile (Skyllingstad & Denbo 1995; Polton et al. 2005; Tamura et al. 2012). The monochromatic profile requires the transport, and the surface Stokes drift velocity is easy to implement and computationally inexpensive, but it leads to an underestimation of the shear near the surface. Breivik et al. (2014) proposed an exponential profile based on the same parameters required to compute the monochromatic profile. This approach gives a more accurate shear with no added numerical cost. It offers a slight improvement over the monochromatic profile, however, it is not a good match for the shear. The Phillips spectrum profile was proposed based on the Phillips spectrum (Breivik et al. 2016). We can compute the profile from the same two parameters as the monochromatic profile, which gives a much better match, especially to represent the shear near the surface. Therefore, the Phillips spectrum profile is the best option.
The wave model WAVEWATCH III, developed mainly by Tolman (1996; 2009), is employed to provide the wave field data in this study. The model domain is set from 70°30'S to 66°N latitude and 0°15'E to 359°45'E longitude. The resolution is 0.5° in latitude and longitude. There are 274 × 720 grid points in the horizontal direction. The model be comes stable after a simulation for one month and the simulated time is from January 2014 to December 2014. The simulate d results include many characteristics of the waves, such as the significant wave height, wave direction, wavelength, etc. The simulated significant wave height was validated, and the Stokes drift velocity was calculated based on the simulated results.
The POM is a three-dimensional baroclinic ocean circulation model based on the original equations developed at Princeton University. The model adopts the σ coordinate in the vertical direction and the
SBPOM does not take sea waves into account. In this study, the Stokes drift is calculated based on the simulated results from WAVEWATCH III and introduced into the momentum equation of SBPOM in the form of velocity to analyze the effect of the Stokes drift.
For a viscous incompressible fluid, the momentum equation is as follows:
where
After the Euler method is used for time averaging,
where
where
where
Thus, the momentum equation with the Stokes drift has the following form:
Figure 2 shows the annual averaged significant wave height (SWH) distribution in the North Pacific Ocean simulated by WAVEWATCH III. The annual averaged SWH gradually decreases with decreasing latitude; SWH at high latitudes is higher than that at low latitudes (Zheng et al. 2015). The SWH on continental coasts is mostly less than 1.5 m. The SWH in high-value areas is greater than 2.5 m and the maximum value is 4.3 m, which is located in the middle of the North Pacific Ocean, corresponding to the prevailing westerly winds in the range of 30°N to 60°N. In the open ocean area corresponding to the northeast trade winds from 10°N to 30°N, the SWH is in the range of 2 to 2.5 m and the SWH is approximately 1.5 m near the Hawaiian Islands due to the island distribution. From the equator to 10°N, winds are weak in most sea areas and the SWH is small, approximately 1.5 m.
Annual averaged significant wave height
The simulation results for January, April, July and October 2014 were selected to represent winter, spring, summer and autumn, respectively. As shown in Figure 3, the reanalyzed data were used to validate the simulated results for these 4 months. The differences between the WAVEWATCH III-simulated results and the ECMWF SWH data demonstrate that the simulated errors are within the range of 0.5 m for most of the
Differences between WAVEWATCH III-simulated results and the ECMWF product
sea area and thus, the simulation exhibits good performance. The simulated results are the best in July and the simulated error of most sea areas is within the range of 0.3 m. In the near-equatorial areas, the simulated error is within the range of 0.2 m. Compared with the validation data, the simulated results are low in most areas at low latitudes, but the western Pacific simulated result is high in January, with an error of 0.5 m. At high latitudes, the western Pacific simulated result is low, except for July, and the simulated error in April is the largest, reaching −0.6 m. The simulated result in the eastern part of the ocean is high and the simulation error in January is the largest, reaching 0.5 m.
The simulated SWH data are validated by data from 18 National Data Buoy Center (NDBC) buoy stations. The locations of the buoys are shown in Figure 4.
Locations of the North Pacific NDBC buoys and TAO buoys
The accuracy of the simulated results is measured by calculating the mean error (ME), the root-mean-square difference (RMSD), the correlation coefficient (R), and the nondimensional skill score (SS) between the simulated and observational values. A negative SS indicates that the simulation is poor, but positive values are common. The larger the value, the better the simulation (Murphy 1995). The statistical measures can be written as follows:
where
The number of buoy values (
buoy | buoy | buoy | |||
---|---|---|---|---|---|
51000 | 8745 | 46006 | 8734 | TAO1 | 8756 |
51003 | 8747 | 46070 | 8740 | TAO2 | 8755 |
51004 | 8750 | 46071 | 8741 | TAO3 | 8753 |
52200 | 8738 | 51101 | 8736 | TAO4 | 8751 |
46001 | 8735 | 21413 | 8739 | TAO5 | 8754 |
46002 | 8742 | 21419 | 8740 | TAO6 | 8755 |
46072 | 8730 | 52404 | 8733 | TAO7 | 8753 |
46075 | 8743 | 52401 | 8739 | ||
46005 | 8737 | 52402 | 8736 |
Table 2 presents the ME, RMSD, R and SS of the simulated values and observational values in 2014, and the results indicate that the simulated results of WAVEWATCH III are consistent with the observational values of the NDBC buoys. The simulated values are low: the ME of 15 buoys is negative, the maximum error is −13.64 cm and the minimum error is 1.03 cm. The RMSD is in the range of 7.52 to 32.19 cm. Except for the value of 46072, the RMSD values are all below 30 cm, with a minimum of 7.52 cm. The lowest R is
Comparison of simulated values and buoy values
buoy | ME (cm) | RMSD (cm) | R | SS |
---|---|---|---|---|
51000 | −13.64 | 21.69 | 0.84 | 0.83 |
51003 | −11.11 | 17.32 | 0.89 | 0.84 |
51004 | −9.71 | 14.39 | 0.92 | 0.86 |
52200 | 1.03 | 7.52 | 0.95 | 0.91 |
46001 | −3.53 | 15.38 | 0.94 | 0.88 |
46002 | −3.59 | 20.22 | 0.88 | 0.81 |
46072 | 12.18 | 32.19 | 0.90 | 0.71 |
46075 | −7.01 | 12.79 | 0.93 | 0.89 |
46005 | −1.39 | 18.37 | 0.91 | 0.80 |
46006 | −6.87 | 23.52 | 0.89 | 0.68 |
46070 | −9.63 | 18.85 | 0.91 | 0.82 |
46071 | 3.75 | 19.33 | 0.88 | 0.77 |
51101 | −8.29 | 21.03 | 0.93 | 0.87 |
21413 | −5.63 | 15.98 | 0.88 | 0.84 |
21419 | −10.65 | 20.81 | 0.93 | 0.91 |
52404 | −7.23 | 23.06 | 0.85 | 0.75 |
52401 | −2.00 | 9.67 | 0.95 | 0.90 |
52402 | −3.68 | 23.52 | 0.92 | 0.86 |
0.8401. Thus, all of the R values are greater than 0.8, indicating that the simulated SWH and the measured SWH are strongly correlated and that the simulation is good. The SS values are positive and greater than 0.7. Except for the values of buoys 46071, 46072, and 52404, the remaining SS values are greater than 0.8, indicating that the simulation is good. The reasons for simulation errors may be that swells are not taken into account or the parameterization scheme of wind stress does not fully reflect the transmission process of wind energy (Shi et al. 2013). Swail & Cox (2000) assessed the 10 m wind speeds from NCEP by using advanced third-generation (3-G) ocean wave prediction models and concluded that this approach performed best. Chelton & Freilich (2005) analyzed the accuracy of the 10 m wind speeds from NCEP and found that there was no evidence of bias in the 10 m wind speeds from NCEP. Therefore, the wave results simulated by WAVEWATCH III driven by wind data from NCEP are acceptable.
The sea surface Stokes drift was calculated with Equation (3) (Fig. 5). The distribution of the monthly averaged sea surface Stokes drift in the North Pacific shows a small spatial difference in July and a large spatial difference in other months. In January, April and October, the surface Stokes drift velocity increases with increasing latitude and the drift velocity of the 30°N to 60°N range is greater than that of the 0° to 30°N range. The highest sea surface Stokes drift velocity is at high latitudes, which is consistent with the conclusions of Tamura (2012), Zhang et al. (2014; 2015) and Wang et al. (2015). Areas with high surface Stokes drift velocity values are distributed north of 30°N and correspond to the position of the prevailing westerly winds. The maximum value is 0.2 m s−1. The winds and waves near the equator (0°–5°N) are weak all year round and the Stokes drift is relatively weak. Therefore, the Stokes drift velocity is mostly less than 0.06 m s−1 and varies slightly depending on the season. At high latitudes, the Stokes drift velocity varies greatly from season to season. In July, the Stokes drift velocity was the lowest and uniform, and the maximum velocity domain was located in the South China Sea, with a maximum value of 0.13 m s−1. The maximum Stokes drift velocities in January, April and October were observed in the high-latitude area, and the area with velocities greater than 0.1 m s−1 was large, especially in January and October, occupying almost all areas north of 40°N. In January and October, local high Stokes drift velocity areas were observed in the Philippine Sea and the South China Sea. In April, there were two areas with high drift velocity, and velocities at low latitudes were much lower than those at high latitudes. The Stokes drift velocity was small in July, mainly due to weak wind speed, small energy input and relatively weak
Distribution of the surface Stokes drift velocity
waves. The Stokes drift velocity is highly seasonal, with higher values in January and lower values in July associated with extratropical cyclones (Bell & Kirtman 2018).
To perform a more intuitive comparison of the sea surface Stokes drift velocity with the surface current velocity, the velocity ratio is defined as follows:
where
Figure 6 shows the zonally averaged curve of the monthly mean Stokes drift velocity and the HYCOM velocity as well as the curve of the velocity ratio
Comparison of the surface Stokes drift velocities and HYCOM velocities and parameter b (the red curve shows the change in the monthly averaged sea surface Stokes drift velocity with latitude; the black curve shows the change in monthly averaged HYCOM surface ocean velocity with latitude; the blue curve represents the change in the velocity ratio b with latitude)
and the surface Stokes drift velocity is approximately equal to the surface current velocity, and the impact is slightly significant. In October,
According to Equation (4), the Stokes drift is vertically integrated to obtain the Stokes transport in the North Pacific Ocean (Fig. 7). The distributions of the Stokes transport and Stokes drift are almost
Distribution of Stokes transport
the same. The high-value area of the Stokes transport coincides with the high-value area of the Stokes drift. The Stokes transport distribution shows a small spatial difference in July and a large spatial difference in other months. In January, April and October, the Stokes transport velocity increases with increasing latitude and the transport velocity of 30°N–55°N is greater than that of 0°N–30°N, except for July. The high Stokes transport velocity is observed at high latitudes, which is consistent with the conclusions of Deng et al. (2007), Zhang et al. (2014; 2015) and Wang et al. (2015). Near the equator (0°–5°N), the Stokes transport velocity is low due to the low Stokes drift and does not vary significantly seasonally. In terms of its size distribution, the transport velocity is higher at high latitudes in January, April and October: the Stokes transport velocity in the area of 30°N–60°N is mostly above 0.35 m2 s−1, reaching its maximum value of more than 0.95 m2 s−1 around 45°N in January. Peaks occur in the area around 45°N. In July, in most areas of the North Pacific, the distribution of Stokes transport velocity was uniform, with values below 0.3 m2 s−1. The seasonal variation in the Stokes transport velocity in the North Pacific is apparent. The Stokes transport velocity in July is low, while the Stokes transport velocity in January, April and October is high, which is consistent with the seasonal variation characteristics of the surface Stokes drift.
Previous research (Deng et al 2007; Zhang et al. 2009; Zhang 2013; Bi 2013) has shown that Stokes transport accounts for a large proportion of ocean transport. McWillianms et al. (1999) and Bi (2013) compared the magnitudes of Stokes transport and Ekman transport and concluded that Stokes transport accounts for up to 50% of the Ekman transport. Polton et al. (2005) defined the Ekman-Stokes number according to the work of Mcwilliams & Restrepo (1999), and this parameter represents the proportion of Stokes transport in the net transport. The Ekman-Stokes number is calculated as follows:
where
where
Figure 8 shows the variation in the latitudinal mean Ekman-Stokes number with latitude obtained according to Equation (39). The Ekman-Stokes number is zero at the equator due to a Coriolis value equal to 0, whereas the Ekman-Stokes number in high-latitude areas is greater than 0.1 and higher than that in low-latitude areas. The Stokes transport in the sea near 50°N in January accounts for more than 20% of the net transport, and the Stokes transport in other months is less than 20%. In July, from 20°N to 55°N, the Stokes transport accounts for approximately 10% of the net transport, and the Stokes transport in other months features a peak near 50°N. North of 60°N, the Ekman-Stokes number in January and October decreases gradually, while the Ekman-Stokes number in April and July increases gradually. The influence of Stokes transport on the mixing layer in the ocean is concentrated in the high-latitude sea area, which occupies a large proportion of the net transport. The divergence of seawater caused by transport changes the structure of the mixed layer in the ocean, leading to the enhancement of the mixing effect in the upper ocean, which accounts for the SST changes. Compared to other items in the standard Ekman balance, the Ekman-Stokes number (Polton & Belcher 2007) also reflects the distribution of the relative magnitudes of the Coriolis-Stokes force in the North Pacific. At high latitudes, the Coriolis-Stokes force is relatively large.
Ekman-Stokes number variation curve (the curve shows the variation in the latitudinal average Ekman-Stokes number with latitude)
Langmuir circulation is generated by the wave-current interaction in the upper ocean and is a longitudinal spiral motion that is axially parallel to the wind direction (Langmuir 1938). The Langmuir circulation strengthens the vertical shear instability of the upper ocean circulation, so the effect of mixing in the upper mixed layer is strengthened and the depth of mixing is increased, and the Langmuir circulation also produces small-scale effects in the upper ocean, e.g. by affecting the turbulent kinetic energy. McWilliams et al. (2004) argued that the increase in turbulent energy caused by the Langmuir circulation is the main mechanism leading to the enhancement of mixing. The Langmuir circulation can enhance the upper turbulence effect and the induced vertical convection can increase the vertical water exchange and play a “rolling” role, which causes the lower layer of cold water to cool the upper layer of water. McWilliams et al. (1997) used the Langmuir number as a measure of the magnitude of the Langmuir turbulence effect and defined it as follows:
where
The monthly averaged Langmuir number in the North Pacific is calculated by using Equation (41) (Fig. 10). The Langmuir number is small at high latitudes
Langmuir number
and large at low latitudes. Therefore, the Langmuir turbulence effect is more evident at high latitudes than at low latitudes. In the range of 0° to 30°N, the Langmuir number in the eastern North Pacific is large and that in the western North Pacific is small, which is consistent with the conclusion of Wang et al. (2015). In the range of 30°N to 60°N, the Langmuir number is less than 0.3. In July, the Langmuir number in the eastern North Pacific was greater than that in the western North Pacific, and in January, April and October, the Langmuir number was less than 0.2 around 30°N. In January and October, the distribution of the Langmuir number is very similar, the differences between areas with Langmuir numbers less than 0.3 are small, and areas with high Langmuir numbers are distributed along the southern coast of North America.
The Stokes drift decreases exponentially with water depth and the maximum influence depth in the vertical direction is referred to as the Stokes drift influence depth. The profile is as follows (Zhang et al. 2014):
where
Stokes drift influence depth
The OISST product was used to validate the simulated SST without the Stokes drift (Fig. 12). When the Stokes drift is not taken into account, the simulated
Differences between SST values simulated without the Stokes drift and OISST data
SST values in most sea areas are high and the error is in the range of ±2°C. In the Japan Sea, the simulated SST is low and the seasonal variation is small. At high latitudes, the simulated SST is high, and the simulation errors in April and October are larger than those in January and July. In the range of 30°N to 45°N, the simulated SST in the middle Pacific is low in January and that in the western Pacific is also low in July, and the simulated SST of the surrounding seas is high. At low latitudes, the simulation error is smaller than that at high latitudes, and the simulation error of most sea areas is within the range of ±0.5°C. In the sea near the equator, the simulated SST in January is high, with a maximum error of 1.5°C, while the simulated SST in October is low, with a maximum error of −1°C.
To study the effect of the Stokes drift on the SST simulation, the Stokes drift calculated using the WAVEWATCH III-simulated results was introduced into the SBPOM model. To analyze the influence of the Stokes drift on the SST simulation in the North Pacific and the difference between the simulated SST with the Stokes drift and that without the Stokes drift, the optimal value α is defined:
Figure 13 shows that the Stokes drift mainly plays the role of cooling in the simulation of SST, especially at high latitudes, and the maximum absolute value of α does not exceed 1.2°C. In the sea area north of 30°N, the Stokes drift has less influence in July than in other months, and the absolute value of α is less than 0.3°C. In January, April and October, the Stokes drift has a greater impact, and the absolute value of α is up to 1.2°C. In January, the area affected by the Stokes drift is mainly in the eastern and western parts of the North Pacific. In April, the area affected by the Stokes drift is mainly in the central and eastern parts of the North Pacific. In October, the area affected by the Stokes drift occupies almost the entire area north of 30°N. South of 30°N, the absolute value of α for the four months is less than 0.6°C and the Stokes drift effect is relatively small. In July, the effect of the Stokes drift on the SST simulation is almost negligible. In January, the simulated SST with the Stokes drift in the eastern North Pacific near the equator is lower by approximately 0.3°C. In April, the Stokes drift causes partial cooling in the central North Pacific, with a temperature reduction of 0.5°C. In October, the simulated SST with the Stokes drift in the eastern North Pacific near the equator rises, with a maximum temperature variation of less than 0.5°C.
Optimal value α
The error between the simulated SST with the Stokes drift and OISST data is calculated (Fig. 14). According to Figures 12, 13 and 14, the simulation effect of SST is improved after taking the Stokes drift into account. The simulation error of most sea areas is in the range of ±0.5°C and the simulated
Differences between SST values simulated with the Stokes drift and OISST data
SST is closer to the OISST data. The simulated results in January, April and October are good. In July, the simulated results are good south of 30°N, whereas north of 30°N, the influence of Stokes drift is small and the simulation error is in the range of ±0.8°C. The Stokes transport is prone to producing a divergence, which causes upwelling (Smith 1998). The upwelling enhances the convection in the upper mixed layer, which affects the surface water, and the enhancement of convection between the upper and lower layers of seawater cools the SST. According to Figures 5, 6 and 7, the Stokes drift and Stokes transport in January, April and October are stronger at high latitudes than at low latitudes of the North Pacific, causing a significant drop in temperature. At low latitudes, the Stokes drift and Stokes transport are weak, causing a weaker temperature drop than at high latitudes. The distribution of this cooling is consistent with the distributions of the Stokes drift and Stokes transport. In July, the Stokes drift and Stokes transport are weak throughout the North Pacific, resulting in a weaker cooling of the SST.
The Coriolis-Stokes force, one of the wave-induced stress terms, changes the energy input, the velocity profile and the sea surface current distribution of the mixed layer (Hasselmann et al. 1970), and changes in the current distribution directly affect the energy conservation equation, further affecting the SST. The Stokes drift changes the temperature and density distribution of the upper ocean and the current velocity shear through the Coriolis-Stokes force. The Enhancement of instability causes the sea surface current to change. Changes in the sea surface current lead to changes in transport. The upwelling generated behind the transport lowers the SST. Figures 8 and 9 show that at high latitudes, the Coriolis-Stokes force is strong and the Stokes drift has a strong cooling effect, while at low latitudes, the Coriolis-Stokes force is weak and the Stokes drift has a weak cooling effect. Figure 10 shows that the smaller the Langmuir number, the more obvious the SST cooling effect. As the Stokes drift becomes stronger, the Langmuir number decreases, the Langmuir turbulence effect becomes dominant and the turbulent energy increases. This increased turbulence enhances the convection of the seawater, bringing water from depth to the surface and lowering the SST. In the eastern North Pacific near the equator in October, after taking the Stokes drift into account, the simulated SST rises. This simulated SST rise may be explained by the effect of wave mixing. Changes in the ocean temperature and salinity resulting from the wave mixing cause the pressure gradient force to change and weaken the upwelling. The lower seawater cannot rise to the sea surface. As a result, the simulated SST increases (Song 2006; Qiao 2010; Song 2012; Liu 2017).
Eighteen NDBC buoys and seven TAO buoys presented in Figure 4 are used as validation points. The RMSD, R and SS between the daily averaged buoy time series of SST and the daily averaged simulated SST with and without the Stokes drift in 2014 were calculated (Table 3). Table 4 shows the RMSD, R and SS between
Comparison of simulated SST and buoy SST
buoy | No Stokes drift | Stokes drift | ||||
---|---|---|---|---|---|---|
RMSD (°C) | R | SS | RMSD (°C) | R | SS | |
51000 | 0.90 | 0.77 | 0.73 | 0.82 | 0.78 | 0.74 |
51003 | 0.89 | 0.80 | 0.77 | 0.83 | 0.83 | 0.81 |
51004 | 1.25 | 0.84 | 0.79 | 1.12 | 0.85 | 0.80 |
52200 | 1.46 | 0.84 | 0.82 | 1.25 | 0.85 | 0.83 |
46001 | 1.23 | 0.87 | 0.86 | 1.13 | 0.88 | 0.88 |
46002 | 1.10 | 0.82 | 0.78 | 0.91 | 0.85 | 0.82 |
46072 | 1.24 | 0.94 | 0.88 | 1.01 | 0.94 | 0.90 |
46075 | 1.41 | 0.88 | 0.85 | 1.13 | 0.88 | 0.87 |
46005 | 2.06 | 0.97 | 0.95 | 1.87 | 0.97 | 0.96 |
46006 | 1.43 | 0.86 | 0.81 | 1.28 | 0.86 | 0.85 |
46070 | 0.77 | 0.87 | 0.82 | 0.60 | 0.89 | 0.83 |
46071 | 0.53 | 0.93 | 0.89 | 0.52 | 0.94 | 0.93 |
51101 | 0.81 | 0.91 | 0.86 | 0.80 | 0.93 | 0.89 |
21413 | 1.22 | 0.88 | 0.83 | 0.99 | 0.89 | 0.84 |
21419 | 1.16 | 0.78 | 0.77 | 1.09 | 0.80 | 0.78 |
52404 | 0.99 | 0.87 | 0.85 | 0.83 | 0.88 | 0.88 |
52401 | 1.07 | 0.79 | 0.77 | 0.90 | 0.80 | 0.77 |
52402 | 0.97 | 0.79 | 0.74 | 0.96 | 0.82 | 0.76 |
TAO1 | 0.65 | 0.91 | 0.87 | 0.52 | 0.92 | 0.89 |
TAO2 | 0.74 | 0.87 | 0.84 | 0.58 | 0.88 | 0.84 |
TAO3 | 0.82 | 0.83 | 0.78 | 0.78 | 0.84 | 0.80 |
TAO4 | 1.12 | 0.85 | 0.80 | 1.02 | 0.85 | 0.81 |
TAO5 | 1.14 | 0.90 | 0.81 | 0.99 | 0.94 | 0.83 |
TAO6 | 0.92 | 0.81 | 0.92 | 0.85 | 0.84 | 0.95 |
TAO7 | 1.21 | 0.87 | 0.82 | 1.03 | 0.88 | 0.82 |
Average | 1.08 | 0.86 | 0.82 | 0.95 | 0.87 | 0.84 |
Comparison of simulated subsurface temperature and buoy data
buoy | No Stokes drift | Stokes drift | ||||
---|---|---|---|---|---|---|
RMSD (°C) | R | SS | RMSD (°C) | R | SS | |
51000 | 0.65 | 0.83 | 0.72 | 0.88 | 0.84 | 0.81 |
51003 | 1.38 | 0.86 | 0.80 | 0.93 | 0.83 | 0.77 |
51004 | 0.94 | 0.96 | 0.83 | 0.92 | 0.92 | 0.92 |
52200 | 0.53 | 0.78 | 0.80 | 0.51 | 0.84 | 0.85 |
46001 | 1.77 | 0.84 | 0.76 | 1.93 | 0.80 | 0.79 |
46002 | 0.86 | 0.91 | 0.91 | 0.84 | 0.85 | 0.87 |
46072 | 1.67 | 0.83 | 0.87 | 1.55 | 0.83 | 0.85 |
46075 | 2.11 | 0.85 | 0.80 | 1.43 | 0.79 | 0.86 |
46005 | 0.91 | 0.84 | 0.74 | 1.27 | 0.92 | 0.83 |
46006 | 0.73 | 0.92 | 0.83 | 0.93 | 0.88 | 0.85 |
46070 | 0.94 | 0.93 | 0.94 | 1.69 | 0.85 | 0.78 |
46071 | 0.82 | 0.86 | 0.84 | 0.95 | 0.90 | 0.86 |
51101 | 0.98 | 0.85 | 0.87 | 0.84 | 0.84 | 0.91 |
21413 | 0.77 | 0.78 | 0.91 | 0.87 | 0.87 | 0.87 |
21419 | 0.74 | 0.91 | 0.88 | 0.82 | 0.88 | 0.90 |
52404 | 1.26 | 0.88 | 0.83 | 1.25 | 0.89 | 0.78 |
52401 | 1.64 | 0.78 | 0.81 | 0.97 | 0.79 | 0.81 |
52402 | 0.77 | 0.87 | 0.93 | 0.95 | 0.77 | 0.89 |
TAO1 | 1.74 | 0.82 | 0.82 | 1.34 | 0.85 | 0.86 |
TAO2 | 0.52 | 0.95 | 0.85 | 0.75 | 0.86 | 0.89 |
TAO3 | 1.09 | 0.82 | 0.87 | 0.89 | 0.91 | 0.91 |
TAO4 | 0.66 | 0.86 | 0.93 | 0.69 | 0.85 | 0.88 |
TAO5 | 0.92 | 0.90 | 0.87 | 0.90 | 0.86 | 0.85 |
TAO6 | 1.37 | 0.83 | 0.78 | 1.13 | 0.84 | 0.81 |
TAO7 | 0.81 | 0.93 | 0.84 | 1.03 | 0.83 | 0.87 |
Average | 1.06 | 0.86 | 0.84 | 1.05 | 0.85 | 0.85 |
the daily averaged buoy time series of sea subsurface temperature and the daily averaged simulated sea subsurface temperature.
As shown in Table 3, there is a large error between the simulated SST without the Stokes drift and the buoy data. The maximum RMSD is 2.06°C and the R values for four buoys (51000, 21419, 52401, and 52402) are less than 0.8. The R value for 46005 is the largest, 0.97, and the values for all buoys are greater than 0.75. After taking the Stokes drift into account, the RMSD is reduced; the maximum RMSD is 1.87°C, the minimum RMSD is only 0.52°C, and the average RMSD is reduced from 1.08°C to 0.95°C. Except for buoy 51000, the buoys have R values greater than 0.8, and the maximum R value is 0.97. Compared with the simulated results without the Stokes drift, the average R of the result with the Stokes drift increased by 0.01, and the average SS increased from 0.82 to 0.84. Therefore, the introduction of the Stokes drift into the ocean circulation model can improve the accuracy of the SST simulation.
In Table 4, the simulated sea subsurface temperature with the Stokes drift is slightly different from that without the Stokes drift. The RMSD of 11 buoys with the Stokes drift is greater than that without the Stokes drift. The maximum RMSD is 2.11°C. The average RMSD with the Stokes drift is smaller than that without the Stokes drift by only 0.1. The R values of these results are all above 0.75, and the average R of the simulated result with the Stokes drift decreases by 0.1. The average SS of the simulated result with the Stokes drift increases by 0.1. For the simulation of subsurface temperature, similar performance can be obtained with the Stokes drift and without the Stokes drift.
In this study, four buoys located in different areas were selected to analyze the simulated results without the Stokes drift, the simulated results with the Stokes drift and the buoy data (Figure 15). The temperature profile simulated with the Stokes drift is similar to the observational temperature profile of the buoys. Based on a comparison of the vertical simulation results with and without the Stokes drift, even though the Stokes drift depth is only 20 m, it indirectly affects the entire mixed layer through the interaction of the Langmuir circulation and the downwelling (Skyllingstad & Denbo 1995; Li et al. 2005). The reason is that the Coriolis-Stokes force and the Stokes-vortex force result in a reduced mean shear and the Langmuir circulation intensifies the downwelling (Polton & Belcher 2007; Tamura et al. 2012; Li et al. 2013). The temperature profiles in Figure 15 show that after taking the Stokes drift into account, the upper boundary of the thermocline can be determined more accurately, and the analysis of the thermocline characteristics can be improved.
Comparison of simulated and buoy-measured temperature profiles
The Stokes drift velocity is calculated by three Stokes drift approximation parameterization schemes, introduced into the SBPOM model in the same way. The four experiments: E-1, E-2, E-3 and E-4 correspond to the spectral profile, the monochromatic profile, the exponential profile and the Phillips spectrum profile, respectively. Using the daily averaged buoy time series from 25 buoys presented in Figure 4 to calculate the RMSD of the simulated results of E-2, E-3, and E-4, the effects of different Stokes drift approximate parameterization schemes are compared and the results are shown in Table 5.
RMSD of experiments E-1, E-2, E-3, and E-4 (°C)
Buoy | E-1 | E-2 | E-3 | E-4 |
---|---|---|---|---|
51000 | 0.82 | 0.86 | 0.85 | 0.85 |
51003 | 0.83 | 0.85 | 0.87 | 0.83 |
51004 | 1.12 | 1.22 | 1.24 | 1.19 |
52200 | 1.25 | 1.35 | 1.33 | 1.32 |
46001 | 1.13 | 1.21 | 1.20 | 1.15 |
46002 | 0.91 | 1.04 | 1.02 | 0.93 |
46072 | 1.01 | 1.19 | 1.17 | 1.12 |
46075 | 1.13 | 1.18 | 1.15 | 1.20 |
46005 | 1.87 | 1.82 | 1.80 | 1.75 |
46006 | 1.28 | 1.39 | 1.35 | 1.34 |
46070 | 0.60 | 0.78 | 0.77 | 0.65 |
46071 | 0.52 | 0.61 | 0.64 | 0.58 |
51101 | 0.80 | 0.90 | 0.89 | 0.83 |
21413 | 0.99 | 1.13 | 1.11 | 1.08 |
21419 | 1.09 | 1.17 | 1.16 | 1.13 |
52404 | 0.83 | 1.03 | 1.01 | 0.98 |
52401 | 0.90 | 1.13 | 1.12 | 1.10 |
52402 | 0.96 | 1.11 | 1.09 | 1.04 |
TAO1 | 0.52 | 0.64 | 0.62 | 0.55 |
TAO2 | 0.58 | 0.69 | 0.67 | 0.70 |
TAO3 | 0.78 | 0.90 | 0.88 | 0.87 |
TAO4 | 1.02 | 1.20 | 1.19 | 1.16 |
TAO5 | 0.99 | 1.10 | 1.09 | 1.04 |
TAO6 | 0.85 | 0.96 | 0.95 | 0.92 |
TAO7 | 1.03 | 1.15 | 1.14 | 1.14 |
Average | 0.95 | 1.06 | 1.05 | 1.02 |
Table 5 shows that there is a difference between the Stokes drift parameterization schemes in the SST simulation. The comparison of the results of the four experiments reveals that the RMSD of E-1 is small, the RMSD of E-2 is greater than that of E-1, and the RMSD for 46005 in E-3 is the only value in E-3 that is smaller than that in E-1. The RMSD for 46005 in E-4 is also smaller than that in E-1, but the RMSD values for the remaining buoys are greater than those in E-1. In terms of the average RMSD, the average RMSD of E-1 is the smallest. The average RMSD of E-4 is the closest to that of E-1, followed by E-3, and finally E-2. The RMSD values for two buoys (46075 and TAO2) in E-4 are greater than those in E-3. The RMSD values for other buoys are lower than those in E-3. The RMSD values of all buoys in E-4 are lower than those in E-2. The RMSD values of E-3 slightly differ from those in E-2, but the RMSD values for 22 buoys in E-3 are less than those in E-2. Of the four experiments, the best simulated results were obtained in E-1, with RMSD of 0.95°C, and the simulated results of E-4 are the closest to those of E-1, with RMSD of 1.02°C. The RMSD of E-3 is 1.05°C, thus following E-4. The RMSD of E-2 is 1.06°C, indicating that the simulation is not as good as that of E-1, E-3, or E-4. In general, the simulated results with the Stokes drift are better than those without the Stokes drift.
The temperature profiles of E-1, E-2, E-3, and E-4 were analyzed by selecting four different buoys from those presented in Figure 15. When analyzing the vertical temperature profiles in the mixed layer, the smallest error between the simulated results and the buoy data is that of E-1, followed by E-4. As illustrated by B, C, and D in Figure 16, compared to E-2 and E-3, E-4 performed best in the vertical temperature simulation in the mixed layer. Below the mixed layer, the simulations of E-2, E-3, E-4 and E-1 are not satisfactory. The improvement of SST simulation with the Stokes drift parameterization schemes is mainly reflected in the mixed layer.
Simulated results of different Stokes drift parameterization schemes
Applying the wave spectrum results of ocean wave models to ocean circulation models is a great challenge. The computational complexity of the spectral parameterization scheme is very high, especially the integration of the wave spectrum on each vertical layer. As evidenced by Table 5 and Figure 16, the Phillips spectrum profile approximate parameterization scheme is the best among the three approximate parameterization schemes. With a lower numerical cost than the spectral profile parameterization scheme, these three approximate parameterization schemes require only the transport and the surface Stokes drift velocity. Therefore, in the case of numerical costs, we can use the Phillips spectrum profile approximate parameterization scheme in the ocean circulation model.
The wave-induced Stokes drift affects the ocean dynamic process and the distribution of ocean elements in the upper ocean. In the numerical simulation of SST, the wave-induced Stokes drift affects the simulation of SST. After taking the Stokes drift into account, the simulated SST decreases and the simulation performance is improved. Compared with the OISST data, the simulation results are still high in most sea areas, but the error is reduced to ±0.5°C. The comparison of the simulated results with and without the Stokes drift reveals that the Stokes drift cools the SST, and similar performance can be obtained in the simulation of sea subsurface temperature simulation. The influence on the SST simulation is mainly reflected in the high-latitude sea areas and is related to the strong Stokes drift in the high-latitude sea areas. The distribution of SST improvement is similar to the distribution of the Stokes drift. The simulation of the temperature profile was also validated. The temperature profile simulated with the Stokes drift is similar to the measured temperature profile of the buoys in the mixed layer, and the Stokes drift can affect the entire mixed layer.
The results of the SST simulation with the spectral profile parameterization scheme are compared with those of the three Stokes drift approximate parameterization schemes. The simulation of SST and the vertical temperature profile in the mixed layer is improved when the Stokes drift is included. Furthermore, the results of the spectral parameterization scheme are the closest to the buoy data. Of the approximate parameterization schemes, the Phillips spectrum profile parameterization scheme is the second best, followed by the exponential profile parameterization scheme and the monochromatic profile parameterization scheme. The results of the SST simulation with the monochromatic profile parameterization scheme are still better than those without the Stokes drift. The wave-induced Stokes drift affects the simulation of SST, reducing the simulation error for SST. Therefore, the effect of the wave-induced Stokes drift should be taken into account in simulations of SST.