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 Direct numerical simulation of surface tension dominated and non-dominated breaking waves A. lafrati, E.F. Campana (INSEAN - Italian Ship Mode] Basin, Italy) ABSTRACT Surface tension effects onto the two-dimensional wave breaking flow produced by a hydrofoil moving beneath the free surface are investigated. The study is carried out numerically by a finite difference approach which solves the Navier-Stokes equations. The air-water in- terface is embedded in the computational domain and it is captured via a Level-Set technique. A heterogeneous unsteady domain decomposition approach is used, al- lowing to focus the computational effort onto the free surface vicinities, while the flow about the body is ap- proximated by a potential flow model. Surface tension effects are investigated by progres- sively reducing the length scale, while keeping Froude and Reynolds number constant. Different flow regimes are recovered, ranging from intense plunging jet, even- tually resulting in large amount of entrapped air, up to a micro-breaker, in which case air entrapment is sup- pressed and the jet is replaced by a bulge growing on the wave crest. At this scale, the surface tension is re- sponsible for the large curvature at the toe of the bulge, and when the bulge slides upon the forward face of the wave, an intense shear layer develops from the toe. In- stabilities of this shear layer are observed, and, when increasing the Reynolds number, the shear layer breaks- up into coherent vortex-structures that interact with the free surface, eventually leading to the formation of large surface fluctuations which propagate downstream. INTRODUCTION In the present paper, attention is drawn to two as- pects of the wave breaking process: the role of the sur- face tension on breakers of progressively smaller length scale and the effects of viscosity on short waves produc- ing micro-breakers. Among the large body of works on breaking waves, many papers have been studying surface tension effects theoretically, experimentally or numerically. Tulin (1996) presented a numerical investigation on the effects of sur- face tension on the jet development. In the study, car- ried out with the aid of a potential flow solver, waves of different length are followed along their evolution to- ward the breaking. The analyzed length scales range from ~ = x (corresponding to the case of negligible surface tension) down to ~ = 25 cm. As the wave- length is decreased, surface tension forces are relatively increased and modifications occur to the jet: initially its tip is rounded and finally the jet is suppressed and re- placed by bulge on the crest of the wave. The growing of the bulge was also shown through some calculations by Longuet-Higgins (19961. Based on potential flow assumptions (with and without surface tension) he also showed the appearance of a train of par- asitic capillary waves upstream the leading edge of the bulge (referred here as the toe). Within the same po- tential flow assumptions, surface tension effects on the bulge-capillary system have been numerically addressed also by Ceniceros & Hou (1999), where results for a wide range of surface tension coefficients are reported. Furthermore, by using a modified boundary integral for- mulation for water waves to include weak viscous ef- fects, they show that, depending on the viscous coef- ficient, the train of these capillary waves may be sup- pressed. That these capillaries might be a source of vortic- ity shedding was initially predicted by Longuet-Higgins (1992~. Numerical and experimental evidence of the vorticity field generated by capillary ripples have been provided by Mui & Dommermuth (1995) and by Lin & Rockwell (1995), respectively. In a successive paper, Longuet-Higgins (1994) sug- gested that the vorticity shed by these capillaries (named Type 1) may also explain the unexpected appearance of longer capillary ripples above the toe (named Type 11), propagating downstream, experimentally observed by Dun- can et al. (1994~. In his theory, Longuet-Higgins as- sumed that these Type II downstream ripples might be primarily caused by instabilities of the shear flow in- duced by the vorticity shed by the parasitic capillaries.
Duncan & Dimas (1996), performing linear stability anal- ysis of a theoretical velocity profile established that the downstream ripples are primarily generated at the breaker. In two more recent papers, Duncan et al. (1999) and Qiao & Duncan (2001), by performing an experimental study on gentle spilling breakers induced by the wave fo- cusing technique, observed that when the bulge becomes fairly steep the toe began to move, sliding down on the forward face of the wave. In these papers, it is also ar- gued that the downstream ripples are very likely induced by an instability of the shear layer generated between the downslope flow of the bulge and the underlying upslope incoming flow. In the present paper, two-dimensional breaking waves, produced by the motion of a submerged hydrofoil, are numerically investigated by solving the unsteady incom- pressible Navier-Stokes equations. The flow close to the free surface is described as a two-phase flow and the Level-Set method is used capture the interface be- tween air and water. The computational effort is re- duced by taking advantage of a heterogeneous domain decomposition approach previously developed (Iafrati et al. 2000, Iafrati & Campana, 2002~. A set of Froude-scaled numerical simulations is pro- duced by varying the speed and the chord of the hydro- foil, so that breakers of different length scales are ob- tained. With this mechanism, the relative importance of viscosity, gravity and surface tension is altered and a substantial range of variation of some of the character- istics phenomena of unsteady breakers is modeled and studied. As to the effects of the surface tension on break- ers of progressively reduced length scale, starting from a scale at which the flow is strongly dominated by the gravity and the perturbed free surface evolve in a plung- ing breaker, different breaking regimes are recovered, up to a very short breaking wave. With this set of data, the modifying action of the surface tension on the jet formation and on the entrapment of air is discussed and differences in the vorticity generation mechanisms are shown. Finally, for the smallest simulated scale, effects of the viscosity are also investigated. At this scale, sur- face tension is highly dominating and it is responsible for the large curvature at the toe of the bulge. Numerical results show the development of an intense shear layer from the toe. At the lowest Reynolds number, weak in- stabilities of the shear layer are observed, responsible for extremely small downstream ripples on the free sur- face. However, being the shear layer still quite stable, it does not break up into separate vortices. When increas- ing the Reynolds number, the instabilities grow, and the shear layer break up into coherent vortex-structures. The strong interaction of these coherent structures with the free surface leads to the formation of large free surface ripples, whose amplitude and wavelength is governed by the vortex-blobs. Details of the evolution of the co- herent structures also shown that the combined action of the vorticity field and of the surface tension produce a cusp-shaped wave pattern, and secondary separations from the cusps are also observed. NUMERICAL MODEL Domain Decomposition Technique Wave breaking is usually responsible for a highly com- plex flow in the free surface vicinities: the air-water in- terface assumes a complicated topology, eventually lead- ing to jets, drops, air entrapment, and an intense vortic- ity field is generated either by the impact processes or by the action of strong velocity gradients. In the attempt of accurately describe complex flows, several numerical models, able to deal with a substan- tial two-phase flow, have been developed and are now at hand (a survey is provided by Scardovelli & Zaleski, 1999~. However, the cost of a detailed computation of the breaking and post-breaking evolution is still very high. On the other hand, in absence of breaking, the wavy flow generated by bodies moving beneath the free sur- face is rather accurately predictable by much simpler potential approaches, in spite of the strong assumptions made about the flow features in the body vicinities (Iafrati & Campana, 19981. Hence, when attention is mainly concerned with the flow details in the free surface re- gion, the use of a domain decomposition approach, that make use of the most suitable assumptions according to the flow region, sounds very attractive. I I r ~~ ~ I T [ I I I I I I I I I I I I T=T 1 ¢1~t 1 $~ N-S domain , . . I I I I , , I , I I I I , I , , I i j I I I I j I I I I j I I I I Matching surface ~ I I 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 _ I ~ I T I r 1 1 i I I 1 I T I I 1 1 ~ r(c~ :~ - BEM domain Figure 1: Sketch of the body and free surface domains used for the domain decomposition. On the basis of the above considerations, an unsteady heterogeneous domain decomposition approach has been developed to tackle the wave breaking flow induced by a submerged hydrofoil moving beneath the free surface (Iafrati & Campana, 20021. In the free surface region, a viscous flow model with an interface capturing tech- nique is adopted, while, in the body region, a potential
flow approximation is used. A coupling procedure is developed, allowing an exchange of information with- out the need of overlapping between the domains. In the body region the flow is governed by a Laplace equation for the velocity potential ~ Van = 0 (1) which satisfies Neumann boundary conditions at the in- flow ~QB, at the outflow PRO, at the bottom of the channel ~QB and all along the body contour Arc (see Fig. 1) while the velocity potential at the matching line ~ is assigned by the integration in time of the unsteady Bernoulli's equation: ~3'P = 1 (63Q~ ~ (2) - ~n Am: on 1 (0Qo ) (3) = 0 (0QB U SAC) (4) ~ = ~ ( 8~) dt (r) (s) where '; = gx2 2 In the above equation Q is the fluid density, 9 is the accel- eration of gravity, x2 is the vertical coordinate oriented upwards with x2 = 0 at the still water level, u is the local fluid velocity and pB iS the pressure acting on the body domain which, by enforcing the continuity of the normal stresses at the matching line, is given by p = p 2,ll In ~ (6) n being the unit vector normal to the matching line, ,u the local fluid viscosity and pF iS the pressure field at the matching line provided by the solution of the Navier- Stokes equations in the free surface domain. Further- more, a steady Kutta condition is applied at the trailing edge of the hydrofoil to properly account for the vortex shedding. More details concerning the coupling strat- egy and its validation are reported in Iafrati & Campana (2002~. Navier-Stokes solver The two phase flow of air and water is approximated as that of a single fluid with density and viscosity vary- ing smoothly through the interface. By assuming the fluids to be incompressible, the continuity equation in generalized coordinates simply reads: Hum 0 (7) where ,~ Hi i (8) Is the volume flux normal to the (m iso-surface and J-i is the inverse of the Jacobian. According to the smooth variation of fluid proper- ties through the interface, surface tension effects are in- cluded by using a continuum model, as suggested by Brackbill et al. (1992) and recently employed also by Sussman & Puckett (2000~. As a result, the momentum equation, in non dimensional form, is: (J suit + ~: (Umui) Q 0(m ( taxi ) 7-i In _ ~ ~ (J i0(m0(d)) Q Re 0(m ( 06t ~~ ) where hi is the ith Cartesian velocity component, d is the Kronecker delta and Or , Re = UrLrQw We = Ur>/~ low ' ~ (10) are the Froude, Reynolds and Weber numbers, respec- tively. Here, Ur and Lr are reference values for veloc- ity and length, cr is the surface tension coefficient while Qw, qw are the values of density and dynamic viscosity in water and are used as reference values. In (9) Gem J_~ (J(m 05t Bm~ji = J_~ 0fm 6341 Xj ~Xj ~Xj Phi (1 1) are metric quantities, ~ is the local curvature and H(d) is the Heaviside function with the distance function d be- ing positive in water and negative in air, so that H(d) = 1 in water and H(d) = 0 in air. The system of the Navier-Stokes equations is solved by using a finite difference method on a non staggered grid. The grid stencil and the numerical approach are similar to those suggested by Zang et al. (19941: carte- sian velocities and pressure are defined at the cell centers whereas volume fluxes are defined at the mid point of the cell faces. A fractional step approach is employed: the momentum equation is advanced in time by neglect- ing pressure terms (Predictor step) whose effects are successively reintroduced by enforcing the continuity of the velocity field (Corrector step). The diagonal part of the dominating diffusive terms, i.e. that originated from Vu, are computed with a Crank-Nicolson scheme, whereas all the other terms are computed explicitly. With respect to Zang et al. (1994), a three-steps Runge-Kutta scheme (Rai & Moin 1991) is adopted here. The grid be- ing fixed in time, the discretized form of the momentum equation at the step n is
· Step l · Step 2 · Step 3 (J-1o~l/\tDI) fur USA= Y1 [C(uin) + DE(tlin) + Ti(dn)] + 2~1 [J-lFr2 +DI(0i ) ~1 ^1 Ri(¢l) HiHi =)lQlJ-1 (J 1 _ or2/\tDI) ~ ~ ~~ t) ~ )2 [C(Uil ) + DE (Gil ) + Ti (d1 )] + (1 [C(uin) + DE (Olin) + Ti (dn)] + 20~2 [J 1Fr2 + DI (pi)] Ui2Ui2 = )2 ~2( -1 + (1 el~-l (J-1cx3 /\tDI) (ut up ~ = )3 [C(Ui ) + DE (ui2) + Ti (d2)] + (2 [C(Uil ) + DE (Ui1 ) + Ti (d1 )] + 2(X3 LJ 1Fr2 + DI (up) un+1 _ ui3 = y3 Hi ~ ~ J_ i) + (2 e2 J- l The coefficients oli, Hi, h are reported in Rai & Moin (1991) and in literature cited therein. In the above equa- tions, for the sake of clarity, a compact notation is used to represent the convective, diffusive and surface tension contributions: C(ui) = d9( (Umui) ~ ( ) e Re i9(m ( t96~ ) DE(ui) = R ~: (/lGml 5fi +~Bmkji'~:j) m 76 | , eWe2 Am (J Taxi H(d)) , while 0(m ( Act ) (12) is the gradient operator in generalized coordinates. In the corrector steps, ~ is the pressure corrector term which is found by enforcing the continuity of the velocity field at the end of the substep, (Rai & Moin, 1991, Kim & Moin, 1985~. The procedure is as follows: once the intermediate velocity field is found, say hi, the fluxes associated to this velocity field (U~ ~ are com- puted by (8) at the mid point of the cell faces through a second order upwind scheme (QUICK). In terms of fluxes, the corrector step can be written as UmUm 'Yt ( pt 0:j ) ( Gm~ i90 - 1 ) (13) so that, by applying the continuity (7) to Ut, the fol- lowing Poisson equation for the pressure corrector is ob- tained: ~ ~ Gmi d9¢~l N\ 1 BUT = _ {3(m ~ Q (j ~ /\t 0(m (-1 ~ {Gmi 04~-~N )~ Gym ~ 0~~1 Aim J (14) When the velocity field is assigned throughout the bound- ary of the computational domain, (13) provides Neu- mann boundary conditions for the solution of the Pois- son equation (14~. The pressure corrector term is related to the pressure field by the following equation: Ri spy ) = (A J 1 _ cot i\tD~) ( pt ,[_1 ) , (15) but, since the solution of this equation is not straightfor- ward, usually, an approximate pressure field is obtained as (Rosenfeld et al. 1991~: Ripply ~ Ri(~) =' pi = 0' + O(~\t) . (16) The system of equation discussed above is spatially discretized by a central finite difference approach, sec- ond order accurate. As already stated, a second order up- wind scheme (QUICK) is used to evaluate volume fluxes at the cell faces. At each substep of the Runge-Kutta scheme, the momentum equation is solved by using an approximate factorization approach of the diffusive part, as suggested by Kim & Moin (1985~. The time step
is chosen so that the Courant number is always smaller than ~ and the stability constraint required by surface tension is satisfied (Brackbill et al. 1992) 1 At < We>/ (i ,,, 0~) /\x3 Since not all the viscous contributions are treated implic- itly, further limitations to the time step can be required. A multigrid technique is adopted for the solution of the Poisson equation for the pressure corrector term, which is the most expensive part of the computational procedure. A corrector scheme is used for restriction and prolongation (Brands 1992) and a LSOR method is employed as high frequency smoother. Since met- ric quantities and fluid properties appear into the coef- ficients of the Poisson equation, a simple average of the metric and of the distance function is used in the restric- tion phase. Free surface capturing via Level-Set technique The air-water interface is captured as the zero level-set of a signed normal distance from the interface dart) which, at t = 0, is initialized by assuming d > 0 in water, d < 0 in air and d = 0 at the interface (Sussman et al. 1994~. The fluid property f is assumed to vary from the air to the water values as follows: ~ fa f (d) = ~ i~ + f 6~ + Em Id sin ( 2~5 ) ~ fw if idl < ~ if dad; (17) Here ~ is the half width of the transition region needed to to evaluate derivatives of fluid properties into the gov- erning equations (Iafrati et al., 2001~. Accordingly, the same smoothing is applied to the Heaviside function, thus obtaining (Sussman & Puckett, 2000) ~ ~ 2 + 2 sin (25) if Ids < ~ (18) with Had) = 0 if d <~ and Had) = 1 if d > 5. During the motion, the distance is transported by the flow, that is the equation Mu Vd=0 (19) is integrated in time to update the distribution of the function d and then to follow the interface as the level- set d = 0. The integration is carried out with the three-steps Runge-Kutta by using the same discretization scheme employed for the convective terms, that is do = dt-i + ~~\tC(<d~-~ ) + (_~ AtC(d~-2) (20) with do - dn, d3 _ d~+i and C(d7) = ~ (Ut dim . When studying free surface waves, disturbances propa- gate towards the lateral boundaries of the computational domain and, to avoid spurious reflections, a numerical beach model is introduced in the transport equation (19) ~1; =u Vdv(<d +x2) . (21) The model is applied onto two beach regions close to the lateral bounds of the computational domain. The damping coefficient v is zero at the inner limit of the beaches and grows quadratically toward the boundaries. To keep constant in time the width of the jump re- gion, the distance function is periodically reinitialized by computing, at each cell center, the minimum dis- tance from the interface. Several efficient procedures have been developed to this aim by Russo & Smereka (2000) and Sussman & Fatemi (1999), among others. Here a direct approach is used which computes the exact distance from the reconstructed interface (Sussman & Dommermuth, 2000~. Attention being mainly focused on two-dimensional applications, the computational ef- fort requested by a direct reinitialization of the distance is not prohibitive (Iafrati & Campana, 20021. For the if d <~ sake of saving the computational effort, the reinitializa- tion of the distance function is carried out only in a nar- row band about the interface and at the boundaries of the computational domain (Sethian, 1999~. APPLICATIONS Parameters of the simulations The numerical model described so far is applied to study the wave breaking flow generated by a NACA 0012 hydrofoil moving, with an angle of attack Al = 5°, be- neath the free surface. By following the experimental results by Duncan (1983), wave breaking conditions are recovered by setting the depth of the hydrofoil to be 0.783 times the chord. In order to investigate the role of surface tension onto the wave breaking dynamics, repeated numerical simulations have been carried out progressively reduc- ing the test scale, starting from a length scale compara- ble to the one used by Duncan (1983) in his experiments. The chord and the speed of the hydrofoil are reduced so that Froude and Reynolds numbers remain constant while the Weber number is decreased. In all cases it has been set Fr = 0.567 and Re = 1000, this latter being considerably smaller than the corresponding experimen- tal value Reedy for the cases A and B. Viscous effects
| case | U(cm/s) | L(cm) | Re | Reedy | We | eventually suppressing, the jet formation, which is fi- A 79.40 20 1000 158720 42 natty replaced by a bulge (Tulin, 1996). Consequently, B 39.70 5 1000 19840 21 the entrapment of air into the water is also reduced. C 19.85 1.25 1000 2480 10.5 D 19.85 1.25 2480 2480 10.5 Table 1: Parameters used for numerical simulations. Re and Reeve denote the Reynolds numbers of the compu- tations and of the experimental value at corresponding scale, respectively. on the wave breaking dynamics have been also inves- tigated at the smallest length scale (case C), by repeat- ing the simulation with the corresponding experimental Reynolds value (case D). A complete description of the set of parameters used for the computations is reported in Table 1. Concerning the set of data employed in the computa- tions, the fluid domain extends from x = - 15 to x = 15 with the matching surface located at y = - 0.2 and the top boundary is at y = 0.4. The grid has 768 x 192 cells, with /\y = 0.0025 and /\x ~ 0.0027 in the breaking region. The density and viscosity jumps are spread on a stripe which half thickness is ~ = 0.02. It is worth to notice that, to reduce the formation of forward propagat- ing waves generated by an impulsive start, a sinusoidal ramp is used to accelerate the hydrofoil up to the final speed which is reached at t = 10. Effects of surface tension on the wave breaking de- o velopment The unsteady process that results into the breaking wave formation and establishment is numerically computed for three different length scale at the same Froude num- ber, thus varying the role of surface tension effects onto the free surface dynamics, compared to inertial and grav- ity terms. To clearly describe how deeply surface tension af- fects the wave breaking establishment in the three cases A,B,C, pictures of the resulting free surface shape and of the vorticity field are shown in Figs. 2,3,4 at three different stages of the breaking process. Since, due to the baroclinic contribution, the behavior of the vortic- ity contours into the transition region can be misleading, in the figures three different density contours are shown denoting the air and the water values together with the mean level. In this way the thickness of the transition region used in the computations can be also argued. The comparison among the three different processes, clearly reveals that, when reducing the length scale, the dominant effects of surface tension progressively reduces, net ._ 0.1 -0.1 on C o. , - , .,, .,,, .., .,,,, .,,,,,, .,,,,, . I,, I, 1.2 1.3 1.4 1.5 1.6 1.7 1.8 x Figure 2: Density (black lines) and vorticity contours at three different stages of the wave breaking establish- ment in case A. The three different density contours denote the air, average and water values, respectively. Beside a significant change in the wave breaking de- velopment, a substantially different vorticity production mechanism is also evident. At the largest scale (case A), when the jet develops, vorticity contours reveal a weak shedding of vorticity in the water. In fact, in the under- side region of the plunging jet, the curvature is rather small so that the flow is able to follow the free surface shape. After the first plunging, air is entrapped and, as a consequence of the momentum exchange due to the
impact, a significant amount of circulation is produced into the water. Successive stages of the motion show that vorticity generated in this manner is then convected with the air bubble. The sequence in Fig. 2 clearly shows the process which leads to the air entrapment and to the formation of a highly rotational region that is convected downstream after wards. Due to the impact, a splash-up is then formed, thus leading to a new, thinner jet and to a second, weaker impact. 0.2 0.1 n o n ; 1.2 1.3 1.4 1.5 1.6 1.7 1.8 x In the case B. due to the reduced length scale, the increased role of surface tension leads to a significant reduction in the jet intensity and in the entrapment of air (Fig. 3~. Water-particle velocities in the jet are also reduced, and the impact is less violent. Furthermore, a large increase of the curvature at the jet root is found, which leads to a flow separation and to the development of a strong shear layer. The mechanisms and the magni- tude of vorticity production are already sensibly differ- ent from those of case A. Figure 3: The same as in Figure 2 but for case B. In Fig. S. a close up view of the local velocity field during the formation of the first splash up is shown. At the beginning, the velocity field suggests that the main portion of the fluid comes from the jet while, in a later stage, most of the fluid arrives from the incoming flow under the forward face of the wave. 0.2 F 0.1 . no 0 1 o -0.1 0.2 n -0 1 Figure 4: The same as in Figure 2 but for case C. To produce test conditions in which surface tension has a dominant effect on the evolution of the breaking wave, a further reduction of the length scale is finally used (case C). Results shown in Fig. 4 reveal that the
formation of the jet is now completely suppressed and replaced by a bulge which forms at the crest. In a succes- sive stage, the toe of the bulge starts its motion sliding on the forward face of the wave. After an initial accel- eration, the toe reaches an almost constant speed. As this motion continues, flow separates due to the sudden change in the free surface slope at the toe, thus originat- ing an intense shear layer that propagates downstream, remaining closely beneath to the free surface. Air en- trapment is completely suppressed. 0.05 , n -0.05 0.05 o -0.05 1.4 1.5 1.6 x Figure 5: Plots of the density contours and of the ve- locity field (every other vector is shown) at two stages of the jet plunge (case A). A careful inspection of the results shown in Fig. 4, reveals that no capillary waves appear on the free surface ahead of the bulge, although, for short waves, the large curvature at the toe is often responsible for the formation of a train of parasitic capillary waves (Longuet-Higgins, 1992~. The lack of these capillary waves in the numer- ical results might have several explanations. Actually, Lin & Rockwell (1995) experimentally show that this pattern only occurs in a rather narrow range of Froude numbers, and that a slight variation of speed deeply af- fects the capillary pattern, reducing the number of crests up to their disappearance. I............... 60 50) 40 1n 20 1R 16 14 12 10 30~ ,:= .J ~......................................................... V . . . . . . O 0.5 1 1.5 2 2.5 3 X 3.5 4 8 1 1.5 2 2.5 x 3 3.5 4 Figure 6: Time sequences of free surface profiles ob- tained for case C. Below a close up view about the time at which the bulge takes its foremost position for the first time is shown. The At between two profiles is 0.5 in the upper figure, 0.1 in the lower one. The same quantity is used to shift vertically two successive pro- files. For the clarity, a vertical scale factor is applied. Another possible explanation for the lack of capil- laries upstream the toe can be related to the damping effect played by viscosity on these small scales (Ceni- ceros & Hou, 19991. Finally, it has to be remarked that, due to the use of a continuum model (Brackbill et al.,
1992), the surface tension contribution to the momen- tum equation is spread on a transition region (18) and this can play an important role on the accurate descrip- tion of capillary waves having amplitude comparable or smaller than 5. 40 35 30 20 10 0 0.5 1 1.5 2 x 2.5 3 3.5 4 Figure 7: Time sequence of the free surface profiles obtained in case D. Here At = 0.1 is the time delay and the vertical shift between two successive profiles. Effects of viscosity on short breaking waves Wave crest profiles history of a long time simulation is shown in Fig. 6 for the case C, in a frame of reference fixed with the hydrofoil. For the sake of clarity, a ver- tical displacement is applied at each profile. The pro- file history shows the initial steepening of the crest and the growing of the bulge. Next, the bulge starts to slide down the forward face of the wave. After the first large downslope movement, the toe reduces its forward speed, reaches its foremost position and begins to retreat. A damped oscillatory motion is then established with reg- ular period. The cause of this back and forth movement can be ascribed to the start from rest, as discussed by Duncan (1981), whilst the damping is due to the action of the numerical beaches posed at both sides of the com- putational domain. A more refined observation of the profile history re- veals some interesting details. Indeed, while the damped motion appears to be quite smooth and the surface pro- files highly regular, a closer inspection of the crest pro- files shows the presence of small downstream propagat- ing surface fluctuations (lower part in Fig. 6~. These fluctuations are produced each time the toe experiences its downslope movement and their amplitude is extremely small. To get more insights about the nature of these ripples and the mechanisms responsible for their propagation, the numerical simulation at the experimental Reynolds number, Reese = 2480 (case D), is performed and the resulting profile history is shown in Fig. 7. Several sig- nificant differences with respect to the case C can be observed: (i) the first downslope motion appears to be much faster; (ii) the bulge reaches an upstream fore- most position; (`iii) much larger downstream propagat- ing ripples appears, which makes much more evident the recurrence of their formation each time the bulge is ad- vancing; (iv) ripples seem to be produced even when the bulge is at rest; (v ~ the oscillatory motion of the bulge in more irregular and its amplitude is damped more rapidly. By carefully looking at the free surface profiles, it can be noticed that, as soon as they appear, ripples ex- hibit a growth of their amplitude and wavelength as they move downstream. Furthermore, initially they propa- gate downstream with a velocity which is about 0.6 times the incoming flow speed. The occurrence of downstream propagating ripples was experimentally observed by Duncan et al. (1994) and theoretically investigated by Longuet-Higgins (1994), where they were denoted as Type II capillaries to distin- guish them from the parasitic capillaries waves (the Type I capillaries), located in front of the toe and propagating upstream. In Longuet-Higgins (1994) it is speculated that Type II capillaries are shear layer instability waves and that their generation is due to the vorticity shed from the Type I capillaries. Duncan et al. (1994), by using a wave focusing technique to induce breaking, argued that the ripples are generated by the instability of the shear layer, produced between the gravity induced downslope flow and the incoming upslope flow, which separates at the toe. Linear stability analysis of measured velocity pro- files (Duncan & Dimas, 1996) provided a model which is consistent with the latter assumption. Present results show the occurrence of Type II rip- ples although Type I are absent, which agrees with the Duncan et al. (1994) hypothesis. To better understand this point, instantaneous distributions of the vorticity are analyzed in the next section.
-0.25 -0.25 -0.25 -0.25 -0.25 1.~ Z ~.5 -0.25 Figure 8: Sequence of the free surface profiles and of the vorticity contours during the first upstream motion of the bulge. They refers to the interval t = 14 to t = 16.2 with time step At = 0.2. The development of shear flow instabilities and the corresponding formation of free surface ripples is clearly evident in this case. Vorticity/free surface interaction The main features of the complex interaction between the vorticity shed from the toe and the free-surface can be drawn by the sequence shown in Fig. 8. There, the later stage of the toe motion between t = 14 to t = 16.2, is shown. During this time, the bulge front moves from x = 1.30 up to x = 0.95 and free surface ripples are established. The instantaneous vorticity fields illustrate the ini- tial growth of the instabilities of the shear layer formed at the toe. Then, instabilities develops into separated co- herent structures which strongly interacts with the free surface. Results shown in Fig. 8 thus reveal that ripples are traces on the free surface of the underlying vortex- blobs. As soon as ripples are formed, their propaga- tion speed is the same at which vortex-blobs are con- vected into the water. Furthermore, the amplitude and the wavelength of the free surface ripples is essentially governed by the growth of the size of the vortex-blobs as they are convected downstream (see the last frames of Fig. 81.
Free surface ripples display a highly non-symmetric shape with respect to the horizontal axis: crests are smooth rounded while, due to surface tension, troughs are cusp- shaped. Under the action of the flow field induced by the vortex-blobs, the vertical symmetry is also lost with the free surface eventually developing an overturning mo- tion in the downstream direction (Fig. 9~. Due to the large curvature at the free surface during this stage, a flow separation from the troughs also takes place, as it is shown in Fig. lO. 0.06 0.04 0.02 o -0.02 -0.04 -0.06 -0.08 1.6 1.65 1.7 1.75 1.8 1.85 1.9 1.95 2 x Figure 9: Free surface profile at t = 15.8 for case D. A significant asymmetry with respect to the vertical axis is induced by the interaction with the vortex-blobs. Shear layer analysis To characterize the shear layer originating at the toe, cuts of the velocity field at different stages of the mo- tion of the bulge are made. In particular four different vertical cuts of the horizontal velocity component, my), are performed when: (i) the bulge is still close to the crest but it is sliding down along the forward face; (ii) the bulge takes its maximum forward position; (iii) the bulge is moving backward; (iv) the bulge takes its max- imum backward position. To make them comparable, cuts are made at the same horizontal distance Ecus from the toe. The corresponding plots are reported in Fig. l 1. To have a rough estimate of the intensity of the shear layer, the total velocity defect, between the underlying flow and the flow at the crest, and the corresponding thickness, along which the shear takes place, are eval- uated and their ratio reported in Table 2. As expected, maximum and minimum intensity of the shear corre- spond to forward and backward motion of the toe, re- spectively. Possible correlations among the values of this ratio and the occurrence of shear layer instabilities is still un- der investigation. However, a comparison between the velocity profiles measured during the forward motion of the bulge at Re = lOOO and at Re = 2480 reveals that, in the former case, the total velocity defect is smaller and the thickness of the shear layer is larger (Fig. 12), thus resulting in a value of the ratio which is 230s-l and this can explains the rather small size of the surface ripples observed in case C. O.1 -0.1 : 1.6 1.8 X 2 Figure Jo: Close up view of the density and vorticity contours of the last frame of the sequence in Fig. 8. Due to the strong curvature, a secondary flow separation is induced at the wave troughs. -0.04 -0.06 -0.08 -0.1 -0.12 -0.14 -0.16 -0.18 -0.2 -0.4-0.2 0 0.2 0.4 0.6 0.8 1 1.2 u Figure ll: Plots of the vertical cuts of the horizon- tal velocity component taken at a horizontal distance recut = 0.1 from the toe (case D). The four different cuts refer to different stages of the bulge movement: red the bulge is moving forward; green the bulge has reached its maximum forward position; blue the bulge is moving backward; pink the bulge has reached its maximum backward position. Phase | i\u/b(s-~ ) | 1. Forward motion 350 2. Maximum forward position 300 3. Backward motion 200 4. Maximum backward position 230 Table 2: Ratio between total velocity defect and the thickness of the shear layer at different stages of the bulge motion for case D (Fig. l l).
o -0.02 -0.04 -0.06 -0.08 -0.1 -0.12 -0.14 -0.16 -0.18 -0.2 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1 1.2 u Figure 12: Comparison of the vertical cuts for case C (red) and case D (blue). The significant increase in the shear layer intensity may explain the much stronger vorticity production and, as a consequence, ripples for- mation observed in case D. CONCLUSIONS In this paper the effects of surface tension on the wave breaking dynamics and the effects of viscosity on micro-breakers have been numerically addressed. By using a domain decomposition approach, an accurate description of the breaking and post-breaking evolution has been obtained by solving the Navier-Stokes equa- tions only in a narrow region encompassing the free sur- face. Results on different length scales are compared, and effects of the increasing relative importance of the sur- face tension are fully recovered. The use of a viscous flow model has also allowed to identify changes in the vorticity production mechanism. In particular, for low surface tension, vorticity into the water is mainly intro- duced by air-entrapment while, for large surface tension, a intense shear layer develops at the toe of the bulge. On the smallest analyzed length scale, the evolu- tion of a micro-breaker has been followed at two differ- ent Reynolds numbers. Downstream propagating rip- ples, scarcely visible at the lowest Reynolds number, have been clearly observed when Re increases. In both cases, shear layer instabilities, generated between the in- coming flow and the downsloping flow induced by the motion of the toe, are responsible for the appearance of these ripples. For the higher Reynolds number, shear layer instabilities develops into intense vortex-blobs; the strong interaction of these coherent structures with the free surface leads to large ripples that propagates down- stream. ACKNOWLEDGMENTS The work by A.I. was supported by the Ounce of Naval Research, under grant N.000140010344, through Dr. Pat Purtell. The work by E.F.C. was supported by the Min- istero Trasporti e Navigazione in the frame of the IN- SEAN research plan 2000-02. REFERENCES Brackbill, J.U., Kothe, D.B., and Zemach, C., "A con- tinuum method for modeling surface tension," Journal Computational Physics, Vol. 100, 1992, pp. 335-354. Brandt, A., "Guide to multigrid development," Multigrid Methods, Hackbush W., Trottenberg U. Eds., Springer- Verlag, Berlin, Germany, 1992. Ceniceros, H.D., and Hou, T.Y., "Dynamic generation of capillary waves,", Physics of Fluids, Vol. 11, 1999, pp. 1042-50. Duncan, J.H., "An experimental investigation of break- ing waves produced by a towed hydrofoil," Proceedings of Royal Society, London, Ser. A, Vol. 377, 1981, pp 331-348. Duncan, J.H., "The breaking and nonbreaking wave re- sistance of a two-dimensional hydrofoil," Journal of Flu- id Mechanics, Vol. 126, 1983, pp. 507-520. Duncan, J.H., Philomin, V., Behres, M., and Kimmel J., "The formation of spilling breaking water waves," Physics of Fluids, Vol.6, 1994, pp. 2558-2560. Duncan, J.H., and Dimas, A.A., "Surface ripples due to steady breaking waves," Journal of Fluid Mechanics, Vol. 329, 1996, pp. 309-339. Duncan, J.H., Qiao, H., Philomin, V., and Wenz, A., "Gentle spilling breakers: crest profile evolution," Journal of Fluid Mechanics, Vol. 379, 1999, pp. 191-222. Iafrati, A., and Campana, E.F., "Unsteady free surface flow around hydrofoils," Boundary Elements XX, Kassab A., Brebbia C.A., Chopra M. Eds., Computational Me- chanics Publications, Southampton, UK, 1998, pp. 451- 460. Iafrati, A., Olivieri, A., Pistani, F., and Campana, E.F., "Numerical and experimental study of the wave break- ing generated by a submerged hydrofoil," Proceedings of the 23r~ Symposium on Naval Hydrodynamics, Val de Reuil, France, 2000. Iafrati, A., Di Mascio, A., and Campana, E.F., "A level- set technique applied to unsteady free surface flows," International Journal for Numerical Methods in Fluids Vol. 35, 2001, pp. 281-297, 2001. Iafrati, A., and Campana, E.F., "A domain decomposi- tion approach to compute wave breaking," submitted for publication to the International Journal for Numerical Methods in Fluids, 2002. Kim, J., and Moin, P., "Application of a fractional-step method to incompressible Navier-Stokes equations," Jour- nalofComputationalPhysics, Vol. 59, 1985, pp. 308- - 323.
Lin, J.C., and Rockwell, D., "Evolution of a quasi-steady breaking wave," Journal of Fluid Mechanics, Vol. 302, 1995, pp. 29-44. Longuet-Higgins, M.S., "Capillary rollers and bores," Journal of Fluid Mechanics, Vol. 240, 1992, pp. 659- 679. Longuet-Higgins, M.S., "Shear instability in spilling break- ers," Proceedings of Royal Society, London, Ser. A, Vol. 446, 1994, pp. 399-409. Longuet-Higgins, M.S., "Progress towards understand- ing how waves break," Proceedings of the 21St Sympo- sium on Naval Hydrodynamics, Trondheim, Norway, 1996, pp. 7-28. Mui, R.C.Y., and Dommermuth, D.G., "The vertical struc- ture of parasitic capillary waves," Journal of Fluids Engi- neering, Vol. 117, 1995, pp. 355-361. Qiao, H., and Duncan, J.H., "Gentle spilling breakers: crest flow field," Journal of Fluid Mechanics, Vol. 439, 2001, pp. 57-85. Rai, M.M., and Moin, P., "Direct simulations of tur- bulent flow using finite-difference schemes," Journal of Computational Physics, Vol. 96, 1991, pp. 15-53. Rosenfeld, M., Kwak, D., and Vinokur, M., "A frac- tional step solution method for the unsteady Incompress- ible Navier-Stokes equations in generalized coordinate system," Journal of Computational Physics, Vol. 94,1991, pp. 102-137. Russo, G., and Smereka, P., "A remark on computing distance functions," Journal of Computational Physics, Vol. 163, 2000, pp. 51-67. Scardovelli, R., and Zaleski, S., "Direct numerical simu- lation of free surface and interracial flow," Annual Review of Fluid Mechanics, Vol. 31, 1999, pp. 567-603. Sethian, J.A., "Level Set Methods and Fast Marching Methods," Cambridge University Press, 1999. Sussman, M., and Dommermuth, D.G., "The numerical simulation of ship waves using Cartesian grid methods," Proceedings of the 23rd Symposium on Naval Hydro- dynamics, Val de Reuil, France, 2000. Sussman, M., and Fatemi, E., "An efficient, interface- preserving level-set redistance algorithm and its applica- tion to interracial incompressible fluid flow," SIAM Jour- nal Scientific Computation, Vol. 20~4), 1999, pp. 1165- 1191. Sussman, M., and Puckett, G., "A coupled level-set and volume-of-fluid method for computing 3D and axisym- metric incompressible two-phase flows," Journal of Com- putational Physics, Vol. 162, 2000, pp. 301-337. Sussman, M., Smereka, P., and Osher, S., "A level set approach for computing solutions to incompressible two- phase flow," Journal of Computational Physics, Vol. 114, 1994, pp. 146-159. Tulin, M.P., "Breaking of ocean waves and downshift- ing," Waves and Nonlinear Processes in Hydrodynamics, Grue J., Gjevik B., and Weber J.E. Eds., Kluwer Acad., Dordrecht, Netherland, 1996. pp. 170-190. Zang, Y., Street, R.L., and Koseff, J.R., "A non-staggered grid, fractional step method for time-dependent incom- pressible Navier-Stokes equations in curvilinear coor- dinates," Journal of Computational Physics, Vol. 114, 1994,pp.18-33.
DISCUSSION D. G. Dommermuth Science Applications International Corporation, USA The authors have made a significant contribution toward developing a numerical procedure for simulating and understanding gravity-capillary waves. For the matching procedure, a Helmholtz decomposition, whereby the upper domain is divided into its potential and vertical components, may give additional insight into the boundary condition that is imposed by the authors. Also, could the authors please comment on what is the next step toward improving our understanding of small-scale breaking waves? AUTHORS' REPLY We thank Dr. Dommermuth for his comments on the paper. Conceming with the use of a Helmholtz decomposition, we actually recognize that it could allow a better understanding of the mechanisms governing the exchange of information between the two sub-domains. Likely, a better way to enforce the matching condition could be also derived. As to the second question the authors are thinking of taking advantage from the two-phase formulation used to solve this problem and move toward the analysis of the small-scale wave-wind interaction, being a very important generation mechanisms for short surface tension dominated breaking waves. This point has been studied experimentally, for instance by Tulin (1996), but little has been done with numerical tools, which can help in the understanding of the fundamental mechanism of this phenomena.
