A Parameterization of Local and Remote Tidal Mixing

Abstract

Vertical mixing is often regarded as the Achilles' heel of ocean models. In particular, few models include a comprehensive and energy‐constrained parameterization of mixing by internal ocean tides. Here, we present an energy‐conserving mixing scheme which accounts for the local breaking of high‐mode internal tides and the distant dissipation of low‐mode internal tides. The scheme relies on four static two‐dimensional maps of internal tide dissipation, constructed using mode‐by‐mode Lagrangian tracking of energy beams from sources to sinks. Each map is associated with a distinct dissipative process and a corresponding vertical structure. Applied to an observational climatology of stratification, the scheme produces a global three‐dimensional map of dissipation which compares well with available microstructure observations and with upper‐ocean finestructure mixing estimates. This relative agreement, both in magnitude and spatial structure across ocean basins, suggests that internal tides underpin most of observed dissipation in the ocean interior at the global scale. The proposed parameterization is therefore expected to improve understanding, mapping, and modeling of ocean mixing.

1 Introduction

When tidal currents flow over uneven seafloor, they generate inertio‐gravity waves called internal tides (Bell, ). Internal tides have long been suspected to play an important role in mixing the deep ocean (Munk, ). Accumulating observations from the past three decades suggest that tides power much of the small‐scale turbulence responsible for irreversible mixing, not only in the deep ocean (Ledwell et al., ; Polzin et al., ) but also in the upper ocean (Hibiya & Nagasawa, ; Kunze, ; Kunze et al., ; Polzin, ; Whalen et al., ). However, mapping tidal mixing at the global scale has proven arduous because of the difficulty in measuring small‐scale turbulence and in attributing the observed turbulence to specific dissipation pathways (MacKinnon et al., ; Waterhouse et al., ). As a result, extant parameterizations of tidal mixing suffer from substantial simplifications and uncertainties in the specified distribution of internal tide energy dissipation. These uncertainties limit our understanding of the drivers, structure, and climatic functions of the overturning circulation (de Lavergne et al., ; Melet et al., ; Sigman et al., ).

Internal tides radiate from the seafloor with a variety of spatial scales. Small‐scale internal tides, referred to as high‐mode internal tides, tend to break into small‐scale turbulence close to generation site (St Laurent & Garrett, ). By contrast, large‐scale or low‐mode internal tides can travel hundreds to thousands of kilometers and fuel dissipation remote from generation site (Dushaw et al., ; Ray & Mitchum, ). Most tidal mixing parameterizations in use in ocean general circulation models (OGCMs) only consider high‐mode internal tides. They rely on a map of internal tide generation and posit that one‐third of this energy source feeds high‐mode waves that dissipate in the local water column (St Laurent et al., ). The remaining two‐thirds of power input are ignored or surmised to participate in sustaining a constant background diffusivity of order 10⁻⁵ m² s⁻¹ (Simmons et al., ). This approach has two principal limitations: (i) the fraction of local dissipation is not uniform and actually depends on resolution and location (Falahat et al., ; St Laurent & Nash, ; Vic et al., ) and (ii) constant background diffusivities disallow energy conservation and do not do justice to observed patterns of mixing rates (de Lavergne et al., ; Eden et al., ; Pollmann et al., ).

Several recent studies tackled the limitation (ii) by explicitly including mixing powered by low‐mode internal tides. Oka and Niwa () employed a static horizontal map of low‐mode dissipation derived from high‐resolution numerical experiments (Niwa & Hibiya, ). Assuming that this dissipation is uniform in the vertical, they assessed the impact of adding remote tidal mixing in an OGCM. Eden and Olbers () proposed an interactive parameterization of low‐mode energy propagation and dissipation. They introduced an equation for the evolution of low‐mode energy within an OGCM. A map of low‐mode internal tide generation, parameterized attenuation rates, and a simple model for the vertical dependence of dissipation then allow solving for the evolving tidal mixing rates. These studies demonstrate the feasibility and importance of replacing background diffusivities by an explicit parameterization of remote tidal mixing. Further advance requires improvements in the realism of the modeled distribution of internal wave energy loss. The vertical structure of dissipation needs particular attention, as it is crucial for ocean ventilation and as it depends on the process causing dissipation (de Lavergne et al., ; Melet et al., ).

To accurately parameterize mixing energized by low‐mode internal tides, it is thus necessary to track where and how they dissipate. Using Lagrangian tracking of internal tide energy beams through observed stratification, de Lavergne et al. () recently estimated column‐integrated internal tide dissipation rates decomposed into contributing processes. Here, we make use of these horizontal maps and of historical microstructure observations to propose a comprehensive and energy‐constrained parameterization of tidal mixing. The parameterization explicitly accounts for the local and remote dissipation of internal tides and obviates assumptions about the fraction of local dissipation. It relies on four static maps of column‐integrated dissipation, each associated with a distinct process and related vertical structure (section ). To gauge the realism of the parameterization, we apply it to an observational hydrographic climatology and compare the obtained three‐dimensional distribution of mixing to a compilation of observational mixing estimates (section ). We then document the inferred global budget of internal tide‐induced dissipation (section ) and conclude with a summary of implications and limitations of the analysis (section ).

2 Process‐Dependent Vertical Structures

The vertical mode of internal tides is a key determinant of their propagation and dissipation characteristics (St Laurent & Garrett, ). The lowest modes have elevated group speeds and relatively slow rates of attenuation by wave‐wave interactions (Olbers et al., ; Onuki & Hibiya, ). They consequently tend to lose much of their energy through interactions with topography (Kelly et al., ). Conversely, high vertical modes tend to lose most of their energy through wave‐wave interactions within the near‐local water column (Nikurashin & Legg, ). These different fates of internal tides were mapped globally for each mode and each of the main three tidal constituents (M₂, S₂, and K₁), using the following ingredients (Figure ; de Lavergne et al., ): (i) an observational climatology of stratification (Gouretski & Koltermann, ); (ii) estimates of internal tide generation projected onto vertical modes (Falahat et al., ; Melet et al., ); (iii) simplified representations of energy sinks (Bühler & Holmes‐Cerfon, ; Hazewinkel & Winters, ; Legg, ; MacKinnon et al., ; Olbers, ); and (iv) a Lagrangian energy tracker.

Diagram summarizing the approach. This study extends the dissipation maps from two to three dimensions by applying a vertical structure appropriate to each process‐specific map of column‐integrated internal tide energy loss. The vertical structures involve the buoyancy frequency N as well as parameters (Hcri, Hbot, rbot) determined by comparison to microstructure observations.

Modes 6 and higher, which cumulate a power of 217 GW, were found to dissipate locally at the half‐degree resolution of the calculation. These locally dissipating modes comprise Modes 6 to 10 (115 GW Falahat et al., ) and much higher modes generated by abyssal hills (102 GW Melet et al., ). The dissipation of Modes 1 to 5 was split into four processes (see de Lavergne et al.,  for detailed methods): wave‐wave interactions (521 GW), incidence on critical topographic slopes (128 GW), shoaling (95 GW), and scattering by abyssal hills (83 GW). Here we organize these contributions to the overall internal tide energy loss (1,044 GW) into four components, based on expectations about induced vertical structures of turbulent kinetic energy production (Figure ):

• E_wwi: attenuation of low modes (1–10) by wave‐wave interactions;

• E_sho: direct breaking of low‐mode waves through shoaling;

• E_cri: low‐mode waves dissipating at critical slopes;

• E_hil: scattering of low‐mode waves by abyssal hills and generation of high‐mode waves by abyssal hills.

The four components are mapped in Figure . They constitute depth‐integrated internal tide energy dissipation rates, with units of W m⁻². Generation and scattering by abyssal hills (E_hil) is most intense along ridges of the Atlantic and Indian basins while nonnegligible throughout most of the open ocean (Figure d). Abyssal hills are dominant features of the ocean floor at horizontal scales smaller than 10 km (Goff, ; Macdonald et al., ). They are thought to be responsible for the bulk of bottom‐intensified mixing above rough ridges (Lefauve et al., ; Muller & Bühler, ; Nikurashin & Legg, ; Polzin, ). The remaining three components encapsulate the dissipation of the first 10 vertical modes (Figures a–c). Modes 6 to 10 are assumed to dissipate through wave‐wave interactions (E_wwi) only, consistent with the overwhelming contribution (99%) of this process to Mode 5 dissipation (de Lavergne et al., ). E_wwi features widespread, strong depth‐integrated dissipation and dominates the overall budget (Figure a). A floor of 10⁻⁵ W m⁻² was imposed on E_wwi to maintain a minimal amount of mixing poleward of the S₂ turning latitude (85.8°); this floor increased total dissipation by 0.1 GW. Shoaling (E_sho) mostly acts at the shelf break and shoreward, so that the induced turbulence is confined to relatively shallow waters (Figure b). Dissipation of low‐mode waves impinging on critical topography (E_cri) occurs primarily at continental slopes (Figure c).

Global maps of the static power fields that enter the parameterization. The color scale shows the base‐10 logarithm of power values. (a) Ewwi: attenuation by wave‐wave interactions of Modes 1–10, comprising 521 GW from Modes 1–5 and 115 GW from Modes 6–10. (b) Esho: direct breaking of low‐mode waves through shoaling (95 GW). (c) Ecri: dissipation of low‐mode waves at critical slopes (128 GW). (d) Ehil: scattering of low‐mode waves by abyssal hills (83 GW) and generation of high‐mode waves by abyssal hills (102 GW). (e) Etid=Ewwi+Esho+Ecri+Ehil.

Power contained in the four static maps must now be distributed in the vertical to obtain a three‐dimensional distribution of the production rate of turbulent kinetic energy. In all the following, this rate will be denoted ϵ and referred to as turbulence production for brevity. Turbulence production (with units of W kg⁻¹) can be related to diffusivity K_ρ, buoyancy frequency N, and frictional heat production ϵ_ν through the simplified balance (Osborn, ): [Image Omitted. See PDF]

Introducing the flux Richardson number or mixing efficiency R_f, the buoyancy flux term can be expressed as K_ρN²=R_fϵ. The constant value R_f=1/6 is generally assumed in stratified oceanic conditions (Gregg et al., ). R_f should be distinguished from the flux coefficient Γ=K_ρN²/ϵ_ν=R_f/(1−R_f), used to deduce K_ρ from microstructure measurements of ϵ_ν, and usually taken to be one‐fifth.

2.1 Low‐Mode Components

Observational and theoretical evidence indicates that turbulence production fueled by the downscale cascade of low‐mode wave energy through wave‐wave interactions scales with the square of the buoyancy frequency (Gregg, ; Henyey et al., ; Müller et al., ; Polzin et al., ; Polzin, ). Accordingly, we define the local power density ϵ_wwi (with units of W kg⁻¹) as [Image Omitted. See PDF] where ρ denotes the in situ density, z height, and $\int dz$ an integral over the full water column.

Shoaling‐induced turbulence often occurs via direct breaking of low‐mode waves, so that the structure of turbulence production roughly matches the profile of the wave energy, proportional to N (D'Asaro & Lien, ; Legg, ; Melet et al., ). We choose [Image Omitted. See PDF]

Incidence on critical slopes triggers boundary turbulence along the slope (Eriksen, ; Legg & Adcroft, ; Moum et al., ; Nash et al., ). The along‐slope concentration of turbulence is not easily mimicked in a coarsely resolved ocean with step‐like topography. To roughly capture the vertical extent and bottom intensification of the turbulence production, we use an exponential decay from the seafloor with an e‐folding scale H_cri equal to the along‐slope height difference (Figure a). The height scale for a given grid column (i,j) is thus defined as the largest local topographic rise, [Image Omitted. See PDF] with H>0 the bathymetry. Where the so defined H_cri is negative, it is set instead as the subgrid‐scale maximum local topographic rise (see Appendix A in de Lavergne et al., ). Turbulence production is then given by [Image Omitted. See PDF] where h_ab denotes the height above bottom. This exponential vertical structure implies that the bulk of turbulent mixing occurs along the slope, but that some reaches higher up, in accord with results of process studies (Legg & Adcroft, ; Legg, ).

Global maps of the two decay scales, (a) Hcri and (b) Hbot, that enter the parameterization. Panel (c) shows Arms normalized by its global average.

2.2 Abyssal Hills

Several vertical structures have been proposed for the dissipation of internal tides excited by abyssal hills (Lefauve et al., ; Melet et al., ; Polzin, , ). To assess these structures, we focus on a densely sampled area of the eastern Brazil Basin: 21–22°S; 16–19°W. During two cruises of the Brazil Basin Trace Release Experiment (BBTRE), in 1996 and 1997, 51 microstructure profiles were collected in this area (Polzin et al., ). Following Polzin (), we first construct a composite profile of turbulence production for the region by averaging microstructure profiles in depth coordinate over the upper 3 km and in height above bottom coordinate over the bottom 1.5 km. The factor (1−R_f)⁻¹=6/5 is used to convert from the measured viscous dissipation to turbulence production. The resulting 4.5‐km deep profile, shown in Figure , implies a depth‐integrated power consumption of 5.5×10⁻³ W m⁻². Somewhat coincidentally, the predicted internal tide energy loss, E_tid=E_wwi+E_sho+E_cri+E_hil, averaged over the same area is 5.6×10⁻³ W m⁻². Only E_hil and E_wwi are important contributors here, amounting to 3.5×10⁻³ and 2.0×10⁻³ W m⁻², respectively.

Parameterized turbulence production in the eastern Brazil Basin (21–22°S; 16–19°W). (a) Composite stratification profile. (b) Composite microstructure profile (black) and parameterized ϵhil profile (blue) for Hbot = 150 m and three values of rbot. An exponential decay with an e‐folding scale of 500 m (St Laurent et al., ) is shown in gray. (c) Same as (b), with the total parameterized profile ϵhil+ϵwwi in blue. (d) Wentzel‐Kramers‐Brillouin (WKB)‐stretched height above bottom as defined by Lefauve et al. () and applied to the composite stratification profile. (e,f) Same as (b,c), with ϵhil distributed according to Lefauve et al. () and three values of HL. (g) WKB‐stretched height above bottom as defined by Melet et al. (). (h,i) Same as (b,c), with ϵhil distributed according to (Melet et al., ) and three values of HM. In all panels, thin black curves delimitate 95% confidence intervals from bootstrapping.

Using a simplified theory for the decay of hill‐generated internal tides within uniform stratification, Polzin () predicted that turbulence production ϵ decreases with h_ab following [Image Omitted. See PDF] where H_bot is a height scale that depends on bottom conditions. Given H_bot = 150 m and ϵ₀=10⁻⁸ W kg⁻¹, this structure provides a good fit to the Brazil Basin composite microstructure profile within 1.5 km of the bottom (Polzin, ). Here, we must express ϵ_hil as a function of the depth‐integrated power E_hil instead of the bottom ϵ₀. To this end, we introduce the parameter r_bot, which represents the fraction of E_hil that dissipates close to the bottom, that is, following equation (6). We obtain a best fit to the bottom‐intensified portion of the composite microstructure profile using r_bot=0.86, together with H_bot = 150 m (Figure b). The remainder of the power, (1−r_bot)E_hil, feeds small‐scale turbulence in the shallower portion of the water column. We choose to distribute this remainder as proportional to N², consistent with turbulence arising from weakly nonlinear wave‐wave interactions (Polzin, ). The complete parameterized vertical structure for ϵ_hil is thus given by [Image Omitted. See PDF]

The factor $\frac{1}{H}+\frac{1}{H_{bot}}$ comes from normalization of the depth‐integrated power. When adding ϵ_wwi, this parameterization tracks the observational profile over the full water column (Figure c). Quantitatively, the percentage power contained in each 750‐m layer of the parameterized profile is within 1% of its observational counterpart (Table ).

Vertical Distribution of Power in the Composite Observational Profile (21–22°S; 16–19°W) and Three Different Parameterizations

Depth	Obs.	Equation (7)	Lefauve et al. ()	Melet et al. ()
range	composite	(rbot = 0.86, Hbot = 150 m)	(HL = 400 m)	(HM = 5 m)
(m)	(%)	(%)	(%)	(%)
0–750	33.1	34.1	29.3	34.6
750–1,500	8.0	7.0	10.4	7.8
1,500–2,250	3.7	2.7	5.8	2.4
2,250–3,000	2.2	2.3	3.5	0.9
3,000–3,750	4.9	5.5	17.9	3.4
3,750–4,500	48.1	48.4	33.1	51.0

Note. The distribution is described by the percentage power contained in each 750‐m‐thick layer. All three parameterizations include ϵ_hil+ϵ_wwi. Only the vertical structure of ϵ_hil changes according to parameterization. The present parameterization of ϵ_hil uses r_bot = 0.86 and H_bot = 150 m. The structure proposed by Lefauve et al. () is applied with H_L = 400 m; that proposed by Melet et al. () is applied with H_M = 5 m.

We assess two alternative parameterizations (Figures d–i; Lefauve et al., ; Melet et al., ). They apply to ϵ_hil only (Figures e and h), so that ϵ_wwi is distributed according to equation (2) in every ϵ_hil+ϵ_wwi profile (Figures f and i). Both parameterizations resort to a stratification‐weighted height above bottom, [Image Omitted. See PDF] also known as Wentzel‐Kramers‐Brillouin stretching. Here m is an exponent equal to 1 (Figure d) or 2 (Figure g). Lefauve et al. () proposed a profile of turbulence production proportional to $N^2\exp (-h_{wkb,1}/H_L)$, where H_L is a height scale. As shown in Figure f and Table , the best fit obtained for H_L = 400 m is unsatisfactory within the bottom portion of the water column: the structure is strongly influenced by the local maximum in stratification near 3.7 km depth (Figure a). Recasting the structure of Polzin () in an energy‐conserving form, Melet et al. () proposed a profile proportional to N²(1+h_wkb,2/H_M)⁻², with H_M an appropriate decay scale. This structure provides a reasonable fit to the composite microstructure profile for H_M = 5 m (Figure i and Table ). However, turbulence production is underestimated at middepth, near the local minimum in stratification. On the whole, parameterization (7) best mimics the measured dissipation across the water column.

Seasonal change of parameterized turbulence production in the eastern Brazil Basin (21–22°S; 16–19°W). (a) Stratification profile representative of (blue) winter and (red) summer, from World Ocean Atlas 2018 (Locarnini et al., ; Zweng et al., ). (b) Parameterized ϵhil profile in (red) summer and (blue) winter. The present parameterization (Hbot = 150 m, rbot=0.86) is shown by thick curves. For comparison, we also show (dashed‐dotted curves) the vertical structure of Lefauve et al. () with HL = 400 m and (thin curves) that of Melet et al. () with HM = 5 m. Blue and red triangles in the bottom right corner highlight the bottom value of blue and red thin curves, respectively.

Parameterization (7) has an additional important characteristic: it ensures that near‐bottom turbulence levels are decoupled from upper‐ocean changes in stratification. This is not the case of structures proposed by Melet et al. () and Lefauve et al. (), according to which an increase in near‐surface stratification implies a reduction of the power input to abyssal mixing. The latter behavior is illustrated in Figure by applying the three parameterizations to climatological summer and winter profiles of stratification in the region 21–22°S; 16–19°W. In summer, larger buoyancy frequencies in the upper 400 m of the water column lead to a decrease of turbulence production in the bottom kilometer using the structure of Lefauve et al. () and in the bottom grid cell using that of Melet et al. (). A long‐term increase in upper‐ocean stratification, as may occur under a climate change scenario, would similarly drain mixing energy from the abyss. In contrast, parameterization (7) predicts unchanged abyssal mixing under unchanged abyssal conditions, as expected for turbulence driven by the rapid decay of bottom‐generated small‐scale internal tides (Polzin, , ).

Properties of parameterization (7) rely on the introduction of two adjustable parameters, H_bot and r_bot. These parameters are expected to vary geographically: more energetic or smaller scale waves should decay more rapidly above the bottom, implying smaller H_bot or higher r_bot or both (Lefauve et al., ). Lefauve et al. () found that internal tide energy dissipation above abyssal hills is largely shaped by the root‐mean‐square amplitude of internal tides at the bottom, A_rms. This characteristic wave amplitude is related to the power input E_hil and the mean local wavenumber of abyssal hills k_hil via the scaling $A_{rms}\propto E_{hil}^{1/2}{k}_{hil}^{3/2}$. (Scaling equation (7) in Lefauve et al. () gives $A_{rms}\propto C{(k_{hil})}^{1/2}J(k_{hil})k_{hil}^2$, where C is the two‐dimensional spectrum for small‐scale bathymetry and J is the Bessel function of the first kind of order one. Scaling equation (7) in Polzin (), and omitting tidal harmonics, gives $E_{hil}(k_{hil})\propto C(k_{hil})J^2(k_{hil})k_{hil}^{-1/2}$. Combining the two scalings gives $A_{rms}\propto E_{hil}^{1/2}{k}_{hil}^{3/2}$.) We mapped A_rms using k_hil estimated by Goff (). In general, A_rms decays away from the crest of the main mid‐ocean ridges (Figure c). The eastern Brazil Basin has relatively short‐scale abyssal hills and relatively large A_rms (Figure ). Hence, on average, turbulence induced by abyssal hills should be less concentrated near the bottom than observed in the 21–22°S; 16–19°W region.

Distribution of the global integrated power, ∫∫Ehildxdy, as a function of (a) Ehil, (b) khil, and (c) Ehil1/2khil3/2∝Arms. The red line shows the eastern Brazil Basin (21–22°S; 16–19°W) value.

To account for this variability, we test a range of scenarios defined by $H_{bot}\propto A_{rms}^l$ and $r_{bot}\propto A_{rms}^p$, with l∈{−2,−1,−0.5,0} and p∈{−0.5,0,0.5,1}. Combined to the values H_bot = 150 m and r_bot=0.86 appropriate to the eastern Brazil Basin, the scenarios define global geographies of the two parameters. We then apply the complete parameterization to the World Ocean Circulation Experiment (WOCE) annual mean climatology of stratification (Gouretski & Koltermann, ) and compare the obtained distribution of ϵ_tid=ϵ_hil+ϵ_cri+ϵ_sho+ϵ_wwi to available microstructure observations across the globe (Figure a; see section for an expanded description). In particular, we sample the parameterized ϵ_tid along microstructure profiles and compare project‐average profiles (Figure ). We find that $H_{bot}\propto A_{rms}^{-1}$ (Figures b and c) and a constant r_bot generate the best overall agreement. Specifically, a steeper decrease of H_bot with A_rms or an increase of r_bot with A_rms promotes underestimation (overestimation) of the near‐bottom dissipation in sampled regions where A_rms is relatively weak (strong). Zonal and vertical patterns of measured dissipation across the Brazil Basin provide the strongest constraint. The scenario $H_{bot}\propto A_{rms}^{-1}$ and r_bot = 0.86 gives a reasonable match between the parameterization and BBTRE across the basin and throughout the water column (Figures and ). In the following, we retain this scenario exclusively and document the resultant climatological distribution of ϵ_tid.

Global (a) microstructure, (b) ship‐based finestructure (Kunze, ), and (c) Argo‐based finestructure (Whalen et al., ) sampling. The number of profiles for each project is shown in parentheses in (a). Finestructure profile locations are shown in blue and black in (b); black sections are those plotted in Figures , S1, and S2. Shading in (c) depicts the maximum depth (in kilometers) of Argo‐based finestructure profiles.

Project‐average microstructure profiles (black) compared to parameterized profiles (blue) of turbulence production. Parameterized profiles are obtained by sampling the global distribution of internal tide energy dissipation—which is based on the WOCE climatology of stratification and the preferred scenario rbot=0.86, Hbot∝Arms−1—at the location and depths of each available microstructure profile. If either the parameterized or measured profile is not defined at some depth, both are set to not‐a‐number at this depth before taking the average; this ensures that differences in coverage do not bias the comparison. Thin lines delimitate 95% confidence intervals from bootstrapping. Project locations are shown in Figure .

Southwestern Atlantic transect of (b) parameterized and (c) measured turbulence production. Selected measurements are BBTRE microstructure profiles identified by red markers in (a). Selected grid points for the parameterization track the observational transect. The chosen scenario is rbot=0.86, Hbot∝Arms−1. Horizontal linear interpolation has been used in seven unsampled grid columns (marked by circles) of Panel (c). In Panel (a), shading shows the “etopo2v2” bathymetry (Smith & Sandwell, ) averaged at 0.5° resolution and white markers locate BBTRE profiles not used in this figure.

3 Comparison to Microstructure and Finestructure Observations

3.1 Microstructure

Our compilation of microstructure data (Figure a) includes field campaigns described and analyzed by Waterhouse et al. () as well as data from eight additional projects: INDOMIX (Bouruet‐Aubertot et al., ); OUTPACE (Bouruet‐Aubertot et al., ); three cruises over the Izu‐Ogasawara ridge, hereafter referred to as IZU (Hibiya et al., ); DoMORE (Thurnherr et al., ); RidgeMix (Vic et al., ); OVIDE (Ferron et al., ); RREX (Petit et al., ); and PROVOLO (Fer et al., ). The compilation encompasses a total of 19 campaigns cumulating 1,171 microstructure profiles. Only two profiles, from Station 5 of the INDOMIX dataset, are excluded from the present analysis. Located in Ombai Strait, these profiles feature dissipation rates in excess of 10⁻⁴ W kg⁻¹, likely related to large internal solitary waves (Bouruet‐Aubertot et al., ), that overwhelm the regional average. Following Waterhouse et al. (), we split the BBTRE dataset into eastern (BBTREe: east of 28°W) and western (BBTREw: west of 28°W) regions (Figure a). Measured viscous dissipation rates are binned into 200‐m depth intervals and multiplied by (1−R_f)⁻¹=6/5 to obtain profiles of turbulence production. Project‐average profiles are shown in black in Figure .

Comparison with the parameterized ϵ_tid profiles (Figure , blue) yields mixed results. Good agreement is found for NATRE, BBTREw, BBTREe, RidgeMix, and DoMORE. Some overestimation of the deepest measurements from the latter three campaigns is nonetheless apparent. A reasonable match in magnitude but departures in shape are obtained for GEOTRACES, OVIDE, DIMES‐DP, HOME, LADDER, and INDOMIX. In several other cases, namely, IZU, RREX, GRAVILUCK, OUTPACE, and TOTO, the parameterization reproduces the shape but substantially overestimates the magnitude of the observational profile below 200 m. In the remaining SOFine, DIMES‐West, PROVOLO, and Fieberling areas, the parameterization predicts lower turbulence production than was measured.

Discrepancies can be explained by (i) biases of the parameterization; (ii) nontidal energy sources; (iii) application of the parameterization to a gridded annual mean stratification field that may depart from the local stratification at the time of measurements; (iv) variability of turbulence that limits representativeness of measurements; and (v) uncertainty in the measured ϵ_ν and assumed R_f. Substantial overestimates of turbulence production by the parameterization most likely stem from (i). Biased internal tide generation rates near Iceland (Lefauve et al., ) may contribute to RREX and OVIDE mismatches. A lack of energy redistribution in the model, possibly linked to overly strong attenuation by wave‐wave interactions (de Lavergne et al., ), could also explain overestimates around generation sites such as IZU. Localized biases in mapped internal tide sources and sinks likely contribute to large magnitude offsets at TOTO, GRAVILUCK, and Fieberling sites. Explanation (ii) is designated for the three Southern Ocean surveys (SOFine, DIMES‐West, and DIMES‐DP): internal waves generated by atmospheric storms or jet‐topography interactions likely dominate turbulence production in these areas (Ledwell et al., ; St Laurent et al., ; Waterman et al., ). In general, nonmodeled processes are expected to contribute to salient mismatches in shape or amplitude.

Microstructure profilers provide instantaneous local measurements of intermittent turbulence, so that spatiotemporal variability and limited sampling will also account for some differences between observations and the ϵ_tid climatology. In particular, seasonal modulation of stratification, the spring‐neap tidal cycle, and differences in seafloor depth between the 0.5° resolution climatology and local observational casts can blur the comparison. For example, analysis of the timing of RREX observations along the Reykjanes ridge crest relative to the spring‐neap cycle indicates that the measurements are biased toward neap tides on average (not shown). This sampling bias could contribute to the discrepancy with the parameterization. Microstructure‐derived dissipations also have intrinsic uncertainty, thought to be about a factor of 2, attributable in part to sensor calibration and to the translation of microscale shear spectra into frictional dissipation rates (Gregg, ; Toole et al., ). We note that a recent reinterpretation of the DoMORE raw data (see Appendix in Thurnherr et al., ) produces dissipation rates (used here) that are larger than those initially estimated (Clément et al., ), leading to better agreement with BBTRE and with the present parameterization. An extra source of uncertainty lies in the assumed mixing efficiency R_f=1/6. Were this mixing efficiency a factor of two underestimate (overestimate), turbulence production rates deduced here from microstructure observations would be underestimated by 25% (overestimated by 10%). The confounding impact of mixing efficiency variations on calibration and validation of the parameterization is therefore expected to be modest.

Viewed overall, the comparison in Figure suggests that the parameterization captures reasonably well contrasts of turbulence production across regions and between ocean layers. For example, the three order‐of‐magnitude difference between Indonesian seas (INDOMIX) and the subtropical northeast Atlantic (NATRE) is correctly reproduced. The predicted vertical structure is also realistic in most cases; exceptions tend to coincide with regions for which there is evidence that internal tides are not the main conduit to small‐scale turbulence. Hence, the constructed two‐dimensional maps (Figures and ) combined with chosen vertical structures (equations 2, 3, 5, and 7) appear to have skill in mimicking turbulence powered by internal tides. This skill is exemplified by the Brazil Basin transect shown in Figure and further endorsed by the two‐dimensional histogram of Figure a. In spite of scatter in microstructure data and variability across four orders of magnitude, turbulence levels and trends are generally comparable in the in situ and climatological datasets. In particular, 85% of the 13,733 values in the microstructure database agree with ϵ_tid within a factor of 10; 41% of values agree within a factor of 2.

Turbulence production estimated from (a) microstructure or (b–d) finestructure observations (y‐axis) versus that predicted by the present parameterization (x‐axis). Panel (b) uses the dataset of Whalen et al. (). Panels (c,d) use the dataset of Kunze (), with depths >2 km excluded in Panel (d). Color scales show the number of values in each dissipation bin; one value corresponds to one 200‐m profile segment. The top color scale applies to Panel (a); the right‐side color scale applies to other panels. Total number of values (n) and correlation coefficients (r) are given in the top left corner of each panel; outside (inside) parentheses is the correlation coefficient computed using raw (using the base‐10 logarithm of) values. All correlations are statistically significant (p value < 0.05). Thin black lines delimitate agreement within a factor of 10; the corresponding percentage of values is given at the top right of each panel.

3.2 Finestructure

Internal wave energy dissipation can be estimated from the finescale strain contained in hydrographic casts (Polzin et al., ; Whalen et al., ; Wijesekera et al., ). Such finestructure dissipation estimates carry large uncertainties because of choices and parameters involved in the inference (Pollmann et al., ; Polzin et al., ). Nonetheless, they allow extensive spatial coverage in comparison to microstructure measurements (Figure ). Here we analyze full‐depth dissipation estimates of Kunze () based on shipboard CTD casts (Figure b) and upper‐ocean dissipation estimates of Whalen et al. () based on Argo CTD profiles (Figure c). Values shallower than 380 m are considered unreliable by Kunze () due to the confounding influence of sharp pycnoclines on calculated strain. For the same reason, Whalen et al. (2015) exclude the mixed layer and mode water from the processing of individual Argo profiles. The Argo‐based dataset employed here spans the 300‐ to 1,900‐m‐depth range and consists of average profiles binned into 1°×1° grid squares.

Comparison of finestructure observational estimates to the present climatology at 400 m depth reveals strong similarities (Figure ), consistent with previous findings (de Lavergne et al., ). The western low‐latitude Pacific stands out as the most dissipative region, followed by western Indian and mid‐Atlantic waters. The imprint of ridges hosting strong internal tide generation and dissipation is visible in all three maps, in spite of the relatively shallow depth (400 m) shown here. Enhanced dissipation above the Izu‐Ogasawara ridge south of Japan is conspicuous in the present (Figure a) and Argo‐finestructure (Figure c) datasets, though less intense in the latter. This local difference in magnitude reflects a more general tendency: basin‐scale horizontal contrasts appear to be amplified in the ϵ_tid climatology relative to finestructure observations. Such a pattern amplification, previously identified in depth‐integrated dissipation rates, could reflect an underestimation of energy redistribution in the two‐dimensional mapping (de Lavergne et al., ).

Internal wave energy dissipation at 400 m depth predicted from (a) the present parameterization, (b) ship‐based finestructure observations (Kunze, ), and (c) Argo‐based finestructure observations (Whalen et al., ).

Other regional discrepancies are apparent. Whalen et al. () and Kunze () infer larger internal wave energy loss near major currents such as the Kuroshio, Gulf Stream, Agulhas Current, and Antarctic Circumpolar Current. These regions host intense mesoscale activity and above‐average near‐inertial energy input from atmospheric storms (e.g., Alford, ; Shum et al., ). It is likely that nontidal processes account for the bulk of the inferred internal wave energy dissipation there (Clément et al., ; Nikurashin et al., ; Pollmann et al., ; Waterman et al., ; Whalen et al., ) and therefore for the discrepancy with the present climatology. The dissipation map of Whalen et al. () also displays a band of enhanced turbulence along the equator (Figure c) that is absent from the ϵ_tid climatology (Figure a). The elevated dissipation rates likely signal equatorial processes unrelated to internal tides and possibly unrelated to inertio‐gravity waves (e.g., Holmes & Thomas, ; Moum et al., ). The inferred dissipation agrees with microstructure measurements at a western Pacific site (154–158°E) on the equator (Whalen et al., ). However, the equatorial band does not stand out in other finestructure estimates (Figure b; Kunze, ; Pollmann et al., ). A dedicated microstructure survey along the equator would be needed to establish whether this mixing pattern is real.

By sampling the parameterized ϵ_tid distribution along each finestructure profile, we can compare mean vertical profiles of turbulence production (Figure ). On average, turbulence production has a very similar vertical distribution in the ϵ_tid climatology and the Argo‐finestructure dataset (Figures d–f). This agreement owes much to the scaling ϵ(z)∝N²(z), valid within most water columns of the open ocean in the 0.3‐ to 2‐km‐depth range. This scaling is built in present and finestructure parameterizations but also endorsed by microstructure data (Gregg, ; Polzin et al., ). Fair agreement of ϵ_tid with the ship‐based mean finestructure profile is also observed above 1.5 km depth (Figure a–c). However, the mean ϵ_tid profile exceeds its finestructure counterpart by almost an order of magnitude deeper than 2 km (Figure a). This excess is most pronounced in regions where abyssal hills and critical slopes are dominant contributors to internal tide dissipation (Figure b). These regions are characterized by elevated deep‐ocean turbulence production catalyzed by rough or steep topography. A transect across the Atlantic near 23°S illustrates the divide (Figure ): although the present parameterization and finestructure estimates compare well in the upper 2 km, finestructure‐inferred dissipation is substantially weaker in the abyss. Indian and Pacific transects exhibit analogous similarities and differences (supporting information, Figures S1 and S2).

Average (black) finestructure‐inferred and (blue) parameterized turbulence production profile. Parameterized profiles are obtained by sampling the global distribution of turbulence production at the location and depths of each available finestructure profile. (a,d) All profiles. (b,e) Profiles where Ecri+Ehil>Ewwi+Esho. (c,f) Profiles where Ecri+Ehil<Ewwi+Esho. Panels (a–c) use the finestructure dataset of Kunze (); Panels (d–f) use that of Whalen et al. (). Thin curves delimitate 95% confidence intervals from bootstrapping. Gray shading in (a–c) marks the 0‐ to 380‐m‐depth range within which finestructure estimates were deemed unreliable by Kunze ().

South Atlantic transect, near 23°S, of (a) climatological stratification, (b) parameterized turbulence production, and (c) finestructure‐inferred turbulence production (Kunze, ). Transect location is shown in Figure b.

The tidal dissipation climatology compares more favorably with deep microstructure data than with deep finestructure inferences (Figures , , , and ). Together with a comparison of neighbor fine structure and microstructure profiles (Figure S3), this suggests that finestructure estimates of Kunze () are biased low in the deep ocean. In particular, the estimates seem to underpredict bottom‐intensified dissipation above rough or steep topography. This assessment concurs with recent studies reporting underestimated dissipation by the finestructure method in regions of rough topography or strong forcing (Bouruet‐Aubertot et al., ; Liang et al., ; Thurnherr et al., ). The identified low bias is not expected to be universal, however, as it depends on implementation choices of the method and regional processes at play (Hibiya et al., ; Kunze & Lien, ; Takahashi & Hibiya, ; Waterman et al., ). Profiles of Figures d–f also hint at a divergence near 1.7 km depth between the present parameterization and the dissipation estimates of Whalen et al. ().

In summary, the ϵ_tid climatology and finestructure observations display similar horizontal patterns and similar vertical structure in the upper ocean (Figures 11–13 and S1–S3). A histogram comparison of the global three‐dimensional distributions (Figures b–d) confirms the strong correlation between ϵ_tid and finestructure‐inferred turbulence production across five orders of magnitude. Agreement within a factor of 10 reaches 93% of values in the Argo‐based dataset and 86% of values in the ship‐based dataset. The latter percentage increases to 89% when excluding depths >2 km. These results support the realism of the present tidal mixing parameterization and suggest that internal tides largely shape the global distribution of turbulence production.

4 Global Distribution of Internal Tide‐Driven Mixing

The parameterization proposed in this study allows visualization and quantification of the global distribution of turbulence production due to internal tides. The zonal sum of ϵ_tid shows that the bulk of turbulence production takes place in the upper kilometer of low and middle latitudes (Figure a). This concentration is largely explained by the influence of stratification on internal tide generation and dissipation rates. In these strongly stratified waters, the parameterized mixing is mainly attributable to low‐mode internal tides dissipating through wave‐wave interactions (Figure b) and to shoaling‐induced wave breaking near continental margins (Figure c). Critical slopes cause turbulence relatively evenly distributed vertically in the upper 2.5 km (Figure d), where most steep continental slopes lie. Deeper than 2.5 km, tidal mixing mostly originates from internal tide generation and scattering by small‐scale topographic roughness (Figure e). The seafloor area distribution of ridges has a distinct footprint in the zonal sum of ϵ_hil, noticeable as a band of relatively high power at abyssal depths (mostly between 2.5 and 4.5 km).

Global zonal sum of (a) the total turbulence production ϵtid and (b–e) its decomposition into the four components: (b) ϵwwi, (c) ϵsho, (d) ϵcri, and (e) ϵhil.

The distribution of the parameterized internal tide energy dissipation as a function of depth or height above bottom is presented in Figures and . The top kilometer of the ocean hosts 70% of the total energy loss (Figure b), with contributions from all four components (Figure a). This leaves only 311 GW, from the total of 1,044 GW, of power input to small‐scale turbulence at depths greater than 1 km. Between 1 and 2 km depth, energy loss amounts to 123 GW, fueled mostly by modes ≤10 dissipating via wave‐wave interactions and critical slopes. Below 2 km depth, power availability drops to 18% of the total and is dominated by the ϵ_hil component.

Global distribution of ϵtid as a function of (a) depth and (b,c) height above bottom. Colors indicate respective contributions of (red) ϵwwi, (orange) ϵsho, (light blue) ϵcri, and (blue) ϵhil. In Panel (c), regions where the seafloor depth is less than 1 km have been excluded.

Global depth distribution of ϵtid in (black) this study, (blue) the standard NEMO parameterization and (orange) the standard CCSM parameterization. (a) Total turbulence production binned into 500‐m‐depth intervals. The gray shading indicates turbulence production implied by a uniform diffusivity of 10−4 m2 s−1. (b) Cumulative turbulence production, given as a percentage of the global total.

Almost 300 GW dissipate in the bottom 500 m of the ocean (Figure b), including 180 GW in the open ocean (Figure c). Components of dissipation linked to topography dwindle rapidly with height above bottom: ϵ_hil+ϵ_cri+ϵ_sho contributes little power at h_ab>2 km. This trend is opposed by the large ϵ_wwi component, so that the overall power distribution decreases only gradually with h_ab beyond the bottom 500 m.

Hence, only about 30% of internal tide‐induced energy dissipation occurs in the immediate vicinity of topography (Figure b), and only about 300 GW contributes to mixing below 1 km depth (Figure a). This budget contrasts with the notion put forth by Munk and Wunsch () that internal tides supply about 1 TW of power to deep‐ocean mixing, primarily along the bottom boundary. In reality, the bulk of the energy of internal tides contributes to mixing the upper ocean, away from topography. This finding highlights the efficient energy redistribution achieved by internal tides and the concentration of vertical gradients in the upper ocean. Notwithstanding, the estimated power distribution does not necessarily imply a sluggish abyssal overturning nor “missing mixing” in the deep ocean: sizeable overturning transports can be maintained below 2.5 km depth with moderate power input to mixing, owing to the weak stratification and large seafloor areas at these depths (de Lavergne et al., , de Lavergne et al., ).

We now compare the predicted distribution of internal tide energy loss to that implied by standard parameterizations of tidal mixing. In the Nucleus for European Modeling of the Ocean (NEMO) model (version 3; Madec & the NEMO team, ), one‐third of the power input to internal tides mapped by Nycander () is distributed in the local water column, using an exponential decay with h_ab and an e‐folding scale of 500 m. A constant diffusivity of 10⁻⁵ m² s⁻¹, reduced to 10⁻⁶ m² s⁻¹ near the equator, is added to represent so‐called background mixing (Mignot et al., ). The Community Climate System Model (CCSM) (version 4; Danabasoglu et al., ) follows an analogous approach but uses a distinct latitude dependence of the background diffusivity (Jochum, ) and the formulation of Jayne and St Laurent () in place of the static map of Nycander (). For a meaningful comparison, the NEMO and CCSM parameterizations are applied here to the WOCE climatology (Gouretski & Koltermann, ). Background diapycnal diffusivities K_ρ are converted into turbulence production rates ϵ using R_f=1/6 and $ϵ=R_f^{-1}{K}_{\rho }{N}^2$ (Osborn, ).

The three parameterizations of tidal mixing produce depth distributions of turbulence production that differ in several ways (Figure ). The CCSM parameterization implies larger power input to mixing in the upper ocean, mostly due to relatively high background diffusivities (averaging 1.68×10⁻⁵ m² s⁻¹ globally compared to 0.83×10⁻⁵ m² s⁻¹ in the NEMO standard). The present parameterization produces the weakest power input at middepths (between 0.5 and 2.5 km) but the largest in the abyss (below 2.5 km). Three main factors underpin this shift of energy toward the abyssal ocean: (i) the present parameterization accounts for variations in the modal distribution of internal tide generation, hence for the stronger local dissipation at abyssal sites relative to middepth sites (de Lavergne et al., ; Vic et al., ); (ii) low‐mode internal tides have been tracked from sources (mostly at steep middepth topography) to sinks (mostly in shallow layers and above abyssal hills); and (iii) the present scheme, unlike the other two, includes internal tide generation by abyssal hills. These effects focus the transition from $\mathcal{O}$(10⁻⁵ m² s⁻¹) to $\mathcal{O}(10^{-4}$ m² s⁻¹) zonal mean diffusivities near 2.5 km depth (Figures a and b).

Diapycnal diffusivity in the present (a,c,e) and NEMO standard (b,d,f) mixing parameterizations. (a,b) Zonal mean diffusivity. (c,d) Diffusivity at 500 m depth. (e,f) Diffusivity at 3,000 m depth. In Panels (a,b) diffusivity averages have been weighted by N2. Note that the standard NEMO parameterization also includes a reduction of diffusivity under sea ice and an elevation in Indonesian seas that are not represented here.

More dramatic differences exist in the horizontal distribution of turbulence production and mixing (Figures c–f). The pronounced lateral heterogeneity of mixing rates mapped here contrasts with the more uniform diffusivities in the NEMO (or CCSM) standard parameterization. Mixing rates are strongly shaped by sources and sinks of internal tide energy, in the abyssal ocean (Figure e) as well as in the pycnocline (Figure c). The resultant complex patterns of mixing are not reproduced when remote tidal mixing is represented by a background diffusivity that varies only with latitude (Figures d and f). Reduced diffusivities in the equatorial band (Gregg et al., ) are not found here; rather, a zonal gradient of diffusivity, between the eastern and western Pacific, is predicted (Figure c) and corroborated by finestructure estimates of Kunze () (Figure b). The more realistic horizontal distribution of mixing produced by the present scheme may have important consequences for the simulated ocean and climate states (Zhu & Zhang, ), including biogeochemical cycles (Tuerena et al., ).

5 Conclusions

Building upon a recent mapping of depth‐integrated internal tide dissipation rates (de Lavergne et al., ), we have proposed a comprehensive and energy‐constrained parameterization of mixing powered by internal tides. This parameterization uses four static two‐dimensional maps of available power (Figure ), each associated with a specific dissipative pathway and a relevant vertical structure of turbulence production (equations 2, 3, 5, and 7). The scheme explicitly accounts for the near‐field dissipation of small‐scale internal tides and the far‐field dissipation of larger scale internal tides, without assuming a constant proportion of near‐field dissipation. Vertical structures incorporate three parameters (H_cri, H_bot, and r_bot) which have been calibrated and mapped (Figure ) with the aid of previous observational and theoretical studies (Legg & Adcroft, ; Lefauve et al., ; Polzin, ; Polzin et al., ).

The proposed parameterization has been applied to an observational climatology of stratification to obtain a global three‐dimensional map of turbulence production. Comparison of this map to a compilation of observational mixing estimates shows that the parameterization has skill in reproducing horizontal and vertical patterns of mixing. The comparison also suggests that, in the ocean interior, internal tides are the principal energy source for mixing and are responsible for the main large‐scale patterns of mixing. This inference is consistent with the large power input to internal tides (∼1 TW) and with the relative temporal stability of basin‐scale dissipation patterns (Ferron et al., ; Kunze, ). The estimated climatology of internal tide energy loss further shows that 70% of the total power lies in the upper kilometer of the ocean. Hence, internal tides contribute first and foremost to turbulence production in the upper ocean—in spite of their generation at the bottom boundary.

The parameterization can also be applied to model oceans. Its implementation in OGCMs merely necessitates a remapping of static maps (Figures and ) onto the model grid and distribution of turbulence production in the vertical according to the simulated stratification. While the available power within each water column is fixed in time, the vertical structure of the energy supply thus evolves with the simulated density profile. Experiments performed with the NEMO model (which will be documented elsewhere) showed that the parameterization obviates the need for a constant background diffusivity above molecular levels, and therefore ensures that explicit diapycnal mixing in the model is energy constrained. This essential property, together with the scheme's low computational cost, motivated its implementation in several climate models participating to Phase 6 of the Coupled Model Intercomparison Project (cf. Voldoire et al., ).

Naturally, the use of static maps representative of modern mean ocean conditions disallows representation of transient changes in the horizontal distribution of internal tide energy loss and limits applicability to other climate states. Nonetheless, we expect only weak sensitivity of this horizontal distribution to climate‐driven changes in ocean stratification (Egbert et al., ). This expectation is backed by an experiment in which we perturbed the climatological stratification entering the two‐dimensional mapping procedure of de Lavergne et al. (). Specifically, we calculated the change in buoyancy frequency between periods 1860–1910 and 2100–2150 of a HadGEM2‐ES simulation (Collins et al., ) forced by historical and RCP8.5 boundary conditions and multiplied the WOCE climatological buoyancy frequency by this change (expressed as the ratio of the later to the earlier period, in h_wkb,1 vertical coordinate). Despite the high climate sensitivity of HadGEM2‐ES (Andrews et al., ), this perturbation led to marginal changes in the geography of internal tide energy sinks (Table ). By contrast, the total energy consumption implied by a constant diapycnal diffusivity K_ρ, calculated as $\int \int \int R_f^{-1}{K}_{\rho }{N}^{2}dM$ with dM the unit mass, increases by 40% between the two periods. Hence, imposing a constant diffusivity can lead to large, spurious, and uncontrolled changes in the magnitude and distribution of turbulence production. This prejudicial behavior and its potential consequences for model drift and simulated climate (Eden et al., ) are avoided by the proposed mixing parameterization.

Process Contributions to M₂ Low‐Mode Internal Tide Dissipation in REF and WARM Experiments

Constituent	M2	M2
Experiment	REF	WARM
Power in Modes 1–5 (GW)	564	564
Shelves (%)	4.8	4.9
Wave‐wave interactions (%)	61.5	60.2
Scattering by hills (%)	9.9	10.7
Critical slopes (%)	14.9	14.9
Shoaling (%)	8.9	9.2

Note. The REF experiment is the reference calculation documented in de Lavergne et al. (). The WARM experiment is identical to REF except for perturbed stratification as described in section . Shelves, defined as areas where the seafloor depth <400 m, have been separated out in the decomposition by process. Only Modes 1 to 5 of the M₂ tidal constituent are included in the total power and percentages given. The impact of the stratification change on the criticality/supercriticality of topographic slopes has been taken into account.

The parameterization shares several features with the IDEMIX framework of Eden and Olbers (). Both approaches are energetically consistent and explicitly include far‐field tidal mixing. But whereas IDEMIX allows low‐mode energy propagation to evolve with the stratification simulated by an OGCM, the present approach relies on static maps of low‐mode energy loss. The choice to map low‐mode dissipation offline is motivated by the aforementioned weak sensitivity of its global geography, by the saved computational cost of online tracking of low‐mode energy, and above all by the greater fidelity of the offline mapping (de Lavergne et al., ). In particular, offline mapping facilitated (i) the use of a Lagrangian (rather than a more diffusive Eulerian) scheme to track low‐mode energy beams, (ii) differential treatment of each of the first five vertical modes, and (iii) a detailed treatment of interactions with topographic slopes. Together with the specific attention given to vertical structures, this results in better congruence with turbulence observations (see de Lavergne et al., ; Pollmann et al., ; and Figure ). Hence, we favored the realism of the geography of internal tide dissipation over its interactivity. This strategy is well founded for mixing due to tides, a relatively stable external forcing, yet it is not for most other sources of mixing. For example, mixing due to lee waves generated by geostrophic flows or to near‐inertial waves energized by atmospheric storms demands fully interactive parameterizations (Jochum et al., ; Melet et al., ).

In spite of its advantages, the proposed tidal mixing parameterization leaves substantial room for improvement. First, its realism is limited by approximations and simplifications of the two‐dimensional mapping procedure and notably by the ad hoc representation of internal tide attenuation by wave‐wave interactions (de Lavergne et al., ). Interactions with balanced flows and with low‐mode near‐inertial waves have not been modeled despite evidence for impacts on propagation and dissipation (Cuypers et al., ; Ponte & Klein, ; Rainville & Pinkel, ). More accurate and comprehensive internal tide generation estimates, accounting for all bathymetric scales and slopes, would also improve the fidelity of the parameterization. Subannual variability, including seasonal and spring‐neap variations of internal tide generation and dissipation, has been ignored in the construction of static maps: its impact deserves further investigation. Second, vertical structures also incorporate important simplifications and uncertainties. For example, the calibration of r_bot and H_bot parameters did not account for the latitudinal sensitivity of triadic wave instabilities (Nikurashin & Legg, ; Richet et al., ) nor for the effects of internal lee waves (Hibiya et al., ). Further work is needed to better constrain these two parameters in order to alleviate biases in the predicted mixing above ridge flanks. Refinement of the vertical structure (2), to account for various scenarios of local or remote dissipation via wave‐wave interactions, should also be pursued. Third, the present parameterization only represents mixing powered by propagating internal tides. Other tidal contributions to mixing, namely, bottom‐trapped internal tides (e.g., Falahat & Nycander, ; Müller, ) and frictional drag of barotropic tides (e.g., Lee et al., ), have not been considered. These contributions could be added provided suitable maps of available power can be obtained.

The degree of agreement between the parameterized distribution of internal tide‐driven mixing and observational estimates of mixing (Figures 8–12) indicates that the parameterization can serve a range of purposes, including forward and inverse modeling, water‐mass transformation estimates, regional to global tracer budgets, and context for field campaigns. Conversely, new field campaigns and numerical studies will help to narrow down uncertainties, expose biases, and identify avenues for improvement of the modeled distribution of turbulence production.

Acknowledgments

We thank all investigators involved in the collection of microstructure measurements analyzed here. Efforts to gather and publish historical microstructure datasets (see https://microstructure.ucsd.edu) are also gratefully acknowledged. We thank E. Kunze, I. Fer, L. Clément, A. Melet, and J. Goff for sharing their published datasets. J. Nycander, L. Clément, E. Kunze, and two anonymous reviewers provided helpful comments on the manuscript. Static maps entering the present parameterization and three‐dimensional fields computed using the WOCE hydrographic climatology are made available at https://doi.org/10.17882/73082. This project has received funding from the European Union's Horizon 2020 research and innovation program under Grant Agreements N° 821001 and N° 824084. A.F.W. acknowledges funding from NSF OCE‐0968721.