An Effective Procedure for Analyzing Molecular Vibrations in Terms of

Article pubs.acs.org/JCTC

An Effective Procedure for Analyzing Molecular Vibrations in Terms of Local Fragment Modes Miquel Huix-Rotllant* and Nicolas Ferré Aix-Marseille Univ, CNRS, ICR, Marseille, 13284, France S Supporting Information *

ABSTRACT: Normal mode analysis is nowadays routinely performed in most quantum chemistry codes as a valuable tool for studying molecular motions. In systems formed by several molecular units, such as noncovalently bound donor−acceptor systems, normal modes are frequently delocalized over the whole structure, leading to a complicated interpretation of vibrations. Here, we propose a simple procedure to localize normal modes, based on the definition of fragment Hessian submatrices that lead to legitimate local fragment modes. These fragment modes can be defined for any arbitrary number of atoms and fragments and form in general a complete basis that can be used to expand the total normal modes. We use this fragment basis to explore the molecular origin of the normal modes of benzene, and several peaks in the infrared spectrum of tetraalanine. Additionally, we propose a criteria to define breathing modes in noncovalently bound systems, based on the fact that motions of the center of mass of fragments are not translational and rotational invariant with respect to the total center of mass of the system. By using the Sayvetz-Eckart conditions of fragments, we can quantify the contribution of rotation and translation of the fragments center of mass on the normal modes of the total system.

1. INTRODUCTION Molecular vibrations are of utmost importance to describe any chemical process. Despite the successful separation of electron and nuclear motions in the so-called Born−Oppenheimer (BO) approximation,1 many systems require the treatment of both nuclear and electronic degrees of freedom at the same time. Paradigmatic cases are conical intersections or Jahn−Teller systems.2,3 In most chemical systems, nuclear motions are naturally expressed in terms of normal-mode analysis of semirigid bodies.4 In its simplest formulation, this analysis is performed by the diagonalization of the second-derivative of the BO energy with respect to the mass-weighted coordinates at a minimum energy molecular configuration. This offers a basis for expanding internal nuclear motions, after separating out rotational and translational motions of the total center of mass.4−6 Normal mode analysis is a very powerful tool to study molecular motions. This facilitates for example the assignment of infrared or Raman spectra in terms of internal vibrations of the molecule. However, normal modes are frequently delocalized over the whole structure, which complicates the assignment of a given absorption band to chemically meaningful vibrations “localized” on molecular subunits. This is specially critical in the assignment of spectra in complex systems such as peptides or molecules absorbed on surfaces, where the peak assignment is often controversial.7,8 Normal mode analysis in terms of local vibrations can be a helpful tool © 2016 American Chemical Society

in this respect, in order to identify the relevant motions or subunits dominating the chemical process. For example, in photochemical processes in noncovalently bound charge transfer complexes, the electron transfer or the charge recombination is mediated by the so-called “breathing modes”, which are delocalized vibrations that increase the electronic couplings between the different subunits. Such vibrations have been detected experimentally to be active after electron−hole separation in important materials such as at the interfaces of conjugated polymers and fullerene derivatives used for organic photovoltaics.9 For such cases, there is a need for a simple methodology to unequivocally identify breathing modes or analyze complex spectra in terms of normal modes localized in fragments. In the literature, several approaches have been proposed that localize the normal modes.10−17 For example, Konkoli and Cremer proposed a method in which a set of internal coordinates is used to describe the local vibrations of the fragment, from which a set of transformation matrices from Cartesian to internal coordinates are constructed.13 This transformation is however not an easy task in general and leads to a complicated form of the vibrational Hamiltonian.4 An easier procedure is that of Jakob and Reiher,17 which is based on the construction of unitary transformations that maximize Received: May 19, 2016 Published: September 19, 2016 4768

DOI: 10.1021/acs.jctc.6b00514 J. Chem. Theory Comput. 2016, 12, 4768−4777

Article

Journal of Chemical Theory and Computation

formed to angular frequencies (in wavenumber units) by applying the transformation

the distance between the center of vibration of the fragments. This procedure is more general in the sense that it does not need to define coordinates or fragment a priori, but the resultant local modes might not always be localized in one single fragment. For this reason, the method of Jakob and Reiher does not always lead to an exact frequency for the localized modes, but only approximate local frequencies are obtained. Additionally, this procedure can be difficult to apply in covalently bound systems in which the centers of vibration might be less clearly defined. Here, we propose a general procedure for localizing normal modes over any arbitrary number of user-defined fragments. The procedure requires only the computation of the total mass-weighted Hessian matrix from electronic structure methods. The local modes are obtained by extracting the fragment Hessians from the total Hessian, which gives a basis of fragment modes, that is, vibrations according to the center of mass of each fragment. The local modes within one fragment are fully orthogonal, whereas the interfragment local modes lead to bilinear couplings. In our method, frequencies for the local modes are exactly obtained. The mode localization method is presented as follows: in section 2 we give the details of the normal mode localization procedure and the analysis of breathing modes, in section 3 we briefly describe the implementation in a Fortran code and the computational details, in section 4 we show some examples of the origin of vibrations of benzene in terms of local oscillators, the analysis of the infrared spectrum of the oligopeptide tetraalanine and the breathing modes in the tetracyanoethylene−mesitylene charge-transfer complex. In section 5 we end with some conclusions.

ωI =

N frag

1 ∂ 2E(R) mI mJ ∂RI ∂RJ

in which ⊕ symbolizes the matrix direct sum (i.e., the construction of a block diagonal matrix from the set of Xf), and f is the index running over the total number of fragments Nfrag. Each fragment has Nfat atoms, ranging from 1 to the total number of atoms in the molecule Nat. Obviously, the constraint frag is that ∑Nf Nfat = Nat. Additionally, we constrain each fragment transformation to be unitary on its own space X†f X f = 1

(7)

This transformation defines a new coordinate system

Qloc = XlocMR

(8)

The new coordinate system is related to the one defined in eq 2 by the transformation Qtot = XtotXloc, †Qloc

(1)

(9)

The local to total normal mode transformation matrix XtotXloc,† can be used to analyze the contribution of fragment modes in the total normal modes. Indeed, the local modes Xloc form a complete basis for describing internal nuclear motions, and thus can be used to analyze the total normal modes in the basis of local modes. One can easily obtain the contribution weight of the local modes in the total modes by extracting the projection of the eigenvectors corresponding to the total modes on the local mode eigenvectors, which is simply obtained by the scalar product of the column j of the fragment part of the total normal mode with the i column of the fragment mode,

(2)

(3)

pjif = X jf ,tot ·X if

where E(R) is the Born−Oppenheimer total energy. The total dimension of the Hessian matrix defined in eq 3 is 3Nat × 3Nat. At a minimum energy structure, ∇E(R) = 0, the total Hessian matrix is positive semidefinite. Its diagonalization thus gives access to the normal modes, λ = Xtot, †HtotXtot , λII ≥ 0

(6)

f =1

constructed from the Cartesian coordinate vector R in Bohrs, the mass-metric MIJ = mI δIJ , where mI is the atomic mass of atom I and the matrix Xtot, which is the eigenvector matrix of the total Hessian in MWC coordinates, is defined as HIJtot =

(5)

Xloc = ⊕ X f

where the mass-weighted Cartesian (MWC) normal coordinates are defined as Qtot = XtotMR

4π 2c 2

which are the frequencies of normal modes entering in the second term of eq 1. The eigenvectors of the total Hessian matrix Xtot give the simplest form to the vibrational Hamiltonian. However, this offers a complexity in the analysis of the normal modes since they are frequently delocalized over the whole molecular system. However, it is more convenient to have local normal modes defined over the constituent fragments of the molecules. Here, we propose an efficient method to define the fragment modes by constructing an appropriate unitary transformation of the normal mode coordinates. The unitary transformation is based on the direct sum of normal modes of the constituent fragments,

2. METHODOLOGY 2.1. Fragment Modes. In semirigid molecular systems, the total vibrational Hamiltonian in the harmonic approximation after separation of rotations and translations is defined (hereafter all equations are in atomic units) as, 1 1 Ĥ (Qtot) = ∇Qtot, †∇Qtot + Qtot, †ω2 Qtot 2 2

λII

(10)

where the center dot (·) symbolizes the dot product between the two vectors, Xf, tot corresponds to the block of the total normal modes corresponding to fragment f, and Xf is the local mode eigenvectors to be defined later in this section. Index j runs over the total internal normal modes, whereas index i runs over the fragment modes. The local mode transformation changes the form of the vibrational Hamiltonian. Inserting the total to local coordinate transformation defined in eq 9 into eq 1 leads to

(4)

where λ and Xtot are the matrices containing the eigenvalues and eigenvectors, respectively. The eigenvalues can be trans4769


Article


Figure 1. Graphical representation of the partitioning of the total Hessian matrix Htot in fragment submatrices Hf (shown in colors), where f is the index for fragments. The fragment submatrices are diagonalized individually to obtain the fragment modes Xf, which are then collected in the matrix Xloc. The latter is used as a basis to expand the real eigenvectors of Htot, thus extracting the information in terms of local modes.

1 1 Ĥ (Qloc) = ∇Qloc, †∇Qloc + Qloc, †XlocΩtotXloc, †Qloc 2 2

HIJf = (11)

where we defined the matrix Ωtot as Ωtot = Xtot, †ω2 Xtot

1 ∂ 2E(R) I, mI mJ ∂RI ∂RJ

J∈f (13)

where f is the index that defines the atom subspace of the fragments (this notation for the index f will be used hereafter). The diagonalization of the local Hessian gives a set of fragment modes,

(12)

which is simply the total Hessian of the system in the units of frequency. It is important to note that the Ω is a nondiagonal matrix in general, with the same structure as eq 3. At this point, we only need to define the construction of the unitary transformation Xf that will lead to the localized fragment modes in the desired form. To define such transformation, we exploit the fact that any submatrix of a positive semidefinite matrix is also positive semidefinite.18 Accordingly, one can extract Nfrag fragment Hessian submatrices from the total Hessian,

λf = Xf , †H f Xf ,

λIIf > 0

(14)

The eigenvalues can be transformed to angular frequencies (in wavenumber units) by applying the transformation ωI f = 4770

λIIf 4π 2c 2

(15) DOI: 10.1021/acs.jctc.6b00514 J. Chem. Theory Comput. 2016, 12, 4768−4777

Article


easiest ways to construct such a transformation is described in the Gaussian white paper on vibrational analysis.20 The translation conditions are simply imposed by constructing the following three vectors

To make the link between the eigenvectors defined in eq 14 with the local unitary transformation defined in eq 6, it is simply necessary to diagonalize the local Hessian matrix N frag

Hloc = ⊕ H f

N at

(16)

f =1

I=1

The main advantage of this choice of unitary transformation is that eq 11 brought to the simple form 1 Ĥ (Qloc) = 2 N

+

N

† † 2 loc ∇Q loc + Q loc, ωf Q f ) ∑ (∇Q loc, f f f f

frag

† f tot f ′ , † loc X Ω ff ′X Q f ′ , ∑ Q loc, f

(17)

Qloc f

in which the local oscillator coordinates are defined as = Xf Mf Rf refers to the mass-weighted coordinates of local tot oscillator f, Ωtot ff ′ is the off-diagonal block of matrix Ω . Note that the center of vibration for the local oscillators is now the center of mass of the fragment as shown in the definition of Qfloc.This equation together with the local unitary transformation are the central results of this article. The local fragment modes can be thus interpreted as local oscillators at frequencies ω f (different from the total normal-mode frequencies ω), bilinearly coupled to the other oscillators localized at different parts of the molecule. To summarize the construction of the local transformation, we graphically represent the partitioning of the total Hessian matrix in submatrices in Figure 1. In general the eigenvalues of each fragment Hessian are all positive, since the rotations and translation vectors of the fragments are not invariant with respect to the whole structure. The number of modes is 3Nfat per fragment. The local modes Xloc are thus defined as a blockdiagonal matrix of dimension 3Nat, collecting the Nfrag blocks with the fragment modes Xf. In this way, when this is used to transform the total Hessian, it leads to diagonal blocks which contain the frequencies representing the local oscillators, which are coupled bilinearly by the transformed off-diagonal blocks (as shown in eq 17.) 2.2. Fragment Sayvetz−Eckart Conditions. For nonlinear molecules, the number of internal vibrations of the total Hessian matrix are 3Nat − 6. The six remaining vectors have an eigenvalue close to 0, corresponding to translations and rotations of the center of mass of the total system. These can be separated out by imposing the so-called Sayvetz−Eckart conditions,4−6 which impose two sets of conditions for each atom I, one for translations

(21)

R0I,1,

I

R0I,2

R0I,3

where and are the x, y, and z coordinates of atom I. The eigenvectors K that diagonalize the moment of inertia tensor IM are then used to construct three rotation vectors of D, N at

D4 = ⊕

I=1 N at

D5 =

⊕ I=1 N at

D6 = ⊕

I=1

1 ((K 2·R0I )K3 − (K3·R0I )K 2) mI 1 ((K3·R0I )K1 − (K1·R0I )K3) mI 1 ((K1·R0I )K 2 − (K 2·R0I )K1) mI

(22)

in which the Ki·RI represents the dot product between the column i of eigenvectors of IM and the Cartesian coordinates of atom I. The local oscillators define semirigid body motions in fragments, in which the center of oscillation is the center of mass of the fragments. Thus, the Sayvetz−Eckart conditions imposed for the total center of mass of the system might not be fully satisfied by the fragment modes. In other words, there might be a contribution of the internal translations and rotations of the center of mass of the fragments in each of the internal total normal modes. To analyze and quantify the amount of translation or rotation of a fragment in the total modes, it suffices to project the total normal modes to the fragment Sayvetz−Eckart conditions,

(18)

dijf = X jf ,tot ·Dif

and another for rotations,

∑ mI R0I (RI − R0I ) = 0

(20)

Nat Nat ⎛ Nat ⎞ 0,2 0 0 ⎜∑ mI (R0,2 − ∑ mI R0I ,1R0I ,3 ⎟ I ,2 + R I ,3) − ∑ mI R I ,1R I ,2 ⎜ I=1 ⎟ I=1 I=1 ⎜ ⎟ N N N at at at ⎜ ⎟ 0,2 0 0 ⎟ + − R ) m R R IM = ⎜− ∑ mI R0I ,1R0I ,2 ∑ mI (R0,2 ∑ I ,1 I ,3 I I ,2 I ,3 ⎜ I=1 ⎟ I=1 I=1 ⎜ ⎟ Nat Nat Nat ⎜ ⎟ 0 0 0,2 − ∑ mI R0I ,2R0I ,3 + R ) ∑ mI (R0,2 ⎜⎜− ∑ mI R I ,1R I ,3 I ,1 I ,2 ⎟ ⎟ ⎝ I=1 ⎠ I=1 I=1

f ,f ′

I

i = 1, 2, 3

where I is the index running over the number of atoms, mI is the atomic mass of atom I, and ei, I is the unit vector of atom I and i = 1,2,3 corresponds to x, y, and z coordinates of atom I, respectively. The transformation imposing the conditions on the rotations of the center of mass involves the construction of the 3 × 3 moment of inertia tensor IM,

frag

∑ mi(RI − R0I ) = 0

mI ei , I

Di = ⊕

(23)

which measures how rotation and translation of the fragment f take part of the total normal modes j. In this equation, i = 1,2,3 corresponds to the x, y, and z translation vectors defined in eqs 20 and i = 4,5,6 corresponds to the rotation vectors defined in eq 22. The label f simply restricts the atoms included in that fragment. Another complementary criteria can be used for quantifying the contribution of breathing modes in a total normal mode. This criteria is based on calculating the displacement of the

(19)

in which R0I is the Cartesian coordinates of the reference position for atom I. These conditions generate a unitary transformation matrix D which sorts out translation and rotation vectors from the internal vibrations. Such unitary transformation is not unique, and several methods have been devised for imposing such conditions.4,19 Perhaps one of the 4771


Article


criteria to extract the breathing modes of the normal modes of noncovalently bound dimers. Local Mode Analysis of Benzene. Here, we analyze the vibrations of a benzene molecule in terms of local oscillators. As described in section 2.1, the local modes are harmonic oscillators coupled bilinearly to the rest of the fragment modes. Benzene can be divided in six equivalent C−H fragments containing two atoms each. Thus, the local oscillators will correspond to one internal vibration and five rotations and translations, as represented in Figure 2. The internal vibration

center of mass of a fragment with respect to the total center of mass of the whole system, that is Nf

f ΔR CM, J =

∑I ∈atf mI (R IJ − R0I ) Nf

∑I ∈atf mI

Nf

=

∑I ∈atf

mI XIJ Nf

∑I ∈atf mI

(24)

where the index J runs over the normal modes, I runs over the atoms belonging to fragment f, R0 is the Cartesian coordinates of the minimum-energy geometry, and R0I is the Cartesian coordinates of atom I, R IJ = R0I + XIJ / mI is the Cartesian coordinates of atom I displaced in the direction of normal mode J. Since the total center of mass is fixed, most normal modes do not involve a change of the center of mass of the constituent fragments. Only a reduced number of them, those in which translation or rotational modes of fragments are mixed with the internal vibrations, the vectors in eq 24, should be different from the zero vector 0. However, these movements are correlated, that is, the movement of the center of mass of one fragment in one direction must be counterbalanced by the movement of another fragment in the opposite direction in order to keep the total center of mass fixed.

3. COMPUTATIONAL DETAILS We have implemented the procedure described in Section 2 in a Fortran code. The program requires only the Cartesian Hessian matrix from Gaussian 09 formatted checkpoint files.21 The fragments can contain any number of atoms in any order. From this list, the program extracts the fragment submatrices as defined in eq 13, which are subsequently diagonalized to obtain the basis of fragment modes and local frequencies (eqs 14 and 15). The total normal modes are then projected in this basis, giving access to the expansion coefficients used for the analysis of spectra and to the couplings between fragment normal modes (see eqs 6 and 10). Additionally, we have implemented the projection of the total normal modes in the fragment Sayvetz−Eckart conditions (eq 23), which gives information on the contribution of each fragment translation and rotation vectors to the total normal modes. We have also implemented the change in the center of mass position of fragments (eq 24), which is output to be graphically represented by typical molecular visualization software. All minimization and Hessian calculations have been performed using the Gaussian 09 package.21 For the benzene molecule, both the geometry optimization and the frequencies have been calculated at the MP2/6-311G** level. The oligopeptide (alanine)4 has been optimized at the B3LYP/631++G** level of theory with a scaling factor of 0.970 for the frequencies to match better the experimental values.22−25 The tetracyanoethylene−mesitylene complex has been optimized at the MP2/6-31++G** level for the ground state minimum and at the CIS/6-31++G** level for the charge-transfer state minimum.22,23,26−28

Figure 2. Classification of the local vibrations of one of the fragments when benzene is partitioned in six C−H fragments. The oscillators are classified according to (a) out-of-plane local oscillators, (b) in-plane local oscillators, and (c) internal vibrations. Internal oscillations conserve the position of the center of mass of the local oscillator, whereas the remainder originate from translations or rotations with respect to the center of mass of the local oscillator.

has a local frequency of 3216.03 cm−1 and corresponds to the C−H stretching (Str). This is the only vibration which keeps the position of the center of mass of the local oscillator. The remainder of the local oscillators correspond to translations and rotations with respect to the center of mass of the local oscillator. These oscillators have nonzero frequencies since they are not invariants with respect to the total center of mass. The lowest frequency oscillators of this type have frequencies of 312.76 and 833.35 cm−1, and correspond respectively to a translation and rotation in the x direction (hereafter named Tx and Rx), which correspond to the direction out-of-plane with respect to the plane defined by the benzene ring. Three inplane translations and rotations occur in the atomic C−H axis at 921.06 cm−1 (hereafter named Tz), and perpendicular to this axis at 1033.74 and 1395.23 cm−1 for the translation and rotation respectively (hereafter named Ty and Ry.) The expansion of the total internal normal modes of benzene on the basis of local oscillators defined in Figure 2 shows interesting features. The total normal modes are either expanded by local vibrations of out-of-plane, in-plane, or internal type, but the three spaces are orthogonal. This has important consequences on the shape and frequencies of the total normal modes of benzene. For example, this feature explains the fact that the frequencies of the total normal modes

4. RESULTS AND DISCUSSION We present three different applications of the analysis of the normal modes with fragments. First, we describe the origin of the normal modes of benzene in terms of local oscillators. Second, we show the decomposition of the infrared spectrum of an oligopeptide into constituent fragments, which facilitates the interpretation of the main spectral peaks. Third, we give a 4772


Article

Journal of Chemical Theory and Computation appear in three regimes of low-frequency (around 0−800 cm−1), midrange (around 800−2000 cm−1), and high frequency (above 3000 cm−1). In the following, we discuss the three spaces in detail. The high frequency total normal modes are represented by the combination of a unique internal local oscillator of frequency 0.399 eV (3218 cm−1). The structure of Ωtot for this space shows a Hückel-like structure, that is, each stretching local oscillator is only coupled to the nearest neighbor with a fixed-valued coupling of 1.288 meV. This leads to the following local oscillator matrix (in eV), ⎡ 0.399 ⎢ ⎢1.288· 10−3 ⎢ ⎢0 ⎢ ⎢0 ⎢ ⎢0 ⎢⎣ 1.288· 10−3

1.288· 10−3 0 0.399

−3

1.288· 10

1.288· 10−3 0.399 0

0

0

0

0

1.288· 10−3 0

1.288· 10−3 0.399

1.288· 10−3 −3

0

0

1.288· 10

0.399

0

0

0

1.288· 10−3

Table 1. Expansion of a Selected Set of Normal Modes of Benzene in Term of the Local Oscillators Depicted in Figure 2, Coming from the Partitioning of Benzene in Six Equivalent C−H Fragments

1.288· 10−3 ⎤ ⎥ ⎥ 0 ⎥ ⎥ 0 ⎥ 0 ⎥ ⎥ 1.288· 10−3 ⎥ ⎥⎦ 0.399

(25)

which is orthogonal to the rest of the local oscillators, and thus the eigenvectors of this matrix (after transformation to cm−1) gives directly the highest frequency spectrum of benzene. In particular, the eigenvectors of such a matrix leads to the spectrum of 3195, 3206, 3206, 3226, 3226, and 3237 cm−1, which matches perfectly with the six highest frequency modes of benzene obtained by diagonalization of the full Hessian matrix, 3198, 3209, 3209, 3224, 3224, 3234 cm−1. A similar construction could be done for the out-of-plane vibrations and the in-plane vibrations, however the structure of the coupling matrix Ωtot is less obvious in these cases. In Table 1, we show several examples of the expansion of out-of-plane and in-plane normal modes in terms of local oscillators. In the table, we selected a restricted number of total normal modes to show the main characteristics of the type of local oscillator combinations that represent the low- and midrange vibrations. The complete list with the analysis of all vibrations is reported in the Supporting Information. As previously mentioned, the out-of-plane total normal modes are represented by combinations of Tx and Rx local oscillators. There is a relation between the frequency of the total modes and the frequency of the local oscillators on which the total mode is expanded. For example, the total mode at frequency 387 cm−1 is very close to the frequency of the local oscillator Tx (313 cm−1), whereas the total mode at 839 cm−1, which involves essentially Rx contributions, oscillates to a similar frequency of the Rx local oscillator (833 cm−1). More complex combinations of local oscillators can occur. For example, the total normal mode at 393 cm−1 involves out-of-plane movements of the carbon atoms only, whereas the hydrogen atoms are essentially fixed. On the basis of local oscillators, such a normal mode is represented by a linear combination of a translation and a rotation in the same direction. As it can be seen in Figure 2, translations involve the movement of both atoms in the same direction, whereas for rotations the atoms move in opposite directions. Thus, a linear combination of a rotation plus a translation decreases the hydrogen atom movement, while it will enhance the carbon atom movement. The in-plane vibrations involve three local vibrations. Some total normal modes follow similar trends to the out-of-plane vibrations. The vibration at 1007 cm−1, for example, involves mainly Tz local oscillators, and its characteristic frequency is close to the frequency of the local oscillator 921 cm−1. More

a

Due to the inversion center of symmetry of benzene, fragments 1 and 4, 2, and 5, and 3 and 6 are equivalent. Therefore, we only specify the contribution sign change, since the numerical value of the coefficient is the same. For the complete list see the Supporting information. Only coefficients larger than ±0.1 are shown.

interestingly, there is the representation of C−C stretching in terms of local oscillators. In this case, the carbon atoms move in an “oblique” direction with respect to the direction defined by the in-plane local oscillators. This would not be possible by combining only one translation and one rotation in the same direction as in the case of out-of-plane vibrations. In contrast, the in-plane space contains Ty and Tz local oscillators, which gives a complete basis to expand vectors in any direction in the plane defined by the benzene ring. This is clearly shown for the vibration at 1641 cm−1 in which two C−C stretching vibrations are involved. Still, the combination of Ty and Ry can lead to some C−C stretching by alternating + and − signs in contiguous fragments, as in the case of the 1446 cm−1 vibration. To conclude, we should point out that it is in general not possible to have an in-plane C−C stretching vibration involving only the movement of carbon atoms with fixed hydrogen atoms (in contrast with the out-of-plane vibrations, where this is possible.) This is because the cancellation of the hydrogen movements involves the combination of a translation and a rotation exclusively in the same direction. This is not possible for C−C stretching vibrations, since these are represented at least by two translations in perpendicular directions. Infrared Spectrum of Alanine Tetramer. The analysis of the infrared spectrum of polypeptides is often complex due to the overlap between vibrations from each of the constituent units. In such a case, the localization of normal modes can be a fundamental tool to simplify the interpretation and determine the contributions of the constituent units to each of the peaks. To illustrate this, we analyze the infrared spectrum of a small oligopeptide, the tetra-alanine (Ala)4. In Figure 3, we show the decomposition of its infrared spectrum in constituent fragments. To construct the fragment spectrum, we have first 4773


Article


expanded each normal mode of the total molecule on the basis of fragment modes. This gives the contribution of the fragment to the normal mode. To produce the spectrum, we multiply the intensity by this coefficient and by a Gaussian function that gives a phenomenological broadening of the peak. In principle, one could recompute the intensities in the new basis, as done by Jacob and Reiher.17 However, the procedure described here is enough for the purpose of offering a simple analysis tool of the total spectrum. We first concentrate on the most intense peaks of the infrared spectrum. The total spectrum is characteristic of short oligopeptides, featuring an intense peak at 1531 cm−1 (peak 1 in the figure) and other less intense peaks at 1275 and 1353 cm−1 forming a shoulder (peak 2), at 1720 and 1742 cm−1 forming the second most intense peak (peak 3), and a smaller peak at 1841 cm−1 (peak 4). The origin of these peaks is apparent by analyzing the normal modes in fragments. Several levels of information can be extracted depending on the chosen fragmentation. To illustrate this, we consider three different cuts in fragments: the central peptide bond to build (Ala)2 dimers, all the peptide bonds to separate the alanine subunits, and all the Cα−CO and OC−N bonds to generate eight fragments. First, we focus on the analysis in the basis of the dimer fragments (see Figure 3,a). In this case, fragment 1 contains the terminal amine group (N terminal), whereas the fragment group contains the terminal carboxylic acid group (C terminal). This is reflected in a nonequivalent contribution to the most intense peaks. This is clearly seen in peak 2 and peak 4, where the main contribution is due to fragment 2, and peak 3, where 66% of the contribution is due to fragment 1. The intense peak 1 has almost the same contribution from both fragments. From the separation in alanine subunits (Figure 3b), we can get a further level of insight on the origin of the peaks. In this case, fragment 1 contains the N terminal group, fragments 2 and 3 are similar intermediate alanine units, and fragment 4 contains the C terminal group. From this partitioning, one can observe more information than in the first fragmentation. On the one hand, peaks 2 and 4 are now mainly due to the C terminal group, whereas peak 3 has mainly a contribution from fragment 2, a small contribution from the N terminal group and fragment 3, and no contribution from the C terminal unit. Peak 1 has contributions from all fragments instead. There is a relation of the present partitioning with the previous dimer partitioning. The sum of alanine units in fragments 1 and 2 is equal to the fragment 1 alanine dimer, whereas the sum of fragments 3 and 4 leads to fragment 2 of the previous partitioning. Equivalently, the sum of the spectra of the corresponding subunits leads to the dimer spectrum. This is an important feature, since the different partitionings give converging results. The fact that the fragments can be of any arbitrary size thus allows us to perform partitionings in chemically meaningful units, which can frequently be different from the monomeric units that form the peptide chain. The most detailed level of information is given by the final partitioning in the peptide bond units, as shown in Figure 3c. In this case, one can readily identify peaks 3 and 4 originating from the carbonyl stretching and the C terminal stretching, respectively. More important are peaks 1 and 2, which mainly originate in the central amine groups. In these peaks, some contribution from the carbonyl stretching is also observed, but the contribution from the N-terminal group is residual. In fact, in this region the N-terminal group is contributing to a small

Figure 3. Contribution of the fragment vibrations to the total infrared spectrum of tetra-alanine (Ala)4. The frequencies are in cm−1, whereas the intensities are in km/mol. Different partitioning of the oligopeptide (Ala)4 in fragments from which a basis of local modes is constructed; (a) partitioning in dimer (Ala)2 units; (b) partitioning in individual (Ala) subunits; (c) partitioning in peptidic bond units. Fragments are indicated by colored squares on the molecular structure. 4774


Article

Journal of Chemical Theory and Computation shoulder between peaks 1 and 3, while its remaining contribution to the spectrum occurs at frequencies lower than 1000 cm−1. The peak at around 371 cm−1 and the peak at 847 cm−1 correspond clearly to the N-terminal vibrations, whereas the peak at 481 cm−1 corresponds to the C-terminal vibrations. Finally, the peak at 615 cm−1 is an amine vibration of fragment 3. The structure of the peaks can be rationalized in terms of the local oscillator frequencies, similar to the previous analysis for benzene. The frequencies for the C−O local stretching are 1703, 1707, and 1709 cm−1 for the internal carbonyl groups, and 1835 cm−1 for the carbonyl in the C-terminal group. The latter is higher in energy due to the coupling with the terminal OH group, which is more flexible than the peptide backbone and couples more strongly to the CO stretching. This is at the origin of the two distinct CO peaks with an intensity ratio 3:1 at the high-frequency part of the spectrum. Similar to the case of benzene, the CO stretching is an internal local oscillator which does not mix strongly with other types of vibrations. This explains the origin of pure CO stretching peaks in the infrared spectrum. This is not the case for the amine vibrations, which are spread over the lower- and midrange frequencies and are mixed. The most intense peak for example is represented mainly by two local oscillators of 1473 and 1500 cm−1, which correspond mainly to in-plane translations and rotations in the y direction of the N−H moiety within the amine fragment. As previously seen for benzene, such vibrations tend to be more mixed, and thus couple strongly to the local oscillators of the same type of the CO fragments, which occur at 1182 cm−1. Breathing Modes. Breathing modes are often discussed in the literature, especially in the context of charge transfer complexes.9,29−31 Such modes are delocalized over weakly interacting fragments and are responsible for increasing the density overlap and thus the electronic coupling. In typical molecular semirigid bodies, normal-mode analysis imposing the Eckart conditions guarantees that the internal vibrations do not change the position of the total center of mass.4−6 However, there is no guarantee that the position of the centers of mass of the fragments remain fixed in the internal modes of the total molecule. In fact, this breaking of the Eckart conditions of the fragments is supposed to explain the origin of such modes such as the breathing modes, which usually occur in the lowfrequency spectrum. Breathing modes can thus favor the crossing of bright exciton states (typically localized on one fragment) to states of charge-transfer character delocalized over two or more molecular entities. In other words, one can roughly define the breathing mode(s) between two noncovalently bound fragments as the normal mode(s) principally built from translations of the fragments’ centers of mass in the direction connecting the two centers of mass. In eqs 23 and 24, we proposed two complementary criteria to analyze the breathing modes in noncovalently bound systems. To illustrate the usefulness of such criteria, we apply them to the analysis of the normal modes of the tetracyanoethylene−mesitylene complex (see Figure 4). First, we concentrate on the second definition (eq 24), which allows us to plot the changes in the position of the fragments center of mass within any global normal mode. The analysis is performed using the lowest charge-transfer excited state minimum energy structure calculated at the CIS level of theory and the 6-31+ +G** basis set. This model is chosen since it features analytic Hessian for excited states and gives reasonable estimates of the frequencies in charge−transfer complexes, even though the

Figure 4. Normal modes changing the position of the center of mass of the fragments in the tetracyanoethylene−mesitylene complex for the charge-transfer minimum CIS/6-31++G**. The left image corresponds to the normal mode, whereas the right image shows how the normal mode variates the center of mass of the fragment (for the definition, see eq 24). Carbons are in gray, hydrogens in white, and nitrogens in blue, whereas the center of mass of each fragment is depicted in pink.

excitation energies are usually too high.32 However, this is usually a constant shift, and the curvature around the excited state minimum is correct.33 For this reason, CIS is sufficient for the present purpose of analyzing excited state vibrations. We start with two fragments, corresponding to the tetracyanoethylene on the one hand and the mesitylene on the other hand. Performing the analysis given in eq 24 allows us to easily identify which modes change substantially the position of the center of mass (CM) (Figure 4). Three types of movements of the center of mass of the fragments are identified in the low frequency part of the vibrational spectrum: the translations along the x-axis and y-axis and the breathing modes corresponding to the stretching parallel to the z-axis. The xand y-axis translations are due to the normal modes with frequencies 52 and 78 cm−1, and move one fragment with respect to the other in the x and y direction, respectively. The breathing modes that move the fragments in the z direction (thus changing the distance between the two fragments) occur at frequencies 116 and 125 cm−1. In this case, the breathing modes are split in two vibrations, which by symmetry correspond in fact to the mirror image of one another. The aforementioned vibrations are the only modes that change non-negligibly the relative positions of the centers of mass of the two fragments. It has been observed experimentally by Raman spectroscopy that in similar systems the 165 cm−1 peak becomes active after the electron transfer occurs.34 This band has been assigned to an intermolecular stretching mode. 4775


Article


in fact, the z-axis translation contribution of the 101 cm−1 mode comes mainly from the movement of the cyano groups, which do not increase the π-stacking interaction between the two fragments. In contrast, the breathing mode at 176 cm−1 has a strong contribution from the z-axis translation of the ethylene group, and thus this mode can be clearly identified to the experimentally observed breathing mode.

However, our previous analysis has identified the breathing modes at much lower frequencies. Actually, looking at the global normal modes of the complex, we find a similar frequency as the one reported experimentally at 176 cm−1. The corresponding normal mode involves an out-of-plane displacement of the ethene double bond in percyanoethylene. More precisely, the ethylene bond is displaced toward the benzene ring, whereas the cyano groups move in the reverse direction. Hence the position of the fragment center of mass remains constant. Therefore, the previous two-fragments analysis could not identify this mode as a breathing mode. Further partitioning in three fragments, namely the CC double bond, the four cyano groups, and the mesitylene is physically motivated by the fact that the movement of the CC double bond in the direction of mesitylene is energetically favored, increasing the π-stacking interaction. Moreover, the cyano groups interact sterically with the methyl substituents of mesitylene, and their displacement in the reverse direction is also favored. Within the three-fragments partitioning, these effects become apparent. In Figure 5, we show the global

5. CONCLUSIONS We have developed a method to localize the normal modes of semirigid molecular systems according to user-defined fragment subunits. The method is based on the mathematical property that the submatrices of a semipositive definite matrix are also semipositive definite, and thus lead to a basis of fragment modes from which the global normal modes can be expanded. This allows us to perform an analysis based on fragments which can be defined arbitrarily, in number and in composition. The vibrational Hamiltonian in the basis of local oscillators shows an interesting structure: oscillators localized in fragments vibrate at characteristic frequencies and are bilinearly coupled to oscillators in other local fragments. This method can be applied to covalently as well as noncovalently bound systems. A strong advantage of the present method with respect to previously reported normal mode localization procedures17 is that the choice of fragments is arbitrary: the basis of fragment modes can be selected by the chemical problem in hand. The analysis of the vibrational spectrum of benzene in six equivalent CH fragments has led to interesting conclusions as to the origin of the vibrational modes of benzene. We have shown that the high frequency modes can be described by a single internal local oscillator with nearest neighbor couplings. We have also shown that out-of-plane modes involving mainly carbon atoms are described by a combination of rotation and translation in the same direction. Finally, we have shown that the in-plane vibrations involve the mixture of three local oscillators, and it is rare to have only carbon atoms moving as in the out-of-plane motions, but rather both hydrogen and carbon atoms must move. We have shown the usefulness of the fragment normal-mode analysis in unraveling the complex infrared spectrum of a small alanine oligopeptide. We have shown that the different fragment partitioning in dimers, alanine units, and elementary chemical units can give an increasing number of details on the nature of the peaks. This allowed us to readily identify the carbonyl stretching peaks appearing around 1700−1800 cm−1 and also quantify the degree of mixing of the amine vibrations with the carbonyl stretching modes occurring around 1500 cm−1. Importantly, since the fragment basis is complete to expand the total internal normal modes, the different partitioning gives convergent results. In other words, we can obtain the same total spectrum from the fragment spectra for each of the chosen partitioning. Finally, we have proposed two criteria to identify the breathing modes in noncovalently bound complexes. We have applied the criteria to analyze the tetracyanoethylene− mesitylene charge-transfer complex, showing that some internal vibrations move the center of mass of the fragments, although the overall center of mass of the system remains unmodified. When choosing the partitioning in two fragments (tetracyanoethylene and mesitylene), we identify two breathing modes around 100 cm−1, which do not match the experimental detected value of 165 cm−1. The analysis of this mode has been achieved through a three-fragment partitioning, motivated by

Figure 5. Identification of the breathing mode in the three fragment partitioning of the complex. The fragments correspond to (1) the cyano groups, (2) the ethylene double bond, and (3) the mesitylene. In the left image, the movement induced by the breathing normal mode is shown, while in the right image the movement of the center of mass of the fragments induced by the breathing mode is seen. For simplicity, only the movements of the ethylene and mesitylene centers of mass are shown.

normal mode at 176 cm−1 together with the movements of the centers of mass of the fragments. While the simple inspection of the normal mode atomic displacements was not sufficient, the analysis based on the three-fragment partitioning demonstrate that this normal mode indeed contributes to the breathing of the complex. From the analysis of the projection of normal modes in the Sayvetz−Eckart conditions of fragments described in eq 23 we arrive to similar conclusions. The partition in two fragments, the mesitylene on the one hand and the tetracyanoethylene on the other hand, one finds that the vibrational modes at 52, 72, 78, 101 and 105, 116, and 125 cm−1 have a major contribution of translation or rotation (larger than 20%), similar to what we found in Figure 4 with the analysis of eq 24. The modes at 52, 72, 105, 116, and 125 cm−1, even though they have an important contribution of translation, are mainly described by rotation contributions. The mode at 101 cm−1 has a clear strong contribution of z-axis translation with minor contributions from rotation. In the two-fragment analysis, the experimental breathing mode of 176 cm−1 has a negligible contribution from translation or rotation modes as we previously found with the analysis of eq 24. When tetracyanoethylene is partitioned in two fragments, it is clearly seen that, 4776


Article

Journal of Chemical Theory and Computation the favorable π-stacking interaction between ethene and mesitylene units and the unfavorable steric interactions between cyano and mesitylene methyl moieties.

■

(21) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; Nakatsuji, H.; Caricato, M.; Li, X.; Hratchian, H. P.; Izmaylov, A. F.; Bloino, J.; Zheng, G.; Sonnenberg, J. L.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Vreven, T.; Montgomery, J. A., Jr.; Peralta, J. E.; Ogliaro, F.; Bearpark, M.; Heyd, J. J.; Brothers, E.; Kudin, K. N.; Staroverov, V. N.; Kobayashi, R.; Normand, J.; Raghavachari, K.; Rendell, A.; Burant, J. C.; Iyengar, S. S.; Tomasi, J.; Cossi, M.; Rega, N.; Millam, J. M.; Klene, M.; Knox, J. E.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Martin, R. L.; Morokuma, K.; Zakrzewski, V. G.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Dapprich, S.; Daniels, A. D.; Farkas, Ã .; Foresman, J. B.; Ortiz, J. V.; Cioslowski, J.; Fox, D. J. Gaussian 09, revision D.01; Gaussian Inc.: Wallingford CT, 2009. (22) Ditchfield, R.; Hehre, W. J.; Pople, J. A. J. Chem. Phys. 1971, 54, 724. (23) Francl, M. M.; Pietro, W. J.; Hehre, W. J.; Binkley, J. S.; DeFrees, D. J.; Pople, J. A.; Gordon, M. S. J. Chem. Phys. 1982, 77, 3654−3665. (24) Becke, A. D. J. Chem. Phys. 1996, 104, 1040−1046. (25) Cao, X.; Fischer, G. Chem. Phys. 2000, 255, 195−204. (26) Frisch, M. J.; Head-Gordon, M.; Pople, J. A. Chem. Phys. Lett. 1990, 166, 281−289. (27) Head-Gordon, M.; Head-Gordon, T. Chem. Phys. Lett. 1994, 220, 122−128. (28) Foresman, J. B.; Head-Gordon, M.; Pople, J. A.; Frisch, M. J. J. Phys. Chem. 1992, 96, 135−149. (29) Spiro, T. G.; Stein, P. Annu. Rev. Phys. Chem. 1977, 28, 501− 521. (30) Lui, C. H.; Malard, L. M.; Kim, S.; Lantz, G.; Laverge, F. E.; Saito, R.; Heinz, T. F. Nano Lett. 2012, 12, 5539−5544. (31) Ling, X.; Liang, L.; Huang, S.; Puretzky, A. A.; Geohegan, D. B.; Sumpter, B. G.; Kong, J.; Meunier, V.; Dresselhaus, M. S. Nano Lett. 2015, 15, 4080−4088. (32) Subotnik, J. E. J. Chem. Phys. 2011, 135, 071104. (33) Liu, X.; Fatehi, S.; Shao, Y.; Veldkamp, B. S.; Subotnik, J. E. J. Chem. Phys. 2012, 136, 161101. (34) Markel, F.; Ferris, N. S.; Gould, I. R.; Myers, A. B. J. Am. Chem. Soc. 1992, 114, 6208−6219.

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jctc.6b00514. Cartesian coordinates for the optimized structures of benzene, tetraalanine and tetracyanoethylene-mesitylene complexes and tables with the expansion of benzene total normal modes in terms of local oscillators (PDF)

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Funding

This work was supported by the project ANR-JCJC BIOMAGNET, and it was granted access to the HPC resources of Aix-Marseille Université financed by the project Equip@ Meso (ANR-10-EQPX-29-01) of the program “Investissements d’Avenir” supervised by the “Agence Nationale de la Recherche”. Notes

The authors declare no competing financial interest.

■ ■

ACKNOWLEDGMENTS M.H.-R. thanks I. Burghardt and H. Tamura for helpful discussions. REFERENCES

(1) Born, M.; Oppenheimer, J. R. Ann. Phys. 1927, 389, 457−484. (2) Domcke, W., Yarkony, D. R., Kö ppel, H., Eds. Conical Intersections: Theory, Computation and Experiment; Advanced Series in Physical Chemistry; World Scientific, 2011; Vol. 17. (3) Köppel, H., Yarkony, D. R., Barentzen, H., Eds. The Jahn-Teller Effect: Fundamentals and Implications for Physics and Chemistry; Springer Series in Chemical Physics; Springer, 2009; Vol. 97. (4) Wilson, E. B.; Decius, J. C.; Cross, P. C. Molecular Vibrations; Dover, 1995. (5) Eckart, C. Phys. Rev. 1935, 47, 552−558. (6) Sayvetz, A. J. Chem. Phys. 1939, 7, 383−389. (7) Mikhonin, A. V.; Ahmed, Z.; Asher, A. I.; Ianoul, S. A. J. Phys. Chem. B 2004, 108, 19020−19028. (8) Breda, S.; Reva, I. D.; Lapinski, L.; Nowak, M. J.; Fausto, R. J. Mol. Struct. 2006, 786, 193−206. (9) Falke, S. M.; Rozzi, C. A.; Brida, D.; Maiuri, M.; Amato, M.; Sommer, E.; De Sio, A.; Rubio, A.; Cerullo, G.; Molinari, E.; Lienau, C. Science 2014, 344, 1001−5. (10) Jaffé, C.; Brumer, P. J. Chem. Phys. 1980, 73, 5646. (11) Child, M. S. Acc. Chem. Res. 1985, 18, 45−50. (12) Duncan, J. L. Chem. Phys. Lett. 1988, 145, 347−353. (13) Konkoli, Z.; Cremer, D. Int. J. Quantum Chem. 1998, 67, 1−9. (14) Konkoli, Z.; Larsson, J. A.; Cremer, D. Int. J. Quantum Chem. 1998, 67, 11−27. (15) Konkoli, Z.; Cremer, D. Int. J. Quantum Chem. 1998, 67, 29−40. (16) Konkoli, Z.; Larsson, J. A.; Cremer, D. Int. J. Quantum Chem. 1998, 67, 41−55. (17) Jacob, C. R.; Reiher, M. J. Chem. Phys. 2009, 130, 084106. (18) Cormen, T. H.; Leiserson, C. E.; Rivest, R. L.; Stein, C. Introduction to Algorithms; The MIT Press, 2003; Chapter 28, p 760. (19) Szalay, V. J. Chem. Phys. 2014, 140, 234107. (20) Ochterski, J. W. White paper: Vibrational Analysis in Gaussian; Gaussian, 1999; p 10 4777


An Effective Procedure for Analyzing Molecular Vibrations in Terms of

Recommend Documents