From bfdb12bee6f933a2d19ef44e099251b8b2d799eb Mon Sep 17 00:00:00 2001 From: Richard West Date: Fri, 15 May 2026 11:01:08 -0400 Subject: [PATCH 1/5] Fix off-by-one in CHEMKIN THERM default-range line I/O read_thermo_block sliced 9 chars per slot ([0:9], [10:19], [20:29]) instead of the documented 10 (cols 1-10, 11-20, 21-30). With a value right-anchored against the right edge of a slot, e.g. ' 1000' (six spaces + 1000) the slice produced ' 100' and the value was misread as 100. Real chemkin files happen not to pack that tightly, so the bug was latent for years. save_chemkin_surface_file also wrote the default-range header with 4 leading spaces instead of 3, shifting every value one column right. The gas-phase save_chemkin_file was corrected in 37dbd53a but the surface variant (added later in 5a65011) inherited a pre-fix copy and never got the same correction, so every surface chemkin file RMG has written since carries the misaligned header. Add a regression test that parses a default-range line packing values to both edges of each 10-column slot. Co-Authored-By: Claude Opus 4.7 (1M context) --- rmgpy/chemkin.pyx | 8 ++++---- test/rmgpy/chemkinTest.py | 22 ++++++++++++++++++++++ 2 files changed, 26 insertions(+), 4 deletions(-) diff --git a/rmgpy/chemkin.pyx b/rmgpy/chemkin.pyx index 9aac0dcabf..1c09b753da 100644 --- a/rmgpy/chemkin.pyx +++ b/rmgpy/chemkin.pyx @@ -1250,9 +1250,9 @@ def read_thermo_block(f, species_dict): meaningfulline, comment = remove_comment_from_line(line) Tmin = Tint = Tmax = None try: - Tmin = float(meaningfulline[0:9].strip()) - Tint = float(meaningfulline[10:19].strip()) - Tmax = float(meaningfulline[20:29].strip()) + Tmin = float(meaningfulline[0:10].strip()) + Tint = float(meaningfulline[10:20].strip()) + Tmax = float(meaningfulline[20:30].strip()) if [Tmin, Tint, Tmax] != [float(i) for i in meaningfulline.split()[0:3]]: logging.warning("Default temperature range line {0!r} may be badly formatted.".format(line)) logging.warning("It should have Tmin in columns 1-10, Tmid in columns 11-20, and Tmax in columns 21-30") @@ -2198,7 +2198,7 @@ def save_chemkin_surface_file(path, species, reactions, verbose=True, check_for_ # Thermodynamics section f.write('THERM ALL\n') - f.write(' 300.000 1000.000 5000.000\n\n') + f.write(' 300.000 1000.000 5000.000\n\n') for spec in sorted_species: f.write(write_thermo_entry(spec, verbose=verbose)) f.write('\n') diff --git a/test/rmgpy/chemkinTest.py b/test/rmgpy/chemkinTest.py index e72d227a92..8c42c434b5 100644 --- a/test/rmgpy/chemkinTest.py +++ b/test/rmgpy/chemkinTest.py @@ -27,6 +27,7 @@ # # ############################################################################### +import io import os from unittest import mock @@ -39,6 +40,7 @@ mark_duplicate_reactions, read_kinetics_entry, read_reaction_comments, + read_thermo_block, read_thermo_entry, save_chemkin_file, save_chemkin_surface_file, @@ -144,6 +146,26 @@ def test_read_thermo_entry_no_temperature_range(self): assert formula == {"H": 6, "C": 2} assert isinstance(thermo, NASA) + @mock.patch("rmgpy.chemkin.logging") + def test_read_thermo_block_temperature_header_columns(self, mock_logging): + # Per the CHEMKIN spec the default-range line is fixed-width: + # Tmin in cols 1-10, Tint in cols 11-20, Tmax in cols 21-30. + # Slot 1: 150 left-anchored (leftmost digit in col 1). + # Slot 2: 1000 right-anchored (rightmost digit in col 20) + # ' 1000' that could be misread as 100. + # Slot 3: 9999 right-anchored (rightmost digit in col 30). + header = "150 1000 9999" + assert len(header) == 30 + f = io.StringIO("THERM ALL\n" + header + "\nEND\n") + + read_thermo_block(f, species_dict={}) + + mock_logging.info.assert_any_call( + "Thermo file has default temperature range 150.0 to 1000.0 and 1000.0 to 9999.0" + ) + for call in mock_logging.warning.call_args_list: + assert "badly formatted" not in call.args[0] + def test_read_and_write_and_read_template_reaction_family_for_minimal_example(self): """ This example tests if family and templates info can be correctly From d1473f828b88cfdf3eb18e309a5f9d821f43f973 Mon Sep 17 00:00:00 2001 From: Richard West Date: Fri, 15 May 2026 11:02:04 -0400 Subject: [PATCH 2/5] Realign THERM default-range header in committed CHEMKIN fixtures Replace ' 300.000 1000.000 5000.000' (4 leading spaces) with ' 300.000 1000.000 5000.000' (3 leading spaces) in every tracked fixture. These files were originally written by the buggy save_chemkin_surface_file (or by a pre-fix copy of the gas-phase writer), so the header was one column wider than the strict cols 1-10/11-20/21-30 spec. The parser's fallback handled the misaligned form, so RMG behaviour is unchanged; the goal is just to stop the "Thermo file has no default temperature ranges" warning from appearing in test logs. mainTest.py contains inline chemkin strings used as test fixtures and gets the same one-space fix. Co-Authored-By: Claude Opus 4.7 (1M context) --- examples/scripts/diffModels/Models/DFT_chem.inp | 2 +- examples/scripts/diffModels/Models/Primary_chem.inp | 2 +- ipython/data/pathway_analysis/minimal/chem.inp | 2 +- ipython/data/regression/new/chem_annotated.inp | 2 +- ipython/data/regression/old/chem_annotated.inp | 2 +- rmgpy/solver/files/collider_model/chem.inp | 2 +- rmgpy/solver/files/listener/chemkin/chem.inp | 2 +- rmgpy/tools/data/diffmodels/chem1.inp | 2 +- rmgpy/tools/data/diffmodels/chem2.inp | 2 +- rmgpy/tools/data/diffmodels/chem3.inp | 2 +- rmgpy/tools/data/diffmodels/surf_model/chem_surface1.inp | 2 +- rmgpy/tools/data/diffmodels/surf_model/chem_surface2.inp | 2 +- rmgpy/tools/data/flux/chemkin/chem.inp | 2 +- rmgpy/tools/data/regression/benchmark/chem_annotated.inp | 2 +- rmgpy/tools/data/regression/tested/chem_annotated.inp | 2 +- rmgpy/tools/data/sim/liquid/chem.inp | 2 +- rmgpy/tools/data/sim/mbSampled/chem.inp | 2 +- rmgpy/tools/data/sim/simple/chem.inp | 2 +- rmgpy/tools/data/various_kinetics/chem_annotated.inp | 2 +- test/rmgpy/rmg/mainTest.py | 6 +++--- test/rmgpy/test_data/chemkin/chemkin_py/NC/chem.inp | 2 +- test/rmgpy/test_data/chemkin/chemkin_py/minimal/chem.inp | 2 +- test/rmgpy/test_data/chemkin/chemkin_py/pdd/chem.inp | 2 +- .../rmgpy/test_data/chemkin/chemkin_py/surface/chem-gas.inp | 2 +- test/rmgpy/test_data/parsing_data/chem_annotated.inp | 2 +- test/rmgpy/test_data/saveOutputHTML/eg6/chem_annotated.inp | 2 +- 26 files changed, 28 insertions(+), 28 deletions(-) diff --git a/examples/scripts/diffModels/Models/DFT_chem.inp b/examples/scripts/diffModels/Models/DFT_chem.inp index c8fa2c3641..86b6651852 100644 --- a/examples/scripts/diffModels/Models/DFT_chem.inp +++ b/examples/scripts/diffModels/Models/DFT_chem.inp @@ -57,7 +57,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ! Thermo library: primaryThermoLibrary Ar Ar 1 G 200.000 6000.000 1000.00 1 diff --git a/examples/scripts/diffModels/Models/Primary_chem.inp b/examples/scripts/diffModels/Models/Primary_chem.inp index 96fdd0677b..33ff2aaa90 100644 --- a/examples/scripts/diffModels/Models/Primary_chem.inp +++ b/examples/scripts/diffModels/Models/Primary_chem.inp @@ -53,7 +53,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ! Thermo library: primaryThermoLibrary Ar Ar 1 G 200.000 6000.000 1000.00 1 diff --git a/ipython/data/pathway_analysis/minimal/chem.inp b/ipython/data/pathway_analysis/minimal/chem.inp index adc95e4d95..9f7180591c 100644 --- a/ipython/data/pathway_analysis/minimal/chem.inp +++ b/ipython/data/pathway_analysis/minimal/chem.inp @@ -20,7 +20,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ! Thermo library: primaryThermoLibrary Ar Ar1 G200.000 6000.000 1000.00 1 diff --git a/ipython/data/regression/new/chem_annotated.inp b/ipython/data/regression/new/chem_annotated.inp index 1716760401..8daf18891c 100644 --- a/ipython/data/regression/new/chem_annotated.inp +++ b/ipython/data/regression/new/chem_annotated.inp @@ -25,7 +25,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ! Thermo library: primaryThermoLibrary Ar Ar1 G200.000 6000.000 1000.00 1 diff --git a/ipython/data/regression/old/chem_annotated.inp b/ipython/data/regression/old/chem_annotated.inp index 1ace47d026..926109b6e2 100644 --- a/ipython/data/regression/old/chem_annotated.inp +++ b/ipython/data/regression/old/chem_annotated.inp @@ -25,7 +25,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ! Thermo library: primaryThermoLibrary Ar Ar1 G200.000 6000.000 1000.00 1 diff --git a/rmgpy/solver/files/collider_model/chem.inp b/rmgpy/solver/files/collider_model/chem.inp index 9ec405ec2b..baa5e30d4a 100644 --- a/rmgpy/solver/files/collider_model/chem.inp +++ b/rmgpy/solver/files/collider_model/chem.inp @@ -20,7 +20,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 Ar Ar1 G200.000 6000.000 1000.00 1 2.50000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 2 diff --git a/rmgpy/solver/files/listener/chemkin/chem.inp b/rmgpy/solver/files/listener/chemkin/chem.inp index 93290297d9..d54bf1f847 100644 --- a/rmgpy/solver/files/listener/chemkin/chem.inp +++ b/rmgpy/solver/files/listener/chemkin/chem.inp @@ -15,7 +15,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 Ar Ar1 G200.000 6000.000 1000.00 1 2.50000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 2 diff --git a/rmgpy/tools/data/diffmodels/chem1.inp b/rmgpy/tools/data/diffmodels/chem1.inp index 7fc3771f7d..35fcd82510 100644 --- a/rmgpy/tools/data/diffmodels/chem1.inp +++ b/rmgpy/tools/data/diffmodels/chem1.inp @@ -19,7 +19,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ethane(1) C 2 H 6 G100.000 5000.000 954.52 1 4.58991157E+00 1.41506360E-02-4.75954111E-06 8.60275134E-10-6.21700682E-14 2 diff --git a/rmgpy/tools/data/diffmodels/chem2.inp b/rmgpy/tools/data/diffmodels/chem2.inp index 7fc3771f7d..35fcd82510 100644 --- a/rmgpy/tools/data/diffmodels/chem2.inp +++ b/rmgpy/tools/data/diffmodels/chem2.inp @@ -19,7 +19,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ethane(1) C 2 H 6 G100.000 5000.000 954.52 1 4.58991157E+00 1.41506360E-02-4.75954111E-06 8.60275134E-10-6.21700682E-14 2 diff --git a/rmgpy/tools/data/diffmodels/chem3.inp b/rmgpy/tools/data/diffmodels/chem3.inp index 12570fd6e4..397571f60b 100644 --- a/rmgpy/tools/data/diffmodels/chem3.inp +++ b/rmgpy/tools/data/diffmodels/chem3.inp @@ -12,7 +12,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ! GRI-Mech3.0 ethane H 6C 2 G 100.000 5000.000 1002.57 1 2.56122991E+00 1.83185519E-02-7.41959133E-06 1.38553537E-09-9.75184837E-14 2 diff --git a/rmgpy/tools/data/diffmodels/surf_model/chem_surface1.inp b/rmgpy/tools/data/diffmodels/surf_model/chem_surface1.inp index 582fcf235f..601f634e86 100644 --- a/rmgpy/tools/data/diffmodels/surf_model/chem_surface1.inp +++ b/rmgpy/tools/data/diffmodels/surf_model/chem_surface1.inp @@ -21,7 +21,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ! Thermo library: surfaceThermoPt111 Binding energy corrected by LSR () from Pt111 X(1) X 1 G 100.000 5000.000 1554.83 1 diff --git a/rmgpy/tools/data/diffmodels/surf_model/chem_surface2.inp b/rmgpy/tools/data/diffmodels/surf_model/chem_surface2.inp index a01200386d..8fbfca2f6a 100644 --- a/rmgpy/tools/data/diffmodels/surf_model/chem_surface2.inp +++ b/rmgpy/tools/data/diffmodels/surf_model/chem_surface2.inp @@ -20,7 +20,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ! Thermo library: surfaceThermoPt111 Binding energy corrected by LSR () from Pt111 X(1) X 1 G 100.000 5000.000 1554.81 1 diff --git a/rmgpy/tools/data/flux/chemkin/chem.inp b/rmgpy/tools/data/flux/chemkin/chem.inp index 4b6d1d5e35..0de3f5962e 100644 --- a/rmgpy/tools/data/flux/chemkin/chem.inp +++ b/rmgpy/tools/data/flux/chemkin/chem.inp @@ -17,7 +17,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 Ar Ar1 G200.000 6000.000 1000.00 1 2.50000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 2 diff --git a/rmgpy/tools/data/regression/benchmark/chem_annotated.inp b/rmgpy/tools/data/regression/benchmark/chem_annotated.inp index 1ace47d026..926109b6e2 100644 --- a/rmgpy/tools/data/regression/benchmark/chem_annotated.inp +++ b/rmgpy/tools/data/regression/benchmark/chem_annotated.inp @@ -25,7 +25,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ! Thermo library: primaryThermoLibrary Ar Ar1 G200.000 6000.000 1000.00 1 diff --git a/rmgpy/tools/data/regression/tested/chem_annotated.inp b/rmgpy/tools/data/regression/tested/chem_annotated.inp index 1716760401..8daf18891c 100644 --- a/rmgpy/tools/data/regression/tested/chem_annotated.inp +++ b/rmgpy/tools/data/regression/tested/chem_annotated.inp @@ -25,7 +25,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ! Thermo library: primaryThermoLibrary Ar Ar1 G200.000 6000.000 1000.00 1 diff --git a/rmgpy/tools/data/sim/liquid/chem.inp b/rmgpy/tools/data/sim/liquid/chem.inp index cdf7f9c1de..9088fb880f 100644 --- a/rmgpy/tools/data/sim/liquid/chem.inp +++ b/rmgpy/tools/data/sim/liquid/chem.inp @@ -49,7 +49,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 Ar Ar 1 G 100.000 5000.000 3142.77 1 2.50000000E+00 3.32026413E-12-1.37161245E-15 2.48094222E-19-1.65799619E-23 2 diff --git a/rmgpy/tools/data/sim/mbSampled/chem.inp b/rmgpy/tools/data/sim/mbSampled/chem.inp index 13a454f23d..637d5d0e5a 100644 --- a/rmgpy/tools/data/sim/mbSampled/chem.inp +++ b/rmgpy/tools/data/sim/mbSampled/chem.inp @@ -37,7 +37,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ! Narayanaswamy C6H5 H 5C 6 G 100.000 5000.000 845.04 1 diff --git a/rmgpy/tools/data/sim/simple/chem.inp b/rmgpy/tools/data/sim/simple/chem.inp index 7fc3771f7d..35fcd82510 100644 --- a/rmgpy/tools/data/sim/simple/chem.inp +++ b/rmgpy/tools/data/sim/simple/chem.inp @@ -19,7 +19,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ethane(1) C 2 H 6 G100.000 5000.000 954.52 1 4.58991157E+00 1.41506360E-02-4.75954111E-06 8.60275134E-10-6.21700682E-14 2 diff --git a/rmgpy/tools/data/various_kinetics/chem_annotated.inp b/rmgpy/tools/data/various_kinetics/chem_annotated.inp index 54acce1ae8..e299354b2c 100644 --- a/rmgpy/tools/data/various_kinetics/chem_annotated.inp +++ b/rmgpy/tools/data/various_kinetics/chem_annotated.inp @@ -19,7 +19,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ! ! Thermo group additivity estimation: group(Cs-CsHHH) + gauche(Cs(CsRRR)) + other(R) + group(Cs-CsHHH) + gauche(Cs(CsRRR)) + other(R) diff --git a/test/rmgpy/rmg/mainTest.py b/test/rmgpy/rmg/mainTest.py index bb5df4bdd8..d3c5ace343 100644 --- a/test/rmgpy/rmg/mainTest.py +++ b/test/rmgpy/rmg/mainTest.py @@ -447,7 +447,7 @@ def setup_class(self): END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ethane(1) H 6 C 2 G100.000 5000.000 954.52 1 4.58987205E+00 1.41507042E-02-4.75958084E-06 8.60284590E-10-6.21708569E-14 2 @@ -483,7 +483,7 @@ def setup_class(self): END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ethane(1) H 6 C 2 G100.000 5000.000 954.52 1 4.58987205E+00 1.41507042E-02-4.75958084E-06 8.60284590E-10-6.21708569E-14 2 @@ -523,7 +523,7 @@ def setup_class(self): END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ethane(1) H 6 C 2 G100.000 5000.000 954.52 1 4.58987205E+00 1.41507042E-02-4.75958084E-06 8.60284590E-10-6.21708569E-14 2 diff --git a/test/rmgpy/test_data/chemkin/chemkin_py/NC/chem.inp b/test/rmgpy/test_data/chemkin/chemkin_py/NC/chem.inp index 769021b250..b2a6b06a7e 100644 --- a/test/rmgpy/test_data/chemkin/chemkin_py/NC/chem.inp +++ b/test/rmgpy/test_data/chemkin/chemkin_py/NC/chem.inp @@ -68,7 +68,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 Ar Ar1 G200.000 6000.000 1000.00 1 2.50000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 2 diff --git a/test/rmgpy/test_data/chemkin/chemkin_py/minimal/chem.inp b/test/rmgpy/test_data/chemkin/chemkin_py/minimal/chem.inp index 9fc1a9b3f6..cf51bc2d75 100644 --- a/test/rmgpy/test_data/chemkin/chemkin_py/minimal/chem.inp +++ b/test/rmgpy/test_data/chemkin/chemkin_py/minimal/chem.inp @@ -14,7 +14,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ! Thermo library: primaryThermoLibrary Ar Ar1 G200.000 6000.000 1000.00 1 diff --git a/test/rmgpy/test_data/chemkin/chemkin_py/pdd/chem.inp b/test/rmgpy/test_data/chemkin/chemkin_py/pdd/chem.inp index b33bc969c8..8bc76faa40 100644 --- a/test/rmgpy/test_data/chemkin/chemkin_py/pdd/chem.inp +++ b/test/rmgpy/test_data/chemkin/chemkin_py/pdd/chem.inp @@ -15,7 +15,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ! Thermo library: primaryThermoLibrary Ar Ar1 G200.000 6000.000 1000.00 1 diff --git a/test/rmgpy/test_data/chemkin/chemkin_py/surface/chem-gas.inp b/test/rmgpy/test_data/chemkin/chemkin_py/surface/chem-gas.inp index 0815daee86..ea252ab0d6 100644 --- a/test/rmgpy/test_data/chemkin/chemkin_py/surface/chem-gas.inp +++ b/test/rmgpy/test_data/chemkin/chemkin_py/surface/chem-gas.inp @@ -11,7 +11,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 Ar Ar1 G200.000 6000.000 1000.00 1 2.50000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 2 diff --git a/test/rmgpy/test_data/parsing_data/chem_annotated.inp b/test/rmgpy/test_data/parsing_data/chem_annotated.inp index 2e81aa0fa7..a85f736a42 100644 --- a/test/rmgpy/test_data/parsing_data/chem_annotated.inp +++ b/test/rmgpy/test_data/parsing_data/chem_annotated.inp @@ -31,7 +31,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ! Thermo library: primaryThermoLibrary O(2) O 1 G100.000 5000.000 4563.27 1 diff --git a/test/rmgpy/test_data/saveOutputHTML/eg6/chem_annotated.inp b/test/rmgpy/test_data/saveOutputHTML/eg6/chem_annotated.inp index 6e7ee71656..b0c721baeb 100644 --- a/test/rmgpy/test_data/saveOutputHTML/eg6/chem_annotated.inp +++ b/test/rmgpy/test_data/saveOutputHTML/eg6/chem_annotated.inp @@ -22,7 +22,7 @@ END THERM ALL - 300.000 1000.000 5000.000 + 300.000 1000.000 5000.000 ! Thermo library: primaryThermoLibrary Ar Ar1 G200.000 6000.000 1000.00 1 From 5fa7bf7393704fc49701457d141d2a4be6dc77c5 Mon Sep 17 00:00:00 2001 From: Richard West Date: Thu, 21 May 2026 13:46:42 -0400 Subject: [PATCH 3/5] Chemkin was not writing isotopes correctly. In this patch, we skip using the Molecule.get_element_count() because that returns just the element.symbol, and deuterium is H. --- rmgpy/chemkin.pyx | 7 ++++++- rmgpy/molecule/element.py | 2 +- test/rmgpy/chemkinTest.py | 14 ++++++++++++++ 3 files changed, 21 insertions(+), 2 deletions(-) diff --git a/rmgpy/chemkin.pyx b/rmgpy/chemkin.pyx index 1c09b753da..5bdfd9e7d2 100644 --- a/rmgpy/chemkin.pyx +++ b/rmgpy/chemkin.pyx @@ -1588,8 +1588,13 @@ def write_thermo_entry(species, element_counts=None, verbose=True): assert thermo.polynomials[1].cm2 == 0 and thermo.polynomials[1].cm1 == 0 # Determine the number of each type of element in the molecule + # Need to use the element's chemkin name, not the element symbol, because of isotopes. + # so we can't just use molecule[0].get_element_count(). if element_counts is None: - element_counts = get_element_count(species.molecule[0]) + element_counts = {} + for atom in species.molecule[0].atoms: + element = atom.element.chemkin_name + element_counts[element] = element_counts.get(element, 0) + 1 # Sort the element_counts dictionary so that it's C, H, Al, B, Cl, D, etc. # if there's any C, else Al, B, Cl, D, H, if not. This is the "Hill" system diff --git a/rmgpy/molecule/element.py b/rmgpy/molecule/element.py index 6ec975fcf9..4f0a6a5704 100644 --- a/rmgpy/molecule/element.py +++ b/rmgpy/molecule/element.py @@ -77,7 +77,7 @@ def __init__(self, number, symbol, name, mass, isotope=-1, chemkin_name=None): self.name = name self.mass = mass self.isotope = isotope - self.chemkin_name = chemkin_name or self.name + self.chemkin_name = chemkin_name or self.symbol if symbol in {'X','L','R','e'}: self.cov_radius = 0 else: diff --git a/test/rmgpy/chemkinTest.py b/test/rmgpy/chemkinTest.py index 8c42c434b5..d2683efdaa 100644 --- a/test/rmgpy/chemkinTest.py +++ b/test/rmgpy/chemkinTest.py @@ -791,6 +791,20 @@ def test_read_thermo_block(self): assert formula == {"H": 6, "C": 2} assert self.nasa.is_identical_to(thermo) + def test_write_thermo_block_for_isotope_uses_chemkin_name(self): + """Isotopic atoms should use their Chemkin names in thermo composition.""" + deuterium = Species().from_adjacency_list( + "1 H u0 p0 c0 i2 {2,S}\n" + "2 H u0 p0 c0 {1,S}" + ) + deuterium.thermo = self.nasa + + result = write_thermo_entry(deuterium, verbose=False) + + first_line = result.splitlines()[0] + assert "D 1" in first_line + assert "H 1" in first_line + def test_write_thermo_block_5_elem(self): """Test that we can write a thermo block for a species with 5 elements""" species = Species().from_adjacency_list( From 9fce4517db93d6c9cd7c51be0562ff961417f67e Mon Sep 17 00:00:00 2001 From: Richard West Date: Fri, 22 May 2026 14:05:30 -0400 Subject: [PATCH 4/5] Optimization in chemkin write_thermo_entry --- rmgpy/chemkin.pyx | 86 +++++++++++++++++++++++++++-------------------- 1 file changed, 50 insertions(+), 36 deletions(-) diff --git a/rmgpy/chemkin.pyx b/rmgpy/chemkin.pyx index 5bdfd9e7d2..e9f5e78df0 100644 --- a/rmgpy/chemkin.pyx +++ b/rmgpy/chemkin.pyx @@ -52,6 +52,7 @@ from rmgpy.reaction import Reaction from rmgpy.rmg.pdep import PDepNetwork, PDepReaction from rmgpy.species import Species from rmgpy.thermo import NASAPolynomial, NASA +from rmgpy.thermo.nasa cimport NASA, NASAPolynomial from rmgpy.transport import TransportData from rmgpy.util import make_output_subdirectory @@ -1568,42 +1569,55 @@ def get_species_identifier(species): ################################################################################ -def write_thermo_entry(species, element_counts=None, verbose=True): +def write_thermo_entry(species, element_counts=None, bint verbose=True): """ Return a string representation of the NASA model readable by Chemkin. To use this method you must have exactly two NASA polynomials in your model, and you must use the seven-coefficient forms for each. """ - - thermo = species.get_thermo_data() - - if not isinstance(thermo, NASA): + cdef NASA thermo + cdef NASAPolynomial poly_low, poly_high + cdef dict counts + cdef list sorted_elements, elements, short_lines + cdef bint extended_syntax + cdef int count, isotope + cdef str string, line, short_line, chemkin_name, symbol, elem_1, elem_2 + cdef object thermo_data + + thermo_data = species.get_thermo_data() + + if not isinstance(thermo_data, NASA): raise ChemkinError('Cannot generate Chemkin string for species "{0}": ' 'Thermodynamics data must be a NASA object.'.format(species)) + thermo = thermo_data assert len(thermo.polynomials) == 2 - assert thermo.polynomials[0].Tmin.value_si < thermo.polynomials[1].Tmin.value_si - assert thermo.polynomials[0].Tmax.value_si == thermo.polynomials[1].Tmin.value_si - assert thermo.polynomials[0].cm2 == 0 and thermo.polynomials[0].cm1 == 0 - assert thermo.polynomials[1].cm2 == 0 and thermo.polynomials[1].cm1 == 0 + poly_low = thermo.polynomials[0] + poly_high = thermo.polynomials[1] + assert poly_low.Tmin.value_si < poly_high.Tmin.value_si + assert poly_low.Tmax.value_si == poly_high.Tmin.value_si + assert poly_low.cm2 == 0 and poly_low.cm1 == 0 + assert poly_high.cm2 == 0 and poly_high.cm1 == 0 # Determine the number of each type of element in the molecule # Need to use the element's chemkin name, not the element symbol, because of isotopes. # so we can't just use molecule[0].get_element_count(). if element_counts is None: - element_counts = {} - for atom in species.molecule[0].atoms: - element = atom.element.chemkin_name - element_counts[element] = element_counts.get(element, 0) + 1 + counts = {} + for atom in species.molecule[0].vertices: + chemkin_name = atom.element.chemkin_name + counts[chemkin_name] = counts.get(chemkin_name, 0) + 1 + else: + counts = element_counts # Sort the element_counts dictionary so that it's C, H, Al, B, Cl, D, etc. # if there's any C, else Al, B, Cl, D, H, if not. This is the "Hill" system # done by Molecule.get_formula - if 'C' in element_counts: - sorted_elements = sorted(element_counts, key = lambda e: {'C':'0','H':'1'}.get(e, e)) + if 'C' in counts: + sorted_elements = sorted(counts, key=lambda e: {'C': '0', 'H': '1'}.get(e, e)) else: - sorted_elements = sorted(element_counts) - element_counts = {e: element_counts[e] for e in sorted_elements} + sorted_elements = sorted(counts) + counts = {e: counts[e] for e in sorted_elements} string = '' # Write thermo comments @@ -1618,9 +1632,9 @@ def write_thermo_entry(species, element_counts=None, verbose=True): string += "! {0}\n".format(line) # Compile element count string - extended_syntax = len(element_counts) > 4 # If there are more than 4 elements, use extended syntax + extended_syntax = len(counts) > 4 # If there are more than 4 elements, use extended syntax elements = [] - for key, count in element_counts.items(): + for key, count in counts.items(): if isinstance(key, tuple): symbol, isotope = key chemkin_name = get_element(symbol, isotope=isotope).chemkin_name @@ -1650,31 +1664,31 @@ def write_thermo_entry(species, element_counts=None, verbose=True): string += '{ident:<16} {elem_1:<20}G{Tmin:>10.3f}{Tint:>10.3f}{Tmax:>8.2f} 1{elem_2}\n'.format( ident=get_species_identifier(species), elem_1=elem_1, - Tmin=thermo.polynomials[0].Tmin.value_si, - Tint=thermo.polynomials[1].Tmax.value_si, - Tmax=thermo.polynomials[0].Tmax.value_si, + Tmin=poly_low.Tmin.value_si, + Tint=poly_high.Tmax.value_si, + Tmax=poly_low.Tmax.value_si, elem_2=elem_2, ) # Line 2 - string += '{0:< 15.8E}{1:< 15.8E}{2:< 15.8E}{3:< 15.8E}{4:< 15.8E} 2\n'.format(thermo.polynomials[1].c0, - thermo.polynomials[1].c1, - thermo.polynomials[1].c2, - thermo.polynomials[1].c3, - thermo.polynomials[1].c4) + string += '{0:< 15.8E}{1:< 15.8E}{2:< 15.8E}{3:< 15.8E}{4:< 15.8E} 2\n'.format(poly_high.c0, + poly_high.c1, + poly_high.c2, + poly_high.c3, + poly_high.c4) # Line 3 - string += '{0:< 15.8E}{1:< 15.8E}{2:< 15.8E}{3:< 15.8E}{4:< 15.8E} 3\n'.format(thermo.polynomials[1].c5, - thermo.polynomials[1].c6, - thermo.polynomials[0].c0, - thermo.polynomials[0].c1, - thermo.polynomials[0].c2) + string += '{0:< 15.8E}{1:< 15.8E}{2:< 15.8E}{3:< 15.8E}{4:< 15.8E} 3\n'.format(poly_high.c5, + poly_high.c6, + poly_low.c0, + poly_low.c1, + poly_low.c2) # Line 4 - string += '{0:< 15.8E}{1:< 15.8E}{2:< 15.8E}{3:< 15.8E} 4\n'.format(thermo.polynomials[0].c3, - thermo.polynomials[0].c4, - thermo.polynomials[0].c5, - thermo.polynomials[0].c6) + string += '{0:< 15.8E}{1:< 15.8E}{2:< 15.8E}{3:< 15.8E} 4\n'.format(poly_low.c3, + poly_low.c4, + poly_low.c5, + poly_low.c6) return string From fc1c6bbae2716810596f6a214e52e98342101d3d Mon Sep 17 00:00:00 2001 From: Richard West Date: Fri, 22 May 2026 14:34:22 -0400 Subject: [PATCH 5/5] Add unit test for warning about badly formatted THERM line --- test/rmgpy/chemkinTest.py | 9 +++++++++ 1 file changed, 9 insertions(+) diff --git a/test/rmgpy/chemkinTest.py b/test/rmgpy/chemkinTest.py index d2683efdaa..138943899b 100644 --- a/test/rmgpy/chemkinTest.py +++ b/test/rmgpy/chemkinTest.py @@ -166,6 +166,15 @@ def test_read_thermo_block_temperature_header_columns(self, mock_logging): for call in mock_logging.warning.call_args_list: assert "badly formatted" not in call.args[0] + @mock.patch("rmgpy.chemkin.logging") + def test_read_thermo_block_warns_for_badly_formatted_temperature_header(self, mock_logging): + header = "15000000001 9999 5000" + f = io.StringIO("THERM ALL\n" + header + "\nEND\n") + + read_thermo_block(f, species_dict={}) + + assert any("badly formatted" in call.args[0] for call in mock_logging.warning.call_args_list) + def test_read_and_write_and_read_template_reaction_family_for_minimal_example(self): """ This example tests if family and templates info can be correctly