Skip to content

Commit c1dfc64

Browse files
committed
Use lowest energy Configuration for energy_correction.
As noted in #2791 the energy_correction was not actually an average of the Configuration energies because some Configurations (reactants and products) are bimolecular. This corrects the calculation, but also switches to using the MINIMUM rather than the average, so that the lowest should end up at zero.
1 parent 97beecc commit c1dfc64

File tree

1 file changed

+3
-7
lines changed

1 file changed

+3
-7
lines changed

rmgpy/rmg/pdep.py

Lines changed: 3 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -788,38 +788,34 @@ def update(self, reaction_model, pdep_settings, requires_rms=False):
788788
# Log the network being updated
789789
logging.info("Updating {0!s}".format(self))
790790

791-
E0 = []
792791
# Generate states data for unimolecular isomers and reactants if necessary
793792
for isomer in self.isomers:
794793
spec = isomer.species[0]
795794
if not spec.has_statmech():
796795
spec.generate_statmech()
797-
E0.append(spec.conformer.E0.value_si)
798796
for reactants in self.reactants:
799797
for spec in reactants.species:
800798
if not spec.has_statmech():
801799
spec.generate_statmech()
802-
E0.append(spec.conformer.E0.value_si)
803800
# Also generate states data for any path reaction reactants, so we can
804801
# always apply the ILT method in the direction the kinetics are known
805802
for reaction in self.path_reactions:
806803
for spec in reaction.reactants:
807804
if not spec.has_statmech():
808805
spec.generate_statmech()
809-
E0.append(spec.conformer.E0.value_si)
810806
# While we don't need the frequencies for product channels, we do need
811807
# the E0, so create a conformer object with the E0 for the product
812808
# channel species if necessary
813809
for products in self.products:
814810
for spec in products.species:
815811
if spec.conformer is None:
816812
spec.conformer = Conformer(E0=spec.get_thermo_data().E0)
817-
E0.append(spec.conformer.E0.value_si)
818813

819-
# Use the average E0 as the reference energy (`energy_correction`) for the network
814+
# Use the lowest E0 as the reference energy (`energy_correction`) for the network
820815
# The `energy_correction` will be added to the free energies and enthalpies for each
821816
# configuration in the network.
822-
energy_correction = -np.array(E0).mean()
817+
energy_correction = -min(sum(spec.conformer.E0.value_si for spec in stationary_point.species)
818+
for stationary_point in self.reactants + self.isomers + self.products)
823819
for spec in self.reactants + self.products + self.isomers:
824820
spec.energy_correction = energy_correction
825821
self.energy_correction = energy_correction

0 commit comments

Comments
 (0)