Cantera Open-Source

Reaction Reduction Mechanism using Cantera Open-Source

Objective-To write a Reaction Reduction Mechanism code for the GRI3.0 mechanism.

The code structure will do the following,

  1. First, it should run a reference Ignition delay and Max temperature Analysis for the unreduced mechanism for a set of state values (Temperature, Pressure and Equivalence ratio) for a time duration of (5 sec).
  2. Then, a sensitivity analysis with respect to temperature must be run on all the reactions in the mechanism and must sort the reactions from the most sensitive to least sensitive.
  3. Then, these sorted reactions should be added to a new mechanism file/object one by one starting from a fixed minimum number of reactions (40).
  4. Now a graph must be plotted plotting the reduced mechanism’s ignition delay and max temperature while having  the reference ignition delay and Temperature as a label for reference (with 5% tolerance).
  5. This graph must have the values of Ignition Delay and Max temperature on the Y-axis, while the number of reactions in the mechanism on the X axis.
  6. Now, repeat the same steps for various state values.

State Values Range –

Temperature – [950 , 1050 , 1150] Kelvin

Pressure – [3 , 4 , 5] bar

Equivalence ratio – [0.5 , 1 , 1.5]

Calculations must be performed with all possible combinations of State Values.

Coding Approach:

Program to read text file:

File parsing technique is used to read the state values from a text file for temperature,pressure and equivalence ratio.

Main Program:

1. Ignition delay and Max temperature Analysis for the unreduced mechanism for a set of state values (Temperature, Pressure and Equivalence ratio).

2. For determining ignition delay and maximum temperature a instance gas object is created by passing the gri30.cti mechanism.

3. A reactor object ‘r’ is created using ‘IdealGasReactor’ and gas object is passed.

4. A reactor network ReactorNet is created by passing reactor to it and solution is advanced.

5. Similarly sensitivity of all the 325 reactions in GRI mechanism file is simulated using add_sensitivity_reaction().The relative and absolute tolerance is set to check convergence.

6. Sensitivities of reactions are sorted from higher to lower values.

7. These sorted reactions should is added to a new mechanism file/object one by one starting from a fixed minimum number of reactions (40).

8. A 5% tolerance is set for both ignition delay and temperature with the reference values.Finally the minimum number of reaction is determined with 5% tolerance level.

Code:

Text File named as new.txt

Temperature 950 1050 1150
Pressure    300000 400000 500000
Phi         0.5  1  1.5

Main Code:

import os
import sys
import cantera as ct
import numpy as np
import matplotlib.pyplot as plt


#program to read text file

class Input_file_reader():
	def __init__(self,file_name):
		self.file_name = file_name
		self.input_dict={}

	def validate_input_file(self):
		if os.path.isfile(self.file_name):
			return True
		else:
			return False

	def read(self):
		if self.validate_input_file():
			for line in open(self.file_name,'r'):
				self.input_dict[line.split()[0]]= line.split()[1:]
		else:
			print('Validation Failed')

a=Input_file_reader('new.txt')
a.read()
state_temp=[float(x) for x in a.input_dict['Temperature']]
state_pres=[float(x) for x in a.input_dict['Pressure']]
state_phi=[float(x) for x in a.input_dict['Phi']]

#Calculating reference ignition delay and Maximum temperature for all state values and sorting sensitivity 

spec =ct.Species.listFromFile('gri30.xml')
rxns =ct.Reaction.listFromFile('gri30.xml')
gas =ct.Solution('gri30.cti')
state_pres=[300000,400000,500000]
state_phi=[0.5,1,1.5]
state_temp=[950,1050,1150]
n_number=gas.n_reactions

for temp in state_temp:
	for pres in state_pres:
		for phi in state_phi:
			gas.TPX= temp,pres,{'CH4':1,'O2':2/phi,'N2':7.52/phi}

			r=ct.IdealGasReactor(gas)
			sim=ct.ReactorNet([r])

			for i in range(n_number):
				r.add_sensitivity_reaction(i)

			sim.rtol=1e-6
			sim.atol=1e-15

			sim.rtol_sensitivity=1e-6
			sim.atol_sensitivity=1e-6
			n=gas.n_reactions
			s_max=[0]*n_number
			states=ct.SolutionArray(gas,extra=['t'])
			for t in np.arange(0,1,1e-3):
				sim.advance(t)
				states.append(r.thermo.state,t=1000*t)
				if (gas.T<=temp+400):
						t_ref=t

				for j in range(n):
					s2=sim.sensitivity('temperature',j)
					if np.abs(s2)>=np.abs(s_max[j]):
						s_max[j]=s2	
			T_maxref=gas.T			
			s_max_abs=np.abs(s_max)
			s_maximum=sorted(s_max_abs,reverse=True)
			sensitivity_final=[]
			indexr=[]
			eq=[]
			for i in range(n_number):
				for j in range(n):
					if np.abs(s_max[j])==s_maximum[i]:
						indexr.append(j)
						eq.append(rxns[j])
			
			#Calculating  ignition delay and Maximum temperature for resuced mechanism
			rxn_reduced=[]
			ign_delay=[]
			T_max=[]
			num=[]
			n_1=1
			tol_check=1
			error_tol=5
			error_time=100
			error_temp=100
			for i in range(0,n_number):
				
				mini_reac=40
				rxn_reduced.append(eq[i])
				gas1=ct.Solution(thermo='IdealGas',kinetics='GasKinetics',species=spec,reactions=rxn_reduced)
				auto_ignition_time=[]
				temp_final=[]
				gas1.TPX= temp,pres,{'CH4':1,'O2':2/phi,'N2':7.52/phi}

				r1=ct.IdealGasReactor(gas1)
				sim1=ct.ReactorNet([r1])
				time=0
				#states=ct.SolutionArray(gas1,extra=['tim'])
				for n in range(1000):
					time +=1e-3
					sim1.advance(time)
					
					#states.append(r1.thermo.state,tim=time)
					if (gas1.T<=temp+400):
						t=time
						#temp=(states.T[n])
				if (n_1>=mini_reac):
					error_time=(abs((t-t_ref)/t_ref))*100
					error_temp=(abs((gas1.T-T_maxref)/gas1.T))*100
					ign_delay.append(t)
					num.append(i)
					T_max.append(gas1.T)
				
				if (error_time<=error_tol and error_temp<=error_tol and tol_check==1):
					print('Minumum no of reactions for corresponding state value is %d'%n_1)
					tol_check=0
				n_1=n_1+1	
			
			plt.subplot(2,1,1)
			plt.plot(num,ign_delay,label='Ignition delay(reduced mechanism)')
			plt.axhline(y=t_ref,color='black',linestyle='--',label='Reference ignition delay=%s s' % t_ref)
			plt.title('Input State Values-Temperature=%s K|Pressure=%d bar|phi=%s' % (temp,(pres/1e5),phi))
			plt.ylabel('Ignition Delay(s)')
			plt.legend()
			plt.grid()
			plt.subplot(2,1,2)
			plt.plot(num,T_max,color='red',label='maximum temperature(reduced mechanism)')
			plt.axhline(y=T_maxref,color='black',linestyle='--',label='Reference max temperature=%d K' % T_maxref)
			plt.title('Input State Values-Temperature=%s K|Pressure=%d bar|phi=%s' % (temp,(pres/1e5),phi))
			plt.xlabel('Number of reactions')
			plt.ylabel('Maximum Temperature(k)')
			plt.legend()
			plt.grid()
			plt.show()						

Results:

Plot for all possible combinations of State Values.

The below plots must shows the Ignition Delay and Max temperature on the Y-axis, while the number of reactions in the mechanism on the X axis.The minimum number of reactions is set to 40.

The plots are generated for all possible combinations of State Values for temperature, pressure and equivalence ratio.Total number of combination is 27.

Figure 1:

Figure 2:

Figure 3:

Figure 4:

Figure 5:

Figure 6:

Figure 7:

Figure 8:

Figure 9:

Figure 10:

Figure 11:

Figure 12:

Figure 13:

Figure 14:

Figure 15:

Figure 16:

Figure 17:

Figure 18:

Figure 19:

Figure 20:

Figure 21:

Figure 22:

Figure 23:

Figure 24:

Figure 25:

Figure 26:

Figure 27:

Conclusion:

1. It is observed that first 50-60 most sensitive reaction are enough to get the results with 5% tolerance in ignition delay and maximum temperature.For higher accuracy we need to compromise on tolerance value which will increase the number of reactions.

2. From the graphs it can be observed that with 180 reactions there is negligible error with the reference values.So adding more reactions above 180 won’t have any impact on results.

3. So reaction reduction mechanism is a very efficient way to discard the reactions which have negligible impact on solution.This technique reduces the computational time.

Related Articles

Leave a Reply

Your email address will not be published. Required fields are marked *