Discontinuity in Klippenstein_Glarborg2016 thermo lib

[sevyharris]

I found a discontinuity in RMG's Klippenstein_Glarborg2016 thermo library for formic acid.

The RMG entry points to this paper https://doi.org/10.1016/j.proci.2014.05.091 and says B Ruscic modified the polynomials.
But the paper in question says the thermo comes from this Argon/Technion report: https://publications.anl.gov/anlpubs/2005/07/53802.pdf

I didn't see any modified version of the thermo for CH2O2 in the paper. It says "The data for HOCHO are well established [29]" which makes me think they just used the Argon/Technion values as is.

Can we just replace the current RMG entry with the original ANL/technion version plotted below?

[sevyharris]

Here's the code to make the plot/get the discontinuous thermo:

import rmgpy.chemkin
import rmgpy
import rmgpy.kinetics
import numpy as np

import matplotlib.pyplot as plt
%matplotlib inline

database = rmgpy.data.rmg.RMGDatabase()

thermo_libraries = [
    'Klippenstein_Glarborg2016',
    'surfaceThermoPt111'
]

database.load(
    path = rmgpy.settings['database.directory'],
    thermo_libraries = thermo_libraries,
    transport_libraries = [],
    reaction_libraries = [],
    seed_mechanisms = [],
    kinetics_families = [],
    kinetics_depositories = [],
    depository = False,
)

original_entry = """HOCHO(71)               C   1H   2O   2     G   200.000  6000.000 1000.00      1
0.46138316E+01 0.64496364E-02 -0.22908251E-05 0.36716047E-09-0.21873675E-13 2
-0.47514850E+05 0.84788383E+00 0.38983616E+01-0.35587795E-02 0.35520538E-04 3
-0.43849959E-07 0.17107769E-10-0.46770609E+05 0.73495397E+01-0.45531246E+05 4"""

species, original_thermo, formula = rmgpy.chemkin.read_thermo_entry(original_entry)
ref_sp = rmgpy.species.Species().from_adjacency_list("""
1 O u0 p2 c0 {3,S} {5,S}
2 O u0 p2 c0 {3,D}
3 C u0 p0 c0 {1,S} {2,D} {4,S}
4 H u0 p0 c0 {3,S}
5 H u0 p0 c0 {1,S}
""")
display(ref_sp)

rmg_data = database.thermo.get_thermo_data(ref_sp)
print(rmg_data)


temps = np.linspace(300, 3000, 1001)
H_RMG = np.zeros_like(temps)
H_original = np.zeros_like(temps)
for i in range(len(temps)):
    H_RMG[i] = rmg_data.get_enthalpy(temps[i])
    H_original[i] = original_thermo.get_enthalpy(temps[i])
plt.plot(temps, H_RMG, label='RMG Klippenstein_Glarborg2016')
plt.plot(temps, H_original, label='Original Argon/Technion', linestyle='dashed')
plt.ylabel('Enthalpy (J/mol)')
plt.xlabel('Temperature (K)')
plt.legend()

[sevyharris]

On closer examination, it looks like the only differences between the RMG entry and the original Technion/ANL version are that the low-temperature polynomials are only off by 8.0 J/mol, and then the high-temperature polynomial is mysteriously shifted 200 J/mol compared to the original.

(this is the diff)

  thermo = NASA(
         polynomials = [
-            NASAPolynomial(coeffs=[3.89836,-0.00355878,3.55205e-05,-4.385e-08,1.71078e-11,-46778.6,7.34954], Tmin=(200,'K'), Tmax=(1000,'K')),
-            NASAPolynomial(coeffs=[4.61383,0.00644964,-2.29083e-06,3.6716e-10,-2.18737e-14,-45330.3,0.847884], Tmin=(1000,'K'), Tmax=(6000,'K')),
+            NASAPolynomial(coeffs=[3.89836,-0.00355878,3.55205e-05,-4.385e-08,1.71078e-11,-46770.6,7.34954], Tmin=(200,'K'), Tmax=(1000,'K')),
+            NASAPolynomial(coeffs=[4.61383,0.00644964,-2.29083e-06,3.6716e-10,-2.18737e-14,-47514.8,0.847884]], Tmin=(1000,'K'), Tmax=(6000,'K')),
         ],
         Tmin = (200,'K'),
         Tmax = (6000,'K'),