Week 3 - Adiabatic Flame Temperature calculation
Introduction:
This report is a summary of the learning outcome of the first three weeks of the course to use python and Cantera for computational combustion applications. The key features of this report will be to demonstrate the knowledge on :
- Fuel, and its types and adiabatic flame temperature of combustion
- Use of python to calculate the adiabatic flame temperature using Newton-raphson method.
- Use of Cantera to calculate chemical equilibrium and adiabatic flame temperature.
The following questions will be addressed for the above-mentioned purpose,
1. Adiabatic flame temperature variation with Equivalence ratio (In constant volume-chamber)
2. Adiabatic flame temperature variation with Heat loss (In constant pressure chamber)
3. Adiabatic flame temperature variation with respect to Alkane, Alkene, Alkyne
4. Adiabatic flame temperature variation with respect to a number of Carbon atom chains
Coding approach & results :
Question 1:
The general equation for stoichiometric combustion of methane was formulated.
For an adiabatic process (at a constant P,V), the total enthalpy of products minus the total enthalpy of the reactants must be zero at T-adiabatic. This enthalpy function has been solved for multiple temperatures starting with an initial guess (T-guess), until the enthalpy function equates to zero. The solution procedure used was newton raphson technique and the GRI 3.0 data has been used to calculate the enthalpy of the species. This process was performed in a loop for different lambda values and the calculated values has been compared with that of Cantera's calculation. It can be seen from the results that the difference between the calculated values and Cantera's output at lambda = 1 is the highest. this is due to the fact that accounts for dissociation and maximum dissociation can only be observed at the highest temperature which in this case is at lambda = 1.
"""
Program to calculate Adiabatic flame temperature of methane combustion
Author : Logesh Paramashivam
Date : 26/05/2021
"""
"""
Stoichiometric reaction
CH4 + 2 (O2 + 3.76 N2) = CO2 + 2 H2O + 7.52 N2
Generic reaction
CH4 + lambda * 2 *(O2 + 3.76 N2) =
Assume for fuel rich products: O2 = 0, CO2 = ...value, CO = ...value
H2O = 2
CO = 4 - (4*lambda)
N2 = lambda * 2 * 3.76
CO2 = (4 * lambda)- 3
for fuel lean products: CO = 0, O2 = ...value
N2 = lambda * 2 * 3.76
CO2 = 1
H2O = 2
O2 = (2 * lambda) - 2
"""
import matplotlib.pyplot as plt
import math
import cantera as ct
def h(T, coeff): # Function to calculate enthalpy of a specific gas (inputs --> temperature and GRI 3.0 coefficients)
R = 8.314 # J/k.mol
a1 = coeff[0]
a2 = coeff[1]
a3 = coeff[2]
a4 = coeff[3]
a5 = coeff[4]
a6 = coeff[5]
# function returns the value of the entalpy for respective co efficients and temperature
return (a1 + a2*T/2 + a3*pow(T,2)/3 + a4*pow(T,3)/4 + a5*pow(T,4)/5 + a6/T) * R * T
def f(T):
# Calculation of enthalpy of each species in the products
h_co2_p = h(T, co2_coeff_h) # calling function h, passing respective (high temp) co efficients and temperature
h_h2o_p = h(T, h2o_coeff_h)
h_n2_p = h(T, n2_coeff_h)
h_co_p = h(T, co_coeff_h)
h_o2_p = h(T, o2_coeff_h)
# Calculation of enthalpy of each species in the reactants
Tstd = 298.15 # kelvin # standard temperature of reactants before start of reaction
h_n2_r = h(Tstd, n2_coeff_l) # calling function h, passing respective (low temp) coefficients and initial temperature
h_o2_r = h(Tstd, o2_coeff_l)
h_ch4_r= h(Tstd, ch4_coeff_l)
# Total enthalpy of reactants
H_rectants = h_ch4_r + 2*lam*h_o2_r + lam*2*3.76*h_n2_r
#print(H_rectants)
# Total enthalpy of products
H_products = (n_co2*h_co2_p) + (n_h2o*h_h2o_p) + (n_n2*h_n2_p) + (n_o2*h_o2_p) + (n_co*h_co_p)
#print(H_products)
# At Tadiabatic, H_products - H_reactants = 0 i.e, f(Tad) = 0
return H_products - H_rectants
def fprime(T):
return (f(T+1e-6) - f(T))/1e-6 # Derivative of f(x)
T_ad = []
iter_ct = []
gas = ct.Solution('gri30.xml')
T_cant = []
for i in range (0, 12):
lam_array = [0.9, 1, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0]
lam = lam_array[i]
n_ch4_r = 1
n_o2_r = 2 * lam
n_n2_r = 7.52 * lam
if lam > 1:
n_h2o = 2
n_co2 = 1
n_n2 = lam * 2 * 3.76
n_o2 = (2* lam)-2
n_co = 0
elif lam == 1:
n_h2o = 2
n_co2 = 1
n_n2 = lam * 2 * 3.76
n_o2 = 0
n_co = 0
else : #lam < 1:
n_h2o = 2
n_co2 = (4*lam - 3)
n_n2 = lam * 2 * 3.76
n_o2 = 0
n_co = 4-4*lam
# Coefficients a1...a7 for enthalpy calculation
# High temperature coefficients T>1000
co2_coeff_h = [3.85746029E+00, 4.41437026E-03, -2.21481404E-06, 5.23490188E-10, -4.72084164E-14, -4.87591660E+04, 2.27163806E+00]
h2o_coeff_h = [3.03399249E+00, 2.17691804E-03, -1.64072518E-07, -9.70419870E-11, 1.68200992E-14, -3.00042971E+04, 4.96677010E+00]
n2_coeff_h = [0.02926640E+02, 0.14879768E-02, -0.05684760E-05, 0.10097038E-09, -0.06753351E-13, -0.09227977E+04, 0.05980528E+02]
o2_coeff_h = [3.28253784E+00, 1.48308754E-03, -7.57966669E-07, 2.09470555E-10, -2.16717794E-14, -1.08845772E+03, 5.45323129E+00]
co_coeff_h = [2.71518561E+00, 2.06252743E-03, -9.98825771E-07, 2.30053008E-10, -2.03647716E-14, -1.41518724E+04, 7.81868772E+00]
# Low temperature coefficients T<= 1000
o2_coeff_l = [3.78245636E+00, -2.99673416E-03, 9.84730201E-06, -9.68129509E-09, 3.24372837E-12, -1.06394356E+03, 3.65767573E+00]
n2_coeff_l = [0.03298677E+02, 0.14082404E-02, -0.03963222E-04, 0.05641515E-07,-0.02444854E-10, -0.10208999E+04, 0.03950372E+02]
ch4_coeff_l = [5.14987613E+00, -1.36709788E-02, 4.91800599E-05, -4.84743026E-08, 1.66693956E-11, -1.02466476E+04, -4.64130376E+00]
T_guess = 1200 # Initial guess for T adiabatic
aplha = 0.6
tol = 1e-6
ct = 1
while (abs(f(T_guess)) > tol):
T_guess = T_guess - aplha*(f(T_guess)/fprime(T_guess))
ct = ct + 1
T_ad.append(T_guess)
iter_ct.append(ct)
gas.TPX = 298.15, 101325,{'CH4':n_ch4_r, 'O2':n_o2_r, 'N2':n_n2_r}
gas.equilibrate('HP','auto')
T_cant.append(gas.T)
print('Adiabatic flame temperature')
print(T_ad)
print('Iterations')
print(iter_ct)
print('cantera')
print(T_cant)
plt.figure()
plt.plot(lam_array, T_ad, label = 'Calculated')
plt.plot(lam_array, T_cant, label = 'Cantera')
plt.legend()
plt.xlabel('Lambda [-]')
plt.ylabel('Temperature [K]')
plt.title('Adiabatic flame temperature as a function of lambda')
plt.show()
Problem 2:
In this problem, the AFT calculation methodology using the newton raphson technique has been used. In this case, an element for heat loss has been introduced to simulate any heat losses in the reactor. This heat loss will be a fraction of the total enthalpy of the used fuel (in this case, methane). As expected, the AFT resulted in lower temperatures with a heat loss factor.
[See code ]
"""
Program to calculate Adiabatic flame temperature of methane combustion with heat loss
Stoichiometric reaction
CH4 + 2 (O2 + 3.76 N2) = CO2 + 2 H2O + 7.52 N2
Author : Logesh Paramashivam
Date : 26/05/2021
"""
import matplotlib.pyplot as plt
import math
import numpy as np
import scipy
def h(T, coeff): # Function to calculate enthalpy of a specific gas (inputs --> temperature and coefficients)
R = 8.314 # J/k.mol
a1 = coeff[0]
a2 = coeff[1]
a3 = coeff[2]
a4 = coeff[3]
a5 = coeff[4]
a6 = coeff[5]
# function returns the value of the entalpy for respective co efficients and temperature
return (a1 + a2*T/2 + a3*pow(T,2)/3 + a4*pow(T,3)/4 + a5*pow(T,4)/5 + a6/T) * R * T
# Co efficients a1...a7 for enthalpy calculation
# High temperature coefficients T>1000
co2_coeff_h = [3.85746029E+00, 4.41437026E-03, -2.21481404E-06, 5.23490188E-10, -4.72084164E-14, -4.87591660E+04, 2.27163806E+00]
h2o_coeff_h = [3.03399249E+00, 2.17691804E-03, -1.64072518E-07, -9.70419870E-11, 1.68200992E-14, -3.00042971E+04, 4.96677010E+00]
n2_coeff_h = [0.02926640E+02, 0.14879768E-02, -0.05684760E-05, 0.10097038E-09, -0.06753351E-13, -0.09227977E+04, 0.05980528E+02]
o2_coeff_h = [3.28253784E+00, 1.48308754E-03, -7.57966669E-07, 2.09470555E-10, -2.16717794E-14, -1.08845772E+03, 5.45323129E+00]
# Low temperature coefficients T<= 1000
o2_coeff_l = [3.78245636E+00, -2.99673416E-03, 9.84730201E-06, -9.68129509E-09, 3.24372837E-12, -1.06394356E+03, 3.65767573E+00]
n2_coeff_l = [0.03298677E+02, 0.14082404E-02, -0.03963222E-04, 0.05641515E-07,-0.02444854E-10, -0.10208999E+04, 0.03950372E+02]
ch4_coeff_l = [5.14987613E+00, -1.36709788E-02, 4.91800599E-05, -4.84743026E-08, 1.66693956E-11, -1.02466476E+04, -4.64130376E+00]
def f(T):
# Calculction of enthalpy of each species in the products
h_co2_p = h(T, co2_coeff_h) # calling function h, passing respective (high temp) co efficients and temperature
h_h2o_p = h(T, h2o_coeff_h)
h_n2_p = h(T, n2_coeff_h)
# Calculation of enthalpy of each species in the reactants
Tstd = 298.15 # kelvin # standard temperature of reactants before start of reaction
h_n2_r = h(Tstd, n2_coeff_l) # calling function h, passing respective (low temp) coefficients and initial temperature
h_o2_r = h(Tstd, o2_coeff_l)
h_ch4_r= h(Tstd, ch4_coeff_l)
# Total enthalpy of reactants
H_rectants = h_ch4_r + 2*h_o2_r + 7.52*h_n2_r
print(H_rectants)
# Total enthalpy of products
H_products = h_co2_p + 2*h_h2o_p + 7.52*h_n2_p
print(H_products)
# At Tadiabatic, H_products - H_reactants = 0 i.e, f(Tad) = 0
return H_products - H_rectants + (H_loss * LHV_CH4)
def fprime(T):
return (f(T+1e-6) - f(T))/1e-6 # Numerical derivative of f(x)
T_guess = 1800 # Initial guess for T adiabatic
aplha = 0.9
tol = 1e-6
ct = 1
H_loss = 0.50
LHV_CH4 = 50000 # J/g
"""
PERFORMING NEWTON RAPHSON LOOP
"""
while (abs(f(T_guess)) > tol):
T_guess = T_guess - aplha*(f(T_guess)/fprime(T_guess))
plt.plot(ct,T_guess,'*',color = 'red')
print (T_guess)
ct = ct + 1
print('Adiabatic flame temperature')
print(T_guess)
print('Iterations')
print(ct)
plt.show()
Problem 3:
In this problem, the adiabatic flame temperature of 2 C atoms of an alkane, alkene, and alkyne has been calculated. Cantera has been used for calculation, the algorithm similar to that one used in problem 1 has been used for this problem but with varying H atoms for each type of fuel based on their generic formula. It can be seen from the results, that the alkane with a single carbon bond requires less energy to break the bond compared to the double bond alkenes and the triple bond alkynes, hence the lower adiabatic flame temperature. It can be noticed that the AFT increases with increasing number of carbon bonds.
"""
Program to calculate Adiabatic flame temperature of Alkanes, alkenes and alkyne combustion
Author : Logesh Paramashivam
Date : 26/05/2021
"""
"""
General structure:
Alkanes: Cn H2n+2
Alkenes: Cn H2n
Alkynes: Cn H2n-2
Stoichiometric reaction
CxHy + (x+(y/4)) (O2 + 3.76 N2) = a.CO2 + b.H2O + c.N2
"""
import math
import matplotlib.pyplot as plt
import cantera as ct
gas = ct.Solution('gri30.xml')
T_alkane = []
T_alkene = []
T_alkyne = []
for i in range (0, 12):
lam_array = [0.9, 1, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0]
lam = lam_array[i]
x = 2
for j in range (0,3):
if j == 0 :
y = (2*x) + 2
n_o2 = lam * (x+ (y/4))
n_n2 = lam * (x+ (y/4)) * 3.76
gas.TPX = 298.15, 101325,{'C2H6':1, 'O2':n_o2, 'N2':n_n2}
gas.equilibrate('HP','auto')
T_alkane.append(gas.T)
elif j == 1 :
y = (2*x)
n_o2 = lam * (x+ (y/4))
n_n2 = lam * (x+ (y/4)) * 3.76
gas.TPX = 298.15, 101325,{'C2H4':1, 'O2':n_o2, 'N2':n_n2}
gas.equilibrate('HP','auto')
T_alkene.append(gas.T)
else :
y = 2
n_o2 = lam * (x+ (y/4))
n_n2 = lam * (x+ (y/4)) * 3.76
gas.TPX = 298.15, 101325,{'C2H2':1, 'O2':n_o2, 'N2':n_n2}
gas.equilibrate('HP','auto')
T_alkyne.append(gas.T)
plt.plot(lam_array, T_alkane,label = 'Alkane - C2H6')
plt.plot(lam_array, T_alkene,label = 'Alkene - C2H4')
plt.plot(lam_array, T_alkyne,label = 'Alkyne - C2H2')
plt.legend()
plt.xlabel('Lambda [-]')
plt.ylabel('Temperature [K]')
plt.title('Adiabatic flame temperature as a function of lambda')
plt.show()Problem 4:
In this problem, an approach similar to problem 3 has been used. Cantera has been used for adiabatic temperature calculation. Here, the AFT has been calculated for increasing the number of carbon atoms for alkane fuel. Based on the number of carbon atoms, the number of monoatomic hydrogen atoms also increases as dictated by the generic formula. It can be seen from the results that the higher number of c atoms release more energy for oxidation and hence the high AFT in increasing order of C atoms.
"""
Program to calculate Adiabatic flame temperature of Alkanes with multiple carbon atoms
Author : Logesh Paramashivam
Date : 26/05/2021
"""
"""
General structure:
Alkanes: Cn H2n+2
Stoichiometric reaction
CxHy + (x+(y/4)) (O2 + 3.76 N2) = a.CO2 + b.H2O + c.N2
"""
import math
import matplotlib.pyplot as plt
import cantera as ct
gas = ct.Solution('gri30.xml')
T_c1 = []
T_c2 = []
T_c3 = []
for i in range (0, 12):
lam_array = [ 0.9, 1, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0]
lam = lam_array[i]
for j in range (0,3):
if j == 0 :
x = 1
y = (2*x) + 2
n_o2 = lam * (x+ (y/4))
n_n2 = lam * (x+ (y/4)) * 3.76
gas.TPX = 298.15, 101325,{'CH4':1, 'O2':n_o2, 'N2':n_n2}
gas.equilibrate('HP','auto')
T_c1.append(gas.T)
elif j == 1 :
x = 2
y = (2*x) + 2
n_o2 = lam * (x+ (y/4))
n_n2 = lam * (x+ (y/4)) * 3.76
gas.TPX = 298.15, 101325,{'C2H6':1, 'O2':n_o2, 'N2':n_n2}
gas.equilibrate('HP','auto')
T_c2.append(gas.T)
else :
x = 3
y = (2*x) + 2
n_o2 = lam * (x+ (y/4))
n_n2 = lam * (x+ (y/4)) * 3.76
gas.TPX = 298.15, 101325,{'C3H8':1, 'O2':n_o2, 'N2':n_n2}
gas.equilibrate('HP','auto')
T_c3.append(gas.T)
plt.plot(lam_array, T_c1,label = '1 Carbon - CH4')
plt.plot(lam_array, T_c2,label = '2 Carbon - C2H6')
plt.plot(lam_array, T_c3,label = '3 Carbon - C3H8')
plt.legend()
plt.xlabel('Lambda [-]')
plt.ylabel('Temperature [K]')
plt.title('Adiabatic flame temperature for Alkanes')
plt.show()