Numerical Simulation of Two-Phase Flow for Wave Propagation and Breaking near Submerged and Vertical Breakwaters R. Bakhtyar1,2, C.E. Kees1, M.W. Farthing1, and C.T. Miller2 1 U.S. Army Engineer Research and Development Center; 2University of North Carolina at Chapel Hill where ν is kinematic viscosity. Introduction Prior Experimental Results Define the sub-grid scale stress tensor as Submerged breakwaters avoid the generation of major reflected waves that affect the beach. They are very useful for erosion control and beach protection by decreasing the wave height[7]. In addition, standing waves can play a significant role in hydrodynamics near the vertical breakwater. However, an accurate understanding of the wave propagation/breaking over a submerged breakwater, and wave overtopping in front of a vertical breakwater has yet to be achieved. and combine with Eqn (8) to give the filtered NS equation of the form Model Equations The Smagorinsky model was used to approximate the sub-grid-scale stress tensor as The Navier-Stokes (NS) equations are used to model Newtonian and incompressible turbulent flow for the time interval [0,Ts] and a physical domain Ω, which are specified as [5]: The general model formulation is specified by an overall conservation of mass equation given as ∂ρ + ∇· (ρv) = 0, Ω ∈ IR2, t ∈ [0, Ts] ∂t (1) and the general NS equation for the case of a constant dynamic viscosity given by 2 ∂(ρv) + ∇· (ρv ⊗ v) + ∇p − 2µ∇·d + µ∇ [(I:d) I] − ρg = 0, ∂t 3 2 Ω ∈ IR , t ∈ [0, Ts] , and the rate of strain tensor is defined as i 1h d = ∇v + (∇v)T 2 (2) and the incompressible NS equation can be written as (5) Turbulence modeling The inner surf and swash zones are complex regions of turbulent flow [8]. Since turbulence has a profound influence on the flow field and momentum transfer, it must be modeled [1]. Because of the complexity of turbulent flows, direct numerical simulation (DNS) to resolve all length scales of the flow field is typically not feasible for the surf and swash zones. A popular alternative to DNS is modeling the conservation of momentum using the Reynolds Averaged NS (RANS) equation, which yields an expression for the average velocity but depends upon a modeling approximation of the dyadic product of fluctuation velocities that arise in the Reynolds stress tensor. The RANS approach thus produces an estimate of the average velocity, while approximating the Reynolds stress effects on this velocity. A key point is that RANS models only yield averaged behavior. A middle ground between DNS and RANS are approaches that attempt to resolve the velocity directly, similar to DNS, but only over a set of dominant energycontaining length scales. The shorter length scale effects are modeled using a closure scheme for the stress contributions due to velocity fluctuations, similar to a RANS approach. Large eddy simulation (LES) is an approach that resolves a range of length scales known to contain the majority of the kinetic energy in the system, while approximating the unresolved scales by modeling the stress contributions using a filtered velocity field. The LES model thus provides a means to simulate important, highly energetic, aspects of the velocity field, while still greatly reducing the computational burden associated with a DNS approach. In the LES model, the velocity is written in filtered form as 0 0 v = v + v = hvi + v , (6) where an implicit spatial filter resulting from the numerical approximation is used, and v0 is the fluctuation velocity. Applying the spatial filter to Eqn (4) yields ∇·v = 0 , (7) and applying the spatial filter to Eqn (5) and simplifying yields the filtered incompressible NS equation given by ∂v 1 + ∇·vv + ∇p − ν∇2v − g = 0 , ∂t ρ τ SGS = −2νT d , (10) (11) where eddy viscosity is approximated as 2 νT = (Csh) 2d:d 1/2 , the filtered rate of strain tensor is i 1h T d = ∇v + (∇v) , 2 (12) 2 (8) h 2d:d ν 38.0m 32.9 2.1 wall wave screens 6.0 3.0 wave paddle (13) 0.6m h is the grid spacing in the numerical approximation, and a dynamic approximation of the Smagorinsky constant is computed as [9] 0 Reh ≤ 1 (14) CS = −0.92 Reh > 1 0.027 × 10−3.23Re Reh = (4) ∂v + ρv·∇v + ∇p − µ∇2v − ρg = 0, ∂t Ω ∈ IR2, t ∈ [0, Ts] . 1 ∂v + v·∇v + ∇·τ SGS + ∇p − ν∇2v − g = 0 . ∂t ρ (3) For the case in which the material derivative of the density vanishes, the so-called incompressible case, the conservation of mass equation reduces to ρ (9) h=0.3 sand bed 1:30 0.15 0.45 1/2 . (15) Figure 3: proteus-mprans numerical reults for [2] Air/Water Flow Validation Test Set Figure 1: Small experimental tank from [13] and the grid Reynolds number is where ρ is the density, t is time, v is the fluid velocity, Ts is the extent of the temporal domain, p is the fluid pressure, µ is the dynamic viscosity, I is the identity tensor, and g is the gravitational acceleration vector. ∇·v = 0, Ω ∈ IR2 τ SGS = vv − v v , The numerical model is validated alongside the two data sets [2, 13]. The numerical results show that a vortex forms when waves spread over a submerged breakwater. The numerical results also indicate that the maximum turbulent kinetic energy, velocities and viscosity occur on the top of breakwater. The results show that for fully standing waves, velocities and wave height are higher in comparison to the comparable partially standing wave case (for the vertical breakwater). These calculations are in overall agreement with earlier observations, while the numerical model describes the water and air phase characteristics in greater detail than current measurements. The model provides a valuable method to advance mechanistic understanding of hydrodynamic characteristics near the breakwater in the nearshore area. Numerical Results The use of turbulent two-phase flow models is still relatively imature for applications in coastal and hydraulic structures. In order to verify and validate the proteusmprans model and place computational modeling on a firmer footing for these applications, we have begun developing a test set spanning dambreak, wave, and flow processes in two- and three-dimensions. The test set can be accessed by contacting the authors. Free surface tracking There are several existing techniques to capture the airwater interfaces. A hybrid level set/Volume Of Fluid (LS-VOF) method denotes one of the approaches for tracking air-water interface that has been used effectively to solve problems in various water resources fields [3, 10, 12]. In the VOF method (Fig. 2a), the filled fraction of computational cells is used for defining the interface [4]. In the LS technique (Fig. 2b) the interface is defined implicitly as the zero level set of scalar function [11] and, consequently, no mesh adjustment is needed to define the air-water interface (i.e., level sets with zero, positive and negative values represent air-water interface, water and air phases). Accurate and suitable approximations of free surface are obtained, by coupling the LS and VOF methods [6]. A hybrid LSVOF method is an Eulerian free surface technique with high-order estimations, and is very compatible with the NS type equations. A precise free surface tracking with a satisfactory computational cost can be obtained using a hybrid LS-VOF technique. This hybrid method reduces the mass conservation errors that exist in the classical level set methods. In air-water two-phase flow, the level set governing equations for determining the free surface can be written as follows: ∂ψ + v·∇ψ = 0 , ∂t (16) k∇ψk = 1 , (17) ρ = ρa [1 − Hf (ψ)] + ρw Hf (ψ) , (18) Figure 4. The air-water-vv project on github. Figure 2: proteus-mprans numerical results for the experiment from [13] References [1] Bakhtyar, R., D. A. Barry, A. Yeganeh-Bakhtiary, and A. Ghaheri, 2009: Numerical simulation of surf-swash zone motions and turbulent flow. Advances in Water Resources, 32 (2), 250–263. [2] Chen, J., C. Jiang, S. Hu, and W. Huang, 2010: Numerical study on the characteristics of flow field and wave propagation near submerged breakwater on slope. Acta Oceanologica Sinica, 29 (1), 88–99. µ = µa [1 − Hf (ψ)] + µw Hf (ψ) (19) and   0 ψ<0 Hf (ψ) = 0.5 ψ = 0 ,  1 ψ>0 [4] Hirt, C. W. and B. D. Nichols, 1981: Volume of fluid (VOF) method for the dynamics of free boundaries. Journal of Computational Physics, 39 (1), 201–225. (20) [5] Kees, C. E., I. Akkerman, M. W. Farthing, and Y. Bazilevs, 2011: A conservative level set method suitable for variable-order approximations and unstructured meshes. Journal of Computational Physics, 230 (12), 4536–4558. [6] Kees, C. E. and M. W. Farthing, 2011: Parallel computational methods and simulation for coastal and hydraulic applications using the Proteus toolkit. PyHPC11 Workshop, Supercomputing 11. where ψ is a level set function that is used to describe the fluid distribution; ρa and µa are the density and dynamic viscosity of air, respectively; Hf is the Heaviside function; and ρw and µw are the density and dynamic viscosity of water, respectively. Kees et al. [5] showed that the LS method, in a complex system like wave breaking processes, cannot normally conserve mass and therefore inadequately predicts the free surface dynamics. In this study, a hybrid LS-VOF method was used in order to further preserve mass conservation, which results in the level set equation ∂ψ + ∇· (ψHf v) = 0 . ∂t [3] Farthing, M. W. and C. E. Kees, 2008: Implementation of discontinuous Galerkin methods for the level set equation on unstructured meshes. Tech. Rep. ERDC/CHL CHETN-XIII-2. [7] Kobayashi, N., L. Meigs, T. Ota, and J. Melby, 2007: Wave transmission over submerged porous breakwaters. ASCE Journal of Waterway, Port, Coastal and Ocean Engineering, 133 (2), 104–116. [8] Longo, S., M. Petti, and I. J. Losada, 2002: Turbulence in the surf and swash zones: A review. Coastal Engineering, 45, 129–147. [9] Mattis, S. A., C. N. Dawson, C. E. Kees, and M. W. Farthing, 2012: Numerical modeling of drag for flow through vegetated domains and porous structures. Advances in Water Resources, 39, 44–59. [10] Osher, S. and R. Fedkiw, 2001: Level set methods: An overview and some recent results. Journal of Computational Physics, 169, 463–502. [11] Sethian, J. A., 1999: Level Set Methods and Fast Marching Methods: Evolving Interfaces in Computational Geometry, Fluid Mechanics, Computer Vision, and Materials Science. Cambridge University Press, New York. (21) A detailed formulation of the equations and numerical methods is given by Kees et al. [5]. [12] Sethian, J. A., 2001: Evolution, implementation, and application of level set and fast marching methods for advancing fronts. Journal of Computational Physics, 169, 503–355. [13] Xie, S., 1981: Scouring pattern in front of vertical breakwaters and their influence on the stability of the foundation of the breakwaters. M.S. thesis, Civil Engineering, TU Delft.