Multi-Sensor Observations Reveal Large-Amplitude Nonlinear Internal Waves in the Kara Gates, Arctic Ocean

1. Introduction

Nonlinear internal waves (NLIWs) are known to have a significant impact on the transport of energy, momentum, sediments, plankton and contaminants in the deep and coastal oceans [1,2,3,4,5]. Recent studies have clearly demonstrated the intensification of vertical mixing and anomalously high levels of turbulent energy dissipation and vertical heat fluxes due to the generation, propagation, and subsequent breaking of short-period NLIWs on the Arctic shelf and continental slope [6,7,8,9]. In recent years, more information has been appearing about new areas of large-amplitude NLIW generation in the Arctic [10,11,12,13,14].

The Kara Gates (KG) is apparently one of such regions [15,16,17]. The strait connects the Barents Sea with the Kara Sea. The region is known for its intense dynamics [18]. The transport of sea water in the strait is about 0.3-0.6 Sv [19]. The average flow velocity is in the range of 0.06-0.26 m/s, and the maximum velocity can be as high as 0.5 m/s [20]. Moreover, the KG is among areas with the most intense marine traffic in the Arctic Ocean. The Northern Sea Route, which is the main transpolar transport route connecting the Pacific and Atlantic Oceans, passes through the KG [21,22]. The ongoing active industrial development of the Russian Arctic, especially in the Kara Sea provides ~2500 ship tracks passing through KG annually [23]. Therefore, the increasing vessel traffic in the KG and oil drilling activities nearby make it essential to study local oceanographic processes that might eventually influence them.

The KG area is known for a very complex bottom topography having a complex sill with many crests unevenly distributed along and across the strait [18]. Interaction of strong tidal currents, dominated by M₂ tide, with undulating topography leads to generation of intense internal tides here [15,18,24,25]. Numerous studies have also reported surface signatures of short-period NLIWs regularly registered here in satellite synthetic aperture radar (SAR) images [16,17,18,26], confirming that strong conversion of barotropic tidal energy into IWs and high levels of tidal energy dissipation occur in the region [27].

Since 1997, the study of IWs in the KG has been carried out using moorings and towed conductivity-temperature-depth (CTD) profilers recording in scanning mode [15,18,20], which made it possible to determine main characteristics of IWs and their relationship with background currents. The results of numerical modeling showed the existence of internal tides with amplitudes of 10-16 m in the strait at the Barents Sea side [25] and formation of a strong hydraulic jump with height >50 m at the Kara Sea side [18]. In these numerical experiments, the generation of short-period NLIW trains was reproduced only qualitatively. A more detailed analysis of generation and evolution of NLIWs in the strait was performed later using a two-dimensional non-hydrostatic MITgcm model that revealed the main kinematic characteristics of NLIWs and some peculiarities of their propagation at both sides of the strait [28].

Despite existing reliable information about the regular formation of NLIW trains in the KG during the interaction of tidal currents with rough bottom topography, so far, no special experiments to study their vertical structure, spatial and kinematic properties have been carried out. It is noteworthy that the strait is characterized by an extremely heterogeneous topography with many underwater ridges and channels [18], which leads to a very complex pattern of the generation and propagation of NLIW packets both along and across the strait [16].

The goal of this paper is to present new results of multi-sensor remote sensing and in situ measurements of NLIWs in the KG in August 2021 during the 58^th cruise of the R/V Akademik Ioffe, documenting, for the first time, the formation and propagation of large-amplitude NLIWs in this Arctic region.

2. Materials and Methods

The KG is a 30 km long and 56 km wide strait between the southern end of the Novaya Zemlya and the northern tip of Vaygach Island connecting the Barents Sea in the west and the Kara Sea in the east (Figure 1). The bottom topography of the strait is characterized by a complex sill with several crests and troughs oriented both along (Figure 1b) and across the strait (see Figures 4, 6 in [18]).

Measurements were collected in the Kara Gates between 07:00 UTC on 12 August and 02:00 UTC on 13 August 2021. Field works were performed at latitudes 70°-71° N, which is close to the critical latitude for the lunar semidiurnal tide M₂ [29,30,31]. Thus, internal tidal waves of main energetic peak at frequency M₂ are forced and must dissipate locally, evolving into packets of short-period NLIWs, which can freely propagate and further transport tidal energy at high latitudes [7,32].

During the expedition, measurements included those by distributed temperature sensors (thermal profilers), frequent CTD casts, observations from an unmanned aerial vehicle (UAV), and quasi-operational analysis of satellite data. Unfortunately, the ship-mounted acoustic Doppler current profiler system was out of service, hence, a standard information on ocean current properties was not available.

Satellite information included optical, infrared and SAR data of medium and high resolution (Sentinel-1, -2, -3, MODIS Aqua/Terra, and Landsat-8). This information was downloaded and processed in near-real time by the Remote Sensing Department of MHI RAS and transmitted to the vessel for planning the measurement sites and assessing the background NLIW field in the strait.

On 12 August 2021, five stations were occupied in the strait, the position of which is shown in Figure 1a with red (start) and green (end) dots and summarized in Table 1. To determine the oceanographic background and estimate the depth of the thermocline, the vertical CTD profiles were performed at the beginning of each station using a Sea-Bird Electronics 19 plus instrument. From these data, the buoyancy frequency was computed as $N_z=\sqrt{-\left(g/{\rho }_0\right)\partial \rho /\partial z}$, where $g=$ 9.8 m/s² is the gravity acceleration, ${\rho }_0$ is the mean water density, and $\rho$ is the vertical density profile. Rapid CTD casts were performed using an autonomous AML BaseX profiler.

In order to record surface manifestations of NLIWs and correct the position of the vessel, photo and video recording was made using a DJI Mavic 2 Pro UAV from heights of 100 to 500 m at all stations where visibility conditions were favorable. We also used the records of marine navigational X-band radar JMA-9122-6XA at one station. The radar was primarily operated to study the properties of wind waves, but it was also useful to detect surface signatures of NLIWs [33]. The radar was located at an altitude of 25 m above sea level and carried out measurements with a rotation frequency of 24 rpm emitting a backscatter signal with a frequency of 9.41 GHz at a wavelength of 3.18 cm with a pulse length of 0.07 μs. The SeaVision system was used to record the radar data.

The key instrument for measuring NLIW parameters was a thermal profiler (TP) TPArctic developed at the MHI RAS with a measuring cable length of 48 m (32 continuous measuring sections each 1.5 m long), a 100 m power/communication cable, and two pressure sensors in the lower and upper sections of the line [34]. It has a digital temperature resolution of 0.0007 °C, a precision in measuring temperature over 1.5 m intervals of 0.1 °C and a sampling rate of 0.5 Hz. The main advantage of this device is its direct connection to a portable computer that provides controlled visualization of vertical profile of seawater temperature enabling to capture NLIW propagation in real-time telemetric mode. In addition, the length of the power/communication cable provides opportunity to vary depth of the instrument by means of casting, depending on the depth of the thermocline. Depth reference in this case was provided using the readings of pressure sensors.

To understand propagation direction of NLIWs, measurements with thermal sensors are often made from three spatially spaced positions in the form of a triangle (e.g., [12,35,36]). At one station (#3915), two more discrete TPs with 10 and 11 Starmon Mini StarOddi autonomous temperature sensors (with an accuracy of ±0.025 °C and sampling rate of 1 s) located at a distance of 2 m from each other in the vertical with pressure sensors DST centi TD StarOddi at the upper and lower depths of the TPs were used. Pressure sensors were used for linear correction of the TPs’ deviation from the vertical position. The total length of these TPs was more than 30 m, the measuring parts were about 20 and 22 m, respectively. Similar sensors were used to measure NLIW properties in other seas [37,38]. Onboard the ship, the TPArctic was located on the port (left) side in the center of the ship, while the discrete TPs were located on the starboard (right) side of the ship (at vessel stern and vessel bow) at a distance of 51 m between them, and from 20 m to 40 m from the TPArctic. Measurements with the TPs were carried out during a free drift of the vessel.

Information about tidal dynamics was derived from the Arctic Ocean Inverse Tide Model on a 5-kilometer grid (hereinafter, Arc5km2018) [39]. The information about the total water depth during the measurements was obtained from ELAC LAZ4700 echosounder. Since not all the bathymetry data were digitally recorded, we also used records from the navigational software OpenCPN. Wind measurements at 22 m height were made using AIRMAR WeatherStation 220 WX with temporal resolution of 1 s.

3. Results

3.1. Background Conditions in the Kara Gates

During the field works in the KG wind conditions were characterized by low to moderate winds with mean velocity of about 5 m/s during the first three stations, rising gradually to about 7 m/s toward the end of observations (Figure 2).

In the KG, the main current is directed from the Barents Sea towards the Kara Sea. Usually, it is located along the northwestern coast of Vaygach Island and may cover about 2/3 of the strait width [40], as seen in a high-resolution (~100 m) sea surface temperature (SST) map of 14 September 2016 (Figure 3a) produced from Landsat-8 data using a two-channel SST algorithm [41]. Landsat-8 perfectly resolves a surface circulation in the strait and clearly depicts a relatively warm northeastward current that typically attains 0.1-0.2 m/s peaking at 0.5 m/s [19]. MODIS Aqua SST map acquired on 9 August 2021, i.e., three days prior to measurements, shows a very similar surface pattern at 1 km resolution with the relatively warm (9-11° C) Barents Sea water covering the southeastern half of the strait (Figure 3b).

The closest to the measurement period and almost cloud-free day on 9 August 2021 provided large set of satellite data, including a number of Sentinel-3 A/B Sea and Land Surface Temperature Radiometer (SLSTR) images. A series of consequent SLSTR Level-1B thermal IR brightness temperature images at 1 km resolution was used to reconstruct the surface velocity field in the study region by applying a 4-Dvar algorithm [42]. The time gap between sequential images acquired at 16:05 and 17:46 GMT was 1 h 41 min. Figure 3c shows the resulting surface velocity vectors plotted over the absolute value of surface current velocity (shown in color). In the strait, velocity vectors were predominantly directed north and northeast. The maximal current speed of 0.8-0.9 m/s was observed in the central and northeastern parts of the strait. In the central and eastern KG, the current direction gradually turned eastward in agreement with the surface thermal expressions demonstrated in Figure 3a,b. The measured maximal surface current velocities were almost twice higher than those reported earlier [19,20], however, the analysis of the ship drift velocity presented below confirmed them.

Tidal Conditions

When studying internal waves, it is important to know the phases of the moon. In our case, the new moon was on 10 August 2021, and 12 August 2021 was the period of the growing moon. Thus, the measurements were carried out during the spring tide, which is characterized by the maximal tidal current velocities along the strait [18]. This is well illustrated in Figure 4a showing tidal variations of sea level during a two-week period on 5-19 August 2021 from Arc5km2018 pointing out the maximum amplitudes of tidal range on 11-13 August 2021.

We used Arc5km2018 model to predict the tidal currents over the period from 12 to 13 August 2021 (48 hours) using four semidiurnal constituents M₂, S₂, N₂, K₂ and four diurnal constituents K₁, O₁, P₁, Q₁. Figure 4b shows the maximum tidal current speed (in color) and the M₂ ellipses in the KG. Three main areas of amplification of tidal currents in the strait can be distinguished: (a) the northwestern part of the strait near the Kusova Land and the Loginov Islands (see Figure 1a for geographic notations); (b) the eastern part of the strait along the Vaygach coast; (c) the southern part of the strait. In the areas (a) and (b) the maximal tidal current velocities expected from the Arc5km2018 model are 0.5‒0.6 m/s. In area (c), the speed of the tidal currents is 0.35‒0.45 m/s. The semidiurnal constituent M₂ makes the main contribution to the total tidal current. The M₂ tidal ellipses are relatively narrow and are aligned along the strait axis, averaging at 0.1-0.15 m/s. Only along the Vaygach coast (areas (b) and (c) the tidal ellipses of the M₂ tide are broader and the current speed is up to 0.3-0.4 m/s.

A very complex bottom topography in the KG and pronounced tidal currents superposed on the mean currents from the Barents to the Kara Seas may cause very high horizontal velocities in the strait of the order of 1 m/s, creating favorable forcing conditions for the generation of large-amplitude NLIWs here.

Surface Signatures of NLIWs in Satellite Data

Satellite data covering the study region were operationally communicated to the expedition team for planning the locations of stations. Figure 5a shows an almost cloud-free Sentinel-2A true color image acquired over the study site on 8 August 2021, i.e., four days prior to field observations. As in situ measurements by thermal profilers were done in free drift, Figure 5 shows start and end points of stations by red and green circles, respectively.

Figure 5b shows an enlarged fragment of Sentinel-2A image with well-seen surface signatures of NLIW trains, most of which in this case were directed south and southwest. However, satellite images of other dates also showed systems of NLIWs travelling northeastward. The locations of stations were selected so that they would begin south of underwater hills and end up after crossing them while the ship freely drifted northeast.

On 11-13 August 2021, the KG region was densely covered by clouds that hindered usage of satellite optical data. In turn, Sentinel-1A/B made two SAR snapshots of the study region on 12 August at 02:43 GMT, i.e., about four hours prior to the ship arrival to the first station, and on 13 August at 02:36 UTC, i.e., just after the completion of the last station (Figure 6). However, due to low winds at the time of image acquisitions, they showed only dark elongated slicks mimicking the patterns of surface circulation (Figure 3c).

3.2. Vertical Thermohaline Measurements at Stations

As distinct thermocline oscillations presumably associated with the NLIW activity were observed only at four stations 3909-3915, below we describe the data only from these stations.

Station #3909

The first station (station #3909) was located in the western part of the KG southwest of one of the sill crests peaking at about 40 m depth (Figure 1b), behind which the formation of NLIW packets was regularly observed earlier [18]. An average water depth during the station was about 170 m (Figure 7).

At this station measurements were carried out only by the TPArctic from 08:25 to 09:28 UTC on 12 August 2021. The beginning of the station coincided with the low water and the beginning of flood (Figure 4a). According to Arc5km2018, an NNE tidal current of 0.08 m/s was present in the beginning of station, gradually rotating clockwise and decreasing to 0.06 m/s (Figure 7).

The vertical CTD cast made prior to thermal profiling revealed stable stratification and thermocline located at the depths of 25-40 m (Figure 7a). The first peak of $N_z$ up to 13 cph was observed at a depth of 18 m. However, the major peak was observed at a depth of about 30-35 m. Pronounced fluctuations of isotherms with heights from a few meters to 15 m and periods from 2 min to 8 min were observed (Figure 7b). Most of the record contained relatively weak oscillations not exceeding 5-6 m. However, starting from the 44^th minute, i.e., closer to the sill crest located northeast, two pronounced and steep displacements of 4°C and 6°C isotherms were seen at 25–45 m depths. The first wave was rather steep, 13–15 m high with a period of 7 min. The second one was about 10 m high and had a front slope much steeper compared to the rear one.

Station #3911

The next station 3911 was planned for the vessel to drift over another sill, the shallowest along the ship route, rising upward from 150-170 m depth up to 35 m depth (Figure 1b). Field works began at 10:59 UTC and lasted for about 2 hours. The measurements were carried out during the flood tide characterized by weak currents (Figure 4a).

The vertical buoyancy profile shows two peaks at 10 m and 33 m depths (Figure 8a). The most intense fluctuations of isotherms were observed between 10 m and 40 m depths (Figure 8b). In the beginning, 5-11 m high oscillations of isotherms at a depth of 25–40 m were seen at a frequency of about 11-12 cph (Figure 8b). After 20 minutes of the free drift northward, the ship passed the sill crest. The profiler first recorded a sharp vertical rise of isotherms up to 8-15 m depth followed by a 30 m high hydraulic jump down to 43 m of the 4°C isotherm with a period of 34 minutes (Figure 8b). The next fluctuation of the 4°C isotherm at the 70^th minute had almost the same period (about 30 min), but a lower height (about 20 m). Another pronounced and steep oscillation of isotherms with a height of about 20 m and a period of about 10 min was also recorded at minute 105 (Figure 8d).

At this station UAV survey was made prior to TP measurements to check whether surface signatures of NLIWs were present nearby. As a result, very characteristic sea surface manifestations of NLIWs in the form of elongated narrow slicks, usually observed under low winds [43], were detected several hundred meters northward from the ship (Figure 9). The curvature of the slick bands suggested that the waves were oriented northeastward.

During this station, the marine radar was in active mode helping to detect surface signatures of NLIWs present in thermal profiling data. The first 30-m high oscillation (centered at 37^th minute in the TP record) has a clear radar manifestation as a curved, about 0.2 km wide and 2-3 km long band of alternating suppressed (from the SW side) and enhanced (from the NE side) radar echo (Figure 9a,b). This feature is similar to standard IW manifestations usually seen in spaceborne SAR data [9,18,44], suggesting that the NLIW propagation direction is toward the northeast.

At 11:16 UTC the IW signature was about 800 m northeast of the ship (Figure 9a), while between 11:31 UTC and 11:35 UTC the ship passed it away (Figure 9b). The mean ship drift velocity was about 0.8-0.9 m/s, suggesting that the observed hydraulic jump and its surface manifestation were rather stationary during this period.

Another pronounced radar manifestation was depicted at 12:37 UTC (Figure 9c). It is characterized by a smaller horizontal extent. The propagation direction of this feature is toward NE (azimuth ~40°). The passage of this feature through the vessel location corresponds well to the deepening of the isotherms centered at 108^th minute of the TP record (Figure 8b). The average propagation velocity of this structure relative to the vessel was ~0.55 m/s. The average ship drift velocity during the interval of 90-105 minutes of the TP record was ~0.6 m/s. From the latter, it follows that the surface manifestation spread in the NE direction at a very low speed of ~0.05 m/s.

Station #3913

Station 3913 was located in the very center of the strait over a smaller sill crest peaking at 60 m surrounded by water depths of 100-120 m. The TPArctic measurements began at 16:18 UTC and lasted for about 80 minutes during the ebb tide. According to Arc5km2018, the northward tidal current reached 0.12 m/s at the beginning of station and raised to 0.19 m/s at the end of measurements (Figure 11).

The vertical profile of buoyancy frequency shows a clear peak in the thermocline layer at 26 m depth, amounting to 20 cph (Figure 11a). At the beginning of the record one can see several oscillations of 10-15 cph with a height of about 5–8 m (Figure 11b). More pronounced fluctuations of isotherms were seen later, after passing the sill crest on the 55^th minute of the record. At this moment, rapid CTD casts down to 60 m depth with a period of about 2-4 minutes were made simultaneously with the thermal profiling. The recorded movement of isopycnals is shown in Figure 11b in green. The oscillations of isopycnals and isotherms coincide in shape and height, but have a phase shift of several minutes since the measurements were made at a distance of 35 m from each other. The maximum heights of oscillations reached 12 m at this station. This result emphasizes that the use of thermal profiling in the strait is quite justified, since the vertical water stratification significantly depends on temperature in this Arctic region.

Station #3915

Station 3915 was made at the end of the ebb tide (Figure 4a) over the deeper region (~ 140 m depth) at some distance from any pronounced bottom topography features (Figure 1b). The previous work showed a formation of a strong hydraulic jump over this part of the strait by the northeastward current over a system of sills [18].

The atmospheric conditions were characterized by 4-7 m/s winds of varying direction and strong fog that hindered UAV measurements. Unfortunately, the recording of the shipborne radar was switched off during the station. The thermal profiling began at 20:01 UTC and lasted for about 70 minutes. Measurements at this station gave the most interesting results in terms of the recorded NLIW amplitudes.

At this station the buoyancy frequency was characterized by a main peak of 21 cph at 18 m depth (Figure 12a). Vertical profiles of temperature and salinity showed a well pronounced two-layer structure with a relatively warm (~8° C) and less saline (33.7 PSU) 15-20 m thick upper layer followed by about 5-m thick pycnocline and much thicker cold (1-3° C) and more saline (~34.2 PSU) bottom layer.

The beginning of the record was characterized by rather calm conditions with only a single 7-8 m high oscillation of nearly all isotherms centered at the 7^th minute of the record (Figure 12b). Yet, the most interesting features started after the 25^th minute of the record. A strong and steep, nearly 30 m high depression of isotherms centered at 32^nd minute with a period of 15 minutes was followed by the more intense 36-40 m high depression centered approximately at the 53^rd minute with a period $T$ of 27 minutes. This second wave had a steeper rear face compared to its front. Normalized wave amplitude of both waves, i.e., the ratio of wave amplitude ${\eta }_0$ to the mean upper-layer depth $h_1$ (~20 m), was about 1.5-2. Note also a relatively flat trough of the second wave that is often observed when the wave amplitude (~40 m) occupies a significant part (~30%) of the total water depth [45,46,47]. This is a characteristic feature, also predicted by weakly nonlinear theory, when an increase in soliton amplitude is accompanied by soliton profile change from a bell-shape to a more rectangular one [48].

These two isothermal depressions, associated with the NLIW propagation, were also well seen in the movement of pycnocline recorded by rapid CTD casts (green lines in Figure 12b). The overall shapes and heights of isothermal and isopycnal depressions were similar, and had a characteristic time delay of about 2-3 minutes as those at the previous station. During the NLIW passing, measurements of the three TPs and rapid CTD casts were made simultaneously. This allowed us to estimate the propagation speed and direction of the NLIWs from characteristic time delays of their appearance in the records of different instruments (Figure 13a).

During the measurements, the ship headed southeast while drifting northeast (azimuth of 52°) with a mean speed of 0.94 m/s (Figure 13b,c). In fact, it was overtaking the wave train moving nearly in the same direction. However, as the ship moved faster than internal waves, they appeared as approaching the ship from the NNE and moving SSW at $C_{rel}~$ 0.1 m/s (see direction of brown arrow in Figure 13b). The summation of the ship drift vector and the NLIW propagation vector relative to the ship shows that the true propagation speed and direction of the NLIW train was $C_{obs}=$ 0.86 m/s to the northeast (azimuth of 56°). Based on the information about the mean period of the NLIWs, $T=\left(15+27\right)/2=21$ minute, and the NLIW speed relative to the vessel, their mean wavelength $\lambda =T{C}_{rel}=21·60·0.1\approx 130$ m, i.e., nearly equal to the water depth.

The phase speed of NLIWs can be also estimated theoretically based on available oceanographic information. The vertical stratification profile taken at this station (Figure 12a) clearly exhibits a two-layer structure allowing to apply an interfacial two-layer model for the assessment of the NLIW phase speed [9,46,49,50].

In our case, we deal with large-amplitude NLIWs for which a fully nonlinear theory should be applied [47,51]. However, several authors showed that a weakly nonlinear Korteweg–de Vries (KdV) model or its extended version (eKdV) are also valid for the moderate and large-amplitude NLIWs [52,53]. The eKdV model, including a cubic nonlinearity term, is even more applicable for such waves, provided the difference between the upper and lower layer is not large and the wave amplitude doesn’t exceed an upper bound close to half-difference between thicknesses of the layers [53,54]. As will be shown below, the latter is about 50 m and the observed 40 m high NLIWs don’t exceed it. Again, the observed broadening of the wave crest with an increase of wave amplitude (Figure 12b) is something that is well caught and predicted by the eKdV [47,55]. Owing to simple use and their wide applicability [46,47,56,57,58,59,60], below we use these models to assess the NLIW phase speed from observations.

Following [61], at first, we have to choose a suitable analytical solitary solution that depends on specific scale parameters, namely, $h_1$ - the thicknesses of the upper layer, $D$ - the full depth, $L$ - the half-width of the soliton, ${\eta }_0$ - the amplitude of NLIWs, and their ratios $L/h_1$, $h_1/D$, $L/D$, ${\eta }_{0}L/h_1^2$. The half-width of the soliton is defined as $L={\left(12\beta /\alpha {\eta }_0\right)}^{1/2}$, where $\alpha$ and $\beta$ are the coefficients of quadratic nonlinearity and dispersion of the KdV equation [62] defined as:

$$\alpha =\frac{3}{2}{c}_0\left(h_1-h_2\right)/h_1{h}_2$$

(1)

$$\beta =\frac{c_0}{6}{h}_1{h}_2$$

(2)

where $c_0$ is the linear IW phase speed defined below using either Eq. (3) or Eq. (5), $h_2$ is the thickness of the bottom layer.

In our case, $h_1=20$ m, $h_2=120$ m, $D=$ 140 m, $L~$50 m, ${\eta }_0=$30-40 m, $c_0=$ 0.42 m/s, yielding $L/h_1=$ 2.5, $h_1/D=$ 0.14, $L/D=$ 0.4, ${\eta }_{0}L/h_1^2=$ 4-5. This result lies somewhere in between the finite depth [$L/h_1\gg 1$, $h_1/D\ll 1$, $L/D\approx$ O(1), ${\eta }_{0}L/h_1^2\approx$O(1)] [63] and the shallow water [$L/D\gg 1$, $h_1/D\approx$ O(1), ${\eta }_{0}L/h_1^2\approx$ O(1)] [64] theories. We, therefore, will apply both theories to estimate the nonlinear phase speed $C$.

In two-layer model, the phase speed of linear internal waves of the lowest internal mode can be defined as [65]:

$${с}_0={\left(\frac{\delta \rho }{{\rho }_0}\frac{g}{k}\right)}^{\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$2$}\right.}{\left(\coth k{h}_1+\coth k\left(D-h_1\right)\right)}^{-\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$2$}\right.}$$

(3)

where $\delta \rho$ is the density difference between the waters below and above the pycnocline, ${\rho }_0$ is the mean water density below the pycnocline, $g$ is the gravity acceleration, $k$ is the wave number defined as $k=2\pi /\lambda$, where $\lambda$ is the characteristic wavelength of NLIWs from observations. In case of the finite depth theory, the nonlinear soliton speed can be found as [63]:

$$C_{fd}=c_0\left\{1+\frac{h_1}{2D}\left[1+\frac{D}{b}\left(1-a^2{b}^2\right)\right]\right\}$$

(4)

where $a$ is the wave number-like parameter satisfying the relationship $a\tan aD=b^{-1}$, and $b=4h_1^2/3{\eta }_0$.

In case of shallow water, the equation for the linear phase speed of the KdV model will take form [66]:

$${с}_0={\left(g\frac{\delta \rho }{{\rho }_0}\frac{h_1{h}_2}{h_1+h_2}\right)}^{\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$2$}\right.}$$

(5)

while that for the nonlinear soliton speed would be equal to:

$$C_{kdv}=c_0+\frac{{\alpha \eta }_0}{3}$$

(6)

with an additional coefficient of cubic nonlinearity, $\kappa$, in the eKdV model:

$$C_{ekdv}=c_0+\frac{{\eta }_0}{3}\left(\alpha +\frac{1}{2}\kappa {\eta }_0\right)$$

(7)

with $\kappa$ defined as:

$$\kappa =\frac{{3c}_0}{{\left(h_1{h}_2\right)}^2}\left[\frac{7}{8}{\left(h_1-h_2\right)}^2-\left(\frac{h_1^3+h_2^3}{h_1+h_2}\right)\right]$$

(8)

Lastly, we also use equation to derive phase speed of strongly nonlinear waves [54]:

$$C_{st}=c_0\left\{1+3\frac{\left(h_2-h_1\right)\left(h_2-h_1+2{\eta }_0\right)}{{\left(h_1+h_2\right)}^2}\left[\sqrt{\frac{\left(h_1-{\eta }_0\right)\left(h_2+{\eta }_0\right)}{h_1{h}_2}}-\frac{h_2-h_1+2{\eta }_0}{h_2-h_1}\right]\right\}$$

(9)

with $c_0$ defined from (5).

At station 3915, the characteristic values are $\delta \rho$ = 1.75 kg/m³, ${\rho }_0$ = 1028 kg/m³, and $\lambda$ = 130 m. The theoretical $C$ values are presented in Table 2. Each of them has a certain range of values because of two values of amplitude ${\eta }_0$ (30 and 40 m) used in the calculations. In most cases, theoretical values are somewhat lower than the observed phase speed of 0.87 m/s. The best correspondence is obtained for the KdV model (eqs. 5, 6) and that for the strongly nonlinear waves (eqs. 5, 9). Notably, the inclusion of cubic nonlinearity term ($\kappa =-$0.001 m^-1s^-1) in the eKdV equation (Eq. 7) significantly decreases the nonlinear phase speed values compared to the KdV model (Table 2).

Following [56,61], one may also assess the total internal wave energy per unit crest length, $E_T$:

$$E_T=\frac{4}{3}\delta \rho g{\eta }_0^{2}L$$

(9)

In our case, the total energy equals to 1.0-1.8 MJ/m for these two large-amplitude waves.

In order to examine the possible NLIW generation mechanism we should look more closely to the flow conditions characterized by Froude number, $Fr=\mathit{u}/{с}_0$, controlling the generation and propagation of NLIWs when the stratified flow passes over undulating topography. The Froude number is defined as the ratio of the total current vector $\mathit{u}\left(u,v\right)$ to the phase speed of long linear internal waves ($c_0$). The latter is defined from the vertical profile shown in Figure 12a using Eq. 3 ($c_0=$ 0.42 m/s) and Eq. 5 ($c_0=$ 0.54 m/s).

Due to the absence of in situ current meter measurements we use ship drift velocity to approximately define $\mathit{u}\left(u,v\right)$. In such case, we assume that the actual ship drift is a sum of ocean current velocity and wind-induced ship drift. As a coarse estimate, we presume that the latter doesn’t exceed 2% (using analogy to Nansen–Ekman ice-drift law [67]) of the maximum wind speed of 7-8 m/s observed prior and during the station (Figure 2).

After subtraction of the wind component, the ship drift velocity vector $\mathit{u}$ is converted into orthogonal $u,v$ components with $u$ oriented cross-strait perpendicular to the NLIW direction (azimuth 146°), and $v$ oriented northeast (azimuth 56°) along the strait and the NLIW propagation direction. In such configuration, we presume that $v$ component dominates in creation of the initial pycnocline disturbance during the ebb tide that might further evolve into the train of nonlinear waves.

During station 3915, the NLIWs were recorded between the 25^th and 65^th minutes when the surface current started to slack and changed direction from northeast to east, presumably due to the ebb tide slackening and reversal, well depicted in Figure 14a and Figure 12c.

In the record, $v$ ranges from 0.1 to 1.4 m/s, being at the maximum prior to and after the passing of NLIWs (Figure 14a). During the NLIW passage, $v$ drops down first by ~0.5 m/s (from 1.2 m/s to 0.7 m/s) and later down to 0.3 m/s at the 50^th minute, followed by a rapid rise to 1-1.4 m/s after the waves’ passing.

The Froude number, calculated for two values of $c_0$ (Table 2), is supercritical for the most of the record (Figure 14b). It passes through the critical value of $Fr=1$ several times. During the NLIW passage, it drops down to near-critical (~ 1) and subcritical (< 1) values, similarly as reported by Fer et al. [8] during the passage of 50 m high NLIWs north of Svalbard.

4. Discussion

The obtained results are solidly grounded on previous observational and numerical studies in the Kata Gates [16,17,18,25,27,28]. However, our record shows much larger internal waves than were observed or predicted before. One of such reasons could be that our measurements were made exactly during the peak spring tide when the tidal forcing was maximal.

Figure 15 shows a generalized 2-d scheme of NLIW locations over the KG topography along the ship route. The recorded NLIWs were detected both immediately over the sill crests or in close proximity to them (stations 3909-3913), as well as at some distance from them (station 3915). The wavelengths of the observed NLIWs (~100-200 m) were of the same scale as narrow sill crests passed by the ship during a free drift.

An intense along strait flow composed of the mean current from the Barents Sea superposed with pronounced barotropic tide results in strong baroclinic response and NLIW generation in the Kara Gates. Vlasenko et al. [32] show that in the vicinity of critical latitude for the M₂ tide (74.5° N) two main options of baroclinic response exist depending on flow regime. In subcritical regime ($Fr$<1), only weak, large wavelength and small-amplitude baroclinic tides are generated. When the flow is supercritical ($Fr$>1), topography-scale short NLIWs are generated, exactly as observed in our case.

Moreover, Vlasenko et al. [32] also showed a possibility for the generation of “double waves” that is caused due to the nonlinearity and the wave energy transfer from the tidal frequency to the higher harmonics. Notably, the shape of high-amplitude NLIWs at station 3915 (Figure 12) is very similar to that predicted by numerical model (see Figure 9 in [32]).

It is also interesting to compare our results with those obtained from nonhydrostatic model by Li et al. [28]. While these authors used much more simplified topography in the strait, the overall structure of NLIWs’ generation was predicted qualitatively well. In particular, they showed that the waves generated on the Kara Sea side were weaker than those on the Barents Sea side, and disintegration of the initial pycnocline disturbances occurred at about 25-30 km from the sill. The model phase speed of the waves travelling to the Kara Sea was 0.72 m/s, somewhat less than in our observations.

However, the major difference between our observations and the model is found for the case when the latter accounted for a steady background current toward the Kara Sea, which is observed in reality. In such case, the model predicted a generation of only a weak internal bore on the Kara Sea side with a small vertical displacement and not disintegrating into the wave train, while our measurements show at least two pronounced large-amplitude waves propagating at about 0.9 m/s. We presume that this discrepancy is most probably related to the lack of the detailed topography in the model simulations. Once included, it might be practical then to clarify what particular sill crest is responsible for the generation of large-amplitude NLIWs at station 3915.

Owing to radial spreading, the intense NLIWs originating from the KG are known to propagate farther offshore to the Barents and Kara Seas [16,17,18] where they could be still enough energetic to be a concern for drilling operations at many offshore oil production sites located nearby. Moreover, one can also expect a possible impact of large-amplitude NLIWs on sediment resuspension both at the strait margins and farther offshore [4].

Our observations evidence that the NLIWs generated in the KG are among the largest ever documented in the Arctic, including those with heights of 36-50 m registered over the northern flanks of the Yermak Plateau [68,69], those of ~50 m registered north of Svalbard [8], and those of 40-50 m registered in Franz Victoria Trough [30,36] and southeast of Hopen Island [13] in the Barents Sea.

While the KG is one of major spots of depth-integrated baroclinic tidal energy dissipation (>10 W/m²) in the Barents Sea [27], it is characterized by rather moderate tidal forcing in the Arctic [70,71]. The latter suggests that similar intense large-amplitude NLIWs might be also forming in other Arctic regions with rough topography and identical flow conditions.

5. Conclusions

In this work we report the results of the experimental observations specifically targeted to resolve short-period NLIWs in the Kara Gates. Synergistic use of remote sensing observations, marine radar, three spaced-apart thermal profilers and rapid CTD casts enabled to record anomalously large NLIWs forming in the KG with amplitudes up to 40 m, i.e., about one third of the total water depth in the strait.

Analysis of high- and medium-resolution satellite optical data clearly revealed main patterns of the surface circulation in the strait, depicting a presence of a relatively warm northeastward current from the Barents Sea toward the Kara Sea attaining velocities of 0.8-0.9 m/s. Prediction of tidal currents using Arc5km2018 model showed the maximal tidal current velocities along the ship route reaching 0.20-0.25 m/s. As a result, the combination of local complex bottom topography and intense currents created favorable forcing conditions for the generation of intense NLIWs here.

Quasi-operational observations of sea surface state in the strait from satellites and UAV was an important component of the experimental work that allowed us to plan the location of stations and trace the surface signatures of NLIW in the strait.

As the vertical density distribution in the study region strongly depended on seawater temperature, shipborne measurements of three spaced-apart thermal profilers allowed us to record numerous pronounced thermocline oscillations in the water layer of 10 to 65 m depth. We associated these oscillations with propagation of NLIWs of depression generated when the stratified flow passed a system of shallow sills in the strait.

The amplitudes of the observed internal waves ranged from 4 to 40 m, revealing the existence of very intense large-amplitude NLIWs in this region. The total internal wave energy per unit crest length for the largest waves was estimated to be equal to 1.0-1.8 MJ/m. The waves were recorded both immediately after passing the shallow sills and at some distance from them. The most intense NLIWs were recorded during the ebb tide slackening and reversal when the background flow was predominantly supercritical.

Measurements of the three spaced-apart TPs enabled to assess the phase speed and propagation direction of the large-amplitude NLIWs. The waves travelled northeastward with phase speeds of about 0.9 m/s. Theoretical predictions of NLIW phase speed made within a two-layer approach showed that the shallow water KdV model and that for strongly nonlinear waves agreed best with observations.

Obviously, the current meter data were lacking in our study to precisely relate the observed thermocline oscillations to the forcing conditions. In this sense, longer fixed-point measurements covering 1-2 complete M₂ tidal cycles would be needed in future to clarify the NLIW generation mechanism and periodicity within a tidal cycle. The future goal is also to study the evolution and possible impacts of NLIWs on water hydrology, turbulent mixing and biogeochemical fluxes using synergistic multi-sensor observations.