Cover Image

PAPERBACK
$203.25



View/Hide Left Panel
Click for next page ( 728


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 727
A Numerical Study of Three Dimensional Viscous Interactions of Vortices with a Free Surface D. Dommermuth (Science Applications International Corporation, USA) D. Yue (Massachusetts Institute of Technology, USA) ABSTRACT We develop semi-implicit and explicit numerical methods for the direct simulation of the three- dimensional Navier-Stokes equations with a free surface. The efficiency of a novel multigrid flow solver permits the simulation of three-dimensional flows with free surfaces at low Reynolds and Froude numbers. The numerical schemes are used to study vortex rings and tubes interacting with walls and free surfaces. In the case of~vortex rings interacting with a no- slip wall our observations of the formation of sec- ondary and tertiary vortex rings agree qualitatively with experimental measurements. When a free surface is present, the results are sensitive to the Froude number. For sufficiently low Froude num- bers, the free surface behaves like a free-slip wall, which agrees qualitatively with experimental ob- servations of vortex rings interacting with clean free surfaces. At intermediate Froude numbers, the normal incidence of a vortex ring with a clean free surface results in the formation of secondary vortex rings. Numerical studies of vortex tubes interacting with free surfaces show two possible mechanisms for the reconnection of vorticity with a free surface in- cluding primary and secondary vorticity reconnec- tions. One type of primary vorticity reconnection should result in a cusp pattern on the free sur- face and secondary vorticity reconnections should manifest themselves as paired dimples on the free surface. The essential stages of the reconnection of secondary vorticity with the free surface are as follows: generation of helical vortex sheets by the primary vortex tube, stripping of the helical vortex sheets due to self-induced straining flows, attach- ment of the helical vortex sheets to the separated 727 free-surface boundary layer, wrapping of U-shaped vortices around the primary vortex tube, feeding of boundary-layer vorticity into the U-vortices, and reconnection of U-vortices with the free surface. We provide evidence which suggests that the stri- ations that may be observed on the free-surface above a pair of trailing tip vortices are caused by helical vorticity. ~ Introduction In the field of free-surface hydrodynamics few things are as complex as the wave and viscous wakes of ships. Consider the region downstream of a moving ship where the steady wave pattern of Kelvin waves is located. This wave pattern is com- posed of two wave systems which move steadily with the ship: the transverse and diverging waves. The transverse waves are normal to the ship track, and the diverging waves obliquely intersect the ship track. These two wave systems are bounded by the cusp line which occurs at an angle of about +19.5 off the ship track. Most of the steady waves dis- perse linearly, but if the conditions are right, rays of the diverging waves that are inside of the cusp lines may evolve into solitons (Brown, et al, 1988~. In addition to the Kelvin waves, there is also a small sector centered near the ship track that con- tains short unsteady waves. Some possible sources of these short unsteady waves include breaking waves and turbulence generation near the ship. Munik, et al, (1987) hypothesize that Bragg scatter from these short unsteady waves is responsible for the narrow V-shaped features that are observed in Synthetic Aperture Radar (SAR) images of ship wakes in very calm seas. The viscous turbulent wake, which is also located near the ship track,

OCR for page 727
has a thrust component that is associated with the propellerts) and a drag component that is due to the ship hull. Axial vortices that are generated in the bow region (Babe, 1969; Mori, 1984; and Takekuma & Eggers, 1984~) and helical vortices that are shed by the propeller (Kerwin, 1986) are embedded in the turbulent wake. Field observa- tions of ships suggests that bubbles rising through the turbulent wake transport surfactant material to the free-surface (Kaiser, et al, 1986~. These bub- bles appear to be entrained by the bow wave, but other sources include propeller cavitation. Slicks consisting of compacted surfactant material form on each side of the ship track possibly due to the action of vortices moving surface water away from the viscous wake's center. These slicks damp cap- illary and gravity-capillary waves that cross the ship track (Dorrestein, 1951) and may contribute to the persistent centerline feature that is observed in SAR images (Lyden, et al, 1988~. Clearly, ship wakes are very complex, and there are many mech- anisms that affect direct and remote observations of ship wakes. Faced with such a diversity of phenomena, we have chosen to focus on a very specialized portion of the problem. Specifically, we are interested in mod- eling the viscous interaction of vorticity with a free surface. Vortices may generate and modu- late short unsteady surface waves and thus affect the observable features. Two sources of vertical flows include turbulent bow-wave breaking (Babe, 1969) and lifting devices (propellers, rudders, & hydrofoils). In addition, coherent vertical struc- tures may also occur in the viscous turbulent wake (Brown & Roshko, 1974, and Roshlco, 1954~. The three-dimensional interactions of these flow fields with the free surface are very complicated, but it is our hope that simple models will provide much needed insight. Some simple models of the large- scale vertical flows that occur in ship wakes include two-dimensional vortex pairs, ring vortices, and three-dimensional vortex tubes. Even these sim- ple models are so complicated that they challenge our basic understanding of vortex interactions with free surfaces. As an initial step toward understanding the inter- actions of ship-scale vortices with a free surface, we consider the direct simulation of the Navier-Stokes equations with a free surface at a finite Froude number. Although the Reynolds numbers that we are able to simulate so far are considerably be low those observed in reality, our three-dimensional method is capable of modeling physical phenom- ena that are simply not possible using inviscid ap- proaches. A notable example is the reconnection of vortex tubes, which is an inherently viscous phe- nomena. Moreover, extensions to our method may eventually allow us to study the effects of turbu- lence, air entrainment, and surfactants. We can partially classify vertical flows near a free- surface in terms of Reynolds, Froude, and Weber numbers. As our velocity scale, we use the cen- terline velocity of a vortex pair, Uc=2I`/,rs, where is the circulation and ~ is the span. In corre- spondence with this velocity, the length scale is the half span, L=~/2, and for ring vortices the ring ra- dius is used. Based on these scales, the Reynolds number is Re=~/~t,, where ~ is the kinematic vis- cosity; the Froude number is Fr=I`/~g~/2~/2~3/2, where 9 is gravity; and the Weber number is We = 2pI`2/7r2sT, where p is the fluid density and T is the surface tension. Although surface contamina- tion plays a very important role in the interaction of vorticity with free surfaces, this effect is difficult to quantify and model, and has not been consid- ered so far in the present simulations. We note that the hydrodynamics of ship-scale surfactant films are discussed by Peltzer, et al, (1990) in this pro- ceedings. In addition, Hirsa, et alms (1990) exper- imental and numerical studies of two-dimensional vortex pairs, which also appears in this proceed- ings, provides a preliminary model of surfactant films. As preparation for our literature review of vortices interacting with free surfaces, we also discuss the related problems of vortex interactions with no- slip walls and the stability and reconnection of vortex tubes. We conclude our introductory sec- tion with a review of numerical methods for simu- lating the incompressible Navier-Stokes equations. In Sec. 2 we provide the mathematical formula- tion for describing the viscous interactions of three- dimensional vortices with free surfaces, and in Sec. 3 we discuss the numerical implementations of our theory. In Sec. 4 we apply our numerical meth- ods to three problems: 1) ring vortex interactions with a wall, 2) ring vortex interactions with a free surface, and 3) vortex tube interactions with a free surface. In Sec. 5 we discuss the long-term direc- tions and promise of this line of research. 728

OCR for page 727
1.1 Interaction of Vortices with WalIs Harvey & Perry's (1971) studies of a wingtip vor- tex near a wall confirm observations of aircraft vor- tices bouncing off the ground. They argue that the rebound effect may be attributed to unsteady boundary-layer separation and the formation of secondary vortices. Their experimental measure- ments show that motion of a vortex moving paral- lel to a wall may be stopped or even reversed due to viscous and inviscid interactions. Inviscid effects alone are ruled out by Saffman (1978), who proves that a vortex pair normally incident to a wall in an inviscid fluid will monotonically approach the wall and then separate apart due to the velocity induced by themselves and their images. Wall~er's (1978) two-dimensional computations show that explosive boundary-layer growth is possible, but his boundary-layer theory is not capable of model- ing the formation of secondary vortices. Peace & Riley's (1983) two-dimensional Navier-Stokes sim- ulations show substantial growth of the boundary layer, but their Reynolds numbers (Re = 0~10~) are too low for secondary vortex formation to oc- cur. WaLker, et ales (1986) studies of vortex rings nor- mally incident to a wall exhibit many of the same phenomena as two-dimensional vortices. An im- portant difference is that ring vortices may stretch. According to inviscid theory, a ring vortex ap- proaching a wall will expand radially outward and the diameter of its core will decrease. However, the experimental results of Walker, et al, and others (see Cerra & Smith, 1983) show much more com- plexity including boundary-layer separation, sec- ondary and tertiary vortex formation, rebounding, azimuthal instabilities, and transition to turbu- lence. For Re < 250, no secondary vortices are observed and the primary vortex rapidly diffuses while being held fixed against the wall. Flow sep- aration and secondary vortex formation occur for Re > 250. For Re > 3000, Walker, et al, have diffi- culties producing laminar vortex rings. They also perform boundary-layer calculations which once again show explosive growth of the boundary layer, but fail to predict the formation of secondary and tertiary vortices. As a prelude to our study of vortex ring interac- tions with a free-surface, we have chosen to simu- late the interactions of a vortex ring with a wall be- cause Wallcer, et alts (1986) documentation permits qualitative validation of our Navier-Stokes com- puter code. In addition, the comparisons between vortex ring interactions with walls and free sur- faces are very useful. Our numerical simulations for the wall case show the formation of secondary and tertiary vortex rings. In fact, the primary vor- tex ring is wrapped in a sheath of counter-sign vor- ticity. Initially, a sharp gradient, which diffuses with time, exists between the primary vortex and the shed vortices. Despite the appearance of rapid dissipation, however, the kinetic energy of these no-slip wall simulations does not decrease signif- icantly more rapidly than our simulations with a free-slip wall. 1.2 Stability and Reconnection of Vor- tex Tubes Spreiter & Sacks (1951) show that the trailing wake shed from an elliptically loaded wing will roll up into a pair of trailing vortices. The span of the trailing vortices in terms of the original span of the foil (L,) is 8 = l/4, and the core ra- dius, assuming constant strength, is rc = 0.0985s. This trailing vortex pair is unstable to sinusoidal disturbances along the length of the tube (Crow, 1970~. The wavelength of the instabilities (Acrou') in terms of the span are Acrou' ~ 8.6s, and the e-folding time of the instability is Crow =l.2lTi, where Ti=2~2/I, is the time it takes the vortices to move one span length due to their own induc- tion. Leonard (1985) indicates that large changes in the curvature of a vortex tube may lead to the rapid formation of helical vortex lines and even vortex breakdown. The exponential growth of helical vor- ticity is related to unstable perturbations of the vortex core in the presence of straining fields. The two-dimensional numerical simulations of Dritschel (1989) suggest that in addition to causing the ini- tial instability, the straining flows may also strip the helical vortices away from the primary vor- tex. According to Leonard's (1985) review, once the helical vortices are formed, the axial flows that are induced by the helical vortices may tend to straighten the primary vortex tube. As Crow (1970) notes, trailing vortex tubes that have gone unstable may eventually join to form crude vortex rings in the late stages of the in- stability. This process of reconnection is also ob- served in the vortex ring experiments of Kambe 729

OCR for page 727
l\kao (1971)1 Pall ~ leaner (19?5), Ostiroa ~ As^ (1977), Schlitz (1987)) and Oslo Izutsu (1988). The bbserstlons in these expert meets seems to locate that vortex reconnection takes place as a rest It of canceH,tlon of vortlc1ty as two eposltely signed vortex tubes lutersect. Ale- l~nder ~ Husband (1988) cut-and-connect rmodel of vortex reconnection (lg88), Stacy is based on spectra shmmdatlons of ant~parsUel vortex tubes, describes vortex reconnection in three stages: 1. I~iscid fnducifon. SelElnduced velocities cause the vortex cores to approach each other and fern a contact zone. Ibe vortex cores flatten and stretch to ferry a Pole cross- sectlon. 2. I. \brt1chy is amUldlated in the con- tact zone and bridges of cross-Dnked vortex lines fbrr1 at both ends of the contact zone. Ibe bridges are orthogonal to the posHlons the vortex tubes bad pclor to contact. 3. ~. Refits of the orlyln~ vortex ties Form between the brldees as the vortex ties pad fit. NO salutes of vortex rlug reco~ectlons Facade W~cke~ans ~ Leonardls (1989) vor- ton ~thod1 Ch~erlaln ~ Lluls fl985' ~te- ~~ence Todd Id Klda1 ~aoka1 ~ Hus- s~nls (1989) t~ee-~menslon~ spectra methods. A notable Bathe of Ch~erla~ ~ gluts ~te- ~rence echelon is tbelr use of m~tEpole expand . - ,. .. . Sons to specify a ramatlon conmbon. According to He~oltz theorems vortex Laments ~ bag ~ ~ -- a Id ~d~ ~ cous) because the vortlcity Held is solenoid. In ad~tlon Kelv1n~s clrc~atlon theorem states that vortex Gents move as maters Dues Id coot be brown ~ ~ lovlscld Had. As a reset vor- tex reco~ctlon is ~ l~erently viscous phenol ena. Slggla ~ Parks (1987) an~ysls of vortex Wants suggests that reco~ectlon shad the place on viscous tl~ scales' as one fight elect. Howlers Mehon, e# sZ1 (1989j provide _Hc~ Id theoretlc~ evidence that the reco~ectlon pro- cess occurs on nearly couvectlve time scales. Thea Wants we s~ported by Kwouls (1989) emery Into etudes of vortex rings co~ectlng fib tree s=~cee ~ch show no time dependence as a ~e don of Reynolds namer. The e~er~nts of Saga ~ Henderson (1984) Id S=p~a (1986) locate that Crow 1nst~1D- ties can develop ~ the trying awes of deeply s-- ~rged by~o~lb. 0= own Navler-Stokes d~a- dons of a sl~so1d~y perturbed vortex the ~v- lng pang to a gee sauce show that reco~ecdon of the parlay vortex the with the bee sauce ~ possible even if the co~utatl~ domain is too short to Perot Crow lastly to occur. Based the experl~nt~ evidence we gay speculate that the preferred ws=kngth of parlay vortex re- co~ecHon corresponds to a Crow ~st~lUt~ but shown by on own ~Hc~ signs Crow st~lUtles me not a necessary conatlon far the co~ectlon of pray vortlcHy with the tree s=- ~ce. 0~ ~erlc~ d~atlons also 1nacate that the generation of bedded vort1chy plays ~ essentl~ role in the reco~ect10n of secondly vortlchy filth the hee-s=~ce. Inst~1Utles in the parlay vor- tex the Open to be the source of Bach vorucH~ Once the beach voracity ~ grad ~ merges with the secondly vort1cHy that is separated Mom the hee~ce bo~d~y layer by the parlay vortex the. Was Dogleg process eventually leads to the rmatlon of U-sbeed vortices that me wrapped goad the parlay vortex the. Ibe D-vortlces gay connect with the gee sauce at their bases or tips. 1.3 Iuteractlou of ~rtlces alto ~ee- Sur~ces Mama ~ Barnes 'p89) e~ed~nt~ gauges of jet Rows interacting filth Gee sauces show that the lnteractlons of vortlchy filth the gee sauce gay generate short instead awes. They Id that sauce waves me generated as vortlc~ struc- t~es in the jet approach the bee sauce. downstream the sauce Titles me caused by l~ge-sc~e vortlc~ structures interacting erectly with gee sauce. They also observe that vortex lines co~ectlng filth the tree sauce me cow On ~ downstream. The early development of the ne=-s=~ce jet UseK is studied experlmen- t~y by Llep_ (1990~. D=lug the tr~dtlon phase Mom Lang to tangent A~1 scourged Id ne=-s=~ce lets convert az~th~ vortlclty (vortlc~ rings) luto stre~lse fortlolt~ His e~- perl~nts show that the presence of a tree sauce accelerates the Talon of stre~lse vortlolt~ 730

OCR for page 727
Sarpkaya, et al, (1988) use experiments and nu- merical simulations to study the normal incidence of a two-dimensional vortex pair with the free sur- face. A novel counter-rotating plate arrangement is used during the experiments to generate vortex pairs. For their numerical simulations the vortex pair is modeled using two point vortices and line vortices are used to model the free surface. Numer- ical results and experiments are presented for two Froude numbers, Fr = 6.4 & 2.8. The Reynolds numbers are not provided. For low Froude num- bers (Fr < 0.8), the experiments indicate that the free-surface acts like a free-slip wall, and the vor- tex centers spread apart. A scar (a dip in the free- surface elevation) appears in front of each vortex and moves with the vortex. At higher Froude num- bers, both experiments and simulations show that a hump of fluid rises above the vortex pair. The vortex pair subsequently spreads apart before dif- fusing due to turbulence. Telste's (1989) numerical studies of the interaction of two point vortices with a free surface suggest that vortices translating beneath a free surface will not generate significant wavetrains at low Froude numbers, Fr < 0.5. For higher Froude numbers Telste's results indicate that wave breaking may inhibit the outward expansion of the vortex pair. For Fr ~ 6.4 the point vortices entrain themselves in a hump of fluid that is located entirely above the mean position of free surface. Telste's numeri- cal simulations indicate that the hump would even- tually collapse on top of itself. For an intermedi- ate Froude number (Fr ~ 2.0), steep wave troughs form on each side of a small hump of fluid. The curvature of the wave troughs is so high that wave breaking appears to be imminent. Willmarth, et al, (1989) also study the normal in- cidence of a two-dimensional vortex pair with a free surface. They compare experimental results to the predictions of an inviscid theory which uses dipole sheets to represent the free surface and the vortex pair. They report their Froude number as Fr = 2.2, and we estimate their Reynolds and We- ber numbers as Re ~ 5, 600 and We ~ 60. As the Kelvin oval (the volume of fluid entrained by the vortex pair) rises upward, the free surface forms an upwelling that is centered above the vortices. Two sharp troughs form at the edges of the upwelling. The vortex pair rises slightly into the upwelling and then sinks as the upwelling springs back. Yu & Tryggvason (1990), call this phenomena, which de pends on the Froude number, an inviscid rebound- ing effect. Although the flow field induced by the vortex pair becomes turbulent during the experi- ment, Willmarth, et al, do not report observing the formation of secondary vortices. Ohring & Lugts' (1989) numerical studies of the viscous interactions of two-dimensional vortex pairs with a free surface indicate that secondary vortices form at low Froude numbers. They use the method of artificial compressibility with a general- ized coordinate system to solve the Navier-Stokes equations. They consider two cases with very low Reynolds numbers (Re < 20) and intermediate Froude numbers (Fr ~ 0.3 & 1.0~. Surface ten- sion effects are not modeled, but exact boundary conditions on the free surface are otherwise speci- fied. For the high Froude number result, the vortex pair entrains itself in the upwelling, and a patch of counter-sign vorticity develops at the base of the upwelling where there is high curvature. At the low Froude number, the height of the upwelling is sig- nificantly reduced and the vortex pair splits apart upon impact with the free surface. Secondary vor- tices form and partially wrap themselves around the primary vortex, but the Reynolds number is 50 low that the vortices lose too much energy for significant inviscid interactions to take place. The normal and oblique incidence of vortex rings with a free surface has been studied experimentally by Kwon (1989~. Kwon generates laminar vortex rings with Reynolds numbers, 600 ~ Re < 3200, and Froude numbers, 0.14 ~ Fr < 0.32. We esti- mate that his Web er numbers are 0.5 < We < 5. Kwon documents single and double reconnections at oblique angles of incidence. According to Kwon the reconnection time depends on the angle of in- cidence but not the Reynolds number. The experiments of Sarpkaya & Henderson (1984) (see also Sarpkaya, 1986) show that trailing vor- tices shed by wings may interact with free surfaces to form striations and scars. The striations oc- cur during the initial phases of the interaction and are normal to the wing's track. The scar features form after the striations, and they are parallel to the wing's track. The scars appear to be localized free-surface disturbances which move with the vor- tices. The early development of the scars is similar to the free-surface upwellings and the steep wave troughs that have been already discussed. In the later stages of the flow, dimples form in the scar region due to normal vorticity connecting with the 73~

OCR for page 727
free surface. They observe that Crow instabilities can occur for deeply submerged wings, which of- ten leads to the connection of vortices with the free surface and a crescent shaped disturbance on the free surface. The nondimensional parameters (in our notation) which characterize their experi- ments are as follows: 7,000 < Re < 21,000 and 0.3 < Or ~ 1. Based on this data, we also estimate that 10 < We < 100. In this parameter regime, Sarpkaya (1986) states that the free-surface signa- tures are not sensitive to changes in the Reynolds and Froude numbers. Sarpkaya claims further that the influence of capillary waves on his free-surface signatures is very small, although no experimental evidence is offered. In contrast to the results of Sarpkaya's experi- ments, Bernal, et al, (1989) and Kwon (1989) re- port that the presence of surface active agents on free surfaces strongly influences the interaction of ring vortices with free surfaces. They find that the motions of a vortex ring normally incident to a contaminated free surface are similar to the mo- tions caused by a no-slip wall. Secondary and tertiary vortices are shed, and the primary vor- tex rebounds. For a cleaner free surface they find that the shed vortices are weaker and the radial expansion of the primary vortex is greater. They report Re ~ 1200 and Fr ~ 0.09, and we esti- mate We ~ 0.9. These physical parameters are less than Sarpkaya's (1986) and may explain why Sarpkaya's observations differ from Bernal, et al, in regard to the influence of surface contamination on surface signatures. Bernal, et al, (1989), also report preliminary experiments with delta wings towed beneath a free surface. They observe what appears to be a Reynolds ridge, which is the flow region beneath a stagnant film on the water sur- face (see Scott, 1982~. According to them, the scar features that are observed by Sarpkaya may in fact be Reynolds ridges. Interestingly, Grosen- baugh & Yueng (1988) report that surface contam- ination strongly influences the onset of turbulence in the bow waves of model ships. They describe a stagnant zone which causes the free surface to separate farther upstream than it would in the ab- sence of surface contamination. This description fits Scott's definition of a Reynolds ridge. As one example of the range of scales that are possible for full-scale ships, consider a lightly- loaded foil traveling slightly below and parallel to the free surface. For negative lift, the trail ing tip vortices will rise upward. We assume a long span, s=O(lOm) and a low circulation, I, = 0~1m2/~. In our notation we have then relatively high Reynolds numbers (Re ~ 105), low Froude numbers (fir ~ 10-2), and high Web er numbers (We ~ 103~. This parameter regime is far out- side the numerical simulations and laboratory ex- periments that we have reviewed in this section. The impact of this mismatch is not clear. We note that our own numerical simulations of vortex interactions with a free surface are near the ap- propriate Froude numbers, but our Reynolds num- bers (Re = 0~103~) are too low. Hopefully, the three-dimensional flow phenomena that we observe at these low Reynolds numbers capture some of the important features of flows at much higher Reynolds numbers. We emphasize that inviscid approximations are inadequate because the salient features of vortex interactions with a free surface such as boundary-layer formation and vortex re- connections are fundamentally viscous. Having noted these limitations, we now proceed to summarize the free-surface phenomena that we observe in our direct numerical simulations of the Navier-Stokes equations. Our numerical studies of vortex rings interacting with a free surface show that secondary vortex rings may be torn away from the free-surface boundary layer at intermedi- ate Froude numbers, Ad at lower Froude numbers the free surface acts like a free-slip wall. Our nu- merical simulations at low Froude numbers agree qualitatively with the vortex ring experiments of Bernal, et al, (1989) and Kwon (1989~. Our stud- ies of vortex tubes illustrate two possible mecha- nisms for the reconnection of vorticity with a free surface including primary and secondary vortic- ity reconnections. The free-surface signatures of these two mechanisms are distinctly different. One type of primary vorticity reconnection should re- sult in a cusp pattern on the free surface, and sec- ondary vorticity reconnection should manifest it- self as paired dimples on the free surface. These numerical experiments show that the generation of helical vortices is an initial stage in the reconnec- tion of the secondary vorticity with the free surface, and our simulations also suggest that the striations that may be observed in the waltes of trailing tip vortices are caused by the interaction of the helical vorticity with the free-surface boundary layer. 732

OCR for page 727
1.4 Numerical Methods A very general method for the direct simulation of the incompressible Navier-Stokes equations with free-surface effects is the MAC (Marker-And-Cell) method (see Hirt, et al, 1975~. MAC methods are distinguished from other numerical methods by their technique for solving for the pressure, treat- ment of the convective terms, formulation on a staggered grid, and capability for particle tracing. MAC methods enforce mass conservation by solv- ing a Poisson equation for the pressure. The MAC equation for the pressure retains temporal and spa- tial derivatives of the dilation term (V u) to inhibit the growth of numerical errors, and the convective terms in the momentum equations are calculated using a combination of upwinding and central dif- ferencing. The MAC time integrator without the upwinded terms is similar to the explicit FTCS (E`orward Time Central Space) method (Roache, 1976). As evidenced by the application of the MAC tech- nique to a wide range of problems including nonlin- ear ship wave problems (see Miyata & Nishimura, 1985), the MAC method is useful for determin- ing qualitative features of very complicated free- surface flows. For simple geometries the high accuracy of spectral methods relative to finite- difference methods has proven to be very de- sirable. A good example is Orszag & Patter- son's (1972) direct simulation of homogeneous and isotropic turbulence using Fourier series. In the- ory Patera's (1984) spectral element method al- lows complex geometric modeling with spectral ac- curacy, but in practice general three-dimensional problems have not been very accommodating. As a result, low-order firLite-difference methods with high resolution appear at this moment to be more robust than spectral methods. The MAC and FTCS methods are subject to stringent stability criteria because of their use of explicit time integrators, and as a result, sev- eral semi-implicit schemes have been proposed for overcoming this deficiency, including SIMPLE (Patankar & Spalding, 1972), SIMPLER (Pan- tankar, 1981), PISO (Issa, 1985), and the Finite- Analytic method (Chen & Chen, 1982~. These four methods solve the momentum equations for the ve- locities and the Poisson equation for the pressure separately, and iterations are used to strongly cou- ple the velocities with the pressures. An interest ing feature of Chen & Chen's method is their use of locally analytic solutions to solve the Navier- Stokes equations. They claim that this technique gives superior modeling of the convective terms without incurring a numerical diffusion penalty. However, the Finite-Analytic method requires sig- nificantly more memory without improving the second-order spatial accuracy of more conventional finite-difference methods. In contrast to these iterative solution procedures, Kim & Moin's (1985) fractional step method uses only one iteration to couple the velocity field with the pressure field to maintain a divergence-free flow. Note that Kim & Main, just like Ghia, et al, (1977), use Neumann boundary conditions in the Poisson equation for the pressure, and as a result, they must satisfy a numerical solvability condition. Kim & Moin are able to achieve second-order accu- racy in time by using an explicit Adams-Bashforth scheme for the convective terms and an implicit Crank-Nicholson scheme for the viscous terms. A1- though Kim & Moin's hybrid time-integrator elim- inates time-step restrictions due to viscosity, their numerical scheme must still meet a Courant con- dition because of their explicit treatment of the convective terms. A fully-implicit scheme has been devised by Hartwich & Hsu (1988) who uses Chorin's (1967) artificial compressibility method as a framework. Chorin's method inserts a time-dependent pres- sure term in the mass-conservation equations, and as the solution nears steady-state, the pressure term drops out. This artifice permits the use of advanced compressible flow techniques (e.g., Warming & Beam, 1978~. In addition, the tech- nique can be generalized to time-dependent prob- lems. Other researchers including Brandt & Di- nar (1977), Fuchs & Zhao (1984), Vanka (1986), and Thompson & Ferziger (1989) solve the incom- pressible equations implicitly. Typically, these im- plicit primit~ve-variable solution schemes use block Gauss-Seidel iterations to drive the divergence to zero, and multigrid methods are used to acceler- ate the solution procedure. The block in lieu of point Gauss-Seidel solvers are necessary because the finite-diRerence approximations to the mass- conservation equation are not diagonally domi- nant. So far we have limited our discussions to the primitive-variable formulation of the Navier-Stokes equations. The widespread use of vorticity formu 733

OCR for page 727
rations also warrant a discussion. More compre- hensive reviews of vortex methods can be found in Saffman & Baker (1979) and Leonard (1980 & 1985~. There are two main vorticity formulations: vortex element methods which use line and sheet singularities and finite-difference approximations of the vorticity-stream function equations. Vor- tex element methods are ideally suited for sim- ulating small regions of rotational flow such as isolated vortex tubes (Leonard, 1985) and sharp density interfaces (Dahm, Scheil, & Tryggvason, 1989~. Generally, large core deformations and vis- cous phenomena are not modeled well by three- dimensional vortex element methods, but Winck- elrnans & Leonard (1989) have developed a vorton (vortex sticks) method to study the initial stages of vortex reconnection. Finite-difference approxi- mations to the vorticity-stream function equations do allow large core deformations with strong vis- cous interactions, but boundary conditions other than free-slip boundary conditions are difficult to implement in three dimensions. We note that vor- tex methods are continually being improved and their range of application will certainly increase with time. Regardless of what formulation is used to simulate the incompressible Navier-Stokes equations, im- plicit schemes are traditionally used for rapidly at- taining steady-state solutions and explicit schemes for simulating time-dependent flows. Explicit schemes are not used to reach steady state because they require too many time steps, and implicit schemes are not used for time-dependent problems because they are computationally more intensive and the time steps are so small based on accuracy considerations that stability is often not a prob- lem for explicit schemes. In some special cases of time-dependent flows, such as creeping flows and capillary waves, implicit schemes are used because the stability limits of explicit schemes are too re- strictive. In this paper, we present a second-order explicit method and a first-order semi-implicit method for solving the three-dimensional Navier-Stokes equa- tions in primitive-variable form with and with- out free surfaces. The unconditional stability of the semi-implicit scheme is useful for obtaining steady-state solutions, and the accuracy of the ex- plicit scheme is desirable for time-domain prob- lems. In both numerical schemes the effects of numerical viscosity are reduced to a minimum by using central differencing. The nonlinear momen- tum equations in the implicit scheme are solved using Newton-Raphson linearization and multigrid iteration. We also use multigrid iteration to solve the Poisson equations for the pressures in both time-stepping procedures. The convergence prop- erties of the numerical schemes are demonstrated for three test problems: (1) the a~cisymmetric stag- nation flow against a wall, (2) the viscous attenua- tion of an axisymmetric standing wave, and (3) the translation and diffusion of a two-dimensional vor- tex. The numerical schemes are used to study sev- eral physical problems including vortex rings and tubes impacting walls and free surfaces. 2 Mathematical Formulation 2.1 Field Equations We consider the unsteady incompressible viscous flow of a Newtonian fluid under a free surface. For an isotropic and homogeneous fluid, the Navier- Stokes equations for conservation of momentum have the form: - + (u V)U = -Up + R-V2U, e (1) where u = u(:e, y, z, I) = (a, v, w) is a three- dimensional vector field. Note that the momen- tum equations have been normalized based on a length scale L and a velocity scale U. Re = UL/L, is the Reynolds number where ~ is the kinematic viscosity. In addition, we have defined p as the hydrodynamic pressures P=P-F2Z, r (2) where P is the total pressure which is equal to the sum of the hydrodynamic plus hydrostatic pres- sures. The pressure terms are normalized by pU2 where p is the density. F2 = U2/gL is the Froude number. The vertical coordinate z is positive up- ward, and the origin is located at the mean free surface. The dyn~nic pressure may be interpreted as a Lagrange multiplier which projects the veloc- ity field onto a divergence-free field as expressed by the mass-conservation equation: V u = 0 . 734 (3)

OCR for page 727
The pressure field itself does not require boundary conditions or initial conditions, but a well-posed initial-boundary-value problem does require these conditions for the velocities. Note, however, that the field equations and boundary conditions for the velocity field can be used to deduce the ini- tial pressure and the behavior of the pressure near the boundaries. For example, the divergence of the momentum equations (eqts. 1) used in combi- nation with the mass-conservation equation (eqt. 3) may be used to derive a useful Poisson equa- tion for the pressure. This equation expressed in indicial notation is as follows: V2p = _ ~Uj UP ant B~j Sirrularly, the momentum equations may also be used to prescribe the normal derivative of the pres- sure on the boundaries of the fluid. Thus, accord- ing to the divergence theorem the pressure is sub- ject to the following solvability condition: /s an /v {3zi ire ' ~ ~ where V is a volume of fluid, S is the surface bounding the volume, and n is the unit outward- pointing normal on that surface. For the velocities, we are required to specify the initial velocity field: u = uO(:z, y, z) at t = 0 . (6) According to Helmholtz's theorem, we can express this initial vector velocity field uO in terms of scalar 0 and vector JO = A?, ye, ?~Az~o velocity poten- tials: Uo=Vo+VX?po. (7) The initial velocity field, like the time-dependent velocity field, is required to be solenoidal. There- fore, we can show that by taking the divergence of equation (7) that TO satisfies Laplace's equation: V2+o = 0 . (8) Here, TO may represent the effects of currents, the disturbances of surface-piercing and submerged bodies, etc. A similar Poisson equation is satisfied by the vector potential To in terms of the initial vorticity field wO V2o =-To ~ (9) where ~ = (wc,~',,wz) may be expressed in terms of the curl of the velocity field at any instant of time: ~ = V x u. (10) Based on this definition, the vorticity field is itself a solenoidal quantity. In addition to posing ini- tial conditions for the velocity field, we must also specify boundary conditions as discussed in the fol- lowing sections. 2.2 Linearized E ree-Surface Boundary Conditions Let the free surface elevation be given by ,1~:r, y, t). We posit the presence of a viscous free-surface boundary layer scaled by fS=O(Re i/2~. We ob- tain a linearization of the free-surface boundary conditions by assuming small surface slopes, i.e., '7~ & t7',=0~) << 1, and consequently ,7 is much smaller than the other length scales of the flow including [. Under these assumptions, the exact (nonlinear) field equations are still valid out- side the boundary layer, but the nonlinear kine- matic and dynamic free-surface boundary condi- tions (Wehausen & Laitone, 1960) reduce to: uz . O vz = 0 (11) 1 -P + F2 ~ ql + use + van = = (12) where We = pU2L/T is the Web er number, p is the water density, and T is the surface tension, and the conditions are applied on z=O. The first set of equations (11) state that the shear stresses and normal stresses are zero on the free surface. Note that the kinematic free-surface condition (12) is still nonlinear in the flow variables and is obtained without otherwise requiring that the velocities be sit. 73s

OCR for page 727
2.3 E`ree-slip Wall Boundary Conditions We next consider the case of a free-slip wall. This type of boundary condition is useful for imposing symmetry boundary conditions and modeling free- surfaces at low Froude numbers. By definition the tangential stresses and the normal flux are zero on a free-slip boundary. Therefore, the boundary conditions for the velocity on a free-slip plane at z = 0 are: uz=0, vz=0, and w=0, (13) where the first two equations give zero tangential stresses and the third equation states that there is no flux across the plane z = 0. These boundary conditions for the velocities en- able us to deduce the behavior on the free-slip wall of the pressure and vorticity, and also the initial conditions for the scalar and vector velocity po- tentials. For example, upon substitution of the boundary conditions (13) into the z-component of (1), we can show that pz=0, where the z-derivative of (3) has been used to eliminate the u~zz=0 term. Similar arguments may be used to show that the vorticity vector is always normal to a free-slip wall, ~ = (O,O,wz) on z=0. Finally, in regard to the initial conditions, we can derive the following rela- tions for the vector velocity potential on the free- slip plane z = 0: (~10 = 0, (~10= 0, B(Z) = 0, and V21,bo = -~0, 0, (~z)O), (14) where we have assumed that the scalar velocity potential is initially zero, .50 = 0. 2.4 No-slip Wall Boundary Conditions On a no-slip wall the velocity is zero: u = 0 . (15) In general, the shear stresses are nonzero and as a result, boundary-layer formation and flow separa- tion are possible. Certainly, the flow field near a no-slip wall is more complex than that for a free- slip wall, but it is still possible to deduce the behav ior of certain flow quantities near the wall. For ex- ample, suppose that a no-slip wall is located on the boundary x = 0, then the no-slip boundary condi- tions (15) and the mass-conservation equation (3) may be used to show u,* = 0. Also, substitution of the no-slip conditions into the z-component of momentum (1), gives the normal derivative of the dynamic pressure: R ZZ (16) e A similar procedure may be used to show that the source term is zero in the Poisson equation for the pressure (4~. Unlike a free-slip wall, the vorticity vector must be parallel to a no-slip wall, w = (~v, O) on z = 0, and as a result a vortex tube cannot terminate normal to a no-slip wall. For a no-slip wall it is difficult to pose unique ini- tial conditions and boundary conditions that are physically meaningful for the vector velocity poten- tial. One method for eliminating this problem is to let the velocity be prescribed by free-slip boundary conditions at time t = 0, and then at time t = 0+ set all velocity components equal to zero. As a re- sult, a sheet of vorticity is instantaneously formed on the no-slip wall. 2.5 Conservation of Energy The vector product of the velocity with the mo- mentum equations ( 1) integrated over the fluid vol- ume gives a formula for the conservation of energy, and the transport theorem in conjunction with di- vergence theorem may be used to simplify the re- sulting equations. In indicial notation the equation is as follows: d ~ Mini . = at JV 2 /s njUit-P6ij + R ~ a0Ui + ~Uj )] Re /v(8zj + ~Zj)0Zj ' (17) where the first term represents the change in ki- netic energy integrated over the material volume of the fluid (V), the second term represents the work done by the stresses on the fluid boundaries (S), and the last term equals the energy that is dissi- pated. (tij denotes the Kroenecker delta function.) 736

OCR for page 727
Note that the stresses do not do any work on slip and no-slip walls. Upon substitution of the exact free-surface boundary conditions into this energy equation, we may derive the following formula: d ~ nisi 1 d ~ 2 ~ - ~ ~ ~ = d! Jv 2 ' 2F2 dt Js' We |S~ ~t( Ri + R2 ) Re JV( Ski + 02'i ) I9~ , (18) where the first term represents the change in ki- netic energy, the second term the change in poten- tial energy, the third term the work done by surface tension, and the last term the dissipated energy. Sf is the projection of the free-surface onto the xy- plane and Rat and R2 are the principal radii of cur- vature. Note that we have assumed that the work done by the stresses on all other boundaries besides the free surface is zero. If the free-surface elevation and slope are small, then the volume integrals are evaluated below the mean waterline (z = 0) and the principal radii of curvature terms simplify to -(~= + ~vu) 3 Numerical Formulation In this section we derive semi-implicit and explicit finite-difference schemes for direct simulation of the three-dimensional Navier-Stokes equations at low Reynolds numbers. (See Appendix A for an outline of the axisymmetric formulation:) We show how the initial velocity field is set up, and demon- strate the implementation of periodic, slip, no-slip, and free-surface boundary conditions. For periodic and/or wall boundary conditions a solvability con- dition for the pressure is also implemented. We then discuss the vectorization of a multigrid solu- tion technique for solving the field equations. We conclude the numerical formulation section with stability analyses and validation studies. The val- idation studies include numerical simulations of axisymmetric stagnation flows, attenuation of ax- isymmetric standing waves, and translation and diffusion of two-dimensional vortices. 3.1 A Semi-Implicit Time-stepping Scheme We derive here a first-order semi-implicit scheme for solving the Navier-Stokes equations. The basic concept of our scheme is based on the PISO scheme (Issa, 1985) which updates the velocities indepen- dently of the pressure. We use Newton-Raphson iteration to solve the nonlinear Navier-Stokes equa- tions, and we assume that the velocities for each implicit time step are equal to some known esti- mates plus some small corrections. Let superscript If" denote the n-th time step, and consider the per- turbations of the velocities at the In + lath time step as follows: Un+l = U + U vn+! = V + V (19) wn+! = W + Hi ~ where (u,v,w~n+~ are the unknown velocities at the next time step, (U,V,W) are the known es- timates, and (u, v, fur are the unknown errors. Note that for a fully-implicit scheme the hydro- dynamic pressure would be expanded in a sim- ilar way. For our serru-implicit scheme the hy- drodynamic pressure is always known in terms of the most recent estimate of (u,v,w~n+~. In (19) we continuously update (u, v, w~n+i until the de- sired accuracy has been achieved. First we es- timate the fluid velocity at the In + Both time step by letting (U. V, W) = (u, v, win, and then we find the error terms, (u, v, w). Our improved es- timate of the velocities is given as (u, v, we+ = (U. V, W) + (u, v, w) and the new pressure can be calculated based on these updated velocities. The Newton-Raphson iterations can stop at this point if the errors are small enough, or they may con- tinue by letting (U. V, W) equal our new estimate of (u, v, we+, etc. Generally, only one to two iter- ations are necessary because the Newton-Raphson procedure has quadratic convergence. Upon substitution of the perturbation expansions (19) into the momentum equations (1) and elimi- nation of quadratic error terms, we derive the fol- lowing set of linearized equations: us + 2(uU)~ + (uV)v + (Uv),, + (uW)z + ( Uw In- -V2u = Ures Rc 737

OCR for page 727
Figure 23e: 1~1 = 1.25 isosurfaces at t=15. 778

OCR for page 727
Figure 24: The interactions of a vortex tube with a no-slip wall. The x~ = 0.5 (helical vor- ticity), ~`,~=0.5 (axial vorticity), and ~ - 0.5 (total vorticity) isosurfaces are plotted at different instants of time: (a) t=5, (b) t=7.5, (c) t=10, (d) t=12.5, (e) t-15, (f) t=17.5, & (g) t=20. The view is from below a no-slip wall, and the primary vor- tex is moving into and to the left of the page due to its images across the centerplane and above the wall. The geometry triad is in a submerged corner on the centerplane of the computational domain. Note that the perspective in Part (g) is slightly different from Parts (a)-(f). The numerical param- eters for this run are provided in Run 3 of Table (43. Figure 24a: x/~ = 0.5 isosurfaces at t=5. Figure 24a: two = 0.5 isosurfaces at t=5. \1 779 Figure 24a: 1~',1 - 0.5 isosurfaces at t=5. - ~1

OCR for page 727
Figure 24b: `/~ = 0.5 isosurfaces at t= 7.5. _~ Figure 24b: ~-0.5 isosurfaces at t=7.5. Figure 24b: 1~9,1 = 0.5 isosurfaces at t=7.5. 780

OCR for page 727
Figure 24c: `~ = 0.5 isosurfaces at t=10. Figure 24c: ~ = 0.5 isosurfaces at t=10. Figure 24c: levy = 0~5 isosurfaces at t-10. \ \1 781

OCR for page 727
- Figure 24d: )~2 + W2 = 0.5 isos~faces at t=12.5. Figure 24d: ~ ~ _ 0.5 isosurfaces at t=12.5. Figure 24d: 1~Vl-0.5 isosurfaces at t=12.5. 782 -

OCR for page 727
1 Figure 24e: I+ = 0.5 isosurfaces at t-15 Figure 24e: iw~ = 0.5 isosurfaces at t=15. Figure 24e: 1~',1 = 0.5 isosurfaces at t=15. 783

OCR for page 727
Figure 24f: ,~ = 0.5 isosurfaces at t_17.5. Figure 24f: 1wl = 0.5 isosurfaces at t=17.5. 784 Figure 24f: levy = 0~5 isosurfaces at t=17.5. y

OCR for page 727
Figure 24g: x/~ = 0.5 isosurfaces at t=20. Figure 24g: ~ = 0.5 isosurfaces at t=20. / 785 Figure 24g: 1~Vl = 0~5 isosurfaces at t=20.

OCR for page 727
Figure 25: The unwrapping of U-shaped vortex tubes around the primary Forte:: tube. The two = 0.5 iso- surface of vorticity is plotted at time t=20. To em- phasize the U-shaped feature this image has been reflected about its midspan relative to the images appearing in Figures (24~. The numerical param- eters for this run are provided in Run 3 of Table (4~. Figure 26: The conservation of energy for a three- dimensional vortex: tube interacting with a no-slip wall as a function of time. The kinetic energy E(t) is compared to the energy that is dissipated D(t) (the Reynolds number term in eqt. (18) integrated over time). The results are plotted as follows: ~ 1 ~ denotes E(t)/E(O) and ~-2-~ de notes (E(O) + D(`t)~/E(O). The numerical param- eters for this run are provided in Run 3 of Table (4~. .00 .95 .90 .85 .80 .75 .69 .65 .60 .55 .50 .45 1 1 1 1 1, 1 ' 1 ' 1 ' 1 ' 1 1 1 ~1 my. ,~ - 1 1 ~1 1 0 2 4 6 8 10 t Am' I , 1 , 1 1 1 1 12 14 16 1 786 8 20 J

OCR for page 727
Figure 27: The helical vorticity evaluated on the wall. The contours of w2 on the wall are plot- ted for different instants of time: (a) t=10, (b) t=15, and (c) t=20. Observe that the features are oriented normal to the axis of the primary vortex tube which is parallel to the Taxis. The numerical parameters for this run are provided in Run 3 of Table (4~. ~9 ~ , 0 .5 1 0 1.5 2.0 2.5 3.0 3 5 4.0 X 787 1 . i ~, . t . ~ |CONTOUR FROM - 7 TO .2 BY ;OC I J 0 .s 1 0 1 5 2.0 2 5 3.0 3.5 4.0 X Figure 27b: t=15.

OCR for page 727
DISCUSSION P. Ananthakrishnan University of California at Berkeley, USA (India) My question/remark is regarding the results corresponding to the interactions due to a vortex ring that is normally incident on the free surface. According to Song, et al (in M. Song, N. Kachman, J.T. Kwon, L.P. Bernal, and G. Tryggvason, Vortex Ring Interaction with a Free Surface,. preprint, 18th Symp. on Naval Hydro.), a dip on the free surface right above the vortex core, which moves radially outward as the vortex ring stretches at the free surface, is observed. They also observe three-dimensional small scale motions and propagating waves at higher Froude numbers. But your results (Fig. 17) show a rise in wave elevation at the axis of symmetry that is comparable in magnitude to the surface depression. Could you comment on the possible reasons for the difference between the observed and computed results? AUTHORS' REPLY Relative to Song, et al's, (1990) experiments, our numerical axisymmetric results are valid for much higher Froude numbers. Even so, we also predict the free-surface depression that occurs above the vortex core as is evident in our Fig. (16). In regard to the small- scale three-dimensional features that are observed in experiments, we believe that some of these features may be explained in terms of U- vortices. For example, Fig, (2b) of Song, et al's (1990) paper shows a U-shaped feature wrapped around the primary vortex tube. The base of the U is on the inside of the ring vortex. Both these features are consistent with the U-vortex phenomena as explained in our paper DISCUSSION Ronald W. Yeung University of California at Berkeley, USA In your results corresponding to Fig. (16) and Fig. (17), you have mentioned the existence of a non-zero mean wave elevation, which you artificially subtract off to obtain presented results. The non-zero mean implies that mass is not conserved; I was wondering if the authors can shed some light on the source of this trouble? AUTHORS' REPLY For our linearized free-surface boundary conditions, the free-surface elevation is much smaller than the boundary-layer thickness (~1 ~ b). The Froude and Reynolds numbers that are illustrated in Figs. (16 & 17) push the limits of this theory so that an interesting physical regime could be investigated. Fully-nonlinear free-surface boundary conditions would eliminate the problem with mass conservation. DISCUSSION Fred Stern The University of Iowa, USA The authors should not use the terminology Direct simulation. to refer to their solutions of the unsteady Navier-Stokes equations for low Reynolds numbers (i.e., laminar flow). Generally, this terminology has become synonymous with the direct simulation of turbulence through extremely high-resolution solutions of the unsteady Navier-Stokes equations for relative high Reynolds numbers as opposed to the use of the Revnolds-avera~ed Navier-Stokes equations. Would the authors please specifically point out the novel aspects of their method since most aspects appear to be familiar. In fact, this should be the focus of the discussion on the numerical formulation which is too long and often confusing. The real value of this work is in the use of such numerical methods for free-surface flows and, in particular, the applications chosen for study. The authors appear to have captured many of the observed phenomena although I have not been able to decipher all that they have from the figures. The free-surface boundary conditions used appear to be identical to those used in my paper. Essentially, these are inviscid approximations in which the viscous-stress conditions are neglected, or if you like, only satisfied to a very low order. It would appear then that many of the observed phenomena are pressure- driven. Also, some of the vorticity in the calculations near the free surface may be erroneous due to these approximations. Would the authors please comment on these points? AUTHORS' REPLY We agree that ~laminar flow simulation" is more appropriate terminology than Direct simulations. Both our semi-implicit and explicit schemes use unique fully-vectorized multi-grid methods to solve the unsteady three-dimensional Navier-Stokes equations at low Reynolds numbers with and without free surfaces. Our linearized free-surface boundary conditions are derived from the exact normal and tangential stress conditions and the exact kinematic free surface condition subject to our assumption of small free-surface slopes. The range of validity of these free-surface boundary conditions are discussed in our paper and in our answer to Professor Yeung's question. A curved free surface, just like the afterbody of a bluff object, may have unfavorable pressure gradients that lead to flow separation. 788