diff --git a/documentation/source/users/rmg/input.rst b/documentation/source/users/rmg/input.rst index f1fdfc5df8..50de5173a4 100644 --- a/documentation/source/users/rmg/input.rst +++ b/documentation/source/users/rmg/input.rst @@ -586,13 +586,15 @@ all of RMG's reaction families. :: generatedSpeciesConstraints( allowed=['input species','seed mechanisms','reaction libraries'], maximumCarbonAtoms=10, - maximumHydrogenAtoms=10, - maximumOxygenAtoms=10, - maximumNitrogenAtoms=10, - maximumSiliconAtoms=10, - maximumSulfurAtoms=10, + maximumOxygenAtoms=2, + maximumNitrogenAtoms=2, + maximumSiliconAtoms=2, + maximumSulfurAtoms=2, maximumHeavyAtoms=10, - maximumRadicalElectrons=10, + maximumRadicalElectrons=2, + maximumSingletCarbenes=1, + maximumCarbeneRadicals=0, + maximumIsotopicAtoms=2, allowSingletO2 = False, ) diff --git a/examples/rmg/commented/input.py b/examples/rmg/commented/input.py index 43ae42be21..b2260c21eb 100644 --- a/examples/rmg/commented/input.py +++ b/examples/rmg/commented/input.py @@ -1,42 +1,42 @@ # Data sources database( - #overrides RMG thermo calculation of RMG with these values. - #libraries found at http://rmg.mit.edu/database/thermo/libraries/ - #if species exist in multiple libraries, the earlier libraries overwrite the - #previous values - thermoLibraries = ['BurkeH2O2','primaryThermoLibrary','DFT_QCI_thermo','CBS_QB3_1dHR'], - #overrides RMG kinetics estimation if needed in the core of RMG. - #list of libraries found at http://rmg.mit.edu/database/kinetics/libraries/ - #libraries can be input as either a string or tuple of form ('library_name',True/False) - #where a `True` indicates that all unused reactions will be automatically added - #to the chemkin file at the end of the simulation. Placing just string values - #defaults the tuple to `False`. The string input is sufficient in almost - #all situations - reactionLibraries = [('C3', False)], - #seed mechanisms are reactionLibraries that are forced into the initial mechanism - #in addition to species listed in this input file. - #This is helpful for reducing run time for species you know will appear in - #the mechanism. - seedMechanisms = ['BurkeH2O2inN2','ERC-FoundationFuelv0.9'], - #this is normally not changed in general RMG runs. Usually used for testing with - #outside kinetics databases - kineticsDepositories = 'default', - #lists specific families used to generate the model. 'default' uses a list of - #families from RMG-Database/input/families/recommended.py - #a visual list of families is available in PDF form at RMG-database/families - kineticsFamilies = 'default', - #specifies how RMG calculates rates. currently, the only option is 'rate rules' - kineticsEstimator = 'rate rules', + # overrides RMG thermo calculation of RMG with these values. + # libraries found at http://rmg.mit.edu/database/thermo/libraries/ + # if species exist in multiple libraries, the earlier libraries overwrite the + # previous values + thermoLibraries=['BurkeH2O2', 'primaryThermoLibrary', 'DFT_QCI_thermo', 'CBS_QB3_1dHR'], + # overrides RMG kinetics estimation if needed in the core of RMG. + # list of libraries found at http://rmg.mit.edu/database/kinetics/libraries/ + # libraries can be input as either a string or tuple of form ('library_name',True/False) + # where a `True` indicates that all unused reactions will be automatically added + # to the chemkin file at the end of the simulation. Placing just string values + # defaults the tuple to `False`. The string input is sufficient in almost + # all situations + reactionLibraries=[('C3', False)], + # seed mechanisms are reactionLibraries that are forced into the initial mechanism + # in addition to species listed in this input file. + # This is helpful for reducing run time for species you know will appear in + # the mechanism. + seedMechanisms=['BurkeH2O2inN2', 'ERC-FoundationFuelv0.9'], + # this is normally not changed in general RMG runs. Usually used for testing with + # outside kinetics databases + kineticsDepositories='default', + # lists specific families used to generate the model. 'default' uses a list of + # families from RMG-Database/input/families/recommended.py + # a visual list of families is available in PDF form at RMG-database/families + kineticsFamilies='default', + # specifies how RMG calculates rates. currently, the only option is 'rate rules' + kineticsEstimator='rate rules', ) # List of species -#list initial and expected species below to automatically put them into the core mechanism. -#'structure' can utilize method of SMILES("put_SMILES_here"), -#adjacencyList("""put_adj_list_here"""), or InChI("put_InChI_here") -#for molecular oxygen, use the smiles string [O][O] so the triplet form is used +# list initial and expected species below to automatically put them into the core mechanism. +# 'structure' can utilize method of SMILES("put_SMILES_here"), +# adjacencyList("""put_adj_list_here"""), or InChI("put_InChI_here") +# for molecular oxygen, use the smiles string [O][O] so the triplet form is used species( label='butane', - reactive=True, #this parameter is optional if true + reactive=True, # this parameter is optional if true structure=SMILES("CCCC"), ) species( @@ -48,8 +48,8 @@ reactive=False, structure=adjacencyList(""" 1 N u0 p1 c0 {2,T} - 2 N u0 p1 c0 {1,T} - """), + 2 N u0 p1 c0 {1,T} + """), ) # You can list species not initially in reactor to make sure RMG includes them in the mechanism species( @@ -63,186 +63,191 @@ structure=SMILES("O=C=O") ) -#Reaction systems -#currently RMG models only constant temperature and pressure as homogeneous batch reactors. -#two options are: simpleReactor for gas phase or liquidReactor for liquid phase -#use can use multiple reactors in an input file for each condition you want to test. +# Reaction systems +# currently RMG models only constant temperature and pressure as homogeneous batch reactors. +# two options are: simpleReactor for gas phase or liquidReactor for liquid phase +# use can use multiple reactors in an input file for each condition you want to test. simpleReactor( - #specifies reaction temperature with units - temperature=(700,'K'), - #specifies reaction pressure with units - pressure=(10.0,'bar'), - #list initial mole fractions of compounds using the label from the 'species' label. - #RMG will normalize if sum/=1 + # specifies reaction temperature with units + temperature=(700, 'K'), + # specifies reaction pressure with units + pressure=(10.0, 'bar'), + # list initial mole fractions of compounds using the label from the 'species' label. + # RMG will normalize if sum/=1 initialMoleFractions={ "N2": 4, "O2": 1, - "butane": 1./6.5, + "butane": 1. / 6.5, }, - #the following two values specify when to determine the final output model - #only one must be specified - #the first condition to be satisfied will terminate the process + # the following two values specify when to determine the final output model + # only one must be specified + # the first condition to be satisfied will terminate the process terminationConversion={ 'butane': .99, }, - terminationTime=(40,'s'), - #the next two optional values specify how RMG computes sensitivities of - #rate coefficients with respect to species concentrations. - #sensitivity contains a list of species' labels to conduct sensitivity analysis on. - #sensitvityThreshold is the required sensitiviy to be recorded in the csv output file -# sensitivity=['CH4'], -# sensitivityThreshold=0.0001, + terminationTime=(40, 's'), + # the next two optional values specify how RMG computes sensitivities of + # rate coefficients with respect to species concentrations. + # sensitivity contains a list of species' labels to conduct sensitivity analysis on. + # sensitvityThreshold is the required sensitiviy to be recorded in the csv output file + # sensitivity=['CH4'], + # sensitivityThreshold=0.0001, ) # liquidReactor( -# temperature=(500,'K'), +# temperature=(500, 'K'), # initialConcentrations={ -# "N2": 4, -# "O2": 1, -# "CO": 1, +# "N2": 4, +# "O2": 1, +# "CO": 1, # }, # terminationConversion=None, -# terminationTime=(3600,'s'), +# terminationTime=(3600, 's'), # sensitivity=None, # sensitivityThreshold=1e-3 # ) - #liquid reactors also have solvents, you can specify one solvent - #list of solvents available at : http://rmg.mit.edu/database/solvation/libraries/solvent/ +# liquid reactors also have solvents, you can specify one solvent +# list of solvents available at : http://rmg.mit.edu/database/solvation/libraries/solvent/ # solvation('water') -#determines absolute and relative tolerances for ODE solver and sensitivities. -#normally this doesn't cause many issues and is modified after other issues are -#ruled out +# determines absolute and relative tolerances for ODE solver and sensitivities. +# normally this doesn't cause many issues and is modified after other issues are +# ruled out simulator( atol=1e-16, rtol=1e-8, -# sens_atol=1e-6, -# sens_rtol=1e-4, + # sens_atol=1e-6, + # sens_rtol=1e-4, ) -#used to add species to the model and to reduce memory usage by removing unimportant additional species. -#all relative values are normalized by a characteristic flux at that time point +# used to add species to the model and to reduce memory usage by removing unimportant additional species. +# all relative values are normalized by a characteristic flux at that time point model( - #determines the relative flux to put a species into the core. - #A smaller value will result in a larger, more complex model - #when running a new model, it is recommended to start with higher values and then decrease to converge on the model + # determines the relative flux to put a species into the core. + # A smaller value will result in a larger, more complex model + # when running a new model, it is recommended to start with higher values and then decrease to converge on the model toleranceMoveToCore=0.1, - #comment out the next three terms to disable pruning - #determines the relative flux needed to not remove species from the model. - #Lower values will keep more species and utilize more memory + # comment out the next three terms to disable pruning + # determines the relative flux needed to not remove species from the model. + # Lower values will keep more species and utilize more memory toleranceKeepInEdge=0.01, - #determines when to stop a ODE run to add a species. - #Lower values will improve speed. - #if it is too low, may never get to the end simulation to prune species. + # determines when to stop a ODE run to add a species. + # Lower values will improve speed. + # if it is too low, may never get to the end simulation to prune species. toleranceInterruptSimulation=1, - #number of edge species needed to accumulate before pruning occurs - #larger values require more memory and will prune less often + # number of edge species needed to accumulate before pruning occurs + # larger values require more memory and will prune less often maximumEdgeSpecies=100000, - #minimum number of core species needed before pruning occurs. - #this prevents pruning when kinetic model is far away from completeness + # minimum number of core species needed before pruning occurs. + # this prevents pruning when kinetic model is far away from completeness minCoreSizeForPrune=50, - #make sure that the pruned edge species have existed for a set number of RMG iterations. - #the user can specify to increase it from the default value of 2 + # make sure that the pruned edge species have existed for a set number of RMG iterations. + # the user can specify to increase it from the default value of 2 minSpeciesExistIterationsForPrune=2, - #filter the reactions during the enlarge step to omit species from reacting if their - #concentration are deemed to be too low + # filter the reactions during the enlarge step to omit species from reacting if their + # concentration are deemed to be too low filterReactions=False, ) options( - #provides a name for the seed mechanism produced at the end of an rmg run default is 'Seed' - name='SeedName', - #if True every iteration it saves the current model as libraries/seeds - #(and deletes the old one) - #Unlike HTML this is inexpensive time-wise - #note a seed mechanism will be generated at the end of a completed run and some incomplete - #runs even if this is set as False + # provides a name for the seed mechanism produced at the end of an rmg run default is 'Seed' + name='SeedName', + # if True every iteration it saves the current model as libraries/seeds + # (and deletes the old one) + # Unlike HTML this is inexpensive time-wise + # note a seed mechanism will be generated at the end of a completed run and some incomplete + # runs even if this is set as False generateSeedEachIteration=True, - #If True the mechanism will also be saved directly as kinetics and thermo libraries in the database + # If True the mechanism will also be saved directly as kinetics and thermo libraries in the database saveSeedToDatabase=False, - #only option is 'si' + # only option is 'si' units='si', - #how often you want to save restart files. - #takes significant amount of time. comment out if you don't want to save + # how often you want to save restart files. + # takes significant amount of time. comment out if you don't want to save saveRestartPeriod=None, - #Draws images of species and reactions and saves the model output to HTML. - #May consume extra memory when running large models. + # Draws images of species and reactions and saves the model output to HTML. + # May consume extra memory when running large models. generateOutputHTML=True, - #generates plots of the RMG's performance statistics. Not helpful if you just want a model. + # generates plots of the RMG's performance statistics. Not helpful if you just want a model. generatePlots=False, - #saves mole fraction of species in 'solver/' to help you create plots + # saves mole fraction of species in 'solver/' to help you create plots saveSimulationProfiles=False, - #gets RMG to output comments on where kinetics were obtained in the chemkin file. - #useful for debugging kinetics but increases memory usage of the chemkin output file + # gets RMG to output comments on where kinetics were obtained in the chemkin file. + # useful for debugging kinetics but increases memory usage of the chemkin output file verboseComments=False, - #gets RMG to generate edge species chemkin files. Uses lots of memory in output. - #Helpful for seeing why some reaction are not appearing in core model. + # gets RMG to generate edge species chemkin files. Uses lots of memory in output. + # Helpful for seeing why some reaction are not appearing in core model. saveEdgeSpecies=False, - #Sets a time limit in the form DD:HH:MM:SS after which the RMG job will stop. Useful for profiling on jobs that - #do not converge. - #wallTime = '00:00:00', + # Sets a time limit in the form DD:HH:MM:SS after which the RMG job will stop. Useful for profiling on jobs that + # do not converge. + # wallTime = '00:00:00', + # Forces RMG to import library reactions as reversible (default). Otherwise, if set to True, RMG will import library + # reactions while keeping the reversibility as as. keepIrreversible=False, - #Forces RMG to import library reactions as reversible (default). Otherwise, if set to True, RMG will import library - #reactions while keeping the reversibility as as. ) # optional module allows for correction to unimolecular reaction rates at low pressures and/or temperatures. pressureDependence( - #two methods available: 'modified strong collision' is faster and less accurate than 'reservoir state' - method='modified strong collision', - #these two categories determine how fine energy is descretized. - #more grains increases accuracy but takes longer - maximumGrainSize=(0.5,'kcal/mol'), - minimumNumberOfGrains=250, - #the conditions for the rate to be output over - #parameter order is: low_value, high_value, units, internal points - temperatures=(300,2200,'K',2), - pressures=(0.01,100,'bar',3), - #The two options for interpolation are 'PDepArrhenius' (no extra arguments) and - #'Chebyshev' which is followed by the number of basis sets in - #Temperature and Pressure. These values must be less than the number of - #internal points specified above - interpolation=('Chebyshev', 6, 4), - #turns off pressure dependence for molecules with number of atoms greater than the number specified below - #this is due to faster internal rate of energy transfer for larger molecules - maximumAtoms=15, - ) + # two methods available: 'modified strong collision' is faster and less accurate than 'reservoir state' + method='modified strong collision', + # these two categories determine how fine energy is descretized. + # more grains increases accuracy but takes longer + maximumGrainSize=(0.5, 'kcal/mol'), + minimumNumberOfGrains=250, + # the conditions for the rate to be output over + # parameter order is: low_value, high_value, units, internal points + temperatures=(300, 2200, 'K', 2), + pressures=(0.01, 100, 'bar', 3), + # The two options for interpolation are 'PDepArrhenius' (no extra arguments) and + # 'Chebyshev' which is followed by the number of basis sets in + # Temperature and Pressure. These values must be less than the number of + # internal points specified above + interpolation=('Chebyshev', 6, 4), + # turns off pressure dependence for molecules with number of atoms greater than the number specified below + # this is due to faster internal rate of energy transfer for larger molecules + maximumAtoms=15, +) -#optional block adds constraints on what RMG can output. -#This is helpful for improving the efficiency of RMG, but wrong inputs can lead to many errors. +# optional block adds constraints on what RMG can output. +# This is helpful for improving the efficiency of RMG, but wrong inputs can lead to many errors. generatedSpeciesConstraints( - #allows exceptions to the following restrictions - allowed=['input species','seed mechanisms','reaction libraries'], - #maximum number of each atom in a molecule + # allows exceptions to the following restrictions + allowed=['input species', 'seed mechanisms', 'reaction libraries'], + # maximum number of each atom in a molecule maximumCarbonAtoms=4, maximumOxygenAtoms=7, maximumNitrogenAtoms=0, maximumSiliconAtoms=0, maximumSulfurAtoms=0, - #max number of non-hydrogen atoms - #maximumHeavyAtoms=20, - #maximum radicals on a molecule + # max number of non-hydrogen atoms + # maximumHeavyAtoms=20, + # maximum radicals on a molecule maximumRadicalElectrons=1, - #If this is false or missing, RMG will throw an error if the more less-stable form of O2 is entered - #which doesn't react in the RMG system. normally input O2 as triplet with SMILES [O][O] - #allowSingletO2=False, + # maximum number of singlet carbenes (lone pair on a carbon atom) in a molecule + maximumSingletCarbenes=1, + # maximum number of radicals on a molecule with a singlet carbene + # should be lower than maximumRadicalElectrons in order to have an effect + maximumCarbeneRadicals=0, + # If this is false or missing, RMG will throw an error if the more less-stable form of O2 is entered + # which doesn't react in the RMG system. normally input O2 as triplet with SMILES [O][O] + # allowSingletO2=False, # maximum allowed number of non-normal isotope atoms: - #maximumIsotopicAtoms=2, + # maximumIsotopicAtoms=2, ) -#optional block allows thermo to be estimated through quantum calculations +# optional block allows thermo to be estimated through quantum calculations # quantumMechanics( -# #the software package for calculations...can use 'mopac' or 'gaussian' if installed -# software='mopac', -# #methods available for calculations. 'pm2' 'pm3' or 'pm7' (last for mopac only) -# method='pm3', -# #where to store calculations -# fileStore='QMfiles', -# #where to store temporary run files -# scratchDirectory = None, -# #onlyCyclics allows linear molecules to be calculated using bensen group addivity....need to verify -# onlyCyclics = True, -# #how many radicals should be utilized in the calculation. -# #If the amount of radicals is more than this, RMG will use hydrogen bond incrementation method -# maxRadicalNumber = 0, +# # the software package for calculations...can use 'mopac' or 'gaussian' if installed +# software='mopac', +# # methods available for calculations. 'pm2' 'pm3' or 'pm7' (last for mopac only) +# method='pm3', +# # where to store calculations +# fileStore='QMfiles', +# # where to store temporary run files +# scratchDirectory=None, +# # onlyCyclics allows linear molecules to be calculated using bensen group addivity....need to verify +# onlyCyclics=True, +# # how many radicals should be utilized in the calculation. +# # If the amount of radicals is more than this, RMG will use hydrogen bond incrementation method +# maxRadicalNumber=0, # ) diff --git a/rmgpy/constraints.py b/rmgpy/constraints.py index f31eb82f80..7fd53ca0b6 100644 --- a/rmgpy/constraints.py +++ b/rmgpy/constraints.py @@ -93,6 +93,16 @@ def failsSpeciesConstraints(species): if (struct.getRadicalCount() > maxRadicals): return True + maxCarbenes = speciesConstraints.get('maximumSingletCarbenes', 1) + if maxRadicals != -1: + if struct.getSingletCarbeneCount() > maxCarbenes: + return True + + maxCarbeneRadicals = speciesConstraints.get('maximumCarbeneRadicals', 0) + if maxCarbeneRadicals != -1: + if struct.getSingletCarbeneCount() > 0 and struct.getRadicalCount() > maxCarbeneRadicals: + return True + maxIsotopes = speciesConstraints.get('maximumIsotopicAtoms', -1) if maxIsotopes != -1: counter = 0 diff --git a/rmgpy/constraintsTest.py b/rmgpy/constraintsTest.py index 829cf8be0c..367427177e 100644 --- a/rmgpy/constraintsTest.py +++ b/rmgpy/constraintsTest.py @@ -33,8 +33,13 @@ """ import unittest +import mock -from rmgpy.constraints import * +from rmgpy.rmg.main import RMG +from rmgpy.constraints import failsSpeciesConstraints +from rmgpy.species import Species +from rmgpy.molecule import Molecule +import rmgpy.rmg.input ################################################################################ @@ -42,14 +47,196 @@ class TestFailsSpeciesConstraints(unittest.TestCase): """ Contains unit tests of the failsSpeciesConstraints function. """ - - def testConstraintsNotLoaded(self): + + @classmethod + def setUpClass(cls): + """ + A function run ONCE before all unit tests in this class. + """ + cls.rmg = RMG() + rmgpy.rmg.input.rmg = cls.rmg + rmgpy.rmg.input.generatedSpeciesConstraints( + maximumCarbonAtoms=2, + maximumOxygenAtoms=1, + maximumNitrogenAtoms=1, + maximumSiliconAtoms=1, + maximumSulfurAtoms=1, + maximumHeavyAtoms=3, + maximumRadicalElectrons=2, + maximumSingletCarbenes=1, + maximumCarbeneRadicals=0, + maximumIsotopicAtoms=2, + ) + + @classmethod + def tearDownClass(cls): + """ + A function run ONCE after all unit tests in this class. + """ + rmgpy.rmg.input.rmg = None + + @mock.patch('rmgpy.constraints.logging') + def testConstraintsNotLoaded(self, mock_logging): """ Test what happens when constraints are not loaded. """ - from rmgpy.species import Species + # Reset module level rmg variable in rmgpy.rmg.input + rmgpy.rmg.input.rmg = None + + mol = Molecule(SMILES='C') + self.assertFalse(failsSpeciesConstraints(mol)) + + mock_logging.debug.assert_called_with('Species constraints could not be found.') + + # Restore module level rmg variable in rmgpy.rmg.input + rmgpy.rmg.input.rmg = self.rmg + + def testSpeciesInput(self): + """ + Test that failsSpeciesConstraints can handle a Species object. + """ spc = Species().fromSMILES('C') - fails = failsSpeciesConstraints(spc) - self.assertFalse(fails) + self.assertFalse(failsSpeciesConstraints(spc)) + + def testExplicitlyAllowedMolecules(self): + """ + Test that we can explicitly allow molecules in species constraints. + """ + mol = Molecule(SMILES='CCCC') + self.assertTrue(failsSpeciesConstraints(mol)) + + self.rmg.speciesConstraints['explicitlyAllowedMolecules'] = [Molecule(SMILES='CCCC')] + self.assertFalse(failsSpeciesConstraints(mol)) + + def testCarbonConstraint(self): + """ + Test that we can constrain the max number of carbon atoms. + """ + mol1 = Molecule(SMILES='CC') + self.assertFalse(failsSpeciesConstraints(mol1)) + + mol2 = Molecule(SMILES='CCC') + self.assertTrue(failsSpeciesConstraints(mol2)) + + def testOxygenConstraint(self): + """ + Test that we can constrain the max number of oxygen atoms. + """ + mol1 = Molecule(SMILES='C=O') + self.assertFalse(failsSpeciesConstraints(mol1)) + + mol2 = Molecule(SMILES='OC=O') + self.assertTrue(failsSpeciesConstraints(mol2)) + + def testNitrogenConstraint(self): + """ + Test that we can constrain the max number of nitrogen atoms. + """ + mol1 = Molecule(SMILES='CN') + self.assertFalse(failsSpeciesConstraints(mol1)) + + mol2 = Molecule(SMILES='NCN') + self.assertTrue(failsSpeciesConstraints(mol2)) + + def testSiliconConstraint(self): + """ + Test that we can constrain the max number of silicon atoms. + """ + mol1 = Molecule(SMILES='[SiH4]') + self.assertFalse(failsSpeciesConstraints(mol1)) + + mol2 = Molecule(SMILES='[SiH3][SiH3]') + self.assertTrue(failsSpeciesConstraints(mol2)) + + def testSulfurConstraint(self): + """ + Test that we can constrain the max number of sulfur atoms. + """ + mol1 = Molecule(SMILES='CS') + self.assertFalse(failsSpeciesConstraints(mol1)) + + mol2 = Molecule(SMILES='SCS') + self.assertTrue(failsSpeciesConstraints(mol2)) + + def testHeavyConstraint(self): + """ + Test that we can constrain the max number of heavy atoms. + """ + mol1 = Molecule(SMILES='CCO') + self.assertFalse(failsSpeciesConstraints(mol1)) + + mol2 = Molecule(SMILES='CCN=O') + self.assertTrue(failsSpeciesConstraints(mol2)) + + def testRadicalConstraint(self): + """ + Test that we can constrain the max number of radical electrons. + """ + mol1 = Molecule(SMILES='[CH2][CH2]') + self.assertFalse(failsSpeciesConstraints(mol1)) + + mol2 = Molecule(SMILES='[CH2][CH][CH2]') + self.assertTrue(failsSpeciesConstraints(mol2)) + + def testCarbeneConstraint(self): + """ + Test that we can constrain the max number of singlet carbenes. + """ + mol1 = Molecule().fromAdjacencyList(""" +1 C u0 p1 c0 {2,S} {3,S} +2 H u0 p0 c0 {1,S} +3 H u0 p0 c0 {1,S} +""") + self.assertFalse(failsSpeciesConstraints(mol1)) + + mol2 = Molecule().fromAdjacencyList(""" +1 C u0 p1 c0 {2,S} {3,S} +2 H u0 p0 c0 {1,S} +3 C u0 p1 c0 {1,S} {4,S} +4 H u0 p0 c0 {3,S} +""") + self.assertTrue(failsSpeciesConstraints(mol2)) + + def testCarbeneRadicalConstraint(self): + """ + Test that we can constrain the max number of radical electrons with a carbene. + """ + mol1 = Molecule().fromAdjacencyList(""" +1 C u0 p1 c0 {2,S} {3,S} +2 H u0 p0 c0 {1,S} +3 H u0 p0 c0 {1,S} +""") + self.assertFalse(failsSpeciesConstraints(mol1)) + + mol2 = Molecule().fromAdjacencyList(""" +1 C u0 p1 c0 {2,S} {3,S} +2 H u0 p0 c0 {1,S} +3 C u1 p0 c0 {1,S} {4,S} {5,S} +4 H u0 p0 c0 {3,S} +5 H u0 p0 c0 {3,S} +""") + self.assertTrue(failsSpeciesConstraints(mol2)) + + def testIsotopeConstraint(self): + """ + Test that we can constrain the max number of isotopic atoms. + """ + mol1 = Molecule().fromAdjacencyList(""" +1 C u0 p0 c0 {2,S} {3,S} {4,S} {5,S} +2 D u0 p0 c0 {1,S} +3 D u0 p0 c0 {1,S} +4 H u0 p0 c0 {1,S} +5 H u0 p0 c0 {1,S} +""") + self.assertFalse(failsSpeciesConstraints(mol1)) + + mol2 = Molecule().fromAdjacencyList(""" +1 C u0 p0 c0 {2,S} {3,S} {4,S} {5,S} +2 D u0 p0 c0 {1,S} +3 D u0 p0 c0 {1,S} +4 D u0 p0 c0 {1,S} +5 H u0 p0 c0 {1,S} +""") + self.assertTrue(failsSpeciesConstraints(mol2)) diff --git a/rmgpy/data/kinetics/family.py b/rmgpy/data/kinetics/family.py index 4731c66940..17b380f727 100644 --- a/rmgpy/data/kinetics/family.py +++ b/rmgpy/data/kinetics/family.py @@ -1329,8 +1329,11 @@ def isMoleculeForbidden(self, molecule): forbidden_structures = getDB('forbidden') + # check family-specific forbidden structures if self.forbidden is not None and self.forbidden.isMoleculeForbidden(molecule): return True + + # check RMG globally forbidden structures if forbidden_structures.isMoleculeForbidden(molecule): return True return False diff --git a/rmgpy/molecule/molecule.pxd b/rmgpy/molecule/molecule.pxd index a307cd0fa5..c719a960b4 100644 --- a/rmgpy/molecule/molecule.pxd +++ b/rmgpy/molecule/molecule.pxd @@ -145,6 +145,8 @@ cdef class Molecule(Graph): cpdef short getRadicalCount(self) + cpdef short getSingletCarbeneCount(self) + cpdef double getMolecularWeight(self) cpdef int getNumAtoms(self, str element=?) diff --git a/rmgpy/molecule/molecule.py b/rmgpy/molecule/molecule.py index 66ca66b680..4ebdead8c6 100644 --- a/rmgpy/molecule/molecule.py +++ b/rmgpy/molecule/molecule.py @@ -884,6 +884,19 @@ def getRadicalCount(self): radicals += atom.radicalElectrons return radicals + def getSingletCarbeneCount(self): + """ + Return the total number of singlet carbenes (lone pair on a carbon atom) + in the molecule. Counts the number of carbon atoms with a lone pair. + In the case of [C] with two lone pairs, this method will return 1. + """ + cython.declare(atom=Atom, carbenes=cython.short) + carbenes = 0 + for atom in self.vertices: + if atom.isCarbon() and atom.lonePairs > 0: + carbenes += 1 + return carbenes + def getNumAtoms(self, element = None): """ Return the number of atoms in molecule. If element is given, ie. "H" or "C", diff --git a/rmgpy/molecule/moleculeTest.py b/rmgpy/molecule/moleculeTest.py index 3ce824c19e..24fbd1e65b 100644 --- a/rmgpy/molecule/moleculeTest.py +++ b/rmgpy/molecule/moleculeTest.py @@ -1089,6 +1089,29 @@ def testRadicalCH2CH2CH2(self): """ molecule = Molecule().fromSMILES('[CH2]C[CH2]') self.assertEqual(molecule.getRadicalCount(), 2) + + def testSingletCarbene(self): + """Test radical and carbene count on singlet carbene.""" + mol = Molecule().fromAdjacencyList(""" +1 C u0 p1 {2,S} +2 C u0 p1 {1,S} +""", saturateH=True) + self.assertEqual(mol.getRadicalCount(), 0) + self.assertEqual(mol.getSingletCarbeneCount(), 2) + + def testTripletCarbene(self): + """Test radical and carbene count on triplet carbene.""" + mol = Molecule().fromAdjacencyList(""" +1 C u2 p0 {2,S} +2 C u0 p1 {1,S} +""", saturateH=True) + self.assertEqual(mol.getRadicalCount(), 2) + self.assertEqual(mol.getSingletCarbeneCount(), 1) + + def testSingletCarbon(self): + """Test that getSingletCarbeneCount returns 1 for singlet carbon atom.""" + mol = Molecule().fromAdjacencyList('1 C u0 p2') + self.assertEqual(mol.getSingletCarbeneCount(), 1) def testSMILES(self): """ diff --git a/rmgpy/rmg/input.py b/rmgpy/rmg/input.py index 18b64210d0..ca9ce1fc77 100644 --- a/rmgpy/rmg/input.py +++ b/rmgpy/rmg/input.py @@ -353,6 +353,8 @@ def generatedSpeciesConstraints(**kwargs): 'maximumSulfurAtoms', 'maximumHeavyAtoms', 'maximumRadicalElectrons', + 'maximumSingletCarbenes', + 'maximumCarbeneRadicals', 'allowSingletO2', 'maximumIsotopicAtoms' ]