1D Flame Speed Analysis

Aim:

To find the flame speed using Cantera for 1D premixed CH4-Air and H2-Air combustions.

Introduction:

The flame speed is the measured rate of expansion of the flame front in a combustion reaction. In an internal combustion engine, the flame speed of a fuel is a property which determines its ability to undergo controlled combustion without detonation. Flame speed is used along with adiabatic flame temperature to help determine the engine's efficiency.

For this report, we assume a 1D premixed flame which is travelling in a constant pressure surrounding, i.e. open to atmosphere. It can be schematically represented as:

The analysis is steady state and the flame is positioned in the chamber such a way that unburned reactants come from left and products go out from right. The calculations have been carried out in a reference frame which has origin at flame front and positive x axis is along burnt products. The flame velocity is assumed zero and other velocities are calculated in relative to this speed. This transformation makes the transient phenomenon, a steady state one. Even if the calculations have been performed in flame front reference, the results are transformed back to original frame of reference with respect to ground.

The analysis is done at initial temperature of 300K and 1 atm constant pressure. The flame speed analysis is done for 3 cases: (For all cases, equivalence ratio is 0.9)

1) CH4 - Air flame (Using GRI 3.0 mechanism)

2) H2 - Air flame (using H2 mechanism)

3) H2 - Air flame (using GRI 3.0 mechanism)

About Cantera's free flame class:

It is a class in cantera which to used to perform analysis of 1D premixed free flowing laminar flame of any fuel-air mixture provided. We can define the grid points in the width of the domain. This flow can be solved by solving governing equations of mass conservation, energy, state equations and species diffusion equations. These all equations are built-in and we just need to ask the solver to run. It does numerical integration using newton-raphson solver to solve coupled differential equations. The solver first solves the domain based on basic grid and then it refines it based on slope and curve parameters of the solution variables. The refinement is done till the tolerance criteria is met. It does the computation using reference frame on flame front as mentioned previously.

How does code work?

1) The 'solution' class is passed to 'gas' object

2) Initial temperature, pressure and species mole fraction are specified.

3) The object 'f' is passed the 'FreeFlame' class. Along with it we pass the 'gas' object and 'width' parameter to 'FreeFlame'

4) The grid refinement parameters are set using 'set_refine_criteria' method. Here, we set 'ratio' as 3, 'slope' as 0.07, 'curve' as 0.14

'Ratio' is the parameter limiting the ratio of lengths of 2 adjacent elements.

'Slope' is the parameter limiting the change in values of solution variables for each node pairs.

'Curve' parameter limits the rate of change of slope.

5) Then we ask the solver to solve, with grid refinement enabled. Then we plot the results.

Code:

Case 1: CH4 - Air flame (Using GRI 3.0 mechanism)

import cantera as ct
import numpy as np
import matplotlib.pyplot as plt

gas = ct.Solution('gri30.xml', 'gri30_mix')
gas.TPX = 300, ct.one_atm, {'CH4':0.9, 'O2':2, 'N2':7.52}

f = ct.FreeFlame(gas,width = 0.03)
f.set_refine_criteria(ratio= 3, slope= 0.07, curve= 0.14)

f.solve(loglevel=1, refine_grid=True, auto=False)

plt.subplots()
plt.suptitle('1D Premixed Laminar Flame of CH4 air mixture (at 300 K and 1 atm initially)')
plt.subplot(2,2,1)
plt.plot(f.grid,f.T)
plt.title('Temperature')
plt.xlabel('Distance from left (m)')
plt.ylabel('Temperature (K)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.subplot(2,2,2)
plt.plot(f.grid,f.u)
plt.title('Velocity')
plt.xlabel('Distance from left (m)')
plt.ylabel('Velocity (m/s)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.subplot(2,2,3)
plt.plot(f.grid,f.int_energy)
plt.title('Internal energy')
plt.xlabel('Distance from left (m)')
plt.ylabel('Internal energy (J/kmol)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.subplot(2,2,4)
plt.plot(f.grid,f.heat_release_rate)
plt.title('Heat Release Rate')
plt.xlabel('Distance from left (m)')
plt.ylabel('Heat Release Rate (J/kmol)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.subplots()
plt.suptitle('Mole Fractions of species in 1D Premixed Laminar Flame of CH4 air mixture')
plt.subplot(2,4,1)
plt.plot(f.grid, f.X[gas.species_index('CH4'),:])
plt.title('CH4')
plt.xlabel('Distance from left (m)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.subplot(2,4,2)
plt.plot(f.grid, f.X[gas.species_index('O2'),:])
plt.title('O2')
plt.xlabel('Distance from left (m)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.subplot(2,4,3)
plt.plot(f.grid, f.X[gas.species_index('N2'),:])
plt.title('N2')
plt.xlabel('Distance from left (m)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.subplot(2,4,4)
plt.plot(f.grid, f.X[gas.species_index('CO2'),:])
plt.title('CO2')
plt.xlabel('Distance from left (m)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.subplot(2,4,5)
plt.plot(f.grid, f.X[gas.species_index('H2O'),:])
plt.title('H2O')
plt.xlabel('Distance from left (m)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.subplot(2,4,6)
plt.plot(f.grid, f.X[gas.species_index('CO'),:])
plt.title('CO')
plt.xlabel('Distance from left (m)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.subplot(2,4,7)
plt.plot(f.grid, f.X[gas.species_index('NO'),:])
plt.title('NO')
plt.xlabel('Distance from left (m)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.subplot(2,4,8)
plt.plot(f.grid, f.X[gas.species_index('OH'),:])
plt.title('OH')
plt.xlabel('Distance from left (m)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)


plt.show()

Case 2: H2 - Air flame (using H2 mechanism)

import cantera as ct
import numpy as np
import matplotlib.pyplot as plt

gas = ct.Solution('H2_mech.cti')
gas.TPX = 300, ct.one_atm, {'H2':1.8, 'O2':1, 'N2':3.76}

f = ct.FreeFlame(gas,width = 0.03)
f.set_refine_criteria(ratio= 3, slope= 0.07, curve= 0.14)

f.solve(loglevel=1, refine_grid=True, auto=False)

plt.subplots()
plt.suptitle('1D Premixed Laminar Flame of H2 air mixture (at 300 K and 1 atm initially) (using H2 mechanism)')
plt.subplot(2,2,1)
plt.plot(f.grid,f.T)
plt.title('Temperature')
plt.xlabel('Distance from left (m)')
plt.ylabel('Temperature (K)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.subplot(2,2,2)
plt.plot(f.grid,f.u)
plt.title('Velocity')
plt.xlabel('Distance from left (m)')
plt.ylabel('Velocity (m/s)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.subplot(2,2,3)
plt.plot(f.grid,f.int_energy)
plt.title('Internal energy')
plt.xlabel('Distance from left (m)')
plt.ylabel('Internal energy (J/kmol)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.subplot(2,2,4)
plt.plot(f.grid,f.heat_release_rate)
plt.title('Heat Release Rate')
plt.xlabel('Distance from left (m)')
plt.ylabel('Heat Release Rate (J/kmol)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.subplots()
plt.suptitle('Mole Fractions of species in 1D Premixed Laminar Flame of H2 air mixture (using H2 mechanism)')
plt.subplot(2,2,1)
plt.plot(f.grid, f.X[gas.species_index('H2'),:])
plt.title('H2')
plt.xlabel('Distance from left (m)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.subplot(2,2,2)
plt.plot(f.grid, f.X[gas.species_index('O2'),:])
plt.title('O2')
plt.xlabel('Distance from left (m)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.subplot(2,2,3)
plt.plot(f.grid, f.X[gas.species_index('N2'),:])
plt.title('N2')
plt.xlabel('Distance from left (m)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.subplot(2,2,4)
plt.plot(f.grid, f.X[gas.species_index('H2O'),:])
plt.title('H2O')
plt.xlabel('Distance from left (m)')
plt.grid(which='major')
plt.minorticks_on()
plt.grid(which='minor', alpha=0.2)

plt.show()

Case 3: H2 - Air flame (using GRI 3.0 mechanism)

The code for case 3 is same as of case 2, except the mechanism file is changed from 'H2_mech.cti' to 'gri30.xml','gri30_mix'

Results:

Case 1: CH4 - Air flame (using GRI 3.0 mechanism)

From all the graphs, it is evident that the flame front is between 0.010 and 0.011 metres. The temperature and velocity plot show unburnt reactants coming in from left of domain, at 0.25 m/s and 300K. When the reactants reach the flame front, the temperature rises up to 1850 K. From internal energy and heat release rate plots, we can say that the main reaction releases most of the heat and some minor reactions release the remaining heat. The products have higher temperature (2150 K) and velocity (3 m/s) due to these remaining minor reactions. This velocity of 0.25 m/s is the flame speed.

From species mole fractions plot, we can see reactants (CH4, O2) having a drop in mole fraction at flame point and products (CO2, H2O) having rise in mole fractions. The CO mole fractions first increase and then decrease. The reason for this first oxidiation of carbon atoms which gives CO and then CO forms CO2. The Ref [1] explains the reason behind curves of N2, NO and OH as Zeldovich mechanism, which talks about formation of NO from OH and N2.

The grid refinement was implemented till their were 155 nodes in the domain. The code stopped till this point because the refinement criteria was fullfilled.

Case 2: H2 - Air flame (using H2 mechanism)

The final temperature of flame is ~2275 K. The flame starts from 0.0105m and has flame speed of ~2.2 m/s. The internal energy plot increases by small amount and then decreases. This must be happening due to some reaction containing O2, because we can see in the plot below that oxygen mole fraction first increases and then decreases. The velocity and heat release rate is higher in H2 flame than CH4 flame.

The fuel H2 is seen to enter from left end at speed of ~2.2 m/s and its mole fraction then decreases to approx 0 when it gets oxidised. Similarly, we can see drop in mole fraction of O2 and rise in mole fraction of H2O. But there is unexplained increase in mole fraction of N2. This can happen due to reduction in total number of moles (of all species) which would mathematically be represented by rise in Mole fraction of N2. By looking at plots we can say that, the since H2 and O2 moles form H2O, there is apparent increase in moles of N2, but in principle the moles of N2 didn't increase at all.

The solver refined the grid till there were 122 points in grid.

Case 3: H2 - Air flame (using GRI 3.0 mechanism)

The results obtained using gri3.0 mechanism is also similar to results obtained from H2 mechanism file. The only visible effect is in velocity of flame. The gri3.0 mechanism file code gives velocity as 2.1 m/s. Also, the gri3.0 mechanism file code gives Carbon oxides as minor species, which are not present in H2 mechanism file code.

Conclusion:

The 1D premixed flame simulation is done and parameters like flame speed and temperature are found.

Answers to some questions:

1) What type of analysis is this, steady-state or transient?

--> It is a steady state analysis. The phenomenon is transient but the equations can be converted to steady state if we change the reference frame from ground to flame front.

2) Why do we need a separate file to simulate Hydrogen Mechanism can't it be done using GRI3.0?

--> The answers obtained using both mechanisms are nearly same but if we use GRI3.0 then we will have to solve for 325 reactions and 53 species, while if we use H2 mechanism, we need to solev 268 reactions and 57 species. The H2 mechanism is optimised for H2 oxidation and hence would give faster and accurate results compared to GRI3.0 mechanism.

3) What do you infer from the graphs, Explain in detail for both the graphs, So if you had to state one value as the 1D Flame which one would you say it is and Why?? What does the sudden increase in velocity or temperature tell us??

--> The sudden increase tells us the location of flame. The flame speed for CH4 flame is ~0.25 m/s while for H2 flame it is ~2.2 m/s. We can say this value as flame speed because the reactant and products can be distinguished based on presense of flame. Also the diagram at the start of report shows the structure of premixed flame. In it the flame is shown at the location where the temperature rises.

4) Substantiate your answer for the previous question by plotting variation in the concentration for "CO2" and explain what you infer from it (to find a function to find the change in concentration of CO2 find it using the help() command)

--> The concentration or mole fraction of CO2 can be seen in the plots, it shows the location where the combustion is complete or in other words the products are present.

5) Does the "width" parameter have an effect on the solution, If so or not so, please explain why?

--> The width parameter has an effect on the solution. The width is actually the length of domain. If we increase the width, we will have long domain and flame will be visible for small part of domain. If we decrease the width, the domain will be smaller and we will see a magnified image of flame. The flame thickness would stay the same in both cases. If the width is higher than required, more time is needed to compute solution. If width is smaller than required, the important flow physics will not be captured. The below images show the temperature plot for two different widths. The above stated points can be verified from these plots.

The flame can always be starting from 1/3 of domain length or width and the ending of flame depends on the flame thickness.