Bonding or Antibonding?

In this example we show how to differentiate between a bonding and antibonding MO for two chosen fragments. As a system we chose the homolytic bond-cleavage of the carbon-chloride bond in chloroethane calculated at the ZORA-OLYP/TZ2P level of theory.

In general, to know if an MO has bonding or antibonding character for two chosen fragments we calculate the bond order

\[\textrm{BO}(\Psi_i) = 2N_i \sum_{\psi_j^A \in \{\psi^A\}} \sum_{\psi_k^B \in \{\psi^B\}} \langle \psi_j^A | \Psi_i \rangle \langle \psi_k^B | \Psi_i \rangle \langle \psi_j^A | \psi_k^B \rangle,\]

where \(\Psi_i\) is the MO of interest with occupation \(N_i\) and \(\{\psi^A\}\) and \(\{\psi^B\}\) are the sets of SFOs from fragment A and B respectively. \(\langle \psi_j^A | \Psi_i \rangle\) is the projection coefficient of \(\psi_j^A\) onto the MO and \(\langle \psi_j^A | \psi_k^B \rangle\) is the overlap between the SFOs. A negative value for \(\textrm{BO}(\Psi_i)\) is antibonding character and a positive value bonding character. The magnitude of \(\textrm{BO}(\Psi_i)\) indicates how strong the bond or antibond is. If the value is close to zero we are likely dealing with a non-bonding combination.

Calculating the bond order for the MOs of chloroethane yields the following bonding and antibonding MOs that have a bond order strength of 0.001 or greater:

BondingAntibonding

MO 4A	MO 8A	MO 11A
MO 12A	MO 13A	MO 14A

Download bonding.adf.rkf

Download bonding.py

bonding.pyOutput

'''
This example shows how to find bonding and anti-bonding orbitals.

This example uses an ADF calculation on chloroethane
calculated at the OLYP/TZ2P level of theory.
The C-Cl bond was homolytically cleaved.
'''

import pyorbb


def bond_order(mo, fmos1, fmos2):
	'''
	This function calculates the bond order for an MO between two sets of FMOS.
	The bonding order is the sum of all overlap populations between the two
	sets of FMOs.
	'''
	total = 0
	for fmo1 in fmos1:
		c1 = fmo1.coefficient(mo)
		for fmo2 in fmos2:
			c2 = fmo2.coefficient(mo)
			S = fmo1.overlap(fmo2)
			# this is the overlap population between fmo1 and fmo2
			total += 2 * mo.occupation * c1 * c2 * S

	return total

# load the orbital data
orbs = pyorbb.Orbitals('bonding.adf.rkf')

# split FMOS between the two fragments
fmos_fragment1 = orbs.fmos.filter(fragment='LeavingGroup')
fmos_fragment2 = orbs.fmos.filter(fragment='Substrate')

# to simplify we only consider one spin type for the MOs
mos = orbs.mos.filter(spin='A')

rows = []
total_bond_order = 0
# go through each MO and calculate the score
for mo in mos:
	score = bond_order(mo, fmos_fragment1, fmos_fragment2)
	total_bond_order += score

	# if the score is too low we don't show it
	if abs(score) < 1e-3:
		continue

	# if the score is high enough we write a new string to be printed later
	rows.append([mo, score > 0, score])

# add a total row
rows.append(['Total', total_bond_order > 0, total_bond_order])

# print the gathered data
print('   MO  Bonding?  Score')
print('──────────────────────')
for mo, bonding, score in rows[:-1]:
	print(f'{str(mo):>5}  {str(bonding):8} {score: .3f}')
print('──────────────────────')
print(f'{str(rows[-1][0]):>5}  {str(rows[-1][1]):8} {rows[-1][2]: .3f}')