proccessing surface charge transfer training reactions

rmgpy/data/kinetics/family.py

@@ -1233,6 +1233,23 @@ def add_rules_from_training(self, thermo_database=None, train_indices=None):
                     Tmin=deepcopy(data.Tmin),
                     Tmax=deepcopy(data.Tmax)
                 )
+            elif isinstance(data, SurfaceChargeTransfer):
+                for reactant in entry.item.reactants:
+                    # Clear atom labels to avoid effects on thermo generation, ok because this is a deepcopy
+                    reactant_copy = reactant.copy(deep=True)
+                    reactant_copy.molecule[0].clear_labeled_atoms()
+                    reactant_copy.generate_resonance_structures()
+                    reactant.thermo = thermo_database.get_thermo_data(reactant_copy, training_set=True)
+                for product in entry.item.products:
+                    product_copy = product.copy(deep=True)
+                    product_copy.molecule[0].clear_labeled_atoms()
+                    product_copy.generate_resonance_structures()
+                    product.thermo = thermo_database.get_thermo_data(product_copy, training_set=True)
+                entry.item.kinetics = data
+                V0 = entry.item.get_reversible_potential(298)
+                data.change_v0(V0)
+                data.Ea = (data.Ea.value_si / 1000, 'kJ/mol')
+                data.V0 = None
             else:
                 raise NotImplementedError("Unexpected training kinetics type {} for {}".format(type(data), entry))
 
@@ -1262,7 +1279,9 @@ def add_rules_from_training(self, thermo_database=None, train_indices=None):
 
             tentries[entry.index].item.is_forward = False
 
-            assert isinstance(entry.data, Arrhenius)
+            if not isinstance(entry.data, Arrhenius):
+                print(self.label)
+                assert False
             data = deepcopy(entry.data)
             data.change_t0(1)
             # Estimate the thermo for the reactants and products
@@ -1291,6 +1310,11 @@ def add_rules_from_training(self, thermo_database=None, train_indices=None):
                 product.thermo = thermo_database.get_thermo_data(product, training_set=True)
             # Now that we have the thermo, we can get the reverse k(T)
             item.kinetics = data
+            if isinstance(item.kinetics, SurfaceChargeTransfer):
+                V0 = item.get_reversible_potential(298)
+                data.change_v0(V0)
+                data.V0 = None
+
             data = item.generate_reverse_rate_coefficient()
 
             item = TemplateReaction(reactants=[m.molecule[0].copy(deep=True) for m in entry.item.products],

rmgpy/data/kinetics/rules.py

@@ -564,12 +564,18 @@ def _get_average_kinetics(self, kinetics_list):
         n = 0.0
         E0 = 0.0
         alpha = 0.0
+        ne = None
         count = len(kinetics_list)
         for kinetics in kinetics_list:
             logA += math.log10(kinetics.A.value_si)
             n += kinetics.n.value_si
             alpha += kinetics.alpha.value_si
-            E0 += kinetics.E0.value_si
+            if isinstance(kinetics, SurfaceChargeTransfer):
+                E0 += kinetics.Ea.value_si
+                if ne is None:
+                    ne = kinetics.ne.value_si
+            else:
+                E0 += kinetics.E0.value_si
         logA /= count
         n /= count
         alpha /= count
@@ -594,15 +600,24 @@ def _get_average_kinetics(self, kinetics_list):
         else:
             raise Exception('Invalid units {0} for averaging kinetics.'.format(Aunits))
 
-        if type(kinetics) not in [ArrheniusEP, SurfaceArrheniusBEP, StickingCoefficientBEP]:
+        if type(kinetics) not in [ArrheniusEP, SurfaceArrheniusBEP, StickingCoefficientBEP, SurfaceChargeTransfer]:
             raise Exception('Invalid kinetics type {0!r} for {1!r}.'.format(type(kinetics), self))
 
-        averaged_kinetics = type(kinetics)(
-            A=(10 ** logA, Aunits),
-            n=n,
-            alpha=alpha,
-            E0=(E0 * 0.001, "kJ/mol"),
-        )
+        if isinstance(kinetics, SurfaceChargeTransfer):
+            averaged_kinetics = SurfaceChargeTransfer(
+                A=(10 ** logA, Aunits),
+                n=n,
+                ne=ne,
+                alpha=alpha,
+                Ea=(E0 * 0.001, "kJ/mol"),
+                )
+        else:
+            averaged_kinetics = type(kinetics)(
+                A=(10 ** logA, Aunits),
+                n=n,
+                alpha=alpha,
+                E0=(E0 * 0.001, "kJ/mol"),
+            )
         return averaged_kinetics
 
     def estimate_kinetics(self, template, degeneracy=1):