UFDC Home  myUFDC Home  Help 



Full Text  
PAGE 1 A NUME RICAL MODEL BASED EXAMINAT ION OF CONDITIONS FOR IGNITIVE TURBIDITY CURRENTS By GOWTHAM KRISHNA A THESIS PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLOR IDA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF MASTER OF SCIENCE UNIVERSITY OF FLORIDA 2009 1 PAGE 2 2009 Gowtham Krishna 2 PAGE 3 To m y Family 3 PAGE 4 ACKNOWL EDEGMENT First, I would like to thank my advisor Dr. Ashish J. Mehta, for all the guidance, support and knowledge that I received from him. Without him, this thesis would not have been completed. Special thanks go to other memb ers of my committee including Dr. James F. Klausner and Dr. Robert G. Dean. Thanks to my colleagues and office mates: Emre Ozdemir, who has been more than a friend to me and Sangdon So, for all the fun, humor and good times we have had. I would also like to thank all my colleagues from the departme nt who have made my two years an enjoyable experience. Lastly, but not the least I would like to thank my family for all the love and support they have given me. 4 PAGE 5 TABLE OF CONTENTS page ACKNOWLEDEGMENT...............................................................................................................4 LIST OF TABLES ................................................................................................................. ..........7 LIST OF FIGURES.........................................................................................................................8 ABSTRACT...................................................................................................................................10 CHAP TER 1 INTRODUCTION..................................................................................................................12 1.1 Study Motivation........................................................................................................... ...12 1.2 Objective and Scope.........................................................................................................14 1.3 Thesis Structure........................................................................................................... .....15 2 TURBIDITY CURRENT:BACKGROUND..........................................................................17 2.1 Definition of Turbidity Current........................................................................................17 2.2 Early Investigations ..........................................................................................................17 2.3 The Structure of Gravity Current......................................................................................19 2.4 KnappBagnold Criterion for Suspension........................................................................23 2.5 Concluding Observation...................................................................................................28 3 TWOPHASE EQUATIONS AND CLOSURES..................................................................31 3.1 Introduction................................................................................................................... ....31 3.2 TwoPhase Flow Equations..............................................................................................32 3.3 Closure of Fluid Stresses..................................................................................................34 3.4 Closure of Sediment Stresses............................................................................................38 3.4.1 Sediment Transport Equation.................................................................................38 3.4.2 Collisional Theory of Granular Flow.....................................................................40 3.4.3 Large Scale Sediment Stress..................................................................................42 3.4.4 Closure at Low Concentrations..............................................................................43 3.4.5 Region of Enduring Contact...................................................................................44 3.5 Initial and Boundary Conditions.......................................................................................45 3.5.1 Initial Condition......................................................................................................45 3.5.2 Bottom Boundary Condition..................................................................................46 3.5.3 Top Boundary Condition........................................................................................47 3.5.4 Lateral Boundary Condition...................................................................................49 5 PAGE 6 4 NUMERICAL IMPLEMENTATION....................................................................................50 4.1 Introduction................................................................................................................... ....50 4.2 Mesh and Grid System......................................................................................................50 4.3 Numerical Discretization..................................................................................................51 4.3.1 Solution of Sediment Phase Equations ..................................................................51 4.3.1.1 Predictor Step...............................................................................................52 4.3.1.2 Corrector Step..............................................................................................58 4.3.2 Fluid Phase Equations ...........................................................................................59 4.4 Stability Analysis............................................................................................................. .61 5 RESULTS AND ANALYSIS.................................................................................................62 5.1 Preamble....................................................................................................................... ....62 5.2 Experimental Data and Model Output..............................................................................64 5.3 Domain of Ignitive Behaviour..........................................................................................65 6 SUMMARY, CONCLU SIONS AND RECOMMENDATIONS ...........................................76 6.1 Summary...........................................................................................................................76 6.2 Conclusions.......................................................................................................................77 6.3 Recommendations for Further Work ...............................................................................79 APPENDIX A Numerical Implementation Scheme.......................................................................................80 LIST OF REFERENCES...............................................................................................................85 BIOGRAPHICAL SKETCH.........................................................................................................88 6 PAGE 7 LIST OF TABLES Table page 31 Summary of numeri cal coefficients adopted by the m odel...............................................38 51 Model runs for gr ain size of 0.21mm.................................................................................68 52 Model runs for gr ain size of 0.28mm.................................................................................68 53 Model runs for gr ain size of 0.35mm.................................................................................69 54 Flow state in experim ent and model run ...........................................................................69 55 Flow states for 0.21mm particles ......................................................................................69 56 Flow states for 0.28mm particles ......................................................................................70 57 Flow states for 0.35mm particles ......................................................................................70 7 PAGE 8 LIST OF FI GURES Figure page 11 Diagrammatic definition of the continental shelf along with the continental slope and related geographic features. ...............................................................................................16 21 Shadow pictures of the profil es of the head of the gravity current.....................................29 22 Frontal region of a gravity current of saltwater advancing along the floor of a freshwater tank.. .............................................................................................................. ...29 23 Inviscid gravity current with m ixing. The lower region (depth h4) is unmixed saline solution, and the top region ( h2) is unmixed ambient fluid. Mixing takes place in the middle region of thickness h3.............................................................................................29 24 Gravitational force components and velo cities associated with th e f luid and the solids....30 25 Schematic drawing of the experim ental setup of Parker et al............................................30 51 Physical interpretation of mode simula tions: (a) downstream gradient in pressure force leads to the development of a boundary layer velocity profile and associated profile in the suspended sediment concen tration. Under nonignitive conditions the concentration decreases due to sediment deposition; (b) conditions are conducive for ignition...............................................................................................................................71 52 Growth of boundary layer velocity and c oncentration profiles dur ing a pretest run with d = 0.21mm u* = 0.06 m/s and = 0o.......................................................................71 53 Velocity and concentration profiles: Experim ental data of Parker et al. (1987) and present model run. The main difference is in the grain size (0.03mm in experiment versus 0.20mm in the model). Both are at a distance of x =1.5m.......................................72 54 Velocity and concentration profiles: Experim ental data of Parker et al. (1987) and present model run. The main difference is in the grain size (0.03mm in experiment versus 0.20mm in the model). Both are at a distance of x =4.5m.......................................72 55 Velocity and concentration profiles: Experim ental data of Parker et al. (1987) and present model run. The main difference is in the grain size (0.03mm in experiment versus 0.20mm in the model). Both are at a distance of x=8.5m.......................................73 56 Shields parameter against bed slope..................................................................................73 57 Shields parameter against bed slope Sedim ent flux per unit width (kg/m3/s) contours for 0.21mm diameter particles...........................................................................................74 8 PAGE 9 58 Shields parameter against bed slope Sediment flux per unit width (kg/m3/s) contours for 0.28mm diameter particles...........................................................................................74 59 Shields parameter against bed slope. Sedi m ent flux per unit widt h (kg/m3/s) contours for 0.35mm diameter particles...........................................................................................75 510 Shields parameter against particle diam eter. Threshold curves for incipient movement and autosuspension curve. Model test run data..............................................75 9 PAGE 10 Abstract of Thesis Presen ted to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree of Master of Science A NUMERICAL MODEL BASED EXAMINAT ION OF CONDITIONS FOR IGNITIVE TURBIDITY CURRENTS By Gowtham Krishna December 2009 Chair: Ashish J. Mehta Major: Coastal and Oceanographic Engineering Turbidity currents form a major mechanism for the movement of sediment in the natural environment. Selfaccelerating turbidity currents over continental slopes are of considerable scientific and engineering inte rest due to their role as agents for submarine sediment transportation from the shelf to the seabed. Su ch currents are called ignitive provided they eventually reach a catastrophic st ate as acceleration results in high sediment loads due to erosion of the sloping bed. In the present study ignition refers to the onset of selfacceleration and erosion characterized by positive values of the rate of change of current velocity and sediment flux down the slope. A numerical model, which treats the fluid and the particles as two separate phases, is applied to investigate the effects of particle size, initial flow friction velocity and mild bed slope on the ignitive condition. Laboratory expe rimental data of Parker et al. (1987) have been included as part of the analysis. Ignition for the smallest of the three select ed sizes (0.21mm) of medium sand typical of Florida beaches was found to depend on the initia l conditions at the head of the slope as determined by the pressure gradient. Bed slope s eemed to be of secondary importance. For the two sands with larger grain si zes (0.28mm and 0.35mm) the bed slope was found to play a more 10 PAGE 11 im portant role when compared to the initial pr essure gradient. For a given pressure gradient, increasing the slope increased the likelihood of se lfacceleration. Thus, in general ignition cannot be defined merely in terms of nontrivial, pos itive values of the velocity gradient and the sediment flux gradient along the slope. Depending on particle size the initial pressure gradient can also play a role. For the selected initial co nditions (grain size, pressure gradient and bed slope), out of the 54 combinati ons tested, all except three satis fied the KnappBagnold criterion for autosuspension irrespective of whether the turbidity current was ignitive or nonignitive. In all 54 cases the current was found to erode th e bed. Further use of the model will require accommodation of wider ranges of sediment size and bed density, and a thorough verification against experimental data. 11 PAGE 12 CHAP TER 1 INTRODUCTION 1.1 Study Motivation A flow resulting from density difference due to spatial variability in the suspended sediment concentration is known as a turbidity current. In nature turbidity currents are witnessed in a variety of environments. The major difference between density currents due to dissolved substances such as salt or due to the effect of temperature on density on one hand, and turbidity currents on the other is the poten tial ability of the turbidity cu rrent to deposit sediment on the bed. Unlike the movement of bedload or suspended load of sediment due to water flow, a turbid current can also move through otherwise st ill water as a diffusioninduced front. Turbidity currents form a major mechanism for the transport of continental shelf sediment into deeper waters. They commonly occur on the shelf due to inflow of turbid waters from riverine sources (Fig. 1.1). Clos er to the river mouth, they ma y persist as underflows over the foreset slope of the clinoform delta. By way of transport of sediment from the shelf to the continental slope, and then into the deeper ocea nic waters, turbidity curre nts are believed to be responsible for breaks in submarine pipelines and telephone cables. The current can erode a submarine channel analogous to a river, and sediment thus transported will contribute to the development of the seabed. Selfaccelerating turbidity currents is a topic of considerable scient ific and engineering interest, because such currents can be responsible for significant trans port of sediment from terrestrial sources to the seab ed. Selfacceleration can occur if the flow velocity and the suspended sediment concentration increase simu ltaneously as bed sediment erodes due to the impelling downslope gravity force and associated turbulence. 12 PAGE 13 Selfacceleration in the present study is defined as the condi tion when the turbid current generates sufficient turbul ence to hold the particles in suspension. It causes the bed to erode, particles are picked up and the suspension becomes de nser. As a result, accelera tion increases, more sediment is picked up and the suspension beco mes even denser. Ignition in this case can be considered to refer to the onset of selfaccele ration without reference to its subsequent fate. Selfacceleration necessarily included a conditi on in which erosion occurs and its rate is greater than that rate of depos ition. Knapp (1938) referred to the process of ignition as autosuspension, a term he is believed to have coin ed. However, the subsequently derived criterion for autosuspension by Bagnold (1962) is based on the condition of equilibrium between eroding bed sediment and depositing suspended (as opposed to bedload) sediment. Thus the criterion for autosuspension is a necessary but not a suffici ent condition for selfacceleration of sediment suspended in a current (Parker, 1982). A restraining influence on an ignitive turbidity current is that the selfsustaining property of the current ceases to exist ei ther due to nonavailability of sediment or due to damping of turbulence cau sed by the excess suspended sediment. When the flow becomes saturated with sediment turbul ence collapse can occur (Winterwerp, 1999). The present study is concerned only with conditions th at lead to ignition. Turbulence collapse is a catastrophic condition that is not considered in the present analysis. Because of the submersed locations and scales of natu ral turbidity currents in the sea, it has been difficult to collect data on selfacceleration in the field. As a result, although selfsustaining currents have been postulated to occur on continental slopes, espe cially close to the seafloor, little is known about their incipient, ignitive dynamics, which is complicated by strong interactions between the sedi ment particles and water flow Therefore, to understand the mechanism of selfacceleration in detail, mathematical models have been developed to simulate 13 PAGE 14 the onset of turbid ity currents in the natural environment. Alt hough such models do not fulfill the need for field data, they are useful for develo ping broad criteria meant to identify conditions under which selfacceleration or ignition occurs. However, most modeling has been based on the use of a single fluid assumed to represent a su spension of sedimentary particles. The modeling approach selected for the present study considers the fluid (water) and sedi ment as two separate phases. 1.2 Objective and Scope In a vertical, onedimensional (1DV) model, tw ophase flow equations for the fluid and the sediment are used to examine the effects of impor tant factors resulting in selfacceleration of turbidity currents. The analysis has been carried out within the framework of the wellknown KnappBagnold criterion for autosuspension of sediment particles in turbulent flows. The mathematical model has been stated wi th the basic equations of continuity and momentum for the two phases, as well as the cl osure schemes and approximations considered. Since the 1DV model is restricted to the coordinate perpendicular to the be d, a vertical column is considered in which essential fluid and se diment properties and dynamic behaviors are described. Initially, a pressure gradient is used as a forcing condition to achieve a fully (or nearly fully) developed boundary layer flow without taki ng the bed slope into account, as the gravity force acts in the directio n perpendicular to the bed. To simula te the subsequent flow conditions over a slope (with a prescribed be d sediment concentration), which is assumed to be mild, the pressure force is switched off and the effects of the gravity components in the pseudostreamwise and vertical directions are invoked. Values of th e velocity and sediment concentration are noted at different elevations and times Based on the depthmean current velocity, the vertical profiles of velocity and concentration can be translated by calculating the approximate distance traveled during each (selected) elapsed time. The sediment flux is then calculated at selected distances 14 PAGE 15 from the start of the slope corres ponding to the times of velocity and concentration sampling in the model. The model is run (numerically) for a range of initi al conditions with respect to particle size, the pressure gradient and bed sl ope in order to identif y the domain within which the system is ignitive, i.e. for which selfacceleration occurs. 1.3 Thesis Structure This thesis is organized into five chapters. Chapter 2 gives a brief introduction to turbidity currents. The main factors affecting gravity currents in general and the place of turbidity currents within the classification of gravity currents are mentioned. Also, a review of studies relevant to the analysis of the model output (in Chapter 5) is presented. Chapter 3 mentions two numerical approaches employed in computing multiphase flows, including the EulerEuler approach and the EulerEu ler approach (used in the present case). Then the governing equations of twopha se flow and their respective closure schemes are described. Chapter 4 presents the numerical method app lied to solve the governing equations along with the initial and boundary conditions us ing the finite difference approximation. In Chapter 5 the results from model applica tions are presented and discussed. Selective data sets from studies mentioned in Chapter 2 ha ve been used to corroborate the findings related to the parametric domain of ignitive turbidity currents within the framework of the KnappBagnold criterion (also desc ribed in Chapter 3). Conclusions based on this analysis and suggestions for further work are given in Chapter 6. In Appendix A, the numerical implementation of the code is given as a flow chart. 15 PAGE 16 Coast Cont inental shelf Continental rise Ocean Figure 1.1 Diagramm atic definition of the continen tal shelf along with the continental slope and related geographic features. 16 PAGE 17 CHAP TER 2 TURBIDITY CURRENTS: BACKGROUND 2.1 Definition of Turbidity Current Information related to th e turbidity current in this and the next two sections (2.2 and 2.3) is meant to provide the reader a ge neral background on turb idity currents. The description of autosuspension in Section 2.4 is di rectly related to the analys is presented in Chapter 5. Turbidity currents are flows of sedimentladen fluid down a slope or along a horizontal surface through water or another fluid. Under cer tain circumstances th e current may not be merely selfsustaining but also selfreinforcing, using the bed sediment as fuel to reach a catastrophic condition. Turbidity currents can be characterized by a di stinctive raised front or head followed by a quasiuniform flow region. Expe riments by Schmidt (1911) showed (Fig. 2.1) the way by which density difference corresponding to a temperature difference influence the shape of the head. The temperature is raised form a few degrees in (a) to 9C in (f) corresponding to a density difference of 1%. The raised nose and the instabilities at th e top surface can be seen from (a) to (f). Figure 2.2 shows the front of a gravity current of salt wa ter advancing in freshwater. The head can be a zone of intense fluid mixing and interfacial wave breaking. Billows on the upper surface of the head may grow and break up into a threedime nsional pattern leaving behind a mixed layer. Another distinct feature of the tu rbidity current is the exchange of the granular particles between the current and the bed. 2.2 Early Investigations An earthquake at the Grand Banks of the s outh coast of Newfoundland (Canada) triggered a turbidity current at estimate d speeds of 60100 km/h (Fine et al., 2005). The earthquake created a large submarine landslide that led to the breakage of 12 submarine telegraph cables and 17 PAGE 18 resulted in a tsunam i. Heezen and Ewing (1952) and Heezen et al. (1954) provided evidence that the break in the telegraph cables was caused by a selfaccelerating turbidity current. As an explanation for the origin of subm arine canyons, Daly (1936) proposed that they were possibly cut by turbidity currents. Kuenen (1937) carried out suppor tive experiments in the laboratory to test Dalys hypothe sis. Theoretical and experime ntal investigations by Knapp (1938) are also among the early works of turbidity and conservative gravity currents. Later, Middleton (1966a, 1966c) carried out a series of laboratory experiments to study the frontal head. He also described laws concerning the uniform flow of density currents and the deposition of sediment by turbidity currents. Plapp and Mitc hell (1960) provided one of the first models for turbidity currents. Most of the early models considered tu rbidity currents to be analogous to conservative gravity flows. Chu et al. (1979) presented an anal ytical model specifically for turbidity current dynamics. They derived the gover ning equations of motions for the continuity of mass and momentum for the solid and the liqui d phases and a diffusion equation for the rate of change of the suspended sediment within the co ntrol volume. Parker (1982) examined turbidity currents applying a slab model and investigated the state of ignition at which the flow entrains sediment and accelerates the current to a catastrophic state. McTigue (1981) used twophase flow equations to study sedi ment transport due to turbidity currents, while Greimann et al. (1999) applied twophase flow equati ons to turbid channel flows. Lin and Mehta (1997) studied the depositiona l turbidity currents and the accompanying sediment transport that takes place in long and shallow closedend ba sins. The experiment simulated a lock exchange setup where the ambient fluid consisted of fresh water and in different test runs sediment laden dense fl uid consisted either of kaolinite, flyash or mud. They show that 18 PAGE 19 along the length of the basin th e se ttling velocity decays exponentially with distance and the rate of accumulation of sediments is highest near the entrance and decreases with distance. 2.3 Structure of Gravity Currents In gravity currents moving on a horizontal surf ace the head remains in a quasisteady state unlike in flow over a slope where the head increases in size with increasing angle of bed inclination. The front advancing over a rigid bed will have its foremost point elevated as a result of the noslip condition acting at the lower boundary. This feature is the nose. The shape of the head/nose does not conform to any unique outline as it changes with the depth of the ambient fluid. Other factors that affect the head shape are the physical pr operties of the two fluids and bed roughness. Flow turbulence and current direction, whether opposing or flowing into the ambient flow, also play a role in determining head shape. Turbulenceinduced mixing is a important featur e of the head of gravity current. There are two main kinds of in stabilities which cause turbulence mixing: Billows: Due to the wellknown KelvinHel mholtz instability, billows are formed on the upper surface of the head and roll up in the region of ve locity shear above the front. Lobes and clefts: Due to the elevated frontal nose, as the current advances some of the lighter fluid in its path is disturbed resulting in instability. Most attempts at understanding the head of the gravity current have been from the standpoint of a conservative current. However, a comparison of the heads of turbidity current with saline gravity current was done by Middleton (1966a). The basic charac teristics influencing the head on a horizontal bed have been discussed in a general se nse and these have also been applied to the heads of a gravity current including a turbidity current on slopes. 19 PAGE 20 The inviscid fluid theory has b een applied to study features of steady grav ity currents such as wave breaking. The propagation of an immiscib le gravity current over a frictionless surface was first proposed by Von Karman (1940). Even in the absence of frictional forces, instabilities leading to the formation of billows can be shown to occur. In the absence of friction the tip of the gravity current would touch the boundary indi cating the absence of a nose. Von Karman considered that the interface would become paralle l to the surface at a large distance behind the head and there the velocity would be constant. Von Karman applied the Bernoulli equation betw een the tip of the current (which in this case touches the ground) where the stagnation po int exists and along the interface far upstream and obtained the relation f fHgU 2 (2.1) Where, is the head height, is the reduced gravity and is the propagation velocity. fHg fUBenjam in (1968) investigated the front of the fr ictionless gravity current in terms of a cavity current, i.e. an empty cavity advancing in a liquid displacing the fluid beneath. He calculated from the continuity equation and from the Bernoulli eq uation applied at the in terface. For frictionless flow the pressure force and the momentum flux per unit mass are conserved. The result is that if the current height h is assumed equal to it will also be equal to onehalf the total water depth H o f the ambient fluid, i.e. h=H/ 2. This implies that the gravity current can progress steadily if it fills one half the height of the ambient fluid. If h>H /2, then external supply of energy would be necessary to sustain the fl ow and this is not possible; whereas if h PAGE 21 Britter and Simpson (1978) carried out a semie mpirical analysis in which the surface is assumed to be smooth but mixing is allowed to take place between the gravity current and the ambient fluid, as for example in the case of largescale gravity currents of saltwater displacing freshwater at tidal or river lo cks, where frontal mixing occurs. The apparatus consisted of a flexible conveyor belt which could move at the speed of the current effectively giving rise to the slip condit ion, thereby simulating a frictionless surface. The experiment represented a gravity current of a de nse fluid (saline solution) miscible with the ambient fluid (freshwater) and moving along the su rface effectively in th e absence of friction. The difference when compared to Von Karman (1940) was that the fluid in Von Karmans analysis was immiscible unlike the present case. In the Britter and Simpson investigation the head was found to be divided into three region s. The bottom region represented the flow of unmixed dense fluid into the gravity current. Th e topmost region consisted of only the unmixed, less dense the ambient fluid. Between these tw o regions was the layer in which mixing took place. This region had nonuniform velocity and salinity profiles. The Britter and Simpson apparatus also made it possible to investigate the properties of billows at the front. It was found that they had the properties of KevinHelmholtz instabilities. Such instabilities could occur for a certain ranges of values of the relative velocity U between the two fluids, distance traveled and the density difference in terms of reduced gravity. Based on the Richardson number Ri = it was observed that the in stab ilities occurred for Ri <1/4. Uhg /'In a flow with nos lip condition (due to th e friction) at the lowe r boundary, the less dense fluid accumulated underneath the nose ascends as the front moves forward, and the associated instabilities occur as unsteady lobes and clefts. B illows still remain as the main feature within which the ambient fluid is mixed with th e gravity current (Bri tter and Simpson, 1979). 21 PAGE 22 Experim ents (Simpson, 1972) have shown that it is feasible to suppress the unstable lobes and clefts by overrunning the less dense fluid with a layer of denser fluid. Si mpson and Britter (1980) conducted laboratory experiments describing the nose height in relation to the head of the current as a function of the Reynolds number, 134Re()/ Uhh where is the speed of the opposing flow, is the height from the floor to the top of the mixed layer and is the kinem atic viscosity of water (Fig 2.3). 1U43hh As m entioned, the turbidity current is a specia l case of a large class of gravity flows. The following classification of gravity flows (Altinak ar 1988) is based on the density difference of the entering fluid and the state of stratification of the ambient fluid. 1. If the entering fluid is denser than the am bient fluid then the current plunges to the bottom and is known as an undercurrent or a bottom gravity current. 2. If the entering fluid is lighter than the ambient fluid, the current will occur on top of the ambient fluid, also called an overflow 3. In case the ambient fluid is a stably st ratified, different scenarios are possible depending on the density of the entering fluid relative to the ambient fluid. For example, if the density of the entering fl uid has an intermediate value between the lowest and the highest densities of the st ratified ambient fluid, an intermediate or intrusion type of current is formed at the appropriate density level. When the gravity current is due to temperatur e or dissolved matter, the buoyancy flux, i.e. the flux of excess body force through a s ection of a current, is defined as (/)a B gq where is the density difference between the fluid comprising the current and the ambient fluid, is the density of the am bient fluid, g is the acceleration due to gravity and is the discharge per unit width. Note that the reduced gravityaq /.agg The buoyancy flux is conserved 22 PAGE 23 throughout the current. In contrast, in a turbid cu rrent this m ay not be the case as the current may erode the bed sediment or permit deposition of the suspended particles. 2.4 KnappBagnold Criterion for Suspension Knapp (1938) was one of the earliest investigat ors to explore the basic criterion for the transport of granular materials in suspension. He hypothesized the concep t of autosuspension stating that sediment partic les in a flowing stream add or subtract energy through two independent processes. If the amount of en ergy added per unit time is less than the amount subtracted per unit time then the net result woul d be that the particles would become a burden on the stream, further implying that the there must a definite limit to the number of such particles that can be kept in suspension. On the other han d, if the amount of energy added is greater than the amount subtracted, then the energy of the stream is increased by the presence of the particle. This would indicate that there w ould be no limit to the number of particles kept in suspension and transported by (a unit volume of) the fluid. The resulting flow is an autosuspension. Bagnold (1962) developed Knapps arguments and suggested a criterion by examining suspended sediment trans port in a laboratory setup. The physical condition considered was a ge ntly sloping bed which enabled a turbid suspension as well as sedimenta tion of sand and silt to occur, w ith turbulence as the cause of resuspension (Fig 2.4). The mass m of the suspended sediment steadily falls but the center of gravity of the suspended sediment stays above the bed, becaus e turbulence acts in th e opposing direction in keeping the sediment suspended .The force du e to turbulence should be acting against the immersed weight of the mass with fall velocity .The power required to do so given by sw 23 PAGE 24 s s fsgmw (2.2) Where is the sediment density and the fluid density. sfThe power expended by the fluid is replenished by the tangential pull acting on the suspension by the downslope gravity com ponent. This is given by sin Ugms fs (2.3) Where Uis the speed of the solids traveling downstream, and is the angle of bed inclination with the horizontal. Therefore the net power expended by the fluid in keeping the sediment suspended is given by )sin/( UwUgms s fs (2.4) When the sediment size is decreased this term becomes zero independent of the magnitude of the mass. Therefore, sin Uws (2.5) As the bed slope increases the power provided by the tangential gravity component is sufficient not only to keep the sediment suspende d but also to contribut e to the power needed against fluid drag at the boundary. The total pow er due to available sediment is given by sin)( UhCgfs (2.6) Where, C is the mean sediment concentration and is the flow height. To maintain the suspension the power expended is h hvCgfs)( (2.7) 24 PAGE 25 Therefore, the excess po wer needed fo r flow against drag is given as )sin()(s fswUhCg (2.8) The power needed to maintain mean velocity uis u0, where is the bed shear stress. From Francis (1957) 0 u0 is given as k h u uo2.13 log 332 10 3 (2.9) From this development Bagnold concluded that the criterion for a sel fsustained turbidity current is given by, k h u u w UhCgs fs2.13 log 33 sin)(2 10 3 (2.10) Middleton (1966a) stated that this theory contains a number of doubtful and overly simplified assumptions. He maintain ed that Bagnold ignored the resi stance at the upper interface, sediment sorting and concentration effects. The form of Bagnolds criterion from equation (2.5) suggests that for larger ratios of settling veloc ity to the horizontal velocity, autosuspension is possible for larger slopes. Middleton contended th at an increase in the slope would increase the Froude number due to which mixing and resistance would rapidly increase at the upper interface. This in turn would indicate th at the optimum conditions for obtai ning autosuspension may be at low slopes. Middleton also questioned the inequality where U is the mean fluid velocity, which should be an equality instead. Middleton did not provide an alternative autosuspension criterion; however, he stated that the ideal conditions in which autosuspension would be favorable would be a combination of the presence of fine grained sediment, low sediment density and a relatively large physical model scale. 1/ swUs 25 PAGE 26 Pantin (1979) argued that the energy budge t considered by Bagnold was incom plete and incorrect. In Bagnolds hypothesis the total availa ble power from equation (2.6) may be denoted by the symbol w The energy needed to keep the sediment suspended is given by equation (2.7). Let this energy be denoted by the symbol and let the work expended against bottom friction, nw 0, u be denoted as The essential power condition given by Bagnold in his autosuspension criterion can be summarized as tw (2.11) tnwww However, Pantin contended that this conditio n is not valid because the power expended in supporting suspended sediment wo uld be independent of the power expended against bottom friction. He noted that the pow er needed to support sediment in suspension would be derived from turbulence within the flui d and the necessary power require d to support turbulence itself derived from the work expended against bottom friction. Therefore is part of and the rem ainder of is dissipated as microturbulence and ev entually heat, which has little influence in keep ing the sediment suspended. Therefore Pantin redefined the power involved by the following relationships nwtwtw (2.12) andtxtwweww n Where is an efficiency factor. A large efficiency would promote autosuspension. xexeParker (198 2) studied the conditions under which ignition of catastrophically erosive turbidity currents is pos sible. He mentioned th at the Bagnold criterion is a necessary but not a sufficient condition for ignition. He presented a threeequation m odel explaining the process involved in ignition. Parker et al. ( 1986) proposed a more refi ned fourequation model 1/ swUs 26 PAGE 27 by adding the equation of m ean turbulence energy to the conservation equations for fluid mass and sediment as well as conservation of moment um of the mixture, and solved the equations numerically to investigate the possibility of se lfaccelerating turbidity current. Parker et al. (1986) observed that the earlier threeequation model does not account for the turbulent energy balance and the solutions thus obt ained are physically not possible. Following the development of the fourequati on model, Parker et al. (1987) conducted laboratory tests to determine the behavior of tu rbidity currents. The experiments consisted of currents laden with noncohesive sediment (silic a flour) on a sloping bed in a channel which was covered with the same sediment (Fig 2.5). Th e velocity and concentration profiles were measured at three locations along the channel. The results obtained were used to calculate tophat functions used in the mathematical model of Parker et al (1986). Let the local mean velocity, volumetric sediment concentration and the mean kinetic energy per unit mass be denoted by u, c and k the layer thickness h and the corresponding layer averages U ,C and K Then the following similarity laws were assumed (,)(,) ();; ()()()uckxz uxzcxz UxCxKx k (2.13) Where h z (2.14) The tophat functions ar e defined such that (2.15) 10 101 if ifkcu 27 PAGE 28 The curren t was considered to be ignitive when th e rate of change of current over distance was positive and the the rate of change of sediment flux over distance was also positive. Based on this definition, two of the 23 runs were show n to generate an ignitive turbidity current. Unfortunately, the ensuing catastr ophic condition could not be r ecorded as the experimental facility was too short. Eidsvik and Brors (1989) studie d selfacceleration of turbidity currents using an algebraic closure model and also the model; they presented results usi ng the latter clos ure model only. They assumed the sediment volume concentration to be very small and predicted the general structure of an ignitive turbidity current by show ing an increase in the velocity, concentration and the turbulent kinetic energy profiles over an extended peri od of time. They simulated conditions applicable to fine grain sediment a nd small slopes. They predicted that ignition can occur at very small angles. 2.5 Concluding Observation Most mathematical models have treated the tu rbidity current as a singlephase conservative current, which makes the analysis approximate and therefore somewhat qualitative. Therefore an attempt has been made in the present study to look in to the phenomenon of ignition through the use of a twophase model in which the role of se diment is included explicitly. Also, the inclusion of the turbulent suspension flux term and particle stresses due to interparticle interactions does not require the imposition of a sediment pick up function, or equations for bedload and suspended load transport. 28 PAGE 29 Figure 2.1 Shadow pictures of th e profiles of the head of the gravity current (Source: Simpson, 1987). Figure 2.2 Frontal region of a gr avity current of saltwater a dvancing along the floor of a freshwater tank (Sour ce: Simpson, 1987). h3 h4 h1 U4 h2 U1 1 2 Figure 2.3 Inviscid gravity current w ith m ixing. The lower region (depth h4) is unmixed saline solution, and the top region ( h2) is unmixed ambient fluid. Mi xing takes place in the middle region of thickness h3 (Adapted from Simpson, 1987). 29 PAGE 30 v sin = sli p v cos vsin Figure 2.4 Gravitational force components and velocitie s associated with the fluid and the solids (Source: Bagnold, 1962). Figure 2.5 S chematic drawing of the experi mental setup of Parker et al. (1987). C(x z ) u(x z ) x z U0, C0 Water Sediment Agitator Sediment m Fluid velo city Solid velocity U uUv sin gm cos gm 30 PAGE 31 CHAP TER 3 TWOPHASE EQUATIONS AND CLOSURES 3.1 Introduction A large number of flow regimes in nature can be interpreted as mixtures of two phases. Different sizes of the same material can be dealt with as separate phases. In multiphase flow a phase can be identified from its interaction wi th the flow fluid. For example gassolid flows (particle laden flows), gasliq uid flows (bubbles in continuous fluid), liquidliquid flows (immiscible fluids separated by a discernible in terphase), liquidsolid flows (slurry flow) and threephase flows (a combination of the above) can be said to be the major classification regimes of a multiphase flow. The numerical approach of com puting multiphase flows can be broadly classified into two categories: EulerLagrange and EulerEuler. EulerLagrange Approach: The fluid is treated as a continuum and the time averaged Navier Stokes equations are used to solve them. Dispersed particles are treated by tracking a large number of them through the flow field. A basic assumption made in the present study is that the second phase occupies a low volumeconcentra tion and has no direct impact either on the generation or dissipation of turbulence. EulerEuler Approach: In this approach both the phases are treated as interpenetrating continua. The derived conservation equations for the phases are consistent with each other and have similar structures. These equations are coupled through pressure or th e interphase exchange coefficient. This approach is computationally more expensive as it evaluates more equations, but is also more accurate. This approach was therefore used. 31 PAGE 32 Sedim ent transport in water is a twophase phenomenon. For the most exhaustive treatment of sediment transport one can assume fluid flow to be governed by the NavierStokes equations and the action of every sediment particle can be considered through its corresponding equation of motion. However, as this is not pragmatic, the sedi ment particles will be considered as part of a continuum. Accordingly, the equations of mo tion for the two phases have the same form. Ensemble averages of the mass and momentum equations are considered As part of this averaging the definition of sediment concentration c is introduced. 3.2 TwoPhase Flow Equations Ensemble averaging is carried out over a size of the order of the grain diameter, but the sediment concentration fluctuates on a scale much larger than the grain diameter, hence another averaging process is required (Hsu et al., 2003). This process, i.e. Favre averaging, is carried out on the already averaged twophase flow equations. The Favre aver aged fluidphase continuity equation is 0)1()1( f f fwc z c t (3.1) and the sedimentphase continuity equation is 0 ssswc z c t (3.2) Where the mass density of water. Sediment particle s are ass umed to have the same diameter d with mass density z is the norm al axis component to the bottom, andfs fw, sw are the zdirection averaged fluid and particle velocities, respectively. Let us consider the fluidphase momentum equations in the x direction: 32 PAGE 33 sf f xz f f fff ffuucSgc z x P cwuc z uc t )1( )1()1( )1( (3.3) and in the zdirection: z c wwcgc z z P cwwc z wc tft sf f zz f f fff ff )1( )1()1( )1( (3.4) Where fu and su are the averaged fluid and sediment horizontal velociti es, respectively, S= slope, f P is the ensemble averaged fluid pressure, and are the average fluidphase stresses including, in both cases the viscous stress and the fl uid phase Reynolds stress. The rem aining terms in the above equations, apart from the first term on the left hand side, which is the unsteady term, and the first term on the right hand side, which is the convective term, the remaining terms indicate Favre averaged drag for ces due to fluidsediment interaction and the relative mean velocity between the two phases, andf xzf zz is the drag coefficient given by n p r fc d U )1( 1 3.0 Re 0.18 (3.5) Where the concentration dependence is taken from the work of Richardson and Zaki (1954) on hindered settling and is the relative velocity between th e fluid and sedim ent phase velocities given as, rU 2 2 sfsf rwwuuU (3.6) If the dynamic viscosity is denoted f the particle Reynolds number is given by 33 PAGE 34 f r f pdU Re (3.7) The exponent n in (3.5), which depends on the partic le Reynolds number, is given by pn Re 45.4 (3.8) The corresponding sedimentphase momentum equations are given as below, in the xdirection as sf s xz sf sss ssuucSgc z x P cwuc z uc t (3.9) and in the zdirection as, z c uucgc z x P cwwc z wc tft sf s zz s f sss ss (3.10) The last term in equations (3.4) and (3.10) is the correlation between concentration fluctuations and the fluid velocity fluctuations. It represents the sediment flux caused by largescale fluid turbulence. The parameters and are the average stresses in the sediment phase. Each stress include s contributions from interparticle coll isions and sedimentfluid interaction and the Reynolds stress due to the velocity fluctuations of the particles. s xzs zz3.3 Closure of Fluid Stresses The tota l stresses in the flui d phase can be written as 0;0 f ff fff x zxzxz zzzzz z R R (3.11) 34 PAGE 35 Where and are the average smallscale stresses c onsisting of the viscous stress and the sm allscale Reynolds stress due to turbulence between the particles or due to particle fluctuations. The Favre aver aging process results in largescale Reynolds stresses and due to turbulent fluctuations in which concentr ation fluctuation is take n in to account. These Reynolds stresses are defined as co rrelations between the concentra tion and velocity fluctuations repres enting the transfer of momentum and are given as 0 f xz0 f zzf xzRf zzR (3.12) (1); (1)ffffffff xzzzRcuwRcw w Turbulent eddy viscosity param eterization is used to solve for the largescale fluid Reynolds stresses. Applying the current clos ure schemes the stresses are given as z uf fft ff xz (3.13) and z w kcf fft f f ff zz 3 4 )1( 3 2 (3.14) It is noted from the continuity equation that the divergence of the fl uid phase velocity is not zero; therefore the second term in the above equations is nontrivial. The kinematic viscosity of the fluid is the fluid phase eddy viscosity is and is the fluid phase turbulent kinetic energy given by fftfk f i f i fuuc c k 1 12 1 (3.15) The fluid phase eddy viscosity is given by, 35 PAGE 36 f f ftck C 12 (3.16) In which is an empirical coefficient and the fluid pha se tur bulent energy dissipation rate is given by C j f i f ij f fx u c c 1 1 1 (3.17) The balance equation for th e turbulent kinetic energy and the turbulent energy dissipation rate is given by fkf 12 )1( )1( )1( )1(f ssf ftf f f f k ft f ft zz f ft xz f f f f fkcww z c c z kc zz w z u wkc z kc t (3.18) The first term on the right hand side is the advection term. The second and the third terms represent the production of kinetic energy or turbul ence. The fourth term indicates the turbulent diffusion term. The next term represents turbul ent or viscous dissipation followed by dissipation due to the turbulent suspension. The last term, which basically involves the corr elation between fluid and sediment velocity fluc tuations, represents the dissip ation of turbul ent energy where is a measure of the degree of correlation. The par ticle response time is given as the time taken for a particle at rest to accelerate to the velo city of the surrounding flui d. From Drew (1976) this time is given by pt s pt (3.19) 36 PAGE 37 The tim e between collisions is calculated based on the mean free path of the colliding pa rticles and the strength of the sediment velocity fluctuations ctcl 2 1 sk 2 1 s c ck l t (3.20) The length scale is defined by cl cgc d lc 024 (3.21) The definition of cg0, the contact value of the radial dist ribution function, is given in the next section. Based on the above definitions, the parameter is found to be dependent on the relative magnitudes of the particle response time and the fluid turbulence time scale which is given by Elghobashi and AbouArab (1983) as ptlt f f lk t165.0 (3.22) Additionally, in equation (3.18) also depends on the time between collisions and is given by ct 1,min 1 cl ptt t (3.23) The balance equation for the turbulent ener gy dissipation is give n from Elghobashi and AbouArab (1983) as 37 PAGE 38 sf ft f f f s f f f f f f f f ft f zz ft f xz ft f f f f f f fww z c k C kc k C c k C z c z z w z u k Cwc z c t 3 3 2 1 12 )1( )1( )1( )1( (3.24) In the absence of appropriate values for the coefficients in the numerical values for t he clear fluid from theffk k taken from Rodi (1984) are used. In the above balance equation, if the averaged sediment concentration is set to zero the standard kequation for a clear fluid is obtained. The set of the numerical co efficients used in this particular equation is given in table 3.1 Table 3.1 Summary of numerical co efficients adopted by the model C 1C 2C 3C k sC s 0.09 1.44 1.92 1.2 1.0 1.3 0.55 1.0 3.4 Closure of Sediment Stress 3.4.1 Sediment Transport Equation Like the fluid stresses, the sediment stre ss can be written down in terms of the sum of the smallscale particles stress and th e largescale Reynolds stress (3.25) s zz s zz s zz s xz s xz s xzR R 0 0 ; The sm allscale stresses are due to interparticl e interactions resulting from collisions or interstitial fluid effects. The largescale Re ynolds stresses are due to sediment velocity fluctuations and associated concen tration fluctuations given as 38 PAGE 39 (3.26) ssss zz ssss xzwwcR wucR ; The sedim ent fluctuation energy, which is an alogous to the fluidphase turbulent kinetic energy, is given as s i s i suuc c k 2 1 (3.27) Since this closure method is more suitable to the largescale stresses, it is assumed that this closure can describe sediment fluctuations even at the smallscale and in the smallscale limit we let = where represents the smallscale velocity. However, at small scales concentratio n fluctuations are abse nt, and as a result the fluctuating part of the concentration is ignored. Accordingly, the concentr ation is the mean value obtaine d from the averaging process. Therefore the above sedimentphase turbul ent kinetic energy can be written as s iu s ius iu s i s isuuk 2 1 (3.28) The concept of granular temperature is intr oduced from Nott and Brady (1994) and Jenkins and Hanes (1998) and the above equation is written as s i s isuuT 3 1 (3.29) Where the brackets denote ensemble averaging. Analogous to the fluid phase turbulent energy, the sediment fluctuation energy, which is the tur bulent kinetic ener gy counterpart for the solids phase, can be derived from sedimentphase momentum equation. The smallscale granular temperature and the sediment fluctuation ener gy are coupled together to provide a single transport equation as follows from Hsu (2002) 39 PAGE 40 sf s s zz s s xz s ss sKkc z Q z w z u z wKc t Kc 2 (3.30) Where denotes the sediment fluctuation energy flux and Q is the rate of energy dissipation. The total sediment stress is used to encompass the smallscale and the largescale stresses. The last term in the above equation is similar to th e last term of the fluid phase turbulent kinetic energy and describes the interaction between the two phases. This can be written down as two terms. One term is the fluid turbulent kinetic energy fkc 2associated with the turbulent eddies acting on the random sediment particles and thereby enhancing th e sediment fluctuation energy. The other term is associat ed with the drag due to the interstitial fluid and additional dissipation mechanism present in the form sKc2. The above equation would need further closures to solve the transport equation in terms of the sediment stresses, flux of the fluctuating energy, and energy dissipation. Let the flux of the sediment fluctuation energy be the sum of the small scale and the large scale com ponents. Therefore, Q0Q1Q (3.31) 10QQQ and energy dissipation is taken to b e collisional dissipation associated with the inelasticity of particles as explained next. 3.4.2 Collisional Theory of Granular Flow Following Jenkins and Hanes (199 8) the closure for the smallscale intergranular stresses is given by the kinetic theory of collisional granul ar flow. This theory is an extension of the classical kinetic theory applie d to dense particle flow and provides explicit closures to dissipation due to the inelastic particleparticl e collisions by introducing a coefficient of 40 PAGE 41 restitution. The assum ption made in th eory is that the normal stress is function of the sedim ent concentration, sediment properties and sediment fluctuation energy. Therefore, s zz z w AEKGcs s s s zz 41 3 20 (3.32) where cgcG0, AE, which is the product of A and E is the sediment viscosity due to collisions and cg0 is a radial distribution function at contact for iden tical spheres, provided by the following expression from Torquato (1995) 0.635c0.49 64.0 8537.0 49.00 22 23 0 pc c c c cg (3.33) In the above equation, p =1. The particle collisiona l shear stress is given by z u AEs s xz 0 (3.34) where, from Jenkins and Hanes (1998) 1 2 25 2 8 ;11 312 s AdcGKE s G 8 (3.35) and, 2 059 11 3sK QAMM z 5 3 2 1 2 G (3.36) Finally the rate of dissipa tion is closed as follows from Jenkins and Savage (1983), s s sKe z w Gc d A 1 4 102 (3.37) 41 PAGE 42 The coefficient of restitu tion e is taken as 0.8. It should be noted that e represents the ratio of the velocities before and after an im pact. For a completely elastic collision e =1and for completely inelastic collision e =0. 3.4.3 Large Scale Sediment Stress The large scale sediment stresses are closed using the eddy viscos ity model. The shear stress is closed as z u Rs st ss xz ~ (3.38) and the normal stress as z w Kc Rs st s s ss xz ~ 3 4 3 2 (3.39) The sediment viscosity is connected to th e sediment fluctuation energy through the sediment mixing length sl sssstKlcC (3.40) where Cs = 0.55.The sediment mixing length is connect ed to the turbulent fluid flow mixing length, as the large scale velocity and concen tration fluctuations are caused by turbulence through the factor as follows (3.41) fsll If the particles are large, then the particle response time is larg e, indicating that the turbulent eddies do not c ontribute to the sedim ent velocity fl uctuation, which in turn would mean that 0 If the particles are fine, then the particle response time will be small, indicating that turbulence does play a role in the sedime nt velocity fluctuation. In that case 1 and the 42 PAGE 43 sedim ent mixing length is nearly equal to th e turbulent mixing length. The turbulent mixing length is given by, f f fk l2/3165.0 (3.42) The energy flux due to largescale se diment fluctuations is given by z K Qs s st s 1 (3.43) where is taken to be 1 for simplicity. s3.4.4 Closure at Low Concentrations The theory of gases, which is used as an analogy for the collisional granular flow theory, would not be a good assumption when the suspension is dilute. Sediment concentration becomes dilute near the flow surface and collisions among particles become infrequent, which violates an important postulate of the kinetic theory for it to remain statically valid. In such a scenario fluid turbulence acts as a major cause of sediment susp ension. As turbulence gives rise to largescale Reynolds stresses, the smallscal e dilute suspension stresses and must be closed appropriately. For this a dam ping parameter 0s xz0 s zz is introduced in terms of the mean free path and the turbulence m ixing length clfl cf fll l (3.44) If the mean free path increases then the sediment is diluting further and the smallscale stresses keep decreasing by the factor Using and from equations ( 3.25), (3.34), (3.38) and (3.44) the total shear stress is 43 PAGE 44 z u AERs st s s xz s xz s xz ~0 (3.45) and the normal stress is z w AEKGc Rs st s s s s zz s zz s zz ~ 3 4 141 3 20 (3.46) Thus the energy flux is given by z K AM QQQs s st s 3 510 (3.47) 3.4.5 Region of Enduring Contact The kinetic theory of granular flow can be applied when the c oncentrations are low so that the time for interparticle collision, i short and in conformity with the kinetic theory of gases. However, when the concentration is higher the coll ision time will not be short, as there will be enduring contacts between the partic les. Therefore kinetic theory can no longer be applied. In between random loose packing c* and random close packing c* sediment behaves in a transitional way from having fluidlike properties to having solidlike behavior. As the concentration increases the collisiona l contribution of the pa rticle normal stress decreases because shearing reduces between the particles and the fluctuations become very small. However, the contribution from enduring contacts increases the norm al stress. Therefore it may be assumed that the particle normal stress is a sum of the collisional normal stress and normal stress due to enduring contact: (3.48) se zz sc zz s zz The collisional norm al stress is closed using equation ( 2.32), whereas for the normal stress due to enduring contacts, a Hertz cont act relation is used and given as 44 PAGE 45 2 3 2 d cK d mse zz (3.49) where is the average compressive volume strain, K is the average nu mber of contacts per particle or coordination number and m is given as follows, 1 39 22d me (3.50) where is the shear modulus and e is Poissons ratio. The ratio of d is related to the concentration difference as follows, 3 2 cc d (3.51) Therefore from the above equations, ccc c cc 02 cccK d mse zz (3.52) and the coordination number K is taken to be a function of concentration, *cc 1 2 2 sin33 c cc cc cK (3.53) 3.5 Initial and Boundary Conditions 3.5.1 Initial Condition The sediment bed, when agitated by turbulent fl ow, is disturbed from the its stationary position and the granules are entrained .They either are kept suspended and move along with turbulent motion or settle down and exchange mo mentum with other stationary particles and 45 PAGE 46 attem pt to dislodge them. This process of sediment transport and its initiation cannot be studied by examining the interaction proce ss of every individual particle. Th erefore, in order to initialize sediment motion an assumed sediment profile zcini is used as the initial condition. The maximum concentration c* is taken at the bed and may simply be assumed to decrease linearly to a concentration of 0.4 at about half the wate r depth. Initially the fluid and sediment phase velocities are considered to be zero, i.e. from Hsu (2002). To accelerate th e computations the sediment vertical velo city is set to zero and the flow is calculated up to an initial time which is taken to be 10 tim es the turnover time 0 ~~ and 0 ~~ sf sfww uu0T 0u h t0TLof the largest eddy with h being the water depth. After time the vertical motion of sediment is calculated using the sedim ent momentum equation. 3.5.2 Bottom Boundary Condition The bottom boundary condition for the sediment pha se is applied at the interface where the porous stationary part of the sediment bed fails, the sediment is sheared and the granular bed is fluidized leading to the formation of a region where the bed moves slowly creating a condition conducive to enduring contacts among particles. If the concentration of particles is larger than the random closepacked concentration then fluidizing the partic les would not be possible ; therefore a failure concentration *cc is defined as the concentration at which the particles are first sheared. The applied external shear can result in erosion or settling and accumulation on the bed. In the porous stationary bed layer, the sediment horizontal velocitysu ~ the vertical velocitysw ~ and the sediment fluctuation energy disappear, giving rise to the condition in the stationary bed sk 0 0 ~ 0 ~ s s sk wu (3.54) 46 PAGE 47 The bed location, which depends on the external sh ear and changes according to the external flow conditions, is dete rmined as part of the solution. Th e failure concentration is given by the Coulomb failure criterion as (3.55) tans zz s xz where is the angle of friction of the sediment. Th erefore, from equations (3.48) and (3.52) and (3.55) sc xz s xz nm d cccK tan2 (3.56) Along with the implementation of the boundary conditions given in (3.54), and are obtained from their clos ures, and solving the above nonlinear equation (3.56) the failure concentration s xzsc xzc is obtained. At this point the bed is set in motion and its location is specified. The stress changes with the external flow condition indicating that s xzc and the bed location also vary along with external flow. 3.5.3 Top Boundary Condition Fluid Phase: The top boundary is cons idered to be a rigid hydraulically smooth lid. In a fully developed turbulent flow the horizontal velocity profile is assumed to follow the logarithmic law of the wall given as z u u uf t t f* *9 ln 1 ~ (3.57) where the von Karman constant in the sedimentladen flow is taken as 0.41and is the friction velocity at the top of th e rigid lid (Hsu, 2002). Given the fl uid velocity calc ulated in the tu* 47 PAGE 48 flow dom ain near the top lid the friction velocity can be f ound and the horizontal velocity gradient at the rigid lid is given as f t fu z u* ~ (3.58) The vertical fluid velocity at the top of the boundary .The noflux boundary condition with respect to the turbulent kinetic energy is used for the top boundary condition 0 ~fw 0fk z (3.59) The turbulent energy dissipation rate f is specified based on the assumption that production equals dissipation in the logarithmic region: z kCf f 2/34/3 (3.60) It is noted from equation 3.18 that there are two more terms that are responsible for dissipation like dissipation due to drag and dissipation due to turbulent suspension flux. Hence it is not steadystate turbulence. The Neumann boundaryconditions are applied to the lateral and bottom fluid pressures, and 0P f is applied to the top boundary as a reference pressure. Sediment phase: At the top the sediment concentration is usually zero or dilute. In the present model the top of the sediment pha se is considered to be the in terface at which the sediment concentration becomes less than a specified concentration of Above this interface the concentratio n is considered to be zero. The calculations (next chapter) are insensitive to as long as it is sm all. In the present case has been taken to be mincmin cmincminc4105The horizontal velocity of the sediment phase has a free slip velocity at the top boundary 48 PAGE 49 0 z us (3.61) indicating that the sediment shear stress is neglig ible for low sediment concentration at the top boundary. It should be noted that this boundary is not the lid but where cmin is defined. Similarly the noflux condition is applied to the sediment fluctuation energy at the top boundary for the sediment phase, 0 z ks (3.62) The location of the top boundary of sediment ch anges according to the external conditions similar to the porous stationary bed. The sediment concentration and the vertical velocity are considered to be zero above the top boundary. 3.5.4 Lateral Boundary Condition The flow is driven by a horizont al pressure gradient. In expe riments that involve steady flow in an open channel, the flow is driv en by gravity through th e horizontal momentum equation. In the present case the gravity term in the fluid momentum equation is converted to the horizontal pressure gradient term .Therefore for the steady gravity driven flow with a given slope or a frictional velocity and hydraulic radius the horizon tal pressure gradient is given as S*ubr b f fr u gS x P2 *1 (3.63) 49 PAGE 50 CHAP TER 4 NUMERICAL IMPLEMENTATION 4.1 Introduction The fluid and sediment phase continuity equations (3.1) and (3.2) along with the momentum equations (3.3), (3.4), (3.9) and (3.10) define the dynamics of the system to be examined. The turbulent kinetic energy balan ce equation (3.18) and the balance equation for turbulent energy dissipation rate (3.24) along with the sediment phase fluctuation energy balance equation (3.30) are solv ed numerically using a finite difference scheme. The implementation scheme is described in this chapter. 4.2 Mesh and Grid System The mesh formed for the computational domai n is implemented usi ng the staggeredgrid approach. In this the fluid horizontal velocity the sediment horizontal velocity fluid pressurefu~su~ f P sediment concentration c, fluid turbulence kinetic energy the turbulent energy dissipa tion ratefkf and the sediment fluctuation energy are calculated at the center of the grid while the sediment phase vertical velocity and the fluid phase ver tical velocity are taken at the top face of the grid. sksw~fw~The quantities do not vary in the downslope ( x ) direc tion; i.e. a uni form flow is assumed. Variations in the z direction are solved for. The computa tional domain is discretized into N+2 grids with single ghost grids at the top and the bottom. As mentioned previously, the flow is driven through a horizontal pressure gradient /fPx and to facilitate these two columns of ghost grids are specified to allow implementation of the pressure gradient. Since only a vertical column in considered for the calculations, th e ghost grids help in se tting up the problem by letting flow develop in the pseudo sense in the horizontal direction. 50 PAGE 51 4.3 Numerical Discretization A m odified twostep projecti on method is used to solve the twophase equations. The numerical scheme proceeds by first discretizing the sediment phase continuity and momentum equations and the necessary sediment phase values like c, and are obtained by the predictorcorrector m ethod .Once these values ar e updated at a new time level, the fluid phase mass and momentum equations are solved agai n using the predictorcorrector method. With these values solved for, after every computational cycle, the and su~fksw~skf equations are solved and updated. 4.3.1 Solution of Sediment Phase Equations The values of and are not directly solved from the governing equations; instead, the following sediment phase vari able definitions are introduced. su~sw~sk The sediment horizontal momentum per unit mass given as ssucU ~ (4.1) Next, the sediment vertical momentum per unit mass is sswcW ~ (4.2) and from the sediment fluctuat ion energy equation, the fluctuat ion energy per unit mass, or specific fluctuation energy, is given as s skcK (4.3) These definitions help avoid any singular ities in the governing equations as the concentration approaches zero. Based on the calculated values of c,U, Wand K the predictorcorrector scheme is described next. 51 PAGE 52 4.3.1.1 Predictor Step The forward tim e difference method is used to di scretize the time derivative. Initially there is a predictor step in which at a new timestep ( n) the sediment phase variables, c,, W, U K are calculated from the known values at the previous time step ( n1). Let j denote the grid index in the z direction. The sediment phase continuity equa tion (3.2) is discretized in the predictor scheme as follows 1 1 n j s n jRHSC t cc (4.4) where c indicates the tentative sediment concentration at timestep n and RHSC indicates the righthand sediment concentration in terms of sediment vertical momentum per unit mass given as j n j n j n jz W W RHSC 1 2/1 1 2/1 1 (4.5) From equation (3.9) the sediment horizontal momentum equation per unit mass is given as follows t t RHSU U Un j n j n j j 1 1 11 (4.6) where is defined in equation 3.5 and 1 n jRHSU is given as 1 1~ n j f s n juccSgSTXCVX RHSU (4.7) and stands for the advection term in the x direction. W hen compared to equation (3.9) it is denoted as CVX 52 PAGE 53 z wU CVXs~ (4.8) This term is discretized with a cen tral difference scheme around the point j as follows j s j j s j j j s jz wUwU z wU CVX 2/1 2/1 2/1 2/1~ ~ ~ (4.9) STXin equation (4.7) stands for th e stress term in the horizonta l direction and is given as z STXs xz (4.10) This term is centered on point j and is discretized as follows: j j s xz j s xz j s xz jz z STX 2/1 2/1 (4.11) 2/1 jU in equation (4.9) is based on a combinati on of central and upwind scheme given as 1 2/1 2/1 2/1~ 1 2 1 ~ 1 2 1 j s j j s j jUwSIG UwSIG U (4.12) which is the sign function which depends on the di rec tion of the vertical velocity and SIG is a numerical coefficient taken to be 1 for upwind and 0 for a central difference scheme in the present model. In equation (4.12) 2 ~~ ~2/12/3 1s j s j s jww w (4.13) 53 PAGE 54 From the vertical sediment m omentum equation (3.10), the vertical specific momentum is given as t t RHSW W Wn j n j n jj 1 2/1 1 2/1 1 2/11 2/1 (4.14) Note the difference in the vertical grid index. This is as a result of the staggered grid where the sediment horizontal and vertical ve locities are calculated at different grid points. In the above equation is given by 1 2/1 n jRHSW 1 2/1 1 2/1~ n j f s n jwccgTSPFSTZCVZ RHSW (4.15) and similar to the horizontal momentum equation, CVZ represents the convection term in the z direction and given as z wW CVZs~ (4.16) Centered on grid index j+ 1/2, the vertical convection term is discretized as follows: 11 1/2 1 1/2(j)/2 s s s j jj j jj jWwWw Ww CVZ zzz (4.17) STZin equation (4.15) is the stre ss in the vertical direction, z STZs zz (4.18) This term is calculated as 54 PAGE 55 11 1/2 1 1/2()ss s zzzz jj zz j jj jSTZ zzz / 2/2 (4.19) PFin equation (4.15) is the pressure term, z P cPFf (4.20) and is discretized as 1 1/21/2 1 1/2(ff f jj jj jj jPP P PFcc zz )/2 z (4.21) 2/1 jcin the above equation is obtaine d as linear interpolation 1 11 2/1 jj jjjj jzz czcz c (4.22) and TS is the turbulent susp ension flux term in the vertical direction z c TSft (4.23) This term is calculated as 1 1/21/2 1/2 1/2 1(jj jftjft j j jjcc c TS zz )/2 z (4.24) With a linear interpolation of the fluid eddy viscosity 1 1 1 2/1 jj j ftj j ftj j ftzz z z (4.25) 55 PAGE 56 For the above equations, is given by a combination of the central and the upwind schem e as follows: 1 jW 2/3 2/1 2/1 2/1 1~ 1 2 1 ~ 1 2 1 j s j j s j jUwSIG WwSIG W (4.26) Similarly, from sediment fluctuation energy e quation (3.30), the specific fluctuation energy is given by t t RHSK K Kn j n j n j j 1 0 1 121 (4.27) 0 being part of the collisional dissipati on term equation (3.37) in the absence of K and where 1 12 n j f n jkcDIF PROD CVK RHSK (4.28) Where is defined in equation (3.23) and CVKbeing the convection term, z wK CVKs~ (4.29) This term is calculated simila rly to convection term in equation (4.9) and is given as j s j j s j j j jz wKwK z UK CVK 2/1 2/1 2/1 2/1~ ~ (4.30) PRODin equation (4.28) is the production term, z w z u PRODs s zz s s xz~~ (4.31) This term is calculated as 56 PAGE 57 2 ~ ~ 2 ~ ~ ~~2/1 2/1 2/1 2/1 j s j s j s zz j s j s j s xz s s zz s s xzjz w z w z u z u z w z u PROD (4.32) The terms() ,s x zj1/2(/)s juz and can be obtained from interpolation as shown previously. 1/2(/)s juz The diffusion term in equation 4.28 is defined as DIF z Q DIF (4.33) This term is calculated similarly to the stress term in the x direction as given in equation (4.11) according to j jj j jz QQ z Q DIF 2/12/1 (4.34) Similar to equations (4.12) and (4.26) is given as 2/1 jK 1 2/1 2/1 2/1~ 1 2 1 ~ 1 2 1 j s j j s j jKwSIG KwSIG K (4.35) In the same way,1/2,jU, j W2/1 jKare obtained similarly to e quations (4.12), (4.26), (4.35). 57 PAGE 58 4.3.1.2 Corrector Step In the corrector step th e predicted value of th e time derivative is obtained at the next time step by substituting the newly obtained predicted values of jc and calculate and These quantities are obtained from jU 2/1jWjK jRHSCjRHSUjRHSWjRHSK j jj jz WW RHSC 2/12/1 (4.36) j n f s juccSgSTXCVX RHSU 1~ (4.37) 2/1 1~ j n f s jwccgTSPFSTZSVZ RHSW (4.38) and j n f jkcDIF PROD CVK RHSK12 (4.39) All the terms on the right hand sides of equatio ns (4.36) to (4.39) with a circumflex (cap) denote the same corresponding quantities as explai ned earlier in the predictor step but are calculated using the newly obtained values of jc from the predicted step. The final corrected value is obtained by averaging the values calculated from the predicted and the corrected steps: jU 2/1jWjK j n j s n j jRHSC RHSC t cc1 12 (4.40) 11 12 1 2nn jj j n jjt URHSCRHSU U t j (4.41) 58 PAGE 59 2/1 1 2/1 2/1 1 2/1 1 2/1 2/1 2 1 2 j n j j n J n j jt RHSW RHSW t W W (4.42) and j n j j n j j n j n j jt RHSK RHSK t K K0 1 0 1 1 1 2 2 1 2 (4.43) The twostep m ethod is employed throughout th e computational domain until convergence is satisfied. The convergence crite rion is reached when the re lative difference between two successive timesteps is less than 0.001. 4.3.2 Fluid Phase Equations Similar to the sediment phase equations, the fluid phase equations are solved using the twophase projection method. Incorporating the flui d phase continuity equation in the horizontal and vertical momentum equations, one obtains sf f f xz f f f f f fucuc c Sg z c x P z u w t u ~~ 1 1 1 1 ~ ~ ~ (4.44) z c c ucuc c g z c z P z w w t wft f sf f f zz f f f f f f 1 ~~ 1 1 1 1 ~ ~ ~ (4.45) 59 PAGE 60 The horizontal fluid pressure gradient term in eq uation (4.44) is a known quan tity as it represents the driving force for the flow. The horizo ntal momentum equation is given as nj f n j j s n nf n n n nf j nf jc t U c Sg SFX CFXtu u 1 1 1 ~ ~1 1 1 (4.46) The vertical fluid velocity is calculated without the vertic al pressure gradient as nj f n j j n n ft s n n f n n n nf j f jc t g z c W c SFZ CFZtw w2/1 2/1 2/1 1 1 1 1 2/1 2/11 1 1 ~ ~ (4.47) The vertical pressure gradient in the mome ntum equation is unknown and it solved along with the fluid continuity equation and the pressure Poisson equation j n f j nf nt c z c tz P c z 1 1 1 (4.48) The above equation results in a tridiagonal matrix and is solv ed using the Thomas algorithm (Hsu, 2002). The fluid phase convection termsCFX, and the stress terms are calculated sim ilarly as in the sediment pha se equations (4.9), (4.17), (4.11) and (4.19). CFZ SFX SFZAfter the pressure is updated at tim estep n the vertical fluid velocity is updated at the same time level using 2/1 2/1 2/11 ~~ j nf f f j nf jz P ww (4.49) 60 PAGE 61 4.4 Stability Analysis The tim estep is dynamically adjusted at ever y time level based on th e numerical stability criterion. From the convective terms of fluid and sediment phase momentum equations, the following Courant condition is obtained sfww z t ~ ~ max1 (4.50) Where the numerical coefficient 1 is taken as 0.3 (Hsu, 2002). The numerical stability due to the diffusion te rm in the fluid and sediment equations is defined by st s ftAE z t /,max2 (4.51) Where st sAE / is the diffusion coeffici ent in equation (3.45). The interaction terms in the momentum equation represent the drag force. Constraining the timestep due to the drag force indicates that the size of the timestep is proportional to the particle response time divided by the partic le sp ecific gravity pt s .This is given as s t tp 2 (4.52) Where 2 =0.1 is a numerical coefficient. The mini mum of criteria (4.50), (4.51), (4.52) is considered as the minimum size of the timestep. 61 PAGE 62 CHAP TER 5 RESULTS AND ANALYSIS 5.1 Preamble apter model simulation results for a se lected set of conditions at the head of the bed w physical interpretation of model use is shown in Fig. 5.1 on a schematic basis. Fig. 5.1a s In this ch ith respect to the initial water level (hel d constant), particle size currentinduced forcing and suspended sediment concentration are presente d. Referring to Fig. 5.1 the model is run first in a pretest mode with a horizontal bottom, which permits the development of the initial conditions with respect boundary la yer velocity and suspended sediment concentration. This run is followed by the test run from which as assessm ent is made of the state of the flow as it develops. The chematizes the pretest mode in which a do wnstream gradient in the pressure force in a horizontal channel leads to the development of the velocity and concentration profiles. During the test run for a sloping channe l, sediment suspended is shown to deposit, which means that autosuspension does not occur in this case. The depthmean suspended sediment flux F ( x0) (per unit width of channel) at any distance x0 along the slope is defined as 00 01 (,)(,)h F uxzcxzdz h(5.1) here c is the concentration in units of kg/m3. Absence of ig w nition would mean that the flow is decelerating, i.e./0xUdUdx where U is the depthmean current velocity. Also the flux gradient/xFdF ative. The opposite conditions prevail when the conditions are conducive to ignition (Fig. 5.1b). T hus in the present study igniti on or the ignitive condition is defined as one in which dU / dx > 0 and dF / dx > 0 and autosuspension is defined as the stage if the Shields parameter is above the Bagnold auto suspension curve as shown in Figure 5.10. It should be pointed out that a condi tion can also exist when these tw o criteria are sa tisfied and yet dx must be neg 62 PAGE 63 dc / dx is n egative, i.e. a depositing turbidity current occurs even when the flow is accelerating so rapidly that dF / dx is also positive. However, such a c ondition cannot last i ndefinitely, because the accelerating current will eventually erode th e bed at such a rate that the gradient dc / dx will become greater than zero. In the present analys is the pretest water de pth is kept constant (at 0.12 m) and conditions for ig d corresponding test ru nition are exam ined by varying the particle size d, the pressure gradient and the bed slope For particle size, three me diumsized beach sands typical ly found in Florida (Dean and Dalrymple, 2002) are selected; these bei ng 0.21 mm, 0.28 mm and 0.35 mm. The pressure gradient is specified in terms of the friction velocity u* per Eq. 3.63. Since the channel is assumed to be hydraulically wide, the hydraulic radius (rb) is replaced by water depth h. The bed slope is varied between 0o and 17o; however, results for slopes lower than 12o indicated that no ignition could occur for the selected range s of diameters and friction velocities. The choice of three values each of d, u* and lead to 54 pretest runs an ns. The selected combina tions of these three variables are given in Tables 5.1, 5.2 and 5.3. The outcome of the pretest mode of model simulation supplies the initial velocity profile and the concentration profile fo r the turbidity current on the slope. As an illustration, model output in terms of the development of the velocity field and concentration field for d = 0.21 mm, u* = 0.06 m/s and = 0o is shown in Fig. 5.2. Note that time is represente d along the abscissa as the model is 1DV. The time profile of the velociti es taken at different he ights is plotted. It is observed that in the absence of any driving fo rce after 150 s, the flow dies down. Another perspective of the same occurr ence can be observed in Figure 5.2b in which the ve locity profile is plotted at different times. It is noticed that from the initia l point at 150 s where the pressure forcing is stopped and the gravitatio nal components are required to take over, in the absence of 63 PAGE 64 any slope, the velocity decreases. In Fig. 5.2c, in which the volum etric concentration profile is plotted, it is discernible that th ere is no suspension; in fact the bed suspension simply dies down. It is also noticed from Fig. 5.2b th at the velocity profile does not quite begin from the origin but at an ordinate above the origi n. This region includes the solid bed below in which the flow was not able to suspend any sediment and also the bed over which these is a flow. This feature can also be noticed in Fig. 5.3 and 5.4 where the results from Parker et al. (1987) are compared against the model. 5.2 Experimental Data and Model Output dicted by numerical models (e.g. Eidsvik and Brs, 1989) elocity and concentration at three Ignitive turbidity currents have been pre ; however it appears that labor atory test results have been unable to fully verify modeling outcomes. As we noted in Chapter 3, the laborato ry tests of Parker et al. (1987) showed ignition in two tests out of the 24 they conducted. Th eir data on the evolution of velocity and concentration profiles offer the po ssibility of testing the perfor mance of the present numerical model. However, since the model has been conf igured by the original developer (Hsu, 2002) for specific types of twophase fl ow conditions (limited to sand grains and very highdensity packing), and because its reconfiguration was be yond the present scope, it was decided to use the model for a qualitative comparison w ith one dataset of Parker et al. Figures 5.3 to 5.5 show measured profiles (dashed lines) of v distances from the upstream end of the la boratory channel (shown in Fig. 2.5) in a single run. Solidline curves indicate the velocity and concentration profiles obtained from the model run. The input data used for this run were the closest possible th at could be used to simulate conditions in Parker et al. test The modeled diameter of sedi ment grain was 0.20mm, compared to 0.03mm in the experiment. The in itial layeraver aged velocity U0 in the model run was 0.13 m/s compared to 0.27m/s in the experiment. This velocity of 0.13 m/s was used to convert time64 PAGE 65 dependent profiles from the model output to distan ce dependent profiles plotted in the figures. The bed slope was 3o in both cases. Finally, the initial volume concentration Cv of the bottom material in the experiment was 0.0041, whereas the model was based on Cv = 0.635, a significantly high value. The velocity and concentrati on profiles at a distance of x = 1.5m in Figs. 5.3, 5.4 and 5.5 show ental evid Behavior Results of model test runs corresponding to th e conditions set in Tables 5.1, 5.2 and 5.3 are given in Tables 5.5, 5.6 and 5.7, respectively. These results are conveniently plotted in terms of the Shields parameter (along the ordinate) expected quantitative differe nc es but at the same time demonstrate that the model results are reasonable. From these plots the quantities/, FdFdt Fx and Ux are calculated and given in Table 5.4. The differences between the experim ence and the mode l output are due to differences in the sediment size. In the experi ments finegrained (but noncohesive) sediment was used. This is coupled with the fact that when compared to the model, the very low volumetric concentration in the experimental run led to an accelerating flow which was sensitive to the initial conditions specified by an initial velocity of 0.27 m/s and a sediment discharge of 0.165 kg/s. Parker et al. mention that in most runs they could obser ve an accelerating yet depositing current. As mentioned, one of the reasons they cite is that the short length of the flume in which the experiments we re carried out possibly did not allow them to witness ignition. In the model larger sediment was used with la rge concentration and a low angle of inclination which led to the damping of the flux. Also, as noted the concentration was very low in the experiments in comparison with the model, in which the stationary bed had random close packing concentration. 5.3 Domain of Ignitive 65 PAGE 66 2 (1) sgd s the specific weight of the particles, is taken as 2.65 in the present study (e.g. Julien, 1995). In Fig. 5.6, is plotted against bed slu (5.2) Where ope and includes the results from the model test runs. s mentioned earlier, in the 24 te sts they performed, o phic development. The points that corre spond to ignition are denoted by the symbol I and th gle where the flux is m any times A nly two were in the regime of accelerating flow and eroding sediment. Howeve r, as also mentioned previous ly they could not observe any catastro ose that correspond to a nonignitive state are denoted by N As an extension of Fig. 5.6, contour s of sediment flux per unit width F (kg/m3/s) for each grain size are plotted in Figs. 5.7, 5.8 and 5.9. It is difficult to identify the presence of an accelerating eroding current based on the intensity of the sediment flux, as erosion occurs albeit at a small rate in comparison to the next condition of an increased an higher. However, it can be said with certainty that selfsustaining capability is encouraged with increasing slope. It is also uncertain w hy, for the smallest diameter of 0.21mm in Fig. 5.7 the flux has its highest value at a slope of 16o but then decreases as the slope is increased. Such a trend is not observed in Figs. 5.8 or 5.9 for the larger diameters. The mass fluxes in Tables 5.5 to 5.7 were calculated at distance x = 20 m from the beginning of the slope. Either erosion or deposition is noticeable at this distance. To calculate the Shields parameter at this point, the depth mean velocity was evaluated and the friction velocity was calculated from duCU (5.3) 66 PAGE 67 The drag co efficient was calibrated to the value 0.018 by taking into account Fig. 5.10 such that the value of Cd did not cause the Shields parameter calculated from the friction velocity obtained from q. 5.3 to fall below the threshold of the It is observed from Fig. 5.7 th at the flux obtained does not increase with increase in the slope above 16o. Comparison with Fig. 5.10 can be shown to indicate that the maximum Shields parameter is not obtained at larg er slopes but rather at the condition where the initial forcing friction velocity at the start was a maximum. In Fig. 5.7 the maximum flux seems to be associated with the flow conditions and not uniquely to the slope. The plot in Fig. 5.8 seems to suggest a trans itory regime where the combined influence of the slope and the initial condition is witnesse d. It is observed that the flux increases with increase in slope but not as cons istently as in Fig. 5.9 where it is clearly obs erved that the flux increases with the slope. In Fig. 5.10 the Shields parameter is plotte d as a function of th e grain size along with Shields threshold curve for incipient grain movement and the Bagnold autosuspension curve. Data points from the model test runs are superimposed. It is observed that all the data points denoted by I (ignitive) and N (nonignitive) in Fig. 5.6 represen t autosuspensi ons. This trend is observed sediment of diameter 0.21mm and 0.28mm.Whereas for sediment diameter 0.35mm three of the cases are between bedload and auto suspension. These three cases are for the slopes of 12o, 13o and 14o and the initial forcing condition in each of them is u = 0.06 m/s (at the onset of the pretest run ). This is in contradicti on to the condition of se lfacceleration where both dU / dx and dF / dx are greater than zero as it includes all the data sets which were mentioned as nonignitive in nature. E incipient m otion. 67 PAGE 68 Table 5.1 Model runs for grain size of 0.21mm Run no. Particle diameter, d (mm) Initial friction velocity, u* (m/s) Bed slope, (degrees) 1 0.21 0.06 12 2 0.21 0.07 12 3 0.21 0.08 12 4 0.21 0.06 13 5 0.21 0.07 13 6 0.21 0.08 13 7 0.21 0.06 14 8 0.21 0.07 14 9 0.21 0.08 14 10 0.21 0.06 15 11 0.21 0.07 15 12 0.21 0.08 15 13 0.21 0.06 16 14 0.21 0.07 16 15 0.21 0.08 16 16 0.21 0.06 17 17 0.21 0.07 17 18 0.21 0.08 17 Table 5.2 Model runs for grain size of 0.28mm Run no. Particle diameter, d (mm) Initial friction velocity, u* (m/s) Bed slope, (degrees) 19 0.28 0.06 12 20 0.28 0.07 12 21 0.28 0.08 12 22 0.28 0.06 13 23 0.28 0.07 13 24 0.28 0.08 13 25 0.28 0.06 14 26 0.28 0.07 14 27 0.28 0.08 14 28 0.28 0.06 15 29 0.28 0.07 15 30 0.28 0.08 15 31 0.28 0.06 16 32 0.28 0.07 16 33 0.28 0.08 16 34 0.28 0.06 17 35 0.28 0.07 17 36 0.28 0.08 17 68 PAGE 69 Table 5.3 Model runs for grain size of 0.35mm Run no. d on Bed slope, Particle In itial fricti iameter, d (mm) velocity, u* (m/s) (degrees) 37 0.35 0.06 12 38 0.35 0.07 12 39 0.35 0.08 12 40 0.35 0.06 13 41 0.35 0.07 13 42 0.35 0.08 13 43 0.35 0.06 14 44 0.35 0.07 14 45 0.35 0.08 14 46 0.35 0.06 15 47 0.35 0.07 15 48 0.35 0.08 15 49 0.35 0.06 16 50 0.35 0.07 16 51 0.35 0.08 16 52 0.35 0.06 17 53 0.35 0.07 17 54 0.35 0.08 17 Table 5.4 Flow state in experiment and model run ( Ux (1/s) Flow state kg/m3/s2) Fx (kg/m4/s) F Reach (m) Data Model Data Mod Data Model el Data Model 1.54.5 8.505 4 3 .710 c./ dep Dec./dep 6.510 .1103 73 1.0031.2Ac 910 4.58.5 5.3103 7.3105 7.4 8.610 0 5 cc./ dep Dec./dep 103 4 46.0 .7103 A Table 5.5 Flowes for 0.21micles Run no. ( kg/m /s2) Fx g/m4/s) Ux (1/s) Flow state stat m partF 3(k 1 0.009 0.021 0.009 ecel./depos. D 2 0.020 0.028 0.008 ecel./depos. D 3 0.034 0.030 0.006 ecel./depos. D 4 0.008 0.016 0.007 ecel./depos. D 5 0.018 0.024 0.006 ecel./depos. D 6 0.032 0.026 0.005 ecel./depos. D 7 0.010 0.018 0.007 ecel./depos. D 8 0.020 0.025 0.005 ecel./depos. D 9 0.024 0.019 0.003 ecel./depos. D 10 0.008 0.014 0.005 ecel./depos. D 11 0.015 0.018 0.004 ecel./depos. D 12 0.025 0.018 0.002 ecel./depos. D 13 0.009 0.015 0.004 ecel./depos. D 14 1.804 0.52 4 0.015 Accel./eros. 15 0.005 0.00 3 0.0004 Accel./eros. 16 0.017 0.02 0 0.002 Accel./eros. 17 0.014 0.01 2 0.001 Accel./eros. 18 0.062 0.03 0 0.005 Accel./eros. 69 PAGE 70 Table 5. Run no. m /s Fx m4/s /s ) state 6 Flow states for 0.28mm particles F 3 ( kg/2) (k g/ ) Ux(1 Flow 19 0.008 0.026 0.013 l./depos. Dece 20 0.015 0.033 0.011 ecel./depos. D 21 0.026 0.037 0.008 ecel./depos. D 22 0.008 0.026 0.011 ecel./depos. D 23 0.016 0.031 0.008 ecel./depos. D 24 0.022 0.029 0.006 ecel./depos. D 25 0.009 0.026 0.009 ecel./depos. D 26 0.016 0.028 0.006 ecel./depos. D 27 0.015 0.017 0.004 ecel./depos. D 28 0.010 0.025 0.006 ecel./depos. D 29 0.015 0.022 0.003 ecel./depos. D 30 0.668 0.46 0 0.010 Accel./eros. 31 1.042 0.39 7 0.015 Accel./eros. 32 1.330 0.39 2 0.014 Accel./eros. 33 0.009 0.00 8 0.002 Accel./eros. 34 1.073 0.37 0 0.012 Accel./eros. 35 1.025 0.33 3 0.011 Accel./eros. 36 1.190 0.32 7 0.010 Accel./eros. Table 5.7 Flow states for 0.35mm particles ( kg/m /s2) F Ux Flow sta Ru nF no. 3x(k g/m4/s) (1/s ) te 37 .007 016 ecel. 0 0.034 0. D ./depos 38 0.01 4ec 4 0.041 0.01 D el./depos. 39 0.02 1ec 3 0.044 0.01 D el./depos. 40 0.009 0.037 0.014 Decel./depos. 41 0.014 0.037 0.011 Decel./depos. 42 0.021 0.035 0.007 Decel./de pos. 43 0.009 0.036 0.011 D ecel./depos. 44 0.014 0.031 0.007 Decel./depos. 45 0.014 0.020 0.003 Decel./depos. 46 0.007 0.022 0.005 Decel./depos. 47 2.063 0.57 4 0.018 Accel./eros. 48 1.725 0.48 0 0.015 Accel./eros. 49 2.524 0.48 4 0.018 Accel./eros. 50 2.401 0.44 5 0.017 Accel./eros. 51 2.226 0.39 7 0.015 Accel./eros. 52 2.364 0.39 7 0.016 Accel./eros. 53 2.266 0.37 6 0.015 Accel./eros. 54 2.164 0.35 7 0.015 Accel./eros. 70 PAGE 71 Figure 5.1 Physical inte rpretation of m odeations: wnstre in pressure force leads to the developmen boundar veloc ity profile and associated profile in the suspended sedim concentrUnder n itive cos the n decreases due to sediment deposition; (b) conditions arucive foion. simul (a) do am gradient t of a a y layer n ent tion. oni g ndition io concentrat e cond r ignit 100 11 0 12 0 130 14 0 150 160 17 0 180 190 200 0 0. 5 1 1.5 time ( s)velocity (m/s) 0 0.2 0.4 0.6 0.8 1 1.2 0 0.05 0.1 Velocity ( m/s)z(m) 0 10 20 30 40 50 60 0 0.01 0.02 0.03 Vol entration(%)z (m) umetric conc z1=0.11886 z2=0.00585 z3=0.00885 z4=0.01155 z5=0.01485 z6=0.01785 z7=0.02085 t=150s t=152s t=154s t=156s t=158s t=160s t=200s t=150s t=152s t=154s t=156s t=158s t=160s t=200s Figure 5.2 Growth of boundary laye r velocity and con centration profiles during a pretest run ith d = 0.21mm, u* = 0.06 m/s and = 0o. wIgnitive Turbidity curren t z u ( z (a) (b) x ) ) c ( z h Foreset slope onignitive turbidity cu rren N t Ignitive turbidity curren t Boundary layer 71 PAGE 72 Figure 5.3 Velocity and concentrat ion profiles: Experim ental data of Parker et al. (1987) and present model run. The main difference is in the grain size (0.03mm in experiment versus 0.20mm in the model). Both are at a distance of x =1.5m. igure 5.4 Velocity and concentrat ion profiles: Experimental data of Parker et al. (1987) and odel run. The main difference is in the grain size (0.03mm in experiment versus 0.20mm in the model). Both are at a distance of x =4.5m. 10 12 z (cm) F present m 0 5 10 15 20 25 0 2 4 6 8 u (cm/s) x =4.5m 0 0.05 0.1 0.15 0.2 0 2 4 6 8 10 12 Cv z (cm) x =4.5m 0 5 10 15 20 25 0 2 4 6 8 10 12 x =1.5m u (cm/s) z (cm) 0 0.05 0.1 0.15 0.2 0 2 4 6 8 10 12 x =1.5m Parker et al. (1987) Cv(cm) z Present m odel Parker et al. (1987) Present model Sediment bed datum Sediment bed datum 72 PAGE 73 Figure 5.5 Velocity and concentrat ion profiles: Experim ental data of Parker et al. (1987) and present model run. The main difference is in the grain size (0.03mm in experiment versus 0.20mm in the model). Both are at a distance of x =8.5m. 0 2 4 6 8 10 12 14 16 18 20 0 5 10 15 20 25 30 Slope in degreesShields Parameter N N N N N N N N N N N N N I I I I I N N N N N N N I I I I I I I I I 0.21mm 0.28mm 0.35mm Figure 5.ieldparater against bed slope. 6 S h s me =8.5m 0 5 10 15 20 25 0 2 4 6 8 10 12 u (cm/s) m) z (c x =8.5m x 12 0 0.05 0.1 0.15 0.2 0 2 4 6 8 10m) c z ( Cv Ignitive nignitive No Approximate boundary 73 PAGE 74 Slope in degreesShields Parameter 12 12.5 13 13.5 14 14.5 15 15.5 16 16.5 17 2 4 6 8 10 12 14 16 2 4 6 8 10 12 14 16 18 20 22 24 Figure 5.7 S nst d s e. Sedimnt flux pe r unit width (kg/m contours 3/s) hields param eter agai be lop e for 0.21mm diameter particles. 70 Slope in degreesShields Parameter 12 12.5 13 13.5 14 14.5 15 15.5 16 16.5 17 5 10 15 20 25 10 20 30 40 50 60 Figure 5.8 Shields parameter against bed slope. Sediment flux pe r unit width (kg/m3/s) contours for 0.28mm diameter particles. 74 PAGE 75 75 Slope in degreesShields Parameter 12 12.5 13 13.5 14 14.5 15 15.5 16 16.5 17 5 10 15 20 10 20 30 40 50 60 70 80 Figure 5.9 Shields parameter against bed slope. Sediment flux pe r unit width (kg/m3/s) contours for 0.35mm diameter particles. 102 2 101 100 101 10 102 101 100 101 Shields Parameter incipient movement Suspension 0.21mm 0.28mm 0.35mm Parker et al diameter in mm Figure 5.10 Shields parameter against particle diameter. Threshold curves for incipient movement and autosuspension curve. Model test run data. Bagnold autosuspension threshold Shields incipient movement threshold PAGE 76 CHAP TER 6 SUMMARY, CONCLUSIONS AND RECOMMENDATIONS 6.1 Summary Turbidity currents form a major mechanism for the transport of sediment from the continental shelf into deeper waters. Selfaccelerati ng turbidity currents is a topic of considerable scientific and engineering interest, because such currents can be selfsustaining and if so be responsible for significant transpor t of sediment derived ultimately from terrestrial sources to the seabed. Selfacceleration can occur if the flow velocity and the suspended sediment concentration increase simultaneously and entrain additional sediment due to the impelling downslope gravity force and associated turbulence. urrent without reference to its s ubsequent fate. Such a fate coul d involve excessive densification of the current as erosion continue s, collapse of turbulen ce, and possible destruction of the system as a stream. A condition meant to determine if a current is able to erode bottom sediment and form a suspension (as opposed to bedload) is the KnappBagnold criterio n for autosuspension. This is a necessary but not a sufficient criterion to ascertain if the ed will be ignitive or depositional, because an autosus pension may not satisfy the conditions for selfacceleration. A selfaccelerating turbid current has also been called a selfsustaining current. Most modeling efforts made to date to dete rmine the conditions fo r ignition are based on the use of singlephase flui spension. In the presen t analysis an available (Hsu 2002) onedimensional (vertical) model, which treats th e fluid and particles as separate phases, is The model investigates the effects of particle size, initial flow friction velocity and bed slope on the relationship between the Shields en trainment parameter and the grain size in the Ignition in the present study refers to the onse t of selfacceleration of an erosive turbidity c turbid current form ds trea ted as a su used. Experimental data of Park er et al. (1987) have been incl uded as part of the analysis. 76 PAGE 77 context of the KnappBagnold curv e for autosuspension. Particle sizes are lim ited to medium sands typical of beaches in Florida. The bed slop e is assumes to be small. Because the model is 1DV, flow velocities and suspended sedim a boundary layer velocity pr ofile is allowed to develop in th e time domain but considered to be locity at the end of the first phase. The bed in this case contains densely pa of the velocity and the concentration profiles are tracked in time. By using the depthmean initial Shields parameter against size to make an assessment of the relationship between the Knapp6.2 Conclusions The following are the main conclu sions reached in this study: 1. Although the twophase flow mode l of Hsu (2002) could only be run for a limited range of grain sizes (medium sands similar to thos e found in Florida beaches) and sediment bed density, comparison with flume experimental data of Parker et al. (1 987) indicated that the model produced reasonable results for th e velocity and concentration profiles in gravitydriven turbidity currents. ent concentration are depe ndent on the vertical coordinate and time but are independent of the ho rizontal coordinates. Accordingly, the model is used in two phases. At first, by applying a horizon tal pressure gradient us ing lateral ghost grids, in the (pseudo) horizontal dir ection. After the boundary layer velocity profile has developed close to a timeindependent struct ure the velocity profile is note d. In the second phase the bottom is tilted, which permits water to flow downslope starting with the ve cked but erodible sandy sediment. The evolutions velocity at the head of the slope time is converted into an approximate distance of travel down the slope. At two selected time s (or equivalent distances), th e depthmean velocity and the suspended sediment flux (per unit flow width) are calculated. These are used to calculate the corresponding flow acceleration and flux gradient, which are then used along with the plot of Bagnold autosuspension curve and the ignitive condition. 77 PAGE 78 2. Selfacceleration of turbidity current over a b ed of sand with the smallest of the three selected sizes (0.21mm) was found to depend on the initial conditions at the head of the slope determined by the imposed pressure grad ient in terms of the bed friction velocity. Bed slope seemed to be of secondary importance. 3. For the two sands with larger grain sizes (0.28mm and 0.35mm) the bed slope was found to play a more important role when compared to the initial pressure gradient. For a given pressure gradient, increasing the slope incr eased the likelihood of selfacceleration. 4. Based on conclusions 2 and 3 it could be concluded that ignition cannot be defined merely in terms of nontrivial, positive values of the velocity grad ient and the sediment flux gradient along the slope. De pending on particle size the in itial pressure gradient can also play a role. 5. For the selected initial conditions (grain size, pressure gradient and bed slope), out of the 54 combinations tested, all except three sati sfied the KnappBagnold criterion for autosuspension irrespective of whether the turbid ity current was ignitive or nonignitive. 6. In all 54 cases the current was found to be erosive. 78 PAGE 79 6.3 Recommendations for Further Work 1. Bearing th e mind that the tw ophase flow model has no re strictions relative to, or imposition of, the sediment pick up function or the need for the use of sediment bedload and suspended load equations, the models evid ent versatility should be demonstrated by reconfiguring it to accommodate a larger ra nge of the sediment sizes. Model results should be carefully tested against experimental data. 2. The model was used to simulate a physical situation in which bounda ry layer flows were generated on a horizontal bed for a certain time and then through numerical manipulation subsequent turbidity currents on inclined beds were simulated. It should be possible to switch the inclined surface after a prescribed time back to a horizontal surface. It should then be possible to witness energy loss thr ough a hydraulic jump, or eventual settling of the suspended sediment. Such a numerical e xperimentation should be verified against experimental data. 79 PAGE 80 APPENDIX A Implementation Scheme Numerical Open files for input and output Read input and problem setup the Setup the problem before computation If restarting, branch around the defaults Numerical p arameters Fluid parameters 1 2 3 Name list inpu t Initialize the constants Initialize the variables 80 PAGE 81 81 Set default array values Writing the run time info Reading the initial input data Calculating the pressure gradient For 1D boundar y la y e r Read sediment and bottom parameters Read the turbulence model Read output forma t Read restart controlling parameters Computing constant terms Minimum time step Generate the grid uniform an nonuniform d PAGE 82 Calculate for va mesh values riable Se te rm t the constant s for plotting Metric coefficients for the pressure solution Wr in ite the mesh formation Term e inat p ro g ra m Initial timestep based on molecular diffusion Initia region quantities lize the Se m consta tup the odel nts Init and ial condition for 82 PAGE 83 Set initial boundary c onditions Left, right and top boundary condition Output calculated data Specif co y boundary nditions Particle velocity using upwind scheme Prepare new c y cle Adjust time steps based on stability criterion Calculate se phase equations diment Calculate intera and large scale stress ction term Convection 83 PAGE 84 Cal p culate roduction Cal Ve culate locity Calculate tentative velocity Calculate fluid stresses Calculate eddy viscosity Advective term Viscous term Production term Strain tensors End 84 PAGE 85 LIST OF REFERE NCES Altinakar, M.S. (1988), Weakly depositing turb idity currents on small slopes, Ph.D. Thesis, cole Polytechnique Fdrale de Lausanne, Switzerland. Bagnold, R.A. (1962), Autosuspension of tr ansported sediment: Turbidity currents, Proceedings of the Royal Society of London, Seri es A, Mathematical and Physical Sciences, 265(1322), A discussion on progress an of Marine Science, 315319. Benjamin, T.B. (1968), Gravity currents and re lated phenomenon, Journal of Fluid Mechanics, 31, 209243. Britter, R.E. and J.E. Simpson (1978), Experiment s on the dynamics of th e gravity current head, Journal of Fluid Mechanics, 88, 223240. Chu, F.H., W.D. Pilkey, and O.H. Pilkey (1979), An analytical study of turbidity current steady flow, Marine Geology, 33, 205220. Daly, R.A (1936), Origin of submarine canyons, Ameican Journal of Science, Series 5, XXXI, 186, 401420. Drew, D.A. (1976), Produc tion and dissipation of energy in a turbulent flow of a particlefluid mixture, with some results on drag reduc tion, Journal of Applied Mechanics, 43, 543547. Eidsvik, K.J. and B. Brs (1989), Selfaccele rated turbidity current prediction based upon ( ) turbulence, Continental Shel f Research, 9(7), 617627. Elghobashi, S.E. and T.W. AbouArab (1983), A twoequation turbulence model for twophase flows, Physics of Fluid, 26, 931938. Fine, I.V., A.B.Rabinovich, B.D Bornhold, R.on, E.A. Kulikov (2005), The Grand Banks landslidegenerated tsunami of November 18, 1929: preliminary analysis and numerical modeling, Marine Geology, 215, 4557. Francis, J.R.D. (1957), London Engineer, 519. Fukushima, Y., G. Parker and H.M. Pantin (1985) Prediction of ignitive turbidity currents in Scripps submarine canyon, Marine Geology, 67, 5581. Greimann, B.P, M. Muste and F.M. Holly (1999), Two phase formulation of suspended sediment transport, Journal of Hydraulic Research, 37(4), 479500. eezen, B.C. and M. Ewing (1952), Turbidity currents and submarine slumps and the 1929 rand Banks earthquake, American Journal of Science, 250, 849873. d need r E Thoms H G 85 PAGE 86 Heezen, B.C., D.B. Ericson and er evid ence for turbidity current llowing the 1929 Grand Banks earthquake, Deep Sea Research. Part I, Oceanographic research su, T.J. (2002), A twophase flow approach fo r sediment transport, Ph.D. Thesis, Cornell ansport: dilute flow, urnal of Geophysical Research, 108, 3057 Jenkins and P.L.F. Liu (2004), On twophase sediment transport: sheet flow of assive particles, Proceedings of Royal Society of London, A, 22232250. llisional sheet flow of sediment driven by a turbulent uid, Journal of Fluid Mechanics, 370, 2952. entation, Cambridge University Press, New York. Union, 19, 501505. ty of Civil Engineers, 107, 659673. iddleton, G.V (1966b), Small scale models of turbidity currents and the criterion for auto6(1), 202208. ology, 315999. M. Ewing ( 1954), Furth fo papers, 193202. H University, Ithaca. Hsu, T.J., J.T. Jenkins and P.L.F. Liu (2003), On twophase sedi ment tr Jo Hsu, T.J., J.T m Jenkins, J.T. and D.M. Hanes (1998), Co fl Julien, P.Y., 1995. Erosion and sedim Knapp, R.T. (1938), Energy balance in streamf lows carrying suspended load, Transactions, American Geophysical Kuenen, P.H. (1937), Experiments in connecti on with Dalys hypothesis on the formation of submarine canyons, Leids geologische Mededel, 8, 327351 Lin, P.C.P. and A.J. Mehta ( 1997), A study of fine sedimentation in an elongated laboratory basin, Journal of Coastal Research, SI(25), 1930. McTigue, D.F. (1981), Mixture theory for susp ended sediment transport, Journal of the Hydraulics Division, American Socie Middleton, G.V. (1966a), Experiments on density a nd turbidity currents, Canadian Journal of Earth Sciences, 3, 523546. M suspension, Journal of Sedime ntary Petrology, 3 Nott, P.R. and J.F. Brady (1994), Pressuredriven flow of suspension: simulation and theory, Journal of Fluid Mechanics, 275, 157199. Pantin, H.M. (1979), Interaction between velocity and effective density in turbidity flow: Phaseplane analysis, with criteria for auto suspension, Marine Ge Parker, G. (1982), Conditions for ignition of ca tastrophically erosive turbidity currents, Marine Geology, 46, 307327. 86 PAGE 87 Parker, G., Y. Fukushim a and H.M. Pantin (1986), Selfaccelerating turbidity currents, Journal of fluid Mechanics, 171, 145181. Parker, G., M. Garcia, Y. Fukushima and W. Yu (1987), Experiments on turb idity currents over an lapp, J.E and J.P Mitchell (1960), A hydrodynamic th eory of turbidity currents, Journal of chmidt, W. (1911), Zur Mechanik der Bo en, Meteorologische Zeitschrift, 28, 355362. ty current, Journal of luid Mechanics, 53, 759768. ents on the dynamics of the front of the gravity urrent, Proceedings of the 2 International Symposium on Stratified Flows, The Norwegian ciety, 106, 485500. elft. iety, 46, 615683. Delft University of Technology, he Netherlands. erodible bed, Journal of Hydr aulic R esearch, 25(1), 123147. P Geophysical Research, 6, 983992. S Simpson, J.E. (1972), Effects of the lower boundary on the head of a gravi F Simpson, J.E. and R.E. Britter (1979), Experimndc Institute of Technology, Trondheim 174183. Simpson, J.E. and R.E. Britter (1980), A labor atory model of an atmospheric mesofront, Quarterly Journal of Royal Me teorological So Simpson, J.E. (1987), Gravity currents: In the environment and the laboratory, Ellis Horwood, England. Richardson, J.F. and W.N. Zaki (1954), Sedimentat ion and fluidization. Part 1, Transation of the Institution of Chemical Engineers, 32, 3553. Rodi, W. (1984), Turbulence models and their a pplication in hydraulics: A state of the art review, International Association of Hydr aulic Engineering and Research. D Von Karman, T. (1940), The engineer grappl es with nonlinear problems, Bulletin of the American Mathematical Soc Winterwerp, H. (1999), On the dynamics of highconcentrated mud suspensions, Communications on Hydraulic an d Geotechnical Engine ering, T 87 PAGE 88 88 Gowtham Krishna was born in 1982, Mysore, India. He completed his undergraduate gical Univer sity, Hyderabad, in f Florida in the Spring term of 2007. BIOGRAPHICAL SKETCH education in 2006 from the Jawaharlal Nehru Technolo Mechanical Engineering in 2006. After briefly work ing in the industry, he joined the University o 