Numerical simulation of water free-surface flows through a front-tracking lattice Boltzmann approach

Hydraulic models for simulating free-surface flows using three-dimensional (3D) schemes have been increasingly used in the last decade for several reasons. The large availability of high-resolution 3D digital topographic information on the one hand and the increasing computer power and efficiency on the other have paved the way for scientists to develop and test always more complex modelling algorithms. Moreover, the need for detailed deterministic hydraulic simulations for urban and land planning projects (Wu & Chau 2006; Manciola et al. 2009) with specific regard to the physically based characterization of the impact of territorial transformations and man-made interventions on the hydraulics of the domain – that means by providing a detailed numerical characterization of the water surface geometry and associated flow depths and velocities – is challenging for hydraulic modellers to accurately simulate the actual conditions with respect to the design (post-intervention) scenario.

At present, standard models for characterizing complex systems are based on the Navier–Stokes continuum equations (NSE). The shallow water equations (SWE) have proven to be computationally efficient and suitable for representing a wide range of free-surface flows (Lai & Khan 2012; Petaccia et al. 2013; Pu et al. 2013). Nevertheless, flows are often fully 3D with strong inertial effects and high free-surface gradients that break the main hypothesis of the SWE (Chau & Jiang 2001; Stoesser et al. 2003; Chau & Jiang 2004; Biscarini 2010; Biscarini et al. 2010; Ferrari et al. 2010; Filonovich et al. 2013).

A possible remarkable alternative to the NSE is represented by the lattice Boltzmann method (LBM), which is based on a discretized Boltzmann kinetic equation (Higuera et al. 1989; Qian et al. 1992; Ubertini et al. 2010). The LBM has been recently adopted for different applications: the simulation of flows in porous media (Ginzburg & D'Humières 2007; Zarghami et al. 2014), and for multicomponent and multiphase flows (Falcucci et al. 2010) among others. The LBM has been also implemented to simulate free-surface flows by employing the shallow water hypothesis (Zhou 2004; Budinski 2012).

In this technical note, a fully 3D front-tracking (FT) formulation of the LBM is proposed. The method, originally presented by Körner et al. (2005) for foam simulations, was further developed by some authors for free-surface hydraulic flows (Biscarini et al. 2011, 2013; Falcucci et al. 2011). The numerical model is here tested against experimental measurements of a partial symmetrical collapse of a vertical dam embankment. The novelty aspect of the present work with respect to previous findings is mainly represented by the detailed comparison of the transient flow field and free-surface evolution, given the availability of literature experimental measurements at different points of the domain and at different times.

The set of discrete speeds must be chosen in order to guarantee mass, momentum and energy conservation, as well as rotational symmetry, in order to recover the correct fluid dynamic equations in the macroscopic limit. Among the possible options proposed for 3D flows, a 19-direction (i.e., b = 19) cubic lattice model is implemented here (namely D3Q19), as a good balance between computational reliability and efficiency (Succi 2001). The resulting numerical mesh is a structured cubic Cartesian lattice where the particles may only reside on the mesh nodes.

In the proposed model, boundary nodes are solved using the half-way bounce-back scheme (Skordos 1993; Inamuro et al. 1995), which allows modelling complex boundaries in terms of elementary mechanical rules (bounce-back, reflections). The bounce-back scheme can be applied to several boundary conditions usually found in hydraulic problems, such as no-sip, free-slip, open flows and moving boundaries. Curvilinear surfaces can be handled through a stair step approximation, which is sufficiently accurate for the typical geometries of the large-scale hydraulic flows discussed herein.

Turbulence is modelled through the large-eddy simulation (LES) approach; in particular, we employ the Smagorinsky–Lilly sub-grid scale with a Smagorinsky constant of 0.15 (Galperin & Orzag 1993).

As the gas dynamics is neglected, the distribution functions coming from the gas cells are never accessed and need to be reconstructed. Here, we use the reconstruction scheme proposed by Körner et al. (2005) by setting that the water pressure equals the air pressure at the interface.

Figure 1(a) shows the error decay with the grid size. Figure 1(b) shows the wave height differences at t = 0.5 s simulated by LBM equation in a test point at position x = 8 m. The resulting convergence order of the scheme is about 2.9.

The proposed LBM is then applied for the simulation of the transient free-surface flow deriving from the partial symmetrical collapse of a vertical dam embankment. A detailed analysis of the transient flow field is performed by comparing numerical results with laboratory experimental measurements on a dry bed (Fraccarollo & Toro 1995; Ferrari et al. 2010). The dam rupture is simulated by a sudden removal of a gate, considered instantaneous. Water in the tank is initially in steady state with a depth equal to 0.6 m. Measurements of pressure, water depth and velocity components for a period of 10 s are available at gauges located in the tank and downstream along the plane bottom at the positions shown in Figure 2. After that time, the water flow slowly varies and reaches the quasi-steady state with a hydrostatic pressure distribution that characterizes almost the entire domain.

Figure 3 shows both numerical and experimental hydrographs at points 4 and 0 along the gate. At the initial stage of the tank emptying (i.e., first 2 seconds), strong three-dimensional inertial effects are clearly visible, as the high free-surface gradient close to the gate produces a local non-negligible acceleration component on the gravity direction. This is more evident at point 0, as point 4 is closer to the wall (i.e., the gate portion that did not collapse).

In Figure 4, the experimental x-components of the velocity, V_x, are compared with our results for point 0. Eight measurements are provided at different heights from the bottom, starting from 50 mm from the bottom up to 0.4 m, every 50 mm. Both experimental and numerical velocity profiles at point 0, shown in Figure 4, are characterized by a typical shape: a fast (i.e., <1.5 s) quasi-linear increase to the maximum value followed by a monotonic velocity decrease in the longitudinal direction. The ending time of the curves reduces with higher z values, as the water surface elevation decreases after the gate break, leaving the higher gauges outside the flow domain.

In this study, a fully 3D FT variant of the LBM is demonstrated to be a valid model in the treatment of instantaneously varying flows interacting with artificial structures. In particular, we demonstrate that the proposed LBM is able to characterize in detail the complex hydraulics of a test case like the symmetrical dam break flows propagating within a rectangular domain. The numerical model is able to predict the experimental hydrographs and velocity profiles at different points of the domain and at different times.

Results are, thus, promising for the use of the proposed LBM for practical engineering applications. Nevertheless, future work should include: the treatment of moving bodies and curvilinear boundaries; the introduction of grid-refinement techniques to further increase computational efficiency; the introduction of additional phases (i.e., solid, landslide material); the development of sensitivity analyses on modelling and geometric parameters (e.g., roughness effects, mesh resolution); the simulation of real case scenarios of unsteady flows through natural topography; and the validation of simulation results with real case measurements.