Determining Electron Affinity of meta-, ortho-, and para-Xylene using DFT

# Introduction

The electron affinity (EA) is a fundamental property of atoms and molecules and is defined as the energy difference between an uncharged species and its negative ion. The negative ion, or anion, results from the addition of an electron from the vacuum to the neutral molecule. [1] There are many areas in the sciences where the ability to absorb an electron and the properties of negative ions or radicals are very important. One such area is in gaseous dielectric materials, such as SF6, in which they can act as electron scavengers and inhibit the possibility of partial discharge (PD). According to the International Electrotechnical Commission, PD is a locally limited discharge that bridges the insulation between the electrodes just partially. [2] Breakdown in insulation materials is initiated by free electrons colliding with other trapped or bound electrons and subsequently produces an avalanche. The ability of electron scavengers to remove these free electrons will ultimately impede or delay the on-set of such an avalanche and improve the breakdown strength of the material. [3]

This paper will try to provide an approach to analyzing an experimental problem encountered by using density functional theory (DFT) to determine the EA of m-, o-, and p-xylene as a function of applied electric field to verify whether this molecule acts as an electron scavenger. A proposed detailed calculation method will be given on determining adiabatic EA of xylene, as well as the possibility of other EA calculation methods. A short discussion will be conducted on the relative accuracy of such calculations, especially in regards to exchange-correlation functional and basis set used.

# Motivation

Our research group at the Electrical Insulation Research Center devised an experiment to study cavity-induced PD within optically clear epoxy. Spherical cavities in the size range of 1.7 to 5 mm were created by injecting a known volume of air into the epoxy with a hypodermic needle on axis in a cylindrical geometry and allowing the system to cure while rotating the sample to keep the air bubble centered. The cured epoxy was then cut to produce samples in the range of 1 cm thick, after which they were placed in between brass electrodes. The applied voltage was ac at 60 Hz, from a 50 kV, 2 kVA transformer. It was assumed that oxygen was consumed during curing of the sample and that nitrogen would be present within the void, which has similar breakdown characteristics to air as can be seen in Figure 1.

Figure 1: Measured and calculated Paschen curves for air, N2, and SF6. [4]

However, no discharges were detected while subjecting the samples to high voltage with electric fields greater than 4 kV/mm, which is sufficient enough for PD to occur in air or nitrogen. As a last attempt a 0.8 mm diameter channel was drilled from the side of the sample into the cavity and flushed with air and the hole was then covered with silicone rubber and allowed to harden. After this procedure, PD was detected within the void. As a result, we came to the conclusion that there must be a gas present within the cavity that was suppressing the PD. Figure 2 shows the chromatogram of the GC/MS analysis of the gas within the void.

Figure 2: GC/MS chromatogram of the gas within the void

As expected, nitrogen was present, yet oddly enough there was also a substantial amount of o-, p-, and m-xylene present. According to the MSDS sheet of the optically clear epoxy, xylene, which encompasses the ortho-, meta-, and para-isomers of dimethylbenzene, is utilized as a thinning agent in order to keep the viscosity of the epoxy resin low enough to ensure optimum mixing with the hardener. The molecular configuration of the three xylene isomers can be seen in Figure 3.

Figure 3: Molecular configuration of ortho-, meta-, and para-xylene

By looking at the xylene molecule, one's chemical intuition would say that it is not electronegative, due to the absence of any electron-withdrawing substituents, and has a negative EA value, meaning the neutral molecule is more stable that the anion. However, nitrogen is considered a pseudohalogen and is definately not electronegative, which would not allow it to suppress discharges. Therefore, we began to ask ourselves a very specific question: under the effect of an electric field do the isomers of xylene act like electron scavengers and suppress PD? Since EA is a fairly difficult parameter to measure experimentally compared to the ionization potential (IP), the use of DFT can be an important and useful tool to determine the EA of certain molecules. [5] I am particularly interested in determining whether the EA of xylene changes under the effect of an electric field.

# Methodology

## Proposed EA Calculation Method for Xylene

The treatment of the EA DFT calculation of xylene will be based upon similar calculations performed by Janet Ho and mentioned in her PhD Thesis. [3] She studied the effect of an electric field on the molecules of Irganox 1010 and 2,6-di-tert-butyl-1,4-benzoquinone by using the local density approximation (LDA) implemented within the SIESTA code. The LDA method treats the electron density as a uniform electron gas. [1] One of the calculation methods within DFT is given in equation 1 and is also known as the adiabatic electron affinity. [5]

(1)
\begin{align} \Delta E=E_{tot}(X)-E_{tot}(X^-) \end{align}

Etot(X) is the total energy of the neutral molecule in its equilibrium geometry
Etot(X-) is the total energy of the anion in its equilibrium geometry

Initially the equilibrium atomic positions of the neutral molecule of all three isomers of xylene will be determined in the absence of an external electric field, so that the forces on each atoms will be at a minimum. Since I am interested in xylene in the gas phase, the calculations will be performed on isolated molecules to reduce any interactions between adjacent molecules. Afterwards, the optimized structure will be given a -1 charge, meaning an electron will be added to the system from the vacuum, and a +1 charge, meaning an electron will be removed to the vacuum, and the equilibrium atomic positions of the two charged states will be determined. At first the EA of the three isomers will be determined using equation 1 in the absence of an electric field. After which the same will be done when a gradually changing electric field will be applied to the molecules under investigation. In Figure 4, p-xylene is used as an example to show possible conjugation of the anion radical after addition of an electron. Stevenson et al. [6] seemed to think that hyperconjugation between the methyl antibonding orbital and the aromatic ring would have a stabilizing effect on the anion radical.

Figure 4: Possible molecular conjugation of p-xylene (1) after addition of an extra electron to the anion radical. Structure (2) depicts the addition of the electron into the LUMO, or the $\pi^\ast$-orbital of the aromatic ring. Structure (3) is one of the valence bond representations given by Stevenson et al. [6]

The total energy of the cation will also be determined in the absence and under the effect of an electric field. The main reason for this calculation method is for comparison purposes, so as to see whether the adiabatic EA of xylene changes in the presence of an electric field and also to see the stabilities of the three different states relative to one another. By calculating the +1 charge state, also known as the IP, one could also determine a further parameter called the Mulliken electronegativity ($\chi$). [7] This parameter is related to both EA and IP and given in equation 2:

(2)
\begin{align} \chi=(EA+IP)/2 \end{align}

In order to also be able to comment on the method's accuracy, I will also calculate the EAs of nitrogen, oxygen, and SF6. These molecules have been rigorously studied before and have well known EAs and electronegativities, which may shed some light on SIESTA's accuracy of determining EA while using the LDA method. Piechota et al. state that LDA can be used to make accurate predictions of the structure and relative energies of covalent as well as ionic systems. However, especially in regard to SF6, they do mention that the use of LDA in conjunction with the Perdew-Wang functional (PWC) yields an EA closer to available experimental data of 1.05 eV. [8]

## Different Methods of Determining EA

### Total Molecular Energies

There are three methods that are particularly useful in determining a molecule's EA and are shown in Figure 5. The proposed calculation method stated above is known as the adiabatic EA, while the other two are known as the vertical detachment energy (VDE) and vertical attachment energy (VAE). VDE is the energy required for the removal of an electron from the anion, whereas the VAE is the energy released from the addition of an electron to a neutral molecule. During each of these processes there is no time for geometry optimization and as a result, the two species are either at the equilibrium geometry of the anion (VDE) or the neutral molecule (VAE). [1] The theoretical definitions are shown mathematically in equations 3 and 4.

(3)
$$VDE=E_{tot}(A)-E_{tot}(X^-)$$
(4)
$$VAE=E_{tot}(X)-E_{tot}(A^-)$$

Etot(X) is the total energy of the neutral molecule in its equilibrium geometry
Etot(X-) is the total energy of the anion in its equilibrium geometry
Etot(A) is the total energy of the neutral molecule at the equilibrium anion geometry
Etot(A-) is the total energy of the anion at the equilibrium neutral geometry

Figure 5: A diagram of potential energy surfaces of an anion (R-) and a neutral molecule (R) representing all diatomic and polyatomic molecules, although not all anions are necessarily energetically more stable. [1]

As can be seen in Figure 5, in addition to electronic effects, molecules also contain vibrational and rotational effects which can contribute to the overall molecular EA. Equations 1, 3, and 4 can be corrected by adding the respective zero-point vibrational energies to the total energies of the neutral and anionic species. However, in most cases the vibrational energies are neglected, as they are usually very similar. [1]

The calculation of VAE as a function of electric field could also be an acceptable method of determining the EA of xylene. However, this process would incorporate the addition of an electron into the equilibrium neutral geometry which, under the effect of an electric field may be unsuitable. Depending on the direction of the applied electric field the addition of an electron from the vacuum may be difficult. If added into the system from the anode side, then electron addition would be easy, whereas the opposite would be difficult due to like charges. Therefore, the adiabatic EA calculation of the equilibrium structures in the absence of an electric field could possibly allow for a more accurate determination of the EA in the presence of an electric field, as the electron has already been added and the factor of the direction of the applied electric field would be eliminated.

### Orbital Energies

De Proft et al. [7] conducted DFT calculations on molecules which had negative EA values, which is very pertinent to xylene as this molecule does exhibit a negative EA. They investigated EAs of neutral, closed-shell molecules according to three different calculation methods based on HOMO and LUMO orbital energies and compared their results to experimental results of electron transmission spectroscopy (ETS). One of their methods which they state as "eqn 6" in Table 1 is based on equation 3 above. The other two methods stated as "eqn 9" and "eqn 11" are given in equations 5 and 6. Table 1 summarizes their results.

(5)
\begin{align} EA \approx -\varepsilon_{LUMO} \end{align}
(6)
\begin{align} EA=-(\varepsilon_{LUMO}+\varepsilon_{HOMO})-IP \end{align}

$\varepsilon_{LUMO}$ is the energy of the LUMO
$\varepsilon_{HOMO}$ is the energy of the HOMO

Table 1: Calculated and experimental VDEs in eV. Calculations performed according to equations 3, 5, and 6 with the PBE functional along with the aug-cc-pVTZ and cc-pVTZ basis sets, respectively. All experimental values were obtained from ETS measurements. [7] Highlighted are the results obtained for m- and o-xylene.

All three of their methods were based upon the assumption that during anion formation, the extra electron enters the LUMO, which allows their results to be directly compared to ETS measurements. Two sets of calculations were done at the PBE/aug-cc-pVTZ and PBE/cc-pVTZ levels. They came to the conclusion that the latter approach in combination with equation 6 yielded the best results. Also, even though results were not given, they state that their methods would be satisfactory for molecules with positive EA values, yet less appropriate for negative EA values due to the fact that the anions are unstable and are strongly resistent towards the uptake of charge. It is interesting to note that from Table 1, the EA of m- and o-xylene were calculated and, as our chemical intuition suggested, revealed negative values. However, I am interested to see whether the EA of xylene will change from a negative to a positive value under the effect of a strong electric field.

## Accuracy of adiabatic EA Calculations using DFT

The computation of EA of atoms and molecules using ab initio methods is an active area of research, especially with respect to the choice of basis set and exchange-correlation functional used. Calculating EA presents challenges because "(1) the correlation energy of anions is usually large and requires sophisticated methods to reproduce accurately, and (2) diffuse functions must usually be included in the basis set to represent the charge density distribution of the anion correctly." [9] DFT is most often defined by its exchange-correlation functional and is usually a combination of an exchange or hybrid exchange functional, which includes the exact exchange energy in Hartree-Fock theory, and a correlation functional. An ab initio method always uses a basis set, which usually consists of atom-centered Gaussian functions. A proper description of the electronic structure of a molecule under investigation requires basis sets with diffuse functions. [1] Naturally, EAs can only be accurately calculated by using the same basis set for the the neutral molecule and anion.

Boesch et al. [9] did extensive research on calculating the gas-phase adiabatic EA of a number of different derivatives of benzoquinone while testing different functionals and basis sets. They tested the exchange functionals of Slater's LDA (S) and Becke's functional (B), combined with different correlation functionals, such as local spin density of Vosko, Wilk, and Nusair (VWN), gradient-corrected form of Lee, Yang, and Parr (LYP), and Perdew's 1986 gradient-corrected functional (P86). [8] Of all methods and basis sets tested, the B3LYP and the BP86 method, both in conjunction with the 6-311G(3d,p) basis set, were determined to be the most accurate. Both methods gave extremely accurate geometeries for the different benzoquinone derivatives and there was excellent agreement between calculated EAs and experimental values. This accuracy can be clearly seen in Table 2, where the B3LYP method was utilized.

Table 2: Calculated and experimental adiabatic EAs for an array of different methylated and halogenated p-benzoquinones using the B3LYP functional along with the 6-311G(3d,p) basis set. [8]

Takahata et al. [5] calculated the EA of certain molecules, grouped into $\pi$-system and $\sigma$-system molecules, using two different functionals (B88-P86 and B3LYP) combined with three different basis sets (6-31++G**, 6-311++G**, and aug-cc-pVTZ). From the different combinations of functionals and basis sets tested, the group came to the conclusion that the combination of B88-P86/6-31++G** (I.A) as well as B3LYP/6-31++G** (II.A) are the recommended approaches. A summary of their results can be seen in Table 3. There is no significant change in the calculated EAs between either two methods, however, a substantial difference in accuracy can be seen for the calculations performed on $\pi$-systems to those performed on $\sigma$-systems. It can also be seen that the calculated VDE values agree better to the experimental adiabatic EA results than the calculated adiabatic EAs, yet the group fails to determine why this could be.

Table 3: Calculated and experimental adiabatic EAsa and VDEsb for different $\pi$- and $\sigma$-system molecules in eV. Methods before and after the symbol // correspond to the energy and the geometry determined, respectively. AAD = average absolute deviation. Experimental adiabatic EA values obtained from the CRC Handbook of Chemistry and Physics, 1996. [5]

Apparently the use of the B3LYP functional is percieved to be the most accurate when performing EA calculations, especially in regard to $\pi$-conjugated systems. Nevertheless, Rienstra-Kiracofe et al. [1] state that because the adiabatic EA is calculated as the difference in energy between the total energies of the neutral molecule and the anion, and not through anion HOMO energy, DFT can be successfully utilized to predict EAs for molecules with inexact functionals and finite basis sets. Studies have been conducted to show that large errors are obtained when calculating individual orbital energies, yet these errors diminish once total energies of the molecules are calculated, as these are not a sum of the individual orbital energies. [1] This would suggest that the LDA method implemented within the SIESTA code would be sufficient to determine the EA of xylene. This could also explain why the group of De Proft et al. had very inaccurate calculation results compared to experimental values. Boesch et al. suggest that the calculation of adiabatic EAs typically coincide more with experimental values when compared to VDE calculations, since these deviate from experimental EAs by an adiabatic correction as much as 0.5 eV [9] Nonetheless, a major reason for additionally calculating nitrogen, oxygen, and SF6 in my proposed calculation method is to see just how accurate the calculations would be with respect to other DFT calculations and experimental results.

# Conclusion

DFT can be considered a very useful technique to investigate molecular properties which are difficult to measure experimentally, especially when such parameters as high voltage and electric field are present. A proposed calculation method has been presented to determine the EA, or the change in EA, of o-, m-, and p-xylene under the effect of an increasing electric field. By utilizing this method, I hope to be able to draw conclusions on why the cavity-induced PD was initially being suppressed in the epoxy samples.

# References

1. Rienstra-Kiracofe J.C., Tschumper G.S., Schaefer III H.F., Atomic and Molecular Electron Affinities: Photoelectron Experiments and Theoretical Computations, Chem. Rev., 2002, Vol. 102(1), 231-282.
2. Lemke, E., Fundamentals of PD Measurements, Lemke Diagnostics GmbH, 2006.
3. Ho J., Effect of Alteration of Antioxidant on the Dielectric Strength of Biaxially-Oriented Polypropylene, PhD Thesis University of Connecticut, 2006, 56-65.
4. Husain E., Nema R.S., Analysis of Paschen Curves for Air, N2, and SF6 using the Townsend Breakdown Equation, IEEE Trans. on Elect. Insul., 1982, EI-17(4), 350-353.
5. Takahata Y., Chong D.P., Density Functional Calculations of Molecular Electron Affinities, J. Braz. Chem. Soc., 1999, 10(5), 354-358.
6. Stevenson C.D., Wagner II E.P., Reiter R.C., An Isotopic and EPR Probe into Hyperconjugative and Induction Effects in p-Xylene, J. Phys. Chem., 1993, 97, 10585-10588.
7. De Proft F., Sablon N., Tozer D.J., Geerlings P., Calculation of negative electron affinity and aqueous anion hardness using Kohn-Sham HOMO and LUMO energies, Faraday Discuss., 2007, 135, 151-159.
8. Bruska M.K., Piechota J., Density functional study of sulphur hexafluoride (SF6) and its hydrogen derivatives, Molecular Simulation, 2008, 34(10-15), 1041-1050.
9. Boesch S.E., Grafton A.K., Wheeler R.A., Electron Affinities of Substituted p-Benzoquinones from Hybrid Hartree-Fock/Density-Functional Calculations, J. Phys. Chem., 1996, 100, 10083-10087.