CFD Simulations AC2-01: Difference between revisions

From KBwiki
Jump to navigation Jump to search
No edit summary
No edit summary
Line 62: Line 62:
The code HybFlow is a finite volume based algorithm for numerical discretization of systems of partial differential equations. It refers to a strong conservative formulation of governing equations which are considered to have the generic structure of the system (1). Following the finite volume method idea the physical domain is subdivided into an ensemble of control volumes (cells) over which the balance laws (1) are locally imposed in integral form. Thanks to the Gauss theorem the integral balance is expressed as surface integrals for the convective and diffusive flux crossing the control volume boundary plus a net production term resulting from the integral of S over the cell extension. Expressing these quantities through approximate expressions, the average time variation of the solution vector is computed in each control volume. Two distinct approaches are used to this aim when convective and viscous fluxes are considered. Given two adjacent control volumes the flux integral over the common face is performed by a mid-point quadrature formula. For convective transport the face mid-point value used in the integration is obtained according to the evolution problem resulting from the interaction of the two adjacent states of the solution vector Q. Therefore the two different values of the solution stored on both cells divided by the common face are considered to interact following the approximate scheme proposed by Roe. The numerical fluxes computed by this scheme satisfy the upwind condition. Second order accuracy is obtained through a linear reconstruction of the solution inside each cell. In this way, the interaction over the common face considers the two states linearly reconstructed from both the control volumes and the numerical accuracy of the discretization method rises to second order.
The code HybFlow is a finite volume based algorithm for numerical discretization of systems of partial differential equations. It refers to a strong conservative formulation of governing equations which are considered to have the generic structure of the system (1). Following the finite volume method idea the physical domain is subdivided into an ensemble of control volumes (cells) over which the balance laws (1) are locally imposed in integral form. Thanks to the Gauss theorem the integral balance is expressed as surface integrals for the convective and diffusive flux crossing the control volume boundary plus a net production term resulting from the integral of S over the cell extension. Expressing these quantities through approximate expressions, the average time variation of the solution vector is computed in each control volume. Two distinct approaches are used to this aim when convective and viscous fluxes are considered. Given two adjacent control volumes the flux integral over the common face is performed by a mid-point quadrature formula. For convective transport the face mid-point value used in the integration is obtained according to the evolution problem resulting from the interaction of the two adjacent states of the solution vector Q. Therefore the two different values of the solution stored on both cells divided by the common face are considered to interact following the approximate scheme proposed by Roe. The numerical fluxes computed by this scheme satisfy the upwind condition. Second order accuracy is obtained through a linear reconstruction of the solution inside each cell. In this way, the interaction over the common face considers the two states linearly reconstructed from both the control volumes and the numerical accuracy of the discretization method rises to second order.


Having no monotonicity restrictions, the diffusive flux computation follows a more straightforward approach, which is equivalent to an unbiased centred scheme. In this case the mean value of to be used with the mid-point quadrature formula on each face is derived considering the average gradients of Q computed on the face itself. The value of is approximated using a finite difference formula between the right and left cell average values of the solution vector.
Having no monotonicity restrictions, the diffusive flux computation follows a more straightforward approach, which is equivalent to an unbiased centred scheme. In this case the mean value of [[Image:Image48.gif]] to be used with the mid-point quadrature formula on each face is derived considering the average gradients of Q computed on the face itself. The value of [[Image:Image71.gif]] is approximated using a finite difference formula between the right and left cell average values of the solution vector.


The steady state flow field is obtained solving the non-linear algebraic system of equations resulting from the spatial discretization. To this aim an implicit iterative Newton method is considered (Adami, 1998). Stability of the numerical algorithm is provided by a time-marching relaxation term resulting from the unsteady approximation in time of the governing equations. The matrix of the implicit method is computed numerically considering finite differences expressions of the residual vector elements with respect to the solution vector components. The resulting linear system is solved at each integration step by the iterative method GMRES (Saad 1994). To obtain an efficient convergence of the linear solution a right preconditioning is coupled with the iterative method. The preconditioning matrix is computed performing an incomplete ILU(0) factorisation of the implicit matrix (Saad, 1994). It is worthwhile to remember that the whole procedure GMRES-ILU(0) makes use of a condensed storage format considering for all the matrices involved only nonzero elements. Concerning the iterative time-marching strategy, the main Navier-Stokes equations are solved implicitly altogether as a 4x4 system while the transport equations for the turbulence model and for the conserved scalar field of the combustion model are solved in a uncoupled fashion. Therefore four separated and consecutively iteration steps are performed after the main flow equation update to march in time both the turbulent and scalar fields. The implicit iteration matrix is build up for every equation using the residual derivatives with respect to the current solution variables while the other components of the solution are kept frozen. The implicit approach neglects therefore any coupling existing between the different scalar fields of .
The steady state flow field is obtained solving the non-linear algebraic system of equations resulting from the spatial discretization. To this aim an implicit iterative Newton method is considered (Adami, 1998). Stability of the numerical algorithm is provided by a time-marching relaxation term resulting from the unsteady approximation in time of the governing equations. The matrix of the implicit method is computed numerically considering finite differences expressions of the residual vector elements with respect to the solution vector components. The resulting linear system is solved at each integration step by the iterative method GMRES (Saad 1994). To obtain an efficient convergence of the linear solution a right preconditioning is coupled with the iterative method. The preconditioning matrix is computed performing an incomplete ILU(0) factorisation of the implicit matrix (Saad, 1994). It is worthwhile to remember that the whole procedure GMRES-ILU(0) makes use of a condensed storage format considering for all the matrices involved only nonzero elements. Concerning the iterative time-marching strategy, the main Navier-Stokes equations are solved implicitly altogether as a 4x4 system while the transport equations for the turbulence model and for the conserved scalar field of the combustion model are solved in a uncoupled fashion. Therefore four separated and consecutively iteration steps are performed after the main flow equation update to march in time both the turbulent and scalar fields. The implicit iteration matrix is build up for every equation using the residual derivatives with respect to the current solution variables while the other components of the solution are kept frozen. The implicit approach neglects therefore any coupling existing between the different scalar fields of [[Image:Image72.gif]] .


Combustion model implementattion
Combustion model implementattion

Revision as of 16:31, 4 September 2008

Bluff body burner for CH4-HE turbulent combustion

Application Challenge 2-01 © copyright ERCOFTAC 2004


CFD Simulations

Conversely from simulation experience on piloted flames, the numerical prediction for the present bluff-body flame are not yet in satisfactory agreement with measurements. Computations of the velocity and turbulence fields are adequate in the up-stream regions of the jet covering the recirculation zone and the necking zone. However, further downstream and starting at about two bluff-body diameter, calculations show increasing deviation from the measurements. This is true regardless of the numerical approach used and of the modifications made to the model constants.

The chemistry model used for this application challenge has been so far based on one of the following assumptions:

1. flamelet combustion regime

2. fast chemistry

3. full or constrained equilibrium

Computations have been presented for temperature and mass fractions of OH and NO in the recirculation zone. The results allows to conclude that detailed chemical kinetics are needed to adequately compute the mass fractions of minor species such as OH or NO even in regions where local extinction phenomena are not so relevant as is the case of the recirculation zone.


Overview of CFD Simulations

The high costs involved in the development of gas turbine components pushed the producers to an increasing use of CFD as design tool. In fact, thanks to the deeper understanding in numerical simulation and to the increasing capability of modern computers, efficient and robust Navier-Stokes solvers are today available for practical application into the design process. At the same time the branch of research concerning the development of physical models for the simulation of realistic industrial conditions followed a less straight evolution and a general and satisfactory approach for turbulent regimes has not fully met yet. In case of turbulent reactive flows the complexity of the problem increases remarkably leading to further uncertainty in the physical model and consequently in the numerical predictions. For an accurate simulation of such flows a large system of differential equations need to be solved to compute, beside the classical conserved quantities, also the burning mixture mass composition. The costs involved in the solution at each grid point of this non-linear system of relationships become even more demanding because of the different scales appearing in the combustion process. Consequently a stiff set of equations has to be solved and the computational costs become extremely high even for a reduced chemistry mechanism. These complications contributed to limit the application of CDF as an effective tool for industrial combustor design. Combined with this also the difficulties in developing an efficient model for turbulent-reactive flows need to be considered responsible for the complex application of CFD in the analysis of industrial combusting systems. For turbulent regimes the commonly used approach refers to a statistical description based on the ensemble averaged formulation of conservation equations. Following the averaging, physical models are required to close the system of equations for turbulent transport of mass, momentum and energy.

Full predictions of the combustor flow field are desirable for a variety of reasons; of particular concern are the reduction of costs involved in the cut and tray testing phase, the increasing in durability or aerothermal performances of the gas turbine and the reduction of pollutants according to the more stringent emission restrictions. Despite the problems posed by the turbulent-combustion phenomena, the great importance of these points has motivated the application of CFD for the analysis of gas combustors. Two different strategies are possible to simplify the numerical simulation: the first one is to consider a simple flow with detailed chemistry while the second one computes complex flows with a reduced chemical mechanism. Considering this last approach, the early models proposed in literature are all derived from the ideas of the code developed at Imperial College by Gosman (1976). These make use of classical pressure correction schemes with a two equations k-ε model using wall functions for the near wall region. An extension of the proposed approach implementing a hybrid upwind scheme with a three level multigrid is described by Shyy et. al.(1987) for an industrial 3D turbofan engine combustor. The combustion model refers to the conserved scalar approach. A similar solver is described by Biswas et al. (1997) using a SIMPLE algorithm with a modified k-ε turbulence model to ensure accuracy for high swirling flows. In this case, a modified EBU model is considered for the combustion. Zheng et al.(1997) describe the numerical development of a 3D implicit multigrid solver for structured grids using a reduced mechanism based on 16 transport equations for the main species of the combustion of CH4. A third order monotone upwind scheme is used for all convection terms within a finite volume approach using a staggered grid technique. The total number of conservation relationship solved for turbulent flows is 23: the turbulence-chemistry interaction is based on an algebraic correlation to correct the reaction rate. Edwards and Roy (1998) suggest the use of a similar structured implicit multigrid solver for 2D problems. The main difference lays in the implementation of a low mach preconditioning of governing equations to yield an efficient approach at all speeds. The turbulent combustion closure considers the EDC (eddy dissipation concept) of Magnussen and Hjertager(1978). Knoll et al. (1996) describe a different numerical procedure which exploits a similar upwind finite volume scheme with a strong implicit Newton time marching method for acceleration of the convergence to the steady state. The matrix free Newton-Krylov method described in the work allows a fully implicit approach, which is supposed to overcome the numerical stiffness of the conservative governing equations for low speed reacting flows. A laminar-burning regime is assumed in the computation reported. The same laminar regime is considered also in the work of Hosangadi et al. (1996) and Ju (1995) describing the application of implicit TVD upwind schemes for unstructured/structured grids and high speed burning conditions.

The key aspect in modelling turbulent reactive flows lays in a realistic representation of the turbulence-chemistry interaction. Considering diffusion flames problems the process can be viewed as a competition between mixing and chemical reaction. Following this rough description, depending on the relative scales associated with the turbulence and the chemistry phenomena, different combustion regimes can be observed. The models used in the above mentioned applications refer generally to simple laminar closure or to semi-empirical closures derived from the EBU concept. A more satisfactory approximation considered for non-premixed flames could be obtained assuming the fast-chemistry burning regime. In this case, the turbulent time scales are presumed to be much higher than the time required by the evolution of chemical reactions. This regime of combustion constitutes a reasonable approximation for practical conditions encountered in gas turbines and making further assumption, such as equal molecular diffusivities and adiabatic flame conditions, it allows the thermodynamical state to be defined as a function of a single conserved scalar. The conserved scalar formalism represents a classical approach used for non-premixed flames in configurations involving the mixing of two inlet streams (Bilger, 1981). It produces a turbulent closure for the combustion through the assumption of a priori PDF distribution that represents the approximation of the conserved scalar statistics. Remarkable example of this procedure are reported by Jones (1982, 1994a, 1994b), Peters (1986) and Liew (1994).

The present work considers the results obtained from a recently developed reactive solver for internal turbomachinery applications (Adami 1998, 1999). The spatial discretization is based on a finite volume method for hybrid unstructured grids having a high geometrical flexibility of great concern in turbine combustors design. The spatial discretization is performed using a TVD unstructured upwind scheme with Roe’s approximate Riemann solver. The time-steady solution is computed through an implicit Newton algorithm coupled with a preconditioned Krylov based linear solver known as ILU(0)-GMRES(Saad, 1994). The turbulent-combustion closure is guaranteed by the conserved scalar approach and a Flamelet chemical model with a presumed b-PDF(Jones, 1994a). The turbulence is a two-equation k-w as proposed by Wilcox(1993). The preconditioning of physical equations described by Adami, 1999 is used to cope with low speed burning regimes allowing the pressure to be computed directly in place of density. This last is computed according to the perfect gas state law with the temperature obtained from the combustion model. Particular attention is here payed to the computation performed in the bluff-body configuration where advantage of the unstructured feature of the solver is considered.

The governing equation

The Navier-Stokes governing equations are considered for the density and velocity components (respectively). The turbulent model of the solver is based on a classical approach using the eddy-viscosity concept and the two-equation model as proposed by Wilcox, 1993. The turbulent quantities are represented by the turbulent kinetic energy and the turbulent dissipation rate (). The combustion model refers to the conserved scalar approach (Jones, 1994). In this regard, considering the fast chemistry assumption, non-premixed flames are viewed as composed by an ensemble of local planar laminar flames. Considering further conditions of adiabatic flow, equal species diffusivity and unit Lewis number, then a single scalar field is needed to describe the thermo-dynamic state within the flow field. In fact, thanks to the proposed simplifications and assuming a proper non-dimensional scaling based on the two inlet conditions for the air and the fuel streams, then the elements mass fraction and enthalpy balance equations describe the same physical process of a conserved scalar transported within the flow field by the convection and diffusion phenomena. A Flamelet approach is then used for computing the burning gas solution through the flame as a function of the scalar value. The turbulent combustion interactions are here accounted by a presumed PDF approach. In this model the gas state obtained from the flamelet approach is averaged using the local PDF of the conserved scalar which is defined using a beta-function with imposed average and variance values as computed from two additional transport equations. The relevant aspect of these two more differential equations, describing the average and variance of conserved scalar field, has to be considered the lack of chemical production/destruction terms. In this approach therefore the combustion actually uncouples from the flow computation and is described by the turbulent mixing of a scalar quantity which is physically conserved within the flow field regardless the combustion process.

The set of all the governing equations to be solved is:

Image46.gif



Here Q is the whole solution vector,Image47.gif and Image48.gif respectively the fluxes components for the convective and diffusive physical transport, while S is the source contribution term to the balance of physical quantities. All variables in the solution represent the mean flow state obtained through a conventional mass averaging of the governing equations. The flux Image47.gif and Image48.gif are expressed as follows:


Image49.gif


The scalar fields Image50.gif represent respectively the turbulence quantities Image51.gif , the conserved scalar average Image52.gif and its variance Image53.gif . Here Image54.gif are the total turbulent-laminar stresses obtained using respectively the effective viscosity Image55.gif while Image56.gif are the scalar Schmidt numbers for turbulent diffusive flux. The pressure field is computed considering the perfect gas state equation Image57.gif. The source terms for the turbulence-scalar equations are given by:

Image58.gif


where Image62.gif and Image63.gif . The constants appearing in the model are:

Image65.gif Image66.gif Image67.gif Image68.gif Image69.gif Image70.gif


The solver Hybflow

The code HybFlow is a finite volume based algorithm for numerical discretization of systems of partial differential equations. It refers to a strong conservative formulation of governing equations which are considered to have the generic structure of the system (1). Following the finite volume method idea the physical domain is subdivided into an ensemble of control volumes (cells) over which the balance laws (1) are locally imposed in integral form. Thanks to the Gauss theorem the integral balance is expressed as surface integrals for the convective and diffusive flux crossing the control volume boundary plus a net production term resulting from the integral of S over the cell extension. Expressing these quantities through approximate expressions, the average time variation of the solution vector is computed in each control volume. Two distinct approaches are used to this aim when convective and viscous fluxes are considered. Given two adjacent control volumes the flux integral over the common face is performed by a mid-point quadrature formula. For convective transport the face mid-point value used in the integration is obtained according to the evolution problem resulting from the interaction of the two adjacent states of the solution vector Q. Therefore the two different values of the solution stored on both cells divided by the common face are considered to interact following the approximate scheme proposed by Roe. The numerical fluxes computed by this scheme satisfy the upwind condition. Second order accuracy is obtained through a linear reconstruction of the solution inside each cell. In this way, the interaction over the common face considers the two states linearly reconstructed from both the control volumes and the numerical accuracy of the discretization method rises to second order.

Having no monotonicity restrictions, the diffusive flux computation follows a more straightforward approach, which is equivalent to an unbiased centred scheme. In this case the mean value of Image48.gif to be used with the mid-point quadrature formula on each face is derived considering the average gradients of Q computed on the face itself. The value of Image71.gif is approximated using a finite difference formula between the right and left cell average values of the solution vector.

The steady state flow field is obtained solving the non-linear algebraic system of equations resulting from the spatial discretization. To this aim an implicit iterative Newton method is considered (Adami, 1998). Stability of the numerical algorithm is provided by a time-marching relaxation term resulting from the unsteady approximation in time of the governing equations. The matrix of the implicit method is computed numerically considering finite differences expressions of the residual vector elements with respect to the solution vector components. The resulting linear system is solved at each integration step by the iterative method GMRES (Saad 1994). To obtain an efficient convergence of the linear solution a right preconditioning is coupled with the iterative method. The preconditioning matrix is computed performing an incomplete ILU(0) factorisation of the implicit matrix (Saad, 1994). It is worthwhile to remember that the whole procedure GMRES-ILU(0) makes use of a condensed storage format considering for all the matrices involved only nonzero elements. Concerning the iterative time-marching strategy, the main Navier-Stokes equations are solved implicitly altogether as a 4x4 system while the transport equations for the turbulence model and for the conserved scalar field of the combustion model are solved in a uncoupled fashion. Therefore four separated and consecutively iteration steps are performed after the main flow equation update to march in time both the turbulent and scalar fields. The implicit iteration matrix is build up for every equation using the residual derivatives with respect to the current solution variables while the other components of the solution are kept frozen. The implicit approach neglects therefore any coupling existing between the different scalar fields of Image72.gif .

Combustion model implementattion

The straightforward application of the solver described for low speed reacting flows proved to be numerically not enough stable for all the investigated conditions. The source for this poor robustness has been attributed to the continuity equation, which is computed directly solving for the density field. This classical compressible approach has therefore been abandoned considering a preconditioning scheme for the governing equations derived from the artificial compressibility method suggested by Chorin (Adami, 1999). The preconditioned formulation refers to the primitive variables and allows the direct computation of the static pressure instead of density from the continuity equation. Owing to this property, the gas state equation is used to compute the density field given the temperature profile resulting from the scalar-Flamelet model. The solving algorithm resembles very closely the approach usually implemented for low speed reacting flows based on pressure correction schemes. The architecture of the solver is drawn in Fig. 1.

The Flamelet model approach requires the computation of the scalar thermodynamic state and as a function of the local fuel mass fraction (the conserved scalar). To provide a computationally efficient scheme the flamelet solution is firstly expressed using a polynomial fit of the curves obtained through the flame front. Given the two inlet reactants composition temperatures and pressure, the laminar diffusion flame is here computed numerically for a 1D-counteflowing model using a commercial package based on the Chemkin routines. The resulting profiles of the solution are then interpolated and integrated analytically using the beta PDF function assuming as parameters the mean and variance of the conserved scalar.

Figure 1: code architecture

Finally the coefficient of the polynomial obtained from the integration are stored into a look-up table which, being accessed during the CFD simulation with the local values of , allows a fast computation of the thermodynamic average quantities inside every grid cell. © ERCOFTAC 2004 CFD results and Comparison against Experiments

The bluff-body burner here considered has been suggested and investigated experimentally by Masri, 1985. It is centred in a coflowing stream of air and consists of a circular bluff-body with an orifice at its centre for the main fuel Fig. 14. The diameters are 3.6 mm for the fuel nozzle and 50mm for the bluff-body. The main fuel jet composition is 50% CH4 and 50% H2 for a bulk Jet velocity of 118 m/s. The outer flowing air has a 40m/s velocity. The jet flame can be extinguished in a region downstream the recirculation zone where turbulence is well developed and the finite rate chemistry effects are significant. The flame may also reignite further downstream where turbulent mixing rates are relaxed. It should therefore be noted that, as observed by Masri, both the piloted and the bluff-body jet flames discussed here consist, generally, of three main zones: stabilization, extinction and reignition zones. The computation of this complex flow structure has been here addressed considering a more flexible unstructured grid. As for the pilot-jet flame, the axial symmetry of the problem is considered solving a single tangential sector of the cylindrical flow.

Figure 9:zoom of the unstructured mesh at inlet section


Figure 10: measured and computed recirculation streaklines

In this case, the elements on the plane section (Fig. 9) are triangles, which allow an easy clustering close to the bluff-body walls and especially at the sharp leading edges with the fuel and air jets. Being D the bluff-body external diameter, a total number of 11000 elements have been used for the discretization grid which has a 2.5xD radial dimension while the outlet plane is placed 5xD downstream the inlet section. The boundary conditions are specified with the same approach described in the previous test. The Flamelet database is obtained from the same reduced mechanism of the previous application.

Figure 11a:axial velocity Y/D=0.6 Figure 11b: axial velocity Y/D=0.9

Figure 11c: axial velocity Y/D=1.3

The steady convergence of the computation has not been obtained for the present case. In fact going on with the computation after the residual level of 10E-6 is firstly reached, then a cyclic behaviour of the convergence is observed. This effect has been associated to a pseudo-unsteady evolution of the flow field behind the bluff-body, which is apparently caused by the numerical discretization. The recirculation zone seems to be not numerically stable and periodically detaches from the bluff-body moving downstream towards the outlet section. Beside the downstream motion of this main recirculation vortex, a new one starts to develop at the bluff-body wall restarting the original flow pattern. It has to be pointed out that the present computation is not time-accurate. Moreover since no mention has been done of a similar phenomenon in the experiments results then this behaviour has to be considered a spurious numerical effect. Despite this fact, the minimum value of 10E-6 for the residual is reached when a single recirculation vortex is attached to the bluff-body wall as observed in the field from the experimental investigation. All the data shown in the following therefore refers to this condition.

In Fig. 11, the velocity profiles are compared with experiments for three cross sections after the inlet. The accordance with measures is reasonably accurate showing that the main recirculation is correctly predicted by the numerical computation. This is qualitatively shown also comparing the flow pattern of Fig. 10. From computation it is also noticeable a small second vortex close to the fuel jet core as mentioned by Masri. The accurate simulation of the flow pattern reflects at the same instance a realistic representation of the turbulence field and therefore of the mixing properties behind the bluff-body which are to be expected in good agreement with experiments. As well known the conserved scalar approach is highly dependent on the turbulent mixing experienced by the reactants in the flow. The realistic simulation of this feature allows the correct prediction of the flame evolution and the accurate placement of the stoichiometric front. This can be observed in the following figures where the conserved scalar and temperature profiles are compared against experiments at several cross sections. Starting from the first section plane at Y/D=0.3 after the inlet, a peak is observed in the T profile inside the outer edge of the hot recirculation zone. This peak seems not present in the experiments and it is probably caused by the entrainment of unburned air flowing close to the bluff-body wall from the coflow stream directly into the main vortex core (Fig. 10).

Figure 12a: T profile Y/D=0.3 Figure 12b:scalar profile Y/D=0.3

Figure 13a: T profile Y/D=0.6 Figure 13b:scalar profile Y/D=0.6

Figure 14a: T profile Y/D=0.9Figure 14b:scalar profile Y/D=0.9

Figure 15a: T profile Y/D=1.3Figure 15b:scalar profile Y/D=1.3

Figure 16a: T profile Y/D=1.8Figure 16b:scalar profile Y/D=1.8

Figure 17a: T profile Y/D=2.4Figure 17b:scalar profile Y/D=2.4

Figure 18a: T axial r/D=0.12 Figure 18b: T axial r/D=0.8

This effect slightly reduces the mixture fraction concentration observed inside the vortex and for r»18 the presence of more air enhance the combustion producing the peak observed in the temperature. At the same instance the vortex outer branch seems to be more stretched by the entrainment of more fresh air in the core and again an overshoot of scalar mixture fraction, this time due to the excess of fuel, is observed from r»20 to r»24. Considering downstream sections the agreement with experiments is generally satisfying although a noticeable lower value of temperature is seen on the flame centreline for Y/D=1.3. Here a higher velocity of the fuel jet (observed from Fig. 11c) still maintains an efficient convective transport for the scalar field, which attains on the centreline a richer fuel-air ratio. This lag in the fuel jet spreading is probably due to the turbulence model response and it is soon recovered by the computation on the following sections as shown from the improved agreement of Fig. 14 and 18. A similar effect was above reported for the pilot-jet flame computation and seems to arise also for the present test concerning the turbulent field behaviour. In fact, although here less pronounced, the computation soon recovers from the under prediction of the scalar diffusion with a strong turbulent production which finally leads to a downstream excess of diffusion. This overdiffusion is clearly shown in Fig. 18 from the axial profile of the flame downstream about Y/D=120. In this case, the increased diffusion slightly reduces the level of the scalar mixture fraction crossing the flame front. Conversely, since quite lower turbulence production is expected within the fuel jet core then an accurate agreement is observed in Fig. 18 close to the flame centreline. This turbulence lag followed by a steep overproduction, although not pronounced as for the previous test, seems to be an intrinsic characteristic of the adopted turbulence model which becomes more critical when the grid is not enough refined and accurate as in the pilot-jet computation and as it starts to be for the section of Fig. 15. It has to be reminded that the large density fluctuations due to the combustion are here accounted in the flow computation simply through the density weighted averaging and no modifications were included into the original turbulent model suggested by Wilcox for inert conditions. Also the application of the eddy-viscosity approach becomes incapable of reproducing some important aspects of the physical flow such as the variable density-mean pressure gradients influence. Despite this evidence, the temperature and scalar profiles resulting from the present computation generally compare accurately with experiments.


Figure 19: radial distribution for species mass fraction, Y/D=0.3

The average levels of main species are reported in Fig. 19, 20 and 21. A good agreement is shown for the O2 and H2 mass fractions in all sections considered. Conversely, a well-pronounced underproduction of CO and CO2 levels is characterising all the numerical prevision. The profiles are reproduced but the mass fraction levels are well below the experiments. This result for CO is to some extend surprising as it is known (Jones, 1994) that the fast chemistry assumption generally overestimates CO emissions especially on the rich side of the flame front.

This chemical mechanism is considered enough accurate for computing CO levels and therefore the source of the inaccuracy observed is more realistically though to be lying in the non-equilibrium effects associated with CO oxidation. In fact, CO is particularly dependent on the local stretch applied by the turbulent field on the laminar flame front.

Figure 20: radial distribution for species mass fraction, Y/D=0.9

Since an accurate prevision has been obtained for both the mixture fraction and the temperature of the flame then this effect has to be attributed to the chemical mechanism and to the Flamelet database computation. As far as the chemical reaction a 22-steps mechanism of CH4 involving 20 species has been considered.

Figure 21: radial distribution for species mass fraction, Y/D=2.4

The turbulent-chemistry interaction is only approximately accounted in the present model through a constant flame stretch for the whole field. This parameter is enforced during the counterflowing Flamelet computation imposing the relative reactants velocity. An increasing value for the flame stretch reduces considerably the amount of equilibrium species found at high temperatures. This explains why for CO the under prediction appears more relevant in the first sections and generally where higher temperatures are found in the hot burning gases. Having no priori information about the level of stretching to be expected in the flow field then it has been guessed using clearly a too high value. Although this effect has a little impact on the temperature profiles, it has a strong influence for CO and CO2 prediction.

Finally, in Fig. 22 the computed fields of temperature, mixture fraction and turbulent kinetic energy are reported showing the overall flame structure.


Figure 22: computed contours of T, and k


The numerical CFD solver HybFlow has been briefly described for simulation of internal turbulent reacting flows. Three application of the numerical procedure have been here described for turbulent non-premixed flames. General agreement with experimental data has been obtained for the turbulent CH4 bluff body flame experimentally studied at SANDIA although some discrepancies where observed in the temperature profiles. Apparently, the reason for these inaccuracies seems to be due to the strong response of the k-w turbulence model in the burning shear layer surrounding the fuel main jet. The low mach preconditioning approach of the solver performed efficiently for this particular flow condition where inlet velocity profiles have a fairly different magnitude that are ranging from about 40 m/s in the jet down to less than 1 m/s in the coflow air.

An accurate comparison against experiments has been performed for the bluff-body burner of Masri, 1985. In this test, the unstructured grid allowed a flexible disposition of the elements especially close the inlet streams and in the main recirculating vortex entraining the hot burning gases. This unstructured mesh features proved to be successfully allowing an accurate computation of the flow pattern developing behind the bluff-body as also of the mixture fraction and temperature profiles. Concerning the species mass fraction the miss prediction observed for CO and CO2 were attributed to an incorrect overestimation of the turbulent flame stretch within the Flamelet database assembly. In this regard, equilibrium phenomena resulted of considerable importance for the computation of CO emissions.

Further developments are still needed for a deeper analysis and understanding of the turbulence model performances for reacting flows. More computation will be then performed extending the application of the unstructured meshes which are promising greater flexibility and efficiency for the analysis of effective internal gas combustor geometries. © ERCOFTAC 2004 References

Adami P., Michelassi, V., Martelli, F., 1998 “Performanches of a Newton-Krylov scheme against implicit and multi-grid solvers for inviscid flows” AIAA paper 98-2429.

Adami, P., 1999 “Numerical Computation of Turbulent Non-Premixed Reacting Flows in Combustion Chambers” 7th IGTC congress, Nov.99, Kobe

Barlow and Frank 1998 “Effects of Turbulence on Species Mass Fractions in Methane/Air Jet Flames” 27th Combustion Symposium, pp. 1087-1095

Bilger, R., W., 1981 “Turbulent Flows with nonpremixed reactants” In Turbulent Reacting Flows, Ed. P.A. Libby and F.A. Williams, Springer Verlag.

Biswas, D., Kawano, K., Iwasaki, H., Ishizuka, M., Yamanaka, S., 1997 "Three-dimensional computation of gas turbine combustors and the validation studies of turbulence and combustion models" ASME 97-GT-362, June.

Edwards, J.R., Roy, C.J., 1998 “ Preconditioned multigrid methods for two-dimensional combustion calculations at all speeds” AIAA J. Vol. 36, No. 2, Feb.

Gosman A.D., Ideriah, F.J.K., 1976 “TEACH-T: A general computer program for two-dimensional, turbulent recirculating flows”, Imperial College Rep., London.

Hosangadi, A., Lee., R.A., York, B.J., Sinha, N., Dash, S.M., 1996 “ Upwind unstructured scheme for three dimensional combusting flows”, J. Of Prop. And Power, Vol. 12, No. 3, May.

Jones W, P., 1994 “Turbulence modelling and numerical solution methods for variable density and combustiong flows” In Turbulent Reactive Flows, Ed. P.A. Libby and F.A. Williams, Academic Press.

Jones W. P., Whitelaw, J.H., 1982 “Calculation methods for reacting turbulent flows: A review” Comb. & Flame, 48, 1.

Jones W.P., Kakhi, M., 1994 “Mathematical modelling of turbulent flames” In M.V. Heitor, F. Cullick and J.H. Whitelaw Ed., Unsteady combustion Kluwer Academic Publishers.

Ju, Y., 1995 “Lower-upper scheme for chemically reacting flow with finite rate chemistry” AIAA J., Vol. 33, No. 8, Aug.

Knoll, D.A., McHugh, P.R., Keyes, D.E., 1996 “Newton-krylov methods for low mach number compressible combustion”, AIAA J. Vol. 34, No. 5, May.

Liew, S.K., Bray, K.N.C., Moss, J.B., 1994 "A stretched Laminar flamelet model of turbulent nonpremixed combustion" Comb. and Flame, 56:199-213.

Magnussen, B.F., Hyertager, B.H., Olsen, J.G., Bhaduri, D., 1978 "Effects of turbulent structure and local concentrations on soot formation and combustion in C2H2 diffusion flames" 17th Symposium on Combustion, The Combustion Institute.

Masri, A.R. and Bilger, R.W., 1985, `Turbulent Diffusion Flames of Hydrocarbon Fuels Stabilised on a Bluff Body', Twentieth Symposium (International) on Combustion. The Combustion Institute, Pittsburgh, pp. 319-326.

Masri, A.R., Dally, B.B., Barlow, R.S. and Carter, C.D., 1994, `The Structure of The Recirculation Zone of a Bluff-Body Combustor', Twenty-fifth Symposium on Combustion, The Combustion Institute, Pittsburgh, pp.1301-1308.

Peters, N., 1986 "Laminar flamelet concepts in tubulent combustion" 21th Symposium on combustion, The combustion Institute, pp. 1231-1250.

Saad, Y., 1994 "Krylov Subspace Thechniques, Coniugate Gradients, Preconditioning and Sparse Matrix Solvers", CFD VKI LS 1994-05 VonKarman Institute for Fluid Dynamics.

Shyy, W., Braaten, M.E., 1987 “A numerical study of flow in gas-turbine combustor” AIAA-87-2132 San Diego California.

Wilcox D.C., 1993 “Turbulence modelling for CFD” DWC Industries, Inc., La Canada.

Zeng, X., Liao, C., Liu, Z., Liu, C., 1997 “Mass-flux-based implicit multigrid method for modelling multidimensional combustion” J. Of Prop. And Power, Vol. 13, No. 1, Jan © copyright ERCOFTAC 2004

Contributors: Elisabetta Belardini - Universita di Firenze

Site Design and Implementation: Atkins and UniS

       Top        	      Next