diff --git a/rmgpy/kinetics/arrhenius.pyx b/rmgpy/kinetics/arrhenius.pyx index 244318a61d..a02806dfb0 100644 --- a/rmgpy/kinetics/arrhenius.pyx +++ b/rmgpy/kinetics/arrhenius.pyx @@ -591,7 +591,7 @@ cdef class ArrheniusBM(KineticsModel): assert w0 is not None or recipe is not None, 'either w0 or recipe must be specified' if Ts is None: - Ts = [300.0, 500.0, 600.0, 700.0, 800.0, 900.0, 1000.0, 1100.0, 1200.0, 1500.0] + Ts = [300.0, 500.0, 600.0, 700.0, 800.0, 900.0, 1000.0, 1100.0, 1200.0, 1500.0, 2000.0] if w0 is None: #estimate w0 w0s = get_w0s(recipe, rxns) diff --git a/rmgpy/kinetics/arrheniusTest.py b/rmgpy/kinetics/arrheniusTest.py index 9888ef1bb2..10449215ad 100644 --- a/rmgpy/kinetics/arrheniusTest.py +++ b/rmgpy/kinetics/arrheniusTest.py @@ -37,7 +37,11 @@ import numpy as np import rmgpy.constants as constants -from rmgpy.kinetics.arrhenius import Arrhenius, ArrheniusEP, PDepArrhenius, MultiArrhenius, MultiPDepArrhenius +from rmgpy.kinetics.arrhenius import Arrhenius, ArrheniusEP, ArrheniusBM, PDepArrhenius, MultiArrhenius, MultiPDepArrhenius +from rmgpy.molecule.molecule import Molecule +from rmgpy.reaction import Reaction +from rmgpy.species import Species +from rmgpy.thermo import NASA, NASAPolynomial ################################################################################ @@ -404,6 +408,129 @@ def test_change_rate(self): self.assertAlmostEqual(2 * kexp, kact, delta=1e-6 * kexp) +################################################################################ + +class TestArrheniusBM(unittest.TestCase): + """ + Contains unit tests of the :class:`ArrheniusBM` class. + """ + + def setUp(self): + """ + A function run before each unit test in this class. + """ + self.A = 8.00037e+12 + self.n = 0.391734 + self.w0 = 798000 + self.E0 = 115905 + self.Tmin = 300. + self.Tmax = 2000. + self.comment = 'rxn001084' + self.arrhenius_bm = ArrheniusBM( + A=(self.A, "s^-1"), + n=self.n, + w0=(self.w0, 'J/mol'), + E0=(self.E0, "J/mol"), + Tmin=(self.Tmin, "K"), + Tmax=(self.Tmax, "K"), + comment=self.comment, + ) + + self.rsmi = 'NC(=NC=O)O' + self.psmi = 'O=CNC(=O)N' + self.arrhenius = Arrhenius(A=(8.00037e+12,'s^-1'), + n=0.391734, + Ea=(94.5149,'kJ/mol'), + T0=(1,'K'), + Tmin=(300,'K'), + Tmax=(2000,'K'), + comment="""Fitted to 50 data points; dA = *|/ 1.18377, dn = +|- 0.0223855, dEa = +|- 0.115431 kJ/mol""" + ) + + self.r_thermo = NASA(polynomials=[ + NASAPolynomial(coeffs=[3.90453,0.0068491,0.000125755,-2.92973e-07,2.12971e-10,-45444.2,10.0669], Tmin=(10,'K'), Tmax=(433.425,'K')), + NASAPolynomial(coeffs=[2.09778,0.0367646,-2.36023e-05,7.24527e-09,-8.51275e-13,-45412,15.8381], Tmin=(433.425,'K'), Tmax=(3000,'K'))], + Tmin=(10,'K'), Tmax=(3000,'K'), E0=(-377.851,'kJ/mol'), Cp0=(33.2579,'J/(mol*K)'), CpInf=(232.805,'J/(mol*K)'), + comment="""Thermo library: Spiekermann_refining_elementary_reactions""" + ) + self.p_thermo = NASA(polynomials=[ + NASAPolynomial(coeffs=[3.88423,0.00825528,0.000133399,-3.31802e-07,2.52823e-10,-51045.1,10.3937], Tmin=(10,'K'), Tmax=(428.701,'K')), + NASAPolynomial(coeffs=[2.89294,0.0351772,-2.26349e-05,7.00331e-09,-8.2982e-13,-51122.5,12.4424], Tmin=(428.701,'K'), Tmax=(3000,'K'))], + Tmin=(10,'K'), Tmax=(3000,'K'), E0=(-424.419,'kJ/mol'), Cp0=(33.2579,'J/(mol*K)'), CpInf=(232.805,'J/(mol*K)'), + comment="""Thermo library: Spiekermann_refining_elementary_reactions""" + ) + + def test_a_factor(self): + """ + Test that the ArrheniusBM A property was properly set. + """ + self.assertAlmostEqual(self.arrhenius_bm.A.value_si, self.A, delta=1e0) + + def test_n(self): + """ + Test that the ArrheniusBM n property was properly set. + """ + self.assertAlmostEqual(self.arrhenius_bm.n.value_si, self.n, 6) + + def test_w0(self): + """ + Test that the ArrheniusBM w0 property was properly set. + """ + self.assertAlmostEqual(self.arrhenius_bm.w0.value_si, self.w0, 6) + + def test_e0(self): + """ + Test that the ArrheniusBM E0 property was properly set. + """ + self.assertAlmostEqual(self.arrhenius_bm.E0.value_si, self.E0, 6) + + def test_temperature_min(self): + """ + Test that the ArrheniusBM Tmin property was properly set. + """ + self.assertAlmostEqual(self.arrhenius_bm.Tmin.value_si, self.Tmin, 6) + + def test_temperature_max(self): + """ + Test that the ArrheniusBM Tmax property was properly set. + """ + self.assertAlmostEqual(self.arrhenius_bm.Tmax.value_si, self.Tmax, 6) + + def test_is_temperature_valid(self): + """ + Test the ArrheniusBM.is_temperature_valid() method. + """ + Tdata = np.array([200, 400, 600, 800, 1000, 1200, 1400, 1600, 1800, 2000]) + validdata = np.array([False, True, True, True, True, True, True, True, True, True], np.bool) + for T, valid in zip(Tdata, validdata): + valid0 = self.arrhenius_bm.is_temperature_valid(T) + self.assertEqual(valid0, valid) + + def test_fit_to_data(self): + """ + Test the ArrheniusBM.fit_to_data() method. + """ + reactant = Molecule(smiles=self.rsmi) + product = Molecule(smiles=self.psmi) + reaction = Reaction(reactants=[Species(molecule=[reactant], thermo=self.r_thermo,)], + products=[Species(molecule=[product], thermo=self.p_thermo)], + kinetics=self.arrhenius, + ) + + arrhenius_bm = ArrheniusBM().fit_to_reactions([reaction], w0=self.w0) + self.assertAlmostEqual(arrhenius_bm.A.value_si, self.arrhenius_bm.A.value_si, delta=1.5e1) + self.assertAlmostEqual(arrhenius_bm.n.value_si, self.arrhenius_bm.n.value_si, 1, 4) + self.assertAlmostEqual(arrhenius_bm.E0.value_si, self.arrhenius_bm.E0.value_si, 1) + + def test_get_activation_energy(self): + """ + Test the ArrheniusBM.get_activation_energy() method. + """ + Hrxn = -44000 # J/mol + Ea = self.arrhenius_bm.get_activation_energy(Hrxn) + self.assertAlmostEqual(Ea, 95074, delta=1e1) + + ################################################################################ class TestPDepArrhenius(unittest.TestCase):