by Vishal Jindal

1. Introduction

Thiophene is a heterocyclic aromatic compound with the formula C₄H₄S, consisting of a planar five-membered ring. The thiophene ring is the most widely used building block for the construction of conjugated polymers mainly because of its high environmental stability [1]. Thiophene-based π-conjugated systems have attracted much attention in developing high performance organic solar cells [2]. Here, in this post, we are using Density Functional Theory (DFT) [3] calculations to determine the vibrational frequencies of an isolated Thiophene molecule.

2. Methodology

Before finding the vibrational frequencies of thiophene, we first optimized the structure of the thiophene ring. A model of the thiophene ring was constructed using Material Studio.

In order to use a plane-wave basis set, we need to make the system periodic, so we used a box with side 15 Å and K-point mesh of 1x1x1 for the isolated thiophene molecule (as shown below).

Calculation Details

Structure optimization calculations were performed using the Vienna Ab initio Simulation Package (VASP), a plane wave basis set pseudo-potential code. The exchange and correlation energies were calculated using the Perdew, Burke, and Ernzerhof (PBE) [4] functional described within the generalized gradient approximation (GGA) [5]. A plane-wave basis set cutoff energy of 450 eV was used. The calculations were considered optimized when the force on every relaxed atom ( ionic convergence limit) was less than 0.02 eV Å^-1. The electron densities were self consistently solved for energy with the convergence limit (electronic convergence) set to 10^-5 eV.

The energy of the optimized structure is -54.237 eV (zero-point energy not corrected). After getting the final structure of the thiophene ring, we calculated the vibrational frequencies by varying the displacements (δb) used in finite-difference calculations for vibrational frequencies [3]. The displacement values (δb) used were 0.005, 0.006, 0.01, 0.025, 0.05, 0.07 and 0.1 Å.

3.  Results

Fig 4. The 5th vibrational mode (f5) for the thiophene molecule (~1500 cm^-1)

Thiophene (C₄H₄S) has 9 atoms. So, it has 3N – 6, i.e. 21 vibrational modes. The rest of the 3 rotational modes and 3 translation modes should be theoretically zero, but they have some small or imaginary value because of the small numerical inaccuracies that inevitably occur in calculating a Hessian via finite difference expressions and a calculation method with finite numerical accuracy [3]. Below are the plots of vibrational frequencies as a function of the inverse of displacement value (δb) used for calculating the Hessian matrix.

Zero-point energy (ZPE) corrected energy for isolated thiophene was also calculated for various displacement values (δb) using the equation given below taking into account 21 vibrational modes and E_o  = – 54.237 eV:

 

Table 1. Corrected zero-point energy and the absolute sum of rotational/ translational modes for different displacement values (δb).

Displacement value used for Hessian Matrix (δb)	Inverse of displacement value (1/δb)	Zero-point Energy (ZPE) Correction (eV)	ZPE Corrected Energy (eV) (Eo = -54.237 eV)	Absolute sum of Rotational/ Translational modes (cm-1)
0.005	200	1.783	-52.454	518.82
0.006	167	1.775	-52.462	467.79
0.01	100	1.771	-52.466	358.96
0.025	40	1.769	-52.468	223.58
0.05	20	1.770	-52.467	135.83
0.07	14	1.772	-52.465	131.04
0.1	10	1.778	-52.459	173.87

4. Conclusion

From the graphs above, we can notice that as the displacement magnitude is lowered, the frequencies of the modes converge to a stable value before they again diverge a little when the displacement value becomes too small. This is expected as harmonic approximation is valid at smaller displacements and can have huge deviations from reality for larger separations. But as the displacement value becomes too small, the energy differences also become too small relative to the SCF convergence which results in slight divergence from the “converged” value.

Also, we know that theoretically, the frequencies of translational/rotational modes should be 0. For DFT vibrational calculations, we can compute the sum of the absolute values of these frequencies (taking sum of last 6 normal modes), and closeness of that sum to 0 can be taken as a measure for better results. We did this for different values of displacement (δb), as shown in the last column of the table in the previous section.

Considering the factors mentioned in the previous two paragraphs, displacement values (δb) of 0.025, 0.01 and 0.006 Å (i.e. 1/δb = 40, 100 and 167 respectively), all give convergent and better results than other values of displacement used to calculate Hessian matrix. Therefore,  0.01 Å can usually be used for getting normal vibrational modes, once we have an optimized structure of the molecule of interest.

5. References

[1] Stergios Logothetidis, Handbook of Flexible Organic Electronics (2015) Materials, Manufacturing, and Applications, Section 4.2.1 (https://doi.org/10.1016/C2013-0-16442-2).

[2] Fan Zhang, Dongqing Wu, Youyong Xu, and Xinliang Feng, Thiophene-based conjugated oligomers for organic solar cells, J. Mater. Chem., 2011,21, 17590-17600.

[3] D. Sholl, J.A. Steckel, Density Functional Theory: A Practical Introduction, Wiley 2009

[4] J.P. Perdew, K. Burke, M. Ernzerhof, Generalized gradient approximation made simple, Phys. Rev. Lett., 77 (1996) 3865-3868.

[5] J.P. Perdew, J.A. Chevary, S.H. Vosko, K.A. Jackson, M.R. Pederson, D.J. Singh, C. Fiolhais, Atoms, Molecules, Solids, And Surfaces – Applications of the Generalized Gradient Approximation for Exchange and Correlation, Phys. Rev. B, 46 (1992) 6671-6687.