How Bond Length Alternation and Thermal Disorder Affect the Optical

Subscriber access provided by UNIV OSNABRUECK

Article

How bond length alternation and thermal disorder affect the optical excitation energies of #-conjugated chains: A combined DFT and MD study Juliana Bois, and Thomas Körzdörfer J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.5b01070 • Publication Date (Web): 09 Mar 2016 Downloaded from http://pubs.acs.org on March 11, 2016

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Journal of Chemical Theory and Computation is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 19

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

How bond length alternation and thermal disorder affect the optical excitation energies of π -conjugated chains: A combined DFT and MD study Juliana Bois and Thomas Körzdörfer∗ Institut für Chemie, Universität Potsdam, Karl-Liebknecht-Straße 24-25, 14476 Potsdam, Germany E-mail: [email protected]

Abstract

1

We dissect the sources of error leading to inaccuracies in the description of the geometry and optical excitation energies of π-conjugated polymers. While the ground-state bond length alternation is shown to be badly reproduced by standard functionals, the recently adapted functionals PBEh* and ωPBE* as well as the double hybrid functional XYGJ-OS manage to replicate results obtained at the CCSD(T) level. By analyzing the bond length alternation in the excited state, a sensitive dependence of the exciton localization on the long range behavior of functionals and the amount of HF exchange present is shown. Introducing thermal disorder through molecular dynamics simulations allows the consideration of a range of thermally accessible configurations of each oligomer, including trans to cis rotations which break the conjugation of the backbone. Thermal disorder has a considerable effect when combined with functionals that overestimate the delocalization of the excitation, such as B3LYP. For functionals with a larger amount of exact exchange such as our PBEh* and ωPBE*, however, the effect is small, as excitations are often localized enough to fit between twists in the chain.

The electronic properties of π-conjugated polymers are of great interest due to their wide use in optoelectronic applications. One of their most defining characteristics is the optical bandgap, which affects with which efficiency they can absorb sunlight. Theoretical investigations into the optical excitations of polymers are essential, as a profound understanding of the electronic processes involved and the relationship between electronic and structural properties is required to allow bottom up material design. 1 Due to their large size, first-principles electronic structure calculations on complete polymer systems are computationally not feasible. A simple solution to this problem is the well established oligomer approach, wherein properties are calculated for oligomers of increasing size and a fitting function based on a suitable conjugation length model is used to extrapolate to the polymer limit. 2–7 The smaller size of the oligomers then allows for the use of first principles methods. A popular method to estimate optical bandgaps along these lines is to calculate the vertical excitation energy for a number of oligomers, then use the Kuhn Fit, 6,8 based on a simple model of coupled harmonic oscillators, to extrapolate to the polymer limit. Using this model, the optical excitation energy EN extrap-

∗

To whom correspondence should be addressed

Introduction

ACS Paragon Plus Environment

1


Page 2 of 19

Vertical excitation energy (eV)

5

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 1: Structure of trans-polyacetylene, where N equals the number of double bonds. olates as

EN = E1

π ), 1 − Dk cos( N +1

4

3

1 0

where N is the number of double bonds along the chain, E1 is the transition energy of an isolated double bond and Dk is the coupling between double bonds by single bonds. 6 In order to yield a reliable extrapolation to the polymer limit even for π-conjugated polymers with a large conjugation length, the Kuhnfit approach still requires excited-state calculations for oligomers with several repeat units. Due to the size of these molecules, density functional theory (DFT) and it’s time-dependent extension, TD-DFT, are often the methods of choice for the calculation of the ground-state geometries and excited-state energies, respectively. However, as will be demonstrated below, this approach has a number of serious limitations. To discuss these limitations, we will here focus on the trans-polyacetylene (t-PA) molecule. t-PA is an established model system for π-conjugated polymers, as it allows to focus on the most fundamental problems with the description of the ground- and excited-state properties of π-conjugated systems. 9–11 At the same time, the relatively small number of electrons per repeat unit limits the computational cost, thus, allowing tractable calculations even for long oligomers. When calculating vertical excitation energies of t-PA oligomers using the popular exchangecorrelation functional B3LYP and comparing these to experimental data, 12 B3LYP poorly reproduces both the absolute vertical excitation energies and the slope of the Kuhn fit, as shown in Figure 2. The finding that B3LYP fails to replicate the experimental slope can be attributed to an overestimation of the effective conjugation length. 13,14 The use of standard exchange-correlation functionals such as

Experimental M06HF B3LYP

2

0.05

0.1

0.15

1/N

Figure 2: Vertical excitation energy of polyacetylene oligomers, with N=number of double bonds, for common DFT functionals compared to experimental data. 12 Kuhn fits are used to extrapolate to the polymer limit. B3LYP leads to a strongly overestimated delocalization of the excitation for long oligomer chains, resulting in a large decrease of excitation energy with increasing chain length. 2,13 In fact, for B3LYP the maximum conjugation length is not reached for the oligomer lengths considered here - the curve of the fit is purely due to the nature of the Kuhn Fit. The HFrich M06HF functional has been suggested as a useful alternative for the prediction of πconjugated polymer bandgaps, 14 as it reproduces the experimental slope more accurately than other functionals and with reduced deviation from experimental results (see Figure 2). In contrast to B3LYP, M06HF strongly localizes the excitation, leading to a much shorter maximum conjugation length. This means M06HF shows the expected saturation of the excitation energy for longer chain lengths, clearly visible in Figure 2. Looking at the drastically overestimated extent of bond length alternation (BLA) for M06HF (see Figure 3), however, there is reason to believe these results may not be as favorable as they seem. While a strong correlation between the BLA of a system and its polarizabilities 16 and absorption properties 17–19 is well established, the correct description of BLA in extended π-systems is a task many standard functionals still struggle with. 20–24 In contrast


2

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

BLA at central carbon (Å)

Page 3 of 19


We aim to capture a range of configurations through MD snapshots, to better replicate experimental conditions for theoretical calculations.

0.12 CCSD(T) M06HF B3LYP

2

0.08

2.1 0.04

DFT/TD-DFT

Ground- and excited-state geometry optimizations were performed using DFT and TDDFT, respectively, with a variety of functionals. M06HF (100% HF), 25 PBEh (25% HF), 26 B3LYP, (20% HF) 27 the range separated CAMB3LYP (19% HF in the short range, 65% HF in the long range), 28 the long range corrected 23 and the global LRC-ωPBE* (ω=0.188 a−1 0 ), hybrid PBEh* (45.73% HF), 23 both with parameters optimized for the BLA of polyenes, and the double hybrid B2PLYP (53% HF and 27% MP2), 29 which is based on B3LYP, were all used with Gaussian. 30 The double hybrid XYGJ-OS (77% HF, 43.6% PT2 correlation) 31 was used with QChem. 32 Grid sizes, accuracies and basis set (6-311G) of the two programs were kept identical. The chosen basis set has been shown to be adequate for the description of BLA in t-PA. 24 All optimizations were performed on linear, planar, all-trans oligomers with no side groups, leading to optimized structures with C2h point group. Vertical excitation energies were obtained using TD-DFT calculations and the functionals described above. The Tamm-Dancoff approximation (TDA) was applied for all TD-DFT calculations, as it has been shown to improve the description of the first singlet excited state of π-conjugated systems, 33 while removing triplet instabilities sometimes observed when using functionals with large amounts of exact exchange. 34,35 Fractional particle curves for the evaluation of the localization/delocalization error were performed using the PSI4 program. 36 Due to the need for diffuse functions on the anion side of the fractional particle curves, a jun-cc-pVTZ basis set 37 was employed for all fractional-electron calculations.

20

10

Computational Details

Number of double bonds

Figure 3: Evolution of the bond length alternation (BLA) with increasing polyacetylene oligomer size using different computational methods. The dotted line shows the experimetal value found for the polymer. 15 to M06HF for example, B3LYP starts of very close to higher level results for the BLA of short oligomers, but then falls off to greatly underestimate the BLA for longer chains. The problem is that without a correct description of the geometry, even excitation energies that seem to replicate experimental results well are unreliable and the methods not easily extendable to other systems. Our aim is to break down the difficulties encountered when describing π-conjugated polymers with DFT/TD-DFT, to facilitate a logical step by step improvement of the methods used to calculate their optical band gaps. Therefore we initially look at BLA and which functionals are able to correctly reproduce experimental and higher level theoretical reference data. We then focus on the localization/delocalization of the excitation, as this is vitally important to correctly reproduce the sloping off of excitation energies for long chain lengths found in experimental fits. Finally, we consider the effects of thermal disorder by means of molecular dynamics (MD) simulations of t-PA in solution. For experimental measurements in solvent, the tPA chains will be present in a wide range of configurations, many of them twisted to break the π-conjugation along the backbone. If enough configurations are twisted or disordered in a way that affects the conjugation length, this should be visible in the absorption spectrum.


3


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

2.2

MD

The OPLSAA force field 38,39 was modified to reproduce bond lengths, angles, dihedral angles, and torsional barriers from PBEh* DFT calculations for an oligomer with N=22. This oligomer length is large enough to have the same torsional barrier heights as the polymer. 40 A detailed description of the parameterization procedure, as well as all parameters used in the final force field can be found in the supporting information. Using these parameters for the oligomer, in combination with the well know TIP3P model for water, 41 a periodic solvent box including one oligomer was created and the energy minimized. After thorough equilibration including an NVT and an NPT step, the system was modelled for 25ns, using timesteps of 0.5 fs, a temperature of 298K (modified Berendsen thermostat) and pressure of 1.0 atm (Berendsen barostat). All steps were completed using Gromacs 5.0.1. 42 For the vertical excitation energies, 100 snapshots were taken at regular time intervals. The oligomer geometry obtained from these snapshots was then directly employed for the TD-DFT calculations, without any intermediate DFT-based geometry relaxation step.

3 3.1

Page 4 of 19

Figure 4: Evolution of the bond length alternation (BLA) with increasing polyacetylene oligomer size for additional functionals, including optimized functionals, long range corrected functionals and double hybrids.The dotted line shows the experimental value of BLA for the polymer solid. 15 underestimate it, it is possible to adjust the amount of exact exchange in any hybrid functional to minimize the error in BLA. By using this empirical procedure on the global hybrid functional PBEh, one finds the adapted variant PBEh* with 45.73% HF exchange. 23 As demonstrated in Figure 4, PBEh* shows excellent agreement with higher level and experimental results. Similarly, it is possible to adjust the range separation parameter in the longrange corrected hybrid functional LRC-ωPBE 47 to ensure a good description of the BLA. 23,45 This leads to the LRC-ωPBE* functional with a range-separation parameter of ω=0.188 a−1 0 . As shown in Figure 4, the BLA obtained from LRC-ωPBE* is in excellent agreement with the CCSD(T) data. A more fundamental approach to improving the description of the BLA is the inclusion of non-local correlation, 46 as present in the double hybrid functionals XYGJ-OS 31 and B2PLYP 29 (Figure 4). The degree of improvement gained, however, varies widely. 24 B2PLYP, which is based on B3LYP with added MP2 correlation, shows improvement over B3LYP, but still does not provide accurate results. XYGJ-OS on the other hand is able to replicate higher level cal-

Results DFT Geometry in the ground state

Most standard functionals fail dramatically at reproducing the ground state BLA of long polyenes. 20,22–24,43–46 This is demonstrated in Figure 4, where we compare the results from several exchange-correlation functionals with reference data obtained from coupled cluster theory at the coupled cluster singles doubles (triples) (CCSD(T)) level using a 6-31G(d) basis set, 23 which is generally accepted to yield reliable values for the BLA in t-PA. 20,22,43 Out of the standard hybrid functionals considered here, CAM-B3LYP comes closest to the required amount of HF exchange. The error of the predicted BLA is, however, still too large. As Hartree Fock (HF) overestimates BLA, and semilocal functionals such as LDAs and GGAs


4

Page 5 of 19

sideration when calculating the excitation energies of π-conjugated polymers. The comparison is undertaken via a two step process, involving geometry optimization using one functional and a subsequent TD-DFT calculation of the vertical excitation energy using another functional. The results are shown in Figure 5. For comparison, both parts of Figure 5 include the three curves obtained when using the same functional (M06HF, PBEh*, and B3LYP) for both the geometry optimization and TD-DFT. In addition to that, the left side of the figure shows PBEh* excitation energies calculated on top of M06HF and B3LYP geometries. The different geometries have a large impact on the final outcome. In comparison to the all PBEh* result, a downshift of -0.3 eV can be seen in the polymer limit when the geometries are calculated with the B3LYP functional. Similarly, an upshift of the excitation energy in the polymer limit by 0.5 eV is found when using M06HF geometries. The right side of Figure 5 additionally includes the Kuhn-Fits for the case where the excitation energies have been calculated using the M06HF or B3LYP functionals on the PBEh* geometries. Here a shift of -0.6 eV and +0.8 eV, for B3LYP and M06HF TD-DFT, respectively, can be observed compared to the all PBEh* curve. Although a larger shift is seen when changing the functionals used for the TD-DFT calculation and the error magnitudes vary depending on the functional used, this analysis reinforces the previously demonstrated dependence of the optical response of a system on its ground state BLA. 16,18 It follows that the errors introduced by using a wrong ground state geometry can be significant, and therefore should not be neglected. 14

culations best out of the functionals not specifically adapted for this problem. However, there are multiple difficulties when employing double hybrid functionals for excited state calculations. Most importantly, the added MP2/PT2 correlation can not be considered using TDDFT, but has to be evaluated using a suitable wavefunction based equivalent for excited states such as configuration singles (doubles) (CIS(D)). 48 However, as MP2 and CIS(D) are based on a different footing, it is not clear a priori whether such an approach would yield trustworthy results. Additionally, double hybrids are comparatively expensive and not widely implemented yet.

3.2

Comparison of error magnitudes Different geometries

Vertical Excitation Energy (eV)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Different TD-DFT functionals

5

4 3 2 M06HF PBEh* on M06HF PBEh* PBEh* on B3LYP B3LYP

1 0 0

0.05

0.1

1/N

0.15 0

M06HF on PBEh* B3LYP on PBEh*

0.05

0.1

0.15

1/N

Figure 5: Vertical excitation energies of polyacetylene oligomers of increasing size (N=number of double bonds) with Kuhn fits to extrapolate to the polymer limit. Right side: The mixed functional curves use B3LYP/M06HF for the geometry optimization step and PBEh* for the TD-DFT. Left side: The mixed functional curves use PBEh* geometries and B3LYP/M06HF for the TD-DFT step.

3.3

A comparison of the magnitudes of errors generated by the incorrect description of the ground-state geometry and the inaccuracies in the TD-DFT calculation of excitation energies further emphasizes that the bond length alternation is a factor which must be taken into con-

Excitation energies from experiment and theory

A more detailed look at the comparison between the theoretical and experimental determination of excitation energies is necessary to understand why M06HF manages to replicate the experimental reference rather well, even though it fails to obtain a correct ground state


5


brationally resolved absorption spectra of the respective oligomer in solution are taken using different solvents. The maxima of the first peaks are then fitted to an appropriate equation, such as the Onsager model, 49 to extrapolate the adiabatic excitation energy to vacuum. 2,50

Adiabatic excitation energy (eV)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 19

Figure 6: Schematic representation of ground state and optically active excited state potential energy modes, highlighting the difference between adiabatic excitation energy (E00 ) and vertical excitation energy (Evert ). geometry. First, from now on we will look at the adiabatic rather than vertical excitation energy. There are a number of reasons this makes sense, including the fact that in this case the experimental data set chosen for comparison merely gives a lower bound for the vertical excitation energy, equal to the adiabatic excitation energy. 12 A detailed description of how the experimental data set was assembled is provided in the Supporting Information of Reference 12. For the adiabatic excitation energy to be a good estimate of the vertical excitation energy, however, the relaxation energy ΔEeq caused by geometry changes between the ground and excited state (see Figure 6) should be negligible. As will be demonstrated below, this is not true for the case of t-PA. The error caused by this assumption can be avoided by using the adiabatic excitation energy instead of the vertical one. Also, it has been argued that the molecular equivalent of the optical band gap is more accurately described by the adiabatic excitation energy than by the vertical excitation energy. 2 This is due to the adiabatic excitation energy including geometry changes between the ground and the excited state, see Figure 6. Computationally, the adiabatic excitation energy is determined by the energy difference between the optimized geometries of the ground and excited state. Alternatively, the adiabatic excitation energy can be estimated from experiment using a two-step procedure. First, vi-

5 4.5 4 3.5 Experimental 1 Experimental 2 Experimental 3

3 2.5 0

0.1

0.2

0.3

1/N

Figure 7: Comparison of three experimental sets, experimental 1 12 (used throughout this paper), experimental 2 51 and experimental 3, 50 of adiabatic excitation energies for polyacetylene oligomers, N= number of double bonds. All three include an extrapolation to the polymer limit using a Kuhn Fit. This procedure is well established 2,50,52,53 and, as we will demonstrate below, the uncertainties are small, leading to minor differences in the adiabatic excitation energies in vacuum found for different measurements. The excitation energies can be split into three regions, dependent on oligomer size. For very short chain lengths, the chain ends dominate, leading to discrepancies in the geometry and excitation energies compared to longer chains. To avoid including this uncharacteristic behavior, we only consider oligomers with N≥ 6 in our work. 18 The second region is the linear region, where excitation energies linearly decrease with respect to 1/N. It has been shown, however, that a linear extrapolation to the polymer limit produces large errors. 6,7,54 This is due to the fact that, once the effective conjugation length is reached, the excitation energy no longer decreases linearly with oligomer size, but slopes off and approaches a limit. 7 A number of fitting procedures have


6

Page 7 of 19

been developed to accurately describe this behavior, thus, allowing to extrapolate excitation energies of oligomers reliably to the polymer limit. 2–7 Compared to other models, the Kuhn Fit introduced above has been shown to provide superior results when used for the extrapolation of the optical band gap, even when large oligomers are not included in the fit. 6 Since t-PA has been studied extensively in the literature and since, as discussed above, the estimation of adiabatic excitation energies from experiment is subject to a number of approximations and extrapolations, several different experimental numbers for the adiabatic excitation energies of t-PA can be found in the literature. To provide an estimate for the error bar on the these numbers, we compared experimental data from three different sources and applied the Kuhn-Fit to all three data sets independently (see Figure 7). Although the individual test sets vary in the number of data points and absolute values found for the individual adiabatic excitation energies, the values obtained for the optical gap in the polymer limit agree within an energy window of only 0.1 eV. This sets the error bar on the experimental numbers of the adiabatic excitation energy to which we compare our calculations.

only provides the adiabatic excitation energy, but additionally allows for visualization of the extent of localization for the optical excitation in the t-PA oligomers. 55 This is demonstrated in Figure 8, which shows the BLA in the ground state as well as in the first singlet excited state in a polyene with N=50 as obtained from B3LYP, CAM-B3LYP, and M06HF. The localization of the excitation is strongly dependent on both the total amount of exact exchange, as well as the long range behaviour of the functional. 55,56 While B3LYP tends to delocalize the excitation over the whole length of the chain, the HF-rich M06HF strongly localizes the excitation. The range-separated CAMB3LYP yields an excition delocalization that lies somewhere between the two extremes.

3.4

Figure 9: Change of energy for fractional electron numbers. DFT functionals are below the desired straight line behavior, causing a favorable energy for fractional electrons and hence delocalizing charges in a system. Exact exchange overestimates the energy of fractional charges, artificially stabilizing whole electron numbers, which leads to very localized charges.

0.4

0.1

M06HF

CAM-B3LYP

-0.1 100

50

100

50

0.2 0 -0.2 113

2

113.5

114

δ

114.5

115

M06HF LRC-ωPBE* PBEh* B3LYP

113

Ground state st 1 excited state

50

4

0

0

0

ΔE (eV)

EN+δ - EN (eV)

6

Localization of adiabatic excitation B3LYP

BLA (Å)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


113.5

114

δ

114.5

115

100

Carbon number Carbon number Carbon number

The fact that semilocal or standard global hybrid functionals such as B3LYP tend to overestimate electron delocalizion whereas HF-rich functionals such as M06HF tend to underestimate it is well known in the literature. 23,55–60 This error of approximate exchange-correlation functionals is often referred to as the localization/delocalization error or the many-electron self-interaction error. 59,61 Frequently, it is discussed in the context of fractional particle numbers, which can be introduced into the DFT

Figure 8: Visualization of the differences in localization of the first optically excited state of polyacetylene resulting from the use of different functionals. The localization is shown through the differences in bond length alternation (BLA) of the ground and first optically excited state for an oligomer with 50 double bonds. An excited state geometry optimization not


7


demonstrated by the BLA-optimized functionals PBEh* and LC-ωPBE* in Figure 9.

formalism by coupling the many-electron system to a bath of particles. 62 The implications of the localization/delocalization error have been discussed in a number of recent review articles. 57,63–65 As demonstrated by Perdew and co-workers, the exact total energy as a function of fractional electron number is a series of straight lines with kinks at integer occupations. 62 This, however, is not the case for the commonly used exchange-correlation functionals. Most DFT functionals show a concave behavior, as demonstrated by the B3LYP curve shown in Figure 9. This error is referred to as the delocalization error, as it leads to a spurious stabilization of fractional particle numbers, thus, encouraging charges (and excitations) to delocalize across the entire chain length. In contrast, large amounts of HF exchange, as present in M06HF, destabilize partial charges compared to the desired straight line behavior, leading to very localized charges becoming more favorable. The area between the curves and the ideal straight line can be used as a measure for the many electron self-interaction error inherent to the respective functional. 59,61 It has been demonstrated that the localization/delocalization error can be minimized by using long-range corrected hybrid functionals in which the range-separation parameter is obtained separately for each molecule using a nonempirical tuning procedure. 57,65–68 However, as the range-separation parameter depends on the size and degree of conjugation of the molecule, the tuning procedure to find the system specific range separation parameter for π-conjugated polymers is problematic. 69,70 For some combinations of functional and polymer a limit is found for oligomers of increasing length, resulting in a range separation parameter which adequately describes larger oligomers and the infinite system. 71 For t-PA and LRC-ωPBE, however, this is not the case, 23,57,70 at least not for the chain lengths considered here. The range separation parameter steadily decreases with increasing chain length, and no single range separation parameter for the system can be found using standard tuning methods. 23,57,70 Even so, including a certain amount of HF exchange in a functional reduces the delocalization error, as is

-1

-1

ω=0.1(a0 )

ω=0.2(a0 )

BLA (Å)

st

1 excited state Ground state

0.1 0 -0.1

-1

-1

ω=0.3(a0 ) BLA (Å)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 19

ω=0.4(a0 )

0.1 0 -0.1 0

25

Carbon number

50

25

50

Carbon number

Figure 10: Change of localization of the first optically active excited state of a polyacetylene oligomer with 26 double bonds for varying amounts of exact exchange using the LRCωPBE functional. The range separation parameter ω is responsible for the amount of exact exchange, with a larger value being equal to a larger amount of exact exchange. Figure 10 shows just how much varying the range separation parameter of a functional such as LRC-ωPBE, and hence the total amount of exact exchange present in the functional, influences the localization of the adiabatic excitation in a polyene with 26 double bonds (N=26). For ω=0.1 a−1 0 the excitation is completely delocalized over the entire chain length. As the range seperation parameter, and hence the total amount of exact exchange, is increased to ω=0.4 a−1 the excitation becomes more and 0 more localized. The range separation parameter of ω=0.188 a−1 0 used in LRC-ωPBE* therefore leads to a semilocalized description of the excitation, spanning about 20 double bonds, while describing the ground state geometry well. In contrast, M06HF leads to a more localized excitation over around 12 double bonds, and greatly overestimates the degree of BLA in the ground state. Comparing the adiabatic excitation energies of M06HF and LRC-ωPBE* in Figure 11, it becomes apparent that a more localized descrip-


8

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Page 9 of 19


obtained for oligomers in solution at room temperature. Although special care is taken to eliminate the solvent effect and extrapolate adiabatic excitation energies in vacuum, these values can not be directly compared to our calculations. Thermal energy in the system allows vibrational and rotational motion of the oligomers, leading to changes away from the minimum energy geometry. Kinks and twists of a certain magnitude in the backbone of the chain disrupt conjugation, which is why it has been suggested that long πconjugated chains behave like a superposition of shorter chains, 72,73 rather than the planar straight chains of C2h symmetry used for our calculations. To capture these different configurations of the system caused by thermal disorder, which are inevitably included in the experimental data, MD calculations are performed on single oligomer chains in solution.

4

3

Experimental M06HF LRC-ωPBE* PBEh* B3LYP

2

1 0

0.05

0.1

0.15

1/N

Figure 11: Adiabatic excitation energies of polyacetylene oligomers (N=number of double bonds) for functionals able to replicate the correct bond length alternation (ωPBE*, PBEh*) compared to commonly used functionals and experimental data. Kuhn fits are used for extrapolation to the polymer limit. tion of the excitation seems to better replicate experimental data. While LRC-ωPBE* and PBEh* show clear improvement over B3LYP, both for the excitation energy values and the slope of the curve, their adiabatic excitation energies still fall off too far for long oligomer chain lengths. This suggests an overestimation of the conjugation length and a too delocalized description of the excitation. An adjustment of the amount of exact exchange in PBEh* or the range separation parameter in LRC-ωPBE* to match the experimental data exactly is not possible, but also not desirable, as it would lead to a wrong description of the BLA in the ground state.

4.1

4

Figure 12: Definition of atom types and bonds used for the OPLSAA force field parameterization.

Parameter fitting

To the best of our knowledge there are no standard force field parameters available that accurately describe the BLA of t-PA. Given the strong influence of the geometry on the excitation energies discussed above, it is vitally important that the correct BLA is incorporated into our force field. To achieve this, a model describing t-PA using a total of 8 atom types for the carbon backbone is parameterized for the OPLS-AA force field.

Results MD

The fact that a functional with full HF exchange manages to replicate the slope of the experimental data, while LRC-ωPBE* and PBEh* fall off too steeply, suggests that the increase of conjugation length with increased chain length is overestimated in our calculations. Considering all calculations are performed on perfectly symetrical, straight chains, this is not surprising. Experimental data is

Four atom types are necessary to define alternating single and double bonds, additional four atom types are used to describe the slightly different BLA found at the start and end of each oligomer chain. The average bond lengths and average angles for all atoms are fitted to PBEh* DFT calculations for an oligomer with 22 double bonds. Mulliken charges from the same cal-


9


culations are used as an estimate for partial charges on each atom. The force constants for previously defined bonding atom interactions 39 adequately describe t-PA oligomers. A correct description of the cis/trans barrier height for rotation around the single bonds is likewise important, as the relative energies of the two conformations determine the abundance at which they can be found in a system of many oligomer chains. The barrier height of the rotation determines the transition probability between the two and therefore the timescales involved. An ab initio reference for the change of energy with rotation around a dihedral angle is calculated by DFT geometry optimizations using PBEh* where the dihedral angle around the central single bond is frozen. PBEh* is a suitable functional for the reference, as it yields excellent BLAs. In addition, it has been demonstrated that, due to a close link between the many electron self interaction error and the error found for torsion energies, functionals with a medium amount of exact exchange and a small localization/delocalization error also describe the torsion potential of t-PA accurately. 40 Relative energy (kJ/mol)

for rotation around the single bond and manually adjust the parameters in the force field until they fit the DFT curve. The resulting fit with the final parameters is provided in Figure 13. Further information regarding the force field parameterization, as well as a comparison between the DFT and MD minimum energy geometries and vibrational frequencies, is provided in the supporting information.

4.2 Percentage of snapshots

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

N=8 N=10 N=20 N=34 N=50

60

40

20

0

1

2

3

4

Number of cis bonds in chain

Figure 14: Percentage of snapshots with cis bonds for md trajectories of several polyacetylene oligomers of increasing length. Except for zero, all chains with at least the desired number of cis bonds were counted.

25 20 15 10

0 0

Thermal results

80

0

30

5

Page 10 of 19

Even though the barrier height for cis/trans rotation is very high (around 30 kJ/mol), MD trajectories of the oligomers show a non negligible percentage of twists about one or more single bonds (see Figure 14). Due to the high barrier, the cis/trans rotations do not occur frequently for the considered timescale of a few nanoseconds, and it is important to note that the MD trajectories only manage to capture a representative fraction of the configurations that would be present in an experimental setup, where the timescales are much larger. For the configurations captured by the MD simulations of very short oligomers the all-trans geometry is still dominant (see Figure 14), but as the oligomer chain length increases, chains with at least one cis bond become more prominent. For chains with 50 double bonds, about 50% of

PBEh* Fitted parameters

50

100

150

200

Dihedral angle

Figure 13: Change in energy with change in the dihedral angle around a single bond of polyacetylene for ab initio calculations using PBEh* and the corresponding fit obtained by parameterization of the OPLSAA force field. A direct fit of the equation used to describe dihedral angles for the OPLS-AA force field 39 to the DFT energy curve is not possible, as changing the dihedral angle also effects other factors contributing to the total energy. It is therefore necessary to optimize the oligomer geometry via an MD step, then calculate single point energies


10

Page 11 of 19

using Gaussian functions with a FWHM of 0.05 eV. This width was chosen empirically to give a smooth final spectrum without artificial fine structure or broadening. This is illustrated in Figure 15, which shows the Gaussians for each excitation energy as well as the final spectrum created through a normalized sum of all overlayed Gaussians. The maximum of this curve gives the vertical excitation energy for the respective oligomer length. 75 Finding an estimate for the relaxation energy is not trivial, as it is sensitive to the oligomer length considered and the choice of functional used for both the geometry optimization and TD-DFT steps. Which data is included in an averaged estimate of the relaxation energy therefore has to be carefully considered. As can be seen from Figure 6 the slope of the excited state potential surface needs to be described correctly to obtain a reliable relaxation energy. Hartree-Fock methods largely overestimate the slope of the excited state potential energy surface and hence the vertical excitation and relaxation energy. 2,76 Data from previous M06HF calculations is hence excluded from the estimate. B3LYP shows very strong dependence of the relaxation energy on the oligomer length, even for very long oligomer chains, making it unsuitable for an averaged correction. The adapted functional PBEh*, however, includes only medium amounts of exact exchange, while reaching a limit for the relaxation energy after a few t-PA repeat units. An averaged relaxation energy of 0.2 eV is obtained using this functional (the relaxation energies for the individual chains are provided in the Supporting Information). While this inherently leads to a small error for the adiabatic excitation energy of short t-PA chains, it gives a well-grounded estimate of the adiabatic excitation energy for medium and long oligomers. The effect of thermal disorder on the adiabatic excitation energies of t-PA oligomers varies greatly with respect to the choice of functional used, see Figure 16. This is a combination of geometry and TD-DFT effects. The oligomers obtained from the MD snapshots have an average BLA and geometry equal to that found using PBEh*, which varies to dif-

the snapshots include at least two cis bonds. This is important, as cis bonds break the πconjugation. In addition to trans to cis rotations, curves in the backbone along multiple carbons or twists smaller than 100◦ in the dihedral angles can also lead to breaks in conjugation. These different types of geometry fluctuations can all be observed during our MD simulation with the adapted OPLS-AA force field. 10

Oscillator strength

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


8 6 4 2 0

2

2.2

2.4

2.6

2.8

Energy (eV)

Figure 15: TD-DFT excitation energies for the first optically active excited state of a polyacetylene oligomer with 22 double bonds for 100 snapshots taken from an md trajectory. Stick spectra broadened with Gaussians, which were summed and normalized to give a final spectrum. While recent studies have developed methods with which it is possible to obtain adiabatic excitation energies directly from MD simulations, 74 this would require an enormous increase in computational cost and is outside of the scope of this work. Instead, vertical excitation energies are calculated using TDDFT on top of MD geometries, after which an approximation for the relaxation energy from the vertical to the adiabatic excitation is made and subtracted for each oligomer. To obtain an average vertical excitation energy for each t-PA oligomer length, several MD snapshots are taken. The geometry of the oligomer is extracted from each one and used as input for a TD-DFT calculation of the vertical excitation energy. The relevant excitation for each snapshot for an oligomer with 22 double bonds is shown by the stick spectrum in Figure 15. To obtain a complete absorption spectrum from these peaks, the vertical excitation energies for each snapshot are broadened


11

4

3

2

x x

1 0

0.05

0.1

Page 12 of 19

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60



Experimental M06HF on MD LRC-ωPBE* on MD LRC-ωPBE* B3LYP on MD Polymer experimental

0.15

1/N

Figure 17: Natural transition orbitals (NTOs) for the optical transition of a polyacetylene oligomer with N=50 using LRC-ωPBE*. On the top the NTOs were calculated using an MD snapshot with multiple twists, on the bottom the NTOs were computed using a geometry optimized by DFT.

Figure 16: Influence of thermal disorder on the adiabatic excitation energies of polyacetylene oligomers (N=number of double bonds) for B3LYP, M06HF and ωPBE*. Excitation energies for ωPBE* without disorder are shown for comparison. Kuhn fits are used for extrapolation to the polymer limit. Additionally, the experimental value for the absorption maximum of the solid state polymer bandgap is provided. 77

very similar, so the average BLA seen in the MD geometries is very close to the BLA found for LRC-ωPBE*. As discussed above, our LRCωPBE* has a maximum effective conjugation length of roughly 20 double bonds. This does not change when adding thermal disorder. For short chains, the percentage of cis-twists is too low to cause a localizing effect. As the oligomer length increases, the number of twists become more significant, however the chains are long enough for most configurations to allow the excitation to fit between twists. Figure 17 shows there is little difference between the natural transition orbitals of a perfectly linear chain with N=50 and those computed for an MD snapshot with several cis bonds for functionals with this degree of localization. Inclusion of thermal disorder hence does not alter the maximum conjugation length, even though it leads to a range of configurations with visually very different geometries. Therefore, no large shift of the band gap estimate for LRC-ωPBE* is found. Figure 16 illustrates the very good agreement of excitation energies for t-PA oligomers calculated using M06HF on top of MD geometries with experimental data. M06HF strongly localizes excitations, as shown previously in Figure

ferent extents from that found by other functionals. Compared to the B3LYP geometry, the MD geometry has greatly increased BLA. As evident from Figure 5, however, the vertical excitation energies calculated using B3LYP are not very dependent on the geometry used. Rather the main change in the excitation energy observed is due to the fact that B3LYP does not reach its maximum conjugation length for the oligomer sizes considered, as evident in Figure 8. Twists and kinks in the chain, caused by thermal motion, lead to a shorter conjugation length, forcing a localization of the excitation. Figure 16 shows that B3LYP demonstrates a sloping behavior, previously not observed, when thermal disorder is considered. This leads to higher excitation energies for long oligomer chains, resulting in a difference of about 0.2 eV for the band gap estimate at the polymer limit. In contrast to this, it can be seen in Figure 16 that almost no change in slope or absolute values of the excitation energies is observed when considering thermal motion for LRC-ωPBE*. The geometries of PBEh* and LRC-ωPBE* are


12

Page 13 of 19

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


ure 16 at 1.9 eV. 77 The onset of absorption for the thin film polymer is found at about 1.5eV, and so the t-PA band gap is quoted in literature as being between 1.5 and 1.9eV. 80,81 These values are in excellent agreement with the extrapolated polymer limit found using LRC-ωPBE*. There are, however, discrepancies between our calculations and this data which make a direct comparison difficult. In the solid state, intermolecular interactions strongly affect the optical properties of a system. 2 Our calculations, however, do not capture intermolecular interactions, as only a single oligomer chain is considered at a time. We hence cannot capture the significant bathochromic shift of the excitation energy in the solid state. 2 In fact, the comparison to the extrapolated solution data suggests that the solid state shift is on the order of 1 eV. Hence, we conclude that the nice agreement of our extrapolated optical band-gap in the polymer limit with the solid-state data is a coincidence, based on a cancelation of errors between the limitation of our TD-DFT calculations and the neglect of the solid-state shift.

8. Consequently, just as with LRC-ωPBE*, geometry fluctuations caused by thermal motion have no significant effect on the effective conjugation length. The difference in adiabatic excitation energies between a DFT/TD-DFT calculation using M06HF, plotted in Figure 11, and the MD/TD-DFT calculation, plotted in Figure 16, therefore originates from the difference in the ground state geometries. The MD simulation enforces a correct description of BLA through its parameterization. This is further supported by the finding that the results shown in Figure 11 for M06HF on MD geometries are very similar to the results of the M06HF TDDFT calculations on PBEh* geometries shown in Figure 5. Intererestinly, the combination of a correct description of BLA in the ground-state geometry (as provided by PBEh* as well as by our MD force field) with a strongly localized description of the excitation through M06HF leads to adiabatic excitation energies very close to experimental data, a Kuhn-Fit with a realistic slope and, consequently, a band gap estimate from the extrapolation to the polymer limit in excellent agreement with experiment. In other words, in agreement with previous findings, 14 the amount of HF exchange needed to describe the delocalization of an excitation in a π-conjugated chain is much larger than what is needed for a realistic description of the ground-state geometry. This also means that the methods available to quantify the localization/delocalization error in a systems ground-state might not be directly applicable to its excited state as well. Given that the mixing of a fraction of HF-exchange with semilocal exchange functionals has been interpreted as a way to mimic non-local correlation effects, 63,78,79 a possible interpretation of this finding is that a correct description of correlation in the ground state does not imply a correct description of the correlation in the excitedstate and vice versa. For electronic applications, the property we are ultimately interested in is the optical bandgap of the polymer in the solid state. An experimental value for the absorption maximum of the bandgap for solid state t-PA is shown in Fig-

5

Conclusions

Our aim was to break down problems faced when calculating the optical band gap of πconjugated polymers using the model system t-PA to allow a stepwise minimization of errors. By comparing the effects of using different functionals for the geometry optimization and TDDFT steps we found that a correct description of the ground-state BLA, as well as the correct choice of functional for the calculation of excitation energies is vital. The previously adapted LRC-ωPBE* and PBEh* functionals, as well as XYGJ-OS are shown to all adequately describe the ground-state BLA in polyenes. A sensitive dependence of the localization of the excitation and hence the excitation energies on the choice of functional, more specifically on the amount of exact exchange present in each functional, is shown. Both the difficulties in describing BLA and the extent of localization of the excitation are linked to the localization/delocalization error found in HF/DFT respectively. While


13


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 19

als or hybrid functionals with an intermediate amount of HF exchange, the latter typically requires much larger fractions of HF exchange.

the functionals adapted to reproduce the BLA LRC-ωPBE* and PBEh* give better adiabatic excitation energies than most standard functionals, a surprisingly high amount of HF exchange as found, e.g., in M06HF, is necessary to replicate the slope of the experimental data. Including thermally accessible configurations with less symmetry and conjugation through MD was shown to have very little effect on the excitation energies for LRC-ωPBE*.

Supporting Information Available: The supporting information provides the details of the force-field fitting procedure, a comparison of the force-field results for the geometry and vibrational frequencies with the DFT reference, as well as a detailed description of the averaging procedure used to determine the optical excitation energy from the MD snapshots. This material is available free of charge via the Internet at http://pubs.acs.org/.

In summary, we find that the prediction of the optical band gap of polymers using the oligomer approach still poses a major challenge to (TD-)DFT approaches. While we can verify the accuracy of functionals for the description of the ground-state geometries by comparison to higher level wavefunction calculations such as CCSD(T) as well as by other means such as analyzing the localization/delocalization error, the validity of a certain functional for the description of the (de-)localization of the excited state in a π-conjugated chain cannot be verified as easily. By using MD simulations, we could show that the discrepency between the TD-DFT calculations of the excitation energy and the experimental reference data is very likely not caused by thermal fluctuations. Instead, we conclude that the poor agreement with experimental data is caused by the failure of TD-DFT to describe the degree of exciton delocalization correctly. In particular, the asumption that a functional that yields a correct description of the ground-state BLA and a small many-electron self-interaction error would also predict a correct delocalization of an excitation could not be verified. Hence, we conclude that our current understanding of the localization/delocalization error inherent to certain approximate exchange-correlation functionals in DFT cannot be straight-forwardly applied to excitations calculated within TD-DFT. A practical workaround for the prediction of optical gaps in the polymer limit using the KuhnFit approach is to employ different exchangecorrelation functionals for the ground-state geometry relaxation and the TD-DFT calculation of excitation energies. While the former can be described well by double hybrid function-

Acknowledgement The authors thank S. Kümmel, R. Baer, J.-L. Brédas and V. Vlcek for helpful discussions. Computer time was provided by the National Energy Research Scientific Computing Center (NERSC), which is supported by the Office of Science of the U.S. Department of Energy under contract DE-AC0205CH11231.

References (1) Risko, C.; McGehee, M. D.; Brédas, J.-L. Chem. Sci. 2011, 2, 1200–1218. (2) Gierschner, J.; Cornil, J.; Egelhaaf, H.-J. Adv. Mater. 2007, 19, 173–191. (3) Kuhn, H. J. Chem. Phys. 1949, 17, 1198– 1212. (4) Hückel, E. Z. Phys. 1931, 70, 204–286. (5) Lewis, G. N.; Calvin, M. Chem. Rev. 1939, 25, 273–328. (6) Torras, J.; Casanovas, J.; Alemán, C. J. Phys. Chem. A 2012, 116, 7571–7583. (7) Meier, H.; Stalmach, U.; Kolshorn, H. Acta Polymer. 1997, 48, 379–384. (8) Kuhn, W. Helv. Chim. Acta 1948, 31, 1780–1799. (9) Su, W. P.; Schrieffer, J. R.; Heeger, A. J. Phys. Rev. Lett. 1979, 42, 1698–1701.


14

Page 15 of 19

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(10) Brédas, J. L.; Dory, M.; André, J. M. J. Chem. Phys. 1985, 83, 5242–5249.

(25) Zhao, Y.; Truhlar, D. G. J. Phys. Chem. A 2006, 110, 13126–13130.

(11) Hsu, C.-P.; Hirata, S.; Head-Gordon, M. J. Phys. Chem. A 2001, 105, 451–458.

(26) Adamo, C.; Barone, V. J. Chem. Phys. 1999, 110, 6158–6170.

(12) Schmidt, M.; Tavan, P. J. Chem. Phys. 2012, 136, 124309.

(27) Becke, A. D. J. Chem. Phys. 1993, 98, 5648–5652.

(13) Milián-Medina, B. n.; Gierschner, J. WIREs Comput. Mol. Sci. 2012, 2, 513– 524.

(28) Yanai, T.; Tew, D. P.; Handy, N. C. Chem. Phys. Lett. 2004, 393, 51 – 57. (29) Grimme, S. J. Chem. Phys. 2006, 124, 034108.

(14) Wykes, M.; Milián-Medina, B. n.; Gierschner, J. Front. Chem. 2013, 1, 35.

(30) 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.; Montgomerywil, J. A.; 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, Ö.; Foresman, J. B.; Ortiz, J. V.; Cioslowski, J.; Fox, D. J. Gaussian09 Revision D.01. Gaussian Inc. Wallingford CT 2009.

(15) Yannoni, C. S.; Clarke, T. C. Phys. Rev. Lett. 1983, 51, 1191–1193. (16) Marder, S. R.; Gorman, C. B.; Meyers, F.; Perry, J. W.; Bourhill, G.; Brédas, J.-L.; Pierce, B. M. Science 1994, 265, 632–635. (17) Baeriswyl, D.; Campbell, D. K.; Mazumdar, S. In Conjugated Conducting Polymers; Kiess, H. G., Ed.; Springer Ser. Solid-State Sci.; Springer Berlin Heidelberg, 1992; Vol. 102; pp 7–133. (18) Tretiak, S.; Chernyak, V.; Mukamel, S. Phys. Rev. Lett. 1996, 77, 4656–4659. (19) Gallagher, F. B.; Mazumdar, S. Phys. Rev. B 1997, 56, 15025–15039. (20) Jacquemin, D.; Femenias, A.; Chermette, H.; Ciofini, I.; Adamo, C.; André, J.-M.; Perpète, E. A. J. Phys. Chem. A 2006, 110, 5952–5959. (21) Sancho-García, J. C.; Pérez-Jiménez, A. J. Phys. Chem. Chem. Phys. 2007, 9, 5874– 5879. (22) Jacquemin, D.; Adamo, C. J. Chem. Theory Comput. 2011, 7, 369–376. (23) Körzdörfer, T.; Parrish, R. M.; Sears, J. S.; Sherrill, C. D.; Brédas, J.-L. J. Chem. Phys. 2012, 137, 124305.

(31) Zhang, I. Y.; Xu, X.; Jung, Y.; Goddard, W. A. Proc. Natl. Acad. Sci. U.S.A. 2011, 108, 19896–19900.

(24) Wykes, M.; Su, N. Q.; Xu, X.; Adamo, C.; Sancho-García, J.-C. J. Chem. Theory Comput. 2015, 11, 832–838.


15


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 19

rill, C. D.; Truhlar, D. G.; Warshel, A.; Xu, X.; Aspuru-Guzik, A.; Baer, R.; Bell, A. T.; Besley, N. A.; Chai, J.D.; Dreuw, A.; Dunietz, B. D.; Furlani, T. R.; Gwaltney, S. R.; Hsu, C.-P.; Jung, Y.; Kong, J.; Lambrecht, D. S.; Liang, W.; Ochsenfeld, C.; Rassolov, V. A.; Slipchenko, L. V.; Subotnik, J. E.; Voorhis, T. V.; Herbert, J. M.; Krylov, A. I.; Gill, P. M.; Head-Gordon, M. Mol. Phys. 2015, 113, 184–215.

(32) Shao, Y.; Gan, Z.; Epifanovsky, E.; Gilbert, A. T.; Wormit, M.; Kussmann, J.; Lange, A. W.; Behn, A.; Deng, J.; Feng, X.; Ghosh, D.; Goldey, M.; Horn, P. R.; Jacobson, L. D.; Kaliman, I.; Khaliullin, R. Z.; Kuś, T.; Landau, A.; Liu, J.; Proynov, E. I.; Rhee, Y. M.; Richard, R. M.; Rohrdanz, M. A.; Steele, R. P.; Sundstrom, E. J.; Woodcock, H. L.; Zimmerman, P. M.; Zuev, D.; Albrecht, B.; Alguire, E.; Austin, B.; Beran, G. J. O.; Bernard, Y. A.; Berquist, E.; Brandhorst, K.; Bravaya, K. B.; Brown, S. T.; Casanova, D.; Chang, C.M.; Chen, Y.; Chien, S. H.; Closser, K. D.; Crittenden, D. L.; Diedenhofen, M.; DiStasio, R. A.; Do, H.; Dutoi, A. D.; Edgar, R. G.; Fatehi, S.; Fusti-Molnar, L.; Ghysels, A.; Golubeva-Zadorozhnaya, A.; Gomes, J.; Hanson-Heine, M. W.; Harbach, P. H.; Hauser, A. W.; Hohenstein, E. G.; Holden, Z. C.; Jagau, T.-C.; Ji, H.; Kaduk, B.; Khistyaev, K.; Kim, J.; Kim, J.; King, R. A.; Klunzinger, P.; Kosenkov, D.; Kowalczyk, T.; Krauter, C. M.; Lao, K. U.; Laurent, A. D.; Lawler, K. V.; Levchenko, S. V.; Lin, C. Y.; Liu, F.; Livshits, E.; Lochan, R. C.; Luenser, A.; Manohar, P.; Manzer, S. F.; Mao, S.P.; Mardirossian, N.; Marenich, A. V.; Maurer, S. A.; Mayhall, N. J.; Neuscamman, E.; Oana, C. M.; OlivaresAmaya, R.; O’Neill, D. P.; Parkhill, J. A.; Perrine, T. M.; Peverati, R.; Prociuk, A.; Rehn, D. R.; Rosta, E.; Russ, N. J.; Sharada, S. M.; Sharma, S.; Small, D. W.; Sodt, A.; Stein, T.; Stück, D.; Su, Y.-C.; Thom, A. J.; Tsuchimochi, T.; Vanovschi, V.; Vogt, L.; Vydrov, O.; Wang, T.; Watson, M. A.; Wenzel, J.; White, A.; Williams, C. F.; Yang, J.; Yeganeh, S.; Yost, S. R.; You, Z.-Q.; Zhang, I. Y.; Zhang, X.; Zhao, Y.; Brooks, B. R.; Chan, G. K.; Chipman, D. M.; Cramer, C. J.; Goddard, W. A.; Gordon, M. S.; Hehre, W. J.; Klamt, A.; Schaefer, H. F.; Schmidt, M. W.; Sher-

(33) Wang, Y.-L.; Wu, G.-S. Int. J. Quant. Chem. 2008, 108, 430–439. (34) Sears, J. S.; Körzdörfer, T.; Zhang, C.-R.; Brédas, J.-L. J. Chem. Phys. 2011, 135, 151103. (35) Peach, M. J. G.; Williamson, M. J.; Tozer, D. J. J. Chem. Theory Comput. 2011, 7, 3578–3585. (36) Turney, J. M.; Simmonett, A. C.; Parrish, R. M.; Hohenstein, E. G.; Evangelista, F. A.; Fermann, J. T.; Mintz, B. J.; Burns, L. A.; Wilke, J. J.; Abrams, M. L.; Russ, N. J.; Leininger, M. L.; Janssen, C. L.; Seidl, E. T.; Allen, W. D.; Schaefer, H. F.; King, R. A.; Valeev, E. F.; Sherrill, C. D.; Crawford, T. D. WIREs Comput. Mol. Sci. 2012, 2, 556. (37) Papajak, E.; Zheng, J.; Xu, X.; Leverentz, H. R.; Truhlar, D. G. J. Chem. Theory Comput. 2011, 7, 3027–3034. (38) Jorgensen, W. L.; Maxwell, D. S.; TiradoRives, J. J. Am. Chem. Soc. 1996, 118, 11225–11236. (39) Siu, S. W. I.; Pluhackova, K.; Böckmann, R. A. J. Chem. Theory Comput. 2012, 8, 1459–1470. (40) Sutton, C.; Körzdörfer, T.; Gray, M. T.; Brunsfeld, M.; Parrish, R. M.; Sherrill, C. D.; Sears, J. S.; Brédas, J.-L. J. Chem. Phys. 2014, 140, 054310.


16

Page 17 of 19

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(54) Egelhaaf, H.-J.; Oelkrug, D.; Gebauer, W.; Sokolowski, M.; Umbach, E.; Fischer, T.; Bäuerle, P. Opt. Mater. 1998, 9, 59 – 64.

(41) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. J. Chem. Phys. 1983, 79, 926–935. (42) Páll, S.; Abraham, M.; Kutzner, C.; Hess, B.; Lindahl, E. In Solving Software Challenges for Exascale; Markidis, S., Laure, E., Eds.; Lect. Notes Comput. Sc.; Springer International Publishing, 2015; Vol. 8759; pp 3–27.

(55) Nayyar, I. H.; Batista, E. R.; Tretiak, S.; Saxena, A.; Smith, D. L.; Martin, R. L. J. Phys. Chem. Lett. 2011, 2, 566–571. (56) Igumenshchev, K. I.; Tretiak, S.; Chernyak, V. Y. J. Chem. Phys. 2007, 127, 114902.

(43) Ho Choi, C.; Kertesz, M.; Karpfen, A. J. Chem. Phys. 1997, 107, 6712–6721.

(57) Körzdörfer, T.; Brédas, J.-L. Acc. Chem. Res. 2014, 47, 3284–3291.

(44) Zhao, Y.; Truhlar, D. G. J. Phys. Chem. A 2006, 110, 10478–10486.

(58) Mori-Sánchez, P.; Cohen, A. J.; Yang, W. Phys. Rev. Lett. 2008, 100, 146401.

(45) Jacquemin, D.; Perpète, E. A.; Scalmani, G.; Frisch, M. J.; Kobayashi, R.; Adamo, C. J. Chem. Phys. 2007, 126, 144105.

(59) Mori-Sánchez, P.; Cohen, A. J.; Yang, W. J. Chem. Phys. 2006, 125, 201102.

(46) Chabbal, S.; Jacquemin, D.; Adamo, C.; Stoll, H.; Leininger, T. J. Chem. Phys. 2010, 133, 151104.

(60) Haunschild, R.; Henderson, T. M.; Jiménez-Hoyos, C. A.; Scuseria, G. E. J. Chem. Phys. 2010, 133, 134116.

(47) Vydrov, O. A.; Scuseria, G. E. J. Chem. Phys. 2006, 125, 234109.

(61) Ruzsinszky, A.; Perdew, J. P.; Csonka, G. I.; Vydrov, O. A.; Scuseria, G. E. J. Chem. Phys. 2006, 125, 194112.

(48) Goerigk, L.; Grimme, S. J. Chem. Theory Comput. 2011, 7, 3272–3277.

(62) Perdew, J. P.; Parr, R. G.; Levy, M.; Balduz, J. L. Phys. Rev. Lett. 1982, 49, 1691– 1694.

(49) Onsager, L. J. Am. Chem. Soc. 1936, 58, 1486–1493. (50) D’Amico, K. L.; Manos, C.; Christensen, R. L. J. Am. Chem. Soc. 1980, 102, 1777–1782.

(63) Kümmel, S.; Kronik, L. Rev. Mod. Phys. 2008, 80, 3–60. (64) Cohen, A. J.; Mori-Sánchez, P.; Yang, W. Chem. Rev. 2012, 112, 289–320.

(51) Brédas, J.-L.; Silbey, R.; Boudreaux, D. S.; Chance, R. R. J. Am. Chem. Soc. 1983, 105, 6555–6559.

(65) Baer, R.; Livshits, E.; Salzner, U. Annu. Rev. Phys. Chem. 2010, 61, 85–109.

(52) Christensen, R. L.; Galinato, M. G. I.; Chu, E. F.; Howard, J. N.; Broene, R. D.; Frank, H. A. J. Phys. Chem. A 2008, 112, 12629–12636.

(66) Stein, T.; Kronik, L.; Baer, R. J. Chem. Phys. 2009, 131, 244119. (67) Stein, T.; Kronik, L.; Baer, R. J. Am. Chem. Soc. 2009, 131, 2818–2820.

(53) Gierschner, J.; Mack, H.-G.; Lüer, L.; Oelkrug, D. J. Chem. Phys. 2002, 116, 8596–8609.

(68) Eisenberg, H. R.; Baer, R. Phys. Chem. Chem. Phys. 2009, 11, 4674–4680.


17


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(69) Refaely-Abramson, S.; Baer, R.; Kronik, L. Phys. Rev. B 2011, 84, 075144. (70) Körzdörfer, T.; Sears, J. S.; Sutton, C.; Brédas, J.-L. J. Chem. Phys. 2011, 135, 204107. (71) Sun, H.; Autschbach, J. J. Chem. Theory Comput. 2014, 10, 1035–1047. (72) Kohler, B. E.; Woehl, J. C. J. Chem. Phys. 1995, 103, 6253–6256. (73) Wood, P.; Samuel, I. D. W.; Schrock, R.; Christensen, R. L. J. Chem. Phys. 2001, 115, 10955–10963. (74) Viani, L.; Corbella, M.; Curutchet, C.; O’Reilly, E. J.; Olaya-Castro, A.; Mennucci, B. Phys. Chem. Chem. Phys. 2014, 16, 16302–16311. (75) de Queiroz, T. B.; Kümmel, S. J. Chem. Phys. 2014, 141, 084303. (76) Lagowski, J. B. J. Mol. Struc.-Theochem 2002, 589-590, 125 – 137. (77) Fincher, C. R.; Ozaki, M.; Tanaka, M.; Peebles, D.; Lauchlan, L.; Heeger, A. J.; MacDiarmid, A. G. Phys. Rev. B 1979, 20, 1589–1602. (78) Gritsenko, O. V.; Schipper, P. R. T.; Baerends, E. J. J. Chem. Phys. 1997, 107, 5007–5015. (79) Cremer, D. Mol. Phys. 2001, 99, 1899– 1940. (80) Yang, S.; Kertesz, M. J. Phys. Chem. A 2006, 110, 9771–9774. (81) Fincher, C. R.; Chen, C. E.; Heeger, A. J.; MacDiarmid, A. G.; Hastings, J. B. Phys. Rev. Lett. 1982, 48, 100–104.


18

Page 18 of 19

Page 19 of 19

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Graphical TOC Entry


19

How Bond Length Alternation and Thermal Disorder Affect the Optical

Recommend Documents