Cover Image

PAPERBACK
$475.00



View/Hide Left Panel
Click for next page ( 763


The National Academies | 500 Fifth St. N.W. | Washington, D.C. 20001
Copyright © National Academy of Sciences. All rights reserved.
Terms of Use and Privacy Statement



Below are the first 10 and last 10 pages of uncorrected machine-read text (when available) of this chapter, followed by the top 30 algorithmically extracted key phrases from the chapter as a whole.
Intended to provide our own search engines and external engines with highly rich, chapter-representative searchable text on the opening pages of each chapter. Because it is UNCORRECTED material, please consider the following text as a useful but insufficient proxy for the authoritative book pages.

Do not use for reproduction, copying, pasting, or reading; exclusively for search engines.

OCR for page 762
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line THE NUMERICAL SIMULATION OF SHIP WAVES USING CARTESIAN GRID METHODS 762 The Numerical Simulation of Ship Waves Using Cartesian Grid Methods M.Sussman (Florida State University, USA) D.Dommermuth (Science Applications International Corporation, USA) ABSTRACT Two different cartesian-grid methods are used to simulate the flow around the DDG 5415. The first technique uses a “coupled level-set and volume-of-fluid” (CLS) technique to model the free-surface interface. The no-flux boundary condition on the hull is imposed using a finite-volume technique. The second technique uses a level-set technique (LS) to model the free-surface interface. A body-force technique is used to impose the hull boundary condition. The predictions of both numerical techniques are compared to whisker-probe measurements of the DDG 5415. The level-set technique is also used to investigate the breakup of a two-dimensional spray sheet. 1 INTRODUCTION At moderate to high speed, the turbulent flow along the hull of a ship and behind the stern is characterized by complex physical processes which involve breaking waves, air entrainment, free-surface turbulence, and the formation of spray. Traditional numerical approaches to these problems, which use boundary-fitted grids, are difficult and time- consuming to implement. Also, as waves steepen, boundary-fitted grids will break down unless ad hoc treatments are implemented to prevent the waves from getting too steep. At the very least, a bridge is required between potential-flow methods, which model limited physics, and more complex boundary-fitted grid methods, which incorporate more physics, albeit with great effort and with limitations on the wave steepness. Cartesian-grid methods are a natural choice because they allow more complex physics than potential-flow methods, and, unlike boundary-fitted methods, cartesian-grid methods require minimal effort with no limitation on the wave steepness. Although cartesian-grid methods (CGM) are presently incapable of resolving the hull boundary-layer, CGM can model wave breaking, free-surface turbulence, air entrainment, spray-sheet formation, and complex interactions between the ship hull and the free surface, such as transom- stern flows and tumblehome bows. The cartesian-grid methods that are described in this paper use the panelized geometry that is used by potential-flow methods to automatically construct a representation of the hull. The hull representation is then immersed inside a cartesian grid that used to track the interface. No additional gridding beyond what is already used by potential-flow methods is required. We note that another variation of this approach is to use cartesian-grid methods to track the free-surface interface and body-fitted grids to model the ship hull. For the calculation of ship waves, VOF and level-set methods have certain advantages and disadvantages. VOF uses the volume fraction (F) to track the interface. F=0 corresponds to gas and F=1 corresponds to liquid. For intermediate values, between zero and one, there exists an interface between the gas and the liquid. The interface between the gas and the liquid is sharp for a pure VOF method. Level-set methods use a level-set function to model the gas-liquid interface. By definition, denotes gas, denotes liquid, and is the interface. For conventional level-set schemes, the interface between the gas and the liquid is given a finite thickness[17], which is unlike conventional VOF schemes[3]. In the case of free-surface flows, where the density ratio between air and water is almost three orders the authoritative version for attribution.

OCR for page 762
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line THE NUMERICAL SIMULATION OF SHIP WAVES USING CARTESIAN GRID METHODS 763 of magnitude, the finite thickness of the interface that characterizes level-set methods has two advantages over VOF. First, the finite thickness tends to smooth jumps in the tangential component of the velocity on the interface. Second, the finite thickness tends to facilitate using multigrid methods to solve various types of elliptic equations that involve the density. The advection algorithm that is used for VOF conserves mass if the flow field is solenoidal. The level-set advection equation tends to accumulate numerical errors. For the level-set method, the level-set function must be periodically reinitialized to maintain a proper thickness for the interface, otherwise the interface would become either too thick or too thin. The reinitialization process is a significant source of errors in the level-set method. Based on accuracy considerations, the calculation of gravity-driven flows tends to favor VOF over level-set methods. The interface is reconstructed from the volume fractions in VOF. During the reconstruction process, the interface normals and curvature are calculated. Typically, the calculation of the interface normal and curvature are less accurate for VOF than for level-set methods. The interface normals and curvature are calculated directly in level-set methods in terms of gradients of the level-set function. As a result, the calculation of the normals and curvature are less costly for level-set methods relative to VOF. The calculation of surface tension effects, which are a function of the curvature of the interface, tends to favor level-set methods over VOF due to considerations of accuracy and efficiency. On highly-stretched, multidimensional grids, VOF methods are less prone to aliasing errors than level-set methods. Level-set methods incur errors as the interface rotates through highly resolved regions into regions that are not resolved well. This type of aliasing error occurs in cartesian-grid methods when the mesh along one coordinate axis is more finely resolved than along another coordinate axis. By definition, the level-set and volume-of-fluid function both allow mixing of gas and liquid. This feature of level- set and volume-of-fluid methods maybe desirable for modeling gas entrainment such as the air that is entrained by a breaking wave. During the reinitialization process, level-set methods and “coupled level set and volume-of-fluid methods” (CLS) use a signed distance function to update the level-set function and the thickness of the interface. Naturally, the distance function could be used to model the intensity of turbulence and amount of gas entrainment as a function of the distance to the interface. Dommermuth, et al., (1998) used a stratified flow formulation to simulate breaking bow waves on the DDG 5415 at a Froude number Fr=0.41. Their numerical results compared well to whisker-probe measurements in the bow region [8]. However, Dommermuth, et al., (1998) identified two issues that required further study. First, their stratified flow formulation allowed the free-surface interface to become too diffuse. Second, the contact-line treatment didnot allow the free surface to rise and fall cleanly along the side of the hull. The two new numerical approaches that are discussed in this paper are attempts to remedy these problems. Both numerical approaches use a signed distance function to represent the hull. The distance of a point to the hull is negative inside the hull and positive outside the hull. The finite-volume approach uses the signed distance to calculate the area and volume fractions for computational cells cut by the hull, whereas the body-force technique uses the signed distance to prescribe a smooth forcing term. The coupled interface-tracking algorithm (CLS) uses level-set to calculate the normals (and curvature if needed) to the free-surface interface that are used in VOF. The advection portion of the algorithm is performed by VOF [16]. The level-set interface-tracking algorithm uses a new isosurface scheme to calculate the zero level-set. Then the minimal distance between the cartesian points and the zero level-set is calculated in a narrow band. The minimal distance is made positive in the water and negative in the air. This signed distance to the free surface is used to reinitialize the thickness of the interface. The two numerical approaches are used to simulate the flow around the DDG 5415. The CLS technique is still under development, so only preliminary results are presented. The level-set technique includes upgrades to the numerical technique that is described in [8]. Those upgrades include a new body-force formulation that is mollified, a new reinitialization procedure, and a new finite-volume treatment of the convective terms. The original numerical procedure is not mollified and does not use reinitialization. In addition, the original central-difference formulation of the convective terms is not as robust as the new treatment using a flux integral formulation. We first review the governing equations and then we discuss the numerical approaches. Finally, we present some preliminary numerical results which illustrate various features of the numerical algorithms. The application of level-set methods to the breakup of spray sheets is also illustrated. 2 FIELD EQUATIONS As in Dommermuth, et al., (1998), consider turbulent flow at the interface between air and water [8]. Let ui denote the three-dimensional velocity field as a function of space (xi) and time (t). For an incompressible flow, the authoritative version for attribution.

OCR for page 762
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line THE NUMERICAL SIMULATION OF SHIP WAVES USING CARTESIAN GRID METHODS 764 the conservation of mass gives (1) ui and xi are normalized by Uo and Lo, which are the characteristic velocity and length scales of the body, respectively. On the surface of the moving body (Sb), the fluid particles move with the body: (2) where Ui is the velocity of the body. Let Vℓ and Vg respectively denote the liquid (water) and gas (air) volumes. Following a procedure that is similar to for x ∈ Vg and for x ∈ Vℓ. The fluid interface [13, 15], we let denote a level-set function. By definition, corresponds to The convection of is expressed as follows: (3) where d/dt=∂/∂t+ui∂/∂xi is a substantial derivative. Q is a sub-grid-scale flux which can model the entrainment of gas into the liquid. Details are provided in [8]. Let ρℓ and µ ℓ respectively denote the density and dynamic viscosity of water. Similarly, ρg and µ g are the corresponding properties of air. The flow in the water and air is governed by the Navier-Stokes equations: (4) where Re=ρℓUoLo/µℓ is the Reynolds number, is the Froude number, and is the Weber number. g is the acceleration of gravity, and σ is the surface tension. Fi is a body force that is used to impose boundary conditions on the surface of the body. P is the pressure. Ti accounts for surface-tension effects. δij is the Kronecker delta symbol. As described in [8], τij is the subgrid-scale stress tensor. Sij is the deformation tensor: (5) ρ and µ are respectively the dimensionless variable densities and viscosities: (6) where λ=ρg/ρℓ and η=µg/µℓ are the density and viscosity ratios between air and water. For a sharp interface, with no mixing of air and water, H is a step function. In practice, a mollified step function is used to provide a smooth transition between air and water. Based on [3, 4], the effects of surface tension are expressed as a singular source term in the Navier-Stokes equations: (7) where is the curvature of the air-water interface expressed in terms of the level-set function: (8) The pressure is reformulated to absorb the hydrostatic term: (9) where Pd is the dynamic pressure and Ph is a hydrostatic pressure term: (10) As discussed in [8], the divergence of the momentum equations (4) in combination with the conservation of mass (1) provides a Poisson equation for the dynamic pressure: (11) where Σ is a source term. Equation 11 is used to project the velocity onto a solenoidal field. 3 ENFORCEMENT OF BODY BOUNDARY CONDITIONS Two different cartesian-grid methods are used to simulate the flow around the DDG 5415. The first technique imposes the no-flux boundary condition on the body using a finite-volume technique. The second technique the authoritative version for attribution.

OCR for page 762
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line THE NUMERICAL SIMULATION OF SHIP WAVES USING CARTESIAN GRID METHODS 765 imposes the no-flux boundary condition via an external force field. Both techniques use a signed distance function ψ to represent the body. ψ is positive outside the body and negative inside the body. The magnitude of ψ is the minimal distance between the position of ψ and the surface of the body. With respect to the volume of fluid that is enclosed by the body (Vb), we define a function : (12) The function can be expressed in terms of a surface distribution of normal dipoles [10]. (13) where n is the outward-pointing unit normal to the body, and R is a Rankine source, R=|x−x′|. ψ is expressed in terms of as follows, (14) where |x−x′|min is the minimal distance between the field point x and the points on the body x′. In practice, the body is discretized using triangular panels. As a result, the calculation of the minimal distance sweeps over all the triangles comprising the body and must account for the possibility that the minimal distance may occur either at the corners of triangle, along the edges of triangle, or inside the triangle. 3.1 Free-slip conditions In the finite volume approach, the irregular boundary (i.e. ship hull) is represented in terms of ψ along with the corresponding area fractions A and volume fractions V. V=1 for computational elements fully outside the body and V=0 for computational elements fully inside the body. The representation of irregular boundaries via area fractions and volume fractions has been used previously in the following work for incompressible flows [1, 19, 5]. Recall the pressure equation, (15) with the following no-flow boundary condition: (16) where nwall is the outward normal drawn from the active flow region into the geometry region. For each discrete computational element Ωi,j,k we define the geometry volume fraction V and area fraction A as Γi+1/2,j,k represents the left face of a computational element; similar definitions apply to Γi−1/2,j,k, Γi,j+1/2,k,…. In order to discretely enforce the boundary conditions (16) at the geometry surface, we use a finite volume approach for discretizing (15). Given an irregular computational element Ωijk (see Figure 1), we have The divergence theorem motivates the following second order approximation of the divergence ∇ · U at the centroid of Ωijk: (17) In terms of geometry volume fractions Vijk and area fractions Ai+1/2,j,k, (17) becomes, (18) the authoritative version for attribution. For a zero flux boundary condition at the wall, the last term in (18), is zero. The finite volume approach, when applied to the divergence operator in (15) becomes:

OCR for page 762
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line THE NUMERICAL SIMULATION OF SHIP WAVES USING CARTESIAN GRID METHODS 766 and Due to the no flow condition (16), the terms and cancel each other. The resulting discretization for p is: where, for example, (px)i+1/2,j,k is discretized as 3.2 No-slip conditions The boundary condition on the body can also be imposed using an external force field. Based on Dommermuth, et al., (1998), the distance function representation of the body (ψ) is used to construct a body force as follows: Figure 1: Diagram of computational element (i, j) that is cut by the embedded geometry. (19) where cf is a friction coefficient. ∆ is used to mollify the body force such that it is gradually applied across the surface of the body. Recall that ψ(x)≤0 corresponds to points within the body. Fi=0 outside of the body. is an adjustment function: (20) To is the adjustment time. The adjustment function smoothly increases to unity from its initial value of zero. The effect of the adjustment function is described in [8]. The adjustment function reduces the generation of nonphysical high- frequency waves. As constructed, the velocities of the points within the body are forced to zero. For a body that is fixed in a free stream, this corresponds to imposing no-slip boundary conditions. 4 INTERFACE TRACKING Two methods are presented in our work for computing ship flows. Both methods use a “front-capturing” type procedure for representing the free surface separating the air and water. The first technique is based on the Coupled the authoritative version for attribution. volume-of-fluid and level set method (CLS) and the second technique is based on the level set method (LS) alone.

OCR for page 762
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line THE NUMERICAL SIMULATION OF SHIP WAVES USING CARTESIAN GRID METHODS 767 4.1 CLS method In this section, we describe the 2d coupled Level Set and Volume of Fluid (CLS) algorithm for representing the free surface. For more details, e.g. axisymmetric and 3d implementations, see [16]. In the CLS algorithm, the position of the interface is updated through the level set equation (level set function denoted by ) and volume of fluid equation (volume fraction of liquid within each cell is denoted by Fij), In order to implement the CLS algorithm, we are given a discretely divergence free velocity field uMAC defined on the cell faces (MAC grid), (21) UMAC, Given and we use a “coupled” second order conservative operator split advection scheme in order to find and The 2d operator split algorithm for a general scalar s follows as (22) (23) where Gi+1/2,j=si+1/2,jui+1/2,j denotes the flux of s across the right edge of the (i, j) th cell and denotes the flux across the top edge of the (i, j)th cell. The operations (22) and (23) represent the case when one has the “x-sweep” followed by the “y-sweep”. After every time step the order is reversed; “y-sweep” (done implicitly) followed by the “x-sweep” (done explicitly). The scalar flux si+1/2,j is computed differently depending on whether s represents the level set function or the volume fraction F. For the case when s represents the level set function we have the following representation for si+1/2,j (ui+1/2,j>0), where The above discretization is motivated by the second order predictor corrector method described in [2] and the references therein. For the case when s represents the volume fraction F we have the following representation for si+1/2,j (ui+1/2,j>0), (24) where The integral in (24) is evaluated by finding the volume cut out of the region of integration by the line represented by the zero level set of The term found in (24) represents the linear reconstruction of the interface in cell (i, j). In other words, has the form (25) A simple choice for the coefficients aij and bij is as follows, (26) (27) The intercept cij is determined so that the line represented by the zero level set of (25) cuts out the same volume in cell (i, j) as specified by In other words, the following equation is solved for cij, the authoritative version for attribution. where and Fn+1 have been updated according to (22) and (23) we “couple” the level set function to the volume After fractions as a part of the level set reinitialization step. The level set reinitialization step replaces the current value of with the exact distance to the VOF reconstructed interface. At the same time, the VOF reconstructed interface uses the current value of to determine the slopes of the piecewise linear reconstructed interface. Remarks:

OCR for page 762
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line THE NUMERICAL SIMULATION OF SHIP WAVES USING CARTESIAN GRID METHODS 768 • The distance is only needed in a tube of K cells wide K=ε/∆x+2, therefore, we can use “brute force” techniques for finding the exact distance. See [16] for details. • During the reinitialization step we truncate the volume fractions to be 0 or 1 if Although we truncate the volume fractions, we still observe that mass is conserved to within a fraction of a percent for our test problems. 4.1.1 CLS Contact angle boundary conditions in general geometries The CLS contact angle boundary conditions are enforced by extending into regions where Vij<1 (i.e. initializing “ghost” values of in the inactive portion of the computational domain). The contact angle boundary condition at solid walls is given by (28) where θ is a user defined contact angle and nwall is the outward normal drawn from the active flow region into the geometry region. In terms of (the free surface level set function) and ψ (the geometry level set function), (28) becomes In figure 2, we show a diagram of how the contact angle θ is defined in terms of how the free surface intersects the geometry surface. Figure 2: Diagram of gas/liquid interface meeting at the solid. The dashed line represents the imaginary interface created thru the level-set extension procedure. The “extension” equation has the form of an advection equation: (29) In regions where is left unchanged. For a 90 degree contact angle (the default for our computations), we have In other words, information propagates normal to the geometry surface. For contact angles different from 90 degrees, the following procedure is taken to find uextend: Remarks: • In 3d, the contact line (CL) is the 2d curve which represents the intersection of the free surface with the geometry surface (ship hull). The vector n2 is orthogonal to the contact line (CL) and lies in the tangent plane of the geometry surface. and ψ are defined within a narrow band of the zero level set of we can also define uextend within • Since both a narrow band of the free surface. • We use a first order upwind procedure for solving (29). The direction of upwinding is determined from the extension velocity uextend. We solve (29) for τ=0…ε. the authoritative version for attribution. • For viscous flows, there is a conflict between the no-slip condition and the idea of a moving contact line. See [6, 8, 9, 12] and the references therein for a discussion of this issue. We have performed numerical studies for axisymmetric oil spreading in water under ice [18] with good agreement with experiments. In the future, we wish to experiment with appropriate slip-boundary conditions near the contact line. 4.2 Level-set method A key part of level-set methods is reinitialization. Without reinitialization, the thickness of the interface between the gas and the liquid can get either too thick or

OCR for page 762
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line THE NUMERICAL SIMULATION OF SHIP WAVES USING CARTESIAN GRID METHODS 769 too thin. Reinitialization is based on the construction of a signed distance function that represents the distance of points from the gas-liquid interface. By definition, the signed distance is positive in the liquid and negative in the gas. At the interface, the distance function is zero. A variety of methods have been utilized for calculating the signed distance function, including a hyperbolic equation [17] and direct methods [16]. The hyperbolic equation methods tend to be less accurate but more efficient than direct methods. Here, we outline a direct method that can be efficiently implemented on parallel computers with second-order accuracy. The numerical scheme can also be generalized to higher order. First, calculate the intersection points (xp) where the zero level-set crosses each of the cartesian axes. At these intersection points calculate the normal to the interface (np). Together, xp and np determine local approximations to the planes that pass through the zero level-set. For points that are within a narrow band of these planes, calculate the minimal distance to the planes. Once the minimal distance is calculated, assign the sign of the distance function based on the sign of the level-set function. For example, consider a zero crossing along the z-axis. Locally, near the zero crossing, the level-set function is fitted with Lagrange polynomials. (30) where zo is offset where interpolated level-set function Lk are Lagrange polynomials and are discrete values of near the zero level-set along the z-axis. K−1 is the degree of the interpolating polynomial. zo is calculated directly for low-order polynomials and iteratively for high-order polynomials. Let xo=(xo, yo, zo), where (xo, yo, zo) is the coordinate of the zero crossing. The unit normal no at the zero crossing is calculated in terms of the level-set function: (31) where the gradient terms are calculated using finite difference formulas of desired order. The minimal distance (s) between a point (xp) and a plane lies along the unit normal to the plane. Denote the position where the point intersection occurs as xs, then (32) where s is expressed in terms of a dot product: (33) Note that higher-order corrections involve curvature terms, etc. As long as |xo−xs|≤∆g, where ∆g is the grid size, then s is potentially the minimal distance to the zero level-set. Other candidates include planes in the neighborhood of xp. On a structured grid, shifts along the cartesian axes can be performed to consider other candidates. Only xp near the zero level- set are required in the reinitialization procedure. A simple procedure for finding points near the zero level-set involves weighted averages. First construct a stair-case approximation (Φ) to the zero level-set: (34) A weighted average along the k-th indice is (35) Similar expressions hold along the i-th and j-th indices. Repeated applications of weighted averages provide a narrow band that encompasses the zero level-set. The narrow band corresponds to the region The signed distance function D is expressed in terms of the level-set function and the minimal distance: (36) Based on [17], is reinitialized as follows: (37) where ∆ is the desired thickness of the interface. 5 FLUX INTEGRAL METHODS We define the temporal and spatial averaging over a time step and a cell as follows: (38) where here, the tilde and overbar symbols respectively denote temporal and spatial averaging. ∆t is the time step, and ∆V is the volume of the cell. the authoritative version for attribution.

OCR for page 762
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line THE NUMERICAL SIMULATION OF SHIP WAVES USING CARTESIAN GRID METHODS 770 As an example, consider the application of the preceding operator to the level-set equation (3): (39) where here superscript n denotes the time level. We focus our attention on the convective term. The convective term accounts for the flux of the level-set function across the faces of the control volume. A second-order approximation for the flux across one face of a cell is provided below: (40) where is the flux across the positive face along the x axis. The limits of integration are provided below: (41) where ∆x, ∆y, and ∆z are the lengths of the cell along the cartesian axes. u+ is the normal component of fluid velocity at the center of positive face along the x-axis. v+ and v− are the normal components of the fluid velocities at the centers of the positive and negative faces along the y-axis. Similar definitions hold for w+ and w−. In a mapped coordinate system, the expression for the flux is (42) where J is the jacobian, and x, y, and z are functions of r, s, and t: (43) For this particular approximation, the jacobian is (44) On any one face the stencil associated with the Lagrangian interpolation of is 3×3×3=27 points, but for the entire cell, the stencil is 5×5×5=125 points. We use a upwind-biased stencil for the momentum equations and a symmetric stencil for the level-set function. The diagonal and cross terms in the momentum equations are treated the same. Generally, we use eight-point Gaussian quadrature to evaluate the flux over each face. Details of the numerical algorithm are described in [7]. Various types of limiters are described in [11]. 6 PRELIMINARY RESULTS In section 6.1, we present preliminary computations of flow past a DDG 5415 ship. In section 6.2, we present preliminary computations of the breakup of a two-dimensional spray sheet. 6.1 Ship Wave Results As a demonstration of the level-set and the coupled level-set and volume-of-fluid formulations, we predict the free- surface disturbance near the bow of the DDG 5415 moving with forward speed. The experiments were performed at the David Taylor Model Basin (DTMB), and are available via the world wide web at http://www50.dt.navy.mil/5415/. This is the same flow that Dommermuth, et al., (1998) originally investigated using their stratified flow formulation [8]. As before, we only consider the high speed case. For this case, a plunging breaker forms near the bow. Air is entrained and the authoritative version for attribution. splash up occurs where the wave reenters the free surface. There is flow separation at the stern, and the transom is dry. A large rooster tail forms just behind the stern. Based on the speed (Uo=6.02Knots) and the length (Lo=5.72m) of the model, the Reynolds and Froude numbers are Re=1.8×107 and Fr2=0.41. The effects of surface tension are not included. The density ratio of air and water is λ=0.0012 and the ratio of the dynamic viscosities is η=0.018.

OCR for page 762
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line THE NUMERICAL SIMULATION OF SHIP WAVES USING CARTESIAN GRID METHODS 771 In regard to the numerical parameters for the level-set formulation, we use a friction coefficient cf=500 in the body- force term (19). The adjustment time is To=0.02. For the level-set formulation, the length and width of the computational domain are L=2.5 and W=1.50. The height of the air above the mean free-surface is h=0.15 and the depth below the mean free-surface is d=1.0. One grid resolution is used with 512×128×129 grid points. Three different levels of grid stretching are used along the y- and z-axes. For the highest grid resolution, the smallest grid spacing is 2.6×10−3 along the y-axis and 3.6×10−4 along the z-axis. For the medium resolution simulation, the smallest grid spacing is 3.8×10−3 along the y-axis and 1.8×10−3 along the z-axis. For the coarsest grid simulation, the smallest grid spacing is 3.8×10−3 along the y-axis and 3.5×10−3 along the z-axis. The grid spacing (4.9×10−3) is constant along the x-axis for all three cases. The thicknesses of the free-surface interfaces for the fine, medium, and coarse simulation are respectively ∆=0.05, 0.025, and 0.0125. The durations of the coarse and medium resolution simulations are t=0.76 and t=0.68, respectively. No special treatment is used for the level-set function inside the ship. These durations correspond to about three quarters of a ship length based on the present normalization. For these durations, the flow is steady near the bow and still evolving near the stern. (The fine resolution simulation is still evolving, and it is not possible at this time to present complete results. More complete results will be provided at the symposium and in the discussion section of this paper1.) The ship is centered in the computational domain with the same fixed sinkage and trim as used in the experiments. In order to construct the body force term, the hull is panelized using approximately 4000 panels. Coarse and medium resolution simulations have been performed using the CLS formulation. The coarse simulation uses 256×64×64 grid points, and the fine resolution uses 512×128×128 grid points. The length, width, and height of the computational domain are L=2, W=0.5, and H=0.5, respectively. The water depth is d=0.25. The grid spacing is constant along all three cartesian axes. In the next phase of our research, we will implement grid stretching, which will allow greater water depths to be simulated. The durations of the CLS simulations are t=0.75. Unlike the level-set results, the CLS results extend the free-surface interface into the hull using the techniques outlined earlier in our paper. The free-surface elevation was measured at DTMB using a whisker probe. Twenty-one transverse cuts were performed near the bow, extending from x=0 to x=0.178 in dimensionless units. The whisker probe measures the highest point of the free surface. In regions where there is wave breaking, the whisker probe measures the top of the breaking wave. Seventeen transverse cuts were performed in the stern, extending from x=1.01 to x=1.22. Figures 3 and 4 compare measurements at the bow and stern to the numerical predictions. The bow measurements include profile and whisker-probe measurements. Comparisons to the bow data are performed at four stations: x=0.0444, x=0.0622, x=0.0800, and x=0.0978. The circular symbol denotes profile measurements. The solid black lines denote the outline of the hull and the whisker-probe measurements. The solid blue line is medium CLS and the dashed blue line is coarse CLS. The solid red line is medium level-set and the dashed red line is coarse level-set. In general, the CLS technique captures the rapid rise up the side of the hull. The level-set technique does less well in this regard. In the outer- flow region the CLS coarse results are slightly better than the CLS fine results. This may be attributed to the shallow depth that is used in the CLS. The level-set results appear to converge better in the outer-flow region, but the results of the fine simulation are required for confirmation. Figure 4 shows the entire flow around the ship for the medium resolution level-set simulation. The stern whisker- probe measurements are overlaid for the purposes of comparison. Although the numerical results are not stationary, the shape of the stern contours show general agreement with laboratory measurements. However, the amplitude of the numerical results are significantly lower than the measurements. Note that the stern is partially dry in the numerical simulations. The outline of the hull is visible in the numerical simulations because the level-set function intersects the hull. 6.2 Spray Sheet Results The Navier-Stokes equations in combination with a level-set formulation are used to study the breakup of two- dimensional sheet of water. The sheet is lo=6mm thick. The length of the sheet is 24mm. The top and bottom of the sheet are bounded by air. The initial mean-velocity of the water is uo=3m/s. The initial rms turbulent velocity of the water is The air is initially quiescent. Based on the sheet thickness (lo) and the mean velocity (uo), the Reynolds number is Re=uolo/µ=18,000 and the Weber number is where µ is the kinematic viscosity of water, ρ is the water density, and σ is the surface tension. The density and viscosity ratios the authoritative version for attribution. 1We would have performed longer simulations, but the NAVO T3E was unexpectedly shutdown for five days of maintenance just before this paper was due.

OCR for page 762
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line THE NUMERICAL SIMULATION OF SHIP WAVES USING CARTESIAN GRID METHODS 772 are λ=ρg/ρℓ=0.0012 and η=µg/µℓ=0.018, which are appropriate for air-water interfaces. This parameter regime roughly corresponds to experiments that were performed by Sarpkaya and Merrill (1998), [14]. Numerical convergence is established using 20482 and 40962 grid points. Second-order accuracy in space is established. A third-order Runge-Kutta scheme is used to integrate the system of equations with respect to time. Mass is conserved to within 0.25% throughout the entire calculation. Figure 5 illustrates the evolution of a two-dimensional spray sheet. The black contour lines indicate the interface between air and water. The water sheet is bounded by air both at the top and the bottom of the sheet. The color contours denote the vorticity. The flow is turbulent within the water sheet and laminar in the air. The mean velocity and rms velocity profiles are initially top-hat functions. The flow is moving from left to right. The turbulent fluctuations in the water are initially immersed below the top of the sheet and above the bottom of the sheet (see Fig. 5: t=0). The turbulence in the water diffuses and interacts with the interfaces (see Fig. 5: t=2.5). The initial interaction is a roughening of the air-water interface. A thin boundary layer forms in the air. The boundary layer is colored blue (negative) at the top of the sheet and colored red (positive) at the bottom of the sheet. As the interface gets rougher and ligaments begin to form, the air separates from the back of the ligaments. The boundary layer thickens, and air is dragged along the top and the bottoms of the sheet. Primary vortex shedding initially occurs behind the ligaments (see lower left of sheet in Fig. 5: t=5). As the primary vortices are shed, their interactions lead to the formation of secondary and tertiary vorticity (see upper middle of sheet in Fig. 5: t=7.5). Vortices are periodically shed from the backs of ligaments (see lower middle of sheet in Fig. 5: t=10). There is evidence of vortex merging both in the air and in the water (see upper left of Fig. 5: t=17.5). Although there is significant flow separation in the air, there is little or no separation in the water. The largest ligaments are formed by eddies impinging on the interface (see upper left of Fig. 5: t=12.5). Cavities form in regions where primary vortices are trapped. The inlets to the cavities shed secondary vorticity, which tends to make the cavities even larger (see middle of sheet in Fig. 5: t=15). At the inlets to the cavities, vortex pairs are formed. Under their own self-induced velocities, the vortex pairs move into the cavities where they diffuse. Note that droplets do not actually form at the tips of the ligaments because 2d flows are not subject to the same instabilities as 3d flows. The turbulent kinetic energy tends to concentrate in the thicker portions of the deformed spray sheet. The flow within the ligaments is relatively benign. In agreement with theory, the pressure at the tips of the longest ligaments roughly scales like P=(Wer)−1 where r is the radius of curvature of the tip. 7 CONCLUSION In this paper, we have outlined the key numerical algorithms for simulating free-surface flows on cartesian grids using level-set and coupled level-set and volume-of-fluid techniques. Preliminary numerical results have been shown for ship waves and spray sheets. The ship wave results indicate that cartesian-grid methods are capable of resolving the flow around a ship if the grid resolution is sufficient. Near the bow and stern, we estimate that the grid spacing along all three cartesian axes should be ∆=0.0005 (based on ship length) in order to resolve breaking waves. On a parallel computer, it is possible to approach this level of grid resolution, but adaptive gridding may also be required to fully resolve the entire flow around a ship [15]. Alternatively, cartesian-grid methods could be embedded in more conventional boundary-fitted methods to capture complex flows near the bow or stern. The spray-sheet results show that cartesian-grid methods are capable of resolving the air and water boundary layer at realistic Reynolds numbers. Acknowledgments. The first author is supported in part by NSF Division of Mathematical Sciences under award number DMS 9996349. The second author is supported by ONR under contract number N00014–97-C-0345. Dr. Edwin P.Rood is the program manager. The numerical simulations have been performed on the T3E computer at the Naval Oceanographic Office using funding provided by a Department of Defense Challenge Project. We are very grateful to Mr. George Innis, Dr. James Rottman, and Mr. Andrew Talcott for assistance with this paper. REFERENCES [1] A.S.Almgren, J.B.Bell, P.Colella, and T.Marthaler. A cartesian grid projection method for the incompressible euler equations in complex geometries. SIAM J. Sci. Comput., 18(5):1289–1309, 1997. [2] J.B.Bell, P.Colella, and H.M.Glaz. A second-order projection method for the incompressible Navier-Stokes equations. J. Comput. Phys., 85:257– 283, December 1989. the authoritative version for attribution.

OCR for page 762
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line THE NUMERICAL SIMULATION OF SHIP WAVES USING CARTESIAN GRID METHODS 773 [3] J.U.Brackbill, D.B.Kothe, and C.Zemach. A continuum method for modeling surface tension. J. Comput. Phys., 100:335–353, 1992. [4] Y.C.Chang, T.Y.Hou, B.Merriman, and S.Osher. Eulerian capturing methods based on a level set formulation for incompressible fluid interfaces. J. Comput. Phys., 124:449–464, 1996. [5] P.Colella, D.T.Graves, D.Modiano, E.G.Puckett, and M.Sussman. An embedded boundary/volume of fluid method for free surface flows in irregular geometries. In proceedings of the 3rd ASME/JSME joint fluids engineering conference, number FEDSM99–7108, San Francisco, CA, 1999. [6] R.G.Cox. The dynamics of the spreading of liquids on a solid surface, part 1. viscous flow. J. Fluid Mech., 168:169–194, 1986. [7] D.G.Dommermuth. Numerical Flow Analysis (nfa) working papers. Technical report, Science Applications International Corporation, 2000. [8] D.G.Dommermuth, G.E.Innis, T.Luth, E.A. Novikov, E.Schlageter, and J.C.Talcott. Numerical simulation of bow waves. In proceedings of the Twenty Second Symposium on Naval Hydro., pages 508–521, Washington, D.C., 1998. [9] L.M.Hocking and A.D.Rivers. The spreading of a drop by capillary action. J. Fluid Mech., 121:425–442, 1982. [10] H.Lamb. Hydrodynamics. Dover Publications, New York, 1932. [11] B.J.Leonard. Bounded higher-order upwind multidimensional finite-volume convection-diffusion algorithms. In W.J.Minkowycz and E.M.Sparrow, editors, Advances in Numerical Heat Transfer, volume I, pages 1–58. Taylor & Francis, 1997. [12] C.G.Ngan and E.B.Dussan V. On the dynamics of liquid spreading on solid surfaces. J. Fluid Mech., 209:191–226, 1989. [13] S.Osher and J.A.Sethian. Fronts propagating with curvature-dependent speed: Algorithms based on hamilton-jacobi formulations. J. Comput. Phys., 79(1):12–49, 1988. [14] T.Sarpkaya and C.Merrill. Spray formation at the free surface of liquid wall jets. In proceedings of the Twenty Second Symposium on Naval Hydro., pages 796–808, Washington, D.C., 1998. [15] M.Sussman, A.Almgren, J.Bell, P.Colella, L.Howell, and M.Welcome. An adaptive level set approach for incompressible two-phase flows. J. Comput. Phys., 148:81–124, 1999. [16] M.Sussman and E.G.Puckett. A coupled level set and volume of fluid method for computing 3d and axisymmetric incompressible two-phase flows. J. Comp. Phys. accepted for publication. [17] M.Sussman, P.Smereka, and S.J.Osher. A level set approach for computing solutions to incompressible two-phase flow. J. Comput. Phys., 114:146– 159, 1994. [18] M.Sussman and S.Uto. Computing oil spreading underneath a sheet of ice. Technical Report CAM Report 98–32, University of California, Los Angles, July 1998. [19] H.S.Udaykumar, H.C Kan, W.Shyy, and R.Tran-Son-Tay. Multiphase dynamics in arbitrary geometries on fixed cartesian grids. J. Comput. Phys., 137(2):366–405, 1997. the authoritative version for attribution.

OCR for page 762
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as the authoritative version for attribution. Figure 3: Flow near bow. THE NUMERICAL SIMULATION OF SHIP WAVES USING CARTESIAN GRID METHODS 774

OCR for page 762
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as the authoritative version for attribution. Figure 4: Flow near stern. THE NUMERICAL SIMULATION OF SHIP WAVES USING CARTESIAN GRID METHODS 775

OCR for page 762
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as the authoritative version for attribution. Figure 5:2d spray sheet. THE NUMERICAL SIMULATION OF SHIP WAVES USING CARTESIAN GRID METHODS 776

OCR for page 762
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as the authoritative version for attribution. Figure 5:2d spray sheet continued. THE NUMERICAL SIMULATION OF SHIP WAVES USING CARTESIAN GRID METHODS 777

OCR for page 762
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line THE NUMERICAL SIMULATION OF SHIP WAVES USING CARTESIAN GRID METHODS 778 DISCUSSION U.Bulgarelli Instituto Nazionale per Studi ed Esperienze di Architettura Navale, Italy In your algorithm do you have already adopted the adaptive grid in the 3D geometry? AUTHOR'S REPLY We are in the process of developing a body-fitted method. The method will be described at the next symposium. DISCUSSION K.Hendrickson Massachusetts Institute of Technology, USA The authors of this paper show an aggressive use of a numerical method which has, to date only, been used for smaller engineering problems. The fact that they are attempting it for this type of problem says much about their patience and ambition. The blending of the volume of fluid and level set methods is quite creative and shows promising results. I believe that the Cartesian Grid Method is wonderfully useful in that many of the gridding difficulties have been removed and/or reduced to panel method that has been dealt with in detail in the literature. I wish them both luck as they push the method further. QUESTIONS 1. At this stage in development of the two methods, it seems that the coupled level-set/volume of fluid technique (CLS) is more accurate/robust in treating the hull boundary conditions mainly because it can better define where the hull lies in the Cartesian grid. Have the authors done any investigation on the effects of the mollified body force term used in the level-set (LS) technique? In using comparisons to the waterline measurements as the benchmark, the exact hull position would likely be a critical point. Is it possible that the mollified body force term is smoothing out the hull to the extent that it is affecting the waterline results? Would less mollification or higher resolution in the region near the body produce better LS results? This can almost be inferred from Figure 3 in the paper. 2. The choice of friction coefficient in the body force term (equation 19) seems to be somewhat arbitrary. Have the authors done any type of parametric study on a range of friction coefficients and their effect on the LS results? 3. Considerable effort has been invested in the LS community to address reinitialization, which is also done in this paper. The reinitialization issue comes about because a Lagrangian thought process has been applied to a Eulerian method. In most LS formulations, the advection of the level-set function, ~, is performed using the velocity of the fluid. This causes the LS function to lose its distance function property and require reinitialization. It is possible to construct a velocity field such that the distance function remains one [1]. Have the authors considered this type of LS formulation? 4. What is the computational cost comparison between the CLS and LS methods at the resolutions submitted in the paper? 5. How do the authors feel CGM compare to other less computationally expensive capabilities such as unRANS or 2D+T methods? 6. What do the authors consider to be the major limitations of the CGM, both CLS and LS, in terms of their applicability to Marine Hydrodynamics and Computational Ship Hydrodynamics? REFERENCES 1. Adalsteinsson, D. and Sethian, J.A. “The Fast Construction of Extension Velocities in Level Set Methods,” J. Comp. Physics, Vol. 148, 1999, pp. 222. AUTHOR'S REPLY 1. The coupled level-set/volume (CLS) of fluid technique has a more accurate treatment of the hull boundary condition than the level-set method (LS). Our tests indicate that mollified body-force terms improve convergence. 2. The body-force term is as large as possible without violating the Courant condition. 3. The advantage of our reinitialization procedure is its accuracy, which can be generalized to any order. Other procedures, the authoritative version for attribution.

OCR for page 762
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line THE NUMERICAL SIMULATION OF SHIP WAVES USING CARTESIAN GRID METHODS 779 such as the one proposed by the discussor, are effective away from the interface. 4. The CLS method is about twice as expensive as the LS method. However, the computational costs associated with both methods are less than ten percent of the Poisson solver. 5. Interface tracking methods are capable of capturing physics that unRANS and 2D+T will never be capable of modeling. Although interface tracking methods are more computationally expensive than unRANS and 2D+T, this will become less of an issue as computers become faster. Ten years from today, interface tracking will be the method of choice for modeling breaking waves and the near-field flow around naval combatants. 6. The treatment of the hull boundary condition is not accurate enough. This issue is currently being addressed by using a body-fitted grid with a level-set treatment of the free-surface elevation in the near field of the hull. In the outer- flow region, the inner solution is matched to a combination of spectral methods and panel methods. This matching procedure reduces the number of grid points and the amount time that is required to generate 3D grids. the authoritative version for attribution.