1 Introduction

In the open ocean, the subsurface pycnocline develops seasonally in mid- and high-latitude regions and is sustained throughout most of the year in low-latitude regions. The subsurface pycnocline divides well-sunlit/nutrient-poor surface mixed layers and poorly lit/nutrient-rich lower layers (Longhurst 2007). The water exchange between these layers is regulated by the vertical eddy diffusivity that causes upward nutrient transport. The order of magnitude of the vertical eddy diffusivity in the pycnocline (including the subsurface pycnocline) of the oligotrophic ocean has been estimated as 10–5 m2 s–1 (Gregg 1987; Ledwell et al. 1993; Lewis et al. 1986). While this is much weaker than that of the surface mixed layer, the emergence of a subsurface chlorophyll maximum at the nutricline indicates that vertical eddy diffusivity plays a substantial role in generating new production (Cullen 2015; Furuya 1990).

In addition to the general consensus on its order of magnitude, at least in the open ocean under calm conditions, the vertical eddy diffusivity is known to have substantial temporal and spatial variability. This is attributed in part to the intermittent nature of turbulence (e.g., Baker and Gibson 1987), although geographical variability related to differences in various forcing and background conditions is also important. Waterhouse et al. (2014) compiled profiles of turbulent dissipation rates at various sites and obtained a global map of vertical diffusivity. While the globally averaged vertical diffusivity above 1000 m was O (10–5) m2 s–1, regional variability was evident, explained in part by the levels of local power input by internal tides and wind-driven near-inertial waves. Using strain-based estimates of mixing from Argo floats, Whalen et al. (2018) suggested that wind-induced mixing in the pycnocline is modulated by mesoscale currents, which also vary geographically in activity and strength. Spatial variability in diffusivity has also been estimated from an inverse method using salinity profiles of Argo floats in the North Pacific (Kouketsu 2018), with results suggesting enhanced diffusivity around prominent topography.

Since the nutrient flux associated with turbulent mixing is likely a major driver of new production in subsurface layers, it is important to refine the estimates of vertical eddy diffusivity in the subsurface pycnocline, especially across a broad geographical range. As this layer is typically above 250 m and often under the influence of wind-driven geostrophic flows and secondary disturbances, a large number of direct microstructure measurements are required. However, few studies that have presented direct estimates of vertical eddy diffusivity have focused on (1) the subsurface pycnocline just below the surface mixed layer, and (2) its geographical distribution, especially in the Pacific. For the first point, Waterhouse et al. (2014) averaged the diffusivity from the bottom of mixed layer to 1000 m, and Whalen et al. (2018) made evaluations for 250–500 m at the shallowest. The inverse method by Kouketsu (2018) for solving tracer equations on the isopycnal surface did not cover the near surface layers above 25.0 σθ.

For the second point, microstructure profiles in the open Pacific have generally been concentrated around turbulence hot spots (Waterhouse et al. 2014, their fig. 1c) or near the equator (e.g., Inoue et al. 2012), prominent topography such as the Hawaiian Ridge (e.g., Klymak et al. 2006; Rudnick et al. 2003), and strong currents such as the Kuroshio (e.g., D’Asaro et al. 2011; Kaneko et al. 2012; Nagai et al. 2009, 2012).

Fig. 1
figure 1

Stations for microstructure measurements during the three research cruises. Color shading indicates seafloor depth, from Smith and Sandwell (1997)

Full size image

Nevertheless, several studies have conducted turbulence measurements in the Pacific with a regional- or basin-scale coverage (Bouruet-Aubertot et al. 2018; Fernández-Castro et al. 2014; Kaneko et al. 2020; Liu et al. 2017; Mori et al. 2008; Moum and Osborn 1986). Some of these studies proposed various mechanisms enhancing the vertical mixing that could cause geographical variability in upper-layer diffusivity. The observations by Fernández-Castro et al. (2018) covered the western parts of the South Pacific subtropical and tropical regions and the central-to-eastern part of the North Pacific tropical region. The mean vertical diffusivity in the upper pycnocline of these regions was 1.2–3.3 × 10–5 m2 s–1, except for the equatorial area (29–49 × 10–5 m2 s–1), which is attributed to the vertical velocity shear of the Equatorial Undercurrent (EUC). Bouruet-Aubertot et al. (2018) conducted microstructure measurements in the tropical South Pacific along ~ 19°S from 160°E to 160°W. In their section, the mean vertical eddy diffusivity in the upper 100–500 m was 0.28–0.60 × 10–5 m2 s–1, stronger in the west and weaker in the east. The enhancement in the west was explained by an energetic baroclinic near-inertial wave caused by a cyclone.

In this study, we estimate the vertical eddy diffusivity in the subsurface pycnocline along meridional and zonal sections of the Pacific along 170°W and ~ 23°N. Our main objective is to quantify the magnitude of the diffusivity in the subsurface pycnocline where nutrients are supplied to the subsurface chlorophyll maximum. The remainder of the paper is organized as follows: In Sect. 2, observations and data processing methods are explained. In Sect. 3, we first describe large-scale hydrographic structure across the Pacific and then present the distributions of the intensity of turbulence in the subsurface layer. Distributions of the internal tide and wind forcing conditions, and fine-scale shear-to-strain ratio are also examined in this section. Finally in Sect. 4, we made comparisons of our results with previous studies and discuss potential mechanisms enhancing the subsurface turbulence and possible impacts on the nutrient transport.

2 Materials and methods

2.1 Observations and data processing

Observations were conducted during research cruises KH-13-7, KH-14-3, and KH-17-4 of the R/V Hakuho Maru. The 170°W section was observed during KH-13-7 (equator to 40°S with a 5° interval, except for 5°S) and KH-14-3 (equator to 50°N with a ~ 5° interval), and a section at ~ 23°N (21°30′N–23°N, with a 10° interval) was conducted during KH-17-4 (Fig. 1 and Table 1). The observations took place during the seasonal stratification seasons in the North and South Pacific, during (from the first station to the last station along the section) December 2013 to January 2014 (KH-13-7), July 2014 (KH-14-3), and August–October 2017 (KH-17-4). The depth of the seafloor at the observation stations exceeds 3000 m, except for ~ 2800 m at 20°N, 170°W near the Mid-Pacific Seamounts (Fig. 1 and Table 1).

Table 1 Summary of VMP observations
Full size table

At each station, conductivity–temperature–depth (CTD; SBE911plus; SeaBird Electronics Inc.) observations were conducted down to the sea floor (several casts during KH-14-3 were down to 1000 m), and expendable CTD (XCTD: Tsurumi Seiki, Co., Ltd.) surveys were carried out between the stations during KH-13-7 and KH-14-3. Horizontal velocity profiles were recorded continuously throughout the cruise using a shipboard acoustic Doppler current profiler (SADCP; Ocean Surveyor 38 kHz; Teledyne RD Instruments). These data from CTD, XCTD, and SADCP were used to describe general hydrographic structure along the transects. Horizontal velocity profiles were also measured at the stations using a lower acoustic Doppler current profiler (LADCP; WorkHorse Monitor 300 kHz; Teledyne RD Instruments) mounted on the CTD frame. The LADCP data with higher vertical and lower horizontal resolutions than the SADCP data were analyzed together with the CTD data to estimate shear-to-strain ratio R. The calculation methods of R are given later in this section.

A vertical microstructure profiler (VMP-250; Rockland Scientific) was deployed 1–9 times at the station same as the CTD (Table 1). The instrument, which was attached with a thin tether, measured and internally recorded microstructure data while free-falling with a nominal speed of 0.4–0.8 m s–1 down to the maximum tether length. The microscale velocity data obtained from two shear probes were processed using the software provided by the manufacturer, which integrates wavenumber spectra of microscale velocity shear within temporal bins (in the present case, 8-s bins) to obtain the turbulent kinetic energy dissipation rate ε (Wolk et al. 2002). If the spectra were not consistent with the Nasmyth universal spectrum, the corresponding data were manually excluded. The two ε values at each bin from the two shear probes were averaged if they were within a factor of 4, or otherwise the smaller one was taken as in Fer et al. (2014). The ε profiles were interpolated at a 1-m interval, and then, multiple profiles were averaged to obtain the mean profile 〈ε〉 for each station. Mean profiles of buoyancy frequency 〈N2〉 at each station were calculated using the conductivity–temperature sensor attached to the VMP-250. The vertical eddy diffusivity was then obtained as 〈Kρ〉 = Γ〈ε〉/〈N2〉 (Osborn 1980), where Γ = 0.2 was assumed.

2.2 Averaging over the subsurface pycnocline

To highlight the geographical distribution, the depth-averaged turbulence properties \(\overline{\langle \varepsilon \rangle }\) and \(\overline{\langle {K}_{\rho }\rangle }\) in the subsurface pycnocline were also estimated with confidence intervals. To do this, we first calculated the vertical average of the ensemble-means (hereafter called ensemble and vertical average) energy dissipation rate \(\overline{\langle \varepsilon \rangle }\) with 95% confidence intervals over the 150-m layer from 10 m below the surface mixed layer, using the bootstrap method. The mixed layer depth was defined as the depth where the ensemble-mean potential density exceeds that at 10 m by 0.125 kg m–3. We adopted this large potential density difference, which is often used for climatological data sets (Levitus 1983), to exclude data within any remnant mixed layer that had developed in the short term (e.g., a nighttime mixed layer). Although the diurnal cycle of turbulence is prominent in the equatorial area even below the mixed layer (e.g., Gregg et al. 1985), we presume that the above criterion selects the layer with less diurnal variation, which was termed the upper core layer (of the Equatorial Undercurrent) by Inoue et al. (2012). In estimating the confidence interval, a fixed decorrelation scale of 10 m was adopted, which was greater than the median values of 6–8 m based on vertical autocorrelation of log (ε).

A 150-m bin was selected so that the subsurface chlorophyll maximum was included within the bin. As the mean profiles at the equator were limited to depths of 190 m and 119 m during austral and boreal summer, respectively, \(\overline{\langle \varepsilon \rangle }\) was calculated using the available data for the former (a 91-m layer) but not calculated for the latter. The ensemble and vertical average of \(\overline{\langle {N}^{2}\rangle }\) were calculated similarly, and \(\overline{\langle {K}_{\rho }\rangle }\) was obtained as for 〈Kρ〉.

2.3 Forcing factors and shear-to-strain ratio

The geographical distribution of \(\overline{\langle \varepsilon \rangle }\) and \(\overline{\langle {K}_{\rho }\rangle }\) is compared with major forcing factors that could enhance the turbulence. We consider the influences of internal tide and wind power on the turbulence in the subsurface pycnocline. This is done so following Waterhouse et al. (2014) that focused on the source terms for the internal wave field, instead of wave-breaking processes.

Energy dissipation caused by the internal tide can occur through multiple pathways, with source waves from either remote or local fields. We use the 0.5°-resolution dataset of temporally and vertically averaged total dissipation rate caused by the internal tide (de Lavergne et al. 2019), which is denoted in the present study as Dtide. From the global map, Dtide values were spatially interpolated onto the positions of our stations, to indicate the potential influence of the internal tide. We selected Dtide instead of simpler parameters such as the barotropic-to-baroclinic conversion rate, in order to consider the effects of far-field internal tide dissipation to diffusivity (de Lavergne et al. 2019). Although Dtide is temporally averaged, our estimations of \(\overline{\langle \varepsilon \rangle }\) and \(\overline{\langle {K}_{\rho }\rangle }\) at each station are derived mostly from snapshot observations. Therefore, some of variability of near-field mixing, included in the observation-based products of \(\overline{\langle \varepsilon \rangle }\) and \(\overline{\langle {K}_{\rho }\rangle }\), cannot be directly compared with Dtide. However, as presented in previous studies such as Waterhouse et al. (2014), dissipation rates from individual profiles were generally related to the averaged power input from tide.

The impact of wind is considered in terms of the local turbulent energy production by the wind ρ0u*3, where u* and ρ0 are the friction velocity and the upper-layer water density, respectively. Hourly wind stress data from the second Modern-Era Retrospective analysis for Research and Applications (MERRA2; Geralo et al. 2017) were extracted at the locations of stations, and the temporal average calculated for 10 days until the mean time of VMP-250 deployments, and then, ρ0u*3 was calculated assuming ρ0 = 1025 kg m–3; herein, this is denoted as Ewind. The duration of 10 days was employed considering the duration of 30–50 days used in Whalen et al. (2018) focusing on deeper layers of 250–500 m, while the geographical pattern of Ewind did not show substantial changes with different durations in ranges of 1–30 days.

To examine the regime of the fine-scale variability influencing the turbulence in the subsurface pycnocline, the shear-to-strain ratio was estimated from the CTD and the LADCP data. For internal waves with a single frequency, this parameter was often denoted as \({R}_{\omega }\) and measures the frequency content as \({R}_{\omega }\simeq ({\omega }^{2}+{f}^{2})/({\omega }^{2}-{f}^{2})\) (e.g., Chinn et al. 2016; Kunze et al. 2006), where ω and f are the internal wave and inertial frequencies, respectively. Although the variations in the subsurface pycnocline just below the mixed layer do not likely follow this simple formulation, the shear-to-strain ratios themselves can be calculated from our data, which in the present study is denoted as R.

In calculating R, the raw CTD and LADCP data were first processed using the Seasave V7 (https://www.seabird.com/software) and LDEO software IX (https://www.ldeo.columbia.edu/~ant/LADCP.html; Thurnherr 2010), respectively. The detrended strain profiles for the 500-m bin below the surface mixed layer (mixed layer + 10 m – mixed layer + 510 m) were obtained from the buoyancy frequency N2 calculated from the 1-m interval CTD data as \(\xi_{z} = {{\left[ {N^{2} } \right]^{^{\prime}} } \mathord{\left/ {\vphantom {{\left[ {N^{2} } \right]^{^{\prime}} } {\left[ {\overline{N}^{2} } \right]}}} \right. \kern-\nulldelimiterspace} {\left[ {\overline{N}^{2} } \right]}}\), where overbar denotes 500-m-bin average and []ʹ indicates vertical detrending that is to eliminate the strong gradient of N2 in the subsurface. The detrended shear profiles \(\left({U}_{z},{V}_{z}\right)\) for the 500-m bin were similarly obtained from the 10-m interval LADCP data. The strain and shear spectra (normalized by the buoyancy frequency) were obtained as \({\phi }_{st}=S\left({\xi }_{z}\right)\) and \({\phi }_{sh}=S\left({U}_{z}\right)/\overline{{N}^{2}}+S\left({V}_{z}\right)/\overline{{N}^{2}}\), where \(S\left(f\right)\) denotes a wavenumber spectrum of a profile f. These spectra were then integrated over the range from \({k}_{min}=1/200\) to \({k}_{max}=1/50\) cpm as

$$\left( {\Phi _{{st}} ,\Phi _{{sh}} } \right) = \int_{{k_{{\min }} }}^{{k_{{\max }} }} {\left( {\phi _{{st}} \left( k \right),\phi _{{sh}} } \right)} dk$$
(1)

and then R was estimated as

$$R = \frac{{\Phi_{sh} }}{{\Phi_{st} }}$$
(2)

The interval of integration [1/200, 1/50] was carefully selected to represent the scale of the bin-averaging (150 m) and not to include the range where shear spectra are attenuating. The confidence intervals of the spectral estimates of R were calculated using the bootstrap method.

3 Results

3.1 General hydrography

The meridional section along 170°W penetrated the subarctic gyre, the two subtropical gyres and the equatorial current systems, and the zonal section ran through the middle of the subtropical gyre in the North Pacific (Fig. 2). The mesoscale variability was relatively low along 170°W and the eastern North Pacific. From central to western North Pacific along the zonal section, the amplitude of the mesoscale variability showed moderate level of ~ 0.2 m, which is attributed to mesoscale eddies associated with the Subtropical Counter Currents.

Fig. 2
figure 2

Station plots on maps of absolute dynamic topography, based on the satellite altimetry data produced by the Copernicus Marine Environment Monitoring Service (http://marine.copernicus.eu/), which is averaged over the observation periods of each cruise: (a) KH-13-7 during austral summer in 2013/2014, (b) KH-14-3 during boreal summer in 2014 and (c) KH-17-4 during boreal summer in 2017

Full size image

The subsurface pycnocline appears either as the seasonal pycnocline in the mid- and high-latitude regions or the permanent pycnocline in the low-latitude regions (Fig. 3). The main pycnocline is shallow in the equatorial area, but generally deeper than the seasonal pycnocline in the higher-latitude regions. The salinity profiles clearly depicted the water mass distributions in these pycnoclines. In the mid- and low-latitude South Pacific, subsurface pycnoclines were mostly occupied by the high salinity water of > 35.0, which is defined as the South Pacific Subtropical Underwater (STUW) in the upper part and the South Pacific Subtropical Mode Water (STMW) in the lower part (Talley et al. 2011; Fig. 3a). In the North Pacific, the occurrence of the North Pacific STUW was confined in the upper layer/lower latitude, and the North Pacific STMW and the North Pacific Central Mode Water (CMW) characterized by the thick layer were more clearly observed in 25–40°N than the South Pacific (Oka and Qiu 2011). On the southern and northern sides of the saline water outcropped to the surface around the equator, the low-salinity pool associated with the Intertropical Convergence Zone (ITCZ) occurred above the pycnocline in 5–12°. In the North Pacific, the prominently fresh subarctic water that was connected to the North Pacific Intermediate Water (NPIW) was observed (Fig. 3b). For the zonal section, the subsurface pycnocline mainly consisted of the North Pacific STUW west of 140°W, while the North Pacific Eastern STMW (Hanawa and Talley 2001) was recognized in 140–130°W (Fig. 3c). Although the mode waters are characterized by pycnostads rather than a pycnocline, we include them in the calculation of the \(\overline{\langle \varepsilon \rangle }\), \(\overline{\langle {N}^{2}\rangle }\), and \(\overline{\langle {K}_{\rho }\rangle }\) terms that are defined as the mean dissipation rate, stratification, and diffusivity in the subsurface pycnocline, respectively, considering the connectivity of the diffusive transport at the scale of O(100) m.

Fig. 3
figure 3

Vertical cross sections of the salinity (color shading) obtained using data from conductivity–temperature–depth (CTD) measurements at the stations and from expendable CTD (XCTD) deployed between the stations: (a and b) meridional sections in the South and North Pacific along 170°W, and (c) zonal section in the subtropical North Pacific along 21°30′N–23°N. The gray contour lines, blue solid lines, and blue dashed lines indicate potential density [kg m–3], mixed layer depth, and the lower bound of the turbulence averaging, respectively

Full size image

Horizontal velocity in the upper layer was moderate along both of the meridional and zonal sections, except in the equatorial area (Fig. 4). The prominent eastward jet of the Equatorial Undercurrent was observed, where maximum eastward speed exceeded 0.7 m s–1 at the core. The surface-intensified eastward flows at slightly higher latitudes are the North and South Equatorial Counter Currents, and the westward currents in 8–16°N and 13–22°S are North and South Equatorial Currents, respectively. The eastward speed of the North and South Pacific Currents originating from the western boundary currents was mostly less than 5 cm s–1 at 170°W (Fig. 4a, b). Being along the middle of the North Pacific subtropical gyre, velocity across the zonal section was not driven by the large-scale current system but by mesoscale disturbances with moderate amplitudes (Fig. 4c; see also Fig. 2c).

Fig. 4
figure 4

Same as Fig. 3 but for the horizontal velocity (color shading) as observed by shipboard acoustic Doppler current profiler: (a and b) eastward velocity along the meridional sections, and (c) northward velocity along the zonal section. Gray color shows that speed was less than 5 cm s–1, and the white areas indicate that the data were not available

Full size image

3.2 Turbulence in the subsurface pycnocline

The turbulent kinetic energy dissipation rate 〈ε〉 in the subsurface pycnocline below the surface mixed layer shows large variability along both meridional and zonal sections (Fig. 5). It is notably enhanced at the equator, exceeding 10–7 W kg–1, while relatively high values of 〈ε〉 (3–30 × 10–9 W kg–1) are also intermittently observed within the pycnocline over wide areas. We found occurrences of narrow layers of elevated turbulence (typically with < 10 m and > 10–8.5 ~ 3 × 10–9 W kg–1) away from the equator at 10°S, 10°N, and 50°N along 170°W (Fig. 5a, b), and 150°E, 170°W, 160°W, and 120°W along the subtropical North Pacific (Fig. 5c). As seen in Fig. 1, these stations, except for 120°W, located around the prominent topographic features (10°S: Samoan Passage, 10°N: Mid-Pacific Seamounts, 50°N: Aleutian Ridge, 150°E: Marcus-Wake Seamount Group, 170°W: Mid-Pacific Seamounts and 160°W: Hawaiian Ridge; Fig. 1).

Fig. 5
figure 5

Same as Fig. 3 but for turbulent kinetic energy dissipation rate (color shading): (a and b) along the meridional sections, and (c) along the zonal section

Full size image

The elevation of 〈Kρ〉 basically followed that of 〈ε〉 (Fig. 6). It exceeds 10–4 m2 s–1 at the equator, and moderately high values (10–5 to 10–4 m2 s–1) are observed in areas mentioned above. The difference between 〈Kρ〉 and 〈ε〉 caused by the 〈N2〉 variation was typically found both at the layers of peak 〈N2〉 (decrease 〈Kρ〉) and those of mode waters (increase 〈Kρ〉) over 30–40°N, 40°S, and 130°W.

Fig. 6
figure 6

Same as Fig. 3 but for vertical eddy diffusivity (color shading): (a and b) along the meridional sections, and (c) along the zonal section

Full size image

The ensemble and vertical mean dissipation rate \(\overline{\langle \varepsilon \rangle }\) and diffusivity \(\overline{\langle {K}_{\rho }\rangle }\) are shown in Fig. 7 and Table 2. We defined stations with elevated turbulence as those in the upper 25% of ε in the subsurface pycnocline of all stations (i.e., values ≥ 3.9 × 10–9 W kg–1, thick dashed lines in Fig. 7a, b), excluding those at the equator where extremely high ε was observed, and we defined stations with values below the 25% level as base-level stations. The estimated \(\overline{\langle \varepsilon \rangle }\) at the base-level stations (open symbols in Fig. 7a, b) range from 0.96 × 10–9 to 3.6 × 10–9 W kg–1, while those at stations of elevated turbulence (solid symbols) are 4.9 × 10–9 to 14 × 10–9 W kg–1 (except for the equator). The \(\overline{\langle \varepsilon \rangle }\) at the equator exceeds 10–7 W kg–1, although the confidence interval is wide because of the patches of extremely high ε (Fig. 5).

Table 2 Mean dissipation rate and diffusivity in the subsurface pycnocline
Full size table
Fig. 7
figure 7

Distributions of mean turbulent energy dissipation rate (a and b) and mean vertical eddy diffusivity (c and d) averaged over the subsurface pycnocline. Error bars indicate 95% confidence intervals. The thick dashed line shows the upper 25% level of ε; stations with \(\overline{\langle \varepsilon \rangle }\) elevated above this level are shown with solid symbols

Full size image

The \(\overline{\langle {K}_{\rho }\rangle }\) at the base-level stations (open symbols in Fig. 7c, d) fell within the range 0.14–0.47 × 10–5 m2 s–1. Even the upper bound of the confidence intervals did not reach 10–5 m2 s–1, except for the station with Central Mode Water at 40°N, 170°W. The elevation of \(\overline{\langle {K}_{\rho }\rangle }\) from the base level was defined in a similar way to that for \(\overline{\langle \varepsilon \rangle }\). The \(\overline{\langle {K}_{\rho }\rangle }\) of the elevated stations away from the equator was 0.51–1.4 × 10–5 m2 s–1, and 13 × 10–5 m2 s–1 at the equator.

3.3 Internal tide, wind, stratification and shear-to-strain ratio

In this subsection, we compare the geographical distributions of the mixing intensity with those of the forcing and background condition. The two forcing factors Dtide and Ewind, and \(\overline{\langle {N}^{2}\rangle }\) scaling velocity from internal waves (e.g., Leaman and Sanford, 1975) that could be forced by tide and/or wind are shown in Fig. 8, and the shear-to-strain ratio R that is indicative of the regime of the fine-scale variability is shown in Fig. 9.

Fig. 8
figure 8

Distributions of total dissipation rate caused by the internal tide Dtide according to the model of de Lavergne (2019) (a and b), turbulent energy production by the wind Ewind calculated from the second Modern-Era Retrospective analysis for Research and Applications (MERRA2; Gelaro et al. 2017) (c and d) and ensemble and vertical average of stratification \(\overline{\langle {N}^{2}\rangle }\) (e and f). Solid symbols indicate elevated turbulence as shown in Fig. 7

Full size image
Fig. 9
figure 9

Shear-to-strain ratio R calculated for the pycnocline below the mixed layer (mixed layer + 10 m – mixed layer + 510 m). Error bars indicate 90% confidence intervals. Solid symbols indicate elevated turbulence as shown in Fig. 7

Full size image

Although the depth of the seafloor at the observation stations exceeds 3000 m, except for ~ 2800 m at 20°N, 170°W near the Mid-Pacific Seamounts (Fig. 1 and Table 1), there is some correspondence between the model-based vertically integrated dissipation rate Dtide and vertical mixing (Fig. 8a, b). The enhanced vertical diffusivity \(\overline{\langle {K}_{\rho }\rangle }\) near the Hawaiian Ridge and the Mid-Pacific Seamounts at 170°W and 160°W along the zonal transect (Fig. 8b) is likely caused by the internal tide, where the model of de Lavergne et al. (2019) predicts high Dtide. The relatively high Dtide is also predicted at 10°S (Fig. 8a) and 150°E (Fig. 8b). However, despite high Dtide at 15°S and over 15–24°N along 170°W (Fig. 8a), and 137°E and 140°E along 23°N (Fig. 8b), turbulence is not elevated (Fig. 7). As seen in Figs. 5, 6, patches of elevated turbulence are observed around the lower bounds of, and in layers deeper than, the averaging bin (10 to 160 m below the mixed layer) at these stations, suggesting weaker penetration of the internal wave energy to the upper layer.

The turbulent energy production by wind (Ewind) averaged over 10 days at the observation stations was moderate with a level comparable to or lower than 10–3 W m–2 (corresponding to a wind stress magnitude of ~ 0.1 N m–2; Fig. 8c, d). While the estimates change slightly for different averaging times (1 day or 30 days; not shown), the pattern is similar to the climatology of wind power on inertial motions during summer (e.g., Alford 2003). Among the stations where turbulence is elevated but not associated with the high Dtide, the wind is strong at 10°N and the equator (during cruise KH-13-7 in austral summer; Fig. 8c), where the mixed layer developed close to 100 m (Fig. 8c). In other areas, the level of Ewind is not apparently linked to the turbulence level.

As seen in Fig. 3, a shallow permanent pycnocline appeared just below the surface mixed layer in the equatorial areas. This resulted in a clear latitudinal pattern of \(\overline{\langle {N}^{2}\rangle }\), highest in the equatorial areas and showing a decreasing trend with latitude (Fig. 8e, f). There was a significant positive relationship between \(\overline{\langle {N}^{2}\rangle }\) and \(\overline{\langle \varepsilon \rangle }\) (r = 0.60, p = 0.00, without the data at the equator showing the exceptionally strong mixing), while the relationship between \(\overline{\langle {N}^{2}\rangle }\) and \(\overline{\langle {K}_{\rho }\rangle }\) was not significant (r = –0.06, p = 0.76, without the data at the equator).

The shear-to-strain ratio R estimated for the 500-m bin below the mixed layer showed a dome-shaped pattern with respect to latitude along 170°W (Fig. 9a) section and an increasing trend from west to east along the subtropical North Pacific section (Fig. 9b). Although the confidence intervals are wide (error bars in Fig. 9 are 5–95% range), significant differences (no overlapping between 2.5–97.5%, not shown in the figure) are observed between those in 15–25°S (high) and 30–40°S (low), and 15–25°N, 35°N (high) and 5–10°N, 45–50°N (low). The meridional and zonal distributions of R have significant positive relationships with that of Dtide (Fig. 8a; log–log linear correlation, r = 0.53, p = 0.02) and Ewind (Fig. 8d; log–log-linear correlation, r = 0.79, p = 0.00), respectively. However, neither R-Ewind relationship along the meridional section and R-Dtide relationship along the zonal section shows significant linear relationship (log–log-linear correlation, r = –0.00, p = 0.99 and r = –0.43, p = 0.18, in Fig. 8b, c, respectively). The R distributions themselves show no correspondence with the stations of the strong turbulence (solid symbols).

4 Discussion and conclusions

Various observations have suggested that the background vertical eddy diffusivity is O(10–5) m2 s–1 (Gregg 1987; Ledwell et al. 1993; Lewis et al. 1986). Our estimates of \(\overline{\langle {K}_{\rho }\rangle }\) at the base-level stations were significantly smaller at 0.14–0.47 × 10–5 m2 s–1, except for one station (Fig. 7c, d and Table 2). We assume that our estimates are relatively low primarily because our observations were conducted only during summer when atmospheric conditions are calm (Alford 2003; Inoue et al. 2017) and subsurface stratification is developed. Whalen et al. (2018) showed that diffusivity is significantly increased from summer to winter, especially in areas with high eddy kinetic energy, where the ratio may be 2–3 at the main pycnocline (250–2000 m). A seasonal difference is also indicated by observations in the eastern subtropical South Pacific. Fernandez-Castro et al. (2014) and Bouruet-Aubertot et al. (2018) made observations during April–May (30–300 m) and February–April (100–500 m), yielding average diffusivities of 3.0 and 0.45 × 10–5 m2 s–1 (ratio of ~ 6.7), respectively. Another possible explanation is the difference in the depth of layers. While we focused on the upper part of the pycnocline where nutrient transport is critical for the production, some studies on the mixing in the pycnocline used the data of slightly deeper layers; e.g., the experiment by Ledwell et al. (1993) injected the tracer into 310 m depth. The stronger stratification in the subsurface layers would cause the lower estimate of the diffusivity even if the energy dissipation rates were at a similar level.

At the equator, the layer we focused on was well below the surface mixing layer, likely corresponding to the upper core layer of the EUC, where diurnal variability is weak (Inoue et al. 2012). Previous studies have suggested that the shear of the EUC enhances turbulence (Inoue et al. 2012; Fernández-Castro et al. 2014). However, as our observations are limited in number (3 profiles) and did not penetrate the core of the EUC, we can neither determine the diffusivity level nor identify the underlying mechanisms at the equator.

While our observations showed elevated turbulence around prominent topography, they were made in the subsurface pycnocline above 300 m at stations with seafloor depth mostly greater than 3000 m. Therefore, the elevation is unlikely to be attributable to the direct breaking of internal waves (St Laurent et al. 2002), but to the interactions of internal waves either generated, reflected, or modulated by the topography. The level of diffusivity near the Hawaiian Ridge (160°W) is consistent with previous studies, in that it was one order of magnitude smaller than that on the ridge but still significantly higher than the background level (Fig. 7d; Klymak et al. 2006; Rudnick et al. 2003). The enhancement of the vertical mixing at 21°30′N, 170°W and 23°N, 150°E (Figs. 7 and 8) is likely related to the internal tide generated by the Mid-Pacific Seamounts and Marcus-Wake Seamount Group (Fig. 1). We assume that the overall effects of these topographic features have significant impacts on the vertical heat and material budget of the subsurface pycnocline in the open Pacific. However, the vertical diffusivity is not always elevated in areas with high Dtide, such as those around 15°S and 15–24°N along 170°W, and 137 and 140°E along 23°N (Figs. 7, 8). One explanation for this apparent discrepancy is that less internal wave energy penetrates into the upper layer from the bottom. Patches of elevated turbulence are mainly distributed in the sharp pycnocline, and the pycnocline at some stations without enhanced turbulence, such as those at 137°E, 15°S and 20°N, was relatively deep. This assumption is consistent with the positive relationship between \(\overline{\langle {N}^{2}\rangle }\) and \(\overline{\langle \varepsilon \rangle }\).

If the fine-scale variances are dominated by linear internal waves with a single frequency, the shear-to-strain ratio indicates the frequency of the waves: high for near-inertial waves and low (close to 1) for higher frequency waves such as semidiurnal internal waves in the low- and mid-latitudes (e.g., Chinn et al. 2016; Kunze et al. 2006). Since our estimates of R fell below 1 in some areas, the above assumption on the frequency content does not hold. Nevertheless, the geographic pattern of R, dome shape with respect to latitude (Fig. 9a) and the increasing trend from west to east (Fig. 9b) suggests that the forcing and/or background conditions formed some kind of regimes of the fine-scale variability. Although the significant linear relationship between R and Ewind along the zonal section suggests the sensitivity of R to the near-inertial internal waves in the middle of the North Pacific subtropical gyre, the suggestion does not hold for the meridional section along 170°W. Moreover, the positive relationship between R and Dtide along the meridional section is opposite to the expectation that increasing internal tide energy lowers the R value. Because the shear-to-strain ratio is under the complex influences of internal tides and wind variability (Chinn et al. 2016) and also latitude (e.g., MacKinnon et al. 2013) and our focus is the subsurface pycnocline just below the surface mixed layer, clarification of the underlying mechanisms needs further accumulation of evidences.

In the present study, we have measured the geographical variability of the vertical eddy diffusivity in the subsurface pycnocline. The diffusivity levels and geographical distributions presented in this study improve estimations of the material transport into oligotrophic euphotic layers and will thus contribute to better understanding biogeochemical cycles. As for the other sources of the nitrogen for the new production, the role of biological nitrogen fixation has increasingly been highlighted (e.g., Horii et al. 2018). An evaluation of the vertical turbulent nutrient transport from the lower layer will be conducted in a future study, which will refine the definition of ocean biogeography based on the source type of the nitrogen for the new production.