National Academies Press: OpenBook

Eighteenth Symposium on Naval Hydrodynamics (1991)

Chapter: A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface

« Previous: Viscous Flow Past a Ship in a Cross Current
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 727
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 728
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 729
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 730
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 731
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 732
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 733
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 734
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 735
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 736
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 737
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 738
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 739
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 740
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 741
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 742
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 743
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 744
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 745
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 746
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 747
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 748
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 749
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 750
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 751
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 752
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 753
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 754
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 755
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 756
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 757
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 758
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 759
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 760
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 761
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 762
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 763
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 764
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 765
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 766
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 767
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 768
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 769
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 770
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 771
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 772
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 773
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 774
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 775
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 776
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 777
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 778
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 779
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 780
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 781
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 782
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 783
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 784
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 785
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 786
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 787
Suggested Citation:"A Numerical study of Three-Dimensional Viscous Interactions of Vortices with a Free Surface." National Research Council. 1991. Eighteenth Symposium on Naval Hydrodynamics. Washington, DC: The National Academies Press. doi: 10.17226/1841.
×
Page 788

Below is the uncorrected machine-read text of this chapter, intended to provide our own search engines and external engines with highly rich, chapter-representative searchable text of each book. Because it is UNCORRECTED material, please consider the following text as a useful but insufficient proxy for the authoritative book pages.

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,

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

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

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

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~

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

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

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)

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=V¢o+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 V2¢o =-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

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

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

At + (TV) + (Uv),` + 2(vV)v + (vW)z +(Vw)z-R V2v = Ores (20) e oft + (uW)I + (UW)~ + (VW)v + (VU))v +2(u)W)z- -V A' = WreJ, Re where the unknown error terms are grouped on the left-hand side and (Ure,,, fares, WreJ) is the known residual error vector given below: Ure`' = -U.-( UU) ~ -(UV)~,-(UW)z-Pc + R V2U Ares = -Vt -( UV )~ -(Viny - (VW)z - Pv + R V2V (21) Wre,' = -Wt-(Uw)s where the first two difference formulas are first- order one-sided formulas for 9/9x, the third for- mula is a second-order centered formula for B/8x, and the last formula is a second-order centered for- mula for the Laplacian operator. Similar formulas are written for the partial derivatives with respect to V and z. We use upwind differencing in (20) to ensure di- agonal dominance, and we use central differencing in (21) to reduce the effects of numerical damping. Since the final converged solution for (u, v, VJ)n+l is central differenced, we are able to realize the efficiency of first-order upwind differencing with- out the numerical damping and the accuracy of second-order central differencing without the slow convergence. We illustrate how to upwind a single convective term: (uV)~ = -(u(V + ~V~g + 2(u(V - ~V~g -(VW)V-(WW)z-Ps + R V2w . (UV)v ~2V (u(V + EVE) + ~2 (u(V-EVE) In these equations for the residual errors, the hy- drodynamic pressure is expressed in terms of the most recent estimate of (u, v, If+. Note that the equations for the error terms (20) and the residual errors (21) are expressed in conservation forms. We solve the set of equations (20) and (21) for (u, v, u') by using finite differencing. Let (Ax, Ay, Liz) be the unstaggered grid spacing along the x,y, and x-axes respectively, and let At be the time step. Also, let subscript ()i,j,k denote the indices of the grid points along the x,y, arid x-axes respectively. We define the following finite- difFerence formulas: 6+F = F`,j,`-Fi-l,j, ~ fix 6- F Fi+l,j,`-Fi~j~k Ax tCF _ Fi+l~j~k-Fi-l,j,k ~ 2Ax V2F _ Fi+l,j,.-2Fi,j,k + Fi-l,j,k (~x)2 + Fi,j+l,6-2Fi~j~k + Fi,j-l,k (~y)2 + Fi,j,`+1-2Fi~j~k + Fi,j,k-1 (~z)2 UV)V ~ (U(V |V|))i,j+l,k + (UlVI)j job (U(V + IVl))i~j-leek _ ~ 2Ay This upwinding scheme always contributes a posi- tive element to the diagonal of the equation system. We use this upwinding procedure to discretize (20) as follows: (22) - ,\t + 6+[u(U + IUI)] + iz [u(u - IUI)] +.7~ [U(V + IVI)] + ~2 [u(V - IVI)] + 2 [v(U + IUI)] + ~2 [v(U - IUI)] + 62Z [u(W + IWI)] + ~2 [u(W - IWI)] + 623 [u'(U + IUI)] + ~2 [(Mu - IUI)] - 1 V2u = Uresij`, e /`t + 2 [u(V + IVI)] + ~2 [u(V-|V| )] +~2 [v(U + IUI)] + ~2 [v(U - IUI)] 738

+~+[v(v + lYl)] + 6V [v(v - lYl)] +~2 [V(W + IWI)] + ~2 [V(W-|W|)] (23) Uijk + 2 [W(V + IVI)] + 2 [W(V-|Vl)] Vijk = -R V V = vresiik At + 2 [u(W + IWI)] + ~2 [u(W - IWI)] +~2 [(Mu + IUI)] + ~2 [(Mu - IUI)] +~2 [V(W + IWI)] + ~2 [V(W - IWI)] s+ `_ + 2 [W(V + IVI)] + 2 [t,,(V - IVI)] +~+[~(W + IWI)] + [x [w(W - IWI)] -R V U' = Wre~is. , where a first-order Euler scheme is used to ap- proximate the time derivative terms. Note that by definition there are no errors in the velocities at the n-th time step, and so only one term is required to approximate each time derivative of (Ui,j,k, Vi,j,i, U)i,j,k). As pointed out earlier, the con- vective terms are upwinded in (23) to ensure the diagonal dominance of the equation system. The discretized forms for (21) are given by: aresiji = - 'i' ·,j~k) _ ~Ctou' -6~tUV)-~(UW)-b~P + R V2U (Vi j k-vnj k) ij. '\ t -~(VV) - [xtVW) - {',P + R V2V (24) wresij. = - 'i' i~j~k) _ tCtUW' 6~(VW) _ [C(WW) _ {cp+ ~ V2W where once again, a first-order Euler scheme is used for the time-derivative terms. The system of equations (23) and (24) are solved using Gauss-Seidel iteration. A good initial guess of the error terms may be obtained by ignoring all of the off-diagonal terms as illustrated below: AtUres;j,,' - di,.k 2\ tVrestj`' diivk 2\ tWres, i di~k (25) di~` = 1 + 2 t~Uu' ~ + At~V`j. ~ + At~W`*~ +2/~( 1 + 1 + 1 RC (~X)2 (^y)2 (^Z)2 d v ~t|Uijk| 2 ~t|Vijil ~t|Wiji ijk = 1 + ^~ + ~V ~z +2~( 1 +_ 1 + Re (^Z)2 (^y)2 (^Z)2 d w ~t|Uijk| ~t|Vijil 2 ~t|WiJk ijk = 1 + + + Re ( (Az)2 + (^y)2 + (Az)2 ) We observe that the FTCS scheme is recovered if all the At terms are eliminated in the denominators of (25) and if only one Newton-Raphson iteration is used. Also, if the time step is small enough such that the diffusion terms are much less than one, 2/~` ~ + ~ + ~ ~ ,2) ~ 1, the system of equations (23) and (24) are parabolic and Gauss-Seidel iteration may be used to effi- ciently solve the equations. For larger time steps, howe~rer, the diffusion terms dominate, the system of equations are elliptic, and multigrid iteration is much more effective. Finally, if the diffusion terms dominate and if the grid spacing along one or two axes is much smaller than the spacings along the other axes, the system of equations are anisotropic and a line Gauss-Seidel sol~rer should be used in conjunction with a multigrid solver. After solving (23) and (24), we first update the ve- locities (u, v, UJ)n+! = (U, V, W) + (U, V, W), then reassign the approximate solutions (U,V,W) = (u, v, 111)n+~, and finally update the pressure. The Poisson equation for the pressure (4) must be dis- cretized carefully because numerical errors may otherwise accumulate and cause the solution to di- verge. In particular, in the analytic derivation of 739

the Poisson equation, we neglected functions of the velocity dilation term (V · u) which may differ from zero in our numerical calculations. Thus, for the numerical implementation of our implicit scheme we retain the time derivative of the dilation term which is a function of the n-th time step. Func- tions of the dilation term at the (n + 1~-th time step are set equal to zero so that mass conserva- tion will be satisfied. Our discrete equation for the pressure IS: v2pn+1 = +_~`cun + scan + scant (26) +2(§Cun+1 [Vvn+l + [Cun+1 licWn+ +6Vvn+l [Zwn+l _ [CUn+1 6cvn+1 -lizun+1 [c wn+1 _ [Cvn+1 ~icwn+ We solve this elliptic equation by using a com- bination of Gauss-Seidel and multigrid iteration. Once the pressure has been updated, the Newton- Raphson iteration may continue to reduce the residual errors or the solution scheme may progress to the next time step. The stability and accuracy of this numerical scheme are discussed later. 3.2 An Explicit Time-stepping Scheme We derive here the second-order extension of the FTCS method. The time integrator uses a two- stage modified Euler method. The first stage is a standard Euler step: Ui,j,k = Ui,jk + ~tfiUjk Vi,j,k = Vi,j,k + ~tfi,j,k Wi, j,k = Al i, j,k + ~ t fi, l;k where the forcing terms are defined by: fi,j,k = -[~(UU) -[V(UV) -[Z(UW)" _~cpn + ~ V2un e fi,j,k = -(UP) -[V(VV)n-[z(VW)" _~cpn + ~ V2vn (28) fi,,;k = -(I) -[V(VW)n-[z(W~)n _~cpn + ~ V2U~n e We observe that the convective terms and the dif- fusion terms are all central differenced. Next we calculate the hydrodynamic pressure as a function of this first estimate of the velocities at the next time step: = +_(6Cun + [Vv" + [Zoo ) +2(6,cubvv + t'~u~u) +`cv~cw-~vubcv (29) _<cu bow - Rev [vu)) ' This equation for the pressure is the same as that for the implicit scheme (Eq. 26). Upon solving this equation, we proceed to the second stage: U' j+k = U',j,k + 2 (li,j,k + l''~'k ) V'n'+k = Vi,j,k + 2 (/i,j,k + fi,~,k ) (30) = Wi,j,k + 2 (li,7;k + fi,,;k ) ~ where the second-stage forcing terms are defined by: f ,j,k = -6x (U U)-[v(U V) _ [2C (u w) _~cp+ 1 V2U , Re (27) fi,j,k = -6~(U V)-[v(v V) _ [zC(v w) -6VP + R-V2v (31) fi,,;k = -b~ (U W)-6~(V W) [Z (w w) _`cp-+ 1 V2W e This two stage process gives a second-order accu- rate estimate of the velocities at the next time step. We conclude the algorithm by updating the pres- sure field: V2pn+1 = +_(6Cun + [Cvn + [Clan) (32) 740

+2(t,l~un+1 [Vvn+l + [cun+1 ~wn+1 +'!ivvn+1[Zwn+1 _ [Cun+1[Cvn+1 _~Cun+l ~cU~n+1 _ ~cvn+1 ~c2L~n+l) . The procedure repeats for the next time step. The stability and accuracy of the above algorithm is discussed later. 3.3 Numerical Solution of the Initial Conditions According to Helmholtz theorem, the initial ve- locity field may be expressed in terms of scalar and vector velocity potentials (Eq. 7~. The vor- tex problems that interest us most in this paper are completely specified by the vector velocity po- tentials. The finite-difference forms of the Poisson equations for the vector velocity potentials (see Eq. 9) are: V2¢ To = -(~:~0 yi, j,k v2 to = -((~0 ji, j,k V2?/JZO = -((`JzO)i,j,k where (~o~i,j,k is the discrete sampling of the initial vorticity field which is prescribed. Upon solving (33) for the potentials, we use a discrete form of Helmholtz theorem to calculate the initial veloci- ties: (Uo~i,j,k = [V¢Zo ~t,b~, (Vo~i,j,k = -Maze + taboo (34) (U)o~i,j,k = {ia:l/)Vo-[V¢~° ~ Formally, the initial vorticity field prescribed in (33) should be solenoidal, but in practice an initial vorticity field which is not necessarily solenoidal may be used. The primitive variable formulation allows this flexibility because any vector field for the initial vorticity will lead to a solenoidal veloc- ity field due to the Helmholtz decomposition. In fact, the vorticity field that is recalculated based on the definition of vorticity (10) is also solenoidal. 3.4 Numerical Implementation of the Boundary Conditions Our numerical simulations use periodic, free-slip, no-slip, arid free-surface boundary conditions. To illustrate how these boundary conditions are im- plemented, we assmne that one of these types of boundary conditions is imposed on the plane z=0. Let indices (i, j,K) denote a grid point on z=0, and (i, j, K + 1) an imaginary set of grid points above the plane z=0. The grid points on the k- axis are ordered k = 1,2,...,K. For Dirichlet boundary conditions the field equations are never evaluated on z = 0 and no special treatment is required for the field equations. For periodic and Neumann boundary conditions, expressions for the imaginary set of grid points in terms of real grid points are formulated and substituted into the field equations. 3.4.1 Periodic Boundary Conditions The field equations for the velocities, pressures, and stream functions are evaluated at the grid level (i, I, K), and the centered difference operators (22) are functions of the imaginary set of grid points. (33) We impose periodic boundary conditions by set ting the imaginary grid points at the grid level (i, j,K + 1) equal to the grid points at the first grid level (i, j, 1~. Similarly, we set the imaginary grid points at (i, j, O) equal to the grid points at (i, j, K). This algorithm gives second-order spatial approximations to the field equations without ad versely affecting vectorization (see Appendix B). 3.4.2 Etree-Slip Boundary Conditions As shown in Sec. 2.3, the horizontal velocities (u, v) and the hydrodynamic pressure p are even functions about z=0, the vertical velocity us is an odd function, and w=0 on the wall. These rela- tions may be used to assign the flow quantities at the imaginary set of grid points in terms of the same quantities at interior grid points: Ui,j,K+t = Ui,j,K-1 Vi,j,K+! = Vi,j,K-1 Eli, j,K+i = -Eli, j,K-~ Pi,j,K+1 = Pi,j,K-1 · (35) 74~

At the (i, j, K) grid level, these relations are back substituted into the horizontal components of the momentum equations and the Poisson equation for the pressure. The resulting finite-difference ap- prox~nations are second-order accurate. The value and behavior of the velocities on the free-slip plane are used to assign the boundary values and be- havior of the vector velocity potentials. For ex- ample, setting ?~b~o=0 and 1,bvo=0 on z=0 ensures that the no-flux condition is met. In addition we assume that ¢~0 and TO are odd functions about the free-slip plane, and assign ?tzo to be even (~`bzo/0z = 0~. At the (i, j, K) level these relations are substituted into the -Poisson equation (33) and the discrete velocity calculations (34~. This form of the vector velocity potentials ensures that the horizontal velocities (u, v) are even functions on the free-slip plane. 3.4.3 No-Slip Boundary Conditions Free-slip boundary conditions are used to set up the initial conditions at time t = 0, and at time t = 0+ we set u = 0 on the no-slip wall. As a re- sult, a sheet of vorticity is instantaneously formed on the no-slip wall. In our simulations of vortex rings impacting the wall, the rings are initially far enough away from the wall such that the strength of the vorticity on wall is less than one percent of the core strength. As shown in Sec. 2.4, the source term in the Poisson equation for the pressure is zero on a no-slip wall. For our numerical simula- tions we also set the time derivative of the dila- tion term equal to zero. The fact that u'* = 0 on the no-slip wall is used together with the Neumann boundary condition for the pressure (16) to assign the pressure at the imaginary grid point (Pi,j,K+iJ: Pi,j,K+! = Pi,j,K-~ + LIZ R Wi,j,K-~ · (36) This equation is then substituted into the 7-point star for the Laplace operator in (26, 29, or 32~. In the multigrid solution of the pressure we find it convenient to store the inhomogeneous velocity term (?l)i,j,K-~) in the vector used for the source strengths (see Appendix B). 3.4.4 Free-Surface Boundary Conditions Again, free-slip boundary conditions are used to assign the initial velocity field. Dirichlet conditions are specified for the hydrodynamic pressure since the pressure is known in terms of the free-surface elevation (see Eq. 11~. The kinematic condition for the free-surface elevation (12) may be evaluated either implicitly or explicitly. As an example, for the first-stage of our explicit scheme, we use: gi,j = q.,j + Atfi,, ~ ~ _ . . . _ ~ i,] wn _ un [c )7n- _ vn; K [ct7n (37) The evaluation of the momentum equations on the free surface depends on the imaginary grid points. Since the horizontal velocities are even functions about the free-surface (see Eq. 11), they are assigned using our free-slip wall method. The mass-conservation equation (3) is used to evaluate the vertical velocity at an imaginary grid point in terms of real quantities: U)',j,K+! = ~i,j,K-1-2^X LOU-2^Z {vv . (38) Finally, the z-derivative of the pressure in the mo- mentum equations is evaluated using a 3-point one- sided difference scheme: PA = 2 ~ z (3Pi, j,K-4Pi, j,K-1 + Pi, j,K-2 ~ . (39) This implementation of the free-surface conditions gives second-order spatial accuracy. In addition, the compact nature of the scheme does not affect the vectorization of either the semi-implicit or the explicit time integrators. 3.5 Solvability Condition for the Pres- sure As noted in Sec. 2, a solvability condition must be satisfied for the pressure if Neumann or peri- odic boundary conditions are specified on all of the boundaries. The numerical implementation of this integral constraint is discussed by Ghia, et al, (1977), which we extend to a multigrid solderer. Let 742

~i,j,k represent the source term in the discretized Poisson equation (see Eqs. 26, 29, or 32), and let ~ be the error of the integral constraint normalized by the fluid volume V. We use trapezoidal rule to evaluate the error: g=~v~vai~i~k . ., .,J,~ - ~ ~& ~ (~'Pn ji,j,k ~ , (~40~) ., j,k where l~v = Az~yAz and lea = Away, Waltz, or AyAz; [v = 1, 1/2, 1/4, or 1/8 if the grid point is respectively in the volume, on a side, edge, or corner of the computational domain; and similarly, I,, = 1, 1/2, or 1/4 depending on whether the grid point is on a side, edge, or corner respectively. Let ai.j,k represent a source strength which satisfies the integral constraint, then one possible solution for ably is simply: amok = ~i,j,k-~ (41) These modified source strengths ~i,j,k are used in- stead of the original source strengths ~i,j,k in the discrete Poisson equations for the pressure. We emphasize that the discrete Laplacian operator with all Neumann or Periodic boundary conditions is still singular, but ~.,j,k is in the same solution space as the generalized inverse. As a result, a point Gauss-Seidel scheme will converge to an op- timum solution for the pressure. Note that trape- zoidal rule is used to evaluate the integral con- straint because the trapezoidal rule is the adjoins of central differencing. Other numerical integra- tion schemes, even more accurate numerical inte- grators, will reduce the accuracy of the pressure solution. We also note that the numerical imple- mentation of the integral constraint has to be ap- plied at each level of a multigrid method, not just the level of the firmest grid. 3.6 Multigrid Solutions of the Field Equations The efficiency of our Navier-Stokes' codes is due to multigrid solution techniques. Multigrid is used to solve the momentum equations in the implicit time-stepping scheme (23), the Poisson equations for the stream functions (33), and the Poisson equations for the pressure (26, 29, and 32~. One obvious advantage of the multigrid method is that the convergence speed does not deteriorate as grid resolution increases in contrast to conventional it- erative techniques such as Jacobi or Gauss-Seidel (lIackbusch, 1985~. Indeed, the multigrid conver- gence rate and memory requirements, which are di- rectly proportional to the number of unknowns, are ideally suited for three-dimensional solutions of the Navier-Stokes equations. Multigrid uses fine grids to reduce the high wavenumber errors and coarse grids to reduce the low wavenumber errors. This interplay of fine and coarse grids is easily coded using recursive languages like PASCAL; however, FORTRAN can mimic recursive algorithms by us- ing array pointers in subroutine calling arguments. We provide a multigrid FORTRAN code for solv- ing the hydrodynamic pressure in Appendix B. 3.7 Convergence Tests Analyses and Error In this section we establish the stability and convergence properties of our implicit and ex- plicit time-stepping procedures. We consider three test cases: the a~cisymmetric stagnation problem, translation and diffusion of two-dimensional vor- tices, and attenuation of axisymmetric standing waves. The axisymmetric stagnation flow problem illustrates the rapid convergence to steady state of our implicit scheme. In addition, we show the spa- tial accuracy of our numerical scheme as a function of grid Reynolds number. The numerical simula- tions of two-dimensional vortices orbiting in a box test the temporal and spatial accuracy of our vor- tex formulations. Finally, numerical simulations of the attenuation of a~cisym~netric standing waves validates the formulation of the free-surface bound- ary conditions for low Reynolds numbers and small wave amplitudes. 3.7.1 Implicit Scheme for the Steady-State Problems Except for the treatment of the pressure, our semi- implicit scheme is directly related to the Briley- McDonald (1973) method for the compressible Navier-Stokes equations. Both schemes are first- order accurate in time and second-order accurate in space. A linearized stability analysis shows that 743

both schemes me ~conatlon~y stale. The an Clarence is that instead of a block-trl~ag~ solver, we use a ~tl~ld solver far the ~men- t~ e~atlons. As a test of on se~-l~llolt schema we con- ~der the ~sy~etrlc stagnation Dow problem far Cab ~ exact solution to the Navler-Stokes Cations ~ known (Sc~chtlug, 1968~. ~ let ad ~ respectively denote the ramp (~) ad ~ (a) velocltles. Ibe problem is nor~- hed by a length I=1 ad the velocity at 1~- lty ~ = e(r = 1~= -~=1. Two sets of re- s~ts, one using ~ ~sy~etrlc version of on scheme (see Appends A) ad one uslug the ~ t~e~mendon~ ~r~l~ (wltbout assuming a~sy~etryj are provided. The co~utatlons me performed with known Dldc~et con~tlons on the boundless ad zero lastly lilted velocity ad the iterations are continued Eta the = change in the velocity between iterations is less than a tolerance of ~=10-e. Par both computer codes the parameters Each control the namer of sweeps over the West grid level me as Buoys: number of Newton-Raphsou lteratlons, al; namer of ~tlgrld lteratlons, ^~'=6; namer ofOauss Seldel lteratlons, =6; ad namer of Jacot1 iterations ~=2. par ~ on n~rlc~ studies, the namer of Jacobi iterations is Amps Wed at Waco = 2.) ~ use a constant grid spac- lng gab A=2-5, where ~ = Ar = 6z ~ the grid size. The co~utatlon~ domain has mat belybt ad mat radius (or b~wldtbj, ad the d~ step is ~=1. To amaze lastly transients ~d to 1~ prove convergence, the viscosity is decreased Un- eaIly hom an lastly We of ~=1 to the prescclbed value after ~amp=5 time steps. Fugues (1) show the =~ Solute error in the veloclt1es ~ ad ~ as a Colon of grid Reynolds namer H~ = /. Ibe plots mea- s~e the resolution of the w~ bo~d~y layer as R~ increases. The ~sy~trlc code performs sUghtly better than the t~ee-~nslon~ code possibly because the t~ee-amenslon~ co~uta- Uon~ domain ~ lager tab the ~sy~etrlc do- m~n by a factor of ~ The errors far both methods Sib relay as ~ decreases below ~2. For R~ ~1 the ~ nor~zed error is 0~10-~]. ~ evaluate the eminency of these codes, ad as a guide to later co~utatlons, we seek the optE ma set of co~utatlon~ patters (as a ~c don of JO far ~ tote eeratlons. In the luterest of economy only the ~sy~etHc code used alto the e~ectatlon that resets far the t~ee-~nslou~ code wlD be sl~l=. ~ ~clD- tate come-on far these ~tl~ld co~utatlons, we Deane a Waste of the Per of operations Fir to be the extent Per of G~ss-Seldel sweeps over the Pest geld level far the entire set gowns (~,e1 p) ~k (1) sbows partly ret sits anew the optl~sj of ~ extensive v~da- don Judy. Note that Ben a ram down of the vlscosHy is not sapped (~omp=°~, the solution may averge even though the l~0cH scheme is If. Fanny' we show the peace of the Build scheme as a Colon of the er of grid poluts. has is shown in Flame (2) far the case ~=1 filth the ~' patters Beef= l l ^~1 =21 "~=l1 Ad= 51 ad ~a=~=41 ad keeping the sax tolerance of ~=10-a. The namer of grid poluts To' is increased by Carry decreasing the geld size ad keeping the sax domain ~men- dons. Note that the log-log plot of Bar versus ~ has ~ average slope of ~ 0.3. has locates that the rate of convergence ~ the m~tlgrld solver only sUghtly dower tab a anew Reckon of the namer of unknowns. 3.~.2 The Expl~lt Scheme for Onsteady Problem Far the time-dependent problem we use a second- order pre~ctor-corrector (R~ge-Kutta) method far Ume 1ntegratlon ad centered Fences far the spittle derlvatlves. 0= prlnclp~ concerns here me fast stablDty ad second accuracy of the time ln- tegratlon. Ibe Une~lsed st~lUty ~ysls of the Navler- Stokes e~atlons ~ usually perjured by assay lag pedo~c bo~d~y con~tlons1 lynorlng the (Un- e~) pressure gradient terms, ad Une~~lng the c~vectlve terms to Stan the anew advectlve ..~ , . mnuslon echelon: ~ + ~ Via= EVE where # is assumed to be constant ad prescclbed. For one-~menslon~ problems, tb6~t~lDty crNerla far the (~st-order) ~rw~d-tl~ centered-space (FIGS) scheme ~ wed known (egg Roam 1976, 744

,B2 < 2a~ < 1, (42) where c'C=~At/~x2 and ,B==UAt/Az. For more spatial dimensions and higher-order time integra- tions, the generalization is in principle straightfor- ward but algebraically quite complicated. For two dimensions, a generalization of the following form has often been suggested and used (Fromm, 1964; Miyata & Nishimura, 1985~: (~z + pg)2 ~ 2(a~ + av) < 1 , (43) which can be derived by essentially letting Ax = Ay and assuming equal phases in the two dimen- sions. Unfortunately, (43) is neither necessary nor sufficient for stability in general. The latter be- ing obvious since (43) does not even guarantee the one-dimensional condition (42~. By not assuming equal phases, the correct analysis is more involved and the final results can be sum- marized as follows. Without loss of generality, eve assume U if 0 and define 7 = V/U. For simplicity we denote the left inequality of (42) applied to z and y by IS and I', respectively, the right inequality of (43) as II, and the phrase 'necessary and suffi- cient for stability' by 'iff'. For 7 < ~-1: II is iff for Rc = U^z/z, < 2~ and In is iff otherwise. The symmetric case is for 7 > (~-1~-~: II is iff for Rat < 2~/, and I', is iff otherwise. For inter- mediate values of ~-1 < 7 < (~-Phi, II is again iff for Rx < 4/~1 + 7~. For Rx > 4/~1 + 7), the results depend on whether 7 is greater or less than 1. For 7 <, > 1, I~,I', eventually become established as iff for large Rat, while for 7 = 1, I=, Iv are never sufficient. The overall solution is rather complicated and is depicted graphically in Figures (3) where the relevant quantities involved in the inequalities (eqts. 42 & 43) evaluated at the maximum At for stability are plotted as a function of Ret, for three values of 7=.7, 1.0, and 1.4. For the second-order Runge-Kutta scheme in time, the resulting (complex) amplification factor is sim- ply: G = 1 + f + f2/2, (44) where f is the amplification factor for the first-order (FTCS) difference scheme. The analysis is some- what more involved but qualitatively unchanged. In practice, the maximum value of At is computed numerically using the best estimates for U. In view of the linearization, the elimination of the (linear) pressure term, and the simplified boundary condi- tions in the stability analysis, more conservative values for the time step are often found necessary in actual computations. We now turn to the accuracy of our second-order (Runge-Kutta) explicit scheme. We consider two problems for which analytic solutions (for the ide- alized cases) exist. The first problem is the decay of a two-dimensional rectilinear Gaussian vortex tube. Since the free decay of this vortex in an unbounded domain is relatively trivial and does not involve the nonlinear convective terms, we con- sider instead the vortex tube placed asymmetri- cally inside a rectangular domain. Free-slip con- ditions are applied on the boundaries of the do- main. Because of the presence of the image vor- tices, the tube orbits within the box. The core vorticity should, however, decay at a rate close to that for unbounded domain, at least for small core sizes compared to the dimension of the domain. Specifically, we consider a 4 x 4 box with a Gaus- sian vortex of core radius rc=0.5 placed initially one unit inboard of a side's midpoint. The ini- tial circulation of the vortex is I' = or. Equal grid spacings with ~x=Ay=4/25 are used, and a grid Reynolds number, R,`=UoAx/l, with UO-1 is again defined. Figures (4) plot, as functions of R^, the relative error in the maximum vortic- ity, Lyman-Wexa I /We ra ~ where (Oman is the maxi- mum calculated value for the vorticity and We:ca = I`/~7rtr2 + but)) is the analytic value in an un- bounded domain. The errors are plotted at two times, t=1 and 10. The corresponding relative er- ror in conservation of total energy are show in Fig- ures (5~. Three time step sizes, At=.01, .02 and .05 are considered for which the minimum R.` for stability are respectively R's=.32, .64, and 1.6. As seen in Figure (4a) the peal' errors in the vortic- ity occur at Rig = 2 for t=1. However, the energy at this point (5a) is well conserved probably due to the conservation-law form we use to formulate the Navier-Stokes equations. Here, we have defined the relative error in the energy conservation formula as ec = Hefty-E(O) + D(t)~/E(O), where E(t) is the kinetic energy and D(t) is the energy that is dissipated (the Reynolds number term in Eq. 17 integrated over time). For grid resolutions lower than Ret 2 2, one expects an increase in error due to oscillations (Roach, 1976), this is indicated in 74s

both our steady-state and unsteady calculations. The errors, however, are not intolerable and may at times even decrease as Rig increases above two. Finally, we consider a time-dependent problem with a free surface. We study the decrease in wave amplitude due to viscous dissipation of an axisymmetric standing wave in a cylindrical basin. For simplicity, we choose the tank radius and the tank depth R`=D,=1 and the standing wave pe- riod T=1. For this axisymmetric problem, we use a 26 x 26 grid, and At=.005. The Froude number is fixed at Fr=.5 and consistent with the free-surface linearization, a relatively small wave slope of A(0)k=.1, where k = TO ~ 3.832 is the wavenumber corresponding to the first-mode standing wave, and A(t) is the (centerline) wave amplitude. The analytic solution for this problem is dAe~a/dt=exp(-2uk2t) (Lamb, 1932~. For com- parison, we integrate for one wave period (t = T) and compute the relative errors, c~ = ~A(T) Ae~a(T)~/A(0) and 62 = ~A(T)-Ae~a(T)~/(A(O)- Ae~a(T)), as functions of the grid Reynolds num- ber Ret = R2/uT. These results are shown in Fig- ure (6~. The relative errors are all less than four percent. For Rig ~ 4 the relative errors increase due to low damping at low spatial resolutions. For Rig < 4 the errors increase as the diffusion stability boundary is neared. 4 Numerical Stuclies The present numerical methods are used to study ring vortices and vortex tubes impinging on walls and free surfaces. The vortex ring studies are performed using an axisymmetric version of the program (see Appendix A). For the vortex tube studies the full three-dimensional code is used. Schematic pictures illustrating the asymmetric ring vortex and the three-dimensional vortex tube simulations are given in Figures (7) and (8) respec- tively. The actual physical and numerical param- eters used in the simulations are provided in three tables: Table (2) for ring vortices interacting with walls, Table (3) for ring vortices interacting with free surfaces, and Table (4) for vortex tubes inter- acting with walls and free surfaces. In later discus- sions of the numerical results the figure headings refer to the detailed information that is provided in the tables. 4.1 Ring Vortices Impinging on Walls and Etree-Surfaces We study the normal incidence of vortex rings with no-slip walls, free-slip walls, and free surfaces. For the no-slip wall cases the radial expansions of the primary vortex rings are stopped by the growth of the boundary layer on the wall. As in the expert meets of Walker, et al, (1987) with a solid wall, we too observe that at sufficiently high Reynolds num- bers, secondary and tertiary vortex rings are sepa- rated from the boundary-layer on the wall. In ad- dition, our numerical simulations show that these secondary and tertiary vortices form their own boundary layers on the wall as they are pinched off from the wall. Once they are formed, the sec- ondary and tertiary vortices are wrapped around the primary vortex. In the later stages of the flow, the secondary vortex is thrust through the center of the primary vortex ring whereupon it merges with the wall boundary layer. Similarly, when the wall is replaced by a free sur- face, our numerical simulations show that sec- ondary vortex rings are shed from the surface at intermediate Froude numbers (F2 > .25~. At lower Froude numbers, free surfaces appear to behave like free-slip walls. The ring expands radially out- ward, the core diameter is reduced due to stretch- ing and the ring cross-section develops a distinctive head and tail shape. These low Froude number re- sults agree qualitatively with the observations of Bernal, et al, (1989) and Kwon (1989) for vortex rings interacting with a clean free surface. These conclusions regarding the interactions of vortex rings with walls and free surfaces are based on the analyses of 12 numerical simulations. For these computations, a Gaussian core distribution is used to specify the initial vorticity field of the vortex ring: silo = `~'c expt_ (r rO) + (x _ zo)2 ~ `45' where we is the peals vorticity, (rO, zO) denotes the center of the core, and 7.c is the core radius. We also assign an image vortex about the centerline in order to make the initial vorticity field solenoidal. The Reynolds and Froude numbers are based on the length scale r0 = 1 and a centerline velocity scale Uc _ I`/~,rrO) = 1. The Reynolds numbers vary from 200 to 400, and the Eroude numbers are 746

less than one. The vortex rings we choose to simu- late are thick with a core to ring radii ratio of 1/2. As indicated in Figure (7) and Tables (2) and (3), the initial conditions are set up so that the vor- tex ring rises up towards the wall or the free sur- face and then spreads out. Free-slip boundary con- ditions are unposed on the remaining boundaries of the cylindrical domain. The numerical parame- ters which control the number of solution iterations over the grid are chosen so that the average diver- gence over the grid points is maintained to within a tolerance of 10-6. Figures (9 & 10) compare the trajectories of ring vortices interacting with free-slip and no-slip walls. The vortex cores are initially circular in both sim- ulations in accordance with the choice of Gaus- sian core distributions. The vortex cores then evolve into elongated shapes reminiscent of Nor- bury's (1973) steadily translating family of vortex- ring solutions. When the vortex rings are within one ring radius of the wall, the two numerical simu- lations begin to differ. The most striking difference is that the vortex ring continues to expand radially in the free-slip wall case, whereas in the no-slip wall case the vortex ring is stopped at about r ~ 2.25 due to a boundary layer forming on the wall. The final positions of the vortex core in the no-slip wall case agree qualitatively with the vortex ring ex periments of Walker, et al, (1987) who observed radial expansions of r ~ 2. Also in agreement with experiments, our numerical simulation shows that the ring vortex rebounds from the no-slip wall. In the free-slip wall case, as the ring expands ra- dially outward, the core diameter is reduced due to stretching and the core's appearance develops a distinctive head and tail shape. This shape of the vortex core agrees qualitatively with the very initial stages of the experiments of Bernal, et al, (1989) for vortex rings interacting with a clean free surface. Sheriff, et alms (1989) inviscid studies also show the development of such head and tail shapes. As an aside, they End that the deformations of the vortex core in the free-slip wall case are primarily responsible for the noise that is generated in the head-on collision of two vortex rings. Figures (11, 12, & 13) show several stages in the interaction of a ring vortex with a no-slip wall. As illustrated in Figure (11) the boundary-layer has just separated from the wall due to the adverse pressure gradient formed by the primary vortex ring. Here, we refer to the original vortex ring as the primary vortex ring, and any vortex rings that are ejected from the boundary layer of the primary vortex as secondary and tertiary vortex rings de- pending on when they are formed. The secondary and tertiary vortex rings have a sense of rotation that is opposite to the primary vortex ring, but as we observe in Figure (11), the separated vor- ticity may itself form a boundary layer which has the same rotation as the primary vortex. This for- mation of multiple boundary layers appears not to have been observed in laboratory studies of vortex rings. We note that the secondary vorticity induces a ve- locity on the primary vortex which causes the re- bounding that is observed in experiments (see Fig. 10~. At this Reynolds number of Rc = 400, the pri- mary vortex ring is sufficiently strong to pinch off the boundary-layer vorticity as illustrated in Fig- ure (12) to form a secondary vortex ring. This secondary vortex ring orbits the primary vortex ring due to the velocities induced by the primary vortex on the secondary vortex ring. In addition, we observe that boundary-layer vorticity is being wrapped around the primary vortex ring such that a tertiary vortex may form as shown in Figure (13~. Here the secondary vortex ring has been thrust through the core of the primary ring where it merges with the boundary-layer vorticity. The formation of secondary and tertiary vortex rings has been observed in the experiments of Walker, et al, (1987~. Unlike the boundary-layer theory they use, however, the present direct simulations of the Navier-Stokes equations permit the actual model- ing of the secondary and tertiary vortex rings. In Figure (14) we show how well energy is con- served by our numerical simulations. Observe that the kinetic energy is reduced to about 10% of its initial value in the no-slip simulation and 30% of its initial value in the free-slip simulation due to diffusion. (Note that the duration of the free-slip simulation is less than that for the no-slip wall be- cause we have stopped the former before the vortex ring hit the outer boundary.) In both numerical simulations the energy is conserved to within 2%. Interestingly, the rates at which the kinetic ener- gies are dissipated are only slightly different be- tween the two simulations. In fact, if the free-slip simulation were continued further, the energies of the two cases may become close again. We hope to investigate this phenomenon in a future paper. In Figure (15) we compare high and low resolu 747

lion simulations of a vortex ring interacting with a no-slip wall. The high resolution simulation is twice as dense as the low resolution run. The agreement that exists between the two simulations indicates that we have adequate resolution of the wall boundary-layer, where the peak vorticity mag- nitude reaches a value five times greater than its initial value in the primary vortex ring. In Figures (16) we report results for the normal incidence of a vortex ring with a clean free-surface at a low Froude number of F2 = .125. Figure (16a) plots the free-surface elevations, and Figure (16b) plots the motion of a single contour of vor- ticity. For this numerical simulation, the mean free-surface elevation is not always correctly pre- dicted to be zero. For time t < 5 the spurious mean free-surface elevation is less than 10% of the peak free-surface elevation, and for t > 5 the mean free-surface elevation is the same magnitude as the peak elevation. We note that the free-surface dis- turbances are quite small, and in fact the attenua- tion of the kinetic and potential energies are con- served to within 5%. The results we present in Figure (16a) have the mean free-surface elevation subtracted out so that we may focus on the fea- tures in the free-surface elevation that we believe are correctly predicted. Observe that the evolu- tion of the vortex core at this low Froude number is very similar to the free-slip case shown in Figure (9~. The most significant difference is that the tail is now slightly longer. The plots of the free-surface elevations indicate that a wave trough is located above the low pres- sure region in the vortex core and a wave crest near the front of the core where the pressures are high- est. The free-surface elevation decreases in am- plitude as the disturbance moves radially outward due to mass conservation, but no waves appear to radiate away from the vortex ring. The wave troughs that travel with the vortex ring are simi- lar to the free-surface depressions observed in the vortex-pair experiments of Willmarth, et al. (19RR) and Sarpkaya, et al, (1988~. In Figures (17, 18, & 19) we report results for the normal incidence of a vortex ring with a clean free- surface at an intermediate Froude number. The magnitudes of the predicted free-surface elevations are beyond the restrictions imposed by the approx- imations used to derive the linearized free-surface boundary conditions (11 & 12), but the gross fea- tures of the flow should be valid. For this numerical simulation, the attenuation of the kinetic and po- tential energies are conserved to within 15%, which is not unacceptable for the present purposes. We observe in Figure (17a) that the free-surface eleva- tion reaches a maxanam height above the center of the vortex ring that is about half the initial core radius, but once again no waves appear to radi- ate away from the vortex ring. From Figure (17b) the motion of the core is similar to the case of the no-slip wall (see Fig. 10~. This is especially clear in Figures (18 & 19), where we observe the sec- ondary and tertiary vortex rings shed by the free surface. These numerical simulations show that vortex shedding can occur at intermediate Eroude numbers even if the free surface is clean. 4.2 Three-dimensional Vortex Tubes Impinging on Walls and Free Sur- faces In this section we study the interaction of vortex tubes with no-slip walls and free surfaces. Specifi- cally, we model two problems: (1) the interaction of a pair of trailing tip vortices with boundaries, and (2) the horizontal translation of a vortex tube slightly submerged below a free surface. As ob- served by Sarpkaya (1986), an interesting aspect of these problems is the reconnection of the vor- tex tubes with the boundaries. Our numerical simulations illustrate two possible mechanisms for the reconnection process: (1) reconnection of the primary vortex tube; and (2) reconnection of the secondary vorticity that is generated in the free- surface boundary layer and swept up by the pri- mary vortex. We hypothesize that reconnection of the primary vortex tubes may occur for deeply submerged trail- ing vortices that have enough time to develop large undulations prior to their impact with the free sur- face. The connection process of the primary vortex with the free surface is similar to the formation of crude vortex rings that is observed in the far walces of aircraft (see Crow, 1970~. Instead of the trail- ing vortex connecting with its neighbor, however, in the present case, the trailing vortex connects with its image above the free surface. A unique characteristic of the primary vortex reconnection is the cusp pattern that is traced by the normal vorticity on the free-surface. If the trailing vortices are very deeply submerged, they may connect with each other prior to contact with the free surface 748

J and the vortex reconnection with the free surface would be different. Two sources of boundary-layer vorticity which may connect with the free surface include shear stresses induced by surfactants and variations in the free- surface elevation at high Froude numbers. Thick sheets of helical vorticity spiral off the primary vor- tex and unevenly pull at the free-surface boundary- layer in the axial direction. The origin of the sheets of helical vorticity appears to be the result of a he- lical instability that is initiated by large changes in curvature along the anus of the primary vortex tube. Self-induced straining flows may help to strip the tightly wound helical vortices off of the primary vortex tube. Figure (20) provides a schematic of the helical vor- ticity spiraling off of the primary vortex tube. As the helical vortices are swept around the primary vortex, they attach themselves to the free-surface boundary layer. The free-surface boundary-layer separates in front of the primary vortex tube and the helical vortex sheets evolve into narrow ridges of cross-axis vorticity that are wrapped around the primary vortex and ride on top of the secondary- vorticity that is shed by the boundary layer. The separated boundary-layer feeds into the undersides of the ridges, and the ridges swell into U-shaped vortices as shown in Figure (21~. The bases of the U-shaped vortices are located close to the free sur- face in the region where the free-surface boundary- layer separates. The bases are either free or an- chored in the secondary vorticity depending on their energy level, and the undersides of the tips of the U-vortices are anchored in the secondary boundary-layer vorticity that is wrapped around the primary vortex. The vortex lines enter the tip of one leg of the U-shape, move down to the base of the U-shape, and exit the tip of the other leg. The U-vortices are large in size at their bases, and the diameters of the legs of the U-vortices be- come slightly smaller near the tips possibly due to stretching around the primary vortex. The ori- entation of the U-vortices permits reconnection of normal vorticity with the free surface at the bases of U-vortices and the tips of the legs. We conjec- ture that reconnection may take place as the bases of the U-vortices open up when they come into contact with the free surface. If our conjecture is correct, the reconnection of the bases of the U- vortices with the free surface should manifest itself as a pair of closely spaced dimples on the free sur face, and these dimples should appear on the out- board side of the pr mary vortex. Similarly, if the primary vortex is sufficiently strong, the two tips of the U-vortices may be wrapped into contact with the free surface. As a result, we may expect more reconnections in the region where the tips of the U- vortices feed into the secondary vorticity shed by the boundary layer. The reconnection of the two tips should also manifest itself as one or two pairs of dimples on the free surface, but the tip dimples will most likely be inboard of the dimples formed by the bases of the U-vortices. In addition, the tip dimples may be smaller than the dimples formed by the thicker bases of the U-vortices. Our conclusions regarding the interactions of vor- tex tubes with free surfaces are based on analy- ses of 7 three-dimensional numerical simulations. The geometry definitions are provided in Figure (8), and the numerical parameters of the simula- tions are provided in Table (4~. A Gaussian core distribution is used to specify the initial vorticity field as follows: (Rio = cen Zcen (e- ICen)2 + (Z-Zcen)2 tic expt- r2 ~ TO + Camp COSt w Y) ZO + Ham p COST W Y) , (46) where tic is the peak vorticity, (zO, zO) denotes the mean position of the core center in the ~z-plane, (~amp, Camp) are sinusoidal perturbations applied to the core position, W. is the width of the tank, and rc is the core radius. We use free-slip bound- ary conditions on the bottom of the computational domain and on the walls in the :tz-planes. For the first six of the simulations we use free-slip bound- ary conditions also on the yx-planes, and for the seventh run we use periodic boundary conditions along the ~-axis. On the top boundary either no- slip or free-surface boundary conditions are em- ployed. As indicated in Figure (8) and Table (4), the ini- tial conditions for the first six vortex tube studies are meant to simulate the rise of a pair of trail- ing vortices up towards either a wall or free sur- face. For these six simulations we assign 3 image vortices in the lower corner of the computational domain where the tubes are initially located. The 749

final computer run is meant to simulate a vortex tube that has risen to the free surface with a large sinusoidal disturbance along its length. For this simulation we assign an image vortex above the free surface, and use periodic boundary conditions along the :~-axis to allow the vortex tube to trans- late parallel to the free surface. For all the computer simulations, the Reynolds and Froude numbers are based on the mean half span a/2 = 1 and centerline velocity scales Uc-- 2It/~,r a)= 1. The Reynolds numbers vary from 200 to 400, and the Froude numbers vary from O to 4. Note that the core radii are very thick for three of the simulations (rc=0.5), but the other computer runs have more realistic core radii of rC=0.25. For these simulations we choose to focus on the recon- nection process itself and not the mechanism by which the vortex tube becomes initially disturbed. As a result, we use a width for the computational domain that is too short to permit Crow instability. As before, the numerical parameters which control the number of solution iterations over the grid have been chosen such that the average divergence over all the grid points is less than a tolerance of 10-6. A typical three-dimensional simulation consisting of 653 grid points and 3,000 time steps takes about 3 hours on a Cray Y-MP. Figure (22) shows the cusp pattern that is formed on the free surface by the reconnection of primary vorticity with the free surface at a low Froude num- ber. A single vortex tube with a sinusoidal pertur- bation is initially submerged slightly below a free surface. The actual physical and numerical param- eters are provided in Run 7 of Table (4~. The por- tion of the vortex tube that is closest to the free surface connects with the free surface. At this low Froude number (F2 = .25), the free surface acts like a free-slip wall, and the image of the vortex tube above the free surface induces high horizon- tal velocities in the region where the primary vor- tex tube is connecting with the free surface. The connecter region moves faster than the deeply sub- merged portions of the primary vortex tube so that as the vortex tube brealts apart, a cusp pattern is formed. Although we are actually plotting the nor- mal vorticity on the free surface in Figure (22), the free-surface elevations are directly related to the normal vorticity according to the solution of Rank- ine's combine vortex (Streeter, 1948~. Thus, we may expect to see a pair of shallow dimples on the free surface that split apart as the primary vortex tube connects with the free surface. In fact, this phenomenon is observed in our numerical results, but the free-surface elevation is contaminated by initial transients. For this numerical simulation the attenuation of the kinetic and potential energies is conserved to within 2% relative to the initial en- ergy. Figures (23a-e) show the evolution of the primary vortex tube as it reconnects with the free surface. Figure (23a) shows the vortex tube a short time af- ter the numerical simulation began. The isolated structures that are not attached to the vortex tube are remnants of helical vorticity that have been torn off of the vortex tube. Typically, these rem- nants have very little kinetic energy relative to the vortex tube. In Figures (23b & c) we see that one end of the vortex tube has developed a head and tail shape that is very similar to the shape formed as a vortex ring approaches a free-slip wall (see Fig. 9~. This phenomenon has also been observed by Stanaway, et al, (1988~. When we consider the vortex tube and its image above the free surface as one unit, we realize that the head and tail fea- ture forms a dipole as pointed out by Melander & Hussain (1988~. The next phase of the recon- nection process is known as bridging in the nota- tion of Melander & Hussain. During this phase, which is illustrated in Figure (23d), vortex lines connect with the free surface in the head region of the dipole. The flat section of the tail in Figure (23d) is composed of vorticity that is normal to the symmetry (free-slip) boundary in the yz-plane and parallel to the free surface in the ~y-plane. Sims ilarly, the bridge, which is the raised cylinder on top of the flat tail, is composed of vorticity normal to the free surface and parallel to the symmetry boundary in the yz-plane. Since vorticity can only be normal to a free-slip plane, the bridge does not touch the yz-plane and the tail does not touch the free surface. This explains the slot above the flat tail Ad the gap on the side of the bridge. In Fig- ure (23e) the reconnection with the free surface is complete, and one end of the vortex tube is con- nected to the free surface and the other end of the vortex tube is connected to the symmetry bound- ary in the back of the computational domain on the yx-plane. Hirsa, et alts (1990) measurements and two- dimensional computations indicate that the inter- actions of a vortex pair with a contaminated free- surface are similar to those for a no-slip wall. This 7sO

provides the motivation for us to consider the sim- ulation of the interaction of a pair of trailing tip vortices with a no-slip wall because such simula- tions should allow us to understand the most im- portant features of vortex tube interactions with contaminated free surfaces. Thus we consider a deeply submerged vortex tube with centerplane symmetry to model the rise of a trailing vortex pair toward a no-slip wall as illustrated in Fig- ure (8~. The initial position of the vortex tube is perturbed sinusoidally according to Equation (46~. The Reynolds number is Re = 400, and we choose a initial core radius he = .25, which according to Spreiter & Sacks (1951) is the right order of mag- nitude for an elliptically loaded wing with span ~ = 2. The other physical and numerical parmne- ters are provided in Run 3 of Table (4~. Figures (24a-g) show the interactions of the vor- tex tube with the no-slip wall. We plot a single isosurface of three relevant quantities at different instants of time: the magnitude of the helical vor- ticity (~/~=0.5), the magnitude of the axial vorticity (~w',~=0.5), and the magnitude of the to- tal vorticity (~=0.5~. Figure (24a) shows the initial unwinding of the he- lical vorticity. An outer sheet of helical vortic- ity covers an inner sheet. At earlier time steps, which correspond to the onset of a helical instabil- ity, the sheets of helical vorticity are very tightly wound around the primary Vortex tube. Dritschel's (1989) two-dimensional simulations suggest that self-induced straining flows may strip the helical vortices away from the primary vortex. The plot of the actual vorticity in Figure (24a) indicates that the initial perturbations in the primary vortex tube have been eliminated. According to Leonard (1985) the asocial flows that are induced by the he- lical vortices may cause this straightening of the primary vortex. We see that the formation of the wall boundary-layer is also visible in the plot of the axial vorticity. The plot of total vorticity in Fi~- ure (24a) shows how the helical vortex sheets wrap around the prunary vortex tube. Figures (24b) show more unraveling of the helical vortex sheets fi am the pranary vortex tube. As the helical sheets unravel they sometimes attach them- selves perpendicularly to the symmetry boundary on the centerplane. Figures (24b) also show that the boundary layer on the wall is being affected by flow fields induced by the helical vortices. Note that the helical vortex sheets have become slightly narrower relative to the preceding figure. In Figure (24c) the helical vortex sheets have merged with the wall boundary layer, which has separated from the wall. The attachment of the helical vortex sheets to the wall boundary layer is a very important mechanism for increasing the low kinetic energy of the helical vortices. We base this conclusion on our tree-surface simulations which show that the helical vortices rapidly diffuse if no energy is fed into them. In contrast, the no-slip wall simulations indicate that the helical vortices gain energy from the wall boundary layer. Figures (24d) show that one sheet of the helical vorticity has evolved into a narrow ridge of cross- axis vorticity that is wrapped around the base of the primary vortex tube. On the other hand, por- tions of the top sheet of helical vorticity appear to be weakening. A slot has opened up between the narrow ridge and the symmetry boundary because the cross-axis vorticity is parallel to the symme- try boundary and cannot connect with the sym- metry boundary. However, the secondary vorticity near the base of the ridge is normal to the symme- try boundary, and connection with the symmetry boundary occurs. The plot of axial vorticity in Fig- ure (24d) shows another boundary layer forming behind the separated secondary vorticity. The sep- arated secondary vorticity also induces a velocity on the primary vortex core that causes a rebound. This phenomenon is more clearly seen from other perspectives (not shown). Figures (24e) show more growth of the ridge at the base of the primary vortex, and the cross-axis vor- ticity at the top of the primary vortex appears to be strengthening. The boundary~layer that is gen- erated by the shed secondary vorticity is clearly visible in Figures (24e). A curved sheath of sec- ondary vorticity is being wrapped around the pri- mary vortex. Figures (24f) show two distinctive ridges of cross- ~ris vorticity wrapped around the primary vortex tube. The plot of asocial vorticity shows that a fat secondary vortex is begirm;ng to form at the top of the vortex tube. At the bottom of the vortex tube, a filament of secondary vorticity has broken away from the sheath of secondary vorticity that is wrapping around the primary vortex. The plot of total vorticity in Figure (24f) illustrates that the vortex lines move down the sheath of secondary vorticity, feed into the ridge of cross-~is vorticity 7s~

from underneath, move parallel along the ridge of cross-ems rorticity, and exit through the filament of secondary vorticity. Figures (24g) show that the ridges have swelled into tubes of cross-a=e vorticity. As illustrated in the plot of axial vorticity, the secondary vortex is clearly attached to the sheath of secondary vortic- ity, whereas the filament of axial vorticity at the base of the tube hangs free. In Figure (25) the image of the ~=0.5 isosurface of vorticity is re- flected about its n~dspan to emphasize the struc- ture of the U-shaped vortex tubes that are wrapped around the primary vortex. Although our no-slip wall simulations do not per- mit the connection of normal vorticity with the wall, we can speculate as to how normal vorticity would connect with a contaminated free surface. As illustrated in Figures (24g & 25), the tips of the U-shaped vortex tubes are wrapped around the primary vortex to the point where contact may oc- cur with the free surface. In addition, as the sec- ondary vorticity is wrapped around the primary vortex core, the velocities induced by the secondary vorticity on the primary vortex can cause the pri- mary vortex to collide with the wall again. This reasoning is supported by the two-dimensional sim- ulations of vortex pairs by Orlandi (1990), which show multiple rebounds of the primary vortex and multiple ejections of secondary vorticity. The same flow field which causes the primary vortex to im- pact the wall a second time should also propel the U-vortices against the free surface. Once the U- vortices are near the free surface, the bases of the U-vortices can open up to form two closely-spaced dimples on the outboard side of the primary vor- tex. Similarly, the tips of the U-vortices, where they connect with sheath of secondary vorticity, can also open up, and they too will form two pairs of closely-spaced dimples on the free surface. In Figure (26) we show how well energy is con- served for the no-slip wall case. We observe that the kinetic energy is reduced to about 50% of its initial value, and the attenuation of the kinetic en- ergy is conserved to within 3% relative to its ini- tial value. As mentioned earlier, the average di- vergence over all the grid points is maintained to within a tolerance of 10-6. This evidence together with earlier validation studies confirm the reliabil- ity of the numerical simulation. As already shown in Figures (24a-c), the longitu dinal structure of the helical vorticity involves in- ner and outer sheets of vorticity that are tucked inside of each other. These sheets of helical vor- ticity would tug at the free-surface boundary layer in both directions along the axis of the primary vortex tube. Figures (27a-c) show the cross-track vorticity evaluated on a no-slip wall at different in- stants of time. These figures show a banded struc- ture, especially when one considers the symmetry conditions that are imposed at the ends of the y- axis. Our numerical simulations of vortex tubes interacting with clean free surfaces at low Froude numbers also develop helical vortices, but the free surface does not form any striations (banded struc- tures). This suggests that the striations that are observed in laboratory experiments (see Sarpkaya & Henderson, 1984; and Sarpkaya, 1986) are prob- ably the result of helical vorticity interacting with a contaminated free surface. We hypothesize that a trailing vortex pair sweeps the free surface clean of surfactants at different rates in the central re- gion of the vortex pair due to the influence of axial flows induced by the helical vorticity. As a result, the concentration levels of the surfactants will have a banded structure that is normal to the track of the trailing vortices. Furthermore, we expect this sweeping phenomenon to be strongest during the initial interactions of the vortex pair with the free surface because the helical vorticity is strongest at that moment. This conjecture regarding the strength is supported by the sequence of contour plots we provide in Figures (27a-c), which show that the cross-track vorticity attenuates with time due to viscous effects. In addition, we may also expect to see dimples at the ends of the striations because the helical vortex sheets ultimately lead to the reconnection of normal vorticity with the free surface. In summary, we believe that the essential stages of the reconnection of secondary vorticity with the free surface are as follows: . Generation of helical vortex sheets by the pri- mary vortex tube due to the onset of a helical instability. · Stripping of the helical vortex sheets due to self-induced straining flows. · Attachment of the helical vortex sheets to the separated free-surface boundary layer. · Wrapping of U-shaped vortices around the pri 752

mary vortex tube. · Feeding of boundary-layer vorticity into the U- vortices . · Reconnection of the tips of the U-vortices with the free surface as the primary vortex wraps the tips of the U-shaped vortices into contact with the free surface. · Reconnection of the bases of the U-vortices with the free surface as the primary vortex tube collides with the wall a second time after the first rebound and carries the U-vortices into contact with the free surface. We have also shown evidence which suggests that striations that may be observed on the free surface above a pair of trailing tip vortices are caused by helical vorticity. 5 Conclusions We have developed numerical methods for the di- rect simulation of the three-dimensional Navier- Stokes equations with a free surface. These time- dependent simulations are feasible due to the small memory requirements and vectorized solution ca- pabilities of a unique application of multigrid finite-difference methods. The numerical schemes have been used to investigate the interactions of ring vortices and vortex tubes with walls and free surfaces. Our numerical studies of vortex tubes interacting with the free surface show that primary and sec- ondary vorticity may connect with the free surface. The free-surface signatures of these two mecha- nisms are distinctly different. One type of primary vortex reconnection forms a cusp pattern on the free surface, and secondary vorticity reconnection should appear as paired dimples on the free sur- face. The essential aspects of the reconnection of secondary vorticity with the free surface include the generation of U-shaped vortices by helical vor- ticity shed from the primary vortex tube, the wrap- ping of the U-shaped vortices around the primary vortex tube, and the reconnection of the bases and tips of the U-shaped vortices with the free surface. We also provide evidence which suggests that free- surface striations observed above a pair of trailing tip vortices are caused by helical vorticity emanat- ing from the primary vortex tube. There are many questions that need answering Ad problems that need resolving as we continue our efforts to develop computational models for three- dimensional viscous flows involving a free surface. Some of the more immediate issues that need to be addressed are: ~ Experimental validations. In the case of vortex rings interacting with walls, our numerical simulations, which show the forma- tion of secondary and tertiary vortex rings, agree qualitatively with experiments. In addition, we show that secondary and tertiary vortex rings gen- erate their own boundary layers, a phenomenon that has apparently not been observed in experi- ments. Similarly, secondary vortex rings are also shed by clean free surfaces at intermediate Eroude numbers. At lower Froude numbers, free surfaces appear to behave like free-slip walls so that as the ring expands outward radially, the core diameter is reduced due to stretching and the core appearance develops a distinctive head and tail shape. These low F,roude number results also agree qualitatively with-experimental measurements of vortex rings interacting with clean free surfaces. Until recently the technology did not exist for quantita- tive comparisons between numerical simula- tions and experimental measurements of time- dependent viscous flows. There have, how- ever, been significant recent advances in both these fields. In the area of flow visualiza- tion, three-dimensional graphical techniques are helping us to understand numerical sim- ulations of complex flows (Schisvone & Pap- athomas, 1990) heretofore not easy to com- prehend. Similar advances have also occurred in the laboratory. A noteworthy example is Willert & Gharib's (1990) Digital Parti- cle Image Velocimetry (DPIV) which allows the analysis of temporal evolutions of two- dimensional cuts of low speed flows. One apt plication of the DPIV technology is the specifi- cation of initial conditions for vortex ring sim- ulations. For example, the use of DPIV mea- surements as input to numerical models would allow detailed comparisons between theories and experiments. Validatioh studies of this sort should result in rapid progress in the un- derstanding of viscous free-surface flows. 7s3

· Radiation conditions. Wave reflections con- ta~n;nate the present simulations of a vortex tube moving parallel to a free surface, and in addition the vortex tube has strong interac- tions with its or images. IN order to alleviate these problems, possible solutions include the implementation of far-field matching bound- aries to eliminate wave reflections (Dommer- muth & Yue, 1987), and the use of multi-pole expansions to eliminate the effects of anages (e.g., Chamberlain & Liu, 1985~. · Nonlinear free-surface boundary conditions. Our present method assumes that the ampli- tude of the free-surface is smaller than the thickness of the boundary layer. This assump- tion is no longer valid at higher Reynolds num- bers and moderate Froude numbers. F`ully- nonlinear free-surface boundary conditions would remove this restriction, but their im- plementation will at least require the use of an adaptive grid (e.g., Thompson & Ferziger, 1989). . Surfaciants. In the open ocean, surfactants inhibit the generation of short waves by wind, and in the absence of wind the damping of short waves increases if surfactants are present (see Dorrestein, 1951~. In ship wakes, slicks consisting of compacted surfactant material form on either side of the ship track possi- bly due to the action of vortices moving sur- face water away from the False center (see Kaiser, et al, 1986~. These slicks damp cap- illary and gravity-capillary waves crossing the ship wake. In addition, experiments of Bernal, et al, (1989) show that the interaction of vor- tices with these slicks may lead to the forma- tion of secondary vorticity. A linear stabil- ity analysis of surface-tension effects indicates that limitations on time-step size would be too restrictive for explicit schemes. The situation does not to improve with the addition of sur- factants since the gradients of surface-tension must be balanced by viscous stresses within the free-surface boundary layer. We anticipate that a fully-implicit time-stepping scheme will be required. · Free-surface turbulence. Even without the ad- ditional complexity of a free surface, Frisch & Orszag (1990) note that the study of tur- bulence continues to pose unique challenges that are no less formidable than those occur- ring in post-Newtonian physics. Some of the unique aspects of free-surface turbulence are illustrated by the entrainment of air due to brealting waves and the exchange of energy at the air-water interface (Bonmarin, 1989~. Solutions to these problems will certainly re- qmre major developments in turbulence theo- ries and models currently being developed for unbounded and wall-bounded flows. In view of these formidable challenges, the present work represents only a modest first step toward what we believe will be a most exciting and fruitful · . Journey. Acknowledgements: This research is financially supported by the Office of Naval Research. Much of this work is based on the multigrid research of DGD which is supported by ONR Code 1215 (con- tract N00014-90-C-0027~. DKPY is supported by ONR Code 11 (contract N00014-90-J-1158~. A1- though the majority of DOD's research in the vis- cous flow area is a personal endeavor, he is grate- ful to Science Applications International Corpora- tion (SAIC) for their understanding and partial fi- nancial support. Most of the reported computa- tions have been performed on the Cray Research, Inc., Cray-2S through a grant Tom the Industry, Science, and Technology Department. Some com- putations have also been performed on the NSF Pittsburgh Supercomputer Center Cray Y-MP. We thank the members of the Naval Hydrodynamics Division at SAIC for their encouragement and sup- port. In particular, we are grateful to R.E. Hall, G.E. Innnis, D. Liepmann, and D.J. Loeser for their fruitful interactions, and we especially thank J.C. Talcott for his graphics expertise. REFERENCES Baba, E. (1969). A new component of Viscous resistance of ships. J. Sac. Naval Arch. Japan, 125, 23-34. Bcrnal, L.P., Hirea, A., Kwon, J.T., tic Willmarth, W.W. (1989). On the interaction of vortex rings and padre with a free surface for varying amounts of surface active agent. Phys. Fluids, A 1 (12), 2001-2004. Bonmarin, P. (1989). Geometric properties of deem water breaking waves. J. Fluid Mech., 209, 405-433. Brandt, A. & Dinar, N. (1977). Multigrid solution to elliptic flow problems. In Numerical methods in Par- tial Differential Equations (ea. S.V. Parter), Academic Press, 53-147. 754

Briley, W.R. & McDonald, H. (1973~. Solution of the three-dimensional compressible Na~rier-Stolces equa- tions by an implicit technique, Proc. Fourth Int. Conf. Num. Methods Fluid Dyn., Boulder, Colorado, liectu Notes in Physics, SS, Springer-Verlag, 105-110. Brown, E.D., Buchsbaum, S.B., Hall, R.E., Penhune, J.P., Schmitt, K.F., Watson, K.M., & Wyatt, D.C. (1988~. Observations of a nonlinear solitary wave paclcet in the Kelvin Wales of a Ship. J. Fluid Mech., 201, 263- 293. Brown, G.L. ~ Roshio, A. (1974~. On density effects and large structures in turbulent mincing layers. J. Fluid Mech., 61, 775-816. Cerra, A.W. & Smith, C.R. (1983~. Experimental oh secretions of vortex ring interaction with the fluid adja- cent to a surface. Dept. of Mech. Eng. and Mechanics, Lehigh University, Report FM-4. Chamberlain, J.P. & Liu, C.H. (1985~. Navier-Stolces calculations for unsteady three-dimensional Cortical flows in unbounded domains. AIAA J., 23: 868-874. Chen, C.-J. & Chen, H.-C. (1982~. The finite analytic method. Iowa Institute of Hydraulic Research, Report No. 232-IV. Chorin, A.J. (1967~. A numerical method for soldering incompressible viscous flow problems. J. Comp. Phys., 2: 12-26. Crow, S. (1970~. Stability theory for a pair of trailing vortices AIAA J., 8: 2172-2179. Dahm, W.J.A., Scheil, C.M., & Tryggvason, G. (1989~. Dynamics of vortex interaction with a density interface. J. Fluid Mech., 205, 1-43. Dommermuth, D.G. & Yue, D.K.P. (1987~. Numerical simulations of nonlinear axisymmetric flows with a free surface. J. Fluid Mech., 178: 195-219. Dorrestcin, R. (1951~. General linearized theory of the effect of surface films on water ripples. I & II. In Proc. Amsterdam Acad. of Sci., ti4, 260-272, 350-356. Dritschel, D.G. (1989~. Strain-induced vortex stripping. In Mathematical Aspects of Vortex Dynamics (ea. R.E. Caflisch), SIAM, 107-119. Fohl, T. ~ Turner, J.S. (1975~. Colliding vortex rings. Phys. Fluids, 18: 433-436. Frisch,U. & Orssag,S.A. (1990~. Turbulence: challenges for theory and experiment. Physics Today, January, 2t 32. Fromm, J. (1964~. The time dependent flow of an in- compressible viscous fluid. Meth. in Comp. Phys., B: 345-383. Fuchs, L. & Zhao,H.-S. (1984~. Solution of thrce- dimensional viscous incompressible flows by a multi- gridmethod. Int. J. Namer. Meth. Fluids, 4: 539-555. Chin, K.N., Hanlcey, W.L., ~ Hodge, J.K. (1977~. Study of incompressible Na~rier-Stolres equations in primitive variables. In Jrd Comp. Fluid Dyn. Conf., AIAA, 156-167. Grosenbaugh, M.A. & Young, R.W. (1988~. Nonlinear bow flows - an experimental and theoretical in~restiga- tion. InProc. 17thSymp. onNavalRydro., The Hague, 49-67. Haclcbusch, W. (1985~. Nulti-grid methods and appli- cations Springer-Verlag. Hartwich, P.-M. & Hsu, C.-H. (1988~. High resolution upwind schemes for the thrce-dimenaonal, incompress- ible Na~rier-Stokes equations. AIAA J., 26 Hanrey, J.K. & Perry,F.J. (1971~. Flowfield produced by trailing vortices in the vicinity of the ground. AIAA J. Eli: 250-260. Hirsa, A., Trygg~on, G., Abdollahi-Alibeil`, J., & Willmarth, W.W. (1990~. Measurement and compu- tations of vortex pair interaction with a clean or con- taminated free surface. In Proc. 18th Symp. on Naval Hvdro., Ann Arbor, To appear. Hirt, C.W., Nichols, B.D., & Romero, N.C. (1975~. SOLA-A numerical solution algorithm for transient fluid flows. Los Alamos Scientific Lab., Report No. LA- 5852. Issa, R.I. (1985~. Solution of the implicitb discretized fluid flow equations by operator- splitting. J. Comp. Phys., 62: 40-65. Kaiser, J.A.C., Garrett, W.D., Ramberg, S.E., Peltzer, R.D., & Andrews, M.D. (1986~. WAKEX 86: A ship wake/films exploratory experiment. Naval Research Laboratory, NRL memorandum Report 6270. Kambe, T. & Takso, T. (1971~. Motion of distorted vortex rings. J. Phys. Sac. Japan, 31 (2), 591-599. Kerwin, J.E. (1986~. Marine propellers. Ann. Rev. Fluid Mech., 18, 367-403. Kida,S., Taksolta, M., & Huesain,F. (1989~. Reconnec- tion of two vortex rings. Phys. Fluid, A 1 (4~: 630-632. Kim, J. & Main, P., (1985~. Application of a fractional- step method to incompressible Na~rier-Stolces equations. J. Comp. Phys., B9: 308-323. Kwon, J.T. (1989~. Experimental study of vortex ring interaction with a free surface. Ph.D. thesis, University of Michigan, Dept. of Aero. Eng. Lamb, H. (1932~. Hydrodynamics, Dover Publications. Leonard, A. (1980~. Vortex methods for flow simula- tion. J. Comp. Phys. 37, 289-335. Leonard, A. (198S). Computing three-dimensional in- compressible flows with vortex elements. Ann. Rev. Fluid Mech. 17, 523-559. Liepmann, D. (1990~. Ph.D. thesis, University of Cali- fornia, San Diego, Dept. of App. Mech. Lyden, J.D., Hammond, R.R., Lysenga, D.R., & Shuch- man, R.A. (1988~. Synthetic aperture radar imaging of surface ship walces. J. of Geophys. Res., fig (C10), 12,293-12,303. Madnia, K. & Bernal, L.P. (1989~. Interaction of a tur- bulent round jet with the free surface. The University 755

of Michigan Prog. in Ship Hydro., Report No. 89-05. Meiron, D.I., Shelly, M.J., Ashurst, W.T., & Orssag, S.A. (1989~. Numerical studies of vortex reconnection. In Mathematical Aspects of Vortex Dynamics (ea. R.E. Caflisch), SIAM, 183-194. Melander, M.V. & Hussain,F. (1988~. Cut-and-connect of two antiparallel vortex tubes. In Studying turbulence using numerical simulation databases - II: proceedings of the 1988 summer program, Center for Turbulence Rc- search, NASA Ames Research Center & Stanford Uni- ~rersity, Report CTR-S88. Miyata, H. & Nishimura, S. (1985~. Finite-difference simulation ofnonlinear strip wares. J.FluidMech., 157: 327-357. Mori, K. (1984~. Neclclace vortex and bow wave around blunt bodies. In Proc. 15th Symp. on Naval Hydro., Hamburg, 303-317. Munlc, W.H., Scully-Power, P., & Zachariasen, F. (1987~. The Bakerian Lecture, 1986. Ships from space. Proc. R. Soc. Land., A 412, 231-254. Norbury, J. (1973~. A family of steady vortex rings. J. Fluid Mech., 57, 417-431. Ohring, S. & Lugt, H.J. (1989~. Two counter-rotating vortices approaching a free surface in a viscous fluid. David Taylor Research Center, Report No. DTRC- 89/013. Ohsima,Y. & Asaka,S. (1977~. Interaction of two vortex rings along parallel axes in air. J. Phys. Soc. Japan, 42 (2), 708-713. Ohsima,Y. & Izutsu,N. (1988~. Cross-linking of two vortex rings. Phys. Fluids, 31 (9), 2401-2403. Orlandi, P. (1990~. Vortex dipole rebound from a wall. Phys. Fluids A, 2 (8), 1429-1436. Ors~ag, S.A. & Patterson, G.S. (1972~. Numerical sim- ulation of three-dimensional homogeneous isotropic tur- bulence. Phys. Rev. Lett. 28, 76-79. Patanlcar, S.V. & Spalding, D.B. (1972~. A calculation procedure for heat, mass, and momentum transfer in three- dimensional parabolic flows. Int. J. Heat Mass Trails., 15: 1787-1806. Patanltar, S.V. (1981~. A calculation procedure for two- dimensional elliptic situations. Namer. Heat Trance., 4: 409-425. Patera, A.T. (1984~. A spectral element method for fluid dynamics: laminar flow in a channel expansion. J. Comp. Phys. 54, 468-488. Peace, A.J. & Riley, N. (1983~. A viscous vortex pair in ground effect. J. Fluid Mech., 129: 409-426. Peltzer, R., J. Kaiser, O. Griffin, W. Barger, & Skop, R. (1990~. Hydrodynamics of ship wake surfactant films. In Proc. 18th Symp. on Naval Hydro., Ann Arbor, To appear. Roache, P.J. (1976~. Computational fluid dynamics. Hermosa. Roshko, A. (1954~. Structure of turbulent shear flows: a new loolc. AIAA J., 14, 1349-14S7. Saffman, P.G. (1978~. The approach of a vortex pair to a plane surface in in~riscid fluid. J. Fluid Mech., f)2, 497-503. Saffman, P.G. & Baker G.R. (1979~. Vortex interactions Ann. Rev. Fluid Mach. 11, 95-122. Sarplraya, T., & Henderson, D.O. (1984~. Surface dis- turbances due to trailing vortices, Naval Postgraduate School, Monterey, CA, Technical Report No. NPS- 6984004. Sarplcaya, T. (1986~. Trailing-vortex wakes on the free surface. In Proc. 16th Symp. on Naval Hydro., Berke- ky, 38-50. Sarpkaya, T., Elnitalcy, J., & Reeler, R.E. (1988~. Wake of a vortex pair on the free surface. In Proc. 1 7th Symp. On Naval Hydro., The IIague, 47-54. Schatzle,P.R. (1987~. An experimental study of fusion of vortex rings. Ph.D. thesis, California Institute of Technology. Scott, J.C. (1982~. Flow beneath a stagnate film on water: The Reynolds ridge. J. Fluid Mech., 11G, 283- 296. Schlichting, H. (1968~. Boundary-layer theory, McGraw-Hill. Schisvone, J.A. & Papathomas, T.V. (1990~. Visualiz- ing meteorological data. Bulletin of the Amer. Meteo- rological Soc. 71 (7), 1012-1020. Sheriff, K., Leonard, A., & Fersiger, J.H. (1989~. Dy- namics of a class of vortex rings. Ames Research Cen- ter, Moffet Field, CA, NASA Technical Memorandum 102257. Siggia, E.D. & Pumir,A. (1987~. Vortex dynamics and the existence of solutions to the Navier-Stokes equa- tions. Phy,. Fluids, 30, 1606-1626. Spreiter, J.R. & Saclts, A.H. (1951~. The rolling up of the trailing vortex sheet~and its effect on the downwash behind wings. J. Aero. Sci., 18: 21-32. Stanaway, S., Shariff, K., & Hussain, F. (1989~. Head- on collision of viscous vortex rings. In Studying turbu- lence using numerical simulation databases - II: proceed- ings of the 1988 summer program, Center for Turbu- fence Research, NASA Ames Research Center & Stan- ford University, Report CTR-S88. Streeter, V.L. (1948). Fluid dynamics, McGraw-Hill. Talcekuma, K. & Eggers, K. (1984). Effect of bow shapes on the tree-surface shear flow. In Proc. 15th Symp. on Naval EIydro., Hamburg. Telste, J.G. (1989). Potential flow about two counter- rotating vortices approaching a free surface. J. Fluid Mech., 201: Thompson, M.C. ~ Fersiger, J.H. (1989). An Adam tiers multigrid technique for the incompressible Navier- Stokes equations. J. Comp. Phys., 82: 9~121. 756

Vanlra, S.P. (1986~. Blocic-implicit multigrid solution of Navier-Stolres equations in primitive ra~iables. J. Camp. Phys., 65: 138-158. Walker, J.D.A. (1978~. The boundary layer due to a rectilinear vortex. Proc. R. Soc. Lond. A 359. 167- 188. Wallrer, J.D.A., Smith,C.R., Cerra,A.W., & Doligalski, T.L. (1987~. The impact of a vortex ring on a wall. J. Fluid Mech. 181. 99-140. Warming, R.F. & Beam, R.M. (1978~. On the con- struction and application of implicit factored schemes for consecration laws. SIA M-A MS proc., 11: 85-129. Wehausen, J.V. & Laitone, E.V. (1960~. Surface warren. Handbuch der Physik, fit: 446-778. Spinger. Willert, C.E. & Gharib, M. (1990~. Digital particle im- age velocimetry applied to an e~rol~ring vortex ring. To appear in Imp. in Fluids. Willmarth, W.W., Tryggvason, G., Hirsa,A., & Yu,D. (1988~. Vortex pair generation and interaction with a free surface. Phys. Fluid`, A 1 (20), 170-172. Winckelmans, G., & Leonard, A. (1989~. Numerical studies of vortex reconnection. In Mathematical Aspects of Vortex Dynamics (ea. R.E. Caflisch), SIAM, 25-35. Yu, D. & Tryggvason, G. (1990~. The surface signature of unsteady, tw~dimensional vortex flows. Submitted for publication. Appendix A: Axisymmetric Flow Equations Our studies of vortex rings interacting with walls and free-surfaces are based on an axisymmetric for- mulation. We summarize the basic equations in this appendix. The finite-difference forms of these equations are directly related to our full three- dimensional formulation. Let v and w respectively denote the radial (T) and axial (z) velocities, then the axisymmetric momentum equations for r > 0 are: Wit + (WW)z + (tower + r = -Pa + R (buzz + burr + r ~ Vc + (VW)z + (VV)r +-= -Pr + R (vzz + Err + r- 2 ~ . (47) For r = 0 the radial velocity is zero, v = 0, and we 771 + Vt7r -P + t2 r7 = i P + Fr2 t7 = use the following limiting form for the axial mo- mentum: I'` + (WU))z + 2(vut~r = -pa +-(wee + 2U)rr) Re Similarly, the equations for mass conservation are: we + or + - = 0 for r. ~ O Ritz + 2Vr = 0 fan r = 0 . (48) The gradient of the momentum equations (47) gives a Poisson equation for the pressure: Pzz + Prr + r = -wound-2vzwr-vrvr- 2 fOT r > 0 Pzz + 2Prr = -want*-2vztl)r-2vrvr for r = 0 . (49) If the thickness of the free-surface boundary layer is greater than the amplitude of the free-surface elevation, then the linearized free-surface boundary conditions for axisy~runetric flow are: vz = 0 for r.>0 = a' for r ~ O (`rlrr + ) for r > 0 W 17rr for r = 0 . (50) We use the following stream function formulation to calculate the initial velocity field for r > 0: jr V = -i ~ To = Vz-U)r ~ and fizz + Herr _ Nor = -red . (51) On r = 0, the limiting forms are v = 0, wO = 0, 1,6 = 0, and u' = Irr. 757

The conservation of energy equation with a lin- earized surface-tension term is as follows: 2 dt /~V r (W + V ~ + 2F2 d! / r q2 = W / r77~77rr +-~ R /V T [2~z + Or + 2 ~ + (VX + ~r)2] . (52) Appendix B: Multigric} solution of the hydrodynamic pressure Consider the following segment of code which is written in quasi-FORTRAN to illustrate a typical multigrid solution of the hydrodynamic pressure: CALL SOURCE IF(WILL PROBLEM) CILL FEASIBLE IF(FREE-SURFICE) CALL SURFICE_PRESSURE DO WHILE MULTIGRID DO WHILE FINE-TO-COARSE GRID CILL SOLVE_PRESSURE CALL RESTRICTION IF(WALL PROBLEM) CALL FEASIBLE ENDDO DO WHILE COARSE-TO-FINE GRID CALL PROLONGATION CALL MINUS ENDDO ENDDO CALL SOLVE_PRESSURE Subroutine SOURCE assigns the source term in the Poisson equation for pressure (see Eqs. 26, 29, and 32), and subroutine FEASIBLE makes the source term compatible with pure Neumann boundary conditions based on (40) and (41~. Since free-surface problems involve mixed Dirichlet and Neumann boundary conditions, subroutine FEA- SIBLE is not called, but instead subroutine SUR- FACE)RESSURE is used to assign the pressure on the free-surface (Eqs. 11 and 12~. Upon completion of subroutines SOURCE, FEASI- BLE, and SURFACELE'RESSURE, we implement a V-cycle multigrid algorithm using an unstaggered grid with uniform mesh coarsening. At the finest grid level the number of grid points along the x-, y-, and z-axes are respectively Imax = 2ipOw + 1, Jmas = 2ipow + 1, and Kmas = 2ipow + 1. The total number of grid points is Nto, = Imas X Jmas X Oman, and the next finest grid level has about 1/8 the number of grid points, and the coarsest grid level has 3x3x3=27 grid points. The total number of grid levels is equal to the maximum of ipoW, jam,, and know. The fine-to-coarse grid iterations start with a call to subroutine SOLVE~RESSURE which uses a vectorized Gauss-Seidel smoother with chequer- board ordering followed by two iterations of vec- torized Jacobi iteration. Upon completion of the Jacobi iterations, SOLVE)RESSURE also cal- culates the residual error. The calling argu- ments of SOLVE~RESSURE include vectors to store at all grid levels the pressure solutions, the source terms, and the residual errors. As a re- sult, each vector requires approximately (1 + 1/8 + 1/64 + ·~)N`o`=8N,o`/7 storage locations. An array pointer is used to position each vector at the correct grid level in the calling arguments of SOLVE)RESSURE. The pressure solver is dis- cussed in more detail in the next paragraph. Subroutine RESTRICTOR uses a 27-point restric- tor to calculate the source terms at the next coarser grid level based on the residual errors at the cur- rent grid level. The 27-point restrictor is calculated in three sweeps over the i, j, and k grid lines using a vectorized 3-point restrictor: 4~1,2,1~. If a won problem is being solved, subroutine FEASIBLE is called to make the source terms realizable at each of the coarse grid levels. The fine-to-coarse grid iterations continue until the coarsest grid level is reached, and then the coarse-to-fine grid iterations begin. Subroutine PROLONGATION uses 27-point pro- longation to interpolate the solutions of the coarse grid equations onto the next finer grid. The 27- point prolongation is calculated in three sweeps over the grid lines using a vectorized piecewise lin- ear interpolation. Subroutine MINUS corrects the solutions at the finer grid levels by using the inter- polated data from coarser grid levels. The coarse- to-fine grid iterations continue until the finest grid level is reached. At this point the V-cycles may stop if the residual errors are small enough; other- wise we start another V-cycle. We conclude the multigrid iterations with a post-smoothing step represented by a call to SOLVE]'RESSURE. 758

The Gauss-Seidel smoother in SOLVE}RESSURE is by far the most computa- tionally intensive portion of the multigrid scheme. In addition, the Gauss-Seidel smoother also has to be very general because it must work on all grid levels. Consider the following FORTRAN code for the Poisson solver with Neumann (wall) boundary conditions: 1. DX2= 1.0/(DI*DI) 2. DY2= 1.0/(DY*DY) 3. DZ2= 1.0/(DZ*DZ) 4. FIC=-0.5/(DI2+DY2+DZ2) 5. DO 10 TGlUS=1,ITER 6. DO 10 ITYPE=1,27 7. CILL INDICE(LEVEL,ITYPE,IBEG,IEND) 8. CALL JUMPER(ITYPE,ISKTP,JSKIP, 1 KSKIP,MOVE) 9. DO 10 IVECT=1,2 10. CILL DOCRlY(IGlUS,IVECT,TBEG,TEND, 1 ISTIR,IFINI,ISTEP) 11. CDIRt IVDEP 12. DO 10 I=ISTAR,IFINI,ISTEP 13. IC=IPOINT(I) 14. PNEW(IC)= ~ -(DX2*(PNEW(IC+MOVE(2)) 2 ~PNEW(IC+MOVE(1))) 3 ~DY2*(PNEW(IC+MOVE(4)) iPNEW(IC+MOVE(3))) +DZ2*(PNEW(TC+MOVE(6)) 6 +PNEW(IC+MOVE(5))) 7 -PSOU(IC))*FIC 15. 10 CONTINUE Lines 1 through 3 calculate off-diagonal elements in the linear system of equations. DX, DY, and DZ are the grid spacings, which depend on grid level, along the x-, y-, and z-axes respectively. FAC is the inverse of the diagonal element. The outer-most DO-loop at line 5 controls the number (ITER) of Gauss-Seidel iterations. The second DO-loop at line 6 sweeps over the 27 different portions of the computational domain: 1 set of interior grid points, 6 sets of side grid points, 12 sets of edge grid points, and 8 sets of corner grid points. Subroutine INDICE returns two pointers (IBEG & IEND) that describe where grid points are stored in a pointer array IPOINT as a function of grid level (LEVEL) and the type of bound- ary (ITYPE). Subroutine JUMPER describes how off-diagonal elements are related to the diago- nal elements as a function of boundary type. The arguments of JUMPER include skip param- eters which depend on the grid lerel. If the pressure is stored i-indices first, j-indices second, and k-indices last, then ISKIP=1, JSKIP=Imac, and KSKIP=ImaS X Jmax at the finest grid level. Subroutine JUMPER returns a 6-integer array n~ned MOVE which stores the skip parame- ters along the i-, j-, and k-a~ces. For the interior points, MOVE(1~=-ISKIP, MOVE(2~= ISKIP, MOVE(3~=-JSKIP, MOVE(4~= JSKIP, MOVE(5~=-KSKIP, aDd MOVE(6~=KSKIP. Sup- pose the grid points are on the i=1 side boundary, and on this boundary we specify free-slip boundary conditions for the pressure. We pretend that a fic- titious set of grid point e exists at i=-l;, and we set these grid points equal to the grid points at i=2. This is a second-order difference formula which can be implemented by setting MOVE(1~=ISKIP (note the positive sign). The pointer array IPOINT stores grid points se- quentially as a function of the boundary type and grid level. For exacnple, for the interior points at the finest grid level, we make the following assign- ments: 10 IBEG=1 IEND=0 DO 10 K=2,KMIl-1 KC=(K-1)*KSKIP DO 10 J=2,JMII-1 JC=(J-1)*JSKIP+KC DO 10 I=2,TMIX-1 IC=I+JC IEND=IEND+1 IPOINT(IEND)=IC CONTINUE Based on this type of ordering of IPOINT, the DO- loop at line 9 can set up a chequer-board ordering of the Gauss-Seidel solver. Subroutine DOCRAY controls sweeps over the odd or even grid points in the forward or backward direction as a function of IGAUS and IVECT. The table below summa- rizes what DOCRAY returns as ~ralues for ISTAR, IFINI, and ISTEP: | IGAUS | IVECT || ISTAR | IFINI | ISTEP | . ~ odd 1 IBEG IEND 2 odd 2 IBEG+1 IEND 2 . ~ even 1 IEND IBEG -2 even 2 ~ IEND-1 IBEG -2 759

Line 11 is a CRAY directive which specifies that the inner-most DO-loop vectorizes. The inner- most DO-loop implements Gauss-Seidel iteration using a 7-point star for the Laplace operator. The diagonal element is IC, and the assignments of the off-diagonal elements are controlled by MOVE. PNEW stores the pressures, and PSOU stores the source terms. Note that PSOU may also be used to store inhomogeneous Neumann boundary con- ditions. The inner DO-loop vectorizes because ISTEP=2 ensures that no diagonal element is a function of another diagonal element until the DO- loop is exited. As a result, half of the interior grid points vectorize in the Gauss-Seidel solver. The residual error that is associated with this vec- torizing scheme has some high wavenumber con- tent due to the alternating sweeps over the grid. As a result, the multigrid scheme will diverge unless measures are taken to smooth the resid- ual error. We perform this additional smoothing by using two Jacobi iterations in addition to the Gauss-Seidel iterations. Jacobi algorithms vector- ize without using chequer-board ordering, but their smoothing rate is slower than vectorized Gauss- Seidel and they require one additional vector for storage. However, as it turns out, we also need an additional vector to store the residual errors which we calculate after the Jacobi iterations are com- plete. We emphasize that the residual errors are calculated after the Jacobi iterations and not dur- ing the Jacobi iterations in order to provide the smoothest residual errors and the fastest conver- gence of the multigrid scheme. 760

I Ret | Nne~ot | NmU't | Nina, | At | Nramp | Ntime | Noper | I.0 1 2 1 .5 4 24 218 1 2 1 .5 2 24 218 1 2 1 .5 1 24 218 1 2 1 .5 0 * * 1 2 1 .25 0 * * = 1 2 1 .2 0 48 -408 1 2 1 .1 0 B5 572 = 1 2 1 .75 4 26 235 1 2 1 1.0 4 28 253 1 2 1 .25 4 33 298 1 2 2 .5 4 24 278 . 1 2 2 .25 4 32 365 1 2 2 1.0 4 26 299 1 2 3 .5 4 24 337 1 2 3 1.0 4 25 350 1 2 4 1.0 4 24 396 1 1 1 1.0 4 225 1175 1 1 2 .5 4 44 288 1 1 2 1.0 4 50 326 1 3 4 1.0 4 23 550 1 3 8 1.0 4 23 893 2 2 2 1.0 4 16 365 2 3 8 1.0 4 1 16 1217 1 2 1 .75 4 32 287 1 2 1 .75 2 32 287 1 2 1 .75 1 36 322 1 2 1 .75 0 * * 1 2 1 .01 0 * * 2 ~.005 0 710 6135 1 2 1 .5 4 33 296 1 2 1 .5 2 32 287 1 2 1 .5 1 38 322 1 2 1 .5 0 * * 1 2 1 1.0 4 34 304 1 2 2 .5 4 33 376 1 2 3 .5 4 34 469 1 3 1 .5 4 34 426 1 1 1 .5 4 258 1347 2 2 1 .5 4 25 442 2 1 1 .5 4 80 630 Table 1: Convergence of implicit scheme to steady-state The grid Reynolds number is denoted by R^, and the time step is At. The number of Newton-Raphson iterations, multigrid iterations, and Gaus~Seidel iterations are respectively NATO ~ NmU`', and N'aU, . The number of time steps and sweeps over the finest grid level are respectively Name and Noper . The number of time steps required to ramp down the initial Viscosity to the prescribed viscosity is Tramp. The initial viscosity is L,o = 1 for all cases. Most of the solutions have converged to a tolerance SO = 10-6. Solutions that have diverged are denoted by an asterisk (*~- 761

| Item || Run 1 | Run 2 | Run 3 | Run 4 | Run 5 | Run 6 | Code? exp. exp. imp. imp. exp. exp. B.C.? no no no no no dip Rc 400 400 400 400 400 400 TC 0.5 0.5 0.5 0.5 0.5 0.5 WC ~4 -4 -4 -4 -4 -4 To 1 1 1 1 1 1 20 -3 -3 -3 -3 -3 -3 Rt 4 ~ 4 4 4 6 Dt 4 4 4 4 B ~ At 0.005 0.0025 0.005 0.0025 0.0025 0.0025 Ntime 3000 8000 3000 6000 6000 6000 Nnew, 1 1 Nmu` ~2 2 2 2 2 2 N9au. 8 5 8 8 8 Imps 129 257 129 257 257 257 Kmar 129 257 129 257 257 257 Table 2: Data for ring vortices impinging on a malt. The 'Code' item indicates whether the explicit (exp.) or semi-implicit (imp.) computer codes had been used to perform the calculation. The 'B.C.' item indicates whether a no-slip (no) or slip (slip) boundary condition is used on the wall at z = 0. Re = I`/~'r ~) is the Reynolds number. The initial core radius and peals ~rorticity are denoted by rc and we. The vortex core is centered at (To' Zo). The radius and depth of the tattle are respectively R. and D,. The time step is At. The number of time steps, Newton-Raphson iterations for the implicit scheme, multigrid iterations, and Gauss-Seidel iterations are respectively denoted by NtjmC ~ Nne,,,' ~ Nmu~t ~ and N9aU, . The number of grid points along the r- and z-axes are respectively Imps and Kmax ~ | Item || Run 1 | Run 2 | Run 3 | Run 4 | Run 5 | Run 6 | Rc 400 400 400 200 200 200 Fr2 0.5 0.25 0.125 0.5 0.25 0.125 tic 0.5 0.5 0.5 0.5 0.5 0.5 . . WC -4 -4 -4 -4 -4 -4 To 1 1 1 1 1 1 Ret 3 6 6 6 6 6 Dt 4 4 4 4 4 4 i\t 0.0025 0.0025 0.0025 0.0025 0.0025 0.0025 Ntime 6000 6000 6000 6000 6000 6000 Nmult 2 2 2 2 2 2 Ima z 257 257 257 257 257 257 . Kma ~257 257 257 257 257 257 Table 3: Data for ring vortices impinging on a free-surface. Most of the entries in this table are defined in Table (2), except for the Froude number Fr = I`/~,r 9~/2 (To)3/2~. 762

| Item || Runt | Run2 | Run3 | Run4 | Runs I Rune I Run7 | B. . C. ? noslp noslp noslp free free flee free Rc 200 200 400 200 200 400 200 Fr 0 0 0 0.25 0.25 0.25 0.25 rc 0.5 0.25 0.25 0.5 0.25 0.25 0.5 sO 4 16 16 4 16 16 2 Do -3 -3 -3 -3 -3 -3 -1 Camp -0.25 -0.25 -0.25 -0.25 -0.25 -0.25 -0.5 be -0.25 -0.26 -0.26 -6.26 -0.26 0 25 -0.5 we 4 4 4 4 4 4 4 Dt 1 4 4 4 4 4 4 4 At 1 o.oo5 0.005 0.005 0.005 0.005 0.005 0.005 Ntimc 1 3000 3000 4000 3000 3000 3000 3000 Nmult 4 4 4 4 4 4 4 NgauJ 1 8 8 8 8 8 8 8 Oman 65 65 65 65 65 65 65 Jma ~65 65 65 65 65 65 65 Oman 65 65 65 65 65 65 65 Table 4: Data for vortex tubes impinging on walls and free surfaces. Runs 1 thru 6 use free-slip boundary conditions on the vertical walls, and Run 7 uses periodic boundary conditions along the :z-axis. The axes of the vortex tubes are aligned with the ~axis. The 'B.C.' item indicates whether a no-slip (noslp) or free-surface (free) boundary condition is used on the wall at z = 0. The Reynolds number is denoted by Re and the Froude number is Fr. The initial core radius and peak vorticity are denoted by TC and we. The mean position of the vortex core is centered at (:rO, GO). The initial amplitudes of the sinusoidal perturbations in the core's position are Camp and Camp The length (art), width (y), and depth (z) of the tattle are respectively Lo, W`, and D`. The time step is /`t. The number of time steps, multigrid iterations, and Gauss-Seidel iterations are respectively denoted by N`,me, NmUl`, and Ngou~ . The number of grid points along the ::-, an, and z-axes are respectively Ima2 ~ Jmae, and Oman · 763

-1 .8 -2 . 0 - 2 . 2 -2 .4 -2 .6 _ - 2 . 8 12 ~ -3.6 - 1 _ IIIII,,I,1,,,,1,,,,1,,,,11,,,I,,,,I,,,,I,,,,I,,,, .5 1 .0 1.5 2.0 2.5 - 3 . 0 -3 .4 -4 .0 -4 .2 -4 .4 -4 .6 0 -~ . ~ -2 .4 -2 .6 -2 .8 -3 . 0 -3 .2 -4 .0 /// Ret Figure la: Radial velocity. 3.0 3.5 4.0 4.5 5.0 Figure 1: The spatial accuracy as a function of grid Reynolds number for a jet impinging a wall. The logic maximum absolute errors in the (a) radial and (b) axial velocities are plotted ver- sus grid Reynolds number (Rig) for axisymmetric ~ 1-~ and three-dimensional ~ 2 ~ re- sults. I, .,,,,,, I I I I I I I I I I, · I, I I I, I I I I I,,, I I, I I I I I ~ - / 1 ~ , , , 1 , 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 Rat Figure lb: Azial velocity. 10.2, I , , , I I I ~I ~I ~I r 10.0 9.8 9.6 9.4 9.2 9 .0 B.e B.6 8.4 8.2 8.0 7.8 - log2 Noper / ': ~ 7 6 ~ / 7 . 2 , 1 1 1 , 1 1 9 10 11 12 13 14 15 16 lOg2 N-tot Figure 2: The performance of multigrid scheme as a function of the number of grid points for a jet im- pinging a mall. The log2 of the normalized number of operations (Noper) is plotted versus the loge of the number of grid points (Ntot). 764

/ 1.4 _ 1.2 _ 1.0 _ .8 _ .6 _ .4 _ .2 _ 0 / 1 ,.~ 1.2 ~ 1 .0 ~ .8: .6 _ .4 _ .2 0 / in' /] /~ _\ ~ ~, .5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 Rat .8 _ 1.6 .2 _ ~1 ~3 it // 1_/\-j Zen it// _ .5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 Rat 3 Figure 3a: ~ = 0.7. / Figure 3c: By = 1.4. /\ it, / / ~, , ~ - .-1 - 1 .0 2 .~ Rat 3.0 Figure fib: ~ = 1.0. 3._ 4.0 Figure 3: The values of stability parameters at the maximum allowable time step as a function of grid Reynolds number for a FTCS scheme in two di- mensions. The parameters are labeled as follows: (-1 ) denotes 2(ct~ + av); (-2 ) de- notes,02t(2Ltr); ( 3 ) den°teS§v/(2~'v); and ( 4-) denotes (by + pV)2/(2(~= + all)). The results are plotted for different velocity ratios (a): (a) ~ = 0.7, (b) ~ = 1.0, and (c) ~ = 1.4. 765 5.0

l 030 028 026 024 02c 020 018 016 0 012 .016 . 008 . 006 004 , 0 .40 .30 .25 20 . 1 0 .05 . _ ,, I I, I, I I I, it,, 1 4 6 8 10 12 14 16 18 20 Rid Figure 4a: T = 1. l 0 0 2 4 6 8 RFigure 4b: T = 10. Figure 4: The relative error in the maximum vor- ticity as a function of grid Reynolds number for a Gaussian vortex core orbiting a boo. The results are plotted for At = 0.01 ( 1 ) ; At = 0.02 ( 2 ) ;andAt=0.05( 3 ) . 766 20

. 024 . 022 .018 .016 .014 .012 .010 . 008 . 006 . 004 . 002 4~ .30 7~ 20 \ 0 i it_ 1 - 8 10 Rip Figure 5a: T = 1. 14 16 l 18 20 I , I , I , I -, I T I · r- i ~ 0 2 4 6 8 10 12 14 16 18 20 Rip Figure fib: T-10. Figure 5: The relative error in an energy conserva- tion lab) as a function of grid Reynolds number for a Gaussian vortex core orbiting a boa. The results are plotted for At = 0.01 ~ 1-~ ; At = 0.02 ~2- ~;andAt=0.05( 3 ~. . 036 . 034 . 032 . 030 . 028 . 026 . 024 . 020 .018 .016 .014 .012 .010 . 008 . 006 . 004 . 002 0 767 ''''l~l~I....,....,....,....,. /- . 022 _ 2\ 2.0 2.5 3.0 3.5 4.0 R:, 4.5 5.0 5.5 6.0 Figure 6: The relative error in the wave amplitude as a function of grid Reynolds number for an atten- uating azisymmetric standing wave. The relative errors c~ ~ 1 ~ and e2 ~ 2 ~ are defined in the text.

- / FREE-SLIP WALL (BACK FACEJ - - - 1 NO-SLIP WALL FREE-SLIP WALL OR FREE-SURFACE Z / \ (TOP FACE) AL l - ~ . ~ ~' ~ 6~^ , FREE-SLIP WALL (BOTTOM FACE) Figure 7: The numerical simulation of a vortex ring impinging a boundary. The initial core radius and peak vorticity are denoted by rc and ~c. The vortex core is centered at (To, Zo). The radius and depth of the tank are respectively R. and D`. FREE SLIP WALLS OR PERK)D C BOUNDARIES {SIDE fACES) ~ \ FREE-SLIP WALL ~ (80 TT0M fACE) Z ~ 'I ~ J ~ _ Rae FREE-SLIP WALL (OUTER RADIUS) * l NO-SLIP WALL OR FREE SURFACE (TOP fACE) ~ ~ r LFREE-SLIP WALL ~ ~ ~(FRONT FACE) ~TO VORTEX TUBE\ ~ Figure 8: The numerical simulation of a vortex tube impinging a boundary. The initial core ra- dius and peak vorticity are denoted by He and '`)c. The mean position of the vortex core is centered at (HO, HO). The length (air), width (y), and depth (z) of the tank are respectively L`, W., and D`. 768

a: l l l l ~ l l l Figure 9: A ring vortex impinging a free-slip wall. The w=-1 contours of the primary vortex core are plotted at different instants of time: t=O, 2, 4, 6, 8, 10, 12, & 14. The ring vortex starts at the bottom of the figure and expands radially outward. The numerical parameters for this run are provided in Run 6 of Table (2~. Figure 10: A ring vortex impinging a no-slip Tall. The w=-1 contours of the primary vortex core are plotted at different instants of time: t=O, 2, 4, 6, 8, 12, 16, & 20. The ring vortex starts at the bottom of the figure and expands radially outward up to r ~ 2.25. The numerical parameters for this run are provided in Run 2 of Table (2~. in in - ~9 r\i ~9 Ln z jets", ~ \: ~ m 6} 1 0 I ~ I,, i I, I, I I,,, I,,, l - z |CONTOUR FROM -1 TO -1 BY O .5 1 .0 1.5 2.0 2.5 3.0 3.5 4 0 4.5 5.0 5 5 6.0 T 1 ; ,~\-; r\i US N l ray . _ ~ 0 .5 at, / \ ~ 1 ~ ~! ~ , 1 0 1.5 2.0 2 5 I CONTOUR FROM -1 TO -1 BY O ,,, i 3.0 3 5 4 0 T 769

z - r~ l: . - -- - ~ - - - - ~l ~ Nl ~2 . 4 . 6 Figure 12: The orbit of a secondary vortex ring. The solid and dashed contour lines represent pos- itive and negative vorticity respectively at time t=15. Note that only a small portion of the compu- tational domain is plotted. The numerical param- eters for this run are provided in Run 2 of Table (2~. .8 1.0 1.2 1.4 1.6 1.8 2.0 r I CONTOUR FROM -5 TO 17 BY .25 Figure 11: The boundary-layer separation due to a ring vortex interacting with a no-slip wall. The solid and dashed contour lines represent positive and negative vorticity respectively at time t=7. Note that only a small portion of the computa- tional domain is plotted. The numerical parame- ters for this run are provided in Run 2 of Table (23. _ r rim '1 .0 1.2 1.4 1 .6 ~ 5 ~ _~\\\~` 1 ~ At\ " ~"2~~' / I/// 770 rid . ., , , , , , , , , , 1.8 2.0 2.2 2.4 2.6 2.9 3.0 r

~ - r~ 1 z - ~ll'\~ --me r 1 CONTOUR FROM -2 .25 TO 2.75 BY . 2 Figure 14: The energy conservation for ring vor- tices impinging on no-slip and free-slip walls 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) inte- grated over time). The results are plotted as fol- lows: ~ 1 ~ denotes Ett)/E(O) for the no- slip case, ~ 2 ~ denotes (E(O)+ D(t)~/E(O) for the no-slip case, ~-3 ~ denotes E(t)/E(O) for the free-slip case, and ~ 4- ~ denotes (E(O) + D(t)~/E(O) for the free-slip case. The nu- merical parameters for these runs are provided in Runs 2 & 6 of Table (2~. - I , 1 ~ ~ I I , 1 ,~ 1 , ili = ~~ High resolution | CONTOUR FROM -1 TO .25 BY 1 .25 _ '1 ~1 .2 1.4 1 .6 1.8 2.0 2.2 2.4 2.6 2.8 3 0 r Figure 13: The formation of a tertiary vortex ring. The solid and dashed contour lines represent pos- itive and negative vorticity respectively at time t=20. Note that only a small portion of the compu- tational domain is plotted. The numerical param- eters for this run are provided in Run 2 of Table (2~. W,'~ ' ' i i 1 ' I \ free-slip wall no-slip wall I\ . 1 1 1 ~ 1 1 1 1 1 1 1 , 1 1 1 0 2 4 6 8 10 12 14 t 1 1 1 1, ,_,, 16 18 20 Figure 15: A comparison of low and high resolution simulations of a ring vortex interacting with a no- slip wall. The solid and dashed contour lines corre- spond to ~ = .25 and the dashed contour lines cor- respond to ~ = -1 at time t=15. Note that only a small portion of the computational domain is plot- ted. The numerical parameters for these runs are provided in Runs 1 & 2 of Table (2~. 771

.05 .04 .03 .02 .01 0 .01 .02 .03 .04 1 1 Ln z l Figure 16a: Free-surface Elevations - - \ ~ 1 ad\ i, \ 1 ~ l ,,,,1,,,,IIII11,,,11,,,,11,I11,,,,1,,,,1,,,,I,,,,1,,,,1,III . . . . . . . . . . . . . . 0 .5 1.0 1.5 rid In rid -1 ~ 2.0 2.5 3.0 3.5 4.0 4.5 5.0 5.5 6.0 ( U) \ / it\ ~Xi , , . . . . . . . . . . . . . . . . \ / Figure 16b: Vorticity Contours ' i ' ~Z g ' ' 3 ~ d ' 3 5 ~ O- ~ 5-5 O 5 5 ~ ~ .5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 5.5 b.0 r Figure 16: A ring vortex interacting with a free surface at a loo' Froude number. Part (a) plots the free-surface elevation as a function of the radial coordinate at times t=2, 6, 10, & 14. Part (b) plots = -1 contours at times t=O, 2, 4, 6, 8, 10, 12, & 14. The ring vortex starts at the bottom of Part (b) and expands radially outward. The numerical parameters for this run are provided in Run 3 of Table (3). 772

.30 .25 .20 . 10 0 - . 05 . 10 n 1 ~ I, i , I,, I, I I I I I I I I I I I I I I I I I I I I I I I ~ ~ I ~ I I ~ ~ ~ I ~ ~ ~ ~ I ~ ~ ~ ~ I 1 1 1 ~1 1 ~ \ ~\\ Figure 17a: Eree-surface Elevations \ ~ 1\~\ ~ '° \ ~ \ ~ \~\'t~' 1.0 1.5 - 7 I I 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 ~ I I I I I I I I 0 .5 2.0 2.5 3.0 3.5 4.0 4.5 r ~ 1 1 1 1 1 1 1 ~1 1 1 1 1 ~ 1 ~ 1 ~ 1 1 ~ 1 1 1 1 ~1 1 1 ~ ~ 5.5 6.0 x, ~ iN~/~ X r~ r ( ~\ Figure 17b: Vorticity Contours . - ~CONTOUR FROM -1 TO -1 BY O : -- , ,,,.,,,,,,,,,,,,,,,,,,,,,,,,, , .51.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 5.5 6.0 T Figure 17: A ring vortez interacting ?I)ith a free surface at an intermediate Froude number. Part (a) plots the free-surface elevation as a function of the radial coordinate at times t=O, 1, 2, 3, 4, 5, 6, & 7. Part (b) plots ~ = -1 contours at times t=O, 2, 4, 6, 8, 1O, 12, & 14. The ring vortex starts at the bottom of Part (b) and expands radially outward up to r ~ 2.25. The numerical parameters for this run are provided in Run 1 of Table (3~. 773

oW N 1 . - 1 ~ l '_ Z I fir z q" - · . , , , ~ ~,~ all 1 ' ~. . . . . ~ ~.5 1.0 1.5 . . . //~ / i, ,, ~ = ~ ~1,!: l 2.0 2.5 T l _ CONTOUR FROM -4 TO 2.75 BY Figure 18: The orbit of a secondary vortex ring that is shed by a free surface at an intermediate Froude number. The solid and dashed contour lines repre- sent positive and negative vorticity respectively at time t=12. Note that only a small portion of the computational domain is plotted. The numerical parameters for this run are provided in Run 1 of Table (3~. . . . . . . . . . . . . . . . . 3.0 3 5 4.0 4.5 _ _ , I, / I I / t\J//4 1~ I CONTOUR FROM -5 .25 TO 2 .5 BY .25 . .., .,, ., ~ ~ .. I ~5 1.0 1 5 2.0 2.5 3.0 3.5 4.0 4.5 r Figure 19: The shedding of a tertiary vortex ring by a free surface at an intermediate Froude num- ber. The solid and dashed contour lines represent positive and negative vorticity respectively at time t=15. Note that only a small portion of the compu- tational domain is plotted. The numerical param- eters for this run are provided in Run 1 of Table (3~. 774

Boundary layer / on the wall - Inner helical vortex sheet Figure 21: The attachment of U-shaped vortex tubes to the separated boundary-layer on a wall. A constant isosurface of vorticity magnitude is shown. The primary vortex tube had been mov- ing from right to left before being stopped by the separated boundary-layer on the wall. To empha- size the U-shaped feature this image has been re- flected about its midspan relative to the image at the top of the page. Although the figures on this page have been adapted from a numerical simula- tion of a vortex tube impinging on a no-slip wall, we expect similar behavior for a contaminated free- surface. vortex sheet Cross-section of attached U-shaped vortex tube Sheath of secondary vorticity wrapping around primary vortex Unattached U-shaped vortex tube Base of U-shaped vortex tube 775 Figure 20: The unwinding of helical vorticity around a primary vortex: tube. A constant isosurface of vorticity magnitude is shown. The view is from below a no-slip wall, and the primary vortex is moving into and to the left of the page due to its images across the centerplane and above the wall. The neighboring vortex that is on the right side of the page is not shown. Primary vortex tube

Figure 22: The free-surface cusp pattern formed by the reconnection of primary vorticity with a free surface. The of = +1 contours of primary vorticity on the free surface are plotted at differ- ent instants of time: t=8, 97 TO, 1l, 12, 13, 14, & 15. The submerged portion of the primary vortex tube is moving from left to right across the page. To emphasize the cusp pattern this image has been reflected across the symmetry boundary at y = 4. The numerical parameters for this run are provided in Run 7 of Table (4~. Figure 23: The evolution of a vortex tube travel- ing slightly submerged below a free-surface. The ~ = 1.25 isosurface of the primary vortex core is plotted at different instants of time: (a) t=2.5, (b) t=5, (c) t=7.5, (d) t=lO, & (e) t=15. The view is looking down on the free surface, and the vortex tubes are moving from right to left due to their images above the free surface. The computational domain has been doubled in length to avoid confu- sion as portions of the vortex tube enter and exit the periodic boundaries at the left and right sides of the page. The geometry triad is in a submerged corner of the computational domain. The numeri- cal parameters for this run are provided in Run 7 of Table (4~. lo to Y ~ ,~'''' '_," Air'?' it, -,,,,,,,,,,....,,,..,!~: 0 5 1 0 1 5 2 0 2.5 3 0 3 5 4.0 x Figure 23a: Iwl = 1.25 isosurfaces at t=2.5. 776

Figure 23b: two = 1.25 isosurfaces at t=5. ., Figure 23c: ALAS-1 25 - . 1sosurfaces at t=7.5. 777

Figure 23e: 1~1 = 1.25 isosurfaces at t=15. 778

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

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

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

- 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 -

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

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

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.

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

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.

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

Next: On the Numerical Solution of the Total Ship Resistance Problem under a Predetermined Free Surface »
Eighteenth Symposium on Naval Hydrodynamics Get This Book
×
Buy Paperback | $275.00
MyNAP members save 10% online.
Login or Register to save!

This volume contains technical papers and discussions covering ship motions, ship hydrodynamics, experimental techniques, free-surface aspects, wave/wake dynamics, propeller/hull/appendage interactions, and viscous effects.

  1. ×

    Welcome to OpenBook!

    You're looking at OpenBook, NAP.edu's online reading room since 1999. Based on feedback from you, our users, we've made some improvements that make it easier than ever to read thousands of publications on our website.

    Do you want to take a quick tour of the OpenBook's features?

    No Thanks Take a Tour »
  2. ×

    Show this book's table of contents, where you can jump to any chapter by name.

    « Back Next »
  3. ×

    ...or use these buttons to go back to the previous chapter or skip to the next one.

    « Back Next »
  4. ×

    Jump up to the previous page or down to the next one. Also, you can type in a page number and press Enter to go directly to that page in the book.

    « Back Next »
  5. ×

    To search the entire text of this book, type in your search term here and press Enter.

    « Back Next »
  6. ×

    Share a link to this book page on your preferred social network or via email.

    « Back Next »
  7. ×

    View our suggested citation for this chapter.

    « Back Next »
  8. ×

    Ready to take your reading offline? Click here to buy this book in print or download it as a free PDF, if available.

    « Back Next »
Stay Connected!