Challenge

Adiabatic Flame Temperature calculation

AAbin Abraham23 Mar 2022

OBJECTIVE

To determine the Adiabatic Flame Temperature for combustion of hydrocarbon fuels and study the variation of equivalence ratio, number of carbon atoms, type of carbon bond and effect of heat loss using Python and Cantera

THEORY

In a chemical reaction the conversion of reactants to products involves heat transfer in the system. Considering the system to be adiabatic, that is, there is no heat transferred into the environment; the temperature attained inside the system is defined as adiabatic flame temperature (AFT). It is an ideal case temperature and the actual flame temperature observed for a particle fuel is always less than the AFT.

Newton-Raphson Method:

It is one of the approximation methods to find the roots of a real-valued function f (x)= 0. The method is based on finding the maximum and minimum values by approximating using a tangent to a continuous and differentiable function. using x_guess as the initial solution approximation, the new value is calculated as follows:

x_new = x_guess - (f(x_guess)/f'(x_guess))

In most of the cases, another ter called the relaxation factor (α) plays a major role, which makes the equation as

x_new = x_guess - α * (f(x_guess)/f'(x_guess))

Enthalpy equation for constant volume combustion

Using the enthalpy definition and first law of thermodynamics, we can write equations as:

h=U+pV

Applying differentials on both sides we get,

dh = q+ Vdp

For an adiabatic process, q=0. A lso, by substituiting pV =nRT, we get the equation as dh = Rd(nT)

i.e, h_products - h_reactants = R*[n_products * T_products - n_reactants*T_reactants]

Enthalpy equation for a constant prrssure system

For the constant pressure combustion, the enthalpy change is 0. This can be explained by the work done by the system to expand the volume of the system. Considering the system to have a heat loss factor we can write the equation as:

dh=0

h_products - h_reactants + h_loss = 0,

where, h_loss= heat loss* Lower Heating Value (LHV)

CODE AND RESULTS

1. Adiabatic flame temperature variation with Equivalence ratio (In constant volume-chamber) :

"""
Adiabatic flame temperature variation with Equivalence ratio (In constant volume-chamber) :

CH4 + 2(O2 + 3.76N2) = CO2 +2H2O + 7.52 N2 +0 O2 +0 CO   STOICHIOMETRIC MIXTURE

CH4 + 2/phi*(O2 + 3.76N2) = CO2 +2H2O + 7.52/phi* N2 +(2/phi -2)*O2 +0*CO   LEAN MIXTURE

CH4 + 2/phi*(O2 + 3.76N2) = (4/phi - 3)*CO2 +2H2O + 7.52/phi* N2 +0*O2+ (4-4/phi )*CO    RICH MIXTURE

"""


import matplotlib.pyplot as plt
import numpy as np

"""define a function to calculate the enthalpy"""

def h (T, coeffs):
    
    R= 8.314 #J/mol-K
    a1= coeffs[0]
    a2= coeffs[1]
    a3= coeffs[2]
    a4= coeffs[3]
    a5= coeffs[4]
    a6= coeffs[5]
    return (a1 + a2*T/2 + a3* pow(T,2)/3 + a4* pow(T,3)/4 + a5* pow(T,4)/5 +a6/T )*T*R

#CH4_coeffs_h =[7.48514950E-02,  1.33909467E-02, -5.73285809E-06,  1.22292535E-09, -1.01815230E-13, -9.46834459E+03,  1.84373180E+01]
#CO2_coeffs_l =[2.35677352E+00,  8.98459677E-03, -7.12356269E-06,  2.45919022E-09, -1.43699548E-13, -4.83719697E+04,  9.90105222E+00]
#H20_coeffs_l =[4.19864056E+00, -2.03643410E-03,  6.52040211E-06, -5.48797062E-09,  1.77197817E-12, -3.02937267E+04, -8.49032208E-01]
CH4_coeffs_l =[5.14987613E+00, -1.36709788E-02,  4.91800599E-05, -4.84743026E-08,  1.66693956E-11, -1.02466476E+04, -4.64130376E+00]
O2_coeffs_h  =[3.28253784E+00,  1.48308754E-03, -7.57966669E-07,  2.09470555E-10, -2.16717794E-14, -1.08845772E+03,  5.45323129E+00]
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_h  =[0.02926640E+02,  0.14879768E-02, -0.05684760E-05,  0.10097038E-09, -0.06753351E-13, -0.09227977E+04,  0.05980528E+02]
N2_coeffs_l  =[0.03298677E+02,  0.14082404E-02, -0.03963222E-04,  0.05641515E-07, -0.02444854E-10, -0.10208999E+04,  0.03950372E+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]
H20_coeffs_h =[3.03399249E+00,  2.17691804E-03, -1.64072518E-07, -9.70419870E-11,  1.68200992E-14, -3.00042971E+04,  4.96677010E+00]
CO_coeffs_h  =[2.71518561E+00,  2.06252743E-03, -9.98825771E-07,  2.30053008E-10,  -2.03647716E-14, -1.41518724E+04, 7.81868772E+00]


"""
function to find the roots
"""

def f(T, mol_r, mol_p):
    #num of moles
    N_mol_r= mol_r[0]+mol_r[1]+mol_r[2]
    N_mol_p= mol_p[0]+mol_p[1]+mol_p[2]+mol_p[3]+mol_p[4]

    #enthalpy of products
    h_CO2_p= h(T, CO2_coeffs_h)
    h_H2O_p= h(T, H20_coeffs_h)
    h_N2_p=  h(T, N2_coeffs_h)
    h_O2_p = h(T, O2_coeffs_h)
    h_CO_p = h(T, CO_coeffs_h)
    h_products= mol_p[0]*h_CO2_p+ mol_p[1]* h_H2O_p +mol_p[2]*h_N2_p+ mol_p[3]*h_O2_p +mol_p[4]*h_CO_p
    #print(h_products)
    
    #enthalpy of reactants
    T_std = 298.15
    h_CH4_r =h(T_std, CH4_coeffs_l)
    h_O2_r  =h(T_std, O2_coeffs_l)
    h_N2_r  =h(T_std, N2_coeffs_l)
    h_reactants = mol_r[0]*h_CH4_r +mol_r[1]*h_O2_r+ mol_r[2]*h_N2_r
    #print(h_reactants)
    h_r = h_products - h_reactants
    #print(h_r)
    print(T)
    return h_r
    return N_mol_r
    return N_mol_p



"""
Employing the Newton-Raphson method
"""

def fprime(T,mol_r, mol_p):
    return (f(T+1e-6, mol_r, mol_p)-f(T, mol_r, mol_p))/1e-6

"""
main code starts here 
"""

T_guess =1500
tol =1e-4
iter=0
alpha= 0.8
phi = np.linspace(0.1,2,50)  #equivalence ratio as array
AFT_python = [] #empty list

for i in range (0,50):
    if phi[i]==1:
        mol_r=[1,2,7.52]
        mol_p=[1,2,7.52,0,0]

    elif phi[i]<1:
        mol_r=[1,2/phi[i],7.52/phi[i]]
        mol_p=[1,2,7.52/phi[i],((2/phi[i])-2),0]
        
    elif phi[i]>1:
        mol_r=[1,2/phi[i],7.52/phi[i]]
        mol_p=[(4/phi[i])-3,2,7.52/phi[i],0,4-(4/phi[i])]
        
    

    while (abs(f(T_guess,mol_r, mol_p))> tol):
            T_guess=T_guess - alpha*(f(T_guess,mol_r, mol_p)/fprime(T_guess,mol_r, mol_p))
    AFT_python.append(T_guess)
    iter=iter+1
    


print(iter)
print(AFT_python)



"""
CANTERA
"""


"""
CH4 + 2(O2 + 3.76N2) = CO2 +2H2O + 7.52 N2 +0 O2 +0 CO   STOICHIOMETRIC MIXTURE

CH4 + 2/phi*(O2 + 3.76N2) = CO2 +2H2O + 7.52/phi* N2 +(2/phi -2)*O2 +0*CO   LEAN MIXTURE

CH4 + 2/phi*(O2 + 3.76N2) = (4/phi - 3)*CO2 +2H2O + 7.52/phi* N2 +0*O2+ (4-4/phi )*CO    RICH MIXTURE
"""

import cantera as ct
import matplotlib.pyplot as plt
import numpy as np
gas=ct.Solution('gri30.xml')

phi=np.linspace(0.1, 2, 50)
AFT_cantera=[]

for i in range(0,50):
    if phi[i]==1:
        mol_r=[1,2,7.52]
        mol_p=[1,2,7.52,0,0]

    elif phi[i]<1:
        mol_r=[1,2/phi[i],7.52/phi[i]]
        mol_p=[1,2,7.52/phi[i],((2/phi[i])-2),0]
        
    elif phi[i]>1:
        mol_r=[1,2/phi[i],7.52/phi[i]]
        mol_p=[(4/phi[i])-3,2,7.52/phi[i],0,4-(4/phi[i])]
        
    gas.TPX=298, 101325, {'CH4':mol_r[0] ,'O2':mol_r[1] , 'N2':mol_r[2] }
    gas.equilibrate('HP')
    AFT_cantera.append(gas.T)
print(gas.T)
plt.plot(phi,AFT_python,'-',color='red',label='python')
plt.plot(phi,AFT_cantera,'-',color='blue',label='cantera')
plt.xlabel('Equivalence ratio')
plt.ylabel('Adiabatic Flame Temperature in K')
plt.title('Effect of equivalence on AFT')
plt.legend()
plt.show()
#print(gas())

It can be observed from the plot that, at lower equivalence ratios, the results obtained from Python and Cantera are very similar. However,at higher equivalence ratios, the results tend to diverge. The values of AFT in the region close to the stoichiometric region are overestimated. This can be attributed to the heat loss in heating the inert nitrogen (N2), resulting in the formation of intermediary NOx compounds. Since the python is unable to capture this, the AFT values observed and more than that calculated using Cantera

2. Flame temperature variation with Heat loss (In constant pressure-chamber) :

"""
Adiabatic flame temperature variation with heat loss (In constant volume-chamber) :

CxHy + (x+y/4)*(O2 + 3.76N2) = xCO2 +y/2*H2O + (x+y/4)*3.76 N2    STOICHIOMETRIC MIXTURE


"""


import matplotlib.pyplot as plt
import numpy as np

"""define a function to calculate the enthalpy"""

def h (T, coeffs):
    
    R= 8.314 #J/mol-K
    a1= coeffs[0]
    a2= coeffs[1]
    a3= coeffs[2]
    a4= coeffs[3]
    a5= coeffs[4]
    a6= coeffs[5]
    return (a1 + a2*T/2 + a3* pow(T,2)/3 + a4* pow(T,3)/4 + a5* pow(T,4)/5 +a6/T )*T*R

#CH4_coeffs_h =[7.48514950E-02,  1.33909467E-02, -5.73285809E-06,  1.22292535E-09, -1.01815230E-13, -9.46834459E+03,  1.84373180E+01]
#CO2_coeffs_l =[2.35677352E+00,  8.98459677E-03, -7.12356269E-06,  2.45919022E-09, -1.43699548E-13, -4.83719697E+04,  9.90105222E+00]
#H20_coeffs_l =[4.19864056E+00, -2.03643410E-03,  6.52040211E-06, -5.48797062E-09,  1.77197817E-12, -3.02937267E+04, -8.49032208E-01]
CH4_coeffs_l =[5.14987613E+00, -1.36709788E-02,  4.91800599E-05, -4.84743026E-08,  1.66693956E-11, -1.02466476E+04, -4.64130376E+00]
O2_coeffs_h  =[3.28253784E+00,  1.48308754E-03, -7.57966669E-07,  2.09470555E-10, -2.16717794E-14, -1.08845772E+03,  5.45323129E+00]
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_h  =[0.02926640E+02,  0.14879768E-02, -0.05684760E-05,  0.10097038E-09, -0.06753351E-13, -0.09227977E+04,  0.05980528E+02]
N2_coeffs_l  =[0.03298677E+02,  0.14082404E-02, -0.03963222E-04,  0.05641515E-07, -0.02444854E-10, -0.10208999E+04,  0.03950372E+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]
H20_coeffs_h =[3.03399249E+00,  2.17691804E-03, -1.64072518E-07, -9.70419870E-11,  1.68200992E-14, -3.00042971E+04,  4.96677010E+00]
CO_coeffs_h  =[2.71518561E+00,  2.06252743E-03, -9.98825771E-07,  2.30053008E-10,  -2.03647716E-14, -1.41518724E+04, 7.81868772E+00]


"""
function to find the roots
"""

def f(T, mol_r, mol_p,h_loss):
    

    #enthalpy of products
    h_CO2_p= h(T, CO2_coeffs_h)
    h_H2O_p= h(T, H20_coeffs_h)
    h_N2_p=  h(T, N2_coeffs_h)
    h_O2_p = h(T, O2_coeffs_h)
    h_CO_p = h(T, CO_coeffs_h)
    h_products= mol_p[0]*h_CO2_p+ mol_p[1]* h_H2O_p +mol_p[2]*h_N2_p
    #print(h_products)
    
    #enthalpy of reactants
    T_std = 298.15
    h_CH4_r =h(T_std, CH4_coeffs_l)
    h_O2_r  =h(T_std, O2_coeffs_l)
    h_N2_r  =h(T_std, N2_coeffs_l)
    h_reactants = mol_r[0]*h_CH4_r +mol_r[1]*h_O2_r+ mol_r[2]*h_N2_r
    #print(h_reactants)
    h_r = h_products - h_reactants + h_loss*LHV
    #print(h_r)
    #print(T)
    return h_r



"""
Employing the Newton-Raphson method
"""

def fprime(T,mol_r, mol_p,h_loss):
    return (f(T+1e-6, mol_r, mol_p,h_loss)-f(T, mol_r, mol_p,h_loss))/1e-6

"""
main code starts here 
"""

T_guess =1500
tol =1e-4
iter=0
alpha= 0.8
h_loss=np.linspace(0.1,1,50)
x=1
y=2*x+2
T_std = 298.15
mol_r=[1,x+y/4,3.76*(x+y/4)]
mol_p=[x,y/2,3.76*(x+y/4)]
LHV=0
LHV =-f(T_std, mol_r, mol_p,0)
AFT = [] #empty list

    
for i in range (0,50):
    while (abs(f(T_guess,mol_r, mol_p,h_loss[i]))> tol):
            T_guess=T_guess - alpha*(f(T_guess,mol_r, mol_p,h_loss[i])/fprime(T_guess,mol_r, mol_p,h_loss[i]))
    AFT.append(T_guess)
    iter=iter+1
    
print(iter)
print(AFT)

plt.plot(h_loss,AFT,'-',color='red')
plt.xlabel('heat loss')
plt.ylabel('Adiabatic Flame Temperature in K')
plt.title('AFT for varying heat loss')
plt.show()

It can be observed that, as the heat loss increases, the corresponding AFT values decreases, and once the heat loss value becomes 1, the AFT becomes the STP temperature, 298.15K.

When the heat loss is 0.35, then the AFT value is 1680K, which is considerably lesser than the AFT value in the previous case (2325K)

3. Adiabatic flame temperature variation with respect to Alkane, Alkene, Alkyne :

"""
Adiabatic flame temperature variation with alkane, alkene, alkyne
C2H6 + 3.5*(O2 + 3.76N2) = 2*CO2+ 3*H2O+ (3.5*3.76)*N2  ALKANE
C2H4 + 3*  (O2 + 3.76N2) = 2*CO2+ 2*H2O+ (3*3.76)*N2  ALKENE
C2H2 + 2.5*(O2 + 3.76N2) = 2*CO2+   H2O+ (2.5*3.76)*N2  ALKYNE


"""


import matplotlib.pyplot as plt
import numpy as np

"""define a function to calculate the enthalpy"""

def h (T, coeffs):
    
    R= 8.314 #J/mol-K
    a1= coeffs[0]
    a2= coeffs[1]
    a3= coeffs[2]
    a4= coeffs[3]
    a5= coeffs[4]
    a6= coeffs[5]
    return (a1 + a2*T/2 + a3* pow(T,2)/3 + a4* pow(T,3)/4 + a5* pow(T,4)/5 +a6/T )*T*R

C2H6_coeffs_l=[4.29142492E+00, -5.50154270E-03,  5.99438288E-05, -7.08466285E-08,  2.68685771E-11, -1.15222055E+04,  2.66682316E+00]
C2H4_coeffs_l=[3.95920148E+00, -7.57052247E-03,  5.70990292E-05, -6.91588753E-08,  2.69884373E-11,  5.08977593E+03,  4.09733096E+00]
#C2H4_coeffs_h=[2.03611116E+00, 1.46454151E-02,  -6.71077915E-06,  1.47222923E-09,  -1.25706061E-13, 4.93988614E+03, 1.03053693E+01]
#C2H2_coeffs_h=[4.14756964E+00, 5.96166664E-03,  -2.37294852E-06,  4.67412171E-10,  -3.61235213E-14, 2.59359992E+04, -1.23028121E+00]
C2H2_coeffs_l=[8.08681094E-01,  2.33615629E-02, -3.55171815E-05,  2.80152437E-08, -8.50072974E-12,  2.64289807E+04,  1.39397051E+01]
O2_coeffs_h  =[3.28253784E+00,  1.48308754E-03, -7.57966669E-07,  2.09470555E-10, -2.16717794E-14, -1.08845772E+03,  5.45323129E+00]
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_h  =[0.02926640E+02,  0.14879768E-02, -0.05684760E-05,  0.10097038E-09, -0.06753351E-13, -0.09227977E+04,  0.05980528E+02]
N2_coeffs_l  =[0.03298677E+02,  0.14082404E-02, -0.03963222E-04,  0.05641515E-07, -0.02444854E-10, -0.10208999E+04,  0.03950372E+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]
H20_coeffs_h =[3.03399249E+00,  2.17691804E-03, -1.64072518E-07, -9.70419870E-11,  1.68200992E-14, -3.00042971E+04,  4.96677010E+00]

"""
function to find the roots
"""

def f(T,x):

    #enthalpy of products
    h_CO2_p= h(T, CO2_coeffs_h)
    h_H2O_p= h(T, H20_coeffs_h)
    h_N2_p=  h(T, N2_coeffs_h)
    h_O2_p = h(T, O2_coeffs_h)
    #print(h_products)
    
    #enthalpy of reactants
    T_std = 298.15
    h_C2H6_r =h(T_std, C2H6_coeffs_l)
    h_C2H4_r =h(T_std, C2H4_coeffs_l)
    h_C2H2_r =h(T_std, C2H2_coeffs_l)
    h_O2_r  =h(T_std, O2_coeffs_l)
    h_N2_r  =h(T_std, N2_coeffs_l)

    if x[i]==1:
        h_fuel=h_C2H6_r
        mol_r=[1, 3.5, 3.5*3.76]
        mol_p=[2, 3, 3.5*3.76]

    elif x[i]==2:
        h_fuel=h_C2H4_r
        mol_r=[1, 3, 3*3.76]
        mol_p=[2, 2, 3*3.76]
            
    elif x[i]==3:
        h_fuel=h_C2H2_r
        mol_r=[1, 2.5, 2.5*3.76]
        mol_p=[2, 1, 2.5*3.76]
        
        
    h_reactants = mol_r[0]*h_fuel +mol_r[1]*h_O2_r+ mol_r[2]*h_N2_r
    h_products= mol_p[0]*h_CO2_p+ mol_p[1]* h_H2O_p +mol_p[2]*h_N2_p
    #print(h_reactants)
    h_r = h_products - h_reactants
    #print(h_r)
    print(T)
    return h_r
    



"""
Employing the Newton-Raphson method
"""

def fprime(T,x):
    return (f(T+1e-6,x)-f(T,x))/1e-6

"""
main code starts here 
"""

T_guess =1500
tol =1e-4
iter=0
alpha= 0.8
x=np.linspace(1,3,3)
AFT = [] #empty list

for i in range (0,3):   

    while (abs(f(T_guess,x))> tol):
        T_guess=T_guess - alpha*(f(T_guess,x)/fprime(T_guess,x))
    AFT.append(T_guess)
    iter=iter+1
    


print(iter)
print(AFT)

plt.bar('C2H6',AFT[0],color='red',width=0.2)
plt.bar('C2H4',AFT[1],color='blue',width=0.2)
plt.bar('C2H2',AFT[2],color='green',width=0.2)
plt.xlabel('Hydrocarbons')
plt.ylabel('Adiabatic Flame Temperature in K')
plt.title('AFT vs saturated and unsaturated hydrocarbons')
plt.show()

AFT _C2H6 =2379.8494589947854

AFT _C2H4 =2564.5336891398315

AFT _C2H2 =2909.354567508851

From the plots it can be observed that the flame temperature increases with increase in the number of carbon bonds. An explanation for this can be the reduction of hydrogen atoms. Due to this, less number of water molecules is formed for alkyne as compared to alkene and alkanes. The result of this is the high flame temperature at thermal equilibrium for combustion of alkyne hydrocarbon compared to alkane and alkenes.

4. Adiabatic flame temperature variation with respect to number of Carbon atom chains :

"""
Adiabatic flame temperature variation with respect to number of Carbon atom chains :
CH4  + 2*(O2+3.76N2)     = CO2+   2*H2O+(2*3.76)* N2 METHANE
C2H6 + 3.5*(O2 + 3.76N2) = 2*CO2+ 3*H2O+ (3.5*3.76)*N2  ETHANE
C3H8 + 5*  (O2 + 3.76N2) = 3*CO2+ 4*H2O+ (5*3.76)*N2  PROPANE

"""


import matplotlib.pyplot as plt
import numpy as np

"""define a function to calculate the enthalpy"""

def h (T, coeffs):
    
    R= 8.314 #J/mol-K
    
    a1= coeffs[0]
    a2= coeffs[1]
    a3= coeffs[2]
    a4= coeffs[3]
    a5= coeffs[4]
    a6= coeffs[5]
    return (a1 + a2*T/2 + a3* pow(T,2)/3 + a4* pow(T,3)/4 + a5* pow(T,4)/5 +a6/T )*T*R

CH4_coeffs_l =[5.14987613E+00, -1.36709788E-02,  4.91800599E-05, -4.84743026E-08,  1.66693956E-11, -1.02466476E+04, -4.64130376E+00]
#C2H6_coeffs_h=[1.07188150E+00,  2.16852677E-02, -1.00256067E-05,  2.21412001E-09, -1.90002890E-13, -1.14263932E+04,  1.51156107E+01]
C2H6_coeffs_l=[4.29142492E+00, -5.50154270E-03,  5.99438288E-05, -7.08466285E-08,  2.68685771E-11, -1.15222055E+04,  2.66682316E+00]
#C3H8_coeffs_h=[0.75341368E+01,  0.18872239E-01, -0.62718491E-05,  0.91475649E-09, -0.47838069E-13, -0.16467516E+05  -0.17892349E+02]   
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_h  =[3.28253784E+00,  1.48308754E-03, -7.57966669E-07,  2.09470555E-10, -2.16717794E-14, -1.08845772E+03,  5.45323129E+00]
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_h  =[0.02926640E+02,  0.14879768E-02, -0.05684760E-05,  0.10097038E-09, -0.06753351E-13, -0.09227977E+04,  0.05980528E+02]
N2_coeffs_l  =[0.03298677E+02,  0.14082404E-02, -0.03963222E-04,  0.05641515E-07, -0.02444854E-10, -0.10208999E+04,  0.03950372E+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]
H20_coeffs_h =[3.03399249E+00,  2.17691804E-03, -1.64072518E-07, -9.70419870E-11,  1.68200992E-14, -3.00042971E+04,  4.96677010E+00]

"""
function to find the roots
"""

def f(T,x):

    #enthalpy of products
    h_CO2_p= h(T, CO2_coeffs_h)
    h_H2O_p= h(T, H20_coeffs_h)
    h_N2_p=  h(T, N2_coeffs_h)
    h_O2_p = h(T, O2_coeffs_h)
    #print(h_products)
    
    #enthalpy of reactants
    T_std = 298.15
    h_CH4_r  =h(T_std, CH4_coeffs_l)
    h_C2H6_r =h(T_std, C2H6_coeffs_l)
    h_C3H8_r =h(T_std, C3H8_coeffs_l)
    h_O2_r  =h(T_std, O2_coeffs_l)
    h_N2_r  =h(T_std, N2_coeffs_l)

    if x[i]==1:
        h_fuel=h_CH4_r
        mol_r=[1, 2, 2*3.76]
        mol_p=[1, 2, 2*3.76]

    elif x[i]==2:
        h_fuel=h_C2H6_r
        mol_r=[1, 3.5, 3.5*3.76]
        mol_p=[2, 3, 3.5*3.76]
            
    elif x[i]==3:
        h_fuel=h_C3H8_r
        mol_r=[1, 5, 5*3.76]
        mol_p=[3, 4, 5*3.76]
        
        
    h_reactants = mol_r[0]*h_fuel +mol_r[1]*h_O2_r+ mol_r[2]*h_N2_r
    h_products= mol_p[0]*h_CO2_p+ mol_p[1]* h_H2O_p +mol_p[2]*h_N2_p
    #print(h_reactants)
    h_r = h_products - h_reactants
    #print(h_r)
    print(T)
    return h_r
    



"""
Employing the Newton-Raphson method
"""

def fprime(T,x):
    return (f(T+1e-6,x)-f(T,x))/1e-6

"""
main code starts here 
"""

T_guess =1500
tol =1e-4
iter=0
alpha= 0.8
x=np.linspace(1,3,3)
AFT = [] #empty list

for i in range (0,3):   

    while (abs(f(T_guess,x))> tol):
        T_guess=T_guess - alpha*(f(T_guess,x)/fprime(T_guess,x))
    AFT.append(T_guess)
    iter=iter+1
    


print(iter)
print(AFT)

plt.bar('CH4',AFT[0],color='red',width=0.2)
plt.bar('C2H6',AFT[1],color='blue',width=0.2)
plt.bar('C3H8',AFT[2],color='green',width=0.2)
plt.xlabel('Hydrocarbons')
plt.ylabel('Adiabatic Flame Temperature in K')
plt.title('AFT vs alkanes')
plt.show()

From the plot, it can be concluded that the AFT increases with an increase in carbon atoms. With an increasing number of carbon atoms, mass of fuel is more which combines with more number of air molecules. i.e, the combustion mixture is enriched for same volume of air, when more number of carbon atoms are present. The energy released due to this increases the product enthalpies and the net enthalpy difference increases. In turn, AFT observed is more than of lesser number of carbon atoms in case of hydrocarbon fuel combustion.

AFT_CH4 = 2325.5981296417945

AFT_C2H6 = 2379.849459070052

AFT_C3H8 = 2392.0975960044666

CONCLUSION

Heat loss from a combustion system results in drop of AFT.
The flame temperature is highest for alkyne, when compared with alkanes and alkenes
The number of carbon atoms in hydrocarbon has a direct relation with AFT.

Adiabatic Flame Temperature calculation

Ship work like this
in your own course.

Talk to a career counsellor

Ship work like thisin your own course.

Talk to a career counsellor

Ship work like this
in your own course.