3. Creating the Reaction Model

The ReactionModel stores information about each elementary reaction step.

Properties

• initial (list[str]): Initial configuration (e.g., ['1 CO* 1', '2 * 1']).

• final (list[str]): Final configuration (e.g., ['1 C* 1', '2 O* 1']).

• activ_eng (float): Activation energy (eV).

• vib_energies_is (list[float]): Vibrational energies for the initial state (meV).

• vib_energies_fs (list[float]): Vibrational energies for the final state (meV).

• vib_energies_ts (list[float], optional): Vibrational energies for the transition state (meV).
  Omit or set to [] for non-activated steps.

• area_site (float, optional): Site area (Ų); required if a gas species participates.

• molecule_is (str, optional): Name of gas-phase molecule in the initial state (only for adsorption).

• molecule_fs (str, optional): Name of gas-phase molecule in the final state (only for desorption).

• site_types (str, optional): Site types required if lattice_type='periodic_cell').

• neighboring (str, optional): Connectivity between sites (e.g., '1-2'; default = None).

• prox_factor (float, optional): Proximity factor.

• angles (str, optional): Angle constraints (e.g., '1-2-3:180'; default = None).

• graph_multiplicity (int or float, optional): Symmetry factor of the step. The computed pre-exponential factor is divided by this value.

• fixed_pre_expon (float, optional): Fixed forward pre-exponential factor (no scaling or symmetry applied).
  Units must follow Zacros conventions:

  • Surface/non-activated desorption → s⁻¹

  • Adsorption (activated/non-activated) → bar⁻¹·s⁻¹

• fixed_pe_ratio (float, optional): Fixed pre-exponential ratio pe_fwd/pe_rev. Must be used with fixed_pre_expon.

Deprecation note: The old molecule column is deprecated and treated as molecule_is with a warning.

How to create it

You can provide a dictionary where each key is a cluster name and each value is a dictionary of elementary step properties:

from zacrostools.reaction_model import ReactionModel

steps_data = {
    'CO_adsorption': {
        'activ_eng': 0.0,
        'area_site': 5.34,
        'initial': ['1 * 1'],
        'final': ['1 CO* 1'],
        'molecule_is': 'CO',
        'prox_factor': 0.0,
        'site_types': 'top',
        'vib_energies_is': [263.427],
        'vib_energies_fs': [240.497, 82.738, 60.132, 60.080, 7.271, 6.553],
    },
    'O2_adsorption': {
        'activ_eng': 0.742,
        'area_site': 5.34,
        'initial': ['1 * 1, 2 * 1'],
        'final': ['1 O* 1, 2 O* 1'],
        'neighboring': '1-2',
        'molecule_is': 'O2',
        'prox_factor': 0.5,
        'site_types': 'top',
        'vib_energies_is': [263.427],
        'vib_energies_ts': [244.692, 81.168, 53.342, 30.084, 2.234],
        'vib_energies_fs': [240.497, 82.738, 60.132, 60.080, 7.271, 6.553],
    }, 
    'CO+O_to_CO2': {
        'activ_eng': 1.693,
        'initial': ['1 CO* 1', '2 O* 1'],
        'final': ['1 CO2* 1', '2 * 1'],
        'neighboring': '1-2',
        'site_types': 'top top',
        'vib_energies_is': [171.188, 145.669, 96.964, 86.255, 56.201, 52.376, 35.933, 24.343, 21.025],
        'vib_energies_ts': [217.941, 81.362, 66.833, 56.918, 50.342, 37.430, 19.074, 12.356],
        'vib_energies_fs': [240.497, 82.738, 60.133, 60.080, 7.272, 6.553, 78.662, 40.796, 40.349],
    },
    'CO2_adsorption': {  # recommended to define the forward direction as the adsorption
        'activ_eng': 0.0,
        'area_site': 5.34,
        'initial': ['1 * 1'],
        'final': ['1 CO2* 1'],
        'molecule_is': 'CO2',
        'prox_factor': 0.0,
        'site_types': 'top',
        'vib_energies_is': [293.279, 163.686, 78.019, 77.959],
        'vib_energies_fs': [171.188, 145.668, 96.961, 86.214, 56.201, 52.375, 35.932, 24.343, 21.022],
    }
}

reaction_model = ReactionModel.from_dict(steps_data)

Alternatively, it can also be created from a CSV file. In this case, the indexes must correspond to the elementary step names:

from zacrostools.reaction_model import ReactionModel

reaction_model = ReactionModel.from_csv('mechanism_data.csv')

Finally, it can also be created from a Pandas dataframe:

import pandas as pd
from zacrostools.reaction_model import ReactionModel

df = pd.read_csv("reactions_data.csv")  # or create dataframe directly
reaction_model = ReactionModel.from_df(df)

---

Adding and removing steps

Use the add_step method to add a new step:

reaction_model.add_step(step_info={
        'step_name': 'CO_diffusion',
        'activ_eng': 1.156,
        'initial': ['1 CO* 1', '2 * 1'],
        'final': ['1 * 1', '2 CO* 1'],
        'neighboring': '1-2',
        'graph_multiplicity': 2.0,
        'site_types': 'top top',
        'vib_energies_fs': [240.497, 82.738, 60.133, 60.080, 7.272, 6.553],
        'vib_energies_is': [240.497, 82.738, 60.133, 60.080, 7.272, 6.553],
        'vib_energies_ts': [218.382, 53.527, 47.612, 28.580, 6.600]
    })

Use the remove_steps method to remove a step:

steps_to_remove = ['O_diffusion']
reaction_model.remove_steps(steps_to_remove)

---

Writing the mechanism_input.dat file

The ReactionModel can generate the mechanism_input.dat file required by Zacros.

Required parameters:

• output_dir (str or Path): Directory where the file will be written.

• temperature (float): Temperature in Kelvin for pre-exponential calculations.

• gas_model (GasModel): Instance of GasModel containing gas-phase species data.

• manual_scaling (dict, optional): Dictionary of manual scaling factors per step.

• stiffness_scalable_steps (list, optional): List of steps that are stiffness scalable.

• stiffness_scalable_symmetric_steps (list, optional): List of steps that are stiffness scalable and symmetric.

• sig_figs_energies (int, optional): Number of significant figures for activation energies.

• sig_figs_pe (int, optional): Number of significant figures for pre-exponential factors.

from zacrostools.gas_model import GasModel

gas_model = GasModel.from_csv('gas_data.csv')
reaction_model.write_mechanism_input(
    output_dir='kmc_simulation',
    temperature=500,
    gas_model=gas_model,
    sig_figs_energies=8,
    sig_figs_pe=8
)

---

Accessing reaction data

The reaction steps data is stored internally as a Pandas DataFrame, accessible via the df attribute.

print(reaction_model.df)

---

Full example

Below is a complete example demonstrating the creation and modification of a ReactionModel:

from zacrostools.reaction_model import ReactionModel
from zacrostools.gas_model import GasModel

# Assume gas_model is already defined
gas_model = GasModel.from_csv('gas_data.csv')

# Initial reaction steps data
steps_data = {
    'CO_adsorption': {
        'activ_eng': 0.0,
        'area_site': 5.34,
        'initial': ['1 * 1'],
        'final': ['1 CO* 1'],
        'molecule_is': 'CO',
        'prox_factor': 0.0,
        'site_types': 'top',
        'vib_energies_is': [263.427],
        'vib_energies_fs': [240.497, 82.738, 60.132, 60.080, 7.271, 6.553],
    }
}

# Create the ReactionModel instance
reaction_model = ReactionModel.from_dict(steps_data)

# Add a new step
reaction_model.add_step(step_info={
        'step_name': 'CO_diffusion',
        'activ_eng': 1.156,
        'initial': ['1 CO* 1', '2 * 1'],
        'final': ['1 * 1', '2 CO* 1'],
        'neighboring': '1-2',
        'graph_multiplicity': 2.0,
        'site_types': 'top top',
        'vib_energies_is': [240.497, 82.738, 60.133, 60.080, 7.272, 6.553],
        'vib_energies_ts': [218.382, 53.527, 47.612, 28.580, 6.600],
        'vib_energies_fs': [240.497, 82.738, 60.133, 60.080, 7.272, 6.553],
    })


# Remove an existing step
reaction_model.remove_steps(['CO_diffusion'])

# Write the mechanism input file
reaction_model.write_mechanism_input(
    output_dir='kmc_simulation',
    temperature=500,
    gas_model=gas_model,
    sig_figs_energies=8,
    sig_figs_pe=8
)

# Access the DataFrame
print(reaction_model.df)

---

Fixing pre-exponential factors

By default, ReactionModel computes the forward pre-exponential factor (pre_expon) and the pre-exponential ratio (pe_ratio) from statistical mechanics.
However, in some cases you may want to fix these values explicitly (for both forward and reverse directions). This can be done by providing two optional parameters:

• fixed_pre_expon (float): User-specified forward pre-exponential factor.

• fixed_pe_ratio (float): User-specified pre-exponential ratio.

⚠️ Important:

• If one of these parameters is provided, the other must also be provided.

• Using fixed_pre_expon/fixed_pe_ratio is incompatible with:

  1. Setting stiffness_scalable_steps="all".

  2. Including a fixed step in stiffness_scalable_steps.

  3. Including a fixed step in stiffness_scalable_symmetric_steps.

If any of these conditions occur, a ReactionModelError is raised.

from zacrostools.reaction_model import ReactionModel

reaction_data = {
    'CO_adsorption': {
        'activ_eng': 0.0,
        'area_site': 5.34,
        'initial': ['1 * 1'],
        'final': ['1 CO* 1'],
        'molecule_is': 'CO',
        'prox_factor': 0.0,
        'site_types': 'top',
        'vib_energies_is': [263.427],
        'vib_energies_fs': [240.497, 82.738, 60.132, 60.080, 7.271, 6.553],
        "fixed_pre_expon": 1.0e13,  # in s^-1 or bar^-1·s^-1 depending on type
        "fixed_pe_ratio": 2.5
    },
}

reaction_model = ReactionModel.from_dict(reaction_data)