Calculating Electron-Transfer Barriers

1. Introduction

Here we describe the workflow for computing thermal activation barriers for outer-sphere electron transfer (ET) processes using the Marcus–Hush model. The procedure outlines how to obtain the free energy difference (\(\Delta G ^\circ\)) and reorganization energy (\(\lambda\)) from density functional theory (DFT) calculations, including the use of non-equilibrium solvation to calculate the solvent contributions.

2. Terminologies and Definitions

\(\boldsymbol{\Delta G ^\circ}\) (Free Energy Difference): Gibbs free energy difference between reactant (pre-SET) and product (post-SET) states.

\(\boldsymbol{\lambda}\) (Reorganization Energy): Sum of the energetic penalties associated with reorganizing. There are two types of reorganization energy:

\(\lambda_i\) : Internal geometry change
\(\lambda_{\circ}\) : Solvent/environmental response

Activation Barrier (Marcus)

\(\Delta G ^ \ddagger\ = \frac{(\lambda + \Delta G ^ \circ )^ 2}{4\lambda}\)

Non-Equilibrium Solvation A solvation model in which the solvent polarization from one electronic state is used to evaluate the energy of another state (using NonEq=Save / NonEq=Read).

3. Brief Theory

Consider the example electron transfer (ET) between a triplet quinone derivative (acceptor) and benzene (donor).

The \(\Delta G ^\circ\) (free energy difference) between the reactant (before ET) and product (after ET) states is the Gibbs free energy difference between the optimized "bottom of the well" geometries of the structures.

The reorganization energy (\(\lambda\)) can be visualized as the energy penalty to take the reactant structures from their optimal geometry and solvation shell to the geometry and solvation shell of the product, while maintaining their initial electronic configuration.

4. Procedure

The Marcus–Hush calculation consists of three major steps:

Compute \(\Delta G ^\circ\) from optimized reactant and product states.
Compute \(\lambda\) using four single-point calculations (Nelsen's four-point method).
Calculate \(\Delta G ^\ddagger\) using the Marcus equation.

Step 1: Compute Free Energy Difference (\(\Delta G ^\circ\))

Optimize the geometries of the reactant and product states separately.
Extract the Gibbs free energies from each calculation.
Compute:

\(\Delta G ^ \circ\ = G_{product} - G_{reactant}\)

Step 2: Compute Reorganization Energy (\(\lambda\))

Reorganization energy is calculated using four single-point calculations: two on reactant surfaces and two on product surfaces under non-equilibrium solvation.

Required Calculations

Reactant-State Energies

\(\boldsymbol{E_1}\) : Donor (D) at its own optimized geometry and solvation Charge/spin = reactant state Save solvation environment using NonEq=Save
\(\boldsymbol{E_2}\) : Acceptor (A) at its own optimized geometry and solvation Charge/spin = reactant state Save solvation environment using NonEq=Save

Product-State Energies Using Reactant Electronic Configuration

\(\boldsymbol{E_3}\) : Donor (D) at the product-optimized geometry Use reactant charge/multiplicity Read product solvation environment with NonEq=Read
\(\boldsymbol{E_4}\) : Acceptor (A) at the product-optimized geometry Use reactant charge/multiplicity Read product solvation environment with NonEq=Read

Computing \(\lambda\)

\(\lambda = (\boldsymbol{E_3} - \boldsymbol{E_1}) + (\boldsymbol{E_4} - \boldsymbol{E_2})\)

Forward and Reverse \(\lambda\)

Compute \(\lambda\) on both reactant and product surfaces.

If the values differ (which is common), take the geometric mean following Nelsen’s recommendation:

\(\lambda_{eff} = \sqrt{\lambda_{forward} \times \lambda_{reverse}}\)

Step 3: Calculate Activation Barrier (\(\Delta G ^ \ddagger\))

Use the standard Marcus equation:

\(\Delta G ^ \ddagger\ = \frac{(\lambda + \Delta G ^ \circ )^ 2}{4\lambda}\)

This value is used as the ET activation barrier.

5. Example Workflow for \(\lambda\) Calculation

We demonstrate the \(\lambda\) calculation workflow using the example shown above.

Four energies are required:

\(\boldsymbol{E_1}\): Optimized D, charge 0, multiplicity 1
\(\boldsymbol{E_2}\): Optimized 3A, charge 0, multiplicity 3
\(\boldsymbol{E_3}\): Cation D⁺ geometry, but electronic configuration of D (charge 0, multiplicity 1); NonEq=Read (from saved solvation from D⁺)
\(\boldsymbol{E_4}\): Anion A⁻ geometry, but electronic configuration of 3A (charge 0, multiplicity 3); NonEq=Read (from saved solvation from A⁻)

Compute \(\lambda\) via:

\(\lambda = (\boldsymbol{E_3} - \boldsymbol{E_1}) + (\boldsymbol{E_4} - \boldsymbol{E_2})\)

Repeat for product-state ET and take the geometric mean.