Strong Correlation in Acene Sheets from the Active-Space Variational

ARTICLE pubs.acs.org/JPCA

Strong Correlation in Acene Sheets from the Active-Space Variational Two-Electron Reduced Density Matrix Method: Effects of Symmetry and Size Kenley Pelzer,† Loren Greenman,†,‡ Gergely Gidofalvi,§ and David A. Mazziotti*,† †

Department of Chemistry and The James Franck Institute, The University of Chicago, Chicago, Illinois 60637, United States Department of Chemistry, University of California, Berkeley, California 94720, United States § Department of Chemistry, Gonzaga University, Spokane, Washington 99258, United States ‡

ABSTRACT:

Polyaromatic hydrocarbons (PAHs) are a class of organic molecules with importance in several branches of science, including medicine, combustion chemistry, and materials science. The delocalized π-orbital systems in PAHs require highly accurate electronic structure methods to capture strong electron correlation. Treating correlation in PAHs has been challenging because (i) traditional wave function methods for strong correlation have not been applicable since they scale exponentially in the number of strongly correlated orbitals, and (ii) alternative methods such as the density-matrix renormalization group and variational twoelectron reduced density matrix (2-RDM) methods have not been applied beyond linear acene chains. In this paper we extend the earlier results from active-space variational 2-RDM theory [Gidofalvi, G.; Mazziotti, D. A. J. Chem. Phys. 2008, 129, 134108] to the more general two-dimensional arrangement of rings—acene sheets—to study the relationship between geometry and electron correlation in PAHs. The acene-sheet calculations, if performed with conventional wave function methods, would require wave function expansions with as many as 1.5 1017 configuration state functions. To measure electron correlation, we employ several RDM-based metrics: (i) natural-orbital occupation numbers, (ii) the 1-RDM von Neumann entropy, (iii) the correlation energy per carbon atom, and (iv) the squared Frobenius norm of the cumulant 2-RDM. The results confirm a trend of increasing polyradical character with increasing molecular size previously observed in linear PAHs and reveal a corresponding trend in two-dimensional (arch-shaped) PAHs. Furthermore, in PAHs of similar size they show significant variations in correlation with geometry. PAHs with the strictly linear geometry (chains) exhibit more electron correlation than PAHs with nonlinear geometries (sheets).

I. INTRODUCTION Polyaromatic hydrocarbons (PAHs) are a class of molecules of significant interest in many areas of science. Because many PAHs are carcinogenic, they present health risks when humans are exposed through dietary sources,1 cigarette smoke,2 or soot.3 PAHs are of interest in astrophysics, where they have been described as “abundant, ubiquitous, and a dominant force in the interstellar medium of galaxies”.4 They are also important as molecular semiconductors with applications as field-effect transistors in electronic displays.5 Two-dimensional PAHs provide a finite approximation to graphene,6 which has recently attracted significant attention for its electronic and optical properties. Accurate electronic structure calculations can contribute significantly to our understanding of the chemistry of PAHs, but recent studies7,8 have shown that the delocalized orbitals in the aromatic system exhibit correlation effects that cannot be captured with most electronic structure methods. Without accurate treatment of correlation effects, calculations r 2011 American Chemical Society

cannot capture important chemical traits of these molecules such as polyradical character or the singlettriplet energy gap. The works of Hachmann et al.,7 Gidofalvi et al.,8 and others913 have studied electron correlation in linear chains of rings, but there have not been any similar studies that investigate strong correlation in PAHs with a two-dimensional layout of the rings. The goal of this study is to apply two-electron reduced density matrix (2-RDM) methods to twodimensional PAHs to, which provide information on the relationship between geometry and strong correlation in PAHs. Although the correlation energy represents a fairly small fraction of the total energy, it can be of great importance in explaining chemical phenomena.14 Strong correlation refers to electron correlation that arises when two or more electron configurations Received: February 21, 2011 Revised: April 28, 2011 Published: May 12, 2011 5632

dx.doi.org/10.1021/jp2017192 | J. Phys. Chem. A 2011, 115, 5632–5640

The Journal of Physical Chemistry A contribute significantly to the wave function at zeroth order of many-body perturbation theory. Multireference wave functionbased methods are able to capture part of this strong correlation by treating the reference wave function as a linear combination of two or more Slater determinants. Because the coefficients in this linear combination are determined by a diagonalization, the cost of these calculations scales exponentially with the number of strongly correlated orbitals.14 When all orbitals are treated as strongly correlated, this method provides a full solution of the Schr€odinger equation in a finite orbital basis set known as a full configuration interaction (FCI). FCI is prohibitively expensive for all but the smallest systems. A multireference method with a more manageable cost than FCI is the configuration interaction complete active-space self-consistent field (CI-CASSCF) technique.15 This method (i) selects an initial active space of the orbitals and electrons that are expected to exhibit the greatest multireference electron correlation, (ii) optimizes the many-electron wave function by a linear combination of determinants constructed from the active orbitals, (iii) minimizes the energy further by performing orbital rotations between the active and inactive spaces to generate an improved set of active orbitals, and (iv) repeats steps (ii) and (iii) until convergence of the energy and wave function.14 This replaces the scaling of FCI rN by rN a , where r is the total number of orbitals and ra is the number of active orbitals. While CICASSCF methods are less expensive than FCI, they are still too computationally costly to manage PAHs with more than four rings.8 Because there are only seven possible geometric arrangements of four rings, CI-CASSCF cannot treat a wide enough variety of structures to permit a thorough examination of geometry and strong correlation. The variational 2-RDM methods offer a cost-effective alternative to multireference wave function-based methods by minimizing the energy as a function of the 2-RDM rather than the wave function8,1644 The two-particle RDM is expressed as Z D2 ð12, 10 20 Þ ¼ Ψð123:::NÞΨð10 20 3:::NÞdð3:::NÞ ð1Þ or equivalently, using second quantization notation and an orbital representation, its elements are defined as 2 i, j Dk, l ¼ Æψjâ†i â†j âlâk jψæ ð2Þ where the operator a†i (ai) creates (annihilates) an electron in spin orbital i. These calculations use N-representability conditions to constrain the reduced density matrices to correspond to an N-electron wave functions.16,2224,4548 Because these methods avoid the expensive diagonalization of wave function-based multireference methods, they are computationally more affordable, and with a matrix factorization to improve efficiency, the cost of 2-RDM methods scales as r6 rather than rN.26 Like CICASSCF calculations, 2-RDM calculations can improve their efficiency further by calculating strong correlation only within the active space.8 In this case, the CI treatment is replaced by a variational 2-RDM calculation in the active space, thus reducing 6 scaling from rN a to ra. Applications of active-space 2-RDM methods have been made to acene and aryne chains,8,43 the metal-to-insulator transition in hydrogen chains and lattices,26,43 and the conical intersections of firefly luciferin.43 The relatively low cost of the variational 2-RDM methods allows them to treat PAHs with up to eight rings. Since the

ARTICLE

number of possible geometric permutations increases with the number of rings (such that eight-ring structures have hundreds of possible permutations), 2-RDM methods can gather data on a wide range of geometries. Density-matrix renormalization group (DMRG) calculations are also a cost-effective alternative to wave function-based methods, and by exploiting spatial locality, they have been used to treat PAHs with up to 12 rings.7 However, some arrangements of rings are problematic for DMRG calculations. Although all molecules are in reality three-dimensional, molecules with a linear arrangement of the rings (a straight chain of rings) can be for the present purposes labeled as “onedimensional.” DMRG, especially when exploiting spatial locality, works best with “one-dimensional” molecules, with greater difficulty treating “two-dimensional” molecules49 such as a chain of rings bent into an arch shape. To date, DMRG has only been used to treat linear PAHs. In contrast, the 2-RDM methods in ref 8, which do not exploit spatial locality, are much less sensitive to the arrangement of the rings in the PAHs. Hence, by applying the 2-RDM methods to two-dimensional arrangements of PAHs, we can elucidate the relationship between geometry— both dimensionality and point-group symmetry—and electron correlation. In this paper we apply active-space variational 2-RDM methods8,32,43 to the study the lowest singlet states of 28 PAHs with sizes ranging from three to eight rings and a wide variety of geometries. The largest structures would require traditional wave function calculations with as many as 1.5 1017 configuration state functions. Our results confirm previous findings on the increase in strong electron correlation with increasing molecular size in linear chains of rings, and we extend this finding by noting a similar trend in arch-shaped chains of rings. When structures with the same number of rings are compared to one another, we find major differences in the degree of electron correlation between molecules with differing geometries. Linear molecules in particular are distinct from the remainder of the molecules, consistently showing greater electron correlation than molecules with more two-dimensional geometries. Our results suggest that electron correlation effects are not only a function of size but also a function of geometry in aromatic systems.

II. THEORY After a discussion of the variational 2-RDM method in section IIA, we describe briefly in section IIB unique features of the method when it is applied to an active space. Finally, in section IIC we define several 1- and 2-RDM-based measures of electron correlation that we will employ in assessing strong correlation in the PAHs. A. 2-RDM Methods. In 2-RDM methods, the energy of a many-electron system can be found by its minimization as a function of the 2-RDM E¼

2 i, j 2 i, j Kk, l Dk, l ∑ i, j, k, l

ð3Þ

where 2Ki,jk,l are the matrix elements associated with the reduced ^ in a finite one-electron basis set Hamiltonian operator 2K ! Zj 1 1 2 1 1 2^ r1 ð4Þ K ¼ þ N 1 2 2 r12 r j 1j

∑

and 2Di,jk,l are the elements of the 2-RDM defined in eq 2. 5633

dx.doi.org/10.1021/jp2017192 |J. Phys. Chem. A 2011, 115, 5632–5640

The Journal of Physical Chemistry A

ARTICLE

Figure 1. Structures of all molecules considered in this work with the number and letter used to identify each molecule in the text. Carbon atoms are shown in red, and hydrogen atoms are shown in blue.

Because not all two-electron density matrices correspond to an N-electron wave function, the 2-RDM must be constrained to represent an N-electron system by nontrivial conditions, known as N-representability conditions. The search for these conditions, called the N-representability problem, can be addressed through p-positivity conditions22,29 where stringency increases rapidly with p. Both 2-positivity and 3-positivity give results that capture a large fraction of the correlation energy,2325,29,3234,45,46 with 2-positivity having a lower computational cost. The 2-positivity conditions restrict three forms of the 2-RDM to be positive semidefinite: 2D, the two-particle RDM whose elements are defined in eq 2, and 2Q and 2G, the two-hole and one-particle/ one-hole reduced density matrices, respectively, whose elements

are given by i, j Qk, l ¼ Æψjâiâjâ†l â†k jψæ

ð5Þ

i, j Gk, l ¼ Æψjâ†i âjâ†l âk jψæ

ð6Þ

2

2

A matrix is positive semidefinite if and only if its eigenvalues are non-negative. Intuitively, the constraint that the eigenvalues must be non-negative corresponds to the fact that probabilities cannot be negative. Constraining 2D to be positive semidefinite implies a non-negative probability of finding two electrons, while constraining 2G to be positive semidefinite implies a non-negative probability of finding one electron and one hole. Applying 5634



ARTICLE

Table 1. All Electron-Correlation Metrics for the Structures in Figure 1a occupation numbers number of rings 3

4

5

6

8

structure

correlation metrics

nHONO

nLUNO

nHONO nLUNO

TVNE

||2Δ||2F

CE/carbon

3a

1.785

0.220

1.565

1.292

1.101

16.6

3c

1.839

0.164

1.674

1.248

1.054

16.1

3b

1.851

0.166

1.685

0.979

0.706

13.7

4a

1.715

0.293

1.422

1.763

1.570

17.1

4c

1.781

0.270

1.511

1.650

1.351

17.2

4g

1.794

0.211

1.583

1.682

1.495

16.6

4f

1.799

0.208

1.591

1.586

1.361

16.8

4e 4b

1.821 1.826

0.183 0.186

1.638 1.640

1.651 1.651

1.447 1.464

16.4 16.4

4d

1.843

0.160

1.683

1.645

1.484

16.2

5a

1.632

0.379

1.253

2.138

2.087

17.6

5f

1.749

0.257

1.492

1.930

1.861

17.0

5d

1.793

0.210

1.584

1.991

1.981

16.8

5b

1.804

0.216

1.588

1.894

1.788

16.9

5h

1.801

0.205

1.595

1.957

1.948

16.6

5c 5e

1.804 1.833

0.202 0.215

1.603 1.619

1.953 1.793

1.929 1.588

16.6 17.0

5g

1.820

0.182

1.638

1.937

1.898

16.6

5i

1.820

0.182

1.639

1.936

1.879

16.6

6b

1.095

0.916

0.178

2.670

2.335

21.0

6a

1.542

0.472

1.070

2.640

2.625

18.0

6c

1.728

0.281

1.447

2.253

2.174

17.5

6d

1.758

0.263

1.495

2.415

2.451

17.0

6f 6e

1.818 1.816

0.204 0.198

1.614 1.617

2.205 2.340

2.297 2.374

16.6 16.7

8a

1.364

0.655

0.709

3.522

3.969

18.6

8b

1.572

0.437

1.136

3.005

3.189

18.2

8c

1.806

0.221

1.586

2.933

3.170

16.4

The molecules are organized first by the number of rings and then the occupation-number gap nHONO nLUNO. We show the occupation numbers nHONO and nLUNO, their difference nHONO nLUNO, the TVNE, the squared Frobenius norm, and the correlation energy (CE) per carbon atom (in mhartrees). a

2-positivity constraints to a 2-RDM calculation results in an energy that is a lower bound to the exact energy in the finiteorbital basis set. B. Active-Space 2-RDM Methods. Improvements in the scaling of the variational 2-RDM methods can be realized by using an active-space approach in which the creation and annihilation operators in the 2-RDM elements treat only certain selected “active-space” orbitals.8 This is the theoretical analogue of CICASSCF methods, but while CI-CASSCF performs a costly full CI calculation in the active space, our calculations treat the active space with the polynomially scaling 2-RDM method. In the active-space variational 2-RDM method, which we will also denote as RDM-CASSCF, the energy and 2-RDM are optimized by the following steps: (i) selection of an initial active space of the orbitals and electrons that are expected to exhibit the greatest multireference electron correlation, (ii) minimization of the energy as a functional of the part of the 2-RDM in the active space, (iii) further minimization of the energy from mixing active and inactive orbitals to generate an improved set of active orbitals, and (iv) repetition of steps (ii) and (iii) until the energy and the 2-RDM converge. In the present study the initial set of active orbitals was selected from the n lowest-energy π-orbitals,

where n is the number of carbons participating in the delocalized π-orbital system. Since previous studies have found polyradical character in the π-orbital system,7,8 we expect significant electron correlation in these orbitals. These orbitals are also of particular interest in extending our results to predict correlation effects in other molecules with delocalized π-orbitals. C. Metrics of Strong Correlation. Strong correlation in the PAHs is measured by the following metrics: (i) the occupation numbers of the natural orbitals, (ii) the von Neumann entropy of the 1-RDM,50 (iii) the correlation energy, and (iv) the squared Frobenius norm of the cumulant part of the 2-RDM.50,51 While the first two metrics are based on the 1-RDM, the latter two metrics are based on the 2-RDM. 1. Natural Occupation Numbers. The 1-RDM can be derived from the 2-RDM and diagonalized to produce eigenvalues and eigenfunctions, known as the natural occupation numbers and natural orbitals, respectively. While single-reference electronic structure methods constrain the occupation number of an orbital to be either 0 or 2, 2-RDM calculations are free from these constraints and can capture fractional occupation numbers, thus assessing the degree of radical character in a molecule. Of greatest interest are the occupation numbers of the (N/2)th and (N/2 þ 1)th 5635



ARTICLE

entropy (TVNE) in which the largest and smallest occupation numbers are removed from the summation to give all molecules with the same number of rings an identical number of entropycontributing occupation numbers. Although the TVNE has the flaw that it disregards information from the highest and lowest natural orbitals, the variation between molecules for the very low or very high occupations is not large. It will also prove useful to assess the von Neumann entropy of the orbitals in the TVNE except the HONO and LUNO, which we call the modified truncated von Neumann entropy (MTVNE). 3. Correlation Energy per Carbon. The correlation energy is the difference between the final energy calculated by the 2-RDM method and the energy from a restricted HartreeFock calculation. To make this metric comparable between any two molecules, the total correlation energy is divided by the number of carbons to obtain the correlation energy per carbon. 4. Frobenius Norm of 2Δ. The magnitude of correlation effects can be measured by the Frobenius norm of the cumulant component 2Δ of the 2-RDM5052 2

Δ ¼ 2 D1 D∧1 D

ð8Þ

Since 1D ∧ 1D represents the mean-field solution of the problem, which does not incorporate strong correlation, subtracting it from 2D results in a metric of the correlation effects in the 2-RDM solution. The squared Frobenius norm assesses the magnitude of the cumulant 2-RDM as follows: Δ

2 F

) )

2

¼ Trð2 Δ† 2 ΔÞ

ð9Þ

The squared Frobenius norm, a size-extensive norm, can be compared between any two molecules.

III. APPLICATIONS Figure 2. Natural-orbital occupation numbers for the (a) linear structures with three, four, five, six, or eight rings and (b) arch-shaped structures with three, four, five, or six rings. The arch-shaped structures are molecules 3c, 4b, 5c, and 6d shown in Figure 1.

orbitals, respectively. Because these orbitals are the naturalorbital analogues of the highest occupied (HOMO) and lowest unoccupied (LUMO) molecular orbitals frequently discussed in chemistry, we will also designate them as highest occupied (HONO) and lowest unoccupied (LUNO) natural orbitals. Occupation numbers can be compared between any two molecules regardless of molecular size. 2. von Neumann Entropy. The information that is contained in the entire set of occupation numbers can be captured in a single metric by calculating the von Neumann entropy (VNE) of the 1-RDM VNE ¼

∑i 2i ln 2i n

n

ð7Þ

where ni is the ith eigenvalue of the spin-averaged 1-RDM lying between 0 and 2. A problem with utilizing the von Neumann entropy as a measure of correlation effects is that it cannot be compared between molecules with a different number of active space orbitals (and hence a different number of elements in the summation), and in our calculations, the number of active space orbitals differs based on the number of carbon atoms in a molecule. We, therefore, also define a truncated von Neumann

A. Methods. Prior to the active-space 2-RDM (RDM-CASSCF) calculations, we obtained an initial set of orbitals from restricted HartreeFock calculations in the DunningHay basis set with the program GAMESS.53 We utilized an [n,n] active space with the n lowest energy π-orbitals serving as our active space. Additional details are provided with the calculations on acene chains by Gidofalvi and Mazziotti.8 All molecules studied are depicted in Figure 1. All possible geometrical arrangements of rings were treated for the three- and four-ring cases. For the five- and six-ring cases, because there were too many possible arrangements of rings to study all of them, we chose nine five-ring and six six-ring molecules for the purpose of studying a range of geometries. In the eight-ring case, we examined three molecules with D2h symmetry to study variations within the same point-group symmetry. Structures 3b, 4c, 5b, and 5e employ a five-membered ring to create a singlet state while all other structures contain only six-membered rings. B. Results. In section IIIB1 we examine the relationship between size and strong correlation in linear and arch-shaped PAHs, and in section IIIB2 we study the relationship between geometry and strong correlation in molecules of a fixed size. All numerical data is presented in Table 1. 1. Size Effects in Linear and Arch-Shaped Structures. Our occupation number data confirms the findings of Gidofalvi and Mazziotti8 and Hachmann et al.7 that linear PAHs develop greater radical character as the number of rings increases. The HONO/ LUNO gaps for the three-, four-, five-, six-, and eight-ring 5636



ARTICLE

Figure 3. Electron density (96% contour) of the HONOs and LUNOs for the three-ring structures. The occupation number of each orbital is given below its image. Structures are ordered from left to right by increasing HONO/LUNO occupation number gap (or decreasing radical character).

Figure 4. Electron density of the HONOs and LUNOs for the six-ring structures. The occupation number of each orbital is given below its image. Structures are ordered from left to right by increasing HONO/LUNO occupation number gap (or decreasing radical character). 99% of the electron density is shown for the six-ring structures. A greater percentage of the density was chosen for the six-ring structures to adjust for a greater spread of the electron density in a larger molecule, making it difficult to compare orbital images for molecules of different sizes.

structures are 1.565, 1.422, 1.253, 1.070, and 0.709, respectively. Figure 2a plots this relationship for the three-, four-, five-, six-,

and eight-ring linear structures, showing the occupation numbers for all of the natural orbitals. The decrease in orbital occupation 5637



ARTICLE

Figure 5. Molecular-orbital diagrams derived from the H€uckel theory analysis for the six-ring triangular, superbenzene, linear, and S-shaped structures. The red line indicates the border between occupied and unoccupied orbitals, assuming that orbitals are filled according to the Aufbau principle. In the case of the triangular structure, the red line crossing two degenerate orbitals indicates that two electrons are present at an energy level containing two degenerate orbitals, hence the HOMO and LUMO are degenerate in the triangular case. This result supports the high radical character of the HONO and LUNO in the 2-RDM calculations.

between the HONO and LUNO orbitals becomes less sharp as the chain lengthens and the degree of strong correlation increases. While the difference in occupation numbers is greatest when comparing HONO and LUNO orbitals, the tendency toward greater polyradical character in larger molecules manifests itself in all of the natural orbital occupation numbers. Furthermore, all of the other metrics for strong electron correlation confirm the increase in correlation effects with an increasing number of rings for linear PAHs. To explore whether this trend also appears in molecules with a nonlinear or two-dimensional arrangement of rings, we examine arch-shaped PAHs. A similar shape with C2v symmetry was calculated for the three-, four-, five-, and six-ring structures, allowing for a meaningful comparison between molecules of different sizes. Figure 2b plots the occupation numbers of all orbitals for the arch-shaped structures, showing, as in the linear case, a trend of increasing polyradical character with increasing molecular size. The differences between the HONO and LUNO occupation numbers for the arch-shaped three-, four-, five-, and six-ring structures are 1.674, 1.640, 1.603, and 1.495, respectively. All other metrics also clearly indicate a significant increase in electron correlation with system size. 2. Geometry Effects. Although both one- and two-dimensional structures in Figure 2a,b show the trend of increasing polyradical character with increasing size, the trend is more pronounced in the linear structures. Furthermore, a comparison of structures with the same number of rings but differing geometries in Table 1 reveals that the linear structures tend to have the greatest polyradical character. Figures 3 and 4 show images of the HONO and LUNO with their occupation numbers for the three- and sixring structures, respectively. Even though in the six-ring case the triangular structure (6b) has significantly more radical character relative to the linear case, with an occupation number gap (nHONO nLUNO) of 0.178 relative to 1.070 for the linear molecule, the six-ring linear acene exhibits a smaller occupation number gap than all of the remaining six-ring structures. To understand the significant differences between the linear, triangular, and superbenzene structures (all of which have notably high symmetry, at D2h, D3h, and D6h, respectively), we used H€uckel theory to obtain an estimate of relative orbital

energies. H€uckel parameters of R = 6.6 eV and β = 2.7 eV were chosen.54 Molecular-orbital diagrams for the six-ring linear, triangular, superbenzene, and S-shaped structures are shown in Figure 5. Each horizontal blue line represents a spatial orbital; the relative position of the lines indicates relative energies with energies increasing along the vertical axis. Multiple lines at the same level reflect degenerate orbitals. The red dotted line divides the occupied orbitals from the unoccupied orbitals, assuming that each spatial orbital contains either zero, one, or two electrons and the orbitals are filled according to the Aufbau principle. The highest symmetry structures, the D3h triangle and the D6h superbenzene, show energy level degeneracies. The triangular structure is noteworthy in that the HOMO and LUMO orbitals are degenerate, as reflected by the red line crossing two degenerate orbitals. Filling the orbitals with electrons by the Aufbau principle and Hund’s rule is consistent with the near radical character of the HONO and LUNO found by the 2-RDM analysis with occupation numbers of 1.095 and 0.916 for the HONO and LUNO, respectively. The superbenzene structure, although it has many orbital degeneracies, does not have a degeneracy between the HOMO and LUMO, and in fact the gap between the HOMO and LUMO is large relative to the gaps between the other energy levels. In the language of organic chemistry, superbenzene is an even alternant aromatic hydrocarbon while the triangular structure is an odd alternant aromatic hydrocarbon. A conjugated hydrocarbon is alternant when the atoms in the conjugation can be divided into two groups where each atom in the first group has only atoms in the second group as its neighbors;55 when the numbers of atoms in each group are equal, the hydrocarbon is even alternant, and otherwise, it is odd alternant. It can be shown that odd alternant hydrocarbons have degenerate HOMO and LUMO orbitals while even alternant hydrocarbons do not. Superbenzene’s H€uckel result is consistent with the relatively large HONO/LUNO occupation number gap (1.614) found in the 2-RDM calculations. The H€uckel molecular orbital diagrams for the linear and S-shaped structures are also in agreement with the occupation numbers from the 2-RDM calculations. The correlation energy per carbon shows the same trend as the HONO/LUNO gap for the six-ring linear, triangular, and superbenzene structures with energies of 21.0, 18.0, and 16.6 mHs, respectively. 5638



ARTICLE

sextet rule.5659 To apply Clar’s sextet rule, the π-electrons of a PAH are grouped into sextets as much as possible; PAHs with a larger number of sextets are expected to be more stable and less reactive. Most of our molecules have only one or two sextets, but we find that the molecules with three sextets (4d, 5c, 5i, 6e, and 8c) have relatively large HONO/LUNO occupation number gaps. This is consistent with the prediction of large band gaps and low reactivity of these molecules. Another trend between geometry and reactivity that has been noted elsewhere is the tendency of PAHs with “zigzag” edges to be less stable than PAHs with “armchair” edges.56,60,61 This also is consistent with our results; the vast majority of our molecules have a structure that is at least partially surrounded by zigzag edges, but 5i, a molecule with pure armchair edges, has the largest HONO/ LUNO occupation number gap of the five-ring structures. Figure 6. Occupation numbers for all orbitals of the six-ring linear, triangular, and superbenzene structures.

In Figure 6 we plot the occupation numbers of all the natural orbitals for the linear, triangular, and superbenzene structures. For most of the orbitals the differences are small with only the HONO and LUNO orbitals showing a striking difference. In assessing polyradical character, we find that the six-ring triangular structure has the greatest TVNE at 2.670, followed by the linear structure (2.640) and the superbenzene structure (2.205). To determine whether the correlation pattern persists in orbitals other than the HONO and LUNO, we can compute the MTVNE defined in section IIC. Of the six-ring structures, the linear structure has the greatest MTVNE (2.099), while the triangular structure has an MTVNE of only 1.982. Structures 6c, 6d, 6e, and 6f have MTVNEs of 1.851, 2.034, 2.023, and 1.885, respectively. Thus it appears that the high TVNE and correlation energy of the triangular structure reflect the unusual situation of a degenerate HONO and LUNO, and, except for this degeneracy, the six-ring structures are not an exception to the rule that linear molecules tend to manifest greater correlation effects. Unlike the von Neumann entropy, the squared Frobenius norm of 2Δ also indicates greater correlation effects in the linear structure, the linear structure having the largest squared Frobenius norm at 2.625, followed by the triangular structure (2.335) and superbenzene (2.297). In interpreting the results of the three-, four-, five-, and six-ring structures, it is difficult to disentangle the effects of varying spatial symmetry and varying dimensionality. To remove symmetry from the problem and elucidate the role of dimensionality, we studied three eight-ring structures of D2h symmetry. The HONO and LUNO occupation numbers of the eight-ring structures are shown in Table 1. In this case, the linear structure (8a) once again has the smallest HONO/LUNO occupation number gap at 0.709, relative to 1.136 and 1.586 for structures 8b and 8c, respectively. The TVNE, Frobenius norm, and correlation energy per carbon follow the same trend. Although the linear structure exhibits the greatest correlation effects, these calculations reveal the interesting result that molecule 8b, with a more two-dimensional arrangement of rings than 8c (meaning that the rings of 8b have a length-to-width ratio closer to one than those of 8c), shows greater correlation effects. It is worth commenting on the consistency of our results with patterns in shape and chemical behavior noted elsewhere. First, our results are consistent with the patterns predicted by Clar’s

IV. DISCUSSION AND CONCLUSIONS Although past electronic structure calculations of linear PAHs indicated an important role of strong correlation in determining the chemistry of these molecules, no data existed on strong correlation effects in PAHs with a nonlinear, two-dimensional arrangement of the rings in a plane. This dearth of information was due to the limitations of most electronic structure methods in treating large PAHs; prohibitive costs and difficulties treating two-dimensional geometries prevent other methods from assessing correlation in these molecules. Applying 2-RDM theory to PAHs of varying geometries has provided a unique opportunity to investigate the relationship between geometry and strong correlation effects. Using 2-RDM methods, this study evaluated electron correlation effects for a sample of 28 PAHs with widely varying sizes and geometries. Several metrics of correlation are calculated in this paper, and the fact that the results find very similar trends among the different metrics suggests that these are robust measures of strong correlation. Our calculations indicate that geometry plays an important role in determining strong correlation in PAHs. Because the PAHs that we studied differ in both symmetry and dimensionality, both of these factors must be considered in interpreting our findings: (i) Quantum theory states that systems with greater symmetry have more energetic degeneracies. Energetic degeneracies can lead to radical character, although the H€uckel molecular orbital diagrams make it clear that degeneracies must be in specific orbitals to have a significant effect on radical character. Irrespective of degeneracies, our results do not find that higher symmetry structures consistently show more electron correlation. The fact that the D6h superbenzene shows far less correlation than the D2h linear molecule suggests that symmetry alone cannot predict correlation effects. (ii) Dimensionality is a more effective predictor when comparing linear PAHs to nonlinear PAHs, in the sense that the linear (“one-dimensional”) structures consistently exhibit more correlation. However, dimensionality alone cannot explain variations in correlation effects observed between differing nonlinear structures. The case of the eight-ring structures demonstrates that a continuum interpretation of the relationship between dimensionality and correlation is not consistent with our results: in this case the longer, narrower (“less two-dimensional”) structure 8c shows less correlation than the roughly square-shaped structure 8b. 5639


The Journal of Physical Chemistry A By applying 2-RDM methods to PAHs of varying geometry, we find strong evidence that geometry plays a role in determining electron correlation. We find that linear PAHs consistently show greater electron correlation effects than more two-dimensional PAHs, except in unusual situations where high-symmetry molecules experience HOMO/LUMO degeneracy. It is also clear that nonlinear PAHs of differing geometries show differing degrees of strong electron correlation. The present results also provide a foundation for exploring the electronic effects of size and geometry in larger two-dimensional PAHs that approach graphene62 as well as three-dimensional PAHs. We conclude that in PAHs, and perhaps in other molecules containing delocalized π-orbitals, strong correlation is influenced by molecular geometry— symmetry and dimensionalityas well as molecular size.

’ AUTHOR INFORMATION Corresponding Author

*E-mail: [email protected].

’ ACKNOWLEDGMENT K.P. acknowledges support from the Krell-DOE Computational Science Graduate Fellowship Program. G.G. is supported in part by a grant to Gonzaga University from the Howard Hughes Medical Institute through the Undergraduate Science Education Program. D.A.M. gratefully acknowledges the NSF, ARO, Microsoft Corporation, the Dreyfus Foundation, and the David-Lucile Packard Foundation for support. ’ REFERENCES (1) Moret, S.; Conte, L. S. J. Chromatogr. A 2000, 882, 245–253. (2) Gray, N.; Kozlowski, L. Ann. Oncol. 2003, 14, 353–357. (3) Indarto, A. Env. Eng. Sci. 2009, 26, 1685–1692. (4) Tielens, A. G. G. M. Annu. Rev. Astron. Astrophys. 2008, 46, 289–337. (5) Geerts, Y.; Kl€arner, G.; M€ullen, K. In Electronic Materials: The Oligomer Approach; M€ullen, K., Kl€arner, G., Eds.; Wiley-VCH: Weinheim, Germany, 1998; p 48. (6) Geim, A. K.; Novoselov, K. S. Nat. Mater. 2007, 6, 183. (7) Hachmann, J.; Dorando, J.; Aviles, M.; Chan, G. K.-L. J. Chem. Phys. 2007, 127, 134309. (8) Gidofalvi, G.; Mazziotti, D. A. Phys. Rev. A 2008, 129, 134108. (9) Bendikov, M.; Duong, H. M.; Starkey, K.; Houk, K. N.; Carter, E. A.; Wudl, F. J. Am. Chem. Soc. 2004, 126, 7416. (10) dos Santos, M. C. Phys. Rev. B 2006, 74, 045426. (11) Qu, Z.; Zhang, D.; Liu, C.; Jiang, Y. J. Phys. Chem. A 2009, 113, 7909. (12) Hajgato, B.; Szieberth, D.; Geerlings, P.; Proft, F. D.; Deleuze, M. S. J. Chem. Phys. 2009, 131, 224321. (13) Ess, D. H.; Johnson, E. R.; Hu, X.; Yang, W. J. Phys. Chem. A 2011, 115, 76. (14) Jensen, F. Introduction to Computational Chemistry, 2nd ed.; Wiley & Sons: West Sussex, U.K., 2007. (15) Roos, B. O.; Taylor, P. R.; Siegbahn, P. E. M. Chem. Phys. 1980, 48, 157–173. (16) Reduced-Density-Matrix Mechanics with Application to ManyElectron Atoms and Molecules. In Advances in Chemical Physics; Mazziotti, D. A., Ed.; Wiley: New York, 2007; Vol. 134. (17) Coleman, A. J.; Yukalov, V. I. Reduced Density Matrices: Coulson’s Challenge; Springer: New York, 2002. (18) Mihailovic, M. V.; Rosina, M. Nucl. Phys. A 1975, 237, 221. (19) Garrod, C.; Mihailovic, V.; Rosina, M. J. Math. Phys. 1975, 10, 1975. (20) Erdahl, R. M. Rep. Math. Phys. 1979, 15, 147. (21) Erdahl, R. M.; Jin, B. Y. On calculating approximate and exact density matrices. In Many-Electron Densities and Reduced Density Matrices; Cioslowski, J., Ed.; Kluwer,: Boston, 2000; pp 5784.

ARTICLE

(22) Mazziotti, D. A.; Erdahl, R. M. Phys. Rev. A 2001, 63, 042113. (23) Nakata, M.; Nakatsuji, H.; Ehara, M.; Fukuda, M.; Nakata, K.; Fujisawa, K. J. Chem. Phys. 2001, 114, 8282. (24) Mazziotti, D. A. Phys. Rev. A 2002, 65, 062511. (25) Zhao, Z.; Braams, B. J.; Fukuda, H.; Overton, M. L.; Percus, J. K. J. Chem. Phys. 2004, 120, 2095. (26) Mazziotti, D. A. Phys. Rev. Lett. 2004, 93, 213001. (27) Gidofalvi, G.; Mazziotti, D. A. J. Chem. Phys. 2005, 122, 194104. (28) Cances, E.; Stoltz, G.; Lewin, M. J. Chem. Phys. 2006, 125, 064101. (29) Mazziotti, D. A. Phys. Rev. A 2006, 74, 032501. (30) Gidofalvi, G.; Mazziotti, D. A. Phys. Rev. A 2006, 74, 012501. (31) Mazziotti, D. A. Acc. Chem. Res. 2006, 39, 207–215. (32) Mazziotti, D. A. In Reduced-Density-Matrix Mechanics: With Application to Many-Electron Atoms and Molecules; Mazziotti, D. A., Ed.; Wiley: New York, 2007; Vol. 134, p 21. (33) Erdahl, R. M. In Reduced-Density-Matrix Mechanics: With Application to Many-Electron Atoms and Molecules; Mazziotti, D. A., Ed.; Wiley: New York, 2007; Vol. 134, p 61. (34) Braams, B. J.; Percus, J. E.; Zhao, Z. In Reduced-Density-Matrix Mechanics: With Application to Many-Electron Atoms and Molecules; Mazziotti, D. A., Ed.; Wiley: New York, 2007; Vol. 134, p 93. (35) Fukuda, M.; Nakata, M.; Yamashita, M. In Reduced-DensityMatrix Mechanics: With Application to Many-Electron Atoms and Molecules; Mazziotti, D. A., Ed.; Wiley: New York, 2007; Vol. 134, p 103. (36) Juhasz, T.; Shenvi, N.; Mazziotti, D. A. Chem. Phys. Lett. 2007, 445, 79. (37) Nakata, M.; Braams, B. J.; Fujisawa, K.; Fukuda, M.; Percus, J. K.; Yamashita, M.; Zhao, Z. J. Chem. Phys. 2008, 128, 164113. (38) Rothman, A. E.; Mazziotti, D. A. Phys. Rev. A 2008, 78, 032510. (39) Kamarchik, E.; Mazziotti, D. A. Phys. Rev. A 2007, 75, 013203. (40) Verstichel, B.; van Aggelen, H.; Van Neck, D.; Ayers, P. W.; Bultinck, P. Phys. Rev. A 2009, 80, 032508. (41) van Aggelen, H.; Verstichel, B.; Bultinck, P.; Van Neck, D.; Ayers, P. W.; Cooper, D. L. J. Chem. Phys. 2010, 132, 114112. (42) Nakata, M.; Yasuda, K. Phys. Rev. A 2009, 80, 5. (43) (a) Greenman, L.; Mazziotti, D. A. J. Chem. Phys. 2009, 130, 184101. (b) Sinitskiy, A.; Greenman, L.; Mazziotti, D. A. J. Chem. Phys. 2010, 133, 014104. (c) Greenman, L.; Mazz-iotti, D. A. J. Chem. Phys. 2010, 133, 164110. (44) Shenvi, N.; Izmaylov, A. F. Phys. Rev. Lett. 2010, 105, 213003. (45) Coleman, A. J. Rev. Mod. Phys. 1963, 35, 668. (46) Garrod, C.; Percus, J. J. Math. Phys. 1964, 5, 1756. (47) Erdahl, R. M. Int. J. Quantum Chem. 1978, 13, 697. (48) Harriman, J. E. Phys. Rev. A 1978, 17, 1257. (49) Hallberg, K. A. Adv. Phys. 2006, 44, 477. (50) Huang, Z.; Kais, S. Chem. Phys. Lett. 2005, 413, 1. (51) Juhasz, T.; Mazziotti, D. A. J. Chem. Phys. 2006, 125, 174105. (52) Mazziotti, D. A. Chem. Phys. Lett. 1998, 289, 419. (53) Schmidt, M. W.; Baldridge, K. K.; Boatz, J. A.; Elbert, S. T.; Gordon, M. S.; Jensen, J. H.; Koseki, S.; Matsunaga, N.; Nguyen, K. A.; Su, S.; Windus, T. L.; Dupuis, M.; Montgomery, J. A., Jr. J. Comput. Chem. 1993, 14, 1347–1363. (54) Brogli, F.; Heilbronner, E. Theor. Chim. Acta. 1972, 26, 289. (55) Coulson, C. A.; Rushbrooke, G. S. Proc. Cambridge Philos. Soc. 1940, 36, 193. (56) Rieger, R.; M€ullen, K. J. Phys. Org. Chem. 2010, 23, 315. (57) Clar, E. The Aromatic Sextet; Wiley: London, 1972. (58) Maksic, Z. B.; Baric, D.; M€uller, T. J. Phys. Chem. A 2006, 110, 10135. (59) Maksic, Z. B.; Baric, D.; M€uller, T. J. Phys. Chem. A 2009, 113, 1151. (60) Stein, S. E.; Brown, R. L. J. Am. Chem. Soc. 1987, 109, 3721. (61) Nakada, K.; Fujita, M.; Dresselhaus, G.; Dresselhaus, M. S. Phys. Rev. B 1996, 54, 17954. (62) Yan, X.; Cui, X.; Li, B.; Li, L. Nano Lett. 2010, 10, 1869.

5640


Strong Correlation in Acene Sheets from the Active-Space Variational

Recommend Documents