Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fix multi pdep arrhenius calculation #1341

Merged
merged 4 commits into from Jul 5, 2018
Merged

Fix multi pdep arrhenius calculation #1341

merged 4 commits into from Jul 5, 2018

Conversation

mjohnson541
Copy link
Contributor

This removes the requirement that PdepArrhenius objects in a MultiPdepArrhenius must have the same pressure lists by changing the getRateCoefficient function to evaluate each PdepArrhenius object independently. Also adds a database test that checks that all library reactions can be evaluated.

@mliu49
Copy link
Contributor

mliu49 commented Apr 6, 2018

Are these two calculation approaches actually the same? I was looking through git history yesterday to understand why we have MultiArrhenius and MultiPDepArrhenius, and discovered that for a brief period, we actually had a MultiKinetics class:

  1. 5821ad3 (7/25/11): MultiArrhenius was replaced by MultiKinetics, which simply sums together rates in the list, allowing any type of kinetics to be used.
  2. d5803cc (9/3/12): MultiArrhenius and MultiPDepArrhenius classes were added, with MultiPDepArrhenius having the fancier calculation method.
  3. d386cb2 (9/3/12): The old classes were removed.

Based on those commits, I don't think that this is an appropriate fix for this issue. At the very least, I think we need to look deeper into the difference between the calculation methods.

@mliu49 mliu49 added the Status: Blocked This PR should not be merged label Apr 6, 2018
@mjohnson541
Copy link
Contributor Author

Suppose they were different, the point of having a MultiPdepArrhenius is to combine multiple PdepArrhenius rates together so they can be the kinetics for a single reaction object. Why should the rate of a MultiPdepArrhenius object be different from the sum of the rates of the PdepArrhenius objects it contains?

I should also clarify that we have a test that checks that a particular MultiPdepArrhenius evaluates to particular values and both versions pass that test.

@codecov
Copy link

codecov bot commented Apr 6, 2018

Codecov Report

Merging #1341 into master will decrease coverage by <.01%.
The diff coverage is 0%.

Impacted file tree graph

@@            Coverage Diff            @@
##           master   #1341      +/-   ##
=========================================
- Coverage   41.31%   41.3%   -0.01%     
=========================================
  Files         178     178              
  Lines       30348   30351       +3     
  Branches     6140    6141       +1     
=========================================
  Hits        12537   12537              
- Misses      16948   16950       +2     
- Partials      863     864       +1
Impacted Files Coverage Δ
rmgpy/data/kinetics/library.py 43.45% <ø> (ø) ⬆️
rmgpy/reaction.py 0% <0%> (ø) ⬆️
rmgpy/data/kinetics/common.py 72.76% <0%> (-0.6%) ⬇️

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 32af08d...d01dc79. Read the comment docs.

@mliu49 mliu49 removed the Status: Blocked This PR should not be merged label Apr 6, 2018
@mliu49
Copy link
Contributor

mliu49 commented Apr 6, 2018

After looking at it some more, I'm now convinced that the previous method wasn't doing anything special. The only difference should be the check to ensure identical plists.

I think this PR can be merged to address ReactionMechanismGenerator/RMG-website#165. Following this change, there is effectively no difference between MultiArrhenius and MultiPDepArrhenius except for the additional toArrhenius method in MultiArrhenius. I propose that we recombine these two classes into MultiKinetics in the near future.

@rwest, I know this is from a long time ago, but do you remember any details regarding the reasoning for the current implementation of MultiArrhenius/MultiPDepArrhenius?

alongd
alongd previously approved these changes Apr 7, 2018
Copy link
Member

@alongd alongd left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The technical changes look good, I agree this should be the solution. Thanks @mjohnson541. If everything else is accepted on everyone, this could be merged.

@rwest
Copy link
Member

rwest commented Apr 7, 2018

I don't recall exactly, @mliu49, but I had enough of a deja-vu to make me search, leading me to this discussion which may be informative:
#147
#329

This test ensures that every library reaction has evaluable kinetics
"""
for entry in library.entries.values():
entry.data.getRateCoefficient(1000.0,1.0)
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Should we also add some additional checks here? For example, check that the rate is non-negative and within a reasonable range?

@mjohnson541
Copy link
Contributor Author

Do we want to proceed with this solution?

@mliu49
Copy link
Contributor

mliu49 commented May 8, 2018

Yes, I think we can proceed with this. Do you want to add additional checks like @cgrambow suggested?

@mjohnson541
Copy link
Contributor Author

mjohnson541 commented May 8, 2018

Yes, we can check that it's nonnegative pretty easily, @cgrambow were you thinking of checking the diffusion and TST limits?

@cgrambow
Copy link

cgrambow commented May 8, 2018

Yeah, that probably makes the most sense.

@mjohnson541
Copy link
Contributor Author

mjohnson541 commented May 29, 2018

So I just pushed positivity, TST and collision limit tests and we have a lot of failures in existing libraries:

======================================================================
ERROR: Kinetics library Sulfur/DMDS: check rates can be evaluated?
----------------------------------------------------------------------
Traceback (most recent call last):
  File "/usr/local/Cellar/numpy/1.12.0/libexec/nose/lib/python2.7/site-packages/nose/case.py", line 197, in runTest
    self.test(*self.arg)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 116, in <lambda>
    test = lambda x: self.kinetics_checkLibraryRatesCanBeEvaluated(library)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 448, in kinetics_checkLibraryRatesCanBeEvaluated
    raise ValueError('library reaction {0} from library {1}, exceeds the collision limit at 1000 K, 1 bar'.format(rxn,library.label))
ValueError: library reaction Sa + Sa <=> S2 from library Sulfur/DMDS, exceeds the collision limit at 1000 K, 1 bar

======================================================================
ERROR: Kinetics library Sulfur/GlarborgH2S: check rates can be evaluated?
----------------------------------------------------------------------
Traceback (most recent call last):
  File "/usr/local/Cellar/numpy/1.12.0/libexec/nose/lib/python2.7/site-packages/nose/case.py", line 197, in runTest
    self.test(*self.arg)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 116, in <lambda>
    test = lambda x: self.kinetics_checkLibraryRatesCanBeEvaluated(library)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 448, in kinetics_checkLibraryRatesCanBeEvaluated
    raise ValueError('library reaction {0} from library {1}, exceeds the collision limit at 1000 K, 1 bar'.format(rxn,library.label))
ValueError: library reaction HSS + H <=> SH + SH from library Sulfur/GlarborgH2S, exceeds the collision limit at 1000 K, 1 bar

======================================================================
ERROR: Kinetics library Sulfur/GlarborgBozzelli: check rates can be evaluated?
----------------------------------------------------------------------
Traceback (most recent call last):
  File "/usr/local/Cellar/numpy/1.12.0/libexec/nose/lib/python2.7/site-packages/nose/case.py", line 197, in runTest
    self.test(*self.arg)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 116, in <lambda>
    test = lambda x: self.kinetics_checkLibraryRatesCanBeEvaluated(library)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 448, in kinetics_checkLibraryRatesCanBeEvaluated
    raise ValueError('library reaction {0} from library {1}, exceeds the collision limit at 1000 K, 1 bar'.format(rxn,library.label))
ValueError: library reaction SO + H <=> HSO from library Sulfur/GlarborgBozzelli, exceeds the collision limit at 1000 K, 1 bar

======================================================================
ERROR: Kinetics library Dooley/methylformate_all_N2bathgas: check rates can be evaluated?
----------------------------------------------------------------------
Traceback (most recent call last):
  File "/usr/local/Cellar/numpy/1.12.0/libexec/nose/lib/python2.7/site-packages/nose/case.py", line 197, in runTest
    self.test(*self.arg)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 116, in <lambda>
    test = lambda x: self.kinetics_checkLibraryRatesCanBeEvaluated(library)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 445, in kinetics_checkLibraryRatesCanBeEvaluated
    raise ValueError('library reaction {0} from library {1}, exceeds the TST limit at 1000 K, 1 bar'.format(rxn,library.label))
ValueError: library reaction C3H6OOH2-2 <=> CH3COCH3 + OH from library Dooley/methylformate_all_N2bathgas, exceeds the TST limit at 1000 K, 1 bar

======================================================================
ERROR: Kinetics library Sulfur/GlarborgNS: check rates can be evaluated?
----------------------------------------------------------------------
Traceback (most recent call last):
  File "/usr/local/Cellar/numpy/1.12.0/libexec/nose/lib/python2.7/site-packages/nose/case.py", line 197, in runTest
    self.test(*self.arg)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 116, in <lambda>
    test = lambda x: self.kinetics_checkLibraryRatesCanBeEvaluated(library)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 448, in kinetics_checkLibraryRatesCanBeEvaluated
    raise ValueError('library reaction {0} from library {1}, exceeds the collision limit at 1000 K, 1 bar'.format(rxn,library.label))
ValueError: library reaction S + NO <=> SNO from library Sulfur/GlarborgNS, exceeds the collision limit at 1000 K, 1 bar

======================================================================
ERROR: Kinetics library Sulfur/DTBS: check rates can be evaluated?
----------------------------------------------------------------------
Traceback (most recent call last):
  File "/usr/local/Cellar/numpy/1.12.0/libexec/nose/lib/python2.7/site-packages/nose/case.py", line 197, in runTest
    self.test(*self.arg)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 116, in <lambda>
    test = lambda x: self.kinetics_checkLibraryRatesCanBeEvaluated(library)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 445, in kinetics_checkLibraryRatesCanBeEvaluated
    raise ValueError('library reaction {0} from library {1}, exceeds the TST limit at 1000 K, 1 bar'.format(rxn,library.label))
ValueError: library reaction S2 <=> S2JJ from library Sulfur/DTBS, exceeds the TST limit at 1000 K, 1 bar

======================================================================
ERROR: Kinetics library Chernov: check rates can be evaluated?
----------------------------------------------------------------------
Traceback (most recent call last):
  File "/usr/local/Cellar/numpy/1.12.0/libexec/nose/lib/python2.7/site-packages/nose/case.py", line 197, in runTest
    self.test(*self.arg)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 116, in <lambda>
    test = lambda x: self.kinetics_checkLibraryRatesCanBeEvaluated(library)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 448, in kinetics_checkLibraryRatesCanBeEvaluated
    raise ValueError('library reaction {0} from library {1}, exceeds the collision limit at 1000 K, 1 bar'.format(rxn,library.label))
ValueError: library reaction H2CCCH + C2H2 <=> C5H5 from library Chernov, exceeds the collision limit at 1000 K, 1 bar

======================================================================
ERROR: Kinetics library CurranPentane: check rates can be evaluated?
----------------------------------------------------------------------
Traceback (most recent call last):
  File "/usr/local/Cellar/numpy/1.12.0/libexec/nose/lib/python2.7/site-packages/nose/case.py", line 197, in runTest
    self.test(*self.arg)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 116, in <lambda>
    test = lambda x: self.kinetics_checkLibraryRatesCanBeEvaluated(library)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 445, in kinetics_checkLibraryRatesCanBeEvaluated
    raise ValueError('library reaction {0} from library {1}, exceeds the TST limit at 1000 K, 1 bar'.format(rxn,library.label))
ValueError: library reaction C4H71-1,2OOH <=> NC4KET12 + OH from library CurranPentane, exceeds the TST limit at 1000 K, 1 bar

======================================================================
ERROR: Kinetics library Nitrogen_Glarborg_Zhang_et_al: check rates can be evaluated?
----------------------------------------------------------------------
Traceback (most recent call last):
  File "/usr/local/Cellar/numpy/1.12.0/libexec/nose/lib/python2.7/site-packages/nose/case.py", line 197, in runTest
    self.test(*self.arg)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 116, in <lambda>
    test = lambda x: self.kinetics_checkLibraryRatesCanBeEvaluated(library)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 445, in kinetics_checkLibraryRatesCanBeEvaluated
    raise ValueError('library reaction {0} from library {1}, exceeds the TST limit at 1000 K, 1 bar'.format(rxn,library.label))
ValueError: library reaction CHCHNO <=> C2H2 + NO from library Nitrogen_Glarborg_Zhang_et_al, exceeds the TST limit at 1000 K, 1 bar

======================================================================
ERROR: Kinetics library Dooley/methylformate_all_ARHEbathgas: check rates can be evaluated?
----------------------------------------------------------------------
Traceback (most recent call last):
  File "/usr/local/Cellar/numpy/1.12.0/libexec/nose/lib/python2.7/site-packages/nose/case.py", line 197, in runTest
    self.test(*self.arg)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 116, in <lambda>
    test = lambda x: self.kinetics_checkLibraryRatesCanBeEvaluated(library)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 445, in kinetics_checkLibraryRatesCanBeEvaluated
    raise ValueError('library reaction {0} from library {1}, exceeds the TST limit at 1000 K, 1 bar'.format(rxn,library.label))
ValueError: library reaction C3H6OOH2-2 <=> CH3COCH3 + OH from library Dooley/methylformate_all_ARHEbathgas, exceeds the TST limit at 1000 K, 1 bar

======================================================================
ERROR: Kinetics library primaryNitrogenLibrary: check rates can be evaluated?
----------------------------------------------------------------------
Traceback (most recent call last):
  File "/usr/local/Cellar/numpy/1.12.0/libexec/nose/lib/python2.7/site-packages/nose/case.py", line 197, in runTest
    self.test(*self.arg)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 116, in <lambda>
    test = lambda x: self.kinetics_checkLibraryRatesCanBeEvaluated(library)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 448, in kinetics_checkLibraryRatesCanBeEvaluated
    raise ValueError('library reaction {0} from library {1}, exceeds the collision limit at 1000 K, 1 bar'.format(rxn,library.label))
ValueError: library reaction NCN + CN <=> NCNCN from library primaryNitrogenLibrary, exceeds the collision limit at 1000 K, 1 bar

======================================================================
ERROR: Kinetics library NOx2018: check rates can be evaluated?
----------------------------------------------------------------------
Traceback (most recent call last):
  File "/usr/local/Cellar/numpy/1.12.0/libexec/nose/lib/python2.7/site-packages/nose/case.py", line 197, in runTest
    self.test(*self.arg)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 116, in <lambda>
    test = lambda x: self.kinetics_checkLibraryRatesCanBeEvaluated(library)
  File "/Users/mattjohnson/RMGCODE/RMG-Py/testing/databaseTest.py", line 445, in kinetics_checkLibraryRatesCanBeEvaluated
    raise ValueError('library reaction {0} from library {1}, exceeds the TST limit at 1000 K, 1 bar'.format(rxn,library.label))
ValueError: library reaction CHCHNO <=> C2H2 + NO from library NOx2018, exceeds the TST limit at 1000 K, 1 bar

----------------------------------------------------------------------
Ran 842 tests in 210.412s

@mjohnson541
Copy link
Contributor Author

mjohnson541 commented May 29, 2018

I tested the rates at 1000 K 1 bar, the collision and TST limit calculations should be in the commit log.

@alongd
Copy link
Member

alongd commented Jun 1, 2018

One of the reactions this reports is:
root: ERROR: library reaction ONONO2 <=> NO2 + NO2 from library primaryNitrogenLibrary, exceeds the TST limit at 1000 K, 1 bar
This is barrierless scission of:
image
I double-checked with the source by M.C. Lin to make sure there are no typos.
Does it make sense that a high-P limit unimolecular decomposition reaction has a TST limit as the error suggests?

@mjohnson541 mjohnson541 self-assigned this Jun 7, 2018
@@ -443,14 +443,14 @@ def kinetics_checkLibraryRatesAreReasonable(self, library):
if k < 0:
boo = True
logging.error('library reaction {0} from library {1}, had a negative rate at 1000 K, 1 bar'.format(rxn,library.label))
if len(rxn.reactants) == 1:
if len(rxn.reactants) == 1 and rxn.allow_max_rate_violation==False:
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The test fails with AttributeError: 'rmgpy.reaction.Reaction' object has no attribute 'allow_max_rate_violation', I think you should first check hasattr

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Sorry, I forgot to actually set the attribute Reaction.init

@mjohnson541 mjohnson541 force-pushed the fixMultiPdepArrhenius branch 2 times, most recently from 2b44210 to ab2430a Compare June 11, 2018 19:35
@alongd
Copy link
Member

alongd commented Jun 12, 2018

@mjohnson541, I still get an error:
TypeError: loadEntry() got an unexpected keyword argument 'allow_max_rate_violation'
ReactionMechanismGenerator/RMG-database#279 fails when I run the tests locally on this -Py branch

@mjohnson541 mjohnson541 force-pushed the fixMultiPdepArrhenius branch 2 times, most recently from 5515206 to 9568307 Compare June 12, 2018 17:15
@mjohnson541 mjohnson541 force-pushed the fixMultiPdepArrhenius branch 2 times, most recently from cf4acda to 27b0e3e Compare June 13, 2018 19:15
@mliu49
Copy link
Contributor

mliu49 commented Jun 22, 2018

Any chance we could separate this from the collision limit check? It would be nice to fix the rate evaluation.

@mjohnson541
Copy link
Contributor Author

The collision limit check stuff is almost fixed I think we're just waiting for @jimchu10 to look at the aromatic ones.

mjohnson541 and others added 4 commits July 3, 2018 18:30
…f the PdepArrhenius objects within it

before there was a requirement that the pressure lists of each PdepArrhenius objects had to be the same
checks that a 1000 K 1 bar all rates are positive and less than the TST/collision limit
@mliu49
Copy link
Contributor

mliu49 commented Jul 4, 2018

I went ahead and rebased this. Any final considerations before merging?

@mjohnson541
Copy link
Contributor Author

mjohnson541 commented Jul 4, 2018 via email

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

5 participants