Estimating pressure and internal-wave flux from laboratory experiments in focusing internal waves

Abstract

Instantaneous measurements of pressure and wave flux in stratified incompressible flows are presented for the first time using combined time-resolved particle image velocimetry (PIV) and synthetic schlieren (SS). Corrections induced by variations of the refractive index in this strongly density-stratified fluid are also considered. The test case investigated here is a three-dimensional geometry consisting of a Gaussian ring-type topography forced by an oscillating tide representative of geophysical applications. Density and pressure are reconstructed from SS or PIV in combination with linear theories and combined SS-PIV. We perform a direct comparison between the experimental results and three-dimensional direct numerical simulations of the same flow conditions and control parameters. In particular, we show that the estimated velocity or density and the hence wave flux from linear theory solely based on SS or PIV can be flawed in regions of focusing internal waves. We also show that combined measurements of SS and PIV are capable of circumventing these limitations and accurately reproduce the results computed from the DNS.

Graphic abstract

Appendices

Appendix 1: Test case: beam propagation

In this test case, internal wave beam propagation in a two-dimensional domain is considered similar to the problem described in Allshouse et al. (2016). Internal wave beams are generated by adding a forcing term to the x-direction of the momentum Eq. (34b), which writes

$$\begin{aligned} F(z, t) = A(z) \cos (\omega t - k_z z), \end{aligned}$$

(35)

where $A(z) = A_0 \times \exp (-(z-0.50625)^2/0.22)$. This synthetic data is then used to validate our approach and make sure that the pressure solver was reproducing the pressure field produced by the internal wave field.

Fig. 16

a Pressure $(p(x, z, t=10))$ from the numerical simulation, b pressure estimated from from synthetic velocity fields and c error measured between a and b, normalized by the maximum value of the pressure

The comparison between the numerical simulation and its prediction is shown in Fig. 16a–c in the case of a fully developed flow, computed after 10 cycles and where the snapshots are measured at $t=10$. The solution was computed interpolating the data onto an equidistributed $[128\times 128]$ grid. the pressure field is generally very well captured, as shown in Fig. 16a, b with only very little discrepancy between SOMAR’s solution and the present pressure solver, which might be explained by the resolution of the interpolation and the use of a sponge layer in SOMAR but not in our approach. It is worth noting that the two fields are nearly identical, with a maximum error of 10%, measured in the sponge layer, far from the internal wave beam. Otherwise, the error was found less than 6% in the rest of the domain and validates our approach.

Appendix 2: Least-squares solution to the Poisson problem

Fig. 17

Comparison between results obtained using the Poisson solver used in the analysis (left) and the least-squares method (right) in the ($x, y=0, z$). Density perturbation $\rho (x, z)$ in (a, b), pressure p(x, z) in (c, d), and horizontal wave flux $\partial _x pu$ in (e, f)

Pan et al. (2016) recently proposed guidelines on the optimal spatial resolution for the Poisson solver is proposed to limit the error propagated. A notable work commenting on the error propagation dynamics of PIV Poisson-based pressure estimation has been carried out by the same author (Pan et al. 2018). The essence of this work emphasizes the importance of boundary conditions in obtaining an accurate pressure estimate. They noted that the Neumann type of boundary conditions is prone to errors and is to be avoided wherever possible. However, it becomes mandatory in PIV domain focusing on rotational flow regions. A recent review by Liu and Moreto (2020) showed that other methods should be emphasized such as omnidirectional integration method of Liu and Katz (2006), Liu et al. (2016) or the least squares approaches based on singular value decomposition as recently applied by Jeon et al. (2018) (see also references therein).

The procedure solves the system of equations (29a) for pressure in a least-square sense and we used the reconstruction from gradient technique proposed by Harker and O’Leary (2008, 2013) to obtain the pressure field p(x, z). We also used the same technique to obtain the density $\rho$ solving for (8).

As shown in Fig. 17a, b, the density field obtained using this method or the Poisson solver are nearly identical. This is to be expected since the boundary conditions are nearly zero on all sides except near the surface at $z=0$. In fact, the pressure field obtained from the reconstruction of gradients and the Poisson problem follows the same observation. The largest source of error in the PIV comes from the acceleration term which is the next largest term after the buoyancy term in the boundary condition (31). Since this term is nearly zero on all sides except for the top boundary, the pressure and waveflux reconstructed from this method are very similar to that of the Poisson problem. In fact, in the very particular case presented here, the Poisson problem performs slightly better as the acceleration term accumulates errors when integrating (29a) from gradients. Nevertheless, we recognized that in other situations, the reconstruction from gradients may outperform the Poisson problem inversion with Neumann boundary conditions. Both methods are available in the supplementary material attached to the manuscript.

Appendix 3: Three-dimensional vertical velocity forcing calculated from experiments and used in the direct numerical simulation

The boundary condition used in the DNS at $z=-0.0125$  m was computed using a fit from the experiment at the same height in the tank. The forcing only considered the vertical velocity $w(x, y, z=-0.0125m, t)$ which was implemented using the expression in a MATLAB code given hereinafter. This expression is a fit of Gaussian modulated waves to the first four POD modes calculated from the time history $w(x, y=0, z=-0.125, t)$ which captured more than 85% of the energy associated with the vertical velocity at this height. A mode-1 symmetry of revolution was assumed following the expression in Eq. (10). The forcing is given in Matlab format as follows: