Week 3 - Adiabatic Flame Temperature calculation
Adiabatic flame temperature variation with Equivalence ratio (In constant volume-chamber) :
AIM: To calculate the Adiabatic Flame Temperature.
OBJECTIVES:
1) Plot the effect of equivalence ratio on the final adiabatic flame temperature and further compare the results with Cantera.
2) Find out the Adiabatic Flame Temperature of Methane Air Combustion in constant pressure with heat loss.
3) To examine the Adiabatic Flame Temperature of Alk-ane, Alk-ene, Alk-yne with no heat loss.
4) Find out the Adiabatic flame temperature variation with respect to a number of carbon atom chains.
import cantera as cnt
import matplotlib.pyplot as plt
import math
import numpy as np
# CH4 + 2(O2 + 3.76N2) = CO2 + 2H2O + 7.52N2
def h(T,co_effs):
a1 = co_effs[0]
a2 = co_effs[1]
a3 = co_effs[2]
a4 = co_effs[3]
a5 = co_effs[4]
a6 = co_effs[5]
return(a1 + a2*T/2 + a3*pow(T,2)/3 + a4*pow(T,3)/4 + a5*pow(T,4)/5 + a6/T)*R*T
R = 8.314 #j/mol-k
CH4_coeff_lw = [5.14987613E+00,-1.36709788E-02,4.91800599E-05,-4.84743026E-08,1.66693956E-11,-1.02466476E+04,-4.64130376E+00]
O2_coeff_lw = [3.78245636E+00,-2.99673416E-03,9.84730201E-06,-9.68129509E-09,3.24372837E-12,-1.06394356E+03,3.65767573E+00]
N2_coeff_lw = [0.03298677E+02,0.14082404E-02,-0.03963222E-04,0.05641515E-07,-0.02444854E-10,-0.10208999E+04,0.03950372E+02]
CO2_coeff_hg = [3.85746029E+00,4.41437026E-03,-2.21481404E-06,5.23490188E-10,-4.72084164E-14,-4.87591660E+04,2.27163806E+00]
N2_coeff_hg = [0.02926640E+02,0.14879768E-02,-0.05684760E-05,0.10097038E-09,-0.06753351E-13,-0.09227977E+04,0.05980528E+02]
H2O_coeff_hg = [3.03399249E+00,2.17691804E-03,-1.64072518E-07,-9.70419870E-11,1.68200992E-14,-3.00042971E+04,4.96677010E+00]
O2_coeff_hg = [3.28253784E+00,1.48308754E-03,-7.57966669E-07,2.09470555E-10,-2.16717794E-14,-1.08845772E+03,5.45323129E+00]
CO_coeff_hg = [2.71518561E+00,2.06252743E-03,-9.98825771E-07,2.30053008E-10,-2.03647716E-14,-1.41518724E+04,7.81868772E+00]
def f(T,psi):
H_CO2_p = h(T,CO2_coeff_hg)
H_H2O_p = h(T,H2O_coeff_hg)
H_N2_p = h(T, N2_coeff_hg)
H_O2_p = h(T, O2_coeff_hg)
H_CO_p = h(T,CO_coeff_hg)
T_std = 298.15 # for reactants
H_CH4_r = h(T_std,CH4_coeff_lw)
H_O2_r = h(T_std,O2_coeff_lw)
H_N2_r = h(T_std,N2_coeff_lw)
# for lean mixture psi < 1
if psi <= 1:
H_reactants = H_CH4_r+(2/psi)*(H_O2_r+3.76*H_N2_r)
H_products = H_CO2_p+2*H_H2O_p+(7.52/psi)*H_N2_p+((2/psi)-2)*H_O2_p
N_reactants = 1+(2/psi)*(1+3.76)
N_productants = 1+2+(7.52/psi)+((2/psi)-2)
else: # for rich mixture
H_reactants = H_CH4_r+(2/psi)*(H_O2_r+(3.76*H_N2_r))
H_products = ((4/psi)-3)*H_CO2_p+2*H_H2O_p+(7.52/psi)*H_N2_p+(4-(4/psi))*H_CO_p
N_reactants = 1+(2/psi)*(1+3.76)
N_productants = ((4/psi)-3)+2+(7.52/psi)+(4-(4/psi))
return(-H_products+H_reactants-R*(-N_productants*T+N_reactants*T_std))
def fprime(T,psi):
return(f(T+1e-6,psi)-f(T,psi))/1e-6
# NRM
alfa = 1
ct = 0
T_guess = 1500
x = np.linspace(0.2,2,50)
T_ad = []
for psi in x:
tol = 1e-6
T_guess = 1500
while (abs(f(T_guess,psi))>tol):
T_guess = T_guess-alfa*((f(T_guess,psi))/(fprime(T_guess,psi)))
#print(T_guess)
ct = ct+1
T_ad.append(T_guess)
plt.plot(x,T_ad,color = 'red',label='Python')
print(ct)
#plt.show()
# cantera
gas = cnt.Solution('gri30.yaml')
gas.TPX = 298.15,101325,{'CH4':1,'O2':2,'N2':7.75}
#print(gas.X)
gas.equilibrate('UV','auto')
print(gas.T)
T_cantera=[]
for phi in x:
gas = cnt.Solution('gri30.yaml')
gas.TPX = 298.15,101325,{'CH4':1,'O2':2/phi,'N2':7.75/phi}
gas.equilibrate('UV','auto')
T_cantera.append(gas.T)
plt.plot(x,T_cantera,color='blue', label='cant')
plt.legend()
plt.show()
AFT with Heat loss
import cantera as cnt
import matplotlib.pyplot as plt
import math
import numpy as np
# CH4 + 2(O2 + 3.76N2) = CO2 + 2H2O + 7.52N2
def h(T,co_effs):
a1 = co_effs[0]
a2 = co_effs[1]
a3 = co_effs[2]
a4 = co_effs[3]
a5 = co_effs[4]
a6 = co_effs[5]
return(a1 + a2*T/2 + a3*pow(T,2)/3 + a4*pow(T,3)/4 + a5*pow(T,4)/5 + a6/T)*R*T
R = 8.314 #j/mol-k
lower_heating_value = 50E6
Mfuel = 16.04E-3
heat_supplied = lower_heating_value*Mfuel
# low tempareures coefficients
CH4_coeff_lw = [5.14987613E+00,-1.36709788E-02,4.91800599E-05,-4.84743026E-08,1.66693956E-11,-1.02466476E+04,-4.64130376E+00]
O2_coeff_lw = [3.78245636E+00,-2.99673416E-03,9.84730201E-06,-9.68129509E-09,3.24372837E-12,-1.06394356E+03,3.65767573E+00]
N2_coeff_lw = [0.03298677E+02,0.14082404E-02,-0.03963222E-04,0.05641515E-07,-0.02444854E-10,-0.10208999E+04,0.03950372E+02]
# high temparetures coefficients
CO2_coeff_hg = [3.85746029E+00,4.41437026E-03,-2.21481404E-06,5.23490188E-10,-4.72084164E-14,-4.87591660E+04,2.27163806E+00]
N2_coeff_hg = [0.02926640E+02,0.14879768E-02,-0.05684760E-05,0.10097038E-09,-0.06753351E-13,-0.09227977E+04,0.05980528E+02]
H2O_coeff_hg = [3.03399249E+00,2.17691804E-03,-1.64072518E-07,-9.70419870E-11,1.68200992E-14,-3.00042971E+04,4.96677010E+00]
O2_coeff_hg = [3.28253784E+00,1.48308754E-03,-7.57966669E-07,2.09470555E-10,-2.16717794E-14,-1.08845772E+03,5.45323129E+00]
CO_coeff_hg = [2.71518561E+00,2.06252743E-03,-9.98825771E-07,2.30053008E-10,-2.03647716E-14,-1.41518724E+04,7.81868772E+00]
def f(T,heat_loss):
H_CO2_p = h(T,CO2_coeff_hg)
H_H2O_p = h(T,H2O_coeff_hg)
H_N2_p = h(T, N2_coeff_hg)
H_O2_p = h(T, O2_coeff_hg)
H_CO_p = h(T,CO_coeff_hg)
T_std = 298.15 # for reactants
H_CH4_r = h(T_std,CH4_coeff_lw)
H_O2_r = h(T_std,O2_coeff_lw)
H_N2_r = h(T_std,N2_coeff_lw)
H_reactants = H_CH4_r+2*(H_O2_r+3.76*H_N2_r)
H_products = H_CO2_p+2*H_H2O_p+7.52*H_N2_p
return(H_products-H_reactants)+(heat_loss*heat_supplied)
def fprime(T,heat_loss):
return(f(T+1e-6,heat_loss)-f(T,heat_loss))/1e-6
# NRM
alfa = 1
ct = 0
T_guess = 1500
x = np.linspace(0,1,21)
T_ad = []
for heat_loss in x:
tol = 1e-6
T_guess = 1500
while (abs(f(T_guess,heat_loss))>tol):
T_guess = T_guess-alfa*((f(T_guess,heat_loss))/(fprime(T_guess,heat_loss)))
ct = ct+1
print(T_guess)
T_ad.append(T_guess)
plt.plot(x,T_ad,color = 'red',label='Python')
plt.xlabel('Heat loss')
plt.ylabel('Adiabatic flame temp')
plt.title('variation AFT with Heat loss')
plt.show()
- AFT varies linearly with heat loss coefficient
Adiabatic flame temperature variation with respect to Alkane, Alkene, Alkyne :
import matplotlib.pyplot as plt
import numpy as np
import cantera as cnt
#For alkane
#C2H6 + 3.5(O2 + 3.76N2) <-> 2CO2 + 3H20 + 13.16N2
gas = cnt.Solution('gri30.xml')
gas.TPX = 298.15, 101325, {'C2H6':1, 'O2':3.5, 'N2':13.16}
gas.equilibrate('HP','auto')
T_alkane = gas.T
print(T_alkane)
#For alkene
#C2H4 + 3(O2 + 3.76N2) <-> 2CO2 + 2H20 + 11.28N2
gas = cnt.Solution('gri30.xml')
gas.TPX = 298.15, 101325, {'C2H4':1, 'O2':3, 'N2':11.28}
gas.equilibrate('HP','auto')
T_alkene = gas.T
print(T_alkene)
#For alkyne
#C2H6 + 2.5(O2 + 3.76N2) <-> 2CO2 + H20 + 9.4N2
gas = cnt.Solution('gri30.xml')
gas.TPX = 298.15, 101325, {'C2H2':1, 'O2':2.5, 'N2':9.4}
gas.equilibrate('HP','auto')
T_alkyne = gas.T
print(T_alkyne)
Temp = [T_alkane, T_alkene, T_alkyne]
plt.figure()
plt.bar(['Ethane', 'Ethene', 'Ethyne'],Temp,color='red')
plt.ylabel('Adiabetic Flame Temperatur (AFT) in K')
plt.title('AFT for various hydrocarbon fuels')
plt.show()
Adiabatic flame temperature variation with respect to a number of Carbon atom chains :
"""
Adiabatic flame temperature variation with respect to Alkane, Alkene and Alkyne
"""
# general hydrocarbon: CxHy + [x + (y/4)](O2 + 3.76N2) => xCO2 + (y/2)H2O + 3.76[x+(y/4)]N2
import matplotlib.pyplot as plt
import cantera as ct
def h(T,co_effs):
R = 8.314 #j/mol-k
a1 = co_effs[0]
a2 = co_effs[1]
a3 = co_effs[2]
a4 = co_effs[3]
a5 = co_effs[4]
a6 = co_effs[5]
return (a1 + a2*T/2 + a3*pow(T,2)/3 + a4*pow(T,3)/4 + a5*pow(T,4)/5 + a6/T)*R*T
#low temperature coefficients
c2h6_coeffs_l = [4.29142492E+00,-5.50154270E-03,5.99438288E-05,-7.08466285E-08,2.68685771E-11,-1.15222055E+04,2.66682316E+00]
ch4_coeffs_l = [5.14987613E+00,-1.36709788E-02,4.91800599E-05,-4.84743026E-08,1.66693956E-11,-1.02466476E+04,-4.64130376E+00]
c3h8_coeffs_l = [0.93355381E+00,0.26424579E-01,0.61059727E-05,-0.21977499E-07,0.95149253E-11,-0.13958520E+05,0.19201691E+02]
o2_coeffs_l = [3.78245636E+00,-2.99673416E-03,9.84730201E-06,-9.68129509E-09,3.24372837E-12,-1.06394356E+03,3.65767573E+00]
n2_coeffs_l = [0.03298677E+02,0.14082404E-02,-0.03963222E-04,0.05641515E-07,-0.02444854E-10,-0.10208999E+04,0.03950372E+02]
#high temperature coefficients
n2_coeffs_h = [0.02926640E+02,0.14879768E-02,-0.05684760E-05,0.10097038E-09,-0.06753351E-13,-0.09227977E+04,0.05980528E+02]
co2_coeffs_h = [3.85746029E+00,4.41437026E-03,-2.21481404E-06,5.23490188E-10,-4.72084164E-14,-4.87591660E+04,2.27163806E+00]
h2o_coeffs_h = [3.03399249E+00,2.17691804E-03,-1.64072518E-07,-9.70419870E-11,1.68200992E-14,-3.00042971E+04,4.96677010E+00]
#standard temperature
Tstd = 298.15
ntype = [1,2,3]
# Part A - Newton Raphson Method
def f(T,ntype):
# Products
h_co2_p = h(T,co2_coeffs_h)
h_n2_p = h(T,n2_coeffs_h)
h_h2o_p = h(T,h2o_coeffs_h)
# Reactants
h_ch4_r = h(Tstd,ch4_coeffs_l)
h_c2h6_r = h(Tstd,c2h6_coeffs_l)
h_c3h8_r = h(Tstd,c3h8_coeffs_l)
h_o2_r = h(Tstd,o2_coeffs_l)
h_n2_r = h(Tstd,n2_coeffs_l)
if ntype == 1:
x = 1
y = 4
h_fuel = h_ch4_r
elif ntype == 2:
x = 2
y = 6
h_fuel = h_c2h6_r
elif ntype == 3:
x = 3
y = 8
h_fuel = h_c3h8_r
H_products = x*h_co2_p + (y/2)*h_h2o_p + 3.76*(x + y/4)*h_n2_p
H_reactants = h_fuel + (x + y/4)*h_o2_r + 3.76*(x + y/4)*h_n2_r
return (H_products - H_reactants)
def fprime(T,ntype):
return (f(T+1e-6,ntype) - f(T,ntype))/(1e-6)
T_guess = 1500
tol = 1e-3
alpha = 0.2
ntype = [1,2,3]
temp = []
for i in range(0,len(ntype)):
while(abs(f(T_guess,ntype[i])) > tol):
T_guess = T_guess - alpha*((f(T_guess,ntype[i])/fprime(T_guess,ntype[i])))
temp.append(T_guess)
print(temp)
# Plots
plt.figure
plt.bar(['Methane','Ethane','Propane'],temp)
#plt.ylabel('Adiabatic Flame Temperature')
#plt.title('Effect of AFT with respect to no carbon atom chain')
plt.show()