Sensitivity analysis of GRI reaction mechanisms

SENSITIVITY ANALYSIS OF GRI MECHANISM REACTIONS

Sensitivity analysis is the method used to identify how much a specific parameter affects the overall property of any mechanism, i.e, how sensitive the mechanism is to the given parameter and its changes. This analysis is performed here to find out how temperature changes in a system of chemical equations will affect its products and rate of reactions accordingly.

Sensitivity analysis is an important parameter that can be used to accurately determine the effect of certain intermediate reactions and how it affects the overall final product obtained from the reaction. The process of recalculating outcomes under alternative assumptions to determine the impact of a variable under sensitivity analysis can be useful for a range of purposes, including:

• Testing the robustness of the results of a model or system in the presence of uncertainty.
• Increased understanding of the relationships between input and output variables in a system or model.
• Uncertainty reduction, through the identification of model inputs that cause significant uncertainty in the output and should therefore be the focus of attention in order to increase robustness (perhaps by further research).

In our project, we will take the stoichiometric combustion equation of methane and using the GRI 3.0 mechanism, we will derive the intermediate equations and calculate its sensitivites to see the top 10 reactions that are the most temperature sensitive using Cantera. We also tabulate the molar concentrations changes of OH, H and CH4 species across time for the given temperature condition.

CODING

\'\'\'
Calculating the sensitivity of individual intermediate reactions
Combustion equation: CH4 + 2(O2 + 3.76N2) -> CO2 + 2H2O + 7.52N2
\'\'\'
import sys
import matplotlib.pyplot as plt 
import cantera as ct 
import numpy as np 
from tabulate import tabulate

# Initializing cantera gas object
gas = ct.Solution(\'gri30.xml\')

# Inputs
T = 1500
P = ct.one_atm
gas.TPX = T,P,\'CH4:1, O2:2, N2:7.52\'
r = ct.IdealGasConstPressureReactor(gas)
sim = ct.ReactorNet([r])

# Initializing variables using gas object
nr = gas.n_reactions # syntax to acquire total number of reactions in the mechanism
eq = gas.reaction_equations() # syntax to acquire all reaction equations in the mechanism
rid = np.arange(1,nr+1) # arranging reaction IDs in an array
n = 10 # Top 10 sensitive reactions to consider

# enable sensitivity with respect to the rate of first 10
for i in range(nr):
	r.add_sensitivity_reaction(i)

# Relative and absolute tolerances for the reaction mechanism PDE\'s
sim.rtol = 1e-6
sim.atol = 1e-15

# Relative and absolute tolerances for the sensitivity mechanism PDE\'s
sim.rtol_sensitivity = 1e-6
sim.atol_sensitivity = 1e-6

# Define array to hold maximum sensitivity values
sens_max = [0]*nr
states = ct.SolutionArray(gas,extra =[\'t_ms\',\'t\'])

# for loop to calculate sensitivities and derive values
for t in np.arange(0,2e-3,5e-6):
	sim.advance(t)
	states.append(r.thermo.state, t_ms = 1e3*t,t = t)
	
	for j in range(nr):
		sens = sim.sensitivity(\'temperature\', j)

		# Loop to derive sensitivity values
		if np.abs(sens)>=np.abs(sens_max[j]):
			sens_max[j] = sens

# syntax to sort given values under a given condition - based upon absolute values of sensitivity
eq=[r for s,r in sorted(zip(np.abs(sens_max),eq),reverse=True)]
rid=[r for s,r in sorted(zip(np.abs(sens_max),rid),reverse=True)]
sens_max=[sens for abs_sens,sens in sorted(zip(np.abs(sens_max),sens_max),reverse=True)]

eqten = eq[0:n]
ridten = rid[0:n]
sensten = sens_max[0:n]

# syntax to formulate a table for the top 10 temperature sensitive reactions
table = tuple(zip(ridten,eqten,sensten))
print(tabulate(table,headers=[\'Reaction ID\',\'Equations\',\'Sensitivity\'],tablefmt=\'orgtbl\'))

# plots for Molar concentrations
plt.figure(1)
plt.subplot(2,2,1)
plt.plot(states.t,states.T)
plt.xlabel(\'Time (ms)\')
plt.ylabel(\'Temperature (K)\')
plt.subplot(2,2,2)
plt.plot(states.t,states(\'OH\').X)
plt.xlabel(\'Time (ms)\')
plt.ylabel(\'OH Mole fraction\')
plt.subplot(2,2,3)
plt.plot(states.t,states(\'H\').X)
plt.xlabel(\'Time (ms)\')
plt.ylabel(\'H Mole fraction\')
plt.subplot(2,2,4)
plt.plot(states.t,states(\'CH4\').X)
plt.xlabel(\'Time (ms)\')
plt.ylabel(\'CH4 Mole fraction\')

eqten = list(reversed(eqten))
sensten = list(reversed(sensten))

# Bar plot for top 10 temperature sensitive reactions in the mechanism
plt.figure(2)
plt.barh(eqten,sensten)
plt.xlabel(\'Sensitivity\')
plt.ylabel(\'Reactions\')
plt.title(\'Top 10 temperature sensitive reactions in the GRI mechanism 3.0\')
plt.tight_layout()
plt.show()

Using the help command for the gas object, we can find how many reactions are taken into account in our program and we can define variables to accomodate the reactions and their equations to calculate their sensitivites.

We use the sorted function to sort the sensitivity values derived using the sensitivity object under the sim function in our Reactor network. After that, the tuple object is used to make a table of top 10 temperature sensitive reactions from the sorted group of new reaction IDs, equations and sensitivities. The barh function is used to make a horizontal plot as demonstrated in the video and use Y-axis as the alignment feature to show our sensitivities in both negative and positive directions.

RESULTS

Sensitivity Bar graph:

The probability of reactions cascading is determined by its positive sensitivity and how sensitive it is to temperature changes. If the sensitivity is negative, then the reaction is less likely to occur within the specified conditions.

Molar concentrations of OH, H and CH4:

Molar concentrations changes denote the point of auto-ignition due to the temperature condition provided and the species concentration changes with respect to the products and reactant composition for our reaction respectively. Methane is being consumed in the reaction for combustion, hence after auto-ignition, the value decreases while the OH and H molar concentrations increase due to being by-products of the combustion reaction after ignition.