More insights into the lock-exchange test

In the past months, we have made significant advances in putting the NUMO prototype, which delivered results reported in this blog so far, back into the GNuMe framework. GNuMe stands for [G]alerkin [Nu]merical [M]odeling [E]nvironment and is a framework that hosts both NUMO and NUMA – Non-hydrostatic Unified Model of the Atmosphere. This merger allows NUMO to leverage a lot of  NUMA technology (both current and future, as the codes share the same infrastructure) and more rapid advancement in NUMO development. The benefits go both ways, as the development of NUMO brings in new features to NUMA as well. Now that we are happy with how the code infrastructure looks like, we can go back to evolving NUMO into a truly unstructured non-hydrostatic model.

Ever since we have obtained the first lock-exchange results, two issued bothered me. Firstly, we experienced some unphysical instabilities and mesh imprinting close to domain boundaries, particularly in the domain corners. Those instabilities eventually led to simulation crash at later times, way past the usual end time of the simulation when the fronts reach the wall. Regardless, we definitely cannot allow this in a production code. Secondly, we compared our results to those available in the literature, but there was quite a bit of spread between reported Froude (Fr) numbers among different researchers. Recently I have revisited those simulations, and I believe I found answers to both issues.

Since we were getting mostly correct physical behavior of the flow, I suspected  that the instabilities are simply due to not enough stabilization of the temperature field. Following Hiester et al. (2011), we used kinematic viscosity ν=1e-6 m^2/s for velocity, and no thermal diffusion. Given no numerical diffusion in the high-order continuous Galerkin method used for those simulations, we always require some stabilization, either through artificial diffusion, filtering, or other methods. I was using Boyd-Vandeven filters, but this was clearly not enough. I decided to introduce some physical diffusion to see if this will stabilize the simulation. Since I was locked in to kinematic viscosity of water at 20C, I decided to use matching thermal diffusion of κ=0.148e-6 m^2/s. The result was a very stable lock-exchange simulation, which I run till 500s with no trace of growing instabilities. Here is the video of the flow:

This led to a further investigation of the lock-exchange results. As mentioned before, Hiester et al. (2011) report the Grashoff number is Gr=1.25e6, with kinematic viscosity ν=1e-6 m^2/s and no thermal diffusivity (κ=0 m^2/s), which leads to Schmidt number Sc =  ν/κ = ∞. A similar set-up is found in the simulation by Fringer et al. (2006), but the reported Froude number is significantly different than Hiester et al. (see the table here). The authors, however, compare the result with the DNS simulation by Härtel et al. (2000), which provides yet another Fr value, and note the possible cause of this discrepancy. The simulation of Härtel et al. was made (apart from significantly higher resolution) with Sc=0.71, which with ν=1e-6 m^2/s gives the thermal diffusivity of κ=1.41e-6 m^2/s. Fringer et al. point out that even though they do not explicitly use thermal diffusivity, the numerical diffusion of the method (first-order upwinding) makes the Sc number lower than this of Härtel et al., and therefore their fronts are slower.

This led me to try two simulations with different values of Sc number, one with Sc=0.71 like in Härtel et al., and the other with Sc=6.74, which is the simulation shown in the video above. I have compiled the results in a table below, taking into account Sc numbers (all simulations with Gr=1.25e6).

Sc no-slip free-slip
NUMO 6.74 0.420 0.482
NUMO 0.71 0.407 0.475
Härtel et al. (2000) 0.71 0.406 0.477
Hiester et al. (2011) 0.417 0.482
Fringer et al. (2006) 0.396 0.428

It turns out that NUMO results match very well the DNS simulation by Härtel et al. when using matching Schmidt numbers, and the simulation of Hiester et al. when using high Sc. Even though Hiester is not using any thermal diffusivity explicitly, it is hard to evaluate the numerical diffusion of their model. The front speed in Fringer et al. is clearly slower than in other simulations, indicating much lower Schmidt number due to high numerical diffusion of the model.

This result gives me another bit of confidence in NUMO results, as it agrees with the DNS simulation very well, while also matching the result of high resolution adaptive finite element code of Hiester et al. (2011) by using only physical diffusivity.

References

Fringer, O.B., Gerritsen, M. & Street, R.L., 2006. An unstructured-grid, finite-volume, nonhydrostatic, parallel coastal ocean simulator. Ocean Modelling, 14(3–4), pp.139–173.

Härtel, C., Meiburg, E. & Necker, F., 2000. Analysis and direct numerical simulation of the flow at a gravity-current head. Part 1. Flow topology and front speed for slip and no-slip boundaries. Journal of Fluid Mechanics, 418, pp.189–212.

Hiester, H.R., Piggott, M.D. & Allison, P.A., 2011. The impact of mesh adaptivity on the gravity current front speed in a two-dimensional lock-exchange. Ocean Modelling, 38(1–2), pp.1–21. Available at: http://dx.doi.org/10.1016/j.ocemod.2011.01.003.