diff --git a/examples/rmg/methylformate/input.py b/examples/rmg/methylformate/input.py new file mode 100644 index 0000000000..5f4fb681ac --- /dev/null +++ b/examples/rmg/methylformate/input.py @@ -0,0 +1,144 @@ +# Data sources +database( + '../RMG-database/input', + thermoLibraries = ['primaryThermoLibrary','DFT_QCI_thermo','GRI-Mech3.0'], + reactionLibraries = ['Methylformate','Glarborg/highP'], + seedMechanisms = ['GRI-Mech3.0'], +) + +# List of species +species( + label='Mfmt', + reactive=True, + structure=SMILES("COC=O"), +) +species( + label='O2', + reactive=True, + structure=SMILES("[O][O]"), +) +species( + label='C2H', + reactive=True, + structure=SMILES("C#[C]"), +) +species( + label='C2H', + reactive=True, + structure=SMILES("C#[C]"), +) +species( + label='CH', + reactive=True, + structure=adjacencyList( + """ + 1 C 3 {2,S} + 2 H 0 {1,S} + """), +) +species( + label='H2O', + reactive=True, + structure=SMILES("O"), +) +species( + label='H2', + reactive=True, + structure=SMILES("[H][H]"), +) +species( + label='CO', + reactive=True, + structure=SMILES("[C]=O"), +) +species( + label='CO2', + reactive=True, + structure=SMILES("C(=O)=O"), +) +species( + label='CH4', + reactive=True, + structure=SMILES("C"), +) +species( + label='CH3', + reactive=True, + structure=SMILES("[CH3]"), +) +species( + label='CH3OH', + reactive=True, + structure=SMILES("CO"), +) +species( + label='C2H4', + reactive=True, + structure=SMILES("C=C"), +) +species( + label='C2H2', + reactive=True, + structure=SMILES("C#C"), +) +species( + label='CH2O', + reactive=True, + structure=SMILES("C=O"), +) +species( + label='CH3CHO', + reactive=True, + structure=SMILES("CC=O"), +) + + +# Bath gas +species( + label='Ar', + reactive=False, + structure=InChI("InChI=1S/Ar"), +) + +# Reaction systems +simpleReactor( + temperature=(1350,'K'), + pressure=(1.0,'bar'), + initialMoleFractions={ + "Mfmt": 0.01, + "O2": 0.02, + "Ar": 0.97, + }, +) + +termination( + time=(0.5,'s'), +) + +simulator( + atol=1e-22, + rtol=1e-8, +) + +model( + toleranceKeepInEdge=0.0, + toleranceMoveToCore=0.005, + toleranceInterruptSimulation=1.0, + maximumEdgeSpecies=100000 +) + +pressureDependence( + method='modified strong collision', # 'reservoir state' + minimumGrainSize=(2.0,'kcal/mol'), + minimumNumberOfGrains=200, + temperatures=(290,'K',3500,'K',8), + pressures=(0.02,'bar',100,'bar',5), + interpolation=('Chebyshev', 6, 4), +) + +options( + units='si', + saveRestart=False, + drawMolecules=False, + generatePlots=True, +) diff --git a/rmg.py b/rmg.py index fd8429a75b..86701c768e 100755 --- a/rmg.py +++ b/rmg.py @@ -251,6 +251,10 @@ def execute(args): # Read input file reactionModel, coreSpecies, reactionSystems, database, seedMechanisms = readInputFile(inputFile) + # Delete previous HTML file + from rmgpy.rmg.output import saveOutputHTML + saveOutputHTML(os.path.join(settings.outputDirectory, 'output.html'), reactionModel) + # Initialize reaction model if args.restart: reactionModel = loadRestartFile(os.path.join(settings.outputDirectory,'restart.pkl.gz')) @@ -335,12 +339,16 @@ def execute(args): # Save the current state of the model core to a pretty HTML file logging.info('Saving latest model core to HTML file...') - from rmgpy.rmg.output import saveOutputHTML saveOutputHTML(os.path.join(settings.outputDirectory, 'output.html'), reactionModel) # Save a Chemkin file containing the current model core logging.info('Saving latest model core to Chemkin file...') - reactionModel.saveChemkinFile(os.path.join(settings.outputDirectory, 'chemkin', 'chem%04i.inp' % len(reactionModel.core.species))) + this_chemkin_path = os.path.join(settings.outputDirectory, 'chemkin', 'chem%04i.inp' % len(reactionModel.core.species)) + latest_chemkin_path = os.path.join(settings.outputDirectory, 'chemkin','chem.inp') + reactionModel.saveChemkinFile(this_chemkin_path) + if os.path.exists(latest_chemkin_path): + os.unlink(latest_chemkin_path) + os.link(this_chemkin_path,latest_chemkin_path) # Save the restart file if desired if settings.saveRestart: diff --git a/rmgpy/data/thermo.py b/rmgpy/data/thermo.py index 4fb9647524..bf13f6fb21 100644 --- a/rmgpy/data/thermo.py +++ b/rmgpy/data/thermo.py @@ -582,7 +582,7 @@ def getThermoDataFromGroups(self, molecule): thermoData = None - if sum([atom.radicalElectrons for atom in molecule.atoms]) > 0: + if sum([atom.radicalElectrons for atom in molecule.atoms]) > 0: # radical species # Make a copy of the structure so we don't change the original saturatedStruct = molecule.copy(deep=True) @@ -646,7 +646,7 @@ def getThermoDataFromGroups(self, molecule): # Correct the entropy for the symmetry number - else: + else: # non-radical species # Generate estimate of thermodynamics for atom in molecule.atoms: # Iterate over heavy (non-hydrogen) atoms diff --git a/rmgpy/measure/output.py b/rmgpy/measure/output.py index 166e559d62..1b51c5f279 100644 --- a/rmgpy/measure/output.py +++ b/rmgpy/measure/output.py @@ -118,10 +118,10 @@ def writeNetworkConfigurations(f, network): """ # Write isomer configurations for isomer in network.isomers: - f.write('isomer("{0}")\n\n'.format(isomer.label)) + f.write('isomer("{0}")\n\n'.format(isomer)) # Write reactant configurations for reactants in network.reactants: - f.write('reactants("{0}", "{1}")\n\n'.format(reactants[0].label, reactants[1].label)) + f.write('reactants("{0}", "{1}")\n\n'.format(reactants[0], reactants[1])) # No need to write product configurations, as these are assumed ################################################################################ diff --git a/rmgpy/rmg/model.py b/rmgpy/rmg/model.py index c328c07b4d..5977e2540d 100755 --- a/rmgpy/rmg/model.py +++ b/rmgpy/rmg/model.py @@ -488,8 +488,8 @@ def update(self, reactionModel, database): self.collisionModel = SingleExponentialDownModel(alpha0=4.86 * 4184) # Save input file - rmgpy.measure.output.writeInput(os.path.join(settings.outputDirectory, 'pdep', 'network{0:d}_{1:d}.py'.format(self.index, len(self.isomers))), - self, Tlist, Plist, (grainSize, numGrains), method, model) + rmgpy.measure.output.writeFile(os.path.join(settings.outputDirectory, 'pdep', 'network{0:d}_{1:d}.py'.format(self.index, len(self.isomers))), + self, Tlist, Plist, (grainSize, numGrains), method, model, Tmin, Tmax, Pmin, Pmax) self.printSummary(level=logging.INFO) diff --git a/unittest/moleculeTest.py b/unittest/moleculeTest.py index 8ab51b33c4..936be75848 100644 --- a/unittest/moleculeTest.py +++ b/unittest/moleculeTest.py @@ -856,7 +856,22 @@ def testPickle(self): self.assertEqual(molecule0.getFormula(), molecule.getFormula()) self.assertTrue(molecule0.isIsomorphic(molecule)) self.assertTrue(molecule.isIsomorphic(molecule0)) - + + def testRadicalCH(self): + """ + Test that the species [CH] has three radical electrons and a spin multiplicity of 4. + """ + molecule = Molecule().fromSMILES('[CH]') + self.assertEqual(molecule.atoms[0].radicalElectrons, 3) + self.assertEqual(molecule.atoms[0].spinMultiplicity, 4) + + def testRadicalCH2(self): + """ + Test that the species [CH2] has two radical electrons and a spin multiplicity of 3. + """ + molecule = Molecule().fromSMILES('[CH2]') + self.assertEqual(molecule.atoms[0].radicalElectrons, 2) + self.assertEqual(molecule.atoms[0].spinMultiplicity, 3) ################################################################################ class TestMoleculeSymmetry(unittest.TestCase):