Enhanced diapycnal mixing with polarity-reversing internal solitary 1 waves revealed by seismic reflection data 2

Survey, Guangzhou 510760, China 9 Corresponding author. hbsong@tongji.edu.cn 10 11 Abstract 12 Shoaling internal solitary waves near the Dongsha Atoll in the South China Sea dissipate their 13 energy and enhance diapycnal mixing, which have an important impact on the oceanic environment 14 and primary productivity. The enhanced diapycnal mixing is patchy and instantaneous. Evaluating 15 its spatiotemporal distribution requires comprehensive observation data. Fortunately, seismic 16 oceanography meets the requirements, thanks to its high spatial resolution and large spatial coverage. 17 In this paper, we studied three internal solitary waves in reversing polarity near the Dongsha Atoll, 18 and calculated their spatial distribution of diapycnal diffusivity. Our results show that the average 19 diffusivities along three survey lines are two orders of magnitude larger than the open-ocean value. 20 The average diffusivity in internal solitary waves with reversing polarity is three times that of the 21 non-polarity-reversal region. The diapycnal diffusivity is higher at the front of one internal solitary 22 wave, and gradually decreases from shallow to deep water in the vertical direction. Our results also 23 indicate that (1) the enhanced diapycnal diffusivity is related to reflection seismic events; (2) 24 convective instability and shear instability may both contribute to the enhanced diapycnal mixing 25 in the polarity-reversing process; and (3) the difference between our results and Richardson-number26 dependent turbulence parameterizations is about 2-3 orders of magnitude, but its vertical distribution 27 is almost the same. 28 29 [Key words] Internal solitary waves, Polarity reversal, Diapycnal mixing, Northeastern South China 30 Sea, Seismic oceanography. 31 32


12
Shoaling internal solitary waves near the Dongsha Atoll in the South China Sea dissipate their 13 energy and enhance diapycnal mixing, which have an important impact on the oceanic environment 14 and primary productivity. The enhanced diapycnal mixing is patchy and instantaneous. Evaluating 15 its spatiotemporal distribution requires comprehensive observation data. Fortunately, seismic 16 oceanography meets the requirements, thanks to its high spatial resolution and large spatial coverage. 17 In this paper, we studied three internal solitary waves in reversing polarity near the Dongsha Atoll, 18 and calculated their spatial distribution of diapycnal diffusivity. Our results show that the average 19 diffusivities along three survey lines are two orders of magnitude larger than the open-ocean value. 20 The average diffusivity in internal solitary waves with reversing polarity is three times that of the 21 non-polarity-reversal region. The diapycnal diffusivity is higher at the front of one internal solitary 22 wave, and gradually decreases from shallow to deep water in the vertical direction. Our results also 23 indicate that (1) the enhanced diapycnal diffusivity is related to reflection seismic events; (2) 24 convective instability and shear instability may both contribute to the enhanced diapycnal mixing 25 in the polarity-reversing process; and (3) the difference between our results and Richardson-number-26 dependent turbulence parameterizations is about 2-3 orders of magnitude, but its vertical distribution 27 is almost the same. 28 29 [Key words] Internal solitary waves, Polarity reversal, Diapycnal mixing, Northeastern South China 30 Sea, Seismic oceanography. 31 32

33
Energy dissipation of internal waves enhances diapycnal mixing. Turbulence in the form of internal 34 wave breaking is the primary mechanism for modifying thermodynamic properties in the ocean (St. balance and provide energy for ocean mixing (Mackinnon and Gregg, 2003). Due to shoaling 41 internal waves and seafloor roughness, turbulent mixing on the continental shelves and slopes is 42 more variable than in the open ocean (Carter et al., 2005). Diapycnal diffusivity observed on 43 continental shelves and slopes can span four orders of magnitude (Gregg and Özsoy, 1999;Nash 44 and Moum, 2001). Internal solitary waves are a kind of nonlinear internal wave, which usually 45 carries a large amount of energy. Numerical simulations indicate that up to 73% of the internal wave 46 energy can be carried by internal solitary waves (Bogucki et al., 1997). Therefore, internal solitary 47 waves propagating to the continental shelf and slope can greatly change the local mixing.  The water is shallow on the continental shelf and slope near the Dongsha Atoll, so internal solitary 105 waves reach the transition point and their polarity changes from depression to elevation. In the 106 summer of 2009, the Guangzhou Marine Geological Survey (GMGS) set up a two-dimensional 107 seismic observation network in the Dongsha area. We found three internal solitary waves during the 108 polarity reversal process on the L1, L2, and L3 survey lines of the seismic data. The survey lines 109 are shown in Figure 1a, b. The streamer used in the acquisition process has a total length of 6 km 110 and 480 channels, the trace interval is 12.5 m, and the sampling interval is 2 ms. The airgun source 111 capacity is 5080 in 3 (1 in=2.54 cm), and the main frequency of the source is 35 Hz. The shot interval 112 is 25 m, and the minimum offset is 250 m. The time interval of shots is about 10 s. Survey lines L1 113 and L2 are the in-lines, which were from the southeast to the northwest. Survey line L3 is a cross-114 line, which was from southwest to northeast. We calculated the mean buoyancy frequency ( Figure  115 1c) of the region around seismic survey lines (latitude range 21.5°-22.5°, longitude range 116°-118°, 116 blue box in Figure 1a) by reanalysis temperature and salinity data with a water depth of 100-350 m. 117 This depth range matches the observation depth of the seismic data. Besides, since the buoyancy 118 frequency changes seasonally, we only selected the buoyancy frequency from July to August in 119 2009, which matches the seismic data observation time. The hydrographic data are provided by 120 Copernicus Marine Environment Monitoring Service (CMEMS). around seismic survey lines (blue box in (a)) and its 95% confidence interval (blue shadow).

128
After a conventional processing of the seismic data, an image of the ocean interior's structure can 129 be obtained. This image can be approximated as a temperature or salinity gradient map of the water 130 column (Ruddick, et al., 2009). The conventional processing of seismic data has 5 main steps, 131 including defining the observation system, noise and direct wave attenuation, velocity analysis, 132 normal moveout (NMO) and horizontal stacking. Then we use a bandpass filter to filter out low-133 frequency noise below 8 Hz and high-frequency noise above 80 Hz. According to the linear 134 characteristics of the direct wave, we use a median filter to extract the direct wave signal, and 135 subtract it from the original signal to achieve the purpose of attenuating the direct wave. 136 Subsequently, we sorted the seismic data from shot gathers into common midpoint gathers (CMPs).

137
Sound speed is a function of depth and obtained through velocity analysis, and then the NMO is 138 applied to CMPs according to the function to flatten the reflection seismic events of the water 139 column. When NMO is applied, the seismic wave with large offset will be stretched, and the 140 stretched seismic waves need to be cut off. Usually, the default method is to use a linear function to 141 remove the stretched seismic waves (Figure 2a). This may lose a lot of shallow reflection signals 142 (Figure 2b). Bai et al. (2017) used the common offset seismic section to supplement the missing 143 information in shallow water, but the low signal-to-noise ratio of the common offset seismic section 144 cannot guarantee the imaging quality. In order to retain more shallow reflection signal, we used a 145 custom function to cut off the NMO stretch (Figure 2c), thereby satisfying the imaging requirement 146 of the shallow water column (Figure 2d). Finally, the seismic section of the water column can be 147 obtained by stacking the processed CMPs. Due to the shallow water depth, the seismic data is 148 seriously affected by swell noise. We filtered out the components of stacked seismic data in wave 149 number range corresponding to swells in the frequency-wave number domain. A detailed description The red dotted line shows the cut off trace, the right part of seismic data is cut off. The unit TWT of (a) 156 and (b) is the two-way travel time of seismic wave from source to receiver. 157 158 2.2. Diapycnal diffusivity estimates from seismic data 159 160 Klymak and Moum (2007b) found that the horizontal wavenumber spectrum of the vertical 161 isopycnal displacement can be interpreted as the internal wave spectrum at low wavenumbers and 162 the turbulence spectrum at high wavenumbers. The high wavenumber components of spectrum are 163 dominated by turbulence, and the spectral energy follows the -5/3 power of the wavenumber. The 164 turbulence part of the horizontal wavenumber spectrum can be expressed by a simplified Batchelor 165 model (Equation 2-1), so the turbulence dissipation ε can be estimated from the observed 166 horizontal wavenumber spectrum. And diapycnal diffusivity can be calculated from Equation 2-2 167 (Osborn, 1980). (2 ) Where T ζ φ represents horizontal wavenumber spectrum, we use the seismic interpretation software to pick up reflection events in the seismic section ( Figure  181 3a). Then we calculate the vertical displacement of the reflection events. The vertical displacement 182 is the distance of the reflection evens deviate from the equilibrium position in the vertical direction. 183 We take the mean water depth of the reflection events as the equilibrium position. Note that the 184 choice of equilibrium position will not affect the calculation result. The spectral energy T ζ φ of the 185 vertical displacement in the horizontal wavenumber domain can be obtained by Fourier transform. 186 In practical applications, we use the slope spectrum (2 ) This conversion changes the wavenumber power law in the turbulence subrange from -5/3 to 1/3, 191 so that it can be distinguished from the internal wave subrange with -1/2 power law (-5/2 in the 192 displacement spectrum). In calculating the turbulence dissipation in the seismic section, it is 193 necessary to grid the section and calculate the dissipation in each grid separately. The horizontal 194 grid is set as 5 km, and the grid step 2.5 km. As the water depth in the seismic data is shallow, the 195 reflection seismic events are less in the vertical direction. In order to ensure more than two events 196 in each grid, we set the vertical grid to be 75 m and the grid step 37.5 m. In each grid, we calculated 197 the spectral slope of each event and took the average as streamlines, it is difficult to record reflection seismic data from those areas with closed streamline 224 at the resolution scale of seismic data. The regional density gradient recorded by the reflection 225 events still exists, and the streamline is parallel to the isopycnal at this time. While areas with closed 226 streamlines are strongly mixed, and the density gradient weakens or even disappears, which cannot 227 be recorded in seismic data. Unfortunately, the internal solitary waves we observed do not satisfy 228 the second assumption. The KdV equation can simulate internal solitary waves with small amplitude 229 and weak nonlinearity, but the polarity reversal of the large-amplitude internal solitary waves we 230 observed cannot be simulated well. Here we did not use theoretical models to fit observations. 231 Although there are studies using theory to successfully simulate the polarity reversal of internal 232 solitary waves (Liu et al., 1998;Zhao et al., 2003), it is difficult to match theories and observations. 233 We used the picked reflection seismic events to calculate the isopycnal displacement ( , ) x z η 234 ( Figure 4b). ( , ) x z η is the distance that reflection seismic events deviate from the equilibrium 235 position, which is determined by the mean depth of two shoulders of one internal solitary wave 236 ( Figure 4a). We smoothed ( , ) x z η with a spline function same as that was used for smoothing 237 turbulence dissipation, so that the resolution of wave-induced velocity is consistent with that of 238 turbulence dissipation. Therefore, the stream function can be expressed as (Holloway et al., 1999): 239 where c is the phase velocity of internal solitary waves. c can be estimated from pre-stack The wave-induced velocity here is on the seismic-resolution scale, which should be taken as its low-264 frequency component only. The results are insufficient to characterize the high-frequency 265 components. But this rough wave-induced velocity is useful, because our purpose of calculating 266 wave-induced velocity is for the vertical mixing scheme. The wave-induced velocity makes the 267 resolution scale of the mixing scheme equal to that of mixing parameters estimated from the seismic 268 data, and the two are comparable. In addition, the error of the wave-induced velocity is mainly  Figure  301 5a is the seismic section of survey line L1. It shows that the water depth becomes shallower from 302 southeast to northwest, and the bottom slope is steeper between 30-60 km. In the deep-water region 303 of 60-100 km, internal waves are developed, and the reflection seismic events fluctuate obviously. 304 Near the seafloor around 80 km, the reflection seismic events are uplifted and discontinuous, 305 forming a fuzzy reflection area. A mode-1 depression internal solitary wave can be identified at 53 306 km, indicating that the transition point has not been reached yet. The internal solitary wave has 307 reversed polarity at 40 km, and a packet of three elevation waves is formed. The reflection seismic 308 events are continuous here, implying no wave breaking. Five elevation waves can be identified 309 around 24-37 km, among which four elevation waves at 24 km may be formed continuously, while 310 the elevation wave at 37 km is formed later. 311 312 length of the head wave becomes wider and the slope becomes gentler. The leading wave is followed 316 by a packet of multiple elevation waves. The reflection seismic events are continuous in the whole 317 section. 318 319 L3 is a cross line whose observation direction is perpendicular to survey line L1 and L2 (Figure 5c). 320 There are multiple depression waves with large amplitudes of 20-35 km, and the reflection seismic 321 events are continuous. The wave polarity is reversing within 10-20 km, and the reflection seismic 322 events are discontinuous in this region. At 10 km, there is a large-amplitude elevation internal 323 solitary wave, and the wave front is almost parallel to seafloor. There is a large-amplitude depression 324 wave at 17 km, and the wave trough has interacted with topography. Most of the reflection seismic 325 events before 10 km are discontinuous and fuzzy, especially in the range of 6-10 km (Figure 5d). It We picked the reflection seismic events in the three sections ( Figure 7) and calculated the horizontal 338 slope spectrum using the method described in section 2.2. Figure 6 shows the average horizontal 339 slope spectrum of the three sections. We calculated the horizontal slope spectrum of all tracked 340 events and averaged in logarithmic space to determine the wavenumber of turbulence subrange. The 341 turbulence subrange of the survey line L1 section is 0.005-0.069 m -1 , as shown by the gray vertical 342 line in Figure 6a. The corresponding wavelength is 15-200 m. The average diapycnal diffusivity is 343 (7.0±1.2)×10 -4 m 2 /s, which is one order of magnitude larger than the open-ocean value (10 -5 m 2 /s).

344
The spectral energy in internal wave subrange is larger than that in turbulence subrange, indicating 345 that the energy is dominated by internal waves. This is confirmed by internal waves in the seismic 346 sections. The horizontal slope spectrum of the L2 section is shown in Figure 6b. The turbulence subrange is 353 0.008-0.068 m -1 , and the corresponding wavelength is 15-133 m. Compared with the survey line L1, 354 the turbulence shifts to a smaller scale. The spectral energy in internal wave subrange has the same 355 order of magnitude as the spectral energy in turbulence subrange, which indicates that the energy is 356 transferring to small-scale turbulence. This process is closely related to the polarity reversal of 357 internal solitary waves. The average diapycnal diffusivity is (1.5±0.1)×10 -3 m 2 /s, which is two 358 orders of magnitude larger than the background value. 359 360 Figure 6c is the horizontal slope spectrum of the L3 section. It can be seen from the spectrum that 361 the turbulence subrange is small, ranging from 0.011-0.07 m -1 . The corresponding wavelength is 362 14-89 m. The internal wave energy is larger and occupies a larger scale range. It can be seen from 363 the seismic section of survey line L3 (Figure 5c) that the wave amplitudes are large. It indicates that 364 the internal waves carry more energy, so the spectral energy in internal wave subrange is larger 365 (Figure 6c). In addition, there are many discontinuous and weak reflections in the seismic section black line is the spectrum, the gray shadow represents the 95% confidence interval, the gray dashed lines 373 represent the diffusivity contour, the black dashed lines represent the spectral slopes in internal wave 374 subrange, turbulence subrange and noise subrange, respectively. The gray vertical lines label the 375 boundaries of turbulence subrange. 376 377 Figure 6 shows that the spectral energy of the L1 section is smaller than that of the other two sections. 378 This may be because the imaging range of the L1 section is different. The observations in the L2 379 and L3 sections are the polarity reversal of internal solitary waves, while the L1 section includes 380 not only the polarity reversal process, but also internal waves in deep water. The spectral energies 381 of these two processes should be different. We calculated the average horizontal slope spectrum of 382 the polarity reversal region and the non-polarity reversal region, respectively (Figure 7). The 383 spectral energy of the polarity reversal region in L1 section is higher than that of the non-polarity 384 reversal region, so does diapycnal diffusivity (Figure 7a). It implies that the wave energy will 385 accelerate to dissipate and transfer to turbulence when its polarity is reversed. Compared with the 386 non-polarity reversal region, the turbulence subrange of the polarity reversal region is smaller. The 387 lower boundary of the turbulence subrange of the polarity reversal region is slightly larger than that 388 of the non-polarity reversal region. It indicates that the turbulence in this region has a smaller scale. 389 The diapycnal diffusivity in the polarity reversal region in L2 section is about 3 times that of the 390 non-polarity reversal region (Figure 7b). The turbulence subrange of the polarity reversal region in 391 L2 section is slightly larger than that of the non-polarity reversal region. From the L2 section, it can 392 be seen that the events are continuous during the polarity reversal process, which indicates that the 393 wave breaking is weak. The internal solitary wave gradually fissions into several tails during the 394 polarity reversal, and energy is dissipated constantly. Therefore, there will be a large turbulence 395 subrange in the lateral direction (Figure 7b). This process can dissipate much more energy compared The diapycnal diffusivity maps of the three survey lines are shown in Figure 8. Figure  km), the diffusivity of the head wave's front is higher, that is, where the slope of the wave front 418 becomes gentle. While the diffusivity of the two elevation waves followed the head wave is small. 419 It indicates that the mixing induced by internal solitary wave polarity reversal is stronger at the 420 beginning, and more energy is dissipated at this time. In the non-polarity reversal region (50-100 421 km), the diffusivity is low. The mode-1 depression internal solitary wave at 52 km increases the 422 diffusivity. There is an abnormal reflection area near the seafloor at 80 km, and the diffusivity is 423 high. In addition, there is also an area with increased diffusivity between 100-250 m at 93 km. This 424 may be related to large-amplitude internal waves.

426
The diffusivity map of the survey line L2 is shown in Figure 8b. The high value is mainly distributed 427 at the front of head wave during polarity reversal process (4 km), which is consistent with the 428 characteristics on L1. The diffusivity after the head wave is low, but it is still higher than that in 429 other regions. The diffusivity in the non-polarity reversal region is almost uniform. The two internal 430 solitary waves at 15-16 km did not increase the diffusivity. There is a low diffusivity area near the 431 seafloor around 16-18 km, which is caused by not tracked reflection seismic events in this area. The 432 diffusivity map of survey line L3 (Figure 8c) is similar to that of L2. The high value is distributed 433 in the polarity reversal region and the diffusivity of head wave is still high. However, unlike the 434 diffusivity map of L2, the high diffusivity is mainly distributed in the shallow part of the head wave 435 (water depth 50-120 m), while the diffusivity of the whole head wave in L2 is high. In the non-436 polarity reversal region, the diffusivity is small and the distribution is uniform too. The diffusivity 437 near the seafloor at 25-27 km is slightly lower than other regions.  salinity, so the reflection seismic events are related to density gradient. The enhanced mixing reflects 451 the structure of density gradient, thereby changing the appearance of the reflection seismic events. 452 Understanding the relation between diffusivity and reflection seismic events can help us analyze 453 the spatial distribution of diapycnal mixing. Figure 8 shows that the reflection seismic event in the 454 high diffusivity region is obviously different from that in the low diffusivity region. In the high 455 diffusivity area (red in Figure 8), the reflection seismic events are fuzzy, discontinuous or bifurcate.

456
While in the low diffusivity area (yellow and blue in Figure 8), the reflection seismic events are 457 clear and continuous. This is because regions with high diffusivity are strongly mixed. The density 458 gradient is smeared by mixing, so that it affects the appearance of reflection seismic events. For 459 example, in the polarity reversal region of three seismic sections, the diffusivity is high, and the 460 reflection seismic events are fuzzy and discontinuous. Especially in the range of 5-10 km in Fig. 8c,  461 the events are obviously broken and weak. The diffusivity is low in areas where the events are clear, 462 such as the region near the seafloor around 45-50 km and the region near the sea floor around 93-463 99 km in Figure 8a, and the region near the sea floor around 24-27 km in Figure 8c.

465
The diffusivity is not only related to the continuity of the reflection events, but also related to the 466 fluctuation intensity of the events. The greater the fluctuation intensity of the events, the higher the 467 spectral energy, and the greater the diffusivity value. There is a mode-1 depression internal solitary 468 wave at 50-58 km in Figure 8a, and the reflection seismic event is clear and continuous at 180 m. 469 But the diffusivity is high, because the reflection events fluctuate more strongly. It can be seen from 470 the figure that, in addition to the amplitude of the internal solitary wave, there are also many high-471 frequency waves at the shoulders of the internal solitary wave. These waves increase the spectral 472 energy and result in a higher diffusivity. In addition, the reflection seismic events before 4 km in 473 Figure 8b is continuous and without obvious fluctuations, but the diffusivity is higher. It can be seen 474 from the figure that the reflection events of this region are thicker than that of other regions. The 475 seismic data processing of the three sections in Figure 8 is the same, so the thicker events in Figure  476 8b do not stem from the low frequency of seismic waves. We think this may be caused by small-477 scale mixing between layers, such as K-H instability. Figure 9 is an enlarged view of 2-3 km in the 478 seismic section of L2 (Figure 5b). The wavelength of the seismic wave (red line) at 80 m is larger 479 than that at the seafloor, which is formed by the overlap of multiple wavelets. It can be seen from The vertical scale of the K-H instability is small and usually appears on the isopycnal. On one hand, 485 K-H instability weakens the density gradient so that the reflected seismic wave energy is reduced. 486 On the other hand, the vertical scale of K-H instability is lower than the seismic wave resolution (a 487 quarter of the seismic wave wavelength), so it causes overlapped wavelets and stretched wavelength 488 (Figure 9). Therefore, the reflection event in this area is thicker. Besides, the horizontal scale of the 489 K-H instability train is large, which may explain the larger turbulence subrange on the horizontal 490 slope spectrum (Figure 7b). 491 492 The high diffusivity is mainly in the leading internal solitary wave during the polarity reversal. We 524 suggest that strong mixing may be caused by internal wave breaking due to convective instability. 525 In Figures 8a and 8c, the reflection seismic events are obviously discontinuous in the high 526 turbulence area, indicating that the density gradient is weakened by internal wave breaking. The 527 trough of the internal solitary wave decelerates first when the polarity is reversed (Shroyer et al.,528 2008), which makes the Froude number (Fr) greater than 1 and causes convective instability. This 529 phenomenon can be found in other observational data. In the high-frequency acoustic section, the 530 backscatter at the top of internal solitary wave is increased when it changes from depression to 531 elevation wave (Orr and Mignerey, 2003), which indicates that the turbulence of the front increased. 532 However, in the seismic section of Figure 8b, we did not find breaking at the front the polarity 533 reversal internal solitary wave. The strong mixing of this internal solitary wave may be induced by 534 shear instability (Figure 9). Therefore, both convective instability and shear instability are 535 responsible for the enhanced mixing in this process. In addition, the non-polarity reversal region in 536 Figure 8a has a higher diffusivity in 50-150 m than other regions. This range is in the thermocline 537 ( Figure 1c). The internal waves usually greatly increase mixing in the thermocline, which is related 538 to shear instability of internal waves (Mackinnon and Gregg, 2003). Shear instability is an important 539 mechanism of internal wave dissipation (Farmer and Smith, 1978), and it more likely occurs in 540 nonlinear internal waves than convective instability (Zhang and Alford, 2014). The results of high-541 frequency acoustic observations show that the enhanced backscatter at the bottom of the thermocline 542 represents higher shear instability when the internal solitary waves are shoaling (Orr and Mignerey,543 2003), which is consistent with the depth range of high diffusivity in our results.

545
What is inconsistent with the observed distribution of mixing is that our results do not show 546 diffusivity in the bottom boundary layer. Because our seismic data was collected in summer, the 547 strong stratification at this time limits the vertical range of the bottom boundary layer (Mackinnon 548 and Gregg, 2003). So that the bottom boundary layer near the Dongsha Atoll is thin and lower than 549 the thickness that can be recorded by seismic data. So, the diffusivity we calculated does not include We compared the vertical distribution of diffusivity with the vertical mixing scheme of internal 561 wave breaking proposed by Vlasenko and Hutter (2002). Although Klymak and Legg (2010) also 562 proposed a mixing scheme for internal wave shoaling and achieved good results in numerical 563 simulation, we cannot use that method to calculate mixing parameters because of lacking high 564 resolution density observation data. Figure 10 shows the vertical distribution of diffusivity from 565 seismic data (solid line) and the diffusivity calculated from mixing scheme (dashed line) at 4 566 positions of the three survey lines (black arrows in Figure 8). The reflection events in the L3 section 567 are broken, and it cannot be guaranteed that the events are parallel to the streamline. Therefore, we 568 did not use the method described in section 2.3 to calculate the wave-induced velocity, and thus did 569 not obtain the diffusivity of the mixing scheme. It can be seen from Figure 10 that the turbulent 570 diffusivity gradually decreases from shallow to deep water. Except for the local low diffusivity value 571 in the deep water at the position D of Figure 10b and 10c, the diffusivity reduction rate at other 572 locations is similar. Figures 10a and 10b show that the parameterized diffusivity is nearly 2--3 orders 573 of magnitude smaller than our result, but they have a similar trend of change. In Figure 10a

601
We have observed the polarity reversal of internal solitary waves by reflection seismic data near the 602 Dongsha Atoll in the South China Sea, and calculated their slope spectra ( Figure 6) and diapycnal 603 diffusivity (Figure 8). The results show that the average diapycnal diffusivities of the three survey 604 lines are about two orders of magnitude greater than the open-ocean value. We calculated the 605 average spectral slope of the polarity reversal and non-polarity-reversal regions (Figure 7), and 606 found that the former is about 3 times larger than the latter. The diffusivity maps reveal that 607 horizontally high diffusivity is mainly in the leading wavefront of an internal solitary wave in 608 reversing polarity, and vertically high diffusivity is mainly in the thermocline (50-100 m).

610
We analyzed the relation between reflection seismic events and diapycnal diffusivity. The result 611 indicates that continuous and clear reflection events correspond to low diffusivity, while 612 discontinuous or fuzzy events correspond to high diffusivity. The strength of the events also affects 613 the magnitude of diffusivity. The stronger the fluctuation, the higher the spectral energy, and the 614 higher the diffusivity. In addition, we observed an area of high diffusivity with a large horizontal 615 scale in L2, and the reflection events did not appear to be discontinuous or fuzzy. We suggest that 616 this enhanced mixing may be induced by the K-H instability (Figure 9). The vertical scale of the K-617 H instability is smaller than the resolution of our seismic data, so we cannot observe clearly in the 618 seismic data. But its high-energy characteristics can be recorded by reflection events. 619 620 Our results show that shoaling internal solitary waves enhance local mixing. The magnitude order 621 of diapycnal diffusivity is consistent with previous studies. We suggest that there are two 622 mechanisms that could account for the enhanced mixing. On one hand, the polarity reversal of 623 internal solitary waves results in convection instability, which induces internal solitary wave 624 breaking. This mechanism appears at the leading edge of one internal solitary wave in the survey 625 lines L1 and L3. The discontinuous reflection events indicate that the internal solitary wave is 626 broken. While in the seismic section of L2, the reflection events are continuous and clear at the 627 leading edge of the internal solitary wave and other strong mixing areas in the three sections. Such 628 strong mixing may be caused by shear instability.

630
We picked four positions from the diffusivity maps to analyze the vertical distribution of diapycnal 631 diffusivity ( Figure 10). Our result shows that the diffusivity gradually decreased from shallow to 632 deep water (excluding the bottom boundary layer). Compared with previous one mixing scheme, 633 the parameterized diffusivity is about 2-3 orders of magnitude smaller. This means that the mixing 634 scheme underestimates mixing induced by internal solitary wave shoaling near the Dongsha Atoll. 635 However, the vertical pattern of the parameterized diffusivity is consistent with our result.

637
Appendix: The uncertainty of diffusivity from seismic data

639
According to formulas 2-1 and 2-2, the parameters used in calculating diffusivity are buoyancy 640 frequency (N), mixing coefficient (Γ), and the Kolmogorov constant (C T ). It can be seen from the 641 formula that the diffusivity is proportional to N. The mean deviation of N we used is about 2% 642 (Figure 1), so the uncertainty of the corresponding diffusion rate is about 0.008 logarithmic units. 643 In addition, diffusivity is proportional to Γ -1/2 , which is 0.1-0.4 (Mashayek et al., 2017), so the 644 corresponding uncertainty of diffusivity is about 0.15 logarithmic units. Similarly, diffusivity is 645 proportional to C T -3/2 . The value of C T is 0.3-0.5 (Sreenivasan, 1996), and the corresponding 646 diffusivity uncertainty is about 0.15.

648
In addition, the key reason for the uncertainty of diffusivity in our calculation is the fitting of the 649 Batchelor model and the slope spectrum (Figure 3b and Figure 6). We evaluated the uncertainty 650 based on the least squares standard deviation between the Batchelor model and the slope spectrum. 651 The uncertainty of the three diffusivity maps in Figure 8 is shown in Figure A1. 652 653 654 Figure A1. The diffusivity uncertainty along survey lines L1 (a), L2 (b) and L3 (c). 655

656
Code and data availability. The bathymetry data were provided by the General Bathymetric