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.
24th Symposium on Naval Hydrodynamics Fukuoka, JAPAN, 8-13 July 2002 Numerical Simulation of Wakes in a Weakly Stratified Fluid James W. Rothmans, Douglas G. Dommermuthl, George E. Innisl, Thomas T. O'Shea1 and Evgeny Novikov2 ([Science Applications International Corooration. USA. ~ , ~Un~vers~ty ot Aorta, San Diego, USA) ABSTRACT This paper describes some preliminary numerical stud- ies using large eddy simulation of full-scale submarine wakes. Submarine wakes are a combination of the wake generated by a smooth slender body and a number of superimposed vortex pairs generated by various control surfaces and other body appendages. For this prelimi- nary study, we attempt to gain some insight into the be- havior of full-scale submarine wakes by computing sep- arately the evolution the self-propelled wake of a slender body and the motion of a single vortex pair in both a non-stratified and a stratified environment. An important aspect of the simulations is the use of an iterative proce- dure to relax the initial turbulence field so that turbulent production and dissipation are in balance. INTRODUCTION The prediction and evolution of the full-scale wakes of submarines has been an objective of naval hydrodynam- ics for many years. Full-scale submarine wakes consist of the wake behind a slender body on which is superim- posed a number of vortex pairs generated by various con- trol surfaces. Coherent vortex systems are major compo- nents of the wakes of submarines. For this preliminary study, we attempt to gain some insight into the behavior of full-scale submarine wakes by using LES to compute separately the evolution of the self-propelled wake of a slender body and the motion of a single vortex pair in a stratified and turbulent environment. An important as- pect of the simulations is the use of an iterative procedure to relax the initial turbulence field so that turbulent pro- duction and dissipation are in balance. In the first section we discuss the studies of the wake of a slender, self-propelled body without any superim- posed vortex pairs. The objective of this research is to be able to simulate accurately a full-scale wake of a self- propelled body using initial conditions on a cross plane a short distance downstream of the submarine. The ini- tial conditions are obtained from laboratory experiments. In the second section we describe the simulations of a control surface vortex pair in both non-stratified and sta- bly stratified and turbulent environments. The ultimate objective of this study is to develop a subgrid scale tur- bulence model to accurately represent the turbulence in the vortex cores and to be able to predict the turbulent entrainment into and detrainment out of the recirculation region of the pair and therefore to be able to accurately compute the motion of the pair at full-scale. WAKE STUDIES Laboratory experiments, in the absence of any back- ground shear, have shown that the far wake behind a self- propelled body consists of very slowly evolving patches of vertical vorticity, of alternating sign, whose horizontal extent is much greater than their vertical extent. These patches are often referred to by the descriptive name of "pancake eddies". The governing parameters for this flow are the Reynolds number Re = UD/zJ and the Froude number Fr = U/(ND), where U is the speed of the body, D is the diameter of the body, zip is the kinematic viscosity of the fluid, and N = ~g/pO ~Q/0X34i/2 is the buoyancy (or Brunt-Vaisala) frequency in which 9 is the accelera- tion due to gravity, pa is the mean density, and ~Q/0X3 is the vertical derivative of the background mean density (assumed to be a constant). In this paper we describe two numerical experi- ments: the first is an attempt to simulate the formation of eddies in late wake due to a turbulent flow for the case with Fr = 2.0 and Re = 105. The second is similar to the first except that the fluid is non-stratified (Fr = x). The numerical method is large eddy simulation (LES). The flow is initialized with a mean wake flow that includes swirl, which is based on experimental mea- surements of the wake behind a propeller-driven slender body, with a superimposed homogeneous turbulent flow field. A relaxation procedure is introduced to establish the proper balance between production and dissipation before the calculation is begun. There are no coherent structures imposed on the turbulence other than the gross characteristics of the mean wake flow. The results are compared with that of a drag wake that is described in a companion paper, Dommermuth, et al. (2002), and it
is found that the two types of wakes evolve quite differ- ently. The mean swirl generates internal waves. Once the mean swirl radiates internal waves, the turbulence pro- duction terms are considerably reduced, especially the shear stresses. As a result, the wake persists much longer than a drag wake. Our results show that coherent vortices appear in the late wake even though the flow is initialized without any coherent structures. PROBLEM FORMULATION A schematic drawing of the flow under considera- tion is shown in Figure 1 in a reference frame in which the body generating the wake is at rest. This figure serves to define much of the nomenclature used here. The wake is considered to be statistically stationary. Since the en- tire wake is too long to compute as a whole, we make the approximation that the flow can be computed within a rectangular box, with axial dimensions much smaller than the total length of the wake, that moves with the mean flow speed U. Within this box the flow is computed by solving numerically the governing equations for an incompressible Boussinesq fluid using large eddy simu- lation (LES). Large eddy formulation of the Boussinesq equations We will assume that the fluid is incompressible and weakly stratified. A large eddy approximation is invoked whereby the large-scale features of the flow are resolved and the small scales are modelled. Let hi denote the filtered three-dimensional velocity field as a function of space xi (i = 1, 2, 3) and time t. Here, the overbar de- notes spatial filtering in a large eddy sense. The origin of the coordinate system is at the centroid of the body, as shown in Figure 1. x is positive downstream, x2 is transverse to the track of the body, and X3 is positive up- ward. The length and velocity scales of the flow are re- spectively normalized by the diameter of the body (D) and the free-stream velocity (U). The numerical method is large eddy simulation (LES). The flow is initialized with a mean wake flow with a superimposed homogeneous turbulent flow field. A re- laxation procedure is introduced to allow the turbulence field to achieve a balance between production and dissi- pation before the calculation is begun, as is described in more detail in Dommermuth, et al. (20021. There is no coherent structure imposed on the turbulence other than the gross characteristics of the mean wake flow. A mixed model (Bardina, et al., 1984) is used to model the subgrid scale stress tensor. The similarity por- tion of the mixed model provides an accurate model of the turbulent stresses, whereas the Smagorinsky portion provides dissipation. The corresponding mixed model for the residual density flux combines a similarity model with an eddy diffusivity model. Figure 1: A schematic diagram of the self-propelled body and the coordinate system in a reference frame in which the body is at rest. The dashed box denotes the slab of fluid that is modelled using LES. but ~ is the wake deficit. Following Orszag & Pao (1974), a Galilean approx- imation is used to relate the spatial development of the wake to the temporal evolution of the LES. In normal- ized variables, this approximation results in the relation x~ = t, where x~ is the distance downstream of the body in the wake and t is the corresponding time in the LES. Based on this Galilean approximation, we further assume that | dto(t, x = Xo) ~ ~ o {L L J dx¢(t = Tome) = (I)) ~ (1) To where o is a physical quantity, hat accents denote time averaging, angle brackets denote spatial averaging, L is the length of the LES computational domain in the axial direction (see Figure 1), T is the duration of time averag- ing, and XO and To are positions in space and time where the wake of the body and the LES correspond. A tilde ac- cent denotes the turbulent fluctuations, which are defined as fib = 0 - (o). As shown in Figure 1, the wake of the body is modelled as a slab of fluid. Initialization The initial velocity field is decomposed into a mean disturbance and a fluctuating disturbance. The magni- tude and distribution of the mean and fluctuating compo- nents are specified based on the measurements of Lin & Pao (1973, 1974a,b) and Lin, Veenhuizen & Liu (1976~. The mean axial velocity is specified as 2
0.1 ( r2 ) ( r2 ) , (2) where aO is the amplitude of the mean wake deficit nor- malized by the free-stream velocity, and rO is the initial characteristic radius of the wake. The cross-stream com- 0.001 portents, (u2) and (U3) are determined from the mean axial vorticity, which is specified (Qua = aw (1 - 2 2 ) exp (_ r ) , (3) where a``, is the amplitude of the mean axial vorticity and rw is the characteristic radius of the axial vorticity field. The propeller swirl consists of an annulus of vorticity of one sign from the tip vortices surrounding a region of opposite sign vorticity from the root vortices, with a net vorticity of zero. Interestingly, in another numerical experiment (not shown) we have found that the swirl is necessary if the axial velocity is to retain its self-similar shape. Without the swirl, the axial velocity rapidly disin- tegrates into small regions of positive and negative flow, and the rate at which they cancel each other out is more rapid than with swirl present. The initial rms velocity fluctuations are also approx- imated using Gaussian distributions. ~ = at exp ( - 2 2 ) (4) where ai are the initial amplitudes of the rms velocity fluctuations and ri are the initial characteristic radii. The fluctuating velocity field is constructed from a realization of fully-developed homogeneous turbulence that is pro- jected onto the rms velocity distribution, as is described in more detail in Dommermuth, et al.~1997~. The rms fluctuations are initially uncorrelated and the turbulent shear stresses are zero. As discussed later, in the Wake relaxation subsection, an iterative procedure is used to relax the wake until the production of turbulent kinetic energy is balanced by dissipation. We assume that the mean and fluctuating portions of the density disturbance are initially zero. Numerical algorithm The governing equations are discretized using second-order finite-differences. A fully-staggered grid is used in the numerical simulations. Periodic bound- ary conditions are used along the sides of the compu- tational domain, and free-slip boundary conditions are 0.0001 105 10 100 1000 Figure 2: The kinetic energy for 7 < t < 1000 for Re = 105: ( ), fluctuations for Fr = oo, ( - - - - - - - - - ), fluctuations for Fr = 2.0, ( - ), mean flow for Fr = x, and ( ), mean flow for Fr = 2.0. The bold solid line represents axial similarity behavior, t-3/2, and the bold dashed line represents swirl similarity behavior, tat. The energy is normalized by U2D3. imposed at the top and bottom. A third-order Runge- Kutta scheme is used to integrate the equations with re- spect to time. The numerical algorithms have been im- plemented using high-performance fortran (PGHPF) on a CRAY T3E. Additional details and convergence studies of a similar numerical algorithm are described in Dom- mermuth, et al.~1997~. RESULTS For the stratified simulation, the Reynolds number is Re = 105 and the Froude number is Fr = 2.0. The initial mean disturbance and the rms fluctuations are based on least-squares fits of laboratory measurements of a cross section of the wake that is seven diameters downstream. For the mean axial velocity, a0 = 0.10 and r0 = 0.25 and for the mean axial vorticity an = 0.80 and ro = 0.20. For the fluctuating portion of the flow, pre- relaxation, ai = 0.40 and ret = 0.07 and r2 = r3 = 0.05. We choose a computational domain that is 12D long, 24D wide, and 12D deep. The horizontal and vertical dimensions of the computational domain are as large as computer resources allow in order to accommodate the propagation of internal waves and the spreading of the wake, which in the stratified case spreads more in the horizontal than in the vertical. In any case, the horizon- tal and vertical dimensions of our computational domain 3
are larger in most cases than the comparable tank sizes 0.1 used in laboratory experiments. Convergence is estab- lished using two different grid resolutions correspond- ing to coarse (128 x 256 x 129 grid points) and fine (256 x 512 x 257 grid points) simulations. The fine-grid results are presented here. For the non-stratified case, we use the same dimensional parameter values as for the stratified case, except that N = 0 rad/sec. Therefore, all the nondimensional parameters are the same except that Fr= x. Wake relaxation A relaxation procedure is used to establish the proper balance between production of turbulent kinetic energy and dissipation at the beginning of the calcula- tion. During the relaxation procedure, the mean portion of the flow is held fixed. The total turbulent kinetic en- ergy is also held fixed, but the spatial distribution of the turbulent fluctuations is free to vary. Once the turbulent production and dissipation are in balance, the relaxation procedure is turned off and the numerical simulation is initiated. This is the same relaxation procedure described in Dommermuth, etal.~2002) and is similar to the pro- cedure used by Orszag & Pao (1974) in their numerical simulations of a self-propelled body. Similarity In Figure 2 the kinetic energy in the mean portion of the flow integrated over the volume of fluid ({v dV (ui) (ui>) and the turbulent kinetic energy ~ ~ (TV dVuini) are plotted versus time for both the strati- fied and non-stratified cases for Re = 105. The decay of the mean and turbulent kinetic energy in the self-propelled body wake proceeds somewhat dif- ferently than for the towed body, which is described in Dommermuth, et al. (2001~. Figure 2 shows the decay of the mean and fluctuating kinetic energy for both the stratified and non-stratified wakes. Initially, for the non- stratified case, the energy in the thrust and drag wakes dominates and the energy decays as x~ 3/2. Eventu- ally, however, the swirl energy (which is much smaller initially) dominates and the entire wake approaches the swirl decay rate of x~ i. For this particular simulation, the non-equilibrium region is much greater than it is for a drag wake. The mean and turbulent energies for the strat- ified simulation decay much less rapidly than the non- stratified simulation. For the stratified simulation, self-similarity is not es- tablished until the very late wake. We believe that once the effects of stratification disrupt the turbulent produc- tion mechanism, the wake persists significantly longer than it would in a non-stratified fluid. This effect is more significant for momentumless wakes than it is for drag 4 0.01 ant= 0.001 _ 0.0001 _ .... . . . .... 10 - t 100 1000 Figure 3: The total potential and kinetic energies for 7 < t < 1000 for Re = 105. The results for stratified fluids (Fr = 2.0) are labeled: ( ), potential energy; ( - - - - - - - - - ), vertical turbulent kinetic energy; ( ), potential energy plus vertical turbulent kinetic energy, ( - ), horizontal turbulent kinetic energy, and ( - - - - - ), total energy. The energy is normalized by U2D3. wakes because momentumless wakes in the absence of stratification decay more rapidly than drag wakes. Energy is redistributed between the kinetic energy and the potential energy and also between the mean and the fluctuating portions of the flow. Figure 3 shows the total energy (E), which includes the turbulent kinetic en- ergy, the kinetic energy in the mean portion of the flow, and the potential energy (~ Jo do iV dVp u31. Fig- ure 3 also shows the total energy in the fluctuating por- ~ tion of the flow (E), which includes the turbulent kinetic energy and the potential energy. The stratified and non- stratified fluids establish self-similarity at the same rate ~ for both E and E. The results show a tendency in the far wake for the energy in the stratified fluid to be higher than the non-stratified fluid. This effect may be attributed to the generation of internal waves, which do not decay as rapidly as turbulence. The formation of pancake eddies Figure 4 show time series of the vertical component of vorticity in the horizontal plane through the wake cen- terline (X3 = 0~. Part (a) illustrates the results for the non-stratified fluid, and part (b) shows the correspond- ing results for a stratified fluid. In this gray-scale fig- ure, white represents positive vorticity with magnitude
As = 4 and black negative vorticity with ~z = - 4. Each frame has the dimensions 24D in both the cross-stream and upstream (to the left) directions. Note that the flow along the streamwise direction (A ~ has been periodically extended. The centers of each frame are located at (from left to right and top to bottom) t ~ 6, 22, 54, 86, 118, 166, 230, 294, 358, 422, 518, 614, 710, 806, 902, and 998. Each frame is scaled by the distance downstream from the initial plane ~ t = 6 ), which is the expected similarity behavior. For the stratified cases, coherent structures, in the form of nearly circular vortex patches begin to appear at t ~ 100. This corresponds to about the time the mean ki- netic energy begins to decay at the self-similar rate. Note that instabilities are evident almost immediately. Further downstream, the size of these patches of vorticity grow and the number of patches in a frame very slowly de- crease. Over the duration of the simulation, the small-scale features that are observed in the bulges of the non- stratified simulation appear to merge to form the large- scale structures that are observed toward the end of the stratified simulation. VORTEX STUDIES Coherent vortex pairs are major components of the wakes of submarines. Recently, comparisons with laboratory and field observations (Delis) and Greene, 1990; Delisi, et al., 1996; Delisi, 1998) indicate that numerical models fail to accurately predict vortex pair migration in strati- fied environments. The reasons for this failure are not fully understood. A possible reason for this failure is that the turbulence models employed in the numerical simu- lations produce inaccurate predictions of the entrainment and detrainment rates of the recirculating region of the vortex pair. Large-eddy-simulation (LES) codes, devel- oped for modelling aircraft trailing vortex pairs in the at- mosphere (Gerz and Ehret, 1997; Han, et al., 2000), have not been carefully compared with detailed measurements of vortex pair motion in stratified environments and may also not accurately predict entrainment and detrainment rates. It is clear that a more complete understanding of the physical mechanisms controlling the motion of vor- tex systems is needed. The governing parameters for this flow are the Reynolds number Re = F0/u, the vortex Froude num- ber Fr = Vo/(Nbo), and the nondimensional turbu- lence intensity ~ = (ebo)~/3/Vo of the background en- vironment, in which F0 is the vortex circulation strength, z' is the kinematic viscosity, ~ is the turbulent dissipa- tion rate, be is the initial vortex separation distance and V0 = Fo/~2~rbo) is the initial vertical speed of the vortex pair. The nondimensional turbulence intensity is the ratio of the turbulent velocity at the length scale of the vortex separation distance to the vortex descent speed. In this paper we describe two preliminary numerical experiments: the motion of a vortex pair for Re = 105, ~ = 0.15 and Fr = x and 4. The values of these pa- rameters are chosen such that the numerical simulations can be compared directly with aircraft measurements in an atmosphere with moderate ambient turbulence. These aircraft measurements are the only data we have access to for which the values of Re are near to that of full-scale submarine wakes. The numerical method is the same large eddy scheme described in the previous section for computing the wake of a self-propelled body. However, for the vor- tex simulation we have modified the sub-grid-scale tur- bulence model in an initial attempt to deal with the ef- fects of strong rotation, such as would be found in the vortex cores, on the small-scale turbulence. The flow is initialized with an approximation to a measured mean ve- locity field of a rolled up aircraft wing vortex, and a su- perimposed homogeneous turbulent velocity field in an unstratified background. A relaxation procedure is intro- duced, as described in the previous section, to allow the turbulence field to establish a balance between produc- tion and dissipation before the calculation is begun. PROBLEM FORMULATION A schematic drawing of the initial flow configura- tion for the vortex pair simulations is shown in Figure 5, showing the locations of the centers of each vortex of the vortex pair. The distance between these centers is be and the circulation strengths areF0 for the vortex on the left and +~0 for the vortex on the right. The self- induced velocity of this vortex pair is directed upward and has magnitude V0. The self-induced velocity is up- wards since we are attempting to simulate a submarine that is slightly buoyant and so its control surfaces need to generate negative lift to keep the ship in level cruise. Superimposed on this vortex pair flow field is a homoge- neous field of turbulence, whose intensity and character- istic length scale will be described below. Large eddy formulation of the Boussinesq equations We will assume that the fluid is incompressible and weakly stratified. A large eddy approximation is invoked whereby the large-scale features of the flow are resolved and the small scales are modelled. As in the previous section, let hi denote the filtered three-dimensional ve- locity field as a function of space xi (i = 1, 2, 3) and time t. The origin of the coordinate system is as shown in Figure 5. x1 is along the axes of the two vortices, x2 is transverse to these two axes, and X3 iS positive upward. The length and velocity scales of the flow are respec- s
(a) (b) Figure 4: A time history of the vertical component of vorticity wz on the horizontal plane through the center of the wake for Re = 105:(a) Fr = oo and (b) Fr = 2.0. In these gray-scale figures, white represents positive vorticity with magnitude wz = 4 and black negative vorticity with wz = - 4. Each frame has the dimensions 24D in both the cross-stream and upstream (to the left) directions. Note that the flow along the streamwise direction (x1) has been periodically extended. The centers of each frame are located at (from left to right and top to bottom) t ~ 6, 22, 54, 86, 118, 166, 230, 294, 358, 422, 518, 614, 710, 806, 902, and 998. 6
- vo 1 rO X3 -rO 5bo ~ ~ ( · , ,X2. ~ x1 1' be .1 5bo -, Figure 5: A schematic diagram of a vortex pair in a stratified fluid, showing the coordinate system and relevant parameters, the vortex separation distance ho, the vortex circulation strength F0 and the self induced vertical velocity VO. lively normalized by the initial vortex separation distance be and the initial self-induced speed of the vortex pair VO. The numerical model used for this study is the same as that described in the previous Large eddy formulation of the Boussinesq equations subsection, except that the mixed model (Bardina, et al., 1984) used to model the subgrid scale stress tensor has been modified in an ini- tial attempt to represent the suppression of sub-grid-scale turbulence due to the strong rotation in the vortex cores. Since turbulence is suppressed in the cores, they tend to diffuse only very slowly. This modified mixed model was chosen as it is a first step towards implementing a dynamic SGS model, such as the one described by Ger- mano, et al.. (1991~. For the modified mixed model, the SGS stress tensor Tij = UjUiUjU in which Hi is the fluid velocity in the ith direction and the overbar symbol denotes spatial filtering, is repre- sented as Tij= (niUiU}Ui) -Cs /\2 ~ Sij - (Sij ~ ~ (Sij where Si; is the strain tensor, S - ~ - (Sib>) 140Ui Air _ + 2 Axe taxi (6) (7) The first term, in parentheses, on the right hand side of (6) is the similarity portion of the mixed model. It pro- vides an accurate representation of the turbulent stresses. The remaining term on the right hand side is Smagorin- sky portion portion of the mixed model. It represents dissipation and is formulated so that there is no turbulent dissipation when the strain tensor is well resolved. cs is the Smagorinsky coefficient and /\ is the width of the spatial filter, which we set equal to the grid spacing. In the simulations described here cs = 0.05. The computational domain is square in the cross plane with sides of length 5bo and rectangular in the axial plane with the axial sides of length 1.5bo. As in previous studies, the smaller length in the axial direction is cho- sen so as to suppress the Crow instability, Crow (1970), (which has an axial wavelength of about 8.6bo) and to re- duce computational costs. Periodic boundary conditions are imposed at all boundaries. In the calculations de- scribed in this paper, the grid resolution is (64, 256, 256) in the (x, y, z) directions. Initialization The initial flow field consists of the combination of a homogeneous turbulent velocity field and the velocity field associated with a vortex pair. The vorticity distri- bution within each vortex is an empirical fit to that mea- sured in large-scale airplane wing vortices. The background homogeneous turbulent velocity field is generated by prescribing an initially random dis- tribution of Fourier modes and then integrating the equa- tions of motion forward in time, holding the total energy constant by adjusting the amplitudes of all the Fourier modes appropriately at each time step, until the velocity field has reached a statistically steady state. The result- ing velocity field has a well developed inertial subrange and is nearly isotropic. After this steady state has been reached, the amplitudes are no longer adjusted and the turbulence is allowed to decay. The density field is not perturbed initially. The strength of the turbulence field is measured by the parameter ~ = cl/3, where ~ is the nondimensional turbulent dissipation rate. This is a mea- sure of the ratio of the turbulent velocity magnitude at the scale of the initial vortex separation to the self-induced speed of the vortex pair. The initial velocity field associated with the vortex pair that we will use is an empirical fit to observed air- plane wake vortices developed by Proctor (19981. In this representation, the lateral components of the vorticity are zero and the axial component ~ associated with each vor- tex is given by two functions, one for inside the vortex core, r < rc, where r is the radial distance from the vortex center and rc is the radial distance to the peak tan- gential velocity, and one for outside, r > rc, the core.
180 1 4() 3 1m 80 60 40 20 10 x3 6 4 2 o 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 ( r Figure 6: The initial radial distribution of vorticity for each vortex of the vortex pair, normalized by their vorticity values at each vortex's center, plotted as a function of the radial distance, r, from the vortex center. The initial vorticity of each vortex is axisymmetric about the vortex center. For r < rc Airs= 2 t1exptlO(rc/B)3/445 x expt1.2527(r/rc)2], (8) and for r > rc ,,_ _ 7.5/ r \-1/4_~r ares ~ /r~\3/43 (9) an</-) = B<B) expL~u~r/~) ~ in which B = 4/~r and to be consistent with airplane data rc = 0.125. A plot of this vorticity distribution as a function of r is shown in Figure 6. The range of typical values of ~ in the atmosphere, as reported by Han, et al. (2000) is about 0.01 to 0.5. As a simulation representative of moderate turbulence, we chose to run the simulations described in this paper for 7' = 0.15. The turbulent and vortex pair velocity fields are su- perimposed and adjusted so that the combined field is divergence free. RESULTS Here we present preliminary results of the rate of as- cent and the energetics of the vortex pair in non-stratified and stratified background environments. Of particular in- terest is the spatial distribution of the turbulent kinetic energy, which tells us where the most intense turbulence activity is going on and is a clue to where the most mix- ing and entrainment/detrainment into the vortex recircu- lating cell is going on. 12 laniinar unstratified - - turbulent unstratified turbulent stratified ~ M'' / J 0 2 -''1 4 6 10 12 Figure 7: The trajectory of the vortex pair for laminar unstrat- ified conditions, turbulent unstratified and turbulent stratified. Vortex ascent rate The time series for the ascent of the vortex pair is shown in Figure 7. An ideal two-dimensional vortex pair in a non-turbulent and non-stratified fluid would ascend, in our nondimensional units, at a rate Z = t, in which Z is the average vertical position of the pair. As shown in Figure 7, the pair in a non-stratified but turbulent back- ground ascends at a slightly slower rate than this. This appears to be due to the diffusion of vorticity by the tur- bulence. In a stratified fluid the vortex pair has to do work against buoyancy forces and this causes the pair to rise at an even slower rate. The rise rate is somewhat de- termined by how much fluid is entrained and detrained by the pair. If this entrainment and detrainment is large enough, then the pair will have a density that is always close to the density of the surrounding fluid, but if this entrainment is weak then the pair will have a density that is heavier than that of its surroundings and will stop ascending, much as what appears to be happening near t = 12 in Figure 7 for the stratified case. The ascent rate for the vortex pair in a non- stratified fluid is consistent with the LES results of Han, et al.~2000) for the same strength of background turbu- lence. Han, etal.~2000) provide a semi-empirical for- mula for the prediction of the vortex ascent rate based on a diffusive decay of the circulation strength of the vor- tices. This formula seems to work well except when the background is strongly turbulent. By the end of the sim- ulation period (t = 12) the center of the vortex pair has drifted to the left about 0.25bo in the non-stratified case and about 0.5bo in the stratified case. Energetics Time series of the kinetic energy of the mean flow, integrated over the whole computational domain, is plot- ted in Figure 8. As can be seen, the mean energy decays 8
1 0.8 0.6 0.4 0.2 .. _ _ larrunar unstratified ' . rbulent unstratified turbulent stratified 2 4 6 t 8 10 12 Figure 8: The time series of the total kinetic energy of the mean flow. steadily with time in both the non-stratified and strati- fied cases. The kinetic energy decays most rapidly in the stratified case since besides dissipation some of the ki- netic energy is being converted to potential energy and some of the mean energy is being converted into internal wave fluctuations. Time series of the kinetic energy of the fluctuating motion in the whole computational domain is shown in Figure 9. In these time series, the fluctuating energy ini- tially increases steadily as energy from the mean flow is converted to fluctuating energy. Up to t ~ 1, which is about a quarter of a buoyancy period, there is no differ- ence between the non-stratified and stratified cases and there is not significant until t = 2 or 3, which is about half a buoyancy period. At this time it appears that the turbulent energy saturates and begins to decline in the non-stratified case and to oscillate (with a slight decay) in the stratified case. The frequency of the oscillations in the stratified case is approximately equal to the buoyancy frequency, suggesting that a good fraction of the energy has gone into internal gravity waves. The evolution of the spatial distribution of the turbu- lent kinetic energy is shown in Figures 10 and 11. Early in the simulation most of the turbulent kinetic energy is is the vortex cores in both the non-stratified and stratified cases. A little later, a significant amount of turbulent pro- duction has occurred at the top of the vortex pair recircu- lation zone and again this appears to be the same in both the non-stratified and stratified cases. After t ~ 4 the two cases start to deviate from each other substantially. In the non-stratified case, the turbulent kinetic energy re- mains concentrating in the cores as the energy decays steadily outside of the cores. In the stratified case, the turbulent kinetic energy actually increases at first in the region between the cores and eventually throughout the entire recirculation zone. 1~ o.s 0.8 0.7 0.6 ~ 0.5 I: 0.4 0.3 0.2 0.1 01 I rl[ :J turbulent unstratified | - - turbulent stratified | , ~ ~ .' 2 4 6 8 10 12 Figure 9: The time series of the total kinetic energy of the turbulent flow. Conclusions We have presented some preliminary results from using large eddy simulation to compute the late wake of a self- propelled body moving at constant speed through a non- stratified and a uniformly stratified fluid at Re = 105 as well as the motion of a pair of counter-rotating vortices in both non-stratified and stably stratified fluids. An important aspect of the simulations is the use of a relaxation procedure to adjust the initial turbulence fields so that turbulent production and dissipation are in balance. Our simulations intentionally have been initi- ated with different turbulent velocity and density fields than would be found in laboratory experiments. We have done this so that the initial conditions would be nearly the same for stratified and non-stratified simulations as well as for simulations with different Reynolds numbers. These similar initial conditions, which would be very dif- ficult to obtain in laboratory experiments, have allowed us to study how the wake behaves without the compli- cations of the initial conditions varying with the overall flow conditions. We have observed that the pancake eddies form without being driven by similarly-sized eddies in the ini- tial disturbance. In fact, these simulations produced pan- cake eddies even though the fluctuating component of the initialization began with random phase. We have seen that the eddies depend strongly on whether the body is towed or propeller driven: the propeller-driven body pro- duces eddies that are smaller and more chaotic than the nearly periodic pancake eddies in the simulated wake of a towed body. This difference is consistent with the much more rapid decay of the mean axial velocity in the wake in the self-propelled case. In the simulations of the motion of a vortex pair in a turbulent background flow in both non-stratified and 9
1 -1 1 o -1 0.4 1 01 0.2 1 -21 o 0.4 1 0.2 - o 2 -1 -21 -2 0 2 -2 0 2 I! ().2 'O Figure 10: Transverse plane (y, z) contour plots of axially-averaged turbulent kinetic energy at (from left to right and top to bottom) t = 0.2, 2.2, 4.2 and 6.2 for Re = 105 and Fr = oo. 10
o -1 1 -1 -2 0 2 0.4 10.2 Jo -2 0 -21 ~0.4 1 0' 0.2 o 0.2 910.4 10.2 Jo Figure 11: Transverse plane (y, z) contour plots of axially-averaged turbulent kinetic energy at (from left to right and top to bottom) t = 0.2, 2.2, 4.2 and 6.2 for Re = 105 and Fr = 4. 11
stratified fluids we found that the ascent rate of the pair is strongly affected by the background turbulence and by the stratification, both effects lead to a decrease in the as- cent rate of the pair. In the non-stratified case this is due to turbulent diffusion of the vorticity and in the the strat- ified case it is additionally due to the conversion of mean kinetic energy into potential energy and to the generation of internal gravity waves. The distribution of turbulent kinetic energy in the vortex pair was studied and found to be quite different in the intermediate and late times of the simulations. The non-stratified simulation shows the turbulent kinetic en- ergy mainly decaying everywhere except in the vortex cores. The stratified case shows that the turbulent kinetic energy actually increases in the recirculation zone at late times. Acknowledgements This research is supported by ONR under contract num- ber N00014-01-C-0191, Dr. L. Patrick Purtell program manager. This work was supported in part by a grant of computer time from the DOD High Performance Com- puting Modernization Program at the Naval Oceano- graphic Office Major Shared Resource Center. We thank Prof. G. R. Spedding at He University of Southern Cali- fornia and Prof. D. D. Stretch at the University of Natal for many helpful discussions. References Bardina, J., Ferziger, J. H. and Reynolds, W. C., "Improved turbulence models based on LES of homogeneous incompressible turbulent flows," Rep. TF-19, 1984, Dept. of Mechanical Engineering, Stanford University. Crow, S. C. "Stability theory for a pair of trailing vortices," American Institute of Aeronautics and Astronautics Journal, Vol. 8, No. 12, December 1970, pp. 2172-2179. Delisi, D.P. and Greene, G. C., "Measurements and implications of vortex motions using two flow visualization techniques", Journal of Aircraft, Vol. 27, 1990, pp. 968-971. Delisi, D.P., Greene, G. C., Robins R. E., and Singh, R., "Recent laboratory and numerical trailing vortex studies", AGARD Conference Proceedings on The Characterization & Modification of Wakes From Lifting Vehicles in Fluids, Vol. 584, 1996, pp. 34-1 - 34-10. Dommermuth, D. G., Gharib, M., Huang, H., Innis, G. E., Maheo, P., Novikov, E., Talcott, J. C., and Wyatt, D. C., "Turbulent free-surface flows: a comparison between numerical simulations and experimental measurements," Proceedings of the Twenty-first Symposium on Naval Hydrodynamics, Office of Naval Research, 1997, pp. 249-265. Dommermuth, D. G., Rottman, J. W., Innis, G. E. and Novikov, E., "Numerical simulation of the wake of a towed sphere in a weakly stratified fluid," submitted to Journal of Fluid Mechanics, 2002. Germano, M., Piomelli, U., Moin, P. and Cabot, W. H., "A dynamic subgrid-scale eddy viscosity model", Physics of Fluids, Vol. As, 1991, pp. 176~1765. Gerz, T. and Ehret, T., "Wing tip vortices and exhaust jets during the jet regime of aircraft wakes," Aerospace, Science and Technology, Vol. 1, 1997, pp. 463-474. Han, J., Lin, Y.-L., Arya, S. P. and Proctor, F. H., "Numerical study of wake vortex decay and descent in homogeneous atmospheric turbulence," American Institute of Aeronautics and Astronautics Journal, Vol. 38, No. 4, April 2000, pp. 643-656. Lin, J.-T. and Pao, Y.-H., "Turbulent wake of a self-propelled slender body in stratified and non-stratified fluids: analysis and flow visualizations," Report No. 11, July, 1973, Flow Research Company. Lin, J.-T. and Pao, Y.-H., "The turbulent wake of a propeller-driven slender body in a nonstratified fluid," Report No. 14, February 1974a, Flow Research Company. Lin, J.-T. and Pao, Y.-H., "Velocity and density measurements in the turbulent wake of a propeller-driven slender body in a stratified fluid," Report No. 36, August 1974b, Flow Research Company. Lin, J. T. and Pao, Y. H., "Wakes in stratified fluids: a review," Annual Reviews of Fluid Mechanics, Vol. 11, 1979, pp. 317-338. Lin, J.-T., Veenhuizen, S.D., Liu, H.-T., "Experimental data on stratified wakes for validation of wake codes," Report No. 73, October 1976, Flow Research Company. Orszag, S. A. and Pao, Y. H., "Numerical computation of turbulent shear flows," Proceedings of the Symposium on Turbulent Diffusion in Environmental Pollution, 1974, pp. 225-236. Academic Press. Proctor, F. H., "The NASA-Langley wake vortex modeling effort in support of an operational aircraft spacing system", Proceedings of The 36th aerospace sciences meeting & Exhibit, AIAA Paper No. 98-0589, 1998. Riley, J. J., Metcalfe, R. W. and Weissman, M. A., "Direct numerical simulations of homogeneous turbulence in 12
density stratified fluids," Nonlinear Properties of Internal Waves, American Institute of Physics, 1981, pp. 79-1 12. 13
DISCUSSION G.R. Spedding University of Southern California, USA Summary The recent development of large eddy simulation (LES) of stratified wakes flows at significant Reynolds number (Re) has been eagerly awaited by the community, not only in ONR-specific applications (wakes of submerged bodies), but also as a significant step forward in resolving some of the conundrums in coherent structure formation from disordered initial conditions, a process that seems quite characteristic of flows in stably-stratified fluids. The paper by Rottman et al. shows some very interesting results from simulations of the wake of a self-propelled, slender body, and for a rising vortex pair. The selection of these configurations presages a time when complex geometries, complete with control surfaces will be simulated at high Re and for long evolution times. The current results can be seen as significant milestones on the way to this goal, and they raise some intriguing and quite general questions. Pancake eddy formation from initially turbulent wakes The most famous consequence of moving a body through a stratified fluid can be seen reproduced in Fig. 4 of this paper, where the lower half of the figure shows large scale coherence and a lack of small-scale features that clearly contrasts with the upper (non-stratified) equivalent. However, careful inspection of the non-stratified result shows that it too contains coherent structure at the same scale as the dominant features in the stratified case. Perhaps then, the principal feature of the stratified result is in the loss of small scales, rather than the generation of the large-scale structures themselves? How does this happen? Perhaps vertical fluctuations are readily transported away as internal wave motions. Then why preferentially the small scales? And what happens to the small scale horizontal fluctuations? Why do they apparently dissipate so fast? Alternatively, as the authors suggest, perhaps there is preferential merging of small scales in a stratified environment. Why would that be? How are the dynamics different? If Fig. 4a is filtered through a low-pass filter, do we get Fig. 4b? The eddies themselves are described as smaller and more chaotic (in some sense) than those generated by the mean defect profile but a direct comparison is not available. If this is the case, why is it so? Does that shed some light on possible formation mechanisms? Which aspect of the initial conditions is important for this result? While some of these questions might seem abstruse and of limited practical significance, they do bear on the generality of the results. It is only through a clear understanding of the important physics that the range of applicability of a set of simulations can be known. Importance of initial conditions A number of quite specific initial conditions are detailed for both the wake and vortex pair. Velocity and length scales for mean and turbulence (and swirl) quantities are taken from laboratory experiments in the self-propelled body case and from airplane wake data for the vortex pair. While the relevant coefficients are given and apparently copied quite faithfully, it is not clear how, or if, these values are important. Some of the external data are quite unusual; the discontinuous vorticity profile in Fig. 6, or the large fluctuation magnitudes in the kinetic energy of the momentumless wake, where fluctuating quantities are 3-4 times the mean value. One of the recurrent questions of stratified fluids research concerns the relative importance of initial conditions (and, quite frequently, the search for initial conditions that are unimportant), so the results can have general application. Having a correctly functioning code, the authors would seem to be uniquely qualified to answer some of these questions, at the very least with respect to this work. What, if anything, changes if the ratio of mean defect:axial:fluctuating velocities is not 1:8:4? It is to be hoped the authors will have an opportunity to examine and report the sensitivity to initial conditions, both for numerical initialization details and for understanding the relative importance of the nhvsical field variables. Comparing averages r ~ - Comparing and interpreting data can be quite sensitive to particular averaging procedures and two examples come to mind. The first involves an unremarked comparison. A general curve describing the evolution of quantities such as the mean wake defect in the case of a towed sphere has been proposed,! whose main features can be summarized in Fig. D1.
log U. u' 3D \ -213 \' -1/4 I ·~. .. NEQ Q2D .................. ! -213 . ·-. log Nt Fig. D1. A general curve showing the low decay rate of mean wake defect (and of kinetic energy) during an intermediate non- equilibrium (NEQ) regime for a decaying stratified flow. The lasting effect is that at late times, when the flow is quasi-two- dimensional (Q2D), the wake defect or kinetic energy remain significantly higher than if a constant 3D power law had been in effect (dotted line). Fig. 2 of the paper demonstrates that stratified momentumless wakes also have a similar region of comparatively low kinetic energy decay rates, lasting even longer than in the towed-sphere case. However, some caution is required because the data from Fig. D1 above concerns averages made only over the wake itself- the averaging box enlarges as the wake itself does, and the outer ambient fluid is excluded. In the current paper, Figs 2 & 3 are averages over the entire domain. The primary effect is to change the interpretation of the contribution of internal waves to these measures. Since the wave motions propagate energy away from the wake, then this effect by itself might be anticipated to increase the decay of a locally-averaged kinetic energy. Since the opposite happens (NEQ decay rates decrease in Fig. D1) a different explanation is required. In this paper, the comparatively low dissipation rates of internal wave modes are offered to explain the low energy decay rates. This is possible, and if true, then the mechanism is different from that of Fig. D1, any similarity being merely coincidental. In any event, the two results are not strictly comparable. They could be if plots like Fig. 2 & 3 were made for local wake averages, instead of global averages. The second possible effect of averaging is in the comparison of the kinetic energy distributions of rising vortex pairs in homogeneous and stratified ambients (Figs 10 & 11), while looking across the vertical axis. Fig. 11, for the stratified case, has the intriguing result that the kinetic energy intensity appears higher in recirculating regions around the vortex cores, rather than being concentrated exclusively in the cores themselves. The physical mechanism for this is not obvious, and one wonders about the role of averaging. In the presence of the density gradient, vertical motions are constrained, eventually, by buoyancy forces. A spanwise, or axial average through some structure thus might reflect an increased uniformity in the vertical direction (z), rather than some change in the structure itself. While the almost two-dimensional motions of Fig. 10 are still free to vary locally in z, this variance might be suppressed in Fig. 11, having the effect of more closely aligning the structures, which then appear to be more concentrated. Comparisons with experimental/field data Of course extensive comparisons are not practicable in a short paper, but it will be very interesting to know how these recent and interesting data compare with other cases, from laboratory experiment to numerical simulation to field observations - if there are any. There are two strains of interest in this: (i) Using comparative studies to deepen the understanding of the important physical mechanisms that govern formation and late-time evolution of these stratified wake flows. It seems possible, for example, that wakes with both zero and nonzero momentum have some period of surprisingly low kinetic energy decay rates. Does this point to some universal feature of flows in a stratified ambient? How is it related (if at all) to the characteristic pancake eddies that are observed? (ii) The second theme is that of obtaining practical, predictive models for real applications. Is there any evidence from field tests that submerged wakes either have long persistence times or higher than expected kinetic
energies close to their source? generally, what are the practical consequences for the main results reported here? These remarks apply equally to the wakes and vortex pair results. The latter, in particular, ought to be comparable to other data, at least in homogeneous fluids. Currently, the trajectories of Fig. 7 show comparisons between a laminar, inviscid solution and turbulent simulations in homogeneous and stratified fluids Missing is the comparison with aircraft wake data for which the simulations have been set up, or with laboratory data on either laminar or turbulent viscous flows. Or, more Concluding remarks It is likely that sophisticated LES calculations such as those appearing this paper will form an increasingly important component in making practical predictions of complex flows. The combined modeling/direct simulation approach in principle is a powerful way of making calculations that do not exceed current computer capacity, but that can properly account for the important physical processes. To test and further evolve the capacities of these programs it will be essential to keep analyzing them and improving them in an environment of companion test and data from direct numerical simulation, from laboratory data and from field studies. Similarly, the performance and reliability of the LES model components themselves will be verifiable if and only if attention is paid to comparative results from varying initial and boundary conditions. The results reported here are intriguing and thought-provoking, and we look forward to more. 'Spedding. G.R. The evolution of initially turbulent bluff-body wakes at high internal Froude number. J. Fluid Mech. 337, 283-301. AUTHORS' REPLY Introduction Prof. Spedding has raised a number of perceptive questions about the physical implications of the results of our simulations. Many of these questions are of such a general nature that we will be unable to answer them in any depth here, either because we have not done the appropriate simulations, we do not have enough space or simply we do not know the answers (although occasionally we will offer some speculations). The questions that we cannot answer are all intriguing and form a good basis for conducting future research. Pancake eddy formation from initially turbulent wakes Prof. Spedding's observations and questions about the differences between the wakes in non- stratified and stratified fluids are the main issue. We cannot fully answer these questions here, much more research is needed, but we will attempt some speculations. We agree that the non-stratified wakes in our simulations develop large-scale structures that are similar in some ways to the large-scale structures in our simulated stratified wakes. We presume that in stratified wakes the vertical fluctuations propagate away in the form of internal waves leaving only the large-scale structures that have nearly horizontal motion only. We suspect that similar large-scale structures exist in the non-stratified wake, but that they have no preference for being horizontal. It does seem that there is some preferential merging in the stratified flow that must explain the lack of small scales in this flow, but we have no quantitative theory for this at the present. The possible formation mechanism for the pancake eddies is still in dispute. We have found that the ultimate structure of the stratified wake is fairly insensitive to the initial conditions. The initial conditions have a strong influence on the walce at intermediate times, but for long enough time we always get the same results. Importance of initial conditions Prof. Spedding's questions about the importance of initial conditions are of very practical relevance. As stated in the previous section we found that the more accurate the initial conditions, the more accurate were our results for small and intermediate times. This is particularly true of the turbulent stresses. As for the eventual or long-time state we found that the initial conditions were not very important. Of course, we have drawn these conclusions based on a rather limited number of simulations and we hope to do more simulations in the future to gain further insight into this issue. Comparing averages We agree with Prof. Shedding that our averaging procedure is different than that used to obtain the A ~ - ~ ~
decay laws in his figure D1 and therefore the comparisons we made with these laws is not completely fair. Our averaging technique would include the energy contained in the internal wave field which would not be included in the local averaging procedure that Prof. Spedding uses to analyze his experimental results. However, we claim that the internal wave energy will be approximately constant beginning sometime early in the NEQ regime. If this is indeed the case, then the slopes of the curves shown in fig D1 should be the same using either kind of averaging although the magnitudes will not. This is speculation, so we are interested in doing some simulations to see how valid this assumption is. Prof. Spedding's speculation about the effect of our axial averaging on the difference of the distributions of turbulent kinetic energy in the non-stratified and stratified cases may well be true and is something that we will investigate in the future. Comparisons with experimental/field data Prof. Spedding first asks if there is any evidence from field tests that submerged wakes have long persistence times and higher than expected kinetic energies. We have no knowledge of such full-scale field tests that are available in the open literature. He then asks about comparisons of our vortex pair simulations with laboratory data and with data from full-scale aircraft generated vortex pairs. We have begun to make these comparisons and will present these comparisons in the near future.
