# 10.6: Semi-Empirical Methods- Extended Hückel

- Page ID
- 1972

An electronic structure calculation from first principles (*ab initio*) presents a number of challenges. Many integrals must be evaluated followed by a self-consistent process for assessing the electron-electron interaction and then electron correlation effects must be taken into account. Semi-empirical methods do not proceed analytically in addressing these issues, but rather uses experimental data to facilitate the process. Several such methods are available. These methods are illustrated here by the approaches built on the work of Hückel.

One of the first semi-empirical methods to be developed was Hückel Molecular Orbital Theory (HMO). HMO was developed to describe molecules containing conjugated double bonds. HMO considered only electrons in \(\pi\) orbitals and ignored all other electrons in a molecule. It was successful because it could address a number of issues associated with a large group of molecules at a time when calculations were done on mechanical calculators.

The Extended Hückel Molecular Orbital Method (EH) grew out of the need to consider all valence electrons in a molecular orbital calculation. By considering all valence electrons, chemists could determine molecular structure, compute energy barriers for rotation about bonds, and even determine energies and structures of transition states for reactions. The computed energies could be used to choose between proposed transitions states to clarify reaction mechanisms.

In the EH method, only the \(n\) valence electrons are considered. The total valence electron wavefunction is described as a product of the one-electron wavefunctions.

\[\Psi _{valence} = \psi _1(1) \psi _2(2) \psi _3(3) \psi _3(3) \dots \psi _j(n) \label {10.34}\]

where \(n\) is the number of electrons and j identifies the molecular orbital. Each molecular orbital is written as an linear combination of atomic orbitals (LCAO).

\[\psi _j = \sum \limits ^N_{r = 1} c_{jr} \varphi j = 1, 2, \dots N \label {10.35}\]

where now the \(\varphi _j\) are the valance atomic orbitals chosen to include the 2s, 2p_{x}, 2p_{y}, and 2p_{z} of the carbons and heteroatoms in the molecule and the 1s orbitals of the hydrogen atoms. These orbitals form the *basis set*. Since this basis set contains only the atomic-like orbitals for the valence shell of the atoms in a molecule, it is called a *minimal basis set*.

Each \(\psi _j\), with \(j = 1…N\), represents a molecular orbital, i.e. a wavefunction for one electron moving in the electrostatic field of the nuclei and the other electrons. Two electrons with different spins are placed in each molecular orbital so that the number of occupied molecular orbitals N is half the number of electrons, n, i.e. N = n/2.

The number of molecular orbitals that one obtains by this procedure is equal to the number of atomic orbitals. Consequently, the indices j and r both run from 1 to N. The c_{jr} are the weighting coefficients for the atomic orbitals in the molecular orbital. These coefficients are not necessarily equal, or in other words, the orbital on each atom is not used to the same extent to form each molecular orbital. Different values for the coefficients give rise to different net charges at different positions in a molecule. This charge distribution is very important when discussing spectroscopy and chemical reactivity.

The energy of the j^{th} molecular orbital is given by a one-electron Schrödinger equation using an effective one electron Hamiltonian, h_{eff}, which expresses the interaction of an electron with the rest of the molecule.

\[h_{eff} \psi _j = \epsilon _j \psi _j \label {10.36}\]

is the energy eigenvalue of the j^{th} molecular orbital, corresponding to the eigenfunction \(\psi _j\). The beauty of this method, as we will see later, is that the exact form of h_{eff} is not needed. The total energy of the molecule is the sum of the single electron energies.

\[E_{\pi} = \sum \limits _{j} n_j \epsilon _j \label {10.37}\]

where \(n_j\) is the number of electrons in orbital \(j\).

The expectation value expression for the energy for each molecular orbital is used to find and then \(E_{\pi}\)

\[\epsilon _j = \dfrac {\int \psi _j \times h_{eff} \psi _j d\tau}{\int \psi _j \times \psi _j d\tau} = \dfrac {\left \langle \psi _j | h_{eff} | \psi _j \right \rangle}{\left \langle \psi _j | \psi _J \right \rangle} \label {10.38}\]

The notation \(\left \langle | | \right \rangle \), which is called a bra-ket, just simplifies writing the expression for the integral. Note that the complex conjugate now is identified by the left-side position and the bra notation \( < | \) and not by an explicit *.

After substituting Equation \(\ref{10.35}\) into \(\ref{10.38}\), we obtain for each molecular orbital

\[ \epsilon _j = \dfrac {\left \langle \sum \limits ^N_{r = 1} c_{jr}\psi _r | h_{eff} | \sum \limits ^N_{s = 1} c_{js} \psi _s\right \rangle}{\left \langle \sum \limits ^N_{r = 1} c_{jr}\psi _r | \sum \limits ^N_{s = 1} c_{js}\psi _s \right \rangle} \label {10.39}\]

which can be rewritten as

\[\epsilon = \dfrac {\sum \limits ^N_{r=1} \sum \limits ^N_{s=1} c^*_r c_s \left \langle \psi _r |h_{eff}| \psi _s \right \rangle}{\sum \limits ^N_{r=1} \sum \limits ^N_{s=1} c^*_r c_s \left \langle \psi _r | \psi _s \right \rangle} \label {10.40}\]

where the index j for the molecular orbital has been dropped because this equation applies to any of the molecular orbitals.

To simplify the notation we use the following definitions. The integrals in the denominator of Equation \(\ref{10.40}\) represent the overlap between two atomic orbitals used in the linear combination. The overlap integral is written as \(S_{rs}\). The integrals in the numerator of Equation \(\ref{10.40}\) are called either resonance integrals or coulomb integrals depending on the atomic orbitals on either side of the operator h_{eff} as described below.

- \(S_{Rs} = \left \langle \psi _r |\psi _s \right \rangle\) is the overlap integral. \(S_{rr} = 1\) because we use normalized atomic orbitals. For atomic orbitals r and s on different atoms, \(S_{rs}\) has some value between 1 and 0: the further apart the two atoms, the smaller the value of \(S_{rs}\).
- \(H_{rr} = \left \langle \psi _r |h_{eff}| \psi _s \right \rangle\) is the
*Coulomb Integral*. It is the kinetic and potential energy of an electron in, or described by, an atomic orbital, \(\varphi _r\), experiencing the electrostatic interactions with all the other electrons and all the positive nuclei. - \(H_{rs} = \left \langle \psi _r |h_{eff} |\psi _s\right \rangle\) is the
*Resonance Integral or Bond Integral*. This integral gives the energy of an electron in the region of space where the functions \(\varphi _r\) and \(\varphi _s\) overlap. This energy sometimes is referred to as the energy of the overlap charge. If \(r\) and \(s\) are on adjacent bonded atoms, this integral has a finite value. If the atoms are not adjacent, the value is smaller, and assumed to be zero in the Hückel model.

In terms of this notation, Equation \(\ref{10.49}\) can be written as

\[\epsilon = \dfrac {\sum ^N_{r=1} \sum ^N_{s=1} c ^*_r c_s H_{rs}}{\sum ^N_{r=1} \sum ^N_{s=1} c ^*_r c_s S_{rs}} \label {10.41}\]

We now must find the coefficients, the c's. One must have a criterion for finding the coefficients. The criterion used is the *Variational Principle*. Since the energy depends linearly on the coefficients in Equation \(\ref{10.41}\), the method we use to find the best set of coefficients is called the Linear Variational Method.

The task is to minimize the energy with respect to all the coefficients by solving the N simultaneous equations produced by differentiating Equation \(\ref{10.41}\) with respect to each coefficient.

\[\dfrac {\partial \epsilon}{\partial c_t} = 0 \label {10.42} \]

for \(t = 1, 2, 3, \dots N\)

Actually we also should differentiate Equation \(\ref{10.41}\) with respect to the \(c^*_t\), but this second set of N equations is just the complex conjugate of the first and produces no new information or constants.

To carry out this task, rewrite Equation \(\ref{10.41}\) to obtain Equation \(\ref{10.43}\) and then take the derivative of Equation \(\ref{10.43}\) with respect to each of the coefficients.

\[\epsilon \sum \limits _r \sum \limits _s c^*_r c_s S_{rs} = \sum \limits _r \sum \limits _s c^*_r c_s H_{rs} \label {10.43}\]

Actually we do not want to do this differentiation N times, so consider the general case where the coefficient is. Here t represents any number between 1 and N.

This differentiation is relatively easy, and the result, which is shown by Equation \(\ref{10.44}\), is relatively simple because some terms in Equation \(\ref{10.43}\) do not involve and others depend linearly on. The derivative of the terms that do not involve c_{t} is zero (e.g.

\[\dfrac {\partial c^*_3 c_4 H_{34}}{\partial c_2} = 0. \]

The derivative of terms that contain is just the constant factor that multiples the, (e.g. \(\dfrac {\partial c^*_3 c_2 H_{32}}{\partial c_2} = c^*_3 H_{32}\) ). Consequently, only terms in Equation \(\ref{10.43}\) that contain contribute to the result, and whenever a term contains, that term appears in Equation \(\ref{10.44}\) without the because we are differentiating with respect to. The result after differentiating is

\[\epsilon \sum \limits _r c^*_r S_{rt} = \sum \limits _r c^*_r H_{rt} \label {10.44}\]

If we take the complex conjugate of both sides, we obtain

Since

\[\epsilon = \epsilon ^*, S^*_{rt} = S_{tr} \]

and

\[ H^*_{rt} = H_{tr},\]

then Equation \(\ref{10.45}\) can be reversed and written as

\[\sum \limits _r c_r H_{tr} = \epsilon \sum \limits _r c_r S_{tr} \label {10.46}\]

or upon rearranging as

\[\sum \limits _r c_r (H_{tr} - S_{tr}\epsilon ) = 0 \label {10.47}\]

There are \(N\) simultaneous equations that look like this general one; N is the number of coefficients in the LCAO. Each equation is obtained by differentiating Equation \(\ref{10.43}\) with respect to one of the coefficients.

This method is called the **linear variational method **because the variable parameters affect the energy linearly unlike the shielding parameter in the wavefunction that was discussed in Chapter 9. The shielding parameter appears in the exponential part of the wavefunction and the effect on the energy is nonlinear. A nonlinear variational calculation is more laborious than a linear variational calculation.

Equations \(\ref{10.46}\) and \(\ref{10.47}\) represent a set of homogeneous linear equations. As we discussed for the case of **normal mode analysis** in Chapter 6, a number of methods can be used for solving these equations to obtain values for the energies, \(\epsilon ' s\), and the coefficients, the \(c'_r s\).

Matrix methods are the most convenient and powerful. First we write more explicitly the set of simultaneous equations that is represented by Equation . The first equation has t = 1, the second t = 2, *etc*. N represents the index of the last atomic orbital in the linear combination.

\[\begin{align} c_1H_{11} + c_2H_{12} + \dots c_nH_{1N} &= c_1S_{11}\epsilon +c_2S_{12}\epsilon + \dots c_NS_{1N}\epsilon \nonumber \\[6pt] c_1H_{21} + c_2H_{22} + \dots c_nH_{2N} &= c_1S_{21}\epsilon +c_2S_{22}\epsilon + \dots c_NS_{2N}\epsilon \nonumber \\[6pt] \vdots \vdots &= \vdots \vdots \nonumber \\[6pt] c_1H_{N1} + c_2H_{N22} + \dots c_nH_{NN} &= c_1S_{N1}\epsilon +c_2S_{N2}\epsilon + \dots c_NS_{NN}\epsilon \nonumber \end{align} \label {10.48}\]

This set of equations can be represented in matrix notation.

\[HC' = SC' \epsilon \label {10.49}\]

Here we have square matrix **H** and **S** multiplying a column vector **C'** and a scalar \(\epsilon\). Rearranging produces

\[HC' - SC' \epsilon = 0\]

\[ (H - S\epsilon )C' = 0 \label {10.50}\]

The problem is to solve these simultaneous equations, or the matrix equation, and find the orbital energies, which are the \(\epsilon ' s\), and the atomic orbital coefficients, the \(c's\), that define the molecular orbitals.

In the EH method we use an effective one electron Hamiltonian, and then proceed to determine the energy of a molecular orbital where \(H_{rs} = \left \langle \psi _r |h_{eff} |\psi _s\right \rangle\) and \(S_{rs} = \left \langle \psi _r |\psi _s\right \rangle\) .

Minimization of the energy with respect to each of the coefficients again yields a set of simultaneous equations just like Equation \(\ref{10.47}\).

\[\sum \limits _r c_r (H_{tr} - S_{tr}\epsilon) =0 \label {10.52} \]

As before, these equations can be written in matrix form

\[HC' = SC'\epsilon \label {10.49}\]

Equation \(\ref{10.49}\) accounts for one molecular orbital. It has energy \(\epsilon \), and it is defined by the elements in the **C' **column vector, which are the coefficients that multiply the atomic orbital basis functions in the linear combination of atomic orbitals.

We can write one matrix equation for all the molecular orbitals.

\[HC = SCE \label {10.53}\]

where **H** is a square matrix containing the H_{rs}, the one electron energy integrals, and **C** is the matrix of coefficients for the atomic orbitals. Each column in **C** is the **C' **that defines one molecular orbital in terms of the basis functions. In extended Hückel theory, the overlap is not neglected, and **S** is the matrix of overlap integrals. **E** is the diagonal matrix of orbital energies. All of these are square matrices with a size that equals the number of atomic orbitals used in the LCAO for the molecule under consideration.

Equation \(\ref{10.53}\) represents an eigenvalue problem. For any extended Hückel calculation, we need to set up these matrices and then find the eigenvalues and eigenvectors. The eigenvalues are the orbital energies, and the eigenvectors are the atomic orbital coefficients that define the molecular orbital in terms of the basis functions.

The elements of the **H** matrix are assigned using experimental data. This approach makes the extended Hückel method a *semi-empirical* molecular orbital method. The basic structure of the method is based on the principles of physics and mathematics while the values of certain integrals are assigned by using educated guesses and experimental data. The H_{rr} are chosen as valence state ionization potentials with a minus sign to indicate binding. The values used by R. Hoffmann when he developed the extended Hückel technique were those of H.A. Skinner and H.O. Pritchard (Trans. Faraday Soc. **49 **(1953), 1254). These values for C and H are listed in Table 10.1. The values for the heteroatoms (N, O, and F) are taken from Pople and Beveridge(*Approximate Molecular Orbital Theory*, McGraw-Hill Book Company, New York, 1970).

Atomic orbital | Ionization potential (eV) |
---|---|

H 1s | 13.6 |

C 2s | 21.4 |

C 2p | 11.4 |

N 2s | 25.58 |

N 2p | 13.9 |

O 2s | 32.38 |

O 2p | 15.85 |

F 2s | 40.20 |

F 2p | 18.66 |

The H_{rs} values are computed from the ionization potentials according to

\[H_{rs} = \dfrac {1}{2} K (H_{rr} + H_{ss})S_{rs} \label {10.54}\]

The rationale for this expression is that the energy should be proportional to the energy of the atomic orbitals, and should be greater when the overlap of the atomic orbitals is greater. The contribution of these effects to the energy is scaled by the parameter K. Hoffmann assigned the value of K after a study of the effect of this parameter on the energies of the occupied orbitals of ethane. The conclusion was that a good value for K is K = 1.75.

The overlap matrix also must be determined. The matrix elements are computed using the definition \(S_{rs} = \left \langle \psi _r |\psi _s\right \rangle\) where \(\varphi _k\) and \(\psi _s\) are the atomic orbitals. *Slater-type orbitals *(STO’s) are used for the atomic orbitals rather than hydrogenic orbitals because integrals involving STO's can be computed more quickly on computers. Slater type orbitals have the form

\[\phi _{1s} (r) = 2\zeta ^{3/2} \text {exp} (- \zeta r)\]

\[\phi _{2s} (r) = \phi _2p (r) = \left (\dfrac {4\zeta ^5}{3} \right )^{1/2} \text {rexp} (- \zeta r) \label {10.55}\]

where zeta, \(\zeta\), is a parameter describing the screened nuclear charge. In the extended Hückel calculations done by Hoffmann, the Slater orbital parameter \(\zeta\) was 1.0 for the H_{1s} and 1.652 for the C_{2s} and C_{2p} orbitals.

Overlap integrals involve two orbitals on two different atoms or centers. Such integrals are called two-center integrals. In such integrals there are two variables to consider, corresponding to the distances from each of the atomic centers, r_{A} and r_{B}. Such integrals can be represented as

\[S_{A_{2s}B_{2s}} = \left (\dfrac {4\zeta ^5}{3}\right ) \int r_A \text {exp} (- \zeta r_A) r_B \text {exp} (- \zeta r_B) d\tau \label {10.56}\]

but elliptical coordinates must be used for the actual integration. Fortunately the software that does extended Hückel calculations contains the programming code to do overlap integrals. The interested reader will find sufficient detail on the evaluation of overlap integrals and the creation of the programmable mathematical form for any pair of Slater orbitals in Appendix B4 (pp. 199 - 200) of the book *Approximate Molecular Orbital Theory* by Pople and Beveridge. The values of the overlap integrals for HF are given in Table 10.2.

Our goal is to find the coefficients in the linear combinations of atomic orbitals and the energies of the molecular orbitals. For these results, we need to transform Equation \(\ref{10.53}\)

\[HC = SCE \label {10.53}\]

into a form that allows us to use matrix diagonalization techniques. We are hampered here by the fact that the overlap matrix is not diagonal because the orbitals are not orthogonal. Mathematical methods do exist that can be used to transform a set of functions into an orthogonal set. Essentially these methods apply a transformation of the coordinates from the local coordinate system describing the molecule into one where the atomic orbitals in the LCAO are all orthogonal. Such a transformation can be accomplished through matrix algebra, and computer algorithms for this procedure are part of all molecular orbital programs. The following paragraph describes how this transformation can be accomplished.

If the matrix \(M\) has an inverse \(M^{-1}\) Then

\[MM^{-1} = 1tag {10.57}\]

and we can place this product in a matrix equation without changing the equation. When this is done for Equation \(\ref{10.53}\), we obtain

\[HMM^{-1}C = SMM^{-1} CE \label {10.58}\]

Next multiply on the left by \(M^{-1}\) and determine \(M\) so the product \(M^{-1}SM\) is the identity matrix, i.e. a matrix that has 1's on the diagonal and 0's off the diagonal is the case for an orthogonal basis set.

\[ M^{-1}HMM^{-1}C = M^{-1}SMM^{-1}CE \label {10.59}\]

which then can be written as

\[H''C'' = C''E'' \label {10.60}\]

where

\[C' = M^{-1}C \label {10.61}\]

The identity matrix is not included because multiplying by the identity matrix is just like multiplying by the number 1. It doesn’t change anything. The \(H''\) matrix can be diagonalized by multiplying on the left by the inverse of \(C''\) to find the energies of the molecular orbitals in the resulting diagonal matrix \(E\).

\[E = C''^{-1}H''C'' \label {10.62}\]

The matrix \(C''\) obtained in the diagonalization step is finally back transformed to the original coordinate system with the \(M\) matrix, \(C = MC''\) since \(C'' = M^{-1}C\).

Fortunately this process is automated in some computer software. For example, in Mathcad, the command genvals(H,S) returns a list of the eigenvalues for Equation \(\ref{10.53}\). These eigenvalues are the diagonal elements of \(E\). The command genvecs(H,S) returns a matrix of the normalized eigenvectors corresponding to the eigenvalues. The i^{th} eigenvalue in the list goes with the i^{th} column in the eigenvector matrix. This problem, where \(S\) is not the identity matrix, is called a general eigenvalue problem, and *gen *in the Mathcad commands refers to *general*.

David M. Hanson, Erica Harvey, Robert Sweeney, Theresa Julia Zielinski ("Quantum States of Atoms and Molecules")