Large-eddy Simulation of the Formation and Evolution of Benthic Ripples Oliver B. Fringer 473 Via Ortega, Room 187 Dept. of Civil and Environmental Engineering Stanford University, Stanford, CA 94305 phone: (650) 725-6878 fax: (650) 725-9720 email: [email protected] Robert L. Street 473 Via Ortega, Room 259 Dept. of Civil and Environmental Engineering Stanford University, Stanford, CA 94305 phone: (650) 723-4969 fax: (650) 725-9720 email: [email protected] Award Number: N00014-05-1-0177 http://ripples.stanford.edu, http://suntans.stanford.edu LONG-TERM GOALS Our goal is to provide a complete simulation code that will represent and predict the sediment transport and bed features on the continental shelf at user-specified resolution by using state-of-the-art algorithms for the physics and numerics of the simulation code. OBJECTIVES Our primary objective is to simulate the ripple climate on the bed of the inner shelf at depths on the order of 20 m and over domains ranging from centimeters to kilometers in support of the Ripples DRI experiments and analyses. Our secondary objective is to use simulation to better understand the physics of ripple formation and sediment transport in this environment. APPROACH In order to simulate the formation and evolution of benthic ripples, we employ a large-eddy simulation code with suspended sediment and moving-bed modules. The suspended-sediment module follows that of Zedler and Street (2001, 2006), and is incorporated into the parallel large-eddy simulation code PCUI (Cui and Street, 2001; 2004) in order to perform high-resolution simulations of sediment-ripple dynamics on parallel computers. We have employed several modifications to the implementation of Zedler and Street, most notably the implementation of the bottom boundary condition for sediment, which does not require knowledge of the near-bed turbulence parameterization. This is outlined in Chou and Fringer (2008a). Our moving bed module responds to the high-resolution shear stress imposed by the fluid at the bed and moves the bed position z according to the formula (Gessler et al., b 1999; Wu et al., 2000) ∂z (1− p') b +∇⋅ qr = D − E , (1) ∂t 1 Report Documentation Page Form Approved OMB No. 0704-0188 Public reporting burden for the collection of information is estimated to average 1 hour per response, including the time for reviewing instructions, searching existing data sources, gathering and maintaining the data needed, and completing and reviewing the collection of information. Send comments regarding this burden estimate or any other aspect of this collection of information, including suggestions for reducing this burden, to Washington Headquarters Services, Directorate for Information Operations and Reports, 1215 Jefferson Davis Highway, Suite 1204, Arlington VA 22202-4302. Respondents should be aware that notwithstanding any other provision of law, no person shall be subject to a penalty for failing to comply with a collection of information if it does not display a currently valid OMB control number. 1. REPORT DATE 3. DATES COVERED 2008 2. REPORT TYPE 00-00-2008 to 00-00-2008 4. TITLE AND SUBTITLE 5a. CONTRACT NUMBER Large-eddy Simulation of the Formation and Evolution of Benthic 5b. GRANT NUMBER Ripples 5c. PROGRAM ELEMENT NUMBER 6. AUTHOR(S) 5d. PROJECT NUMBER 5e. TASK NUMBER 5f. WORK UNIT NUMBER 7. PERFORMING ORGANIZATION NAME(S) AND ADDRESS(ES) 8. PERFORMING ORGANIZATION Stanford University,Dept. of Civil and Environmental Engineering,473 REPORT NUMBER Via Ortega, Room 187,Stanford,CA,94305 9. SPONSORING/MONITORING AGENCY NAME(S) AND ADDRESS(ES) 10. SPONSOR/MONITOR’S ACRONYM(S) 11. SPONSOR/MONITOR’S REPORT NUMBER(S) 12. DISTRIBUTION/AVAILABILITY STATEMENT Approved for public release; distribution unlimited 13. SUPPLEMENTARY NOTES 14. ABSTRACT 15. SUBJECT TERMS 16. SECURITY CLASSIFICATION OF: 17. LIMITATION OF 18. NUMBER 19a. NAME OF ABSTRACT OF PAGES RESPONSIBLE PERSON a. REPORT b. ABSTRACT c. THIS PAGE Same as 7 unclassified unclassified unclassified Report (SAR) Standard Form 298 (Rev. 8-98) Prescribed by ANSI Std Z39-18 r where p' is the porosity of the bed and q is the bed-load transport, which we determine through an empirical formula similar to that of Meyer-Peter and Muller (1948), and depends on the near-bed shear stress. Because we employ terrain-following coordinates, movement of the bed requires implementation of numerical algorithms on moving grids in a manner similar to the method of Hodges and Street (1999), in which the boundary-following grid is updated in response to the moving free surface. WORK COMPLETED We have developed a moving curvilinear grid formulation that enables simulation of the Navier-Stokes equations, along with the transport equations for the suspended sediment, on a grid that moves in response to moving boundaries. The method guarantees consistency with continuity, in that if the suspended sediment concentration is uniform, it remains uniform in response to flow and transport in the absence of erosion or settling. RESULTS We outline the methodology behind the moving grid formulation and present results that demonstrate the potential ramifications of employing methods that are not consistent with continuity (Gross et al., 2002). Details of the formulation can be found in Chou and Fringer (2008b). On a two-dimensional moving, generalized curvilinear-coordinate finite-volume grid, as shown in Figure 1, the equation governing transport of a scalar field C in grid cell i,j is given by d ( ) ( ) ( ) ( ) ( ) J C + U −U C − U −U C + V −V C − V − V C = 0, (2) dt i, j i, j g e e g w w g n n g s s where e, w, n, and s subscripts denote the east, west, north, and south faces of the two-dimensional computational cell with volume J , U and V are the contravariant volume fluxes (in units of m2 s-1) in i,j the transformed curvilinear coordinates, and the g subscript denotes the contravariant volume flux associated with grid movement. The time-rate of change of the cell volume is governed by the continuity equation in the moving curvilinear coordinate system, dJ i, j +(U −U ) −(U −U ) +(V −V ) −(V −V ) = 0 , (3) dt g e g w g n g s where, after assuming continuity due to incompressibility, i.e. U −U +V −V = 0 , (4) e w n s equation (3) leads to the equation governing discrete conservation of space, dJ i, j + U − U + V −V = 0 . (5) dt g ,e g ,w g ,n g,s 2 In order for equation (2) to be consistent with continuity, or CWC, substitution of a uniform scalar concentration field into its time-discrete form must yield the time-discrete form of (5). In general, this is not a straightforward task and special care must be taken when discretizing equations (2) and (5) in order to ensure CWC. A simple example of violation of CWC can be demonstrated by discretizing equation (2) with the second-order Adams-Bashforth method, such that Jn+1C n+1 − Jn C n 3 [( ) ( ) ( ) ( ) ] i , j i , j i , j i , j + Un −Un Cn − Un −Un Cn + Vn −Vn Cn − Vn −Vn Cn Δt 2 g e e g w w g n n g s s , (6) 1 [( ) ( ) ( ) ( ) ] − U n−1 −U n−1 C n−1 − U n−1 −U n−1 Cn−1 + V n−1 −V n−1 Cn−1 − V n−1 −V n−1 C n−1 = 0 2 g e e g w w g n n g s s and discretizing equation (5) with the first-order forward Euler method, such that J n+1 − J n i , j i , j +U n −U n +V n −V n = 0 . (7) Δt g ,e g ,w g ,n g ,s The discretization (6) is not consistent with (7) because, upon substitution of the constant scalar field C=C into equation (6) and assuming continuity of fluid (4), one obtains 0 J n+1 − J n 3 ( ) 1 ( ) i , j i , j + U n −U n +V n −V n − U n−1 −U n−1 +V n−1 −V n−1 = 0 , (8) Δt 2 g ,e g ,w g ,n g ,s 2 g ,e g ,w g ,n g ,s which violates CWC because it is not identical to equation (7). While CWC would be satisfied if both equations (2) and (5) were discretized with the forward Euler method, discretization (7) violates conservation of space unless special care is taken in computing the contravariant volume fluxes associated with grid movement. This is demonstrated in Figure 2, which shows that in one dimension, equation (7) would yield the correct new volume at the new time step. However, the two-dimensional form yields a value for the new volume that is smaller than the true value as would be computed from the cell vertices which are updated independently of equation (7). Chou and Fringer (2008b) developed a technique to avoid this problem by requiring that fluxes associated with grid movement are always calculated at time step n+1/2, rather than time step n. As an example of the effects of not ensuring CWC in a simulation of transport over a moving bottom, we employed a simulation with a bed that moves according to a pre-specified function (rather than due to a dynamic bed model). With an initially uniform concentration field and in the absence of settling or erosion, Figure 3(a) shows that a CWC method guarantees that the concentration remains uniform as the bed oscillates in time. However, Figure 3(b) shows that a non-CWC method leads to exponential growth of the near-bed concentration field because the non-CWC method effectively introduces a source of mass in the near-bed cells whose volume changes considerably due to the grid motion. While this case shows pronounced effects of violating CWC, the effects in a real suspended sediment simulation are not as evident because of the transient nature of the flow in response to settling and erosion. However, as shown in Figure 4, the suspended sediment concentration in a CWC simulation differs significantly from that in a non-CWC simulation, particularly in the near-bed suspended sediment concentration field, which is much higher for the non-CWC case. 3 Figure 1: Depiction of a two-dimensional finite volume grid cell in a generalized curvilinear coordinate system. Figure 2: Example of grid movement when space is conserved in one dimension (a) and when space is not conserved in two dimensions (b) using the discretization (7). Figure (b) shows that the shaded area incorrectly underestimates the actual area of the new cell. 4 Figure 3: The influence of a moving grid on a uniform concentration field for a CWC method (left panels) and a non-CWC method (right panels). T is the period of oscillation of the bed, and the velocity field results only due to bed movement. 5 Figure 4: Suspended sediment concentration over a growing ripple field in a wave-induced current that oscillates with period T. The left column uses a CWC method while the right column does not. REFERENCES Chou, Y.J., and O. B. Fringer. 2008a. Consistent discretization for simulation of flows with moving generalized curvilinear coordinates, Int. J. Numer. Methods Fluids, submitted. Chou, Y.J., and O. B. Fringer, 2008b, Modeling dilute sediment suspension using LES with a dynamic mixed model, Physics of Fluids, in press. Cui, A., and R. L. Street, 2001, Large-eddy simulation of turbulent rotation convective flow development. J. Fluid Mech., 447, 53–84. 6 Cui, A., and R. L. Street, 2004, Large-eddy simulation of coastal upwelling flow. Environmental Fluid Mechanics, 4, 197–223. Gessler, D., Hall, B., Spasojevic, M., Holly, F., Pourtaheri, H., and N. Raphelt, 1999, Application of 3D mobile bed, hydrodynamic model. J. Hydr. Engrg., 125(7), 737-749. Gross, E. S., Bonaventura, L., and G. Rosatti, 2002, Consistency with continuity in conservative advection schemes for free-surface models, Int. J. Numer. Meth. Fluids, 38 (4), 307-327. Hodges, B. R., and R. L. Street, 1999, On simulation of turbulent nonlinear free-surface flows. J. Comput. Phys., 151, 425–457. Meyer-Peter, E., and R. Muller, 1948, Formulas for bed-load transport. Int. Assoc. Hydraulic Research, 2nd Mtg., Stockholm. Zedler, E. A., and R. L. Street, 2001, Large-eddy simulation of sediment transport: currents over ripples. J. Hydraul. Eng., 127, 442–452. Zedler, E. A., and R. L. Street, 2006, Sediment transport over ripples in oscillatory flow. J. Hydr. Engrg. ASCE, 132(2), 1–14. PUBLICATIONS Chou, Y.J., and O. B. Fringer, 2007, Modeling Sediment Suspension in High Reynolds Number Flow Using Large Eddy Simulation, Proc. Fifth International Symposium on Environmental Hydraulics. Chou, Y.J., and O. B. Fringer, 2008b, Modeling dilute sediment suspension using LES with a dynamic mixed model, Physics of Fluids, in press. HONORS/AWARDS/PRIZES O. B. Fringer, Office of Naval Research Young Investigator Award, 2008. 7