Week 9 - Senstivity Analysis Assignment

Objective:

To write a code to calculate the most sensitive rections among all reactions from GRI mechanism considering the sensitivity parameter as Temperature.

Introduction:

1.Sensitivity Analysis:

Sensitivity analysis is a process to get the relationship between the parameters and output of a model.The main aim is to calculate the ignition delay of all reactions in GRI mechanisms and calculate the 10 most sensitive reactions.The sensitvity analysis will indicate which reactions are sensitive to ignition delay.

This will be helpful in tuning the reactions to reduce the Ignition delay by altering the rate coefficient of the reactions .The reaction mechanism contains a lots of reactions by performing sensitvity analysis we can find out the reactions which are most sensitive and tune the reactions based on the sensitive parameters like Temperature,species.

To find:

• To find the 'n' number of sensitive reactions from the Gri mechanism based on the sensitive parameter as temperature.
• Plot the reaction strings in order of sensitivity.
• To plot the sensitivities vs reactions.

Coding Approach:

• The GRI mechanism consits of 325 reactions and 53 species.
• The sensitivties of all the 325 reaxtions must be calculated using the program.
• Import the necessary modules like cantera,numpy and matplotlib.
• Instantiate the solution class and create the gas object and specify the inputs like temperature,pressure and mole fractions.
• Create the IdealgasReactor and pass the gas object .
• Create a reactor network to hold the reactor, the reactors are passed as a list to the reactor network.
• The number of reactions are obtained by using the gas.n_reactions attribute.
• The equations are obtained by the gas.reaction equations.
• The sensitivity parameter temperature is added for all the reactions.
• The most sensitive reactions that affect the auto ignition are plotted.

Code:

"""
Program to calculate the most sensitvity reactions that affect the auto ignition
"""
import sys
import cantera as ct 
import numpy as np 
import matplotlib.pyplot as plt 
#Inputs
gas=ct.Solution('gri30.xml')
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])
#Number of reactions in GRI mechanism
nr=gas.n_reactions
print(nr)
#Number of sensitive reactions to be considered
N=10
max_s=[0]*nr
plot=[0]*N
for i in range(nr):
    r.add_sensitivity_reaction(i)
# tolerance for the solution 
sim.rtol = 1.0e-6
sim.atol = 1.0e-15
# tolerance for the sensitivity coefficients
sim.rtol_sensitivity = 1.0e-6
sim.atol_sensitivity = 1.0e-6
#solution array that contains thermodynamic values
states = ct.SolutionArray(gas, extra='t')
#Advancing the state of reactornetwork
for t in np.arange(0,2e-3,5e-6):
    sim.advance(t)
    for j in range(nr):
    	s=sim.sensitivity('temperature',j)
    	states.append(r.thermo.state,time_in_ms=1000*t,t=t)
    	if (np.abs(s)>np.abs(max_s[j])):
    		max_s[j]=s
#list to store the sensitive reactions
reac_lab=[]
#loop to extract the 10 most sensitive reactions
for k in range(N):
	maximum=np.max(max_s)
	minimum=np.min(max_s)
	if np.abs(maximum)>np.abs(minimum):
		plot[k]=maximum
		reac_index=max_s.index(maximum)
		reac_lab.append(sim.sensitivity_parameter_name(reac_index))
		max_s.remove(maximum)
	else:
		plot[k]=minimum
		reac_index=max_s.index(minimum)
		reac_lab.append(sim.sensitivity_parameter_name(reac_index))
		max_s.remove(minimum)

for r in range(N):
	reac_lab[k]=reac_lab[k].replace('IdealGasConstPressureReactor_1:','')
reactions=reac_lab
sensitivity=plot
plt.figure(1)
plt.barh(reactions,sensitivity)
plt.xlabel("sensitivity")

plt.figure(2)
plt.subplot(2,2,1)
plt.plot(states.t,states.T)
plt.xlabel('Time in ms')
plt.ylabel('temperature in K')
plt.subplot(2,2,2)
plt.plot(states.t,states('OH').X)
plt.xlabel('Time in ms')
plt.ylabel('OH concentration')
plt.subplot(2,2,3)
plt.plot(states.t,states('H').X)
plt.xlabel('Time in ms')
plt.ylabel('H concentration')
plt.subplot(2,2,4)
plt.plot(states.t,states('CH4').X)
plt.xlabel('Time in ms')
plt.ylabel('CH4 concentration')
plt.show()

Plots:

From the above plots `CH_3+O_2<=> CH_3O+O` has highest sensitivity(positive) which indicates that increase in this rate of reaction results in decrease in ignition delay.

`CH_3+H_2O_2<=> CH_4+HO2` has has highest sensitivity(negative) which indicates that increase in this rate of reaction results in increase in ignition delay.

The plots show the variation of ignition delay with respect to temperature,concentration of OH,H,CH4.