Skip to content

Reject zero rate coefficients and improve collision limit checks - #2992

Open
calvinp0 wants to merge 9 commits into
mainfrom
nan_fix
Open

Reject zero rate coefficients and improve collision limit checks#2992
calvinp0 wants to merge 9 commits into
mainfrom
nan_fix

Conversation

@calvinp0

@calvinp0 calvinp0 commented Jul 30, 2026

Copy link
Copy Markdown
Member

Motivation or Problem

An RMG job I ran on an ammonia oxidation model crashed during the final model check:

Making seed mechanism...
Performing final model checks...
/home/supervisor/Code/RMG-Py/rmgpy/rmg/main.py:1488: RuntimeWarning: divide by zero encountered in log
  violator_list = rxn.check_collision_limit_violation(t_min=self.Tmin, t_max=self.Tmax, p_min=self.Pmin, p_max=self.Pmax)
Traceback (most recent call last):
  File "/home/supervisor/Code/RMG-Py/rmg.py", line 4, in <module>
    __main__.main()
  File "/home/supervisor/Code/RMG-Py/rmgpy/__main__.py", line 102, in main
    rmg.execute(**kwargs)
  File "/home/supervisor/Code/RMG-Py/rmgpy/rmg/main.py", line 1244, in execute
    self.check_model()
  File "/home/supervisor/Code/RMG-Py/rmgpy/rmg/main.py", line 1488, in check_model
    violator_list = rxn.check_collision_limit_violation(t_min=self.Tmin, t_max=self.Tmax, p_min=self.Pmin, p_max=self.Pmax)
  File "rmgpy/reaction.py", line 1704, in rmgpy.reaction.Reaction.check_collision_limit_violation
  File "rmgpy/reaction.py", line 1737, in rmgpy.reaction.Reaction.check_collision_limit_violation
  File "rmgpy/kinetics/arrhenius.pyx", line 129, in rmgpy.kinetics.arrhenius.Arrhenius.get_rate_coefficient
  File "rmgpy/kinetics/arrhenius.pyx", line 139, in rmgpy.kinetics.arrhenius.Arrhenius.get_rate_coefficient
TypeError: Cannot convert 'complex' with non-zero imaginary component to 'double' (this most likely comes from the '**' operator; use 'cython.cpow(True)' to return 'nan' instead of a complex number).

It appears the issue chain was:

  1. generate_reverse_rate_coefficient() buils a Tlist/klist and fits an Arrhenius to it
  2. One forward rate coefficient in that list underflows to exactly 0.0 at the low end of the range
  3. Arrhenius.it_to_data validates its input against inf, NaN and mixed signs but not against zero. Zero passes every check and reaches the least-squares step, which is performed in log space, hence divide by zero encountered in log
  4. lstsq returns NaN for A, n and Ea. It does not raise. The degenerate expression is returned to the caller as though the fit succeeded.
  5. Evaluating that expression raises the TypeError above, far from the cause

Description of Changes

  1. Reject zero rate coefficients when fitting Arrhenius in log space
  2. Keep collisioon limit checks alive when a rate cannot be evaluated. check_collision_limit_violation now evaluates each direction independently and skips a condition it cannot evaluate, with a warning, instead of propogating out of check_model and discarding an otherwise-complete run.
  3. Derive the reverse rate derivative from Keq instead of kb/kf. compute_rate_derivative formed the reverse contribution as (kb/kf)*[products]. Every reactor builds kb as kf/Keq so when kf underflows to exactly zero that ratio is 0.0/0.0, which in C evaluates to NaN without raising and spreads through the sensitivity matrix. Since kb≡kf/Keq, the ratio is identically 1/Keq so the stored Keq is used directly. This is the same value whenever kf > 0 and stays correct when kf underflows, where the true limist is [products]/Keq rather than zero.
  4. Checks the reverse collision limit when reactants lack transport data. calculate_col_limit raises ValueError when a species in the requested direction has no transport data; in the forward direction this was handled with continue, but then skipps the reverse direction check for that condition even though it depends only on the products
  5. Fix argument documentation and tables markup in canteramodel - authored by @rwest

Testing

7 unit tests.

· test/rmgpy/kinetics/arrheniusTest.py::TestArrhenius::test_fit_to_data_with_zero_rate
· test/rmgpy/reactionTest.py::TestCollisionLimitViolation - 4 tests
· test/rmgpy/solver/simpleTest.py::SimpleReactorTest::test_compute_rate_derivative_kf_underflow
· test/rmgpy/rmg/mainTest.py::TestCheckModelCollisionLimits::test_failing_reaction_does_not_abort_the_check

calvinp0 and others added 5 commits July 30, 2026 13:58
Arrhenius.fit_to_data and ArrheniusChargeTransfer.fit_to_data validated
their input rate coefficients against inf and NaN, and against mixed
signs, but not against zero. Zero passes both checks and reaches the
least-squares step, which is performed in log space: log(0) is -inf, so
numpy emits "divide by zero encountered in log" and lstsq then returns
NaN for every fitted parameter without raising.

The NaN expression is returned to the caller as if the fit had
succeeded. It surfaces much later and far from the cause -- evaluating
it raises "TypeError: Cannot convert 'complex' with non-zero imaginary
component to 'double'" from the '**' operator, which reads like a
negative temperature rather than a degenerate fit. It is also possible
for such an expression to reach the Chemkin output instead, in which
case the run completes and writes NaN kinetics with no error at all.

A rate coefficient of exactly zero typically comes from underflow at the
low end of a fitted temperature range, which is common when generating a
reverse rate coefficient for a strongly endothermic reaction. Reject it
at the same point inf and NaN are rejected, so the failure is reported
where it originates.
check_collision_limit_violation evaluated a rate coefficient for every
extreme T/P condition and compared it against the collision limit for
that condition. Any reaction whose rate could not be evaluated -- an
unsupported kinetics type raising ReactionError from
generate_reverse_rate_coefficient, or a degenerate reverse fit -- raised
out of check_model and aborted the run at the very end, after the model
had otherwise completed.

Evaluate each direction independently and record the rate, its collision
limit and the condition together, so that a condition which cannot be
evaluated is skipped with a warning instead of propagating. Storing the
three values as one tuple rather than in parallel lists also removes the
possibility of the comparison loop pairing a rate with another
condition's limit: the collision limit was previously appended before
the rate was evaluated, so a rate that failed left an orphan limit
behind and shifted every later comparison. Since the collision limit
grows with sqrt(T), that skew compared a high-temperature rate against a
smaller low-temperature limit and reported violations that do not exist.

In check_model, keep the broad catch as a backstop -- nothing detected
by this end-of-run diagnostic justifies discarding a completed run --
but log the traceback so a genuine bug behind the failure stays
diagnosable rather than being reduced to a one-line message.
compute_rate_derivative formed the reverse contribution to d(dy/dt)/dk
as (kb/kf) * [products]. Every reactor builds kb as kf/Keq, so when a
large activation energy drives kf to exactly zero, kb is zero too and
that ratio is the indeterminate 0.0/0.0. In C this does not raise: it
evaluates to NaN, which then spreads through flux, whole columns of
rate_deriv and the sensitivity right-hand side. The function only runs
under sensitivity analysis, so the integration itself was unaffected,
but any job with sensitivity enabled could produce a NaN sensitivity
matrix with no indication of where it came from.

Since kb is always kf/Keq, the ratio kb/kf is identically 1/Keq, so
divide by the stored Keq directly. This gives the same value whenever
kf > 0 and remains correct when kf underflows, where the true limit is
[products]/Keq rather than zero -- for a strongly reverse-dominated
reaction that value is large, so treating it as zero misreported the
reaction as forward-only.

Guard with `not (Keq[j] > 0)` rather than a test against zero. Keq is
inf for an irreversible reaction under SimpleReactor and zero under the
liquid and surface reactors, so both sentinels must be handled, and
writing the test this way also rejects a NaN Keq arising from bad
thermochemistry.
calculate_coll_limit raises ValueError when a species involved in the
requested direction has no transport data. In the forward direction that
was handled with `continue`, which skips the rest of the loop body for
that condition and therefore also skips the reverse-direction check --
even though the reverse check depends only on the products, which may
have perfectly good transport data.

Fall through instead, so each direction is evaluated on its own merits.
A reaction whose reactants lack transport data is now still checked in
the reverse direction rather than being silently exempted from the
collision limit report.
The docstrings for generate_cantera_conditions and
Cantera.generate_conditions documented arguments named T0List, P0List
and V0List, which do not exist; the parameters are Tlist, Plist and
Vlist. They also used a mixture of backticks and single quotes for the
same argument names, and listed the supported reactor types as
unmarked continuation lines.

Correct the names, mark the reactor types up as a list, and give both
tables the closing border row that reStructuredText simple tables
require -- without it docutils rejects the table with "Malformed table.
No bottom table border found." Separate the second table from the
sentence above it with a blank line, since a table has to begin a new
block to be recognised at all.
Copilot AI review requested due to automatic review settings July 30, 2026 11:19

Copilot AI left a comment

Copy link
Copy Markdown

Choose a reason for hiding this comment

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

Copilot was unable to review this pull request because the user who requested the review has reached their quota limit.

@calvinp0
calvinp0 requested review from alongd and rwest July 30, 2026 11:20
@github-actions

Copy link
Copy Markdown

Regression Testing Results

⚠️ One or more regression tests failed.
Please download the failed results and run the tests locally or check the log to see why.

Detailed regression test results.

Regression test aromatics:

Reference: Execution time (DD:HH:MM:SS): 00:00:00:53
Current: Execution time (DD:HH:MM:SS): 00:00:00:58
Reference: Memory used: 830.35 MB
Current: Memory used: 829.80 MB

aromatics Passed Core Comparison ✅

Original model has 15 species.
Test model has 15 species. ✅
Original model has 11 reactions.
Test model has 11 reactions. ✅

aromatics Failed Edge Comparison ❌

Original model has 106 species.
Test model has 106 species. ✅
Original model has 358 reactions.
Test model has 358 reactions. ✅

Non-identical thermo! ❌
original: [CH]1C2=CC3C1C=CC23
tested: [CH]1C2=CC3C1C=CC23

Hf(300K) S(300K) Cp(300K) Cp(400K) Cp(500K) Cp(600K) Cp(800K) Cp(1000K) Cp(1500K)
167.21 73.60 28.78 36.79 44.00 50.25 59.65 65.52 74.04
169.15 73.17 31.27 38.45 44.76 50.28 59.14 65.47 72.92

thermo: Thermo group additivity estimation: group(Cs-(Cds-Cds)CsCsH) + group(Cs-(Cds-Cds)CsCsH) + group(Cs-(Cds-Cds)(Cds-Cds)CsH) + group(Cs-(Cds-Cds)CsHH) + group(Cds-CdsCsCs) + group(Cds-CdsCsH) + group(Cds-CdsCsH) + group(Cds-CdsCsH) + polycyclic(s2_4_5_diene_1_5) + polycyclic(s3_4_5_ene_3) + polycyclic(s2_5_5_diene_1_5) - ring(Cyclobutene) - ring(Cyclopentene) - ring(Cyclopentene) + radical(cyclopentene-allyl)
thermo: Thermo group additivity estimation: group(Cs-(Cds-Cds)CsCsH) + group(Cs-(Cds-Cds)CsCsH) + group(Cs-(Cds-Cds)(Cds-Cds)CsH) + group(Cs-(Cds-Cds)CsHH) + group(Cds-CdsCsCs) + group(Cds-CdsCsH) + group(Cds-CdsCsH) + group(Cds-CdsCsH) + polycyclic(s2_4_5_diene_1_5) + polycyclic(s3_4_5_ene_3) + polycyclic(s3_5_5_ene_1) - ring(Cyclobutene) - ring(Cyclopentene) - ring(Cyclopentane) + radical(cyclopentene-allyl)

Non-identical thermo! ❌
original: [CH]1C2C=CC13C=CC23
tested: [CH]1C2C=CC13C=CC23

Hf(300K) S(300K) Cp(300K) Cp(400K) Cp(500K) Cp(600K) Cp(800K) Cp(1000K) Cp(1500K)
172.50 80.61 27.59 35.91 43.44 49.89 59.09 64.21 71.65
174.31 74.05 26.66 34.03 40.90 47.10 57.16 64.03 72.57

thermo: Thermo group additivity estimation: group(Cs-(Cds-Cds)(Cds-Cds)CsCs) + group(Cs-(Cds-Cds)CsCsH) + group(Cs-(Cds-Cds)CsCsH) + group(Cs-CsCsHH) + group(Cds-CdsCsH) + group(Cds-CdsCsH) + group(Cds-CdsCsH) + group(Cds-CdsCsH) + polycyclic(s2_4_4_ene_1) + polycyclic(s2_4_5_diene_1_5) + polycyclic(s3_4_5_ene_1) - ring(Cyclobutene) - ring(Cyclobutane) - ring(Cyclopentene) + radical(bicyclo[2.1.1]hex-2-ene-C5)
thermo: Thermo group additivity estimation: group(Cs-(Cds-Cds)(Cds-Cds)CsCs) + group(Cs-(Cds-Cds)CsCsH) + group(Cs-(Cds-Cds)CsCsH) + group(Cs-CsCsHH) + group(Cds-CdsCsH) + group(Cds-CdsCsH) + group(Cds-CdsCsH) + group(Cds-CdsCsH) + polycyclic(s2_4_4_ene_1) + polycyclic(s1_4_5_diene_1_6) + polycyclic(s3_4_5_ene_1) - ring(Cyclobutene) - ring(Cyclobutane) - ring(Cyclopentene) + radical(bicyclo[2.1.1]hex-2-ene-C5)

Non-identical thermo! ❌
original: [CH]1C2=CC3C1C3C=C2
tested: [CH]1C2=CC3C1C3C=C2

Hf(300K) S(300K) Cp(300K) Cp(400K) Cp(500K) Cp(600K) Cp(800K) Cp(1000K) Cp(1500K)
100.48 61.70 25.50 33.41 40.70 47.02 56.22 61.78 71.32
98.15 66.21 25.82 33.30 40.19 46.24 55.47 61.34 70.49

thermo: Thermo group additivity estimation: group(Cs-CsCsCsH) + group(Cs-(Cds-Cds)CsCsH) + group(Cs-(Cds-Cds)CsCsH) + group(Cs-(Cds-Cds)CsHH) + group(Cds- Cds(Cds-Cds)Cs) + group(Cds-CdsCsH) + group(Cds-CdsCsH) + group(Cds-Cds(Cds-Cds)H) + polycyclic(s2_3_5_ene_1) + polycyclic(s2_3_6_ene_1) + Estimated bicyclic component: polycyclic(s3_5_6_ane) - ring(Cyclohexane) - ring(Cyclopentane) + ring(Cyclohexene) + ring(Cyclopentene) - ring(Cyclopropane) - ring(Cyclopentene) - ring(Cyclohexene) + radical(cyclopentene-allyl)
thermo: Thermo group additivity estimation: group(Cs-CsCsCsH) + group(Cs-(Cds-Cds)CsCsH) + group(Cs-(Cds-Cds)CsCsH) + group(Cs-(Cds-Cds)CsHH) + group(Cds- Cds(Cds-Cds)Cs) + group(Cds-CdsCsH) + group(Cds-CdsCsH) + group(Cds-Cds(Cds-Cds)H) + polycyclic(s2_3_5_ene_1) + polycyclic(s2_3_6_diene_1_3) + Estimated bicyclic component: polycyclic(s3_5_6_ane) - ring(Cyclohexane) - ring(Cyclopentane) + ring(1,3-Cyclohexadiene) + ring(Cyclopentene) - ring(Cyclopropane) - ring(Cyclopentene) - ring(1,3-Cyclohexadiene) + radical(cyclopentene-allyl)

Non-identical thermo! ❌
original: C=CC1C=CC2=CC1C=C2
tested: C=CC1C=CC2=CC1C=C2

Hf(300K) S(300K) Cp(300K) Cp(400K) Cp(500K) Cp(600K) Cp(800K) Cp(1000K) Cp(1500K)
83.22 84.16 35.48 45.14 53.78 61.40 73.58 82.20 95.08
83.22 82.78 35.48 45.14 53.78 61.40 73.58 82.20 95.08

Identical thermo comments:
thermo: Thermo group additivity estimation: group(Cs-(Cds-Cds)(Cds-Cds)CsH) + group(Cs-(Cds-Cds)(Cds-Cds)CsH) + group(Cds-Cds(Cds-Cds)(Cds-Cds)) + group(Cds- CdsCsH) + group(Cds-CdsCsH) + group(Cds-CdsCsH) + group(Cds-CdsCsH) + group(Cds-Cds(Cds-Cds)H) + group(Cds-Cds(Cds-Cds)H) + group(Cds-CdsHH) + Estimated bicyclic component: polycyclic(s3_5_6_ane) - ring(Cyclohexane) - ring(Cyclopentane) + ring(1,3-Cyclohexadiene) + ring(Cyclopentadiene)

Non-identical thermo! ❌
original: [CH]1C2C=CC3=CC2C13
tested: [CH]1C2C=CC3=CC2C13

Hf(300K) S(300K) Cp(300K) Cp(400K) Cp(500K) Cp(600K) Cp(800K) Cp(1000K) Cp(1500K)
200.28 76.10 25.91 33.35 40.28 46.48 56.24 62.60 71.29
144.26 70.81 25.26 32.45 39.06 44.87 53.78 59.67 69.60

thermo: Thermo group additivity estimation: group(Cs-(Cds-Cds)CsCsH) + group(Cs-(Cds-Cds)CsCsH) + group(Cs-(Cds-Cds)CsCsH) + group(Cs-CsCsHH) + group(Cds- Cds(Cds-Cds)Cs) + group(Cds-CdsCsH) + group(Cds-CdsCsH) + group(Cds-Cds(Cds-Cds)H) + polycyclic(s2_4_4_ene_1) + polycyclic(s3_4_6_diene_1_5) + polycyclic(s3_4_6_ene_1) - ring(Cyclobutene) - ring(Cyclobutane) - ring(Cyclohexene) + radical(cyclobutane)
thermo: Thermo group additivity estimation: group(Cs-(Cds-Cds)CsCsH) + group(Cs-(Cds-Cds)CsCsH) + group(Cs-(Cds-Cds)CsCsH) + group(Cs-CsCsHH) + group(Cds- Cds(Cds-Cds)Cs) + group(Cds-CdsCsH) + group(Cds-CdsCsH) + group(Cds-Cds(Cds-Cds)H) + polycyclic(s2_4_4_ene_1) + polycyclic(s3_4_6_ene_1) + Estimated bicyclic component: polycyclic(s2_4_6_ane) - ring(Cyclohexane) - ring(Cyclobutane) + ring(Cyclohexene) + ring(Cyclobutene) - ring(Cyclobutane) - ring(Cyclobutene) - ring(Cyclohexene) + radical(cyclobutane)

Non-identical kinetics! ❌
original:
rxn: [CH]1C2=CC=CC1C=C2(48) <=> [CH]1C2=CC3C1C=CC23(62) origin: Intra_R_Add_Endocyclic
tested:
rxn: [CH]1C2=CC=CC1C=C2(48) <=> [CH]1C2=CC3C1C=CC23(62) origin: Intra_R_Add_Endocyclic

k(1bar) 300K 400K 500K 600K 800K 1000K 1500K 2000K
k(T): -46.27 -30.58 -21.19 -14.94 -7.15 -2.49 3.67 6.72
k(T): -47.51 -31.51 -21.94 -15.56 -7.62 -2.87 3.42 6.54

kinetics: Arrhenius(A=(1.08454e+19,'s^-1'), n=-0.859, Ea=(86.724,'kcal/mol'), T0=(1,'K'), comment="""Estimated from node Backbone1_2R!H-inRing_1R!H-inRing in family Intra_R_Add_Endocyclic.""")
kinetics: Arrhenius(A=(1.08454e+19,'s^-1'), n=-0.859, Ea=(88.43,'kcal/mol'), T0=(1,'K'), comment="""Estimated from node Backbone1_2R!H-inRing_1R!H-inRing in family Intra_R_Add_Endocyclic.""")
Identical kinetics comments:
kinetics: Estimated from node Backbone1_2R!H-inRing_1R!H-inRing in family Intra_R_Add_Endocyclic.

Non-identical kinetics! ❌
original:
rxn: [CH]1C2=CC=CC1C=C2(48) <=> [CH]1C2C=CC13C=CC23(65) origin: Intra_R_Add_Endocyclic
tested:
rxn: [CH]1C2=CC=CC1C=C2(48) <=> [CH]1C2C=CC13C=CC23(65) origin: Intra_R_Add_Endocyclic

k(1bar) 300K 400K 500K 600K 800K 1000K 1500K 2000K
k(T): -49.69 -33.15 -23.24 -16.65 -8.43 -3.52 2.99 6.21
k(T): -50.88 -34.04 -23.95 -17.24 -8.88 -3.88 2.75 6.03

kinetics: Arrhenius(A=(1.08454e+19,'s^-1'), n=-0.859, Ea=(91.423,'kcal/mol'), T0=(1,'K'), comment="""Estimated from node Backbone1_2R!H-inRing_1R!H-inRing in family Intra_R_Add_Endocyclic.""")
kinetics: Arrhenius(A=(1.08454e+19,'s^-1'), n=-0.859, Ea=(93.051,'kcal/mol'), T0=(1,'K'), comment="""Estimated from node Backbone1_2R!H-inRing_1R!H-inRing in family Intra_R_Add_Endocyclic.""")
Identical kinetics comments:
kinetics: Estimated from node Backbone1_2R!H-inRing_1R!H-inRing in family Intra_R_Add_Endocyclic.

Non-identical kinetics! ❌
original:
rxn: [CH]1C2=CC=CC1C=C2(48) <=> [CH]1C2=CC3C1C3C=C2(67) origin: Intra_R_Add_Endocyclic
tested:
rxn: [CH]1C2=CC=CC1C=C2(48) <=> [CH]1C2=CC3C1C3C=C2(67) origin: Intra_R_Add_Endocyclic

k(1bar) 300K 400K 500K 600K 800K 1000K 1500K 2000K
k(T): -15.17 -8.42 -4.36 -1.66 1.73 3.77 6.50 7.88
k(T): -14.18 -7.68 -3.77 -1.16 2.10 4.07 6.70 8.03

kinetics: Arrhenius(A=(1.12e+11,'s^-1'), n=0.26, Ea=(36.869,'kcal/mol'), T0=(1,'K'), comment="""Estimated from node Backbone3_Sp-4R!H=1R!H_Sp-3R!H-2R!H_Sp-2R!H-1R!H_Ext-2R!H-R in family Intra_R_Add_Endocyclic.""")
kinetics: Arrhenius(A=(1.12e+11,'s^-1'), n=0.26, Ea=(35.513,'kcal/mol'), T0=(1,'K'), comment="""Estimated from node Backbone3_Sp-4R!H=1R!H_Sp-3R!H-2R!H_Sp-2R!H-1R!H_Ext-2R!H-R in family Intra_R_Add_Endocyclic.""")
Identical kinetics comments:
kinetics: Estimated from node Backbone3_Sp-4R!H=1R!H_Sp-3R!H-2R!H_Sp-2R!H-1R!H_Ext-2R!H-R in family Intra_R_Add_Endocyclic.

Non-identical kinetics! ❌
original:
rxn: C1=CC2C=C[C]1C=C2(49) <=> [CH]1C2C=CC3=CC2C13(80) origin: Intra_R_Add_Endocyclic
tested:
rxn: C1=CC2C=C[C]1C=C2(49) <=> [CH]1C2C=CC3=CC2C13(80) origin: Intra_R_Add_Endocyclic

k(1bar) 300K 400K 500K 600K 800K 1000K 1500K 2000K
k(T): -59.93 -40.64 -29.08 -21.38 -11.79 -6.05 1.56 5.34
k(T): -25.05 -14.47 -8.15 -3.94 1.29 4.42 8.54 10.57

kinetics: Arrhenius(A=(6.50724e+19,'s^-1'), n=-0.859, Ea=(106.547,'kcal/mol'), T0=(1,'K'), comment="""Estimated from node Backbone1_2R!H-inRing_1R!H-inRing in family Intra_R_Add_Endocyclic. Multiplied by reaction path degeneracy 6.0""")
kinetics: Arrhenius(A=(6.50724e+19,'s^-1'), n=-0.859, Ea=(58.664,'kcal/mol'), T0=(1,'K'), comment="""Estimated from node Backbone1_2R!H-inRing_1R!H-inRing in family Intra_R_Add_Endocyclic. Multiplied by reaction path degeneracy 6.0""")
Identical kinetics comments:
kinetics: Estimated from node Backbone1_2R!H-inRing_1R!H-inRing in family Intra_R_Add_Endocyclic.
Multiplied by reaction path degeneracy 6.0

Non-identical kinetics! ❌
original:
rxn: C1=CC2C=C[C]1C=C2(49) <=> [CH]1C2=CC3C1C3C=C2(67) origin: Intra_R_Add_Endocyclic
tested:
rxn: C1=CC2C=C[C]1C=C2(49) <=> [CH]1C2=CC3C1C3C=C2(67) origin: Intra_R_Add_Endocyclic

k(1bar) 300K 400K 500K 600K 800K 1000K 1500K 2000K
k(T): -8.89 -3.16 0.28 2.58 5.46 7.19 9.52 10.69
k(T): -8.00 -2.50 0.81 3.02 5.79 7.46 9.70 10.83

kinetics: Arrhenius(A=(1.49409e+13,'s^-1'), n=0.283, Ea=(31.249,'kcal/mol'), T0=(1,'K'), comment="""Estimated from node Backbone2_Sp-3R!H=1R!H_N-4R!H->S_2R!H-inRing_5R!H-inRing_Ext-5R!H-R_Ext-6R!H-R_Ext-7R!H-R_1R!H-inRing in family Intra_R_Add_Endocyclic. Multiplied by reaction path degeneracy 3.0""")
kinetics: Arrhenius(A=(1.49409e+13,'s^-1'), n=0.283, Ea=(30.033,'kcal/mol'), T0=(1,'K'), comment="""Estimated from node Backbone2_Sp-3R!H=1R!H_N-4R!H->S_2R!H-inRing_5R!H-inRing_Ext-5R!H-R_Ext-6R!H-R_Ext-7R!H-R_1R!H-inRing in family Intra_R_Add_Endocyclic. Multiplied by reaction path degeneracy 3.0""")
Identical kinetics comments:
kinetics: Estimated from node Backbone2_Sp-3R!H=1R!H_N-4R!H->S_2R!H-inRing_5R!H-inRing_Ext-5R!H-R_Ext-6R!H-R_Ext-7R!H-R_1R!H-inRing in family Intra_R_Add_Endocyclic.
Multiplied by reaction path degeneracy 3.0

Errors occurred during edge comparison ⚠️ ERROR conda.cli.main_run:execute(148): `conda run python scripts/checkModels.py aromatics-edge stable_regression_results/aromatics/chemkin/chem_edge_annotated.inp stable_regression_results/aromatics/chemkin/species_edge_dictionary.txt test/regression/aromatics/chemkin/chem_edge_annotated.inp test/regression/aromatics/chemkin/species_edge_dictionary.txt` failed. (See above for error)
✅ All Observables varied by less than 0.500 on average between old model and new model in all conditions! aromatics Passed Observable Testing ✅

Regression test liquid_oxidation:

Reference: Execution time (DD:HH:MM:SS): 00:00:01:53
Current: Execution time (DD:HH:MM:SS): 00:00:02:00
Reference: Memory used: 912.38 MB
Current: Memory used: 915.01 MB

liquid_oxidation Passed Core Comparison ✅

Original model has 37 species.
Test model has 37 species. ✅
Original model has 239 reactions.
Test model has 239 reactions. ✅

liquid_oxidation Failed Edge Comparison ❌

Original model has 214 species.
Test model has 214 species. ✅
Original model has 1591 reactions.
Test model has 1588 reactions. ❌
The original model has 7 reactions that the tested model does not have. ❌
rxn: C[CH]CC(C)OO(31) <=> [OH](22) + CCCC(C)=O(28) origin: intra_H_migration
rxn: C[CH]CCCOO(63) <=> [OH](22) + CCCCC=O(60) origin: intra_H_migration
rxn: CC(C[CH]COO)OO(118) <=> CC(CC[CH]OO)OO(133) origin: intra_H_migration
rxn: CC(C[CH]COO)OO(118) <=> C[C](CCCOO)OO(132) origin: intra_H_migration
rxn: CC(CC(C)OO)O[O](90) + CC(CCCOO)O[O](108) <=> oxygen(1) + CC([O])CC(C)OO(110) + CC([O])CCCOO(123) origin: Peroxyl_Disproportionation
rxn: CC(CC(C)OO)O[O](90) + CC(CCCOO)O[O](108) <=> oxygen(1) + CC(=O)CC(C)OO(95) + CC(O)CCCOO(152) origin: Peroxyl_Termination
rxn: CC(CC(C)OO)O[O](90) + CC(CCCOO)O[O](108) <=> oxygen(1) + CC(=O)CCCOO(115) + CC(O)CC(C)OO(143) origin: Peroxyl_Termination
The tested model has 4 reactions that the original model does not have. ❌
rxn: C[CH]CC(C)OO(33) <=> CCC[C](C)OO(78) origin: intra_H_migration
rxn: CCCC[CH]OO(102) <=> C[CH]CCCOO(48) origin: intra_H_migration
rxn: CC(C[CH]COO)OO(115) <=> [OH](22) + CC(CCC=O)OO(116) origin: intra_H_migration
rxn: CC(C[CH]COO)OO(115) <=> [OH](22) + CC(=O)CCCOO(112) origin: intra_H_migration

Non-identical kinetics! ❌
original:
rxn: CCC(CC)O[O](37) + CCCCCO[O](36) <=> oxygen(1) + CCC([O])CC(67) + CCCCC[O](69) origin: Peroxyl_Disproportionation
tested:
rxn: CCC(CC)O[O](35) + CCCCCO[O](36) <=> oxygen(1) + CCC([O])CC(69) + CCCCC[O](67) origin: Peroxyl_Disproportionation

k(1bar) 300K 400K 500K 600K 800K 1000K 1500K 2000K
k(T): 3.54 4.28 4.73 5.02 5.39 5.62 5.91 6.06
k(T): 8.02 7.64 7.35 7.11 6.75 6.48 5.99 5.64

kinetics: Arrhenius(A=(3.2e+12,'cm^3/(mol*s)'), n=0, Ea=(4.064,'kcal/mol'), T0=(1,'K'), comment="""Estimated from node Root_Ext-5R-R_7R!H->C_N-7C-inRing_Ext-5R-R in family Peroxyl_Disproportionation.""")
kinetics: Arrhenius(A=(3.18266e+20,'cm^3/(mol*s)'), n=-2.694, Ea=(-0.265,'kcal/mol'), T0=(1,'K'), comment="""Estimated from node Root_Ext-5R-R_7R!H->C_N-7C-inRing in family Peroxyl_Disproportionation.""")
kinetics: Estimated from node Root_Ext-5R-R_7R!H->C_N-7C-inRing_Ext-5R-R in family Peroxyl_Disproportionation.
kinetics: Estimated from node Root_Ext-5R-R_7R!H->C_N-7C-inRing in family Peroxyl_Disproportionation.

Errors occurred during edge comparison ⚠️ ERROR conda.cli.main_run:execute(148): `conda run python scripts/checkModels.py liquid_oxidation-edge stable_regression_results/liquid_oxidation/chemkin/chem_edge_annotated.inp stable_regression_results/liquid_oxidation/chemkin/species_edge_dictionary.txt test/regression/liquid_oxidation/chemkin/chem_edge_annotated.inp test/regression/liquid_oxidation/chemkin/species_edge_dictionary.txt` failed. (See above for error)
✅ All Observables varied by less than 0.100 on average between old model and new model in all conditions! liquid_oxidation Passed Observable Testing ✅

Regression test nitrogen:

Reference: Execution time (DD:HH:MM:SS): 00:00:00:57
Current: Execution time (DD:HH:MM:SS): 00:00:01:00
Reference: Memory used: 914.76 MB
Current: Memory used: 913.49 MB

nitrogen Passed Core Comparison ✅

Original model has 41 species.
Test model has 41 species. ✅
Original model has 360 reactions.
Test model has 360 reactions. ✅

nitrogen Failed Edge Comparison ❌

Original model has 133 species.
Test model has 133 species. ✅
Original model has 983 reactions.
Test model has 983 reactions. ✅

Non-identical thermo! ❌
original: O1[C]=N1
tested: O1[C]=N1

Hf(300K) S(300K) Cp(300K) Cp(400K) Cp(500K) Cp(600K) Cp(800K) Cp(1000K) Cp(1500K)
116.46 53.90 11.62 12.71 13.49 13.96 14.14 13.85 13.58
141.64 58.66 12.26 12.27 12.09 11.96 12.26 12.72 12.15

thermo: Thermo group additivity estimation: group(O2s-CdN3d) + group(N3d-OCd) + group(Cd-HN3dO) + ring(Cyclopropene) + radical(CdJ-NdO)
thermo: Thermo group additivity estimation: group(O2s-CdN3d) + group(N3d-OCd) + group(Cd-HN3dO) + ring(oxirene) + radical(CdJ-NdO)

Non-identical kinetics! ❌
original:
rxn: NCO(66) <=> O1[C]=N1(126) origin: Intra_R_Add_Endocyclic
tested:
rxn: NCO(66) <=> O1[C]=N1(126) origin: Intra_R_Add_Endocyclic

k(1bar) 300K 400K 500K 600K 800K 1000K 1500K 2000K
k(T): -49.54 -33.65 -24.16 -17.85 -10.01 -5.35 0.80 3.82
k(T): -66.25 -46.19 -34.19 -26.21 -16.28 -10.36 -2.54 1.31

kinetics: Arrhenius(A=(6.95187e+18,'s^-1'), n=-1.628, Ea=(88.327,'kcal/mol'), T0=(1,'K'), comment="""Estimated from node Backbone0_N-2R!H-inRing_N-1R!H-inRing_Sp-2R!H-1R!H in family Intra_R_Add_Endocyclic.""")
kinetics: Arrhenius(A=(6.95187e+18,'s^-1'), n=-1.628, Ea=(111.271,'kcal/mol'), T0=(1,'K'), comment="""Estimated from node Backbone0_N-2R!H-inRing_N-1R!H-inRing_Sp-2R!H-1R!H in family Intra_R_Add_Endocyclic.""")
Identical kinetics comments:
kinetics: Estimated from node Backbone0_N-2R!H-inRing_N-1R!H-inRing_Sp-2R!H-1R!H in family Intra_R_Add_Endocyclic.

Errors occurred during edge comparison ⚠️ ERROR conda.cli.main_run:execute(148): `conda run python scripts/checkModels.py nitrogen-edge stable_regression_results/nitrogen/chemkin/chem_edge_annotated.inp stable_regression_results/nitrogen/chemkin/species_edge_dictionary.txt test/regression/nitrogen/chemkin/chem_edge_annotated.inp test/regression/nitrogen/chemkin/species_edge_dictionary.txt` failed. (See above for error)
✅ All Observables varied by less than 0.200 on average between old model and new model in all conditions! nitrogen Passed Observable Testing ✅

Regression test oxidation:

Reference: Execution time (DD:HH:MM:SS): 00:00:01:31
Current: Execution time (DD:HH:MM:SS): 00:00:01:36
Reference: Memory used: 794.63 MB
Current: Memory used: 796.74 MB

oxidation Passed Core Comparison ✅

Original model has 59 species.
Test model has 59 species. ✅
Original model has 694 reactions.
Test model has 694 reactions. ✅

oxidation Passed Edge Comparison ✅

Original model has 230 species.
Test model has 230 species. ✅
Original model has 1524 reactions.
Test model has 1524 reactions. ✅

✅ All Observables varied by less than 0.500 on average between old model and new model in all conditions! oxidation Passed Observable Testing ✅
Errors occurred during observable testing ⚠️ WARNING: Initial mole fractions do not sum to one; normalizing.

Regression test sulfur:

Reference: Execution time (DD:HH:MM:SS): 00:00:00:35
Current: Execution time (DD:HH:MM:SS): 00:00:00:40
Reference: Memory used: 910.28 MB
Current: Memory used: 914.19 MB

sulfur Passed Core Comparison ✅

Original model has 27 species.
Test model has 27 species. ✅
Original model has 74 reactions.
Test model has 74 reactions. ✅

sulfur Failed Edge Comparison ❌

Original model has 89 species.
Test model has 89 species. ✅
Original model has 227 reactions.
Test model has 227 reactions. ✅
The original model has 1 reactions that the tested model does not have. ❌
rxn: O(4) + SO2(15) (+N2) <=> SO3(16) (+N2) origin: primarySulfurLibrary
The tested model has 1 reactions that the original model does not have. ❌
rxn: O(4) + SO2(15) (+N2) <=> SO3(16) (+N2) origin: primarySulfurLibrary

Errors occurred during edge comparison ⚠️ ERROR conda.cli.main_run:execute(148): `conda run python scripts/checkModels.py sulfur-edge stable_regression_results/sulfur/chemkin/chem_edge_annotated.inp stable_regression_results/sulfur/chemkin/species_edge_dictionary.txt test/regression/sulfur/chemkin/chem_edge_annotated.inp test/regression/sulfur/chemkin/species_edge_dictionary.txt` failed. (See above for error)
✅ All Observables varied by less than 0.100 on average between old model and new model in all conditions! sulfur Passed Observable Testing ✅

Regression test superminimal:

Reference: Execution time (DD:HH:MM:SS): 00:00:00:24
Current: Execution time (DD:HH:MM:SS): 00:00:00:25
Reference: Memory used: 958.82 MB
Current: Memory used: 955.40 MB

superminimal Passed Core Comparison ✅

Original model has 13 species.
Test model has 13 species. ✅
Original model has 21 reactions.
Test model has 21 reactions. ✅

superminimal Passed Edge Comparison ✅

Original model has 18 species.
Test model has 18 species. ✅
Original model has 28 reactions.
Test model has 28 reactions. ✅

Regression test RMS_constantVIdealGasReactor_superminimal:

Reference: Execution time (DD:HH:MM:SS): 00:00:02:58
Current: Execution time (DD:HH:MM:SS): 00:00:02:23
Reference: Memory used: 2479.17 MB
Current: Memory used: 2401.59 MB

RMS_constantVIdealGasReactor_superminimal Passed Core Comparison ✅

Original model has 13 species.
Test model has 13 species. ✅
Original model has 19 reactions.
Test model has 19 reactions. ✅

RMS_constantVIdealGasReactor_superminimal Passed Edge Comparison ✅

Original model has 13 species.
Test model has 13 species. ✅
Original model has 19 reactions.
Test model has 19 reactions. ✅

✅ All Observables varied by less than 0.100 on average between old model and new model in all conditions! RMS_constantVIdealGasReactor_superminimal Passed Observable Testing ✅

Regression test RMS_CSTR_liquid_oxidation:

Reference: Execution time (DD:HH:MM:SS): 00:00:28:16
Current: Execution time (DD:HH:MM:SS): 00:00:15:41
Reference: Memory used: 3544.81 MB
Current: Memory used: 2741.31 MB

RMS_CSTR_liquid_oxidation Failed Core Comparison ❌

Original model has 35 species.
Test model has 35 species. ✅
Original model has 142 reactions.
Test model has 141 reactions. ❌
The original model has 3 species that the tested model does not have. ❌
spc: CH3
spc: CCCC(C)O(46)
spc: CCCC=O(94)
The tested model has 3 species that the original model does not have. ❌
spc: CCH2
spc: [CH2]CCCCOO(66)
spc: CC1CC(C)O1(96)
The original model has 10 reactions that the tested model does not have. ❌
rxn: [CH3](10) + CCCC=O(94) <=> CCCC(C)[O](41) origin: R_Addition_MultipleBond
rxn: CCCC(C)[O](41) + pentane(2) <=> CC[CH]CC(7) + CCCC(C)O(46) origin: H_Abstraction
rxn: CCCC(C)[O](41) + pentane(2) <=> C[CH]CCC(11) + CCCC(C)O(46) origin: H_Abstraction
rxn: CCCC(C)[O](41) + pentane(2) <=> [CH2]CCCC(12) + CCCC(C)O(46) origin: H_Abstraction
rxn: CCCC(C)[O](41) + CCCC(C)OO(24) <=> CCCC(C)O[O](20) + CCCC(C)O(46) origin: H_Abstraction
rxn: CCCC(C)[O](41) + CCCC(C)OO(24) <=> C[CH]CC(C)OO(34) + CCCC(C)O(46) origin: H_Abstraction
rxn: CCCC(C)[O](41) + CCC(CC)OO(27) <=> CCC(CC)O[O](22) + CCCC(C)O(46) origin: H_Abstraction
rxn: CCCC(C)[O](41) + CCCCCOO(78) <=> CCCCCO[O](61) + CCCC(C)O(46) origin: H_Abstraction
rxn: CCCC(C)[O](41) + CCC(CC)OO(27) <=> [CH2]CC(CC)OO(38) + CCCC(C)O(46) origin: H_Abstraction
rxn: CCCC(C)[O](41) + CCCCCOO(78) <=> C[CH]CCCOO(65) + CCCC(C)O(46) origin: H_Abstraction
The tested model has 9 reactions that the original model does not have. ❌
rxn: C[CH2](6) + [CH2]CC(5) <=> pentane(2) origin: R_Recombination
rxn: C[CH]CC(C)OO(34) <=> [OH](26) + CC1CC(C)O1(96) origin: Cyclic_Ether_Formation
rxn: CCCCCO[O](61) <=> [CH2]CCCCOO(66) origin: intra_H_migration
rxn: [O]O(13) + [CH2]CCCCOO(66) <=> oxygen(1) + CCCCCOO(78) origin: H_Abstraction
rxn: OO(23) + [CH2]CCCCOO(66) <=> [O]O(13) + CCCCCOO(78) origin: H_Abstraction
rxn: [CH2]CCCCOO(66) + CCCC(C)OO(24) <=> CCCC(C)O[O](20) + CCCCCOO(78) origin: H_Abstraction
rxn: [CH2]CCCC(12) + [CH2]CCCCOO(66) <=> C=CCCC(17) + CCCCCOO(78) origin: Disproportionation
rxn: C[CH]CCC(11) + [CH2]CCCCOO(66) <=> C=CCCC(17) + CCCCCOO(78) origin: Disproportionation
rxn: [OH](26) + CCCCCOO(78) <=> H2O(42) + [CH2]CCCCOO(66) origin: H_Abstraction

Errors occurred during core comparison ⚠️ ERROR conda.cli.main_run:execute(148): `conda run python scripts/checkModels.py RMS_CSTR_liquid_oxidation-core stable_regression_results/RMS_CSTR_liquid_oxidation/chemkin/chem_annotated.inp stable_regression_results/RMS_CSTR_liquid_oxidation/chemkin/species_dictionary.txt test/regression/RMS_CSTR_liquid_oxidation/chemkin/chem_annotated.inp test/regression/RMS_CSTR_liquid_oxidation/chemkin/species_dictionary.txt` failed. (See above for error)
RMS_CSTR_liquid_oxidation Failed Edge Comparison ❌

Original model has 102 species.
Test model has 99 species. ❌
Original model has 387 reactions.
Test model has 384 reactions. ❌
The original model has 3 species that the tested model does not have. ❌
spc: [O]OOO(21)
spc: CCC(CC)OOOO
spc: CCCCCOOOO
The original model has 3 reactions that the tested model does not have. ❌
rxn: oxygen(1) + [O]O(13) <=> [O]OOO(21) origin: R_Recombination
rxn: oxygen(1) + CCC(CC)O[O](22) <=> CCC(CC)OOO[O](40) origin: R_Recombination
rxn: oxygen(1) + CCCCCO[O](61) <=> CCCCCOOO[O](77) origin: R_Recombination

Errors occurred during edge comparison ⚠️ ERROR conda.cli.main_run:execute(148): `conda run python scripts/checkModels.py RMS_CSTR_liquid_oxidation-edge stable_regression_results/RMS_CSTR_liquid_oxidation/chemkin/chem_edge_annotated.inp stable_regression_results/RMS_CSTR_liquid_oxidation/chemkin/species_edge_dictionary.txt test/regression/RMS_CSTR_liquid_oxidation/chemkin/chem_edge_annotated.inp test/regression/RMS_CSTR_liquid_oxidation/chemkin/species_edge_dictionary.txt` failed. (See above for error)
✅ All Observables varied by less than 0.100 on average between old model and new model in all conditions! RMS_CSTR_liquid_oxidation Passed Observable Testing ✅

Regression test fragment:

Reference: Execution time (DD:HH:MM:SS): 00:00:00:31
Current: Execution time (DD:HH:MM:SS): 00:00:00:36
Reference: Memory used: 758.49 MB
Current: Memory used: 758.44 MB

fragment Passed Core Comparison ✅

Original model has 10 species.
Test model has 10 species. ✅
Original model has 2 reactions.
Test model has 2 reactions. ✅

fragment Passed Edge Comparison ✅

Original model has 33 species.
Test model has 33 species. ✅
Original model has 47 reactions.
Test model has 47 reactions. ✅

✅ All Observables varied by less than 0.100 on average between old model and new model in all conditions! fragment Passed Observable Testing ✅
Errors occurred during observable testing ⚠️ WARNING: Initial mole fractions do not sum to one; normalizing.

Regression test RMS_constantVIdealGasReactor_fragment:

Reference: Execution time (DD:HH:MM:SS): 00:00:03:19
Current: Execution time (DD:HH:MM:SS): 00:00:03:01
Reference: Memory used: 2520.80 MB
Current: Memory used: 2477.75 MB

RMS_constantVIdealGasReactor_fragment Passed Core Comparison ✅

Original model has 10 species.
Test model has 10 species. ✅
Original model has 2 reactions.
Test model has 2 reactions. ✅

RMS_constantVIdealGasReactor_fragment Passed Edge Comparison ✅

Original model has 27 species.
Test model has 27 species. ✅
Original model has 24 reactions.
Test model has 24 reactions. ✅

✅ All Observables varied by less than 0.100 on average between old model and new model in all conditions! RMS_constantVIdealGasReactor_fragment Passed Observable Testing ✅
Errors occurred during observable testing ⚠️ WARNING: Initial mole fractions do not sum to one; normalizing.

Regression test minimal_surface:

Reference: Execution time (DD:HH:MM:SS): 00:00:00:29
Current: Execution time (DD:HH:MM:SS): 00:00:00:33
Reference: Memory used: 917.19 MB
Current: Memory used: 918.16 MB

minimal_surface Passed Core Comparison ✅

Original model has 11 species.
Test model has 11 species. ✅
Original model has 3 reactions.
Test model has 3 reactions. ✅

minimal_surface Passed Edge Comparison ✅

Original model has 38 species.
Test model has 38 species. ✅
Original model has 38 reactions.
Test model has 38 reactions. ✅

✅ All Observables varied by less than 0.500 on average between old model and new model in all conditions! minimal_surface Passed Observable Testing ✅

beep boop this comment was written by a bot 🤖

@alongd alongd left a comment

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

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

Thanks, @calvinp0! Added a few comments

import scipy.stats
if not all(np.isfinite(klist)):
raise ValueError("Rates must all be finite, not inf or NaN")
if any(klist==0):

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

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

I'd rather we filtered the bad samples where klist is built than raise here. The five reverse_*_rate helpers (reaction.py:1127, 1150, 1174, 1201, 1227) each build klist[i].

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

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

import scipy.stats
if not all(np.isfinite(klist)):
raise ValueError("Rates must all be finite, not inf or NaN")
if any(klist==0):

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

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

The same hole is left open in two places: SurfaceChargeTransfer.fit_to_data and StickingCoefficient.fit_to_data. Can you add this raise there as well?

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

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

added

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

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

also, StickingCoefficient.fit_to_data has no isfinite or sign validation at all unlike the other three, so inf/NaN rates reach the log fit too. i can make it a separate PR? or fold it into this?

Comment thread rmgpy/reaction.py Outdated
kf_list.append(self.get_rate_coefficient(condition[0], condition[1]))
try:
kf = self.get_rate_coefficient(temp, pressure)
except (ReactionError, KineticsError, TypeError, ValueError,

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

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

I'd remove TypeError and ValueError from here, they catch programmatic bugs, not scientific bugs. If you saw them here in the wild, then add a separate exception?

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

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

Agreed - narrowed to (ReactionError, KineticsError, ZeroDivisionError, OverflowError)`

@calvinp0 calvinp0 Aug 8, 2026

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

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

I want to raise though, with this comment and the previous one to extend the zero rate check - together if that check ever fires, we will skip the whole reaction with a traceback rather than one condition with a warning. and the proceeding comments were about preserving. i can instead raise Kinetics Error from that guard instead of Value Error?

Comment thread rmgpy/reaction.py Outdated
condition = '{0} K, {1:.1f} bar'.format(conditions[i][0], conditions[i][1] / 1e5)
violator_list.append([self, 'reverse', ratio, condition])
try:
kr = self.generate_reverse_rate_coefficient().get_rate_coefficient(temp, pressure)

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

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

generate_reverse_rate_coefficient() is inside the condition loop, so it gets rebuilt up to 4x per reaction. The new ValueError makes the same failure now recomputed and re-logged once per condition. Put it above?

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

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

true, fixed

Comment thread rmgpy/reaction.py
)
else:
reverse_checks.append((kr, limit_r, condition))
for direction, checks in (('forward', forward_checks), ('reverse', reverse_checks)):

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

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

If every condition gets skipped, check_model still prints "No collision rate violators found in the model's core." (main.py:1695), which isn't true. Could you return a skipped count alongside violator_list so the summary can say how much was actually checked?

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

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

true. it now returns (violator_list, skipped) and check_model tracks skipped check and skipped reactions separately, reporting when either is non-zero. but know this changed cpdef return type from list to tuple

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

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

since its a public API on the defintion change, people will need to re-compile again but i guess that is expected for every RMG update.

Comment thread rmgpy/solver/base.pyx
rderiv = kr[j] / kf[j] * C[ip[j, 0]] * C[ip[j, 1]]
else: # three reactants
rderiv = kr[j] / kf[j] * C[ip[j, 0]] * C[ip[j, 1]] * C[ip[j, 2]]
if not (Keq[j] > 0.0):

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

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

For a reversible reaction Keq can't be 0 (reaction.py:863 raises), so the only thing that trips this branch is NaN thermo. Worth a comment saying so.
The bigger thing: the irreversible sentinel isn't consistent. simple.pyx:277 sets Keq = np.inf, which passes the check and gives C/inf = 0.0. liquid.pyx:168, surface.pyx:315 and mbSampled.pyx:239 have no else, so Keq stays at the 0.0 from np.zeros_like and trips the check instead. Could we add else: self.kb[j] = 0.0; self.Keq[j] = np.inf to those three so they agree by design?

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

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

done - kb = 0.0; Keq = np.inf added to liquid, surface and mbSampled and kept the guard rather than assuming Keq > 0

calvinp0 added 3 commits July 31, 2026 18:02
Raising on a zero rate coefficient reports the problem at its source, but
it also fails a fit that could have succeeded: a reverse rate coefficient
that underflows at the cold end of the fitting range still has a perfectly
good fit over the rest of the range, and the underflowed samples carry no
information for a fit performed in log space.

Filter those samples out where klist is built, in the five reverse_*_rate
helpers, so the common case produces a usable fit instead of an error. The
validation in fit_to_data stays as a backstop for callers that construct
klist themselves, and is extended to StickingCoefficient.fit_to_data and
SurfaceChargeTransfer.fit_to_data, which perform the same log-space fit and
had the same gap.

Addresses review comments on rmgpy/kinetics/arrhenius.pyx.
SimpleReactor marks an irreversible reaction with kb = 0 and Keq = inf,
but the liquid, surface and mass-flow-batch reactors left Keq at the 0.0
from np.zeros_like. Consumers therefore had to accept either sentinel:
compute_rate_derivative's `not (Keq[j] > 0)` guard was written to catch
the zero case for that reason.

Set kb = 0 and Keq = inf explicitly in the other three reactors so the
sentinel agrees by design. The guard is unchanged and still correct, but
now the only thing that can trip it is NaN thermochemistry -- a reversible
reaction cannot have Keq == 0, since get_equilibrium_constant raises in
that case. Say so in the comment.

Addresses review comment on rmgpy/solver/base.pyx.
Three changes to check_collision_limit_violation, all from review:

Narrow the caught exceptions to ReactionError, KineticsError,
ZeroDivisionError and OverflowError. TypeError and ValueError indicate
programmatic bugs rather than a reaction whose kinetics cannot be
evaluated, and swallowing them here would hide real defects.

Generate the reverse rate coefficient once per reaction rather than once
per condition. It does not depend on T or P, so it was being rebuilt up to
four times, and a reaction whose reverse coefficient cannot be generated
was reported once per condition instead of once.

Return the number of skipped direction/condition pairs alongside the
violator list, and count skipped reactions in check_model, so that a model
where nothing could be evaluated no longer reports "No collision rate
violators found in the model's core" as though it had been fully checked.
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.

4 participants