This script aims to extract model sources in a clear and informative format. The script first shows what all the kinetic and thermo sources are in a model. Then it goes through each reaction and species to show their source and what the assigned uncertainties are. This can be used with any RMG-generated CHEMKIN file that is annotated.
In [1]:
from rmgpy.tools.uncertainty import Uncertainty
from IPython.display import display
import copy
In [2]:
chemFile = '/home/mjliu/Documents/DTBS/new_runs/pyro/chem_annotated.inp'
dictFile = '/home/mjliu/Documents/DTBS/new_runs/pyro/species_dictionary.txt'
In [3]:
uncertainty = Uncertainty(outputDirectory='testUncertainty')
In [4]:
uncertainty.loadModel(chemFile, dictFile)
In [5]:
uncertainty.loadDatabase()
In [6]:
uncertainty.extractSourcesFromModel()
In [7]:
print 'All Kinetic Sources'
for sourceType in uncertainty.allKineticSources.keys():
if sourceType == 'Library':
print '============'
print 'Library kinetics'
print ''
print '\tReactions: ', uncertainty.allKineticSources['Library']
elif sourceType == 'PDep':
print '============'
print 'PDep kinetics'
print ''
print '\tReactions: ', uncertainty.allKineticSources['PDep']
elif sourceType == 'Rate Rules':
print '============'
print 'Rate rule kinetics'
print ''
for familyLabel, entries in uncertainty.allKineticSources['Rate Rules'].iteritems():
print '\t', familyLabel
for entry in entries:
print '\t\t', entry
elif sourceType == 'Training':
print '============'
print 'Training reaction kinetics'
print ''
for familyLabel, entries in uncertainty.allKineticSources['Training'].iteritems():
print '\t', familyLabel
for entry in entries:
print '\t\t', entry
else:
print sourceType
raise Exception('Kinetics source mut be Library, PDep, Rate Rules, or Training')
In [8]:
print 'All Thermo Sources'
for sourceType in uncertainty.allThermoSources.keys():
if sourceType == 'Library':
print '============'
print 'Library thermo'
print ''
print '\tSpecies: ', uncertainty.allThermoSources['Library']
elif sourceType == 'QM':
print '============'
print 'QM thermo'
print ''
print '\tSpecies: ', uncertainty.allThermoSources['QM']
elif sourceType == 'GAV':
print '============'
print 'Group additivity thermo'
print ''
for groupType, entries in uncertainty.allThermoSources['GAV'].iteritems():
print '\t', groupType
for entry in entries:
print '\t\t', entry
else:
raise Exception('Thermo source must be GAV, QM, or Library')
In [9]:
# Assign all the uncertainties using the Uncertainty() class function
# ThermoParameterUncertainty and KineticParameterUncertainty classes may be customized and passed into this function
# if non-default constants for constructing the uncertainties are desired
uncertainty.assignParameterUncertainties()
In [10]:
T = 623 # temperature in Kelvin for which to evaluate kinetics
P = 1e5 # Pa
In [11]:
for rxn, source in uncertainty.reactionSourcesDict.iteritems():
print '======'
print rxn
display(rxn)
if 'Library' in source:
print 'Library reaction'
elif 'PDep' in source:
print 'PDep reaction'
elif 'Rate Rules' in source:
print 'Rate rule estimate'
family = source['Rate Rules'][0]
sourceDict = source['Rate Rules'][1]
originalTemplate = sourceDict['template']
print '\tFamily = ', family
print '\tOriginal Template = ', [group.label for group in originalTemplate]
print '\tExact = ', sourceDict['exact']
rules = sourceDict['rules']
training = sourceDict['training']
if rules:
print '\tRate rule sources:'
for ruleEntry, weight in rules:
print '\t\t', ruleEntry, '=', weight
if training:
print '\tTraining sources:'
for ruleEntry, trainingEntry, weight in training:
print '\t\t', ruleEntry , 'mapped to', trainingEntry , '=', weight
elif 'Training' in source:
print 'Training reaction'
family = source['Training'][0]
training = source['Training'][1]
print '\t Family = ', family
print '\t\t', trainingEntry, '=', weight
print ''
print 'Rate coefficient at {} K = {:.2e}'.format(T, rxn.kinetics.getRateCoefficient(T,P))
# Uncomment the following lines if you want to verify that the parsing has been performed correctly by
# checking the values for both the original and reconstructed kinetics
# print '---------'
# print 'Original kinetics:'
# print rxn.kinetics
# print ''
# print 'Reconstructed kinetics from parsing:'
# reconstructedKinetics=uncertainty.database.kinetics.reconstructKineticsFromSource(rxn,source,fixBarrierHeight=True)
# print reconstructedKinetics
rxnIndex = uncertainty.reactionList.index(rxn)
print 'Uncertainty dln(k) = ', uncertainty.kineticInputUncertainties[rxnIndex]
In [12]:
for species, source in uncertainty.speciesSourcesDict.iteritems():
print '=========='
print species
display(species)
if 'Library' in source:
print 'Thermo Library: ', source['Library']
if 'QM' in source:
print 'QM: ', source['QM']
if 'GAV' in source:
print 'Group additivity:'
for groupType, groupList in source['GAV'].iteritems():
print '\t', groupType
for group, weight in groupList:
print '\t\t', group, '=', weight
spcIndex = uncertainty.speciesList.index(species)
print ''
print 'Uncertainty dG = ', uncertainty.thermoInputUncertainties[spcIndex], ' kcal/mol'