Week 3 - Adiabatic Flame Temperature calculation

The Evaluation of Adiabatic Flame Temperature using python and Cantera

1. Write a program to analyze the effect of equivalence ratio on the final adiabatic flame temperature of methane contained in a constant volume chamber and also compare the results obtained from cantera.
2. Write a program to analyze the effect of frictional heat loss on the final adiabatic flame temperature of methane contained in a constant pressure chamber.
3. Write a program to calculate the AFT of the combustion process for a given equivalence ratio and heat loss fraction.
4. Write a program to analyze the variation of Adiabatic Flame Temperature for different hydrocarbons (Ethane, Ethene and Ethyne) contained in a constant pressure chamber.
5. Write a program to analyze the variation of Adiabatic Flame Temperature for different alkanes (Methane,Ethane, and Propane) contained in a constant pressure chamber.

• Adiabatic Flame Temperature

When a combustion reaction takes place energy is released to the combustion products. If no heat is lost in this process, the temperature of the combustion products is known as the "Adiabatic Flame Temperature." For methane combustion in air at 1 atmosphere the Adiabatic Flame Temperature is 2328K or 2055C. For hydrogen burning in air at 1 atmosphere the Adiabatic Flame temperature is 2400K or 2127C.

In other words, the adiabatic flame temperature is also defined as the temperature attained when all of the chemical reaction heat released heats combustion products. The adiabatic flame temperature is one that occurs when the combustion chamber is well insulated with no heat losses (adiabatic conditions). The peak adiabatic flame temperature occurs at around Φ = 1 in an ideally insulated combustion chamber.

• Processes involved in combustion

In the study of combustion, there are two types of adiabatic flame temperature depending on how the process is completed:

• constant volume
• constant pressure

both of which describe the temperature that combustion products theoretically can reach if no energy is lost to the outside environment.

Constant Volume

The constant volume adiabatic flame temperature is the temperature that results from a complete combustion process that occurs without any work, heat transfer or changes in kinetic or potenrial energy. Its temperature is higher than the constant pressure process because no energy is utilized to change the volume of the system. This is the maximum temperature that can be achieved for given reactants because any heat transfer from the reacting substances and/or any incomplete combustion would tend to lower the temperature of the products.

In constant volume combustion, the equation can also be written as

`U_(Reactant)(T_(Initial), P_(Initial)) = U_(Product)(T_(ad), P_f)`

We can rearrange the above eqution in the following form

`H(Reactant) - H(Product) - V(P_(Init) - P_(f))` = 0

According to ideal gas relation PV = nRT

`P_(Initial) . V = N_(reactants) . R_(U) .T_(Initial)`

`P_(f) . V = N_(reactants) . R_(U) .T_(ad)`

The equation(ii) can be written as

`H_(Reactant) - H_(Product) - R_(U) .(N_(react).T_(Init) - N_(Prod).T_(ad))`

Consant Pressure

In constant pressure combustion, if air and fuel mixture burns adiabatically at constant pressure, the absolute enthalpy of the reactant at the initial state will be equal to the enthalpy of products at the final state.

`H_(Reactant)(T_i, P) = H_(Product)(T_(ad), P)`

The above equation may also be written as

`h_(Reactant)(T_i, P) = h_(Product)(T_(ad), P)` (per mass of mixture basis)

The above equation will be called it a constant-pressure adiabatic flame temperature.

Note:

In constant volume, the heat that is generated from combustion process going to be used to increase the internal energy. In constant pressure, the heat that is generated from combustion process going to be used to increase the enthaly.

• Stoichiometric equation

The general stoichiometric equation for a combustion of hydrocarbon can be represented as follows:

`C_xH_y+ar(O_2+3.76N_2)rightarrow xCO_2+(y/2)H_2O+(3.76ar)N2`

The term stoichiometric refers to the presence of just oxygen to completly burn all the fuel. The AFT is obtained at the stoichiometric ratio is maximum. However such ideal conditions are rarely observed in combustion process and therefore there is always some additional species leftover at the end of the process.

• Fuel-air equivalence ratio

The equivalence ratio , `phi`, is commonly used to indicate quantitatively whether a fuel-oxidizer is rich, lean, or stoichiometric. The equivalence is defined as follows

`phi = frac{frac{Fuel}{Air}}((frac{Fuel}{Air})_(st)`

From this definiton, we see that for fuel-rich mixture, `phi>1`, and for fuel lean mixture, `phi<1`. For stoichiometric combution, `phi` equals unity. In many combustion applications, the equivalence ratio is the single most important factor in determining the system's performance. Other parameters frequently used to define relative stoichiometric are present stoichiometric air, which relates to the equivalence ratio as

% stoichiometric air = `frac{100%}{phi}`

and present excess air, or

% excess air = ``frac{(1-phi)}{phi} .100``

• General equation of combustion process

The general equation for the combustion of a hydrocarbon can be represented as follows:

`C_xH_y+g(O_2+3.76N_2)rightarrow aCO_2+bH_2O+cN2+dO2+eCO`

• For stoichiometric Mixture (ϕ=1)
• For Lean Mixture (ϕ<1)
• For Rich Mixture (ϕ>1)
• The Enthalpy of the system

The enthalpy of the system depends on the species involved in the process and the temperature of the system. The enthalpy of the system can be calculated using the formula

`H = (a_1+a_2frac{T}{2}+a_3frac{T^2}{2}+a_4frac{T^3}{3}+a_5frac{T^4}{2}+frac{a_6}{T}).R.T`

Where a_1, a_2, a_3, a_4, a_5, a_6 `rightarrow` The numerical coefficients supplied in "http://combustion.berkeley.edu/gri-mech/version30/files30/thermo30.dat".

• Newton Raphson Method:

It is a root-finding algorithm which produces successively better approximations to the roots of a real-valued function.

The general formula for Newton Raphson method can be given as-

`x_(n+1) = x_n -alphafrac{f(x)}{fprime(x)`

where f(x) `rightarrow` Real-valued function: `fprime(x)``rightarrow` Derivative of the function; `alpha rightarrow` Multiplicity of the root

Case 1 - Effect of equivalence ratio on adiabatic flame temperature

1.Process: Combustion of hydrocarbon in a constant volume chamber

2. Hydrocarbon:

• Methane: x=1; y=4

3. Steps involved in python Code:

NumPy: To perform numerical analysis.

Matplotlib.pyplot: To create charts for better visualization of the results.

Cantera: It is a standard chemical kinetics library that contains standard chemical mechanisms GRI mechanisms, thermodynamic properties of several species that are included in a reaction.

1. First of all, i have created the range of equivalence ratios that will be participating in the reaction using linspace in-build function ranging from 0.1 till 2 where 200 points will be taken into calculations.
2. Next we have defned the functions, which will help to get enthalpy of reactants and products by using the eqauations from NASA polynomials form `H = (a_1+a_2frac{T}{2}+a_3frac{T^2}{2}+a_4frac{T^3}{3}+a_5frac{T^4}{2}+frac{a_6}{T}).R.T`.
3. Defined the low and high temperatures coefficients for participating species.
4. Next we defined the function which finding the root iteratively. Based on the range of equivalence ratio, this function calculate the enthalpies. As the reactants are at initial state, enthalpy of the reactants will be calculated using low temperature co-efficients.
   
5. In the next step, we defined the fprime by using the forward differencing scheme.
6. In the next step we setup the loop for performing the newton method for root finding by considering all values of equivalence ratio from which we will get AFT of methane at all values of phi.

Code 1:

"""
In this program we are going to calculate the AFT of Methane at different equivalence ratios

CH4 + 2(O2 + 3.76N2) = CO2 + 2H2O + 7.52N2

Program by - Anam | Skill Lync 2021
General Equation:
CxHy + g(O2 + 3.76N2) = aCO2 + bH2O + cN2 + dO2 + eCO
x = 1; y = 4; ar = x+y/4; g = ar*(1/phi)
Stoichiometric equation (phi=1):
CxHy + ar(O2 + 3.76N2) = xCO2 + (y/2)H2O + 3.76(x+y/4)N2
CH4 + 2(O2 + 3.76N2) = CO2 + 2H2O + 7.52N2
Lean Mixture Equation (phi<1):
A lean mixture of combustion of methane has an extra oxygen atom formed in the product side:
CH4 + (2/phi)(O2 + 3.76N2) = CO2 + 2H2O + (7.52/phi)N2 + ((2/phi)-2)O2
Rich Mixture Equation (phi>1)
Similarly for a rich mixture CO(carbon monoxide) is formed
CH4 + (2/phi)(O2 + 3.76N2) = ((4/phi)-3)CO2 + 2H2O + (7.52/phi)N2 + ((-4/phi)+4)CO
"""
import matplotlib.pyplot as plt
import cantera as ct
import math
import numpy as np
#import matplotlib.pyplot as plt

# To define the range of the equivalence ratio (0.1,0.2,.....2.0)
phi = np.linspace(0.5,2,101)
def h(T, co_effs):
	"""
		Function to evaluate enthalpy of all the reactant and products

	"""
	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 (Reactants) at STP condition
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_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 (Products)
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]
o2_coeffs_h = [3.28253784E+00,1.48308754E-03,-7.57966669E-07,2.09470555E-10,-2.16717794E-14,-1.08845772E+03,5.45323129E+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]

def f(T,phi):
	R = 8.314   #J/mol-K
	x = 1
	y = 4
	"""
		To Calculate the enthalpy of the reactants and products and to perform root finding problem
		Function that represents the Root finding problem
		For calculation of enthaly of products in Jouls

	"""
	if phi<=1:  # Conditional loop based on equivalence ratio enthalpy of products will be find out.
		h_co2_p = h(T,co2_coeffs_h)    # Enthalpy in joul basis
		h_h2o_p = h(T,h2o_coeffs_h)
		h_o2_p = h(T,o2_coeffs_h)
		h_n2_p = h(T,n2_coeffs_h)

		H_products = h_co2_p + 2*h_h2o_p  + ((2/phi)-2)*h_o2_p + (7.52/phi)*h_n2_p
		N_prod = 3 + (7.52/phi) + ((2/phi)-2)

        # For calculation of enthaly of reactant in Jouls
		Tstd = 298.15
		h_ch4_r = h(Tstd,ch4_coeffs_l)
		h_o2_r = h(Tstd,o2_coeffs_l)
		h_n2_r = h(Tstd,n2_coeffs_l)

		H_reactants = h_ch4_r + ((2/phi)*(h_o2_r + 3.76*h_n2_r))
		N_reactant = 1 + 4.76*(2/phi)

	else:
		h_co2_p = h(T,co2_coeffs_h)
		h_co_p = h(T,co_coeffs_h)
		h_h2o_p = h(T,h2o_coeffs_h)
		h_n2_p = h(T,n2_coeffs_h)
		h_o2_p = h(T,o2_coeffs_h)

		H_products = ((4/phi)-3)*h_co2_p  + (4-(4/phi))*h_co_p + 2*h_h2o_p + (7.52/phi)*h_n2_p    # Enthalpy in Joule basis
		N_prod = ((4/phi)-3) +2 + (4-(4/phi)) + (7.52/phi)

		Tstd = 298.15
		h_ch4_r = h(Tstd,ch4_coeffs_l)
		h_o2_r = h(Tstd,o2_coeffs_l)
		h_n2_r = h(Tstd,n2_coeffs_l)

		H_reactants = h_ch4_r + ((2/phi)*(h_o2_r + 3.76*h_n2_r))
		N_reactant = 1 + 4.76*(2/phi)
	return (H_products-H_reactants - R*((N_prod*T) - (N_reactant*Tstd)))

	#Newton Raphson Method to solve for AFT iteratively
def fprime(T,phi):
	# Define the numerical derivative of the function by using the forward differencing
	return (f(T+1e-6,phi)-f(T,phi))/1e-06

#Method solve for AFT till satisfying the criteria of tolerance
tol = 1e-03
#iter = 0
alpha = 0.9
Temp = []   # Define arrays to store the values of the temperature


for i in phi:
	T_guess = 1500
	R = 8.314 #J/mol-k
	while(abs(f(T_guess, i)) > tol):
		T_guess = T_guess - alpha*((f(T_guess, i)/fprime(T_guess, i)))
		#iter = iter+1

	Temp.append(T_guess)
print(max(Temp))
print(Temp)

plt.figure()
plt.plot(phi,Temp,color='red')
plt.xlabel('Equivalence ratio')
plt.ylabel('Adiabatic Flame Temperature')
plt.title('AFT comparison python vs cantera')
#plt.show()
#####################################################################################
#import cantera as ct
#import numpy as np
#import matplotlib.pyplot as plt
gas = ct.Solution('gri30.xml')

npoints = 101
phi = np.linspace(0.5,2,101)
tad = np.zeros(npoints)
fuel_species = 'CH4'
for i in range(npoints):
	gas.set_equivalence_ratio(phi[i], 'CH4', 'O2:2.0, N2:7.52')
	gas.TP = 298.15,101325
	gas.equilibrate('UV','auto')     # UV where we have fixed the volume
	tad[i] = gas.T

#tad.append(gas.T)
print('At phi = {0:10.4g}, Tad = {1:10.6g}'. format(phi[i], tad[i]))
plt.plot(phi, tad)
plt.legend(('Python', 'Cantera'))
plt.grid()
plt.show()

Results Discussions:

• In the above graph, we can easily observe that we got the maximum value of AFT at `phi = 1`in the python curve. But in the cantera curve we got the maximum value of AFT at `phi = 1.08`.
• The reason behind we got the high value in python case is that no of participating species at the product side using python is less compared to the cantera. That's why we got the high value in python case. Cantera has a standard chemical kinetic library, the no of species will be more at the product side.

2. Effect of Heat loss on adiabatic flame temperature

Case 2.1:

Process:

• Combustion of hydrocarbon in a constant pressure chamber.
• In this case, we are taken the equivalence ratio to 1.

We find out the effect of heat loss coefficient on AFT of methane(CH4) for constant pressure reactor because this is very close to the practical situation.

`Q = (m_f)*LHV`

Where Q is the heat loss

`m_f` is the mass of fuel. and `m_f = (LHV).(N_(fuel)). (M_(fuel)) = LHV.(m_(fuel))`

LHV is the lower heating value.

`N_(fuel)`is the no of moles of fuel.

`M_(fuel)` is the molecular weight of the fuel.

`m_(fuel)` is the mass of fuel.

The Heat loss co-efficient is vary from 0 to 1. When heat loss co-efficient is one then Heat loss is maximum. When heat loss co-efficient is 0 the Heat loss is minimum.

In the following code we study the effect of heat loss co-efficient.

Code 2:

"""
In this program we are going to calculate the AFT of Methane at different equivalence ratios

CH4 + 2(O2 + 3.76N2) = CO2 + 2H2O + 7.52N2

Program by - Anam | Skill Lync 2021
General Equation:
CxHy + g(O2 + 3.76N2) = aCO2 + bH2O + cN2 + dO2 + eCO
x = 1; y = 4; ar = x+y/4; g = ar*(1/phi)
Stoichiometric equation (phi=1):
CxHy + ar(O2 + 3.76N2) = xCO2 + (y/2)H2O + 3.76(x+y/4)N2
CH4 + 2(O2 + 3.76N2) = CO2 + 2H2O + 7.52N2
Lean Mixture Equation (phi<1):
A lean mixture of combustion of methane has an extra oxygen atom formed in the product side:
CH4 + (2/phi)(O2 + 3.76N2) = CO2 + 2H2O + (7.52/phi)N2 + ((2/phi)-2)O2
Rich Mixture Equation (phi>1)
Similarly for a rich mixture CO(carbon monoxide) is formed
CH4 + (2/phi)(O2 + 3.76N2) = ((4/phi)-3)CO2 + 2H2O + (7.52/phi)N2 + ((-4/phi)+4)CO
"""

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

"""
Fuel: CH4
Reactor: Constant pressure
phi = 1 stoichiometric condition
Maximum possible heat loss coefficient will be calculated as
LHV = 50 MJ/kg
HHV = 55 MJ/kg
MW = 16.04 gram/mol
Nfule = 1 mole
mfuel = 16.04 gram/mol
"""

x = 1 
y = 4
R = 8.314  #J/kmol-k
Q = 802000      #Q = (LHV)*mf = 50*e*6*16.04*e-03
H_loss = np.linspace(0,1,21)
H_product = []
H_reactants = []


def h(T, co_effs):
	"""
		Function to evaluate enthalpy of all the reactant and products

	"""
	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 (Reactants) at STP condition
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_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 (Products)
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]
o2_coeffs_h = [3.28253784E+00,1.48308754E-03,-7.57966669E-07,2.09470555E-10,-2.16717794E-14,-1.08845772E+03,5.45323129E+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]

def f(T,H_loss):
	"""
		Function to solve for heat balance and AFT
		To Calculate the enthalpy of the reactants and products and to perform root finding problem
		Function that represents the Root finding problem
		For calculation of enthaly of products in Jouls

	"""
	h_co2_p = h(T,co2_coeffs_h)    # Enthalpy in joul basis
	h_n2_p = h(T,n2_coeffs_h)
	h_h2o_p = h(T,h2o_coeffs_h)

	H_products = h_co2_p + 2*h_h2o_p + 7.52*h_n2_p 

        # For calculation of enthaly of reactant in Jouls
	Tstd = 298.15
	h_ch4_r = h(Tstd,ch4_coeffs_l)	 
	h_o2_r = h(Tstd,o2_coeffs_l)
	h_n2_r = h(Tstd,n2_coeffs_l)

	H_reactants = h_ch4_r + 2*h_o2_r + 7.52*h_n2_r

	return H_products-H_reactants + (H_loss*Q)

	#Newton Raphson Method to solve for AFT iteratively
def fprime(T,H_loss):
	# Define the numerical derivative of the function by using the forward differencing
	return (f(T+1e-6,H_loss)-f(T,H_loss))/1e-06

#Method solve for AFT till satisfying the criteria of tolerance
T_guess = 1500
tol = 1e-03
alpha = 0.9
Temp = []   # Define arrays to store the values of the temperature


for i in H_loss:
	T_guess = 1500
	R = 8.314 #J/mol-k
	while(abs(f(T_guess, i)) > tol):
		T_guess = T_guess - alpha*((f(T_guess, i)/fprime(T_guess, i)))

	Temp.append(T_guess)
print(max(Temp))
print(Temp)

plt.figure()
plt.plot(H_loss,Temp,color='red')
plt.xlabel('H_loss coefficients')
plt.ylabel('Adiabatic Flame Temperature')
plt.title('Variation of AFT with Heat loss coefficients')
plt.show()

• In the above graph we can easily observe that when heat loss is zero, adiabatic flame temperature is maximum. And when heat loss is maximum, adiabatic flame temperature is minmum.
• AFT varies linearly with the heat loss cofficient.
• For the heat loss coefficient of 0.35, the AFT is 1750K approximately.

Case 2.2: Effect of fuel type (Alkanes, Alkynes, and Alkenes) on adiabatic flame temperature

The number of carbon atoms, we conidered here is 2. For Alkane(C2H6), Alkenes(C2H4), Alkynes(C2H2).

Code 3:

"""
In this program we are going to calculate the AFT of Methane at different equivalence ratios

CH4 + 2(O2 + 3.76N2) = CO2 + 2H2O + 7.52N2

Program by - Anam | Skill Lync 2021
General Equation:
CxHy + g(O2 + 3.76N2) = aCO2 + bH2O + cN2 + dO2 + eCO
x = 1; y = 4; ar = x+y/4; g = ar*(1/phi)
Stoichiometric equation (phi=1):
CxHy + ar(O2 + 3.76N2) = xCO2 + (y/2)H2O + 3.76(x+y/4)N2
CH4 + 2(O2 + 3.76N2) = CO2 + 2H2O + 7.52N2
Lean Mixture Equation (phi<1):
A lean mixture of combustion of methane has an extra oxygen atom formed in the product side:
CH4 + (2/phi)(O2 + 3.76N2) = CO2 + 2H2O + (7.52/phi)N2 + ((2/phi)-2)O2
Rich Mixture Equation (phi>1)
Similarly for a rich mixture CO(carbon monoxide) is formed
CH4 + (2/phi)(O2 + 3.76N2) = ((4/phi)-3)CO2 + 2H2O + (7.52/phi)N2 + ((-4/phi)+4)CO
"""

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


def h(T, co_effs):
	"""
		Function to evaluate enthalpy of all the reactant and products

	"""
	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 for fuel Alkynes(C2h2), Alkenes(C2h4), Alkane(c2h6)
c2h2_coeffs_l = [8.08681094E-01,2.33615629E-02,-3.55171815E-05,2.80152437E-08,-8.50072974E-12,2.64289807E+04,1.39397051E+01]
c2h4_coeffs_l = [3.95920148E+00,-7.57052247E-03,5.70990292E-05,-6.91588753E-08,2.69884373E-11,5.08977593E+03,4.09733096E+00]
c2h6_coeffs_l = [4.29142492E+00,-5.50154270E-03,5.99438288E-05,-7.08466285E-08,2.68685771E-11,-1.15222055E+04,2.66682316E+00]
 #Low temperature coefficients (Reactants) at STP condition
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 (Products)
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]

Tstd = 298.15
h_C2H2_r = h(Tstd,c2h2_coeffs_l)	 
h_C2H4_r = h(Tstd,c2h4_coeffs_l)
h_C2H6_r = h(Tstd,c2h6_coeffs_l)
H_alkyne = h_C2H2_r
H_alkene = h_C2H4_r
H_alkane = h_C2H6_r
H = [H_alkyne,H_alkene,H_alkane]
H = np.array(H)
def f(T,H):
	"""
		To Calculate the enthalpy of the reactants and products and to perform root finding problem
		Function that represents the Root finding problem
		For calculation of enthaly of products in Jouls

	"""
	h_co2_p = h(T,co2_coeffs_h)    # Enthalpy in joul basis
	h_n2_p = h(T,n2_coeffs_h)
	h_h2o_p = h(T,h2o_coeffs_h)
 
	h_o2_r = h(Tstd,o2_coeffs_l)
	h_n2_r = h(Tstd,n2_coeffs_l)

	if H == H_alkyne:
		x = 2
		y = 2
		ar = x + (y/4)
		H_products = x*h_co2_p + (y/2)*h_h2o_p + (3.76*ar)*h_n2_p     #C2H2
		H_reactants = H_alkyne + ar*h_o2_r + (3.76*ar)*h_n2_r

	elif h == H_alkene:
		x = 2
		y = 4
		ar = x + (y/4)
		H_products = x*h_co2_p + (y/2)*h_h2o_p + (3.76*ar)*h_n2_p     #C2H4
		H_reactants = H_alkene + ar*h_o2_r + (3.76*ar)*h_n2_r 

	else :
		x = 2
		y = 6
		ar = x + (y/4)
		H_products = x*h_co2_p + (y/2)*h_h2o_p + (3.76*ar)*h_n2_p     #C2H6
		H_reactants = H_alkane + ar*h_o2_r + (3.76*ar)*h_n2_r

	return H_products-H_reactants 

	#Newton Raphson Method to solve for AFT iteratively
def fprime(T,H):
	# Define the numerical derivative of the function by using the forward differencing
	return (f(T+1e-6,H)-f(T,H))/1e-06

#Method solve for AFT till satisfying the criteria of tolerance
T_guess = 1500
tol = 1e-03
#iter = 0
alpha = 0.9
Temp = []   # Define arrays to store the values of the temperature

for i in H:
	T_guess = 1500
	R = 8.314 #J/mol-k
	while(abs(f(T_guess, i)) > tol):
		T_guess = T_guess - alpha*((f(T_guess, i)/fprime(T_guess, i)))
		#iter = iter+1

	Temp.append(T_guess)
print(max(Temp))
print(Temp)

plt.figure()	 
plt.bar(['C2H2', 'C2H4', 'C2H6'], Temp)
plt.ylabel('Adiabatic Flame Temperature')
plt.title('AFT vs Fuel type (python)')
plt.show()
########################################################################
import cantera as ct
import numpy as np
import matplotlib.pyplot as plt
#CH4 + 2(O2 + 3.76N2) = CO2 + 2H2O + 7.52N2

gas = ct.Solution('gri30.xml')

gas.TPX = 298.15,101325,{'C2H2':1, 'O2':2.5, 'N2':9.4}
gas.equilibrate('HP','auto')
T_c2h2 = gas.T

gas.TPX = 298.15,101325,{'C2H4':1, 'O2':3, 'N2':11.28}
gas.equilibrate('HP','auto')
T_c2h4 = gas.T

gas.TPX = 298.15,101325,{'C2H6':1, 'O2':3.5, 'N2':13.16}
gas.equilibrate('HP','auto')
T_c2h6 = gas.T

temp = (T_c2h2, T_c2h4, T_c2h6)
temp = np.array(temp)
print(temp)

plt.figure()	 
plt.bar(['C2H2', 'C2H4', 'C2H6'], temp)
plt.ylabel('Adiabatic Flame Temperature')
plt.title('AFT vs Fuel type(cantera)')
plt.show()

• In this section, we are observing the adiabatic flame temperature versus different fuels.
• in the above figure we note that Alkyne(C2H2) fuel got a high adiabatic flame temperature which is approximately 2900K.
• The reason behind Alkyne gave the maximum temperature is that it contains the less number of hydrogen atoms. The fuel which contains the less number of hydrogen atoms gives the maximum temperature.
• The trend with the type of fuel is the same for python as well as cantera. But in cantera case temperature value differs by atleast 400K which is because cantera solves for all species.

Case:2.3: Effect of carbon atoms(1 to 3) with alkane fuel on Adiabatic flame temperature:

Code 4:

"""
In this program we are going to calculate the AFT of Methane at different equivalence ratios

CH4 + 2(O2 + 3.76N2) = CO2 + 2H2O + 7.52N2

Program by - Anam | Skill Lync 2021
General Equation:
CxHy + g(O2 + 3.76N2) = aCO2 + bH2O + cN2 + dO2 + eCO
x = 1; y = 4; ar = x+y/4; g = ar*(1/phi)
Stoichiometric equation (phi=1):
CxHy + ar(O2 + 3.76N2) = xCO2 + (y/2)H2O + 3.76(x+y/4)N2
CH4 + 2(O2 + 3.76N2) = CO2 + 2H2O + 7.52N2
Lean Mixture Equation (phi<1):
A lean mixture of combustion of methane has an extra oxygen atom formed in the product side:
CH4 + (2/phi)(O2 + 3.76N2) = CO2 + 2H2O + (7.52/phi)N2 + ((2/phi)-2)O2
Rich Mixture Equation (phi>1)
Similarly for a rich mixture CO(carbon monoxide) is formed
CH4 + (2/phi)(O2 + 3.76N2) = ((4/phi)-3)CO2 + 2H2O + (7.52/phi)N2 + ((-4/phi)+4)CO
"""

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


def h(T, co_effs):
	"""
		Function to evaluate enthalpy of all the reactant and products

	"""
	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 for fuel Alkynes(C2h2), Alkenes(C2h4), Alkane(c2h6)
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_l = [4.29142492E+00,-5.50154270E-03,5.99438288E-05,-7.08466285E-08,2.68685771E-11,-1.15222055E+04,2.66682316E+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]
 #Low temperature coefficients (Reactants) at STP condition
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 (Products)
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]

Tstd = 298.15
h_CH4_r = h(Tstd,ch4_coeffs_l)	 
h_C2H6_r = h(Tstd,c2h6_coeffs_l)
h_C3H8_r = h(Tstd,c3h8_coeffs_l)
C_1 = h_CH4_r
C_2 = h_C2H6_r
C_3 = h_C3H8_r
C = [C_1,C_2,C_3]
C = np.array(C)
def f(T,C):
	"""
		To Calculate the enthalpy of the reactants and products and to perform root finding problem
		Function that represents the Root finding problem
		For calculation of enthaly of products in Jouls

	"""
	h_co2_p = h(T,co2_coeffs_h)    # Enthalpy in joul basis
	h_n2_p = h(T,n2_coeffs_h)
	h_h2o_p = h(T,h2o_coeffs_h)
 
	h_o2_r = h(Tstd,o2_coeffs_l)
	h_n2_r = h(Tstd,n2_coeffs_l)

	if C == C_1:
		x = 1
		y = 4
		ar = x + (y/4)
		H_products = x*h_co2_p + (y/2)*h_h2o_p + (3.76*ar)*h_n2_p     #C2H2
		H_reactants = C_1 + ar*h_o2_r + (3.76*ar)*h_n2_r

	elif C == C_2:
		x = 2
		y = 6
		ar = x + (y/4)
		H_products = x*h_co2_p + (y/2)*h_h2o_p + (3.76*ar)*h_n2_p     #C2H4
		H_reactants = C_2 + ar*h_o2_r + (3.76*ar)*h_n2_r

	elif C == C_3:
		x = 3
		y = 8
		ar = x + (y/4)
		H_products = x*h_co2_p + (y/2)*h_h2o_p + (3.76*ar)*h_n2_p     #C2H6
		H_reactants = C_3 + ar*h_o2_r + (3.76*ar)*h_n2_r

	return H_products-H_reactants 

	#Newton Raphson Method to solve for AFT iteratively
def fprime(T,H):
	# Define the numerical derivative of the function by using the forward differencing
	return (f(T+1e-6,H)-f(T,H))/1e-06

#Method solve for AFT till satisfying the criteria of tolerance
T_guess = 1500
tol = 1e-03
#iter = 0
alpha = 0.9
Temp = []   # Define arrays to store the values of the temperature


for i in C:
	T_guess = 1500
	R = 8.314 #J/mol-k
	while(abs(f(T_guess, i)) > tol):
		T_guess = T_guess - alpha*((f(T_guess, i)/fprime(T_guess, i)))
		#iter = iter+1

	Temp.append(T_guess)
print(max(Temp))
print(Temp)

plt.figure()
plt.bar(['CH4', 'C2H6', 'C3H8'], Temp)
plt.ylabel('Adiabatic Flame Temperature')
plt.title('AFT vs Fuel type (Effect of carbon atom)')
plt.show()

• The fuel which contains more carbon atom gives the maximum value of adiabatic flame temperature.

Conclusion:

• In this challenge, i have studied about Adiabatic flame temperature and its python and cantera coding.
• I have observed that AFT will be impacted with changing the euivalence ratios (constant volume chamber), type of fuel, number of carbon atom of alkane.
• I have also studied the effect of heat loss coefficient on AFT of constant prssure chamber where combustion is being done using methane.