Conditional Moment Closure for Methane Oxidation Using Two Conditional Variables and Stochastic Processes by Jorge R. Lozada-Ramirez B . A . S c , Universidad de las Americas-Puebla, Mexico, 1998 A THESIS S U B M I T T E D I N P A R T I A L F U L F I L M E N T OF T H E REQUIREMENTS FOR T H E DEGREE OF M A S T E R OF A P P L I E D SCIENCES in T h e Faculty of Graduate Studies (Department of Mechanical Engineering) We accept this thesis as conforming t o the required standard T H E U N I V E R S I T Y OF B R I T I S H C O L U M B I A June 16, 2004 © Jorge R. Lozada-Ramirez, 2004 Library Authorization In presenting this thesis in partial fulfillment o f the requirements for an advanced degree at the University o f British Columbia, I agree that t h e Library shall make it freely available for reference a n d study. I further agree that permission for extensive c o p y i n g o f this thesis for scholarly purposes m a y b e granted by the head o f m y d e p a r t m e n t or by his or her representatives. It is understood that copying or publication o f this thesis for financial gain shall not be allowed without my written permission. N a m e o f A u t h o r (please print) Title o f Thesis: CW.'-rViwf \)<k C\ck.M<?J OtorX Degree: Date (dd/mm/yyyy) /Wa^^-r- STT>cU«xj-Uy M<Lc Canada ^ Y e a r : 20 U**.lr<J T h e University o f British Columbia Vancouver, B C f (?, fi/[A Sr. Department of ClcS ore £ ruuJrtee.r/r>« i) . " OH Abstract ii Abstract 'Conditional Moment Closure for Methane Oxidation Using Two Conditional Variables and Stochastic Processes' b y Jorge R . Lozada-Rajiurez T h e c o n d i t i o n a l m o m e n t closure m e t h o d using t w o c o n d i t i o n i n g scalar variables is applied t o derive t h e t r a n s p o r t equation o f species mass f r a c t i o n , temperature, a n d scalar dissipat i o n i n a decaying, isotropic, homogeneous t u r b u l e n t methane-air flow. T h e s t r a i n tensor i n t h e t r a n s p o r t equation o f scalar dissipation of the c o n d i t i o n i n g variables is simulated using stochastic processes. T h e results o f this m o d e l are t h e n compared t o D N S a n d conditional moment closure w i t h one variable for the same test case. Contents iii Contents Abstract ii Contents iii List of Figures vi Nomenclature vii Acknowledgements x Introduction 1 1.1 Objective 1 1.2 D o u b l y C o n d i t i o n a l M o m e n t Closure ( D C M C ) 1 1.3 Stochastic Processes 2 1.4 Thesis O u t l i n e 2 1 2 Combustion, Turbulence, and Stochastic Processes 4 2.1 Premixed C o m b u s t i o n 4 2.2 Non-premixed Combustion 5 2.3 T r a n s p o r t Equations 5 2.3.1 T r a n s p o r t o f Mass, M o m e n t u m , a n d Energy 5 2.3.2 T r a n s p o r t o f M i x t u r e Fraction 6 2.3.3 T r a n s p o r t o f Species Mass Fraction 6 2.4 Chemical Kinetics a n d Reaction Mechanisms 7 2.5 Scales o f Turbulence a n d Kolmogorov Scales 9 2.6 Cascade o f T u r b u l e n t K i n e t i c Energy 10 2.7 Turbulence S i m u l a t i o n and M o d e l l i n g 10 Contents iv 2.7.1 Direct N u m e r i c a l S i m u l a t i o n (DNS) 11 2.7.2 Large-Eddy S i m u l a t i o n ( L E S ) 11 2.7.3 T i m e or Ensemble Averaging (Reynolds Decomposition) 11 2.7.4 Turbulence Models 12 2.8 Turbulent Combustion 13 2.9 Turbulent Combustion Simulation and Modelling 13 2.9.1 DNS 14 2.9.2 Favre Averaging 14 2.9.3 LES 14 2.9.4 L a m i n a r Flamelets 15 2.9.5 P r o b a b i l i t y Density Function ( P D F ) M o d e l l i n g 16 2.9.6 C o n d i t i o n a l M o m e n t Closure 16 2.10 Stochastic Processes a n d Monte C a r l o Methods 3 2.10.1 Monte C a r l o Simulations o f Turbulence 17 2.10.2 Monte C a r l o Simulations o f T u r b u l e n t Reacting F l o w s 19 2.11 S u m m a r y 21 The D C M C Method with Stochastic Processes 22 3.1 C o n d i t i o n a l M o m e n t Closure w i t h T w o C o n d i t i o n a l Variables ( D C M C ) 3.2 D e r i v a t i o n of t h e D C M C Equations 24 3.3 M o d e l l i n g Scalar Dissipation 26 3.4 Scalar T r a n s p o r t Closure Hypotheses 29 3.5 D C M C Closure Hypotheses 30 3.6 M a t h e m a t i c a l solution o f t h e Z-gradient squared t e r m 31 3.7 Models for the T r a n s p o r t E q u a t i o n o f Scalar D i s s i p a t i o n 32 3.7.1 Periodic Forcing w i t h R a n d o m Phase S h i f t i n g ( P F ) 33 3.7.2 Coupled M a p Lattices ( C M L ) 34 3.8 4 17 Summary . . . . 24 36 Simulation of Methane Oxidation with D C M C - S P 37 4.1 T h e Reference D N S and C M C Databases 37 4.1.1 37 Flow F i e l d Contents 4.2 4.3 v 4.1.2 Chemical Kinetics 38 4.1.3 I n i t i a l Conditions 39 S i m u l a t i o n o f CH 4 Oxidation Using D C M C - S P . . 40 4.2.1 Simulation o f t h e S t r a i n Field 42 4.2.2 Predictions o f Reactive Scalars 49 Summary 55 5 Conclusions and Recommendations 56 Bibliography 58 List of Figures vi List of Figures 2.1 T u r b u l e n t kinetic energy spectrum, f r o m Peters [1] 10 3.1 C o n d i t i o n a l average of t e m p e r a t u r e 23 3.2 D N S calculation o f temperature, f r o m Bushe [8] 23 4.1 DNS computational domain 38 4.2 D N S i n i t i a l velocities 4.3 I n i t i a l Conditions of the Scalar Fields 40 4.4 D C M C - S P code structure 42 4.5 Strain fields simulated w i t h P F 44 4.6 u-velocity fields simulated w i t h C M L 46 4.7 Velocity fields simulated w i t h C M L 47 4.8 Strain fields simulated w i t h C M L 48 4.9 D N S , C M C , a n d D C M C Favre averages, of species mass fractions 50 38 4.10 C o n d i t i o n a l averages at t = 30.0 51 4.11 a-variations 52 i n t h e intermediates field a t t = 30.0 4.12 a-variations i n t h e N O field at t = 30.0 53 4.13 a-variations i n t h e temperature field at t = 30.0 54 Nomenclature Nomenclature A Deterministic coefficient, amplitude of periodic force a Conditioning variable a B Deterministic coefficient C Constant c Speed of sound, concentration Cp Constant pressure specific heat D Molecular diffusivity, dimension dW Wiener process E Activation energy F Fuel, a stochastic process 9 Gravitational force h Specific enthalpy I Intermediates J Diffusion of energy, molecular diffusion k Arrhenius constant, reaction rate coefficient, energy level in coupled map lattice L Length I Length Le Lewis number m Mass O Oxidiser P Products P Pressure Pr Prandtl number R° Ideal gas constant Re Reynolds number vii Nomenclature sSc Surface stretch rate Sij Strain tensor t Time u Velocity V Volume V Velocity W Molecular mass X Mole fraction, stochastic variable X Coordinate axis Y Mass fraction, stochastic variable Z Mixture fraction Schmidt number viii Nomenclature Greek Letters a Value o f the c o n d i t i o n i n g variable a r Diffusion 7 Parameter i n t h e coupled map lattice, ratio o f specific heats Sij Kronecker d e l t a e Dissipation o f t u r b u l e n t kinetic energy £ Value o f the c o n d i t i o n i n g variable Z T) Characteristic length microscale 0 P r o d u c t o f t h e gradients o f scalar dissipation 0 R a n d o m phase shift K T u r b u l e n t kinetic energy A Parameter i n t h e coupled map l a t t i c e H Viscosity v K i n e m a t i c viscosity, stoichiometric coefficients II Product p Density S Summation T Characteristic t i m e microscale, parameter i n t h e coupled m a p l a t t i c e Tij Viscous tensor v Characteristic velocity microscale $ Gradient squared o f the c o n d i t i o n i n g variables X Scalar dissipation A scalar field, periodic force w i t h r a n d o m phase shifting CJ Chemical source t e r m w P e r i o d o f force ix Acknowledgements x Acknowledgements T h e w o r k presented i n this thesis was made possible thanks t o the s u p p o r t and encouragement of m a n y people. I would like t o express m y utmost appreciation t o m y supervisor, D r . K e n d a l Bushe, f o r his patience and timely, expert advice. Thanks are also due t o D r . A n d r e a Frisque for her constant encouragement and f o r those long discussions regarding stochastic processes. T h e financial s u p p o r t o f Mexico's N a t i o n a l C o u n c i l of Science and Tecnology (CONACyT) and o f Westport Innovations, I n c . are also gratefully acknowledged. T h a n k s are also given t o C o l i n B l a i r , A n d r e a , a n d Ray G r o u t f o r proof-reading preeliminary drafts o f t h i s thesis a n d f o r their always i n s i g h t f u l discussions on combustion, turbulence, m o d e l l i n g , p r o g r a m m i n g , a n d politics. T h a n k s also go t o C o l i n and J e n Bruce f o r t h e i r u n conditional s u p p o r t and for becoming m y f a m i l y away f r o m home. Special thanks are due t o M i n a H a l p e r n for her continuous support, encouragement, a n d m o t i v a t i o n i n reaching this goal. T h a n k y o u , M i n a , for all those smiles i n t h e times o f hardship. F i n a l l y , I w o u l d like t o dedicate t h i s thesis t o m y parents Jorge a n d Tere, a n d t o m y siblings Daniel a n d Raquel f o r always believing i n me. Papa, M a m a , D a n y y Raque; muchas gracias p o r su apoyo y p o r sus oraciones. Este documento esta dedicado a ustedes como u n pequeno t r i b u t o a su confianza y amor. Gracias p o r creer e n m i . Jorge R. Lozada-J^imirez June, 2004 i Chapter 1 Introduction T h e c o m p l e x i t y of t h e t r a n s p o r t equations t h a t describe t u r b u l e n t combustion makes necessary t h e use o f averaging procedures so t h a t practical problems can b e analyzed n u m e r i c a l l y i n a n efficient, t i m e l y fashion. T h i s averaging process creates a closure p r o b l e m by generating terms w h i c h are n o t explicit functions o f t h e averaged variables [2]. Generally speaking, closure problems are solved by m o d e l l i n g t h e unclosed terms t h r o u g h t h e use o f theoretical or empirical hypotheses or experimental results. I n t h e following chapters o f this thesis, a closure m e t h o d for t h e t r a n s p o r t equations o f t u r b u l e n t combustion w i l l be f o r m u l a t e d . T h e results o b t a i n e d from t h e c o m p u t a t i o n a l i m p l e m e n t a t i o n o f the resulting m o d e l o n a decaying, isotropic, homegeneous, non-premixed t u r b u l e n t flow o f methane a n d air w i l l be s h o w n a n d compared t o a reference D N S database o f a similar flow. 1.1 Objective T h e m a i n objective o f t h i s research project is t o explore t h e a p p l i c a t i o n of stochastic processes t o the m o d e l l i n g o f t h e strain tensor i n t h e scalar dissipation t r a n s p o r t equation t h a t is der i v e d using t h e D o u b l y C o n d i t i o n a l M o m e n t Closure ( D C M C ) m e t h o d . T h e resulting system o f equations is used t o simulate t h e o x i d a t i o n of methane by a i r i n decaying, homogeneous, isotropic turbulence. I t is expected t h a t t h e predictions of species mass fractions a n d temperature w i l l compare well w i t h those obtained from a direct numerical s i m u l a t i o n database a n d show a n improvement over t h e results provided b y a singly c o n d i t i o n a l moment closure s i m u l a t i o n o f the same test case. 1.2 Doubly Conditional Moment Closure (DCMC) T h e C o n d i t i o n a l M o m e n t Closure ( C M C ) o f K l i m e n k o [3] and B i l g e r [4] was developed w i t h the purpose o f f i n d i n g closure for t h e chemical source t e r m o f the species mass f r a c t i o n t r a n s p o r t Chapter 1. Introduction 2 e q u a t i o n . This m e t h o d is based o n t h e assumption that t h e fluctuations o f scalars, such as the chemical species mass f r a c t i o n a n d temperature, are related t o variations o f other scalars. I n t h i s thesis, t h e C M C m e t h o d has been expanded by using t w o different c o n d i t i o n i n g scalar variables. One o f those variables, m i x t u r e fraction, is used t o express t h e degree o f mixedness of t w o o r more elements i n a non-premixed, reacting flow. variable named a has been defined t o capture t h e effects o f A second non-reactive scalar fluctuations due t o s t r a i n along isosurfaces o f m i x t u r e fraction. T h r o u g h the use o f the d o u b l y c o n d i t i o n a l m o m e n t closure ( D C M C ) methodology, m a c r o - m i x i n g phenomena can be decoupled from t h e chemistry, w h i l e conserving the effects o f m i c r o - m i x i n g v i a the scalar dissipation t e r m . 1.3 Stochastic Processes T h e s t r o n g coupling between fluid dynamics a n d chemistry i n t u r b u l e n t reacting flows is frequently modelled using t h e scalar dissipation t e r m for w h i c h a n u m b e r o f different models have been proposed i n the past. I n t h i s thesis, a t r a n s p o r t e q u a t i o n is used t o calculate scalar dissipation. T h i s t r a n s p o r t equation includes t h e strain tensor w h i c h is often modelled i n t u r b u l e n c e research. Recognizing the chaotic behaviour of turbulence, t h e use o f t w o stochastic processes t o m o d e l t h e s t r a i n tensor i n t h e t r a n s p o r t equation o f scalar dissipation t e r m is proposed and analyzed. Stochastic processes are events t h a t depend o n stochastic variables t h a t e x h i b i t a n unpredictable behaviour. I n t h e first m o d e l , a periodic force w i t h r a n d o m phase shift ( P F ) is used t o emulate t h e chaotic behaviour o f t h e s t r a i n tensor. T h e second m o d e l uses the Coupled M a p L a t t i c e ( C M L ) o f Beck [5] a n d Hilgers a n d Beck [6], [7] t h a t generates realisations o f t h e fluctuating velocity increments as a f u n c t i o n o f t i m e f o r a f u l l y developed t u r b u l e n t flow. T h e strain tensor is t h e n calculated using those velocity increments. T h e results from b o t h codes are then compared t o t h e Direct N u m e r i c a l S i m u l a t i o n database of B u s h e et al [8] a n d the C M C results o f Bushe a n d Bilger [9]. 1.4 Thesis Outline I n t h i s chapter, t h e f u n d a m e n t a l grounds and m o t i v a t i o n f o r t h i s research project have been discussed. I n Chapter 2, t h e fundamentals o f combustion a n d turbulence, as well as t u r b u l e n t c o m - Chapter 1. Introduction 3 b u s t i o n , axe addressed. T h e t r a n s p o r t equations that describe t u r b u l e n t combustion, as well as current t u r b u l e n t combustion models, are also described. I n Chapter 3, the D C M C equations are derived. T h e closure assumptions and m o d e l l i n g o f t h e s t r a i n tensor a n d o f scalar dissipation are also discussed i n t h i s chapter. Chapter 4 pertains t o the i m p l e m e n t a t i o n o f the equations a n d models discussed i n Chapt e r 3 t o a decaying, homogeneous, isotropic t u r b u l e n t reacting flow. T h e turbulence m o d e l and chemical k i n e t i c mechanism used f o r this project are addressed a n d compared w i t h a direct numerical simulation. Finally, i n Chapter 5, the conclusions a b o u t the research are a r t i c u l a t e d and recommendations made f o r f u t u r e research. 4 Chapter 2 Combustion, Turbulence, and Stochastic Processes C o m b u s t i o n can be denned as the self-sustaining process t h r o u g h w h i c h two o r m o r e compounds, fuels and oxidizer, are modified t h r o u g h a chemical reaction w i t h a n associated release of energy. T h e s t u d y of combustion phenomena implies a c o m b i n a t i o n o f t h e r m o d y n a m i c s , fluid mechanics/dynamics, chemical kinetics, a n d transport phenomena. C o m b u s t i o n pro- cesses c a n be classified according to t h e m i x i n g behavior a n d t h e flow velocity regime of t h e system. Premixed systems a n d non-premixed systems can exist i n b o t h t h e l a m i n a r a n d t u r bulent regimes o f fluid flow. T h i s research utilizes m a t h e m a t i c a l m o d e l l i n g a n d computer s i m u l a t i o n o f combustion i n a deacying, isotropic, homogeneous, non-premixed t u r b u l e n t flow o f air a n d methane. T h e D o u b l y C o n d i t i o n a l M o m e n t Closure ( D C M C ) m e t h o d has been used t o solve t h e t r a n s p o r t equations o f mass fraction o f species. I t is proposed t o use a stochastic process t o model t h e strain tensor i n t h e t r a n s p o r t equation o f scalar dissipation responsible for the m i c r o - m i x i n g o f species. Different options f o r the stochastic process are discussed. 2.1 Premixed Combustion P r e m i x e d combustion occurs i f the fuel a n d oxidizer are completely m i x e d before t h e chemical reaction takes place. Examples o f p r e m i x e d combustion are Bunsen burners and s p a r k - i g n i t i o n engines. Premixed combustion is often characterized by t h e b u r n i n g velocity, w h i c h is a f u n c t i o n o f the chemical composition, i n i t i a l temperature, a n d pressure o f the p r e m i x e d flow. Chapter 2. Combustion, Turbulence, and Stochastic Processes 2.2 5 Non-premixed Combustion I n non-premixed combustion, t h e fuel and oxidizer are fed into a combustion chamber o r reactor f r o m different sources. P r a c t i c a l examples of non-premixed combustion include diesel engines, gas turbines a n d gas stoves. Unlike premixed combustion, where m i x i n g has already t a k e n place before the chemical reaction occurs, and the f u e l / a i r r a t i o is homogeneous i n t h e d o m a i n , i n non-premixed combustion fuel and oxidizer have t o m i x i n order f o r t h e reaction t o proceed. T h e state o f mixedness o f t h e reactants is characterized b y m i x t u r e f r a c t i o n , w h i c h is defined as: Here, Yp a n d Y o axe defined as t h e mass fractions o f t h e fuel a n d oxidizer, respectively. M i x t u r e f r a c t i o n for the fuel stream is, Z(xj,t) — 1, a n d Z(xi,t) = 0 for t h e oxidizer stream. 2.3 Transport Equations I n general, t h e t r a n s p o r t equations used i n combustion express t h e coupling between t h e fluid mechanical a n d chemical reaction phenomena t a k i n g place i n t h e combustion event. When dealing w i t h reactive m i x t u r e s o f gases, for example, t h e t r a n s p o r t equations are usually derived m a k i n g use o f t h e kinetic t h e o r y of gases i n order t o capture the specific behavior o f the system under study. Unless n o t e d otherwise, the following equations axe derived f o r t h e i r a p p l i c a t i o n t o problems i n v o l v i n g ideal gases. 2.3.1 Transport of Mass, Momentum, and Energy T h e t r a n s p o r t o f mass balances t h e amount o f mass entering a n d leaving a n y given c o n t r o l volume. dp dt dpuj _ dx{ Q T h e equation o f t r a n s p o r t o f m o m e n t u m indicates the change i n m o m e n t u m caused b y t h e e x t e r n a l forces (e.g. i n e r t i a l , pressure, viscous, a n d body forces) acting o n t h e fluid as, dt dxj dxj dxj Chapter 2. Combustion, Turbulence, and Stochastic Processes 6 where t h e viscous tensor for a N e w t o n i a n fluid, using the molecular viscosity, pi, a n d t h e Kronecker delta 8%j, is described by: (dui duj\ 2 du \ / k "3 The equation o f transport o f energy balances the t h e r m a l , chemical, a n d kinetic energy o f t h e system. Following Veynante and Vervisch [10], t h e t r a n s p o r t equation o f energy can be u expressed i n terms o f the t o t a l specific enthalpy, ht = h + -_-, o f t h e system as follows: 2 dpht dpuiht opuiht dp op dd dt + ~dx~ di = + ( dx- V \ ) rh i+ UjTij i + UiP9i where the diffusion o f energy, i n this case o f specific enthalpy, is described b y the Fourier Law, as: dh dxi i ~ p J r N ^fPr _—' \Sci \ dYj J dxi Here: 1. T h e P r a n d t l number, Pr = ^j^-, is t h e non-dimensional group t h a t measures t h e r a t i o of t r a n s p o r t o f m o m e n t u m due t o diffusion t o t h e temperature. C is the specific heat p at constant pressure a n d A is the t h e r m a l conductivity. 2 . Scj, t h e S d u n i d t number o f the P h diffusivity of the P h species, is defined as S c j = species relative t o a m a j o r o r reference species. 3. Yj is defined as the mass fraction o f t h e P h 2.3.2 where Dj is the species, as mi m Transport of Mixture Fraction M i x t u r e f r a c t i o n is a conserved scalar. I t is only used t o keep track o f t h e m i x i n g process as i t evolves i n t i m e a n d i n space. dpZ dt 2.3.3 | dpZ dxi u 1 d f dxi \ \ _o dxij d Z p D z ( 2 i ) Transport of Species Mass Fraction T h e t r a n s p o r t equation of the mass fractions for J = 1...N chemical species consists o f accum u l a t i o n , advection, diffusion, a n d p r o d u c t i o n terms: dpYj ^ + doYj dJ\ -dx- = -dx- - Ui +0}I M Chapter 2. Combustion, Turbulence, and Stochastic Processes 7 U s i n g Fick's L a w of diffusion [2], t h e molecular diffusion o f species is defined as: (2.3, J! = - D , ^ : A 2.4 Chemical Kinetics and Reaction Mechanisms T h e chemical source t e r m , CJJ, i n E q . (2.2) describes the n e t rate o f p r o d u c t i o n o f species due t o chemical reactions. T h e chemical source t e r m for M reactions is defined as [11], M /Y N K * / = w £ K * - 4 , * ) * n ( ^ ) ' » J=l K=l where Wj is defined as t h e molecular mass of t h e I th m e t r i c coefficients i n the K th w i t h Aj (-) 2 4 ' V species. i/ T K a n d v" K are t h e stoichio- chemical reaction, N N 7=1 7=1 being t h e chemical symbol o f the I species i n t h e chemical reaction, th k is t h e A r r h e n i u s t y p e constant o f t h e f o r m , (2.5) k = BT exp(-^y a Here, J3/fT corresponds t o the frequency factor for t h e K Q rif th reaction step, define t h e temperature dependence o f t h e frequency factor for the K TH represents the activation energy required by t h e K th OLK is used t o reaction, while EK reaction t o proceed, a n d R° being t h e x,= f'W universal gas constant. F i n a l l y , the mole fraction o f the P h species Xj is defined as, m w i t h t h e ideal-gas law, pV = mR°T, (2.7) and, N m = p^(y /Wi). r (2.8) 7=1 S u b s t i t u t i n g E q . (2.8) i n t o E q . (2.7) leads t o , N V = R TpY,(Yj/Wr). i=i 0 P (2.9) Chapter 2. Combustion, Turbulence, and Stochastic Processes 8 U s i n g Eqs. (2.5), (2.6), and (2.9), E q . (2.4) can now be w r i t t e n as */ = W, £ ) K « - B T«« K exp ( - J | ) JJ (^y^ • (2-10) T h e most basic analyses of combustion often rely on the assumption t h a t the chemical reactions occur at faster rates t h a n those associated t o the t r a n s p o r t phenomena o c c u r r i n g i n the flow. This assumption is often referred t o as 'fast chemistry' a n d usually denotes e q u i l i b r i u m states for chemical systems. B y invoking fast chemistry, systems can be analyzed b y using the basic laws of thermodynamics, which greatly simplifies the p r o b l e m i n h a n d . However, some chemical reactions occur at rates similar t o those o f t r a n s p o r t processes; most transient phenomena like i g n i t i o n , e x t i n c t i o n , and re-ignition cannot be predicted w i t h the use of fast chemistry. K i n e t i c mechanisms express t h e rates at w h i c h species are created or consumed d u r i n g a chemical reaction. For a system of reactions r = 1...R a n d species s = 1...S, s=l where k T s=l is the reaction rate coefficient, E q . (2.5) a n d the stoichiometric coefficients o f reac- t a n t s , frl\ a n d products uf). T h e rate of f o r m a t i o n for the I species as a f u n c t i o n o f the th concentration of S species, cs, as given by W a r n a t z [2]: *"=(#) V / chem,r =»{#-#)n<#s = \ Reaction rates are strongly non-linear functions of temperature v i a the reaction rate coefficient k. Pressure is also k n o w n t o affect the reaction rates by m o d i f y i n g t h e concentration o f species. Complex reaction mechanisms, e.g. for the o x i d a t i o n of methane, can consist o f m a n y different global reactions. Researchers have proposed methods t o reduce t h e n u m b e r o f reactions t o include o n l y those t h a t are deemed most relevant for t h e study o f the phenomena o f interest. Most r e d u c t i o n processes are based on assumptions of p a r t i a l e q u i h b r i u m a n d quasi-steady state, w h i c h implies t h a t reduced mechanisms are applicable o n l y under certain l i m i t i n g conditions. Peters [12], identified t h e r e d u c t i o n i n c o m p u t a t i o n t i m e for the n u m e r i cal s i m u l a t i o n of combustion events, as well as the applicability o f a s y m p t o t i c methods t o the analysis o f flame structures as t h e t w o most i m p o r t a n t applications of reduced kinetic mechanisms. A s an example, Bilger et al [13] proposed a m e t h o d t o derive a four-step mechanism Chapter 2. Combustion, Turbulence, and Stochastic Processes 9 for a methane-air flame w h i c h is i n good agreement w i t h t h e more complex mechanisms f r o m w h i c h i t originated for a l i m i t e d class o f combustion problems. T h e reduced chemical kinetic mechanism used i n this research project is described i n C h a p t e r 4. 2.5 Scales of Turbulence and Kolmogorov Scales T u r b u l e n t flows c a n be f o u n d i n different scales o f m a g n i t u d e ranging f r o m large scales such as those found i n oceans o r atmospherically relevant scales, t o the very small scales such as those studied i n combustion. Turbulence is a high-energy phenomenon. T h i s k i n e t i c energy is later dissipated, m a i n l y t h r o u g h viscosity. T h e effects o f molecular viscosity, however, are relevant only at very small scales. T h i s implies t h a t , i n order for the energy t o b e dissipated at those scales, t h i s energy must be transferred f r o m the largest scales d o w n t o t h e smaller scales. T h e vortices, also k n o w n as eddies, present i n a t u r b u l e n t flow v a r y i n scale a n d t h e y can b e correlated t o the development o f velocity gradients. I t has been f o u n d t h a t t h e largest vortices are responsible f o r most o f t h e transport o f m o m e n t u m . However, viscosity acting a t molecular scales is responsible for t h e dissipation o f energy t h r o u g h a conversion o f k i n e t i c energy i n t o heat a t the s m a l l scales o f m o t i o n [14]. Since chemical reactions t y p i c a l l y take place a t this same molecular level, i t is i m p o r t a n t t o consider these s m a l l turbulence scales. T h e K 4 1 theory, proposed b y Kolmogorov [15], states t h a t , at the smallest scales o f t u r b u lence, t h e fluid m o t i o n is self-similar a n d isotropic, giving o r i g i n t o the concept o f universality of turbulence at such small scales. Following Kolmogorov, t h e use of t h e kinematic viscosity o f the fluid, v, a n d t h e energy dissipation rate normalised per u n i t mass, e, p e r m i t s t h e d e f i n i t i o n of t h e following microscales, also k n o w n as the Kolmogorov scales. • Characteristic length microscale, • Characteristic time microscale, • Characteristic velocity microscale, Chapter 2. Combustion, Turbulence, and Stochastic Processes 10 A t these micro-scales, t h e resulting local Kolmogorov Reynolds n u m b e r is equal t o unity, w h i c h means t h a t there is no turbulence at smaller scales. 2.6 Cascade of Turbulent Kinetic Energy Based on Kolmogorov's ideas, researchers have developed the concept of t h e energy cascade, w h i c h states t h a t t u r b u l e n t kinetic energy transfers f r o m large t o small scales v i a vortex stretching a n d once t h e energy reaches the smaller scales of turbulence, viscous dissipation occurs. Tennekes a n d L u m l e y [14] provide a complete discussion o n t h e vortex-stretching mechanism. F i g . 2.1 ( f r o m Peters [1]) shows the t u r b u l e n t k i n e t i c energy s p e c t r u m as a f u n c t i o n o f wave number k, where k is the inverse of t h e vortex o r eddy size. F r o m t h i s figure, i t can be seen t h a t t h e kinetic energy exhibits a monotonic decay i n t h e l o g a r i t h m i c p l o t , w i t h a slope of —5/3 along t h e i n e r t i a l subrange. T h i s is followed by a n a b r u p t decay i n t h e viscous subrange, where viscous dissipation affects t h e smallest scales o f turbulence. large scales log E(k) I 0 energy containing integral inertial subrange 1 r viscous subrange 1 1 7,-' log k Figure 2 . 1 : T u r b u l e n t kinetic energy s p e c t r u m , from Peters [1] 2.7 Turbulence Simulation and Modelling T h e Navier-Stokes equations t h a t describe turbulence can, i n p r i n c i p l e , be solved d i r e c t l y using numerical techniques. I n reality, t h e strong coupling a n d non-linearity o f the equations, along w i t h the w i d e range o f length a n d t i m e scales involved i n the s o l u t i o n a n d t h e chaotic n a t u r e o f turbulence, make such direct calculations possible o n l y i n very l i m i t e d cases. T i m e a n d mass averaging techniques have been developed t o circumvent this p r o b l e m . T h e i n t r o d u c t i o n o f such techniques, however, produces t e r m s for w h i c h the solution is either n o t k n o w n or Chapter 2. Combustion, Turbulence, and Stochastic Processes 11 cannot be expressed as an explicit function of the averaged variables. Oftentimes, these terms have to be modelled. The following sections describe the simulation and modelling techniques most widely used i n turbulence. 2.7.1 Direct Numerical Simulation (DNS) In direct numerical simulation (DNS), the Navier-Stokes equations are solved numerically. This requires very fine spatial and temporal discretisations so that the transport equations can be solved from the large scales down to the Kolmogorov microscales. In general, it is estimated that the computation time for a DNS simulation increases approximately at a rate of Re* [2]. Consequently, at present only low Reynolds number flows can be solved directly. This limitation prevents the use of D N S for most practical applications. Currently, most research using D N S is aimed at the analysis of small turbulent structures, and the development and validation of closure models. 2.7.2 Large-Eddy Simulation (LES) W i t h Large-Eddy Simulation (LES), the large scales of turbulence are calculated directly in the same way as i n DNS. The smaller scales are modelled using turbulence models and through processes of numerical filtering. L E S provides high-quality results at a much lower computational expense than D N S . 2.7.3 Time or Ensemble Averaging (Reynolds Decomposition) The Reynolds decomposition [14], is frequently used to obtain statistical information from a turbulent flow. The time-average velocity of the flow can be defined as, which can also be expressed as: Ui{x,t) N so that instantaneous fluctuations of velocity around the time-average mean can be expressed u'iixjt) = Ui{x,t) — Ui(x). Chapter 2. Combustion, Turbulence, and Stochastic Processes 12 A n i m p o r t a n t , a n d w i d e l y used parameter i n turbulence characterization is t h e t u r b u l e n t k i n e t i c energy K , which is obtained f r o m t a k i n g the root mean square and time-averaging t h e velocity fluctuations as follows [16], - 1 2 — w h i c h , for a three-dimensional flow field becomes, '2 . _ _ Uj + K ~ *2+ M'2 i 3 2 • I f n o a d d i t i o n a l energy is transferred t o the system, the kinetic energy dissipates i n t i m e a t a r a t e given by, DK m- e= 2.7A Turbulence Models F r o m p r o b a b i l i t y theory, i t can be shown t h a t t h e ensemble average o f the velocity fluctua- t i o n s around t h e mean value is equal t o zero [14], or i n m a t h e m a t i c a l f o r m (u[) = 0. T h e second moment (u^u'^j t h a t represents the statistical correlation between u,- a n d Uj becomes a n e x t r a u n k n o w n i n the system o f equations f o r which a n equation as an e x p l i c i t f u n c t i o n o f t h e averaged variables is n o t available. T h i s is commonly k n o w n as a closure p r o b l e m a n d the unclosed terms are often modelled using physical or m a t h e m a t i c a l tools. Launder a n d Spalding [17], a n d Pope [18] summarise different m a t h e m a t i c a l turbulence models. Those models range f r o m one-equation t o multi-equation m a t h e m a t i c a l expressions t o characterise t u r b u l e n t phenomena a n d were derived as a n alternative t o previous models such as t h e m i x i n g - l e n g t h model. T h e numerical designation for the models is based o n the n u m b e r of p a r t i a l differential equations w h i c h are needed t o calculate the modelled parameters. G a t s k i a n d Rumsey [19], for example, identified variations of previous models t h a t have originated 'zero-equation' models, w h i c h are algebraic relations t o define the v o r t i c a l o r eddy viscosity, a n d 'half-equation' models where a n ordinary differential equation is solved t o calculate t h e eddy viscosity. Results obtained f r o m zero- t o m u l t i - e q u a t i o n models may vary considerably, due t o the different degrees o f complexity a n d detail t h a t they may include i n their analyses. M o s t models are derived f o r specific applications, e.g. flow around a cylinder, boundary-layer flows, and so f o r t h and t h e y are usually t u n e d t o provide acceptable t o good a p p r o x i m a t i o n s t o experimental results. Chapter 2. Combustion, Turbulence, and Stochastic Processes 2.8 13 Turbulent Combustion Most practical applications involving combustion occur in turbulent regimes. Directly or indirectly, a l l of the characteristics of turbulence have a major effect on the evolution of a combustion event. Vorticity, for example, increases the interface area between fuel and oxidiser in non-premixed combustion, or between reactants, intermediates, and products i n premixed flames, which results i n an increased rate of reaction. In premixed combustion, turbulence has several important effects. It increases the area of the reactant-intermediate-products interface, which i n turn increases the rate of reaction. This increase i n area is due to the appearance of 'wrinkles' in the flame front which are produced by the action of turbulent vortices. As a consequence of this stretching, the reaction zone becomes thinner, promoting higher flame-propagation velocities. If the accelerated mixing rate produced by the stretch effect reaches a point where mixing takes place faster than reaction, extinction of the flame occurs [20]. In non-premixed combustion, the most significant role of turbulence is to enhance mixing. This is achieved through the increase in the surface area of the fuel-oxidiser interface due to the strain and stretch created by vorticity and the high diffusion velocity present i n turbulent flows which is a result of an increase in gradients. Since non-premixed combustion is a primarily mixing-limited phenomenon, the result of these effects is of paramount importance. 2.9 Turbulent Combustion Simulation and Modelling There are a number of different techniques and approaches to solve the equations that describe the evolution of a turbulent reacting flow. Well established methods have been used for the study and analysis of turbulence. These methods have been modified and adapted to account for the interactions between fluid mechanics/dynamics and chemistry. Basically, most of the research work is dedicated to closure methods for the chemical source term. The works by Mell et al [21], Swaminathan and Bilger [22], Bilger [23], and Peters [1], amongst others, offer a variety of options for the modelling of the chemical source term. A brief description of the most widely used techniques and their applications is given below. Chapter 2. Combustion, Turbulence, and Stochastic Processes 14 2.9.1 D N S E c h e k k i et al [24] report t h e use of D N S , coupled w i t h a multi-step chemical kinetics mechanism for t h e s i m u l a t i o n o f a 2-D p r e m i x e d t u r b u l e n t reacting flow. Papers b y Bushe et a l [8] and S w a m i n a t h a n a n d Bilger [22] describe t h e use o f D N S w i t h reduced chemical k i n e t i c mechanisms i n t h e study o f non-premixed combustion. T h e D N S database of t h e o x i d a t i o n o f methane i n a t u r b u l e n t flow of Bushe et al [8] is used as t h e baseline f o r the comparison w i t h the results o f t h e model developed i n t h i s thesis. 2.9.2 Favre Averaging C o m b u s t i o n processes are characterized by considerable density fluctuations, w h i c h has p r o m p t e d the development o f the density-weighted, or Favre averaging. T h e Favre average o f a property, such as mass f r a c t i o n for example, c a n be defined as y = ^ P (2.11) T h e instantaneous local value o f mass fraction can be defined i n terms o f i t s average a n d fluctuation as (2.12) Y{x,t)=Y{x,t)+Y'{x,t) Similarly, the instantaneous local value o f mass f r a c t i o n is defined i n t e r m s of its Favre average a n d fluctuation as (2.13) Y(x,t)=Y(x,t)+Y"(x,t) C o m b i n i n g the d e f i n i t i o n o f E q . (2.11) w i t h E q . (2.12) leads t o a n expression f o r the Favre average o f mass f r a c t i o n i n terms o f t h e average mass f r a c t i o n as ~ Y (p + d) (7 + Y') = \PTl>)\ —I — = P Y + pfY 1 P^_ , 2 M ) P T h e t e r m p'Y' i n E q . (2.14) denotes t h e correlation between the fluctuations o f density a n d mass f r a c t i o n . T h i s is a n unclosed t e r m which can be computed using a t r a n s p o r t equation, or modelled using theoretical or empirical relations. 2.9.3 L E S C h e m i c a l reactions occur at small scales where t h e flow is modelled i n L E S , n o t resolved; t h i s implies t h e need f o r models for the chemical source t e r m i n the t r a n s p o r t equations of mass Chapter 2. Combustion, Turbulence, and Stochastic Processes 15 f r a c t i o n o f species [16]. LES results have been used as benchmarks for a priori t e s t i n g of o t h e r simulations. T h e piloted flame of Pitsch and Steiner [25], as well as t h a t of Steiner and Bushe [26], have been used i n comparative studies o f computer simulations. I n t h e case of t h e work by Steiner a n d Bushe, L E S was incorporated w i t h the newly-developed c o n d i t i o n a l source-term estimation. Laurence [27] discussed applications of L E S i n i n d u s t r i a l settings, and identified t h e ongoing work on numerical methods for finite volumes i n u n s t r u c t u r e d grids as one of the m a i n activities aimed at p r o v i d i n g more o p p o r t u n i t i e s t o such applications o f the method. 2.9.4 Laminar Flamelets F l o w visualisation of reacting flows w i t h low-intensity turbulence shows t h a t curved flame f r o n t s propagate outwards from the source o f i g n i t i o n t o the boundaries o f the system. In t h i s fashion, the t u r b u l e n t flame can be conceptualised as an aggregation o f such flame fronts. T h i s is also k n o w n as t h e flamelet regime. Based on t h i s concept, the L a m i n a r Flamelets, as described b y Peters [28], [1], are t h i n (their w i d t h being smaller t h a n the K o l m o g o r o v lengthscale), reactive sheets t h a t react a n d diffuse w i t h i n a non-reactive field. T h i s m e t h o d assumes t h a t once i g n i t i o n has occurred, t h e chemical time scale is greatly accelerated, reducing the reaction zone i n t o a t h i n sheet. T h i s concept was f u r t h e r developed into assuming t h a t i f the t h i n reaction zone is smaller t h a n the Kolmogorov length-scale, t h e n the flow i n t h a t zone is l a m i n a r , due t o the fact that there exists no turbulence below t h e K o l m o g o r o v microscales. M u c h research has been dedicated t o t h e development of t h e flamelet m o d e l . Peters [1] p r o v i d e d a complete description o f the m o d e l , as well as its f o r m a l a p p l i c a t i o n t o p r e m i x e d a n d non-premixed t u r b u l e n t combustion. L e n t i n i [29] explored t h e a p p l i c a t i o n o f the stretched l a m i n a r flamelet m e t h o d combined w i t h t h e K — e m o d e l for turbulence i n a t u r b u l e n t nonp r e m i x e d flow a n d reported good agreement w i t h results obtained f r o m experiments p e r f o r m e d at the flamelet regime. Cook et al [30] presented the large-eddy laminar flamelet m o d e l i n w h i c h the chemistry at the smaller scales was modelled using the flamelet theory. T h e flamelet m o d e l offers good predictions o f t h e m a j o r chemical species t h a t originate i n combustion of hydrocarbon-based fuels at w h a t is k n o w n as the 'flamelet regime'. However, t h e flamelet theo r y applies under very specific conditions of t h e flow, e.g. the flamelet thickness is smaller t h a n Chapter 2. Combustion, Turbulence, and Stochastic Processes 16 the smaller length-scales of turbulence, which limits the effective application of the flamelet method. Swaminathan and Bilger [22] pointed out that the flamelet model falls short of providing good predictions of the minor species in combustion of hydrocarbon fuels. This point was also discussed by Mell et al [21] when comparing the flamelet and the conditional moment closure ( C M C ) models with a DNS database for different flow configurations to study the effects of variations of scalar dissipation on flamelet and C M C results. A n improvement over the laminar flamelet model is the unsteady laminar flamelet theory, which introduces a Langrangian time frame that records the time dependence i n the flame structure. Mauss et al [31] presented results that seem to indicate that the transient variations of the flamelet structure of a steady turbulent non-premixed flame are responsible for the scattering of experimental data points obtained with one-point Raman measurements. Pitsch [32] obtained good agreement between experiments and the results of the unsteady flamelet model for a turbulent non-premixed CH^ffyf^-air 2.9.5 flame. Probability Density Function (PDF) Modelling In the P D F method, a transport equation for the joint P D F of relevant scalars is solved to model the evolution of a combustion event. These scalars can include mass fraction of chemical species, flow velocities, temperature, mixture fraction, and so forth. A n important advantage of this method is that the chemical source term is solved for directly, without the need of a model. Modelling is necessary, however, for the molecular mixing terms i n scalar space, as well as for the turbulent transport terms in physical space [33]. Due to the high-dimensionality inherent to this method, traditional computing techniques are insufficient. This has prompted researchers to apply numerical methods capable of dealing with problems with a large number of variables and dimensions, such as Monte Carlo methods. These applications will be reviewed with more detail in the following sections of this thesis. 2.9.6 Conditional Moment Closure The Conditional Moment Closure ( C M C ) of Klimenko [3] and Bilger [4] determines closure for the chemical source term of the transport equation of the mass fraction of species. This method is based on the assumption that the fluctuations of scalars, such as mass fraction Chapter 2. Combustion, Turbulence, and Stochastic Processes 17 of chemical species and temperature, are associated with variations of other scalars. The advantage of using the C M C method is that macromixing phenomena can be decoupled from the chemical kinetics, while conserving the effects of micromixing via the scalar dissipation term. The method is described and its equations derived in Chapter 3. 2.10 Stochastic Processes and Monte Carlo Methods The theory behind stochastic processes originated from the observations of the trajectories of pollen particles moving on a water surface made by Brown in 1827 and later independently developed into an elegant mathematical theory by Einstein [34], [35], Smoluchowski [36], and Langevin [37]. The experiments by Perrin [38], [39], [40], [41], [42] confirmed the previous theoretical analyses and his work was awarded the Nobel Prize of Physics i n 1926. The initial applications of the method concentrated on physics, mathematics and chemistry. Today, stochastic processes are extensively applied i n fields such as biology, finance, economics, demography, atmospheric and ocean sciences, amongst many others. Different numerical methods, such as the Monte Carlo ( M C ) method, have been developed as tools for the study of stochastic processes. The M C method is intended to solve multi-dimensional problems with a large number of variables. This can be a very inefficient method when used to solve simple problems given the large number of computations required. 2.10.1 Monte Carlo Simulations of Turbulence Kramer [43] reviewed the application of the M C method i n generic turbulence problems. In this work, a thorough description of the two different (Eulerian and Lagrangian) types of M C simulation is presented. The M C method is remarkably useful for the simulation of turbulent flows given its ability to treat multi-dimensional problems with a large number of degrees of freedom. Another attractive feature of this method is its ability to generate chaotic particle trajectories and velocity field structures. In the Eulerian approach, a stochastic velocity field model replaces the turbulent velocity field to be simulated. This field is then solved i n a spatial domain, which is computationally less demanding when compared to the solution of the Navier-Stokes partial differential equations. The stochastic field is tailored to the specific application and is frequently used to simulate trajectories of particles suspended in a flow. Chapter 2. Combustion, Turbulence, and Stochastic Processes 18 A number of methods exist to numerically solve the equations derived using the Eulerian approach. The Fourier method is based on a spectral formulation for a homogeneous Gaussian stochastic field. A major drawback of this method is the generation of a spatial period that has detrimental effects on turbulent diffusion simulation. Using a randomization method, the stochastic Fourier integral is discretised by stochastic wave numbers. Despite the fact that this method provides useful data to simulate velocity fields with strong long-range correlations, those fields have non-Gaussian statistics. According to Elliot et al [44] and Sabelfeld [45] this limits the applications of this method for ideal mathematical or physical turbulence models. More physically meaningful models have also been developed. Among them are the moving-average and the wavelet methods. McCoy [46] introduced a Gaussian stochastic field represented in physical space using white noise convolved against a kernel. Long-range correlations i n the velocity field can be appropriately simulated using the wavelet method. There exists a combination of Fourier-wavelet methods that present improvements over its precursors. This method has been reported to have shown excellent results when simulating velocity fields with the condition of self-similarity scaling of their structures. A Lagrangian M C simulation utilises particles evolving i n a stochastic model. A l l the effects of turbulence are incorporated within the properties of each particle. Despite the fact that the Lagrangian M C simulation does not make reference to a fluid velocity field, Minier and Pozorsky [47] and Welton [48] have used it to simulate the motion of a turbulent flow by analysing the stochastic evolution of a group of particles. The average properties of this group of particles (and of the ones in the vicinity of the target field) define the properties of the fluid, and of the flow, at that specific location. The methods developed to deal with M C simulations using the Lagrangian approach include the simulation of fluid motion and particles, the simulation of immersed bodies and hybrid L E S / M C schemes. In the simulation of fluid motion and particles, a set of fluid particles is followed according to a stochastic formulation. Physical, or pseudo-physical properties are defined for the particles. The evolution of the dynamics of the flow is assumed to be the average of the simulated particles. This method offers flexibility i n the way that it can provide some insight into the mixing process of the unresolved scales of the turbulent flow by the use of a Chapter 2. Combustion, Turbulence, and Stochastic Processes 19 simple Langevin formulation. If more details are desirable, the stochastic formulation must be changed. Stochastic differential equations i n the form of Markov processes are often used with this purpose in mind [49]. The simulation of immersed bodies follows a similar procedure. Pozorsky and Minier [50] reported an approach i n which two sets of stochastic particles are defined. The first set represents the immersed particles, while the second is intended to simulate the flow. The physical and pseudo-physical properties of both sets of particles are different from one another. A major drawback of these methods is that they can only simulate homogeneous turbulent flows. In the hybrid L E S / M C schemes the flow field is resolved using a L E S formulation, while the mixing (and chemical reaction, if reacting flow is being considered) simulations are resolved using a M C formulation. This scheme calls for the continuous feedback from one formulation to the other. Jaberi et al [51] and Obliego et al [52] have reported an increase in accuracy in their L E S calculations after the addition of the M C component. The apparent reason for this is that the M C simulation provides some information on the effects of the small scales otherwise neglected by the L E S model. This model has also been used successfully i n nonhomogeneous turbulent flows, which extends its potential for applicability to more realistic problems. Le Maitre et al [53] proposed a method in which a stochastic spectral finite element method undertakes the formalism of a M C simulation and offered numerous results claiming improvements i n computational expenses, amongst other advantages over the M C method. Unfortunately, this research group d i d not provide a comparison between the results of the two formulations. 2.10.2 Monte Carlo Simulations of Turbulent Reacting Flows Most applications of the M C method i n turbulent reacting flows follow one of the Lagrangian formulations described in the preceding section. Pope [54] described what appears to be the first application of the M C simulation to a turbulent reacting flow. This method has shown great usefulness given its capabilities of dealing with problems with a large number of variables and dimensions. Its first and foremost application is intended to numerically solve the P D F transport equations derived for a trubulent reacting flow. Chapter 2. Combustion, Turbulence, and Stochastic Processes 20 T h e advantages o f the m e t h o d are numerous. A s described i n different papers [54], [55], t h e P D F m e t h o d provides closure t o the chemical source t e r m i n the equation w i t h o u t the need of any modelling and can be applied t o either p r e m i x e d or n o n - p r e m i x e d flows. I t is also k n o w n t h a t this m e t h o d incurs large c o m p u t a t i o n a l demands, has no a b i l i t y t o model wall-type boundaries, and is difficult t o a p p l y i n complex flows. According t o Pope [54], the c o m p u t a t i o n a l requirements o f this m e t h o d increase linearly w i t h the n u m b e r of variables used (in this case dimensions, as w e l l ) , compared t o t h e exponential increase o f computations required by traditional finite-difference formulations. T h i s fact provides a way of dealing w i t h problems i n v o l v i n g m a n y species. As discussed i n the previous section, t h e M C m e t h o d can deal w i t h t h e s i m u l a t i o n of walls and complex, inhomogeneous flows, hence solving the more relevant drawbacks o f the P D F m e t h o d . As described by Pope [54], stochastic particles were used t o m i m i c t h e chemical and t h e r m o d y n a m i c properties of the flow. These particles, however, were defined as a numerical t o o l a n d have no physical meaning. A plug-flow reactor p r o b l e m w i t h imperfect m i x i n g was simulated using the M C m e t h o d , using a single scalar t o m o d e l the concentration of reaction products. Pope's numerical results were i n good agreement t o those obtained i n experiments performed by B a t t [56]. Valino [57] proposed an E u l e r i a n M C f o r m u l a t i o n for the s i m u l a t i o n of t h e P D F o f a single scalar i n a t u r b u l e n t flow. His proposal was t o use stochastic particles ' j u m p i n g ' f r o m node t o node i n t h e spatial d o m a i n according t o specific rules (a stochastic equation). Hulek a n d L i n d stedt [58] offered a comparison between d a t a obtained from D N S and f r o m M C simulations. T h e i r f o r m u l a t i o n used the M C m e t h o d t o solve the terms modelled w i t h t h e P D F m e t h o d i n t h e flamelet equation. T h i s m e t h o d , called flamelet acceleration, was also tested a n d i t appeared t o be an improvement w h e n compared w i t h t h e standard flamelet s o l u t i o n . H u l e k and L i n d s t e d t [59] later modelled p r e m i x e d t u r b u l e n t flames using P D F and M C methods. T h e M C s i m u l a t i o n was used to solve the mass density f u n c t i o n evolution equation. Stochastic particles evolved following deterministic processes representing the closed terms of the equation a n d stochastic processes designed t o model the effect o f the unclosed terms. T h e y showed the comparison between several different numerical formulations a n d experimental results, where i t can be seen that t h e P D F - M C simulations were i n good agreement w i t h experiments. Chapter 2. Combustion, Turbulence, and Stochastic Processes 21 The special ability of the M C method to mimic turbulent mixing was explored by Kawanabe et al [60]. The methodology was applied to a plug-flow flame i n unsteady turbulent combustion. The turbulence model used was a standard K — e and the mixture is forced to avoid homogeneity. Cannon et al [61] and Kraft and Fey [62] have explored stochastic reactor models using the P D F and M C methods. Cannon's group offered some practical results, while Kraft and Fey mainly developed a purely analytical solution to the stochastic reactor model. Cannon's research was aimed at the prediction of C O and N O in premixed combustion of methane, using a partially stirred reactor (PaSR). Conceptually, the P a S R could be thought of as a M C process itself, since it also randomly selects particles that are then forced to interact amongst themselves. This method is an improvement over the perfectly stirred reactor ( P S R ) since finite-rate mixing effects can be included i n the solution. The mixing term was modelled using a deterministic method rather than a M C method. The good agreement between the 9-step simulation and the detailed kinetic solution demonstrated the robustness and accuracy of this particular model. 2.11 Summary In this Chapter, the fundamentals of turbulent combustion, have been reviewed. Turbulence, scales of turbulence, and turbulence simulation and modelling were also treated, along with the interactions between turbulence and combustion. Descriptions of the regimes and equations that define turbulent combustion, as well as some of the most important analytical tools available for the solution of those equations were also presented in this chapter. A detailed description of the applications of the Monte Carlo simulations in turbulence and turbulent combustion was provided. It has been observed that there exist a number of different options i n which the M C method can be used to solve numerically demanding problems i n turbulent combustion. Future possible applications include particulate matter and pollutant formation mechanisms and mass transport in permeable membranes. O f particular interest to this research is the capability of the method to be included to solve problems expressed in terms of different modelling theories such as conditional moment closure, L E S , or flamelet models. A n increasing number of research opportunities exist in this field, with promises of continuous and prolific growth in the near future. 22 Chapter 3 The D C M C Method with Stochastic Processes C o n d i t i o n a l moments are averages which are calculated f o r those members o r realisations t h a t satisfy some predefined c o n d i t i o n . T h e conditional moment closure m e t h o d ( C M C ) , inde- pendently developed b y K l i m e n k o [3] and Bilger [4], provides closure for t h e chemical source t e r m i n t h e t r a n s p o r t equation o f species a n d enthalpy o r temperature b y assuming t h a t the v a r i a b i l i t y i n temperature a n d t h e mass fractions of species can b e linked t o t h e fluctuations o f a given scalar variable. T h e reaction progress variable and m i x t u r e f r a c t i o n are c o m m o n choices o f scalar variables for p r e m i x e d a n d non-premixed combustion analysis, respectively. I f t h e mass fractions o f species are expressed as functions of n o t only space a n d t i m e , b u t also as functions of t h e c o n d i t i o n i n g variable, their conditional averages c a n then be used t o a p p r o x i m a t e t h e conditional average of t h e source t e r m as d\Z = £ « u \T\Z = £, Yj\Z = £ j . T h e overbar i n this n o t a t i o n indicates t h a t t h e mean value has been calculated f r o m averaging over n = 1...N realisations. T h e vertical b a r means 'given', and i n t h e case o f Yj denotes the value of t h e mass f r a c t i o n given t h a t m i x t u r e fraction has a value o f £. As an i l l u s t r a t i o n o f t h e theoretical basis of the C M C m e t h o d , F i g . 3.1 shows a plot o f i n d i v i d u a l n = 1..JV realisations of a simulation o f temperature (blue fines). W h e r e the u n c o n d i t i o n a l average o f temperature, T = YliLiTi/N, is represented b y t h e h o r i z o n t a l line at around T « 500, a n d t h e red curve is t h e calculated conditional average o f t e m p e r a t u r e , T = (T\Z = Q, where t h e angle brackets indicate t h a t only those realisations t h a t comply w i t h the c o n d i t i o n o f Z = £ are considered i n the averaging process. I t is clear from F i g . 3.1 t h a t the i n f o r m a t i o n o b t a i n e d from conditional averages is more representative of t h e actual relationship between T a n d Z t h a n that obtained from u n c o n d i t i o n a l averages. Chapter 3. The DCMC Method with Stochastic Processes 23 Z Figure 3.1: Conditional average of temperature; blue(.) D N S , red(-) (T\Z), black(-) T Researchers have previously explored the application of the C M C method and have found that using only one conditional variable offers good results i n terms of species mass fraction and temperature predictions. However, the presence of triple flames and prediction of ignition, extinction, and re-ignition phenomena are not well reproduced by just one conditioning variable. Figure 3.2, from Bushe et al [8], shows a slice of a 3-dimensional D N S calculation of temperature. The contour lines represent isopleths of mixture fraction, which were formed by the intersection of the slicing plane and the corresponding isosurfaces of mixture fraction. The white isopleth corresponds to the stoichiometric mixture fraction. Figure 3.2: D N S calculation of temperature, from Bushe [8] It is clear that this temperature field exhibits strong fluctuations in regions of constant mixture fraction. This is evident in the 'cold pockets' shown in regions A and B in the figure. These cold zones are typical of extinguishing flames. From this figure, it can also be Chapter 3. The DCMC Method with Stochastic Processes 24 observed that there is a transition between hot zones, located to the right of region B. This transition could correspond to re-ignition or bifurcation of the flame. It has been proposed that using a second conditional variable would solve those issues [9], [16], [63]. This result is of paramount importance, given the close relationship between these phenomena and the formation of pollutants. In this thesis, the use of two conditioning variables has been explored. One variable is mixture fraction Z, which characterises the state of mixedness of the flow for a given instant in time and space. The other conditional variable is named a and is introduced to capture fluctuations along an isosurface of mixture fraction. This is achieved through the definition of the scalar a as a conserved scalar with a spatial gradient perpendicular to that of mixture fraction. 3.1 Conditional Moment Closure with Two Conditional Variables ( D C M C ) The use of a second conditional variable has been proposed as a way to solve the shortcomings of the original C M C method. Bushe [16] used a second conditioning variable that filtered between regions of a surface defined by a constant mixture fraction, hence detecting differences in the stretch rate history of the flow. The results of a simulation of auto-ignition of n-Heptane showed an improvement over the single-conditional method i n detecting autoignition and double flames when compared to DNS databases. However, the same results failed to show the presence of triple flames present in the reference data. Cha et al [63] used the scalar dissipation rate as the second conditioning variable for the simulation of a methane-air turbulent reacting flow. After comparing their results with those obtained from an experimental jet and a DNS simulation, it was found that there exists good agreement in the prediction of extinction phenomena, but they obtained an early prediction of the onset of re-ignition. High fluctuations around the conditional means, specifically for low values of the scalar dissipation rate, were identified as a possible source for the discrepancies with the reference data. 3.2 Derivation of the D C M C Equations During the development of this project, two conditioning conditional variables have been used. • Mixture fraction Z — Z (x{,t) , 25 Chapter 3. The DCMC Method with Stochastic Processes • and a = a (xi,t). These two random variables are related to the non-random variables £ and a using: w|(Z = C,o = a) * co (r|(Z = C,« = a ) , Y / | ( Z = C,a = a ) ) , (3.1) such that the result is a function of what values are chosen for the conditioning variables £ and a . The definition of the conditional average of the mass fraction of the P h species is C,a(xi,t)=a), Y[ = Yi{t,a;Xi,t) = (Yi(xi,t)\Z{xi,t) = so that the instantaneous value of species mass fraction using conditional average and fluctuations around that mean can be expressed as + I- Yi(xi,t;n) = 17U( t;n)a(x,t;n) Xj) Eq. ) Y i (2.3) is substituted i n E q . (2.2) to obtain a new version of the transport equation of species mass fractions as dpYi dt dpuiYj dxi + = _d_ ( 9YI\ dxi \ dxi Di ^ + 2 ) J Using the experimental results from Vargaftik [64], it has been proposed [16] that pD\ = constant, and that Dz = Di, implying that Le = 1. Using these assumptions, E q . (3.2) can be written as dYr | dYr u dt 1 dxi D b^Yj Co! = dxidxi p To express E q . (3.3) in terms of conditional averages, it is required to derive partial derivatives, which were obtained by applying the chain rule of differential calculus. • The partial derivative of species mass fraction with respect to time is dt dt dC dt da dt dt ' • The first partial derivative of species mass fraction with respect to space is m dx{ = w + dx{ mM + dC, dxi m ^ + da dxi dY[ dxi ( ' ' 26 Chapter 3. The DCMC Method with Stochastic Processes The second partial derivative of species mass fraction with respect to space is &Yi dWjdZdZ c^l?dada d Yj dZ da dC, dxi dxi da dC,da dxi dxi dxidxi 2 2 dxi dxt 2 3 Y } dZ d Yj da dYi d Z +2 1-2 1'dxidC,dx{ ' ~dxidadxi ' d( dxtdxi 2 2 2 n dY J^ F*_ &Y L + + + L da dxidxi ( 3 6 ) dxidxi dxidxi Equations (3.4), (3.5), and (3.6) were substituted in E q . (3.3) producing, *' *• + *r dC, + dY Yj F da ° &Y'MO dC, dx\ dxi d T l ^Yj da z n 2 dC,da dxi dx{ 2 dx^dC, dxi da da dxi dxi dxida dxi 1 where the following definitions are used: • The transport equation of the conditional average of the species mass fraction is = dy? ! i dy E dt Y l D dxi dY 2 dxidxi • The transport equation of mixture fraction, E q . (2.1), is written here as dt Z * dx{ dxidxi ^ ^ • The transport equation of the scalar a is „ E da a = m + da U i d^- ^ da 2 D d^dx- n = - 0 ,„ „„. ( 3 - 1 0 ) • The transport equation of the fluctuations around the conditional average of species mass fraction is * Y 3.3 dt * dxi dxidxi ^ ^ Modelling Scalar Dissipation Scalar dissipation, Xz = M !^' 7 ^ used to represent the coupling that exists between fluid dynamics and the chemical reactions taking place i n a turbulent reacting flow [23]. Previous research [1] indicated that large values of scalar dissipation rates lead to extinction of flames. Conversely, ignition phenomena are associated to small scalar dissipation rates. Chapter 3. The DCMC Method with Stochastic Processes 27 Peters [1] proposed a model in which the scalar dissipation rate diminished as a function of mixture fraction with time i n a one-dimensional mixing layer. This model is frequently used i n the flamelet formulation discussed i n section 2.6.4. However, the model is not capable of capturing extinction or ignition. A n improvement to better reproduce extinction phenomena was proposed by Peters [1] i n which a conditional stoichiometric scalar dissipation rate, assumed to have a lognormal P D F , was used to calculate the right proportion of burning flamelets during an extinction event. Using C M C and the stretched diffusion concept, Bushe [16] proposed a model for the conditional scalar dissipation rate, which was found to provide good approximations of values of the scalar dissipation rate when compared to DNS of a non-reactive flow of Eswaran and Pope [65]. This same model, however, did not provide good estimations of ignition due to long decreasing times of the scalar dissipation [16]. A n improvement of this model was to incorporate a second conditioning variable to the system of transport equations of species mass fractions, which was also a passive scalar. This model was found to provide good predictions of ignition and double-flame formation of an n-Heptane non-premixed flame when compared to data obtained by other researchers [66], [67], [68]. However, this formulation d i d not predict the presence of triple flames found i n other research papers [67]. In this thesis, the D C M C method has been applied to define conditionally averaged transport equations for the scalar dissipation of Z and a. The conditional average of scalar dissipation is defined as Xl = Xl(t,(x;xi,t) = (xz{xi,t)\Z(xi,t) = t,a{xi,t) = a) leading to an instantaneous value of Xz{%k,t) = Xz\c=Z(xi,t),a=a(xi,t) (3.12) + x' z Ruetsch and Maxey [69] proposed a transport equation for the scalar dissipation, with the terms on the rhs of the equation representing production, dissipation, and diffusion: DXz Dt OXi OXj d fdZ\ dxj \dxi J + D ______ dxjdxj (3.13) The Sjj term i n E q . (3.13) corresponds to the strain tensor that denotes the production of scalar dissipation through vortex stretching mechanisms. In premixed, turbulent reacting Chapter 3. The DCMC Method with Stochastic Processes 28 flows, strain is known to have a strong effect over the stretch history of the flame, which consequently affects the flame structure [11]. The definition of Sij [14] is, (3-14) In turbulent, non-premixed flames, high levels of strain imply high levels of scalar dissipation, which i n turn generate local extinction of the flame, since the chemical reactions cannot take place given the extremely high rate at which reactants are diffused into and through the reaction zone. Williams [11] provided a flamelet-based analysis of this phenomenon. Thevenin and Candel [67] studied the effect of variable strain on the ignition of a non-premixed flow of air and hydrogen through the use of computer simulations. While their research predicted self-ignition and ignition times i n agreement with the experimental results of Cheng et a l [70], the triple-flame phenomenon was observed only for the case where there was no strain present. Using the definition of a material derivative, E q . (3.13) can be expanded to obtain z dZ 2 dxz dt = ~Ui-dxj uj dxi dxj 1 OXjOXj -2D dxjdxii J (3.15) \dxjdxij In order to express E q . (3.15) i n terms of the conditioning variables, it is required to obtain derivatives analogues to Eqs. (3.4), (3.5), and (3.6). The definitions of Eqs. (3.9) and (3.10) are also used to produce dxi at ^ &n HiZ AC da ^ x l d Z dZ dC dxi dxi 2 TP dxi X z ° l ™ dxidxj S lJ i + D - ^ dxidxi d xl da da ^^d ^ dZ da da dxi dxj dCda dxi dxi 2 2 2 d z\ f d z \ 2 dxida dxi 2 2 ( \ dxj dxij \ dxj dxi J d Xz dZ 2 dxidC dxi (3.16) where the transport equation for the conditionally averaged fluctuation of scalar dissipation is Ezr = X z dx'z , + dt dx dxi Ui- D d x!z dxidxi 2 (3.17) The conditionally averaged transport equation of the scalar dissipation Xa follows analogously to Eqs. (3.16) and (3.17). Chapter 3. The DCMC Method with Stochastic Processes 3.4 29 Scalar Transport Closure Hypotheses Equations (3.7) and (3.16) contain unclosed terms for which closure hypotheses are formulated in this thesis. W i t h the following scalar transport hypotheses, some of the terms i n those equations are closed: • Z and a are conserved scalars, hence its lack of a source term as shown i n Eqs. (2.1) and (3.10); • The spatial gradient of a is perpendicular to the spatial gradient of Z to better capture the effects of fluctuations along isosurfaces of mixture fraction. Using these hypotheses, E q . (3.7) is reduced to — = &Yi dZ dZ dC dxi dxi &Y, ^Yj da D 2 2D d 2 T l d Z dxidC, dxi 2D 2 da da dx{ dxi 1E dxida dxi d W l d , Q Yj (3.18) E q . (3.16) is also reduced to dx~z dt _ dxz * dx{ E * x ™ 2 B*S„+ ' dxi dxj dC dxi dxi %3 da 2 dxidxi dxi dxi 2 dxidC, dxi ( &z \ f d z \ 2 +2DS^^-2D dxida dxi \dxjdxiJ (3.19) \dxjdxi J Multiplying Eq.(3.18) by p and taking its conditional average yields dW DC" t LOI = pEyT + pEyJ „ <?Y, dxidC, dZ dxj 0 dZ_dZ_ &Yj dxi dxi da d Yj da da dxi dxi D 2 2 -pD— dxida (3.20) dxi Conditionally-averaging E q . (3.19) yields dt = dxi dZ dZ -E^r-u ^-2——S dxi dxi dxj i J J 2 X ij d Tz dZ dZ d Tzda da d Xz~ dZ r,>o dxi ~—z= - » -—— dC dxi h U da dxi dx{ V IU-dxidC, dxi 2 + dx + D dxidxi 2 2 2 2 / dZ \ ( dZ \ \ dxj dxi J \ dxjdxi J 2 dxida oxi 2 (3.21) 30 Chapter 3. The DCMC Method with Stochastic Processes 3.5 D C M C Closure Hypotheses Equations (3.20) and (3.21) were reduced to a conditionally-averaged form i n the previous section. Whilst closure was determined for some unclosed terms, both equations still contain terms that can be further closed using the D C M C closure hypotheses described in the following modelling assumptions. • D C M C Modelling Assumption 1: the fluctuations around the conditional mean of species mass fraction and around the conditional mean of scalar dissipation are only functions of the scalar variables Z and a, independent of time and space. This produces tW BY 1 cPY' 7 ^+m^-D^ut OX, = 0. (3.23) OXiOXi • D C M C Modelling Assumption 2: assuming that the flow is homogeneous and isotropic, the ensemble of conditional averages of species mass fractions and scalar dissipation are equal at each point i n space for the same instant i n time, and the spatial gradient of the conditional mean of the species mass fraction is equal to zero [16]. This yields _ 3 * i ^dx- = 0 (3.24) = 0 (3.25) = 0 (3.26) = 0 (3.27) = 0 (3.28) = 0 (3.29) d Xl dZ = dxidC, dxi 0 (3.30) d x~z da dxida dxi 0. (3.31) ~D d W l dxidxi d Yj dZ 2 D dxidC, cPYj dxjda p dxi D da dxi dx z OXi D > dxidxi d2T 2 D 2 - Substituting Eqs. (3.22) and (3.24) - (3.27) in E q . (3.20) results in _dYj -d7- T " I = p PTi W P dZdZ d Yi dx^dx- - da 2 n n da L^^dx-dx-- ( 3 - 3 2 ) 31 Chapter 3. The DCMC Method with Stochastic Processes Here, t h e t e r m o n t h e lhs represents t h e conditional average o f chemical source t e r m , w h i l e the second a n d t h i r d terms o n t h e rhs represent t h e conditional average o f diffusion o f t h e species mass fractions i n t h e scalar field. S u b s t i t u t i n g Eqs. (3.17) a n d (3.28) - (3.31) i n E q . (3.21) results i n dxTz dZ dZ = - 2 Sij + D—^xl2 dt ~dxidxj" ' dC oc lJ 2 8Z 2 +D Xa da ' da " 22 /v 2D- dZ 2 dxjdxidxjdxi (3.33) Here, t h e first t e r m o n t h e r h s represents t h e conditional average o f p r o d u c t i o n o f scalar dissipation, the second and t h i r d terms represent t h e conditional average o f diffusion o f scalar dissipation, and t h e t h i r d t e r m represents t h e conditional average o f dissipation o f scalar dissipation. 3.6 Mathematical solution of the Z-gradient squared term We define dz d z 2 9 = 2 (3.34) dxjdxi dxjdxi dXz dxz dxj dxj (3.35) E x p a n d i n g 6 yields, dXz dxz dxj dxj (dxj) d (OZ^dZY dxj \dxi dxi J d_ (MdZ\ dx* \dxi dxi „ d Z dZ " 2 dxjdxi dxi J , d Z dZ 2 \dxjdxidxjdxi) dxjdxi dxi C o m p a r i n g Eq. \dxidxi) (3.36) t o E q . (3.34) i t can be seen t h a t 0 = 4$;vv > d an s 1 4x. 0. (3.37) 32 Chapter 3. The DCMC Method with Stochastic Processes Using E q . (3.12) a n d t a k i n g p a r t i a l derivatives o f Xz w i t h respect t o space, E q . (3.35) can b e expressed as \ dxj J \ dC ) 9XJ \da dxj ) \ dxj ) .^dx^dx^dZ dxidxi da ^dxl dZ dxi da dxj dC, dxj ' " dxi da dxj ' " d£ dxj da dxj | x | ^dXzdx'z ^ dx~z dZ dx! dxj dxj c?C dxj dxj 2 z | jTXz da d ' da dxj dxj X z ( 3 3 g ) Using t h e scalar t r a n s p o r t closure hypotheses and t h e m o d e l l i n g assumptions, E q . (3.38) reduces t o S u b s t i t u t i n g E q . (3.39) i n E q . (3.37) a n d conditionally-averaging, yields 4 V d( J A\da) I t is f u r t h e r assumed t h a t / Xa \ Xa \Xz) (3.41) xl' Using t h i s assumption i n E q . (3.40) and s u b s t i t u t i n g t h i s result i n E q . (3.33) produces the final version o f t h e conditionally-averaged t r a n s p o r t equation o f scalar dissipation o f m i x t u r e f r a c t i o n and o f t h e scalar a dxi dt = -2-^L^S dxi dxj +D^Xz dxZ dt = r>da da -^(^V-^/^Vxi 2 V dC j 2\da)x~z 13 + D ^ x dxi dxj +D^Xz (3.42) a _ P_ (dx^\ i3 + D ^ 2 V T dC J a . _ & (d)g\ 2 V da Xz J xi (3-43) T h e first t e r m o n t h e r h s o f these transport equations describes t h e c o n d i t i o n a l average o f p r o d u c t i o n o f scalar dissipation. T h i s t e r m represents t h e m a p p i n g o f t h e s t r a i n tensor Sij o n t o t h e corresponding diffusing scalar's surface. T h e second a n d t h i r d t e r m s represent t h e c o n d i t i o n a l average o f dissipation o f scalar dissipation, while t h e last t w o terms represent t h e c o n d i t i o n a l average o f diffusion o f scalar dissipation. A l l t h e terms i n these equations are closed, w i t h the exception o f t h e t e r m containing t h e s t r a i n tensor Sij, w h i c h has t o be modelled. Chapter 3. The DCMC Method with Stochastic Processes 3.7 33 Models for the Transport Equation of Scalar Dissipation Two different models using stochastic processes are proposed to model the strain tensor Sij. The first model uses a periodic force with random phase shifting to simulate the chaotic behaviour of the strain tensor. A second model uses the Coupled M a p Lattices ( C M L ) developed by Beck and Hilgers [5], [6], [7] to simulate velocity increments i n an isotropic, homogeneous, fully-developed turbulent flow, from which the corresponding strain tensor is calculated. Pitsch and Fedotov [71] simulated the fluctuations of the scalar dissipation rate in nonpremixed combustion. They used this approach i n combination with a flamelet formulation to simulate a turbulent, non-premixed reactive flow of methane and air. They reported that the extinction process occured at smaller time scales than those relevant to the overall phenomenon, like turbulent or chemical time scales. This seems to indicate that the deviations are due to Gaussian white noise, resulting in a Wiener process. This research group also found delimiting regimes that drive the burning and extinguishing states and that even small-amplitude scalar dissipation fluctuations have considerable effects on the state of the flame, taking it from the burning to the extinguished state. This finding i n particular, further reinforced the need to study the influence of the scalar dissipation rate on a turbulent, non-premixed reacting flow, which is one objective of this thesis. 3.7.1 Periodic Forcing with Random Phase Shifting (PF) The use of forcing schemes involving harmonic functions has been explored i n DNS simulations of isotropic, homogeneous turbulence [72], [73]. The general idea behind this method is to generate an influx of energy at low Fourier wavenumbers that compensates for the dissipative effects at the small scales of turbulence. Eswaran and Pope [73], used a stochastic forcing scheme with an Ornstein-Uhlenbeck formulation to study the validity and capabilities of the method i n an isotropic, homogeneous turbulent flow. Hilgers et al [74] explored the use of stochastic resonance to the study of motion of an overdamped particle in a bistable potential, subject to a external periodic perturbations and coloured noise. Chapter 3. The DCMC Method with Stochastic Processes 34 Similar to the work by Eswaran and Pope [73], the first model proposed i n this thesis uses a combination of harmonic signals to mimic the chaotic influence of the turbulent flow field on the strain tensor by using the product of the signals i n the x and y directions: = A sin (w x + 6 ) sin (u> y + 0 ). n nx nx ny ny (3-44) Here, • The magnitude of the amplitudes is A i n+ = A /2, and is scaled by the Reynolds number n of the flow. • The period of the signals is w + i = 2w and is used as an emulator of the viscosity of the n n flow. Larger periodicity implies a higher number of transitions, which are interpreted as lower viscosity. • The random phase shift 0 is produced with Gaussian random numbers. This shifting n is used to produce random fluctuations i n the strain field. • n = 1...N is the number of signals that are then combined to obtain the assumed strain tensor field in the x — y space: n Stf(x,y) = _T* . i=l n (3.45) The coupling between the strain field and the flow is obtained by forcing proportionality between the amplitudes A n of the signals and the magnitude of the characteristic Reynolds number of the flow. This implies that this scheme can only simulate strain fields for one energy level at a time and must be calibrated for each energy level independently. For decaying, isotropic, homogeneous, turbulence, characteristic parameters can be used to scale the magnitude of the strain field accordingly with turbulence theory [14] as, ISyl = ^ (3-46) The structure of the strain fields simulated with the P F model can be modified by varying the number of signals used in the calculation of E q . (3.45), by varying the scaling of the amplitudes between signals, or by varying the periodicity of Eq. (3.44). The P F scheme was used in this thesis to simulate strain fields with fluctuations and magnitudes relevant to the combustion process studied. The magnitudes were compared to Chapter 3. The DCMC Method with Stochastic Processes 35 those obtained from the D N S of Bushe [8] and the simulated fields in the x — y space were then transformed to the a — £ space using a simple mapping convention. The mapping is performed by assuming a linear correspondence between a and x, and £ = erf(y). 3.7.2 Coupled Map Lattices ( C M L ) The second model that has been implemented uses the turbulent flow simulation proposed in the works by Beck [5] and Hilgers and Beck [6], [7]. This model is a coupled map lattice that mimics the energy cascade in a turbulent flow. The system of equations is as follows Xd,n{i) F(x , -i{i)) d n +<«- i(0(i->»-i)«iJr-}(o. 1 ) (3.47) In this model, • F(x) is a map that models the energy input from the large scales and is defined as F(x) = 1 - 2a: , where x is the position of stochastic particles in D spatial dimensions. 2 • The parameter k corresponds to the energy level in the cascade. High-values of A; imply higher energy-dissipation processes, which leads to a lower energy level. • The coupling parameter g is associated with, and should be scaled to, the molecular viscosity, p, of the fluid. • The damping parameter A is defined as A = exp(—jrk), where the product JT is proportional to the inverse Reynolds number. • The product C £ is used to determine the random fraction of the driving momentum that each of the 'daughter' eddies receive from the previous higher scale. The calculated temporal velocity difference u is later used to calculate the strain rate tensor using E q . (3.46), or, under the assumption of decaying, isotropic, homogeneous turbulence, it can be further demonstrated [14] that, (3.48) Chapter 3. The DCMC Method with Stochastic Processes 36 Using the C M L model, the strain fields were simulated i n the x — y plane and later mapped onto the a — £ plane using the same mapping convention as used i n the P F model: (a, C) —> (z,erf(y)). It is clear that while this model is capable of simulating strain fields for a wide range of energy levels, the main focus of this research is on the simulation of strain fields at the energy levels relevant to the combustion process studied. This means that the C M L had to be capable of simulating velocity increment fields at relatively small scales. This calibration was made possible by using the variable k that denotes the energy level at which the velocity fields are calculated. The target value of k was obtained by running a number of simulations and comparing the results with statistics from the DNS of Bushe et al [8]. A consequence of this is that all the previous energy levels still have to be calculated, which inevitably increased the computation time. 3.8 Summary The D C M C equations that describe turbulent combustion events using two conditioning variables have been derived i n this chapter, resulting i n the following system of equations ZJ W l - ~ d * / dZ DYL dZ 2 D ^ dt d Y dC, ^ dxi dxi dt dxi dxj dXa~ dC, daaa dt dxi dxj 13 2 X z 2 X 7 D 2 2 \da ) x~z a _ D fdxA _ D fdxZ\ xl 2 V d£ J 2 V da J x~a 2 13 da da ^ dx% dxi 2 V #C / da da 2 2 +D^Xz+D^x-a. 2 (3.49) The properties and characteristics of the second conditional variable a have also been described and applied during the derivation of the equations. In this chapter, closure for the D C M C equations has been proposed using scalar transport and C M C closure hypotheses. The resulting equation contains a strain rate tensor term that has been modelled using two different stochastic processes. The results obtained from the use of these models are presented in the following chapter. 37 Chapter 4 Simulation of Methane Oxidation with D C M C - S P T h e system of equations developed in the previous chapter, E q . (3.49), is used to simulate the combustion of methane and air. The present chapter describes the results obtained from this computational simulation. The chemical kinetic mechanism, as well as the initial and boundary conditions developed for the generation of both the D N S and C M C databases have been used i n this research. The D N S database has also been used as a source of information for the evolution of the turbulent flow. This means that the results presented in this thesis have been produced i n an a priori sense. A comparison between the D C M C - S P results of this thesis and the C M C [9] and D N S [8] for a similar test case shows that the D C M C - S P implementation offers an improvement over the singly-conditional moment closure method. 4.1 The Reference D N S and C M C Databases T h e D N S database of Bushe et al [8] described the results of a number of independent runs of a shear-free, decaying, turbulent mixing layer. A C M C simulation of the same test case [9] provided good predictions of temperature and mass fractions of the major species while falling short of providing good predictions of NO and intermediate mass fractions. 4.1.1 Flow Field The flow field simulated in the D N S database is a three-dimensional shear-free, temporal mixing layer, discretised by a 240 x 120 x 120 grid, in which the governing equations for incompressible flow were resolved. The domain is periodic in two directions, while allowing outflow in the third direction. The domain is shown i n Fig. 4.1. Here, region A denotes the approximate spatial domain where the mixing between fuel and oxidiser occurs. Chapter 4. Simulation of Methane Oxidation with DCMC-SP 38 q Oxidiser Fuel Figure 4.1: D N S computational domain The initial turbulent flow field was forced losing the pseudo-spectral scheme of Ruetsch and Maxey [69] that simulates a statistically stationary, incompressible turbulent flow. The initial Taylor Reynolds number was 59 and the Prandtl and Schmidt numbers used were both set to a constant value of 0.75. 4.1.2 Chemical Kinetics The chemical kinetic mechanism developed for the DNS database [8] has been used in this research project to determine the effectiveness of the model proposed i n Chapter 3. This mechanism involves three steps and is a modification of the two-step reduced mechanism originally proposed by Williams [75], and later modified by Swaminathan and Bilger [76] for the oxidation of methane, with the addition of a third step for NO formation with the Zel'dovich mechanism. The resulting three-step mechanism is written here as F + Oxi-> I+P I-\-Oxi 2P N + Oxi 2 -> -> 2NO where F is CH4, Oxi is O2, and / - P = (\* + lco) \^H 0+\c0 2 2 A n important feature of this mechanism is the presence of a steady-state approximation for the mass fraction of the Hydrogen radical. This expression has been designed to emulate the Chapter 4. Simulation of Methane Oxidation with DCMC-SP 39 reduction i n the concentration of the Hydrogen radical at low temperatures. T h e mechanism was f u r t h e r simplified by Bushe et al [8] by assuming peak temperatures i n the 1200-2230 K range, so t h a t the c o m p u t a t i o n a l expense is reduced, while keeping acceptable levels of accuracy. T h e mechanism was non-dimensionalised by assuming t h a t the constant pressure specific heat C , the r a t i o of specific heats 7, a n d the speed of sound c remain constant. I t was p also f o u n d t h a t the region where the chemical reaction takes place needed t o be broadened t o avoid the c o m p u t a t i o n a l expense associated w i t h the calculation of sharp gradients. Initial Conditions 4.1.3 Since the chemical kinetic mechanism used i n the D N S is not capable of a u t o - i g n i t i o n , t h e i n i t i a l conditions of the scalars were taken f r o m a planar, l a m i n a r flame, w h i c h also served as the i n i t i a l conditions for the C M C simulations. T h e i n i t i a l conditions as a f u n c t i o n of m i x t u r e f r a c t i o n are shown i n F i g . 4.2. T h i s laminar flame was located i n the plane formed by points pqrs i n F i g . 4 . 1 . 0.35 CH 0.3 -*- 0.25 4 INT PROD NOx10 7 0.15 P 0.1 I 0.05 0, 4 \ y $ \ 0.2 0.8 (a) Initial Conditions of Species Mass Fractions (b) Initial Conditions of Temperature Figure 4.2: I n i t i a l Conditions of the Scalar Fields 4.2 Simulation of CH 4 Oxidation Using D C M C - S P Using the chemical kinetic mechanism described i n the previous section, E q . (3.49) is expanded t o include the species present i n the chemical kinetics mechanism, as well as the c o n d i t i o n i n g Chapter 4. Simulation of Methane Oxidation with DCMC-SP 40 scalar variables, resulting i n the following system of equations _dY 02 _d¥ P-^- F ~W , d Yp -- dY 2 F -- d Yj ' ^ = » _3YNO -r= ~dt dXa -. dt ^ P F 2 + -^ pD F r 2 Xa F / D n x > - M + d Y -fo?-P x* 2 P D 2 dY 2 , NO NO X d dT a Xa „ dT 2 pD 2 2 "T+dcJ^X'+da*^ = 2 ftr* = 2 dY + —^-pD z+ , UNO 2 ^ ^r^*° dY = d< ___ + d Yp + P _&T p n 2 p _ 2 _ _ 9 2 P P—QJ- n Xz w + _d¥ ~dT p n 2 r = _ , 2 = "F + -^ pD _dYj P -- da da ? ft^ _ dZdZ_2 _ __ 2 V aC 7 V __V___T\ __ V d£ J 2\da) _P_(9xi\ tJ 2 __ z 2 V da ) Xa~ 2 dxi dxj +D^Xz _ __ + D^x-a. 2 xl (4.1) This system of equations was solved for a 50 x 50, a — £ grid, with the initial calculation of the strain fields i n the x — y space and later mapped onto the a — £ space using the simple mapping convention described in Sections 3.7.1 and 3.7.2. The boundary conditions for the strain fields were set to be periodic for both directions. The boundary conditions for all the reacting scalars were set to be periodic in the direction of a and constant i n the direction of C The initial conditions of the reacting scalars and of scalar dissipation of mixture fraction were taken directly from the laminar flame described i n Section 4.1.3. The results obtained from the numerical implementation of the system i n E q . (4.1) are presented i n this section. Once the strain field was simulated using either of the stochastic processes, the O D E solver D V O D E [77] was used to solve E q . (4.1). F i g . 4.3 serves as an illustration of the structure of the code. Chapter 4. Simulation of Methane Oxidation with DCMC-SP Initial conditions DNS database Calculate strain tensor with stochastic process |Sy| = ( 2 u ( t ) / E ( t ) ) 1/2 41 ODE solver dt df _ ~dt~ <?X . dt ' Write output fields. Compare to DNS and CMC data. Figure 4.3: D C M C - S P code structure 4.2.1 Simulation of the Strain Field As discussed in Chapter 3, the strain field has been simulated using two different stochastic processes. This section summarises the results obtained for both cases. Unless expressed otherwise, the figures i n this section show results obtained after a calculation time of t — 5 time units. In this computational implementation, the following model has been used: dZ n dZ „ -2-—-—Sij = dxi dxj n — -2xzniUjbij -2xzniUjSij. Where re, and rij axe unity vectors in the direction of the Xi and Xj directions of the gradients of Z. T h e value of niUj is then [—1 + 1]. Bushe and Cant [78] showed that positive values of niTij are statistically more frequent; as a result, it is assumed i n this research that niUj is negative i f the strain Sij is negative and positive otherwise. This model applies analogously to the corresponding term i n the transport equation of Xa- Chapter 4. Simulation of Methane Oxidation with DCMC-SP 42 Periodic Forcing with Random Phase Shifting ( P F ) This model was validated using the stochastic variables X(x, y; t) and Y(x, y; t) i n the solution of Eqs. (3.44) and (3.45). Using statistics from the DNS database of Bushe et al [8], it was found that the appropriate values to simulate a strain field with similar characteristics using the P F model were, • Ai = 1.585 this parameter is proportional to the Taylor Reynolds number of the D N S database • &ix = uiy = 2n which indicates the periodicity of the signals. This parameter was set empirically by comparing the structure of the strain fields from the D N S database and the strain fields obtained using the P F model. The magnitude of the strain tensor was calculated using the values of the parameters e and v from the D N S database. In all cases, the strain field was first calculated in the x — y plane and later transformed to the a — C plane following the mapping rules described in Section 3.6. Figure 4.4(a) shows the results of the calculation of the strain field i n the x — y plane using ^ = 3 periodic signals. Figure 4.4(b) shows the same strain field after being transformed to the a — £ plane. Comparing Figs. 4.4(a) and 4.4(b) it can be seen that the main effect of the coordinate transformation is an apparent reduction of the area of influence of the structures close to mixture fraction values of 0 and 1, followed by an elongation of the area of influence of those structures between the aforementioned regions of Z. These two figures also show that the overall magnitudes of the strain field are maintained, while the reduction - elongation - reduction effect makes some of the structures disappear, especially for regions in the vicinity of Z = 0. T h i s finding is of relevance since the areas where the elongation effect is prevalent, is where values of mixture fraction indicate a more vigorous mixing process, where chemical processes are taking place with more frequency. As expected, the reduction - elongation - reduction effect does not take place along the a direction. The strain fields are also shown to be periodic in the x, y, a, and £ directions. Chapter 4. Simulation of Methane Oxidation with DCMC-SP 43 Figures 4.4(c) a n d 4.4(d), show the s t r a i n fields i n the a — ( space for \I> = 4 a n d \& = 5 p e r i o d i c signals, respectively. F r o m Figs. 4.4(b), 4.4(c), a n d 4.4(d), it c a n be seen that the number of signals used has a significant effect o n the structure a n d propagation of the calculated field. (c) * = 4, Q- C (d) * = 5, a - C F i g u r e 4.4: S t r a i n fields simulated w i t h P F T h e s t r a i n field calculated w i t h three periodic signals shows well-defined regions of positive a n d negative values of strain, w h i c h have the physical meaning of stretch and compression regions i n the flow. These structures tend to be very large for this p a r t i c u l a r case a n d do not show m u c h intermittency. T h e structures c a n be considered to be eddies, w h i c h i m p l i e s that large structures correspond to the energy-containing subrange or the i n e r t i a l subrange Chapter 4. Simulation of Methane Oxidation with DCMC-SP 44 i n the cascade model. Consequently, smaller, finer structures would represent the behaviour of a flow at the dissipation subrange, at which combustion is assumed to take place. Figures 4.4(c) and 4.4(d) show a much finer structure with more intermittency than Fig. 4.4(b) while yielding similar strain magnitudes to those obtained from the D N S database. Coupled M a p Lattice ( C M L ) In Chapter 3, the C M L equation of Hilgers and Beck [6], [7] that was used to calculate the strain field i n this thesis was described. Opposite to the P F model, the C M L model simulates velocities from which the strain field has to be calculated. T w o independent velocity fields were simulated and assigned arbitrary directions u and v, which are orthogonal to each other. The simulations were run i n the x — y space, with parameters 7T = 0.02, g = 0.0194, and C = 1.411 i n Eq. (3.47) and transformed to the a — C, space using the mapping rules described i n Section 3.6. Figure 4.5(a) shows the the u-velocity field simulated using the C M L model at a level of fc = 3, which is assumed to be i n the integral sub-range of the energy cascade of F i g . 2.1. Figures 4.5(b) and 4.5(c) show the results of the simulation of velocity for levels offc= 9 and fc = 15, respectively. The structures i n these figures do not correspond to eddies, since eddies represent a vortical motion and these figures are only a two-dimensional representation of one component of velocity. A comparison between Figs. 4.5(a), 4.5(b) and 4.5(c) shows that the velocity magnitude decreases about 4 orders of magnitude from levelsfc= 3 tofc= 15. This proves the capabilities of the model to represent the dissipation of energy from one level to the next. It can also be seen that the structure of the velocity field is more uniform as the fc level is increased. Figures 4.6(a), 4.6(b), and 4.6(c) show the vector form of the two-dimensional velocity fields at different levels of energy for the area delimited by 0.2 < a < 0.6 and 0.2 < £ < 0.6. The fields were calculated using the results of the simulations of velocities i n two orthogonal directions for levels offc= 3,fc= 9, and fc = 15, respectively. Chapter 4. Simulation of Methane Oxidation with DCMC-SP 45 C o m p a r i n g these figures, i t can be seen t h a t some o f t h e structures present i n F i g . 4.6(a) have been completely dissipated, or have been dissipated into smaller structures i n F i g . 4.6(c), g i v i n g rise t o a more orderly field. 0.2 0.4 a 0.6 0.8 0.2 (a) k = 3 0.4 a 0.6 (b) it = 9 lo.02 0.8 0.6 3 0.4 0.02 0, 0.2 0.4 a 0.6 • H 0.8 ___-0.04 (c) k = 15 Figure 4.5: u-velocity fields simulated w i t h C M L Figs. 4.7(a), 4.7(b), and 4.7(c) show the strain fields obtained using the C M L model a n d t h e d e f i n i t i o n o f E q . (3.48) f o r energy levels o f k = 3, k = 9, a n d k = 15, respectively. T h e effects o f t h e decaying energy o n t h e s t r a i n fields are evident i n these figures. I t can be seen t h a t t h e magnitude o f the s t r a i n field decays considerably as k advances f u r t h e r i n t o t h e viscous subrange o f the cascade model. T h i s leads t o a strain field w i t h fewer peaks a n d more neutral, or zero-strain areas as shown i n F i g . 4.7(c). Using the values o f t h e turbulence Chapter 4. Simulation of Methane Oxidation with DCMC-SP 46 parameters e and v from the reference D N S database, it was determined that the level of energy in the C M L model that would give the appropriate values corresponds to k = 9, which is shown in Fig. 4.7(b). 0.6 i \ 0.6i , „ \ / 1 \ f t \ * 0.5 - * * / m *• 0.5 ** * • . • » - *• f f * - •#* • - y * -* • 0.4 — to • •*? f / P • - • ' * * » . . 0.3 m 0.3 ft., ft., * * - •» .* r » - ^ ^ / -» »+ / / m * + / / v / - 5 \ I —- \ • • 0.3 N ' t « i r - > 0.4 / # _ * » ! t / / / / _ , \ . / - . _ I V 0.5 0.6 a 0.3 (a) k 0.4 a (b) k = 9 0.6 r 0.5h 0.4 h » 0.3K °-8s \ » 0.3 0.4 a * - M 0.5 0.6 (c) k = 15 Figure 4.6: Velocity fields simulated with C M L 0.5 0.6 Chapter 4. Simulation of Methane Oxidation with 0.2 0.4 _ 0.6 (c) DCMC-SP 47 0.8 k = 15 Figure 4.7: S t r a i n fields simulated w i t h C M L C o m p a r i n g F i g . 4.7(b) w i t h F i g . 4.4(d), it can be seen t h a t there are i m p o r t a n t differences a n d similarities between the structure of the fields produced by the P F a n d the C M L models. F i r s t , the P F model is based on spectral methods which are inherently filtering smooth due t o the t h a t occurs as a consequence of Fourier-like transformations. T h i s is proven by the fact t h a t the structures produced by the P F are larger t h a n those simulated w i t h the C M L model. I t can also be observed t h a t , i n general, the fields produced by the P F behave very m u c h as an amplified p o r t i o n o f a selected area of the C M L - p r o d u c e d field. T h i s discrepancy could be eliminated by using a larger number of signals i n the P F model t h a t would induce more v a r i a b i l i t y along the d o m a i n . One of the similarities between the models is their capability t o Chapter 4. Simulation of Methane Oxidation with DCMC-SP 48 provide the right orders of magnitude for the field when compared to the calculated value of the strain tensor using v and e parameters from the DNS database. A t t = 5 time units, the DNS data yields Sij = y^e/2i/ = ±0.5656, while the both the C M L and P F models provide values i n the range —0.8 < Sij < 0.8 for the case of k = 9, and ^ = 5, respectively. 4.2.2 Predictions of Reactive Scalars Two different sets of calculations using the P F and C M L models proposed in this thesis were performed and compared to the DNS [8] and C M C [9] reference databases. In addition to these, a third set of computations was performed using a generator of Gaussian random numbers that was used i n place of the P F and C M L models to mimic the chaotic variations of strain. This was done with the objective of determine whether or not a simple set of random numbers would provide similar results to those obtained with the periodic forcing and coupled map lattice models. The results obtained with this scheme are identified as D C M C - R A N in the following figures. Figure 4.8 expresses results i n terms of Favre averages which, for the case of species mass fractions, were calculated in a similar fashion as in the C M C formulation as f ^ lopWKpjQdc: fi \CP(Q<K P where P(Q is the P D F of mixture fraction which is represented by the /3-PDF taken from the DNS reference database [8]. In Fig. 4.8(a) it can be seen that the D C M C - P F and D C M C - C M L models provide predictions of the Favre average of oxidiser mass fraction i n excellent agreement with the D N S results with a noticeable improvement over the C M C results, while the D C M C - R A N model provides only a marginal improvement over the C M C predictions. The Favre-averaged methane mass fraction is also well predicted by both the D C M C - P F and D C M C - C M L models, as shown i n Fig. 4.8(b), with a tendency to under-predict, especially after 35 time units. Figures 4.8(c) and 4.9(a) show that the D C M C - P F and D C M C - C M L models provide good predictions of the Favre average and conditionally-averaged intermediate mass fraction after 30.0 time units, respectively. This result has important implications given that the intermediate includes CO. It is expected that i f the intermediate mass fraction is properly simulated, the ability to capture its effect in the overall generation of pollutants w i l l be positively affected. 49 Chapter 4. Simulation of Methane Oxidation with DCMC-SP 0.125 0.12 0.073 0.115 0.072 'fuel 0.11 0.105 0.071 -a+ 0 0.1 — DNS CMC DCMC-PF DCMC-CML DCMC-RAN 10 DNS 0.07 20 Time 30 40 0.069 50 -e+ 0 • CMC DCMC-PF DCMC-CML DCMC-RAN 10 (a) Oxidiser 20 Time 30 40 50 30 40 50 (b) Fuel , x 10 0.045 -e0.04 + 0 - DNS CMC DCMC-PF DCMC-CML DCMC-RAN 0.035 'int 2 'prod 0.03 4^"* 0.025 0.02 (c) Intermediates 10 20 Time (d) Products 650i 600 f 550 500 10 20 30 Time (e) NO 40 50 450 20 30 Time (f) Temperature Figure 4.8: D N S , C M C , and D C M C Favre averages of species mass fractions Chapter 4. Simulation of Methane Oxidation with DCMC-SP 50 T h e Favre-averaged mass fraction o f products is over-predicted for all times by b o t h D C M C - b a s e d simulations as shown i n F i g . 4.8(d). T h i s over-prediction, along w i t h the underp r e d i c t i o n o f the Favre average o f fuel mass fraction, is merely a result o f t h e overprediction of the reaction rates. (c) (T\Z) Figure 4.9: C o n d i t i o n a l averages at t = 30.0 P r o b a b l y the most significant achievement of this project is the good agreement between the D N S and D C M C - P F and D C M C - C M L predictions o f NO. Figure 4.8(e) shows a n i m - p o r t a n t improvement i n t h e predictions o f t h e Favre average o f NO mass f r a c t i o n over t h e C M C results. T h i s improvement is also evident i n t h e conditionally-averaged NO mass fract i o n predictions shown after 30.0 t i m e u n i t s i n F i g . 4.9(b). T h e good predictions of NO Chapter 4. Simulation of Methane Oxidation with DCMC-SP 51 mass f r a c t i o n are closely related t o the quality of the predictions o f temperature, since the Zel'dovich NO mechanism used t o calculate the NO mass fraction is strongly dependent on t e m p e r a t u r e . Figure 4.8(f) shows the Favre average of temperature, where i t can be seen t h a t the D C M C - P F a n d D C M C - C M L results show a considerably improvement i n the predicted results over the C M C results. T h i s result is f u r t h e r confirmed by F i g . 4.9(c). _x1CT 0.2 0.4 (c) n i n t __„x10^ 6 0.6 0.8 , DCMC-RAN Figure 4.10: a-variations i n the intermediates field at t = 30.0 Figures 4.10, 4.11, a n d 4.12 show comparisons between the variations around the condit i o n a l mean i n t h e direction o f a for different scalar fields at 30.0 t i m e u n i t s calculated w i t h the different models proposed i n this research. Here, Hi = Y} (a, £) — ( Y / | a ) . Chapter 4. Simulation of Methane Oxidation with 0.2 0.4 0.6 0.8 0.2 0.4 (b) U , (a) IIJVO, D C M C - P F NO 0.2 0.4 0.6 DCMC-SP 0.6 52 0.8 DCMC-CML 0.8 (c) I W , D C M C - R A N F i g u r e 4.11: a-variations i n the N O field at t = 30.0 From these figures, it can be seen that the variations simulated by the P F and C M L models have similar structures and magnitudes for all cases shown, while the D C M C - R A N implementation simulates variations one order of magnitude smaller in the same fields. For all cases, it is clear that the maximum and minimum variations occur in the vicinity of ( « 0.35, which is where the chemical reaction is assumed to occur with more intensity. While the P F and C M L show consistent results with this hypothesis, the R A N implementation shows large variations in the temperature field at a very different location, as shown in Fig. 4.12(c). This is due to the lack of structure in the tensor field calculated with the D C M C - R A N model. Chapter 4. Simulation of Methane Oxidation with DCMC-SP 53 T h e variations shown i n F i g . 4.12(a) a n d 4.12(b) seem t o indicate t h e presence o f ext i n c t i o n a n d re-ignition phenomena, characterised by 'cold' regions between ' h o t ' regions. I n t h e case o f t h e figures shown, these variations are rather small, b u t large enough t o produce changes i n the simulated values o f N O and intermediates, as shown i n F i g . 4.11(a), 4.11(b), 4.10(a), a n d 4.10(b). B o t h the P F a n d R A N simulations took approximately 2.7 hours t o complete i n a P 4 X e o n 2.4 G H z cluster using a single processor, while the C M L simulations took approximately 30.0 hours t o complete using t h e same hardware. xKf xlO" 3 (a) n , D C M C - P F (b) n , D C M C - C M L T T xlCf 0.2 0.4 „ 0.6 4 0.8 (c) U , D C M C - R A N T Figure 4.12: a-variations i n t h e temperature field at t = 30.0 3 Chapter 4. Simulation of Methane Oxidation with DCMC-SP 4.3 54 Summary T h e results obtained from the n u m e r i c a l i m p l e m e n t a t i o n of the models p r o p o s e d i n t h i s thesis have been presented i n this chapter. T h e most i m p o r t a n t observations are s u m m a r i s e d as follows: • T h e stochastic processes used to simulate the s t r a i n tensor proved t o offer good results i n terms of magnitude a n d structure o f the strain fields. T h e P F m o d e l requires at least 5 signals to offer consistent results, w h i l e the C M L m o d e l provides t h e structure a n d m a g n i t u d e r e q u i r e d for the test combustion process at a d i s s i p a t i o n l e v e l of fc = 9 • It was observed t h a t by u s i n g a second c o n d i t i o n a l variable a d d i t i o n a l variations c a n be s i m u l a t e d w i t h the D C M C m o d e l . These variations w o u l d have been ignored u s i n g the single c o n d i t i o n a l moment closure. • I n general, the D C M C - P F a n d D C M C - C M L models p r o v i d e d i m p r o v e d predictions of species mass fractions a n d temperature when compared w i t h the C M C m o d e l . This is especially evident i n the Favre-averaged predictions of intermediates a n d NO mass fractions and temperature. • W h i l e the c o m p u t a t i o n a l expense o f the D C M C - R A N is comparable to that o f the D C M C - P F m o d e l , it was observed t h a t the lack of structure inherent t o the R A N i m plementation results i n inconsistent simulations o f the scalar fields. 55 Chapter 5 Conclusions and Recommendations T h e a p p l i c a t i o n o f the D C M C m e t h o d w i t h stochastic processes to s i m u l a t e t u r b u l e n t comb u s t i o n has been explored throughout the development of this thesis. findings I n this chapter the p r o d u c e d w i t h this research are summarised a n d recommendations for future work are a r t i c u l a t e d . T h e a p p l i c a t i o n of stochastic processes and M o n t e C a r l o methods i n c o m b u s t i o n were described i n C h a p t e r 2, along w i t h the general fundamentals of turbulent c o m b u s t i o n . E x a m p l e s of relevant applications to t h i s research were also discussed i n C h a p t e r 2. T h e d e r i v a t i o n of the m a t h e m a t i c a l models proposed i n this thesis t o o k place i n C h a p t e r 3. D C M C transport equations of species mass fractions, temperature a n d scalar d i s s i p a t i o n were formulated a n d the unclosed terms were m o d e l l e d u s i n g two tools: scalar t r a n s p o r t a n d c o n d i t i o n a l moment closure hypotheses. T w o different stochastic processes were also proposed i n C h a p t e r 3 of t h i s thesis t o model the s t r a i n tensor i n the transport equation o f scalar dissipation. T h e n u m e r i c a l i m p l e m e n t a t i o n of the models developed i n C h a p t e r 3 showed encouragi n g results i n terms of simulations of the strain tensor, as w e l l as predictions of t e m p e r a t u r e a n d intermediates a n d NO mass fraction. T h e a d d i t i o n of the second scalar v a r i a b l e induced variations that cannot be s i m u l a t e d b y the single c o n d i t i o n a l moment closure m e t h o d . These results are discussed and summarised i n C h a p t e r 4. A n e x a m p l e o f the a p p l i c a b i l i t y a n d advantages of the D C M C m e t h o d w i t h stochastic processes i n s i m u l a t i o n of turbulent c o m b u s t i o n has been demonstrated w i t h the results obt a i n e d from this research. A s w i t h m a n y other projects, however, areas o f improvement have Chapter 5. Conclusions and Recommendations 56 also been identified that would make the model more robust and applicable to a more general range of problems. It is then proposed to: • Explore ways in which this or a similar model could work linked to a C F D program. W i t h this model it is possible to use the fluid flow data provided by a C F D program as a real-time input. The code would then calculate the strain field and the solution for Eq. 4.1, returning updated values of density, temperature, viscosity, etc. to the C F D code. • Investigate the effects of increased variations around the conditional mean of a. These variations could be induced by modifying the initial conditions and coupling between the scalar dissipation of a and Z. • Investigate the use of the models presented in this thesis including variations in space. This implies a different approach for the closure of Eqs. (3.20) and (3.21) to that presented in Chapter 3. It is possible that the differences between the results presented in this research and those found i n the D N S reference database could be due to the lack of a spatial variable. The use of such as a conditional variable will provide the model with a capability to discriminate among isosurfaces of mixture fraction and o. • Test a different case. As with most numerical models, it is important to verify the validity and applicability of the hypotheses used i n the development of the model. The combustion of Hydrogen-Oxygen and Hydrogen-Air flames has been studied extensively and could present itself as a good candidate for comparison purposes with the models presented in this thesis. 57 Bibliography [1] N . Peters, Turbulent Combustion, Cambridge University Press, 2000. [2] J . Warnatz, U . Maas, R. W . Dibble, Combustion. Physical and Chemical Fundamentals, Modeling and Simulation, Experiments, Pollutant Formation, Springer, 2001. [3] A . Y . Klimenko, Multicomponent diffusion of various admixtures in turbulent flow, Fluid Dynamics 25 (1990) 327-334. [4] R. W . Bilger, Conditional moment closure for turbulent reacting flow, Phys. Fluids 5 (2) (1993) 436-444. [5] C . Beck, Chaotic cascade model for turbulent velocity distributions, Physical Review E 49 (5) (1994) 3641-3652. [6] A . Hilgers, C . Beck, Coupled map lattices simulating fully developed turbulent flows, A I P Conference Proceedings 411 (1997) 11-17. [7] A . Hilgers, C . Beck, Hierarchical coupled map lattices as cascade models for hydrodynamical turbulence, Europhysics Letters 45 (5) (1999) 552-557. [8] W . K . Bushe, R . W . Bilger, G . R. Ruetsch, Incorporating realistic chemistry into direct numerical simulations of turbulent non-premixed combustion, in: Annual Research Briefs. Center for Turbulence Research, N A S A Ames/Stanford University, 1997, pp. 195-211. [9] W . K . Bushe, R . W . Bilger, Conditional moment closure modeling of turbulent nonpremixed combustion with reduced chemistry, in: Annual Research Briefs. Center for Turbulence Research, N A S A Ames/Stanford University, 1999, pp. 139-154. [10] D . Veynaute, L . Vervisch, Turbulent combustion modeling, Progress i n Energy and Combustion Science 28 (2002) 193-266. Bibliography 58 [11] F . A . Williams, Combustion Theory, Addison-Wesley Publishing Company, Inc., 1985. [12] N . Peters, Reduced kinetic mechanisms and assymptotic approximations for methane-air flames, Springer, 1991, C h . Reducing mechanisms, pp. 48-67. [13] R . W . Bilger, M . B . Esler, S. H . Starner, Reduced kinetic mechanisms and assymptotic approximations for methane-air flames, Springer, 1991, C h . O n reduced mechanisms for methane-air combustion, pp. 86-110. [14] H . Tennekes, J . L . Lumley, A First Course in Turbulence, The M I T Press, 1992. [15] A . N . Kolmogorov, The local structure of turbulence in incompressible viscous fluid for very large reynolds numbers, Dokl. Akad. Nauk. SSSR 30 (4) (1941) 3201. [16] W . K . Bushe, Conditional Moment Closure Methods for Autoignition Problems. P h D . Thesis, University of Cambridge, 1996. [17] B . E . Launder, D . B . Spalding, Lectures in Mathematical Models of Turbulence, Academic Press, 1972. [18] S. B . Pope, Turbulent Flows, Cambridge University Press, 2000. [19] T . B . Gatski, C . L . Rumsey, Closure Strategies for Turbulent and Transitional Flows, Cambridge University Press, 2002, C h . Linear and Nonlinear E d d y Viscosity Models, pp. 9-46. [20] W . K . Bushe, Combustion Lecture Notes, University of British Columbia, 2003. [21] W . E . M e l l , V . Nilsen, G . Kosaly, J . J . Riley, Investigation of closure models for nonpremixed turbulent reacting flows, Phys. Fluids 6 (3) (1994) 1331-1356. [22] N . Swaminathan, R . W . Bilger, Assessment of combustion submodels for turbulent nonpremixed hydrocarbon flames, Combustion and Flame 116 (1999) 519-545. [23] R. W . Bilger, Future progress i n turbulent combustion research, Prog. Energy Combust. Sci. 26 (2000) 367-380. [24] T . Echekki, C . J . H . , Unsteady strain rate and curvature effects in turbulent premixed methane-air flames, Combustion and Flame 106 (1995) 184-202. Bibliography 59 [25] H . Pitsch, H . Steiner, Large-eddy simulation of a piloted methane/air diffusion flame (sandia flame d), Phys. Fluids 12 (10) (2000) 2541-2554. [26] H . Steiner, W . K . Bushe, Large-eddy simulation of a turbulent reacting jet w i t h conditional source-term estimation, Phys. Fluids 13 (3) (2001) 754-769. [27] D . Laurence, Closure Strategies for Turbulent and Transitional Flows, Cambridge University Press, 2002, C h . Large Eddy Simulation of Industrial Flows?, pp. 392-406. [28] N . Peters, Laminar diffusion flamelet models in non-premixed turbulent combustion, E n ergy Combust. Sci. 10 (1984) 319-339. [29] D . Lentini, Assessment of the stretched laminar flamelet approach for nonpremixed turbulent combustion, Combust. Sci. and Tech. 100 (1994) 95-122. [30] A . W . Cook, J . J . Riley, G . Kosaly, A laminar flamelet approach to subgrid-scale chemistry in turbulent flows, Combustion and Flame 109 (1997) 332-341. [31] F . Mauss, D . Keller, N . Peters, A lagrangian simulation of flamelet extinction and reignition in turbulent jet diffusion flames, in: Twenty-Third Symposium (International) on Combustion. The Combustion Institute, 1990, pp. 693-698. [32] H . Pitsch, Unsteady flamelet modeling of differential diffusion in turbulent jet diffusion flames, Combustion and Flame 123 (2000) 358-374. [33] W . P . Jones, Closure Strategies for Turbulent and Transitional Flows, Cambridge University Press, 2002, C h . The Joint Scalar Probability Density Function Method, pp. 582-625. [34] A . Einstein, Uber die von der molekularkinetischen theorie der warme geforderte bewegung von in ruhenden fliissigkeiten suspendierten teilchen, Annalen der Physik 4 (17) (1905) 549-561. [35] A . Einstein, Investigations on the theory of the Brownian movement. E d . R . F u r t h , Tr. A . D . Cowper, Methuen, 1926. [36] M . von Smoluchowski, Zur kinetischen theorie der brownschen molekularbewegung und der suspensionen, Annalen der Physik 4 (21) (1906) 756-781. Bibliography 60 [37] P. Langevin, Sur la theorie du mouvement brownien, Comptes Rendus de l'Academie des Sciences (Paris) 146 (1908) 530. [38] J . B . Perrin, L a discontinuity de l a matiere, Revue du Mois 1 (1909) 323-344. [39] J . B . Perrin, Mouvement brownien et realite moleculaire, Annales de chimie et de physiqe V I I I (18) (1909) 5-114. [40] J . B . Perrin, Brownian Movement and Molecular Reality, Taylor and Francis, 1909. [41] J . B . Perrin, Les Atomes, Alcan, 1913. [42] J . B . Perrin, Atoms. T r . D . L . Hammick, Constable, 1970. [43] P. R . Kramer, A review of some monte carlo simulation methods for turbulent systems, Lecture Notes, Department of Mathematical Sciences, Rensselaer Polytechnic Institute . [44] F . W . J . Elliot, D . Horntrop, A . J . Majda, A fourier-wavelet monte carlo method for fractal random fields, J . Comp. Phys. 132 (2) (1997) 384-408. [45] K . K . Sabelfeld, Monte Carlo methods i n boundary value problems, Springer-Verlag, 1991. [46] A . McCoy, P h D Thesis, Department of Mathematics, University of California at Berkeley, 1975. [47] J . P. Minier, J . Pozorski, Wall-boundary conditions i n probability density function methods and applications to a turbulent channel flow, Phys. Fluids 11 (9) (1999) 2632-2644. [48] W . C . Welton, Two-dimensional pdf/sph simulations of compressible turbulent flows, J . Comput. Phys. 139 (2) (1998) 410-443. [49] D . J . Thomson, Criteria for the selection of stochastic models of particle trajectories i n turbulent flows, J . Fluid Mech. 180 (1987) 529-556. [50] J . Pozorski, J . P. Minier, Probability density function modeling of dispersed two-phase turbulent flows, Phys. Rev. E 59 (1) (1999) 855. [51] F . A . Jaberi, P. J . Colucci, S. James, P. G i v i , S. B . Pope, Filtered mass density function for large^eddy simulation of turbulent reacting flows, J . Fluid Mech. 401 (1999) 85-121. Bibliography 61 [52] A . Obieglo, J . Gass, D . Poulikakos, Comparative study of modeling a hydrogen nonpremixed turbulent flame, Combust. Flame 122 (2000) 176-194. [53] O. Le Maitre, O. M . Knio, H . N . Najm, R . G . Ghanem, A stochastic projection method for fluid flow: I. basic formulation, J . Comp. Phys. 173 (2001) 481-511. [54] S. B . Pope, A monte carlo method for the pdf equations of turbulent reactive flow, Comb. Sci. Tech. 25 (1981) 159-174. [55] E . E . O'Brien, The Probability Density Function (pdf) Approach to Reacting Turbulent Flows, Springer-Verlag, 1980. [56] R. G . Batt, Turbulent mixing i n a low speed shear layer., J . Fluid Mech. 82 (1977) 53. [57] L . Valifio, A field monte carlo formulation for calculating the probability denstity function of a single scalar i n a turbulent flow, Flow, Turb. Comb. 60 (1998) 157-172. [58] T . Hulek, R . P. Lindstedt, Modelling unclosed nonlinear terms i n a p d f closure for turbulent flames, Math. Comp. Modelling 24 (8) (1996) 137-147. [59] T. Hulek, L . R. P., Computation of steady-state and transient premixed turbulent flames using pdf methods, Comb. Flame 104 (1996) 481-504. [60] H . Kawanabe, M . Shioji, M . Ikegami, Cfd simulation for turbulent mixing w i t h the stochastic approach, Heat Transf. A . Res. 30 (6) (2001) 503-511. [61] S. M . Cannon, B . S. Brewster, L . D . Smoot, Stochastic modeling of co and no in premixed methane combustion, Comb. Flame 113 (1998) 135-146. [62] M . Kraft, F . Harald, Some analytic solutions for stochastic reactor models based on the joint composition pdf, Comb. Theory M o d . 3 (1999) 343-358. [63] C . C h a , G . Kosaly, H . Pitsch, Modeling extinction and reignition i n turbulent nonpremixed combustion using a doubly-conditional closure approach., Physics of Fluids 13 (12) (2001) 3824-3834. [64] N . B . Vargaftik, Tables on the Thermophysical Properties o f Liquids and Gases, John Wiley and Sons, 1975. Bibliography 62 [65] V . Eswaran, S. B . Pope, Direct numerical simulations of the turbulent mixing of a passive scalar, Phys. Fluids 31 (3) (1988) 506-520. [66] U . C . Miiller, N . Peters, Development of reduced reaction schemes for the ignition of diesel fuels i n a non-premixed turbulent flow field, IDEA-Programme Internal Report . [67] D . Thevenin, S. Candel, Effect of variable strain on the dynamics of diffusion flame ignition, Combust. Sci. and Tech. 91 (1993) 73-94. [68] E . Mastorakos, T . A . Baritaud, T . J . Poinsot, Numerical simulations of autoignition in turbulent mixing flows, Combust. Sci. and Technol. 125 (1997) 243-282. [69] G . R . Ruetsch, M . R. Maxey, Small-scale features of vorticity and passive scalar fields i n homogeneous isotropic turbulence, Phys. Fluids 3 (6) (1991) 1587-1597. [70] T . S. Cheng, J . A . Wehrmeyer, R . W . Pitz, O. Jarrett, G . B . Northam, Finite-rate chemistry effects in a mach 2 reacting flow, A I A A / A S M E / S A E 27th Joint Propulsion Conference . [71] H . Pitsch, S. Fedotov, Investigation of scalar dissipation rate fluctuations in non-premixed turbulent combustion using a stochastic approach, Combustion Theory and Modelling 5 (2001) 41-57. [72] S. A . Orszag, G . S. Patterson, Statistical models and turbulence, lecture notes in Physics, Springer, 1972, C h . Numerical simulation of turbulence. [73] V . Eswaran, S. B . Pope, A n examination of forcing in direct numerical simulations of turbulence, Comput. Fluids 16 (3) (1988) 257-278. [74] A . Hilgers, M . Gremm, J . Schnakenberg, A criterion for stochastic resonance, Physics Letters A 209 (1995) 313-316. [75] F . A . Williams, Reduced kinetic mechanisms and assymptotic approximations for methane-air flames, Springer, 1991, C h . Overview of asymptotics for methane flames. [76] N . Swaminathan, R . W . Bilger, Direct numerical simulation of turbulent nonpremixed hyrdocarbon reaction zones using a two-step reduced mechanism, Combust. Sci. Tech. 127 (1998) 167. Bibliography 63 [77] P. N . B r o w n , G. D . Byrne, A . C . H i n d m a r s h , Vode: A variable coefficient ode solver, S I A M J . Sci. Stat. C o m p u t . 10 (1989) 1038-1051. [78] W . K . Bushe, R. S. Cant, Results o f direct numerical simulations o f p r e m i x e d comb u s t i o n , i n : Cambridge University Engineering D e p a r t m e n t I n t e r n a l R e p o r t C U E D / A T H E R M O / T R 6 3 , 1996.
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Conditional moment closure for methane oxidation using...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Conditional moment closure for methane oxidation using two conditional variables and stochastic processes Lozada-Ramfrez, Jorges R. 2004
pdf
Page Metadata
Item Metadata
Title | Conditional moment closure for methane oxidation using two conditional variables and stochastic processes |
Creator |
Lozada-Ramfrez, Jorges R. |
Date Issued | 2004 |
Description | The conditional moment closure method using two conditioning scalar variables is applied to derive the transport equation of species mass fraction, temperature, and scalar dissipation in a decaying, isotropic, homogeneous turbulent methane-air flow. The strain tensor in the transport equation of scalar dissipation of the conditioning variables is simulated using stochastic processes. The results of this model are then compared to DNS and conditional moment closure with one variable for the same test case. |
Extent | 7020563 bytes |
Genre |
Thesis/Dissertation |
Type |
Text |
File Format | application/pdf |
Language | eng |
Date Available | 2009-11-24 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
DOI | 10.14288/1.0103622 |
URI | http://hdl.handle.net/2429/15624 |
Degree |
Master of Applied Science - MASc |
Program |
Mechanical Engineering |
Affiliation |
Applied Science, Faculty of Mechanical Engineering, Department of |
Degree Grantor | University of British Columbia |
Graduation Date | 2004-11 |
Campus |
UBCV |
Scholarly Level | Graduate |
Aggregated Source Repository | DSpace |
Download
- Media
- 831-ubc_2004-0541.pdf [ 6.7MB ]
- Metadata
- JSON: 831-1.0103622.json
- JSON-LD: 831-1.0103622-ld.json
- RDF/XML (Pretty): 831-1.0103622-rdf.xml
- RDF/JSON: 831-1.0103622-rdf.json
- Turtle: 831-1.0103622-turtle.txt
- N-Triples: 831-1.0103622-rdf-ntriples.txt
- Original Record: 831-1.0103622-source.json
- Full Text
- 831-1.0103622-fulltext.txt
- Citation
- 831-1.0103622.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.831.1-0103622/manifest