Week 3 - Adiabatic Flame Temperature calculation
#equivalence ratio of AFT and methane
#aft of methane
# ch4+2(2O2+3.76N2)=CO2+2H2O+7.52N2
import math
import matplotlib.pyplot as plt
def h(T,coeffs):
# function to evaluate enthalphy
R = 8.134 #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)*R*T
# defining the coefficient
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
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
N2_coeffs_h = 0.02926640E+02,0.14879768E-02,-0.05684760E-05,0.10097038E-09,-0.06753351E-13,-0.09227977E+04,0.05980528E+02
def f(T):
#function that represents the root finding pb
h_co2_p = h(T,co2_coeffs_h)
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
print(H_products)
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
print(H_reactants)
return H_products-H_reactants
def fprime(T):
return(f(T+1e-6) - f(T))/1e-6
T_guess = 1500
tol = 1e-3
iter = 0
alpha = 0.2
while (abs(f(T_guess))>tol):
T_guess = T_guess-alpha*((f(T_guess)/fprime(T_guess)))
plt.plot(iter,T_guess,'*',color = 'red')
iter = iter+1
plt.title('equ ratio on aft amd methane')
plt.xlabel('equ ratio')
plt.ylabel('T')
print(iter)
print(T_guess)
plt.show()
----------------------------------------------------------------
#aft heat loss from a reaction with constant pressure
#function to enthalphy at constant presuure
# ch4+2(2O2+3.76N2)=CO2+2H2O+7.52N2
import matplotlib.pyplot as plt
import math
import numpy as np
def h(T,coeffs):
R = 8.134 #J/Mol-k
a1 = coeffs[1]
a2 = coeffs[2]
a3 = coeffs[3]
a4 = coeffs[4]
a5 = coeffs[5]
a6 = coeffs[6]
return (a1+a2*T/2+a3*pow(T,2)/3+a4*pow(T,3)/4+a5*pow(T,4)/5+a6/T)*R*T
# defining the coefficients
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
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
N2_coeffs_h = 0.02926640E+02,0.14879768E-02,-0.05684760E-05,0.10097038E-09,-0.06753351E-13,-0.09227977E+04,0.05980528E+02
def f(T,Heat_loss):
# fun to evaluate enthalphy at constant pressure
Tstd = 298.15 #k
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
print(H_reactants)
h_co2_p = h(T,co2_coeffs_h)
h_H2o_p = h(T,H2o_coeffs_h)
h_N2_p = h(T,N2_coeffs_h)
H_products = h_co2_p+2*h_H2o_p+7.52*h_N2_p
return H_reactants - H_products
def fprime(T,Heat_loss):
return (f(T+1e-6,Heat_loss)-f(T,Heat_loss))/1e-6
H_loss = 0
T_aft = []
H_l = []
while (H_loss<1):
T_guess = 1500
tol = 17691804E-03
iter = 0
alpha = 1
while (abs(f(T_guess,Heat_loss[i]))>tol):
T_guess = T_guess - alpha*((f(T_guess,Heat_loss[i])/fprime(T_guess,Heat_loss[i])))
plt.plot(H_loss,T_guess,'*',color = 'red')
iter = iter + 1
T_ft.append(T_guess)
plt.plot(H_l,T_aft)
plt.xlabel('Heat_loss')
plt.ylabel('Flame_temparature')
print(iter)
print(T_guess)
plt.show()
result:
from this program we can see the effect of heat loss on AFT.
heat loss was taken a fraction of the max possible heat loss,as the heat loss is high flame temparature is low.
--------------------------------------------------------------------------------------------
#aft of methane
#C2Hy+(2+y/4)(O2+3.76N2) ----> 2CO2+(y/2)H20+3.76*(2+y/4)N2
import math
import matplotlib.pyplot as plt
import numpy as np
def h(T,coeffs):
# function to evaluate enthalphy
R = 8.134 #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)*R*T
# defining the coefficient
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
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
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
N2_coeffs_h = 0.02926640E+02,0.14879768E-02,-0.05684760E-05,0.10097038E-09,-0.06753351E-13,-0.09227977E+04,0.05980528E+02
def f(T,n):
Tstp = 298.15 #k
R = 8.134 #J/Mol-k
y = n
h_C2H2_r = h(Tstp,C2H2_coeffs_l)
h_C2H4_r = h(Tstp,C2H4_coeffs_l)
h_C2H6_r = h(Tstp,C2H6_coeffs_l)
h_O2_r = h(Tstp,O2_coeffs_l)
h_N2_r = h(Tstp,N2_coeffs_l)
H_loss = 0.35
if n==2: #ethylene LHV = 1258020 J/MOL
H_reactants = h_C2H2_r+(2+y/4)*h_O2_r+3.76*(2+y/4)*h_N2_r - (1258020*H_loss)
if n==4:
H_reactants = h_C2H4_r+(2+y/4)*h_O2_r+3.76*(2+y/4)*h_N2_r - (1323474*H_loss)
if n==6:
H_reactants = h_C2H6_r+(2+y/4)*h_O2_r+3.76*(2+y/4)*h_N2_r - (1428275*H_loss)
h_CO2_p = h(T,CO2_coeffs_h)
h_H2O_p = h(T,H2O_coeffs_h)
h_N2_p = h(T,N2_coeffs_h)
H_produts = 2*h_CO2_p+(y/2)*h_H2O_p+3.76*(2+y/4)*h_N2_p
return H_produts - H_reactants
def fprime(T,n):
return (f(T+1e-6,n)-f(T,n))/1e-6
T_guess = [1500 for i in range(3)] #k
n = [2,4,6]
tol = 1e-3
iter = 0
alpha = 0.2
i = 0
for i in range(0,3):
while (abs(f(T_guess[i],n[i]))>tol):
T_guess[i] = T_guess[i]-alpha*(f(T_guess[i],n[i])/fprime(T_guess[i],n[i]))
iter = iter + 1
plt.plot(n[i],T_guess[i] ,'o',color = 'red')
plt.title('heat loss of 0.35')
plt.xlabel('C2HX')
plt.ylabel('T')
print(T_guess)
print(iter)
plt.show()
result:
we can see that no of hydrogen atoms is double of hydrocarbons.while AFT decreases no of hydrogen atoms increases .At hydrogen atoms less water molecules is formed which means heat loss is less
----------------------------------------------------------------------------------------------------
#aft of methane
# ch4+2(2O2+3.76N2)=CO2+2H2O+7.52N2
#CxHy+(x+y/4)(O2+3.76N2) ----> xCO2+(y/2)H20+3.76*(X+y/4)N2
import math
import matplotlib.pyplot as plt
import numpy as np
def h(T,coeffs):
# function to evaluate enthalphy
R = 8.134 #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)*R*T
# defining the coefficient
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 = 2.13583630E+00,1.81188721E-02,-1.73947474E-05,9.34397568E-09,-2.01457615E-12,-7.04291804E+03,1.22156480E+01
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
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
N2_coeffs_h = 0.02926640E+02,0.14879768E-02,-0.05684760E-05,0.10097038E-09,-0.06753351E-13,-0.09227977E+04,0.05980528E+02
def f(T,n):
Tstp = 298.15 #k
R = 8.134 #J/Mol-k
x = n
y = 2*n+2
h_CH4_r = h(Tstp,CH4_coeffs_l)
h_C2H6_r = h(Tstp,C2H6_coeffs_l)
h_C3H8_r = h(Tstp,C2H6_coeffs_l)
h_O2_r = h(Tstp,O2_coeffs_l)
h_N2_r = h(Tstp,N2_coeffs_l)
H_loss = 0.35
if n==1: #ethylene LHV = 1258020 J/MOL
H_reactants = h_CH4_r+(x+y/4)*h_O2_r+3.76*(x+y/4)*h_N2_r - (1258020*H_loss)
if n==2:
H_reactants = h_C2H6_r+(x+y/4)*h_O2_r+3.76*(x+y/4)*h_N2_r - (1323474*H_loss)
if n==3:
H_reactants = h_C3H8_r+(x+y/4)*h_O2_r+3.76*(x+y/4)*h_N2_r - (1428275*H_loss)
h_CO2_p = h(T,CO2_coeffs_h)
h_H2O_p = h(T,H2O_coeffs_h)
h_N2_p = h(T,N2_coeffs_h)
H_produts = x*h_CO2_p+(y/2)*h_H2O_p+3.76*(x+y/4)*h_N2_p
return H_produts - H_reactants
def fprime(T,n):
return (f(T+1e-6,n)-f(T,n))/1e-6
T_guess = [1500 for i in range(3)] #k
n = [1,2,3]
tol = 1e-3
iter = 0
alpha = 0.2
i = 0
for i in range(0,3):
while (abs(f(T_guess[i],n[i]))>tol):
T_guess[i] = T_guess[i]-alpha*(f(T_guess[i],n[i])/fprime(T_guess[i],n[i]))
iter = iter + 1
plt.plot(n[i],T_guess[i] ,'o',color = 'red')
plt.title('effect of carbon atoms on aft')
plt.xlabel('no of carbons in CxHy')
plt.ylabel('T')
print(T_guess)
print(iter)
plt.show()
result :
from this program we can see that with the carbon atoms increases AFT also increases