Theoretical Investigation of Static Characterization on Nonlinear

Nov 30, 2006 - The structural distortion of a charged soliton (half-width calculated as L = 13) is shown to be more extended than that of a neutral so...
0 downloads 0 Views 243KB Size
26488

J. Phys. Chem. B 2006, 110, 26488-26496

Theoretical Investigation of Static Characterization on Nonlinear Elementary Excitations in trans-Polyacetylene Haibo Ma, Chungen Liu,* and Yuansheng Jiang† Institute of Theoretical and Computational Chemistry, Key Laboratory of Mesoscopic Chemistry of MOE, School of Chemistry and Chemical Engineering, Nanjing UniVersity, Nanjing 210093, China ReceiVed: January 13, 2006; In Final Form: September 19, 2006

Semiempirical quantum chemical studies on neutral and positively charged H(CH)nH homologues have been performed for systems with n up to 101, where different kinds of nonlinear excitation are found with increasing chain length. The Pariser-Parr-Pople (PPP) model has been employed and solved with the density matrix renormalization group (DMRG) method. The geometrical and electronic distortions induced by defects are obtained and compared with previous theoretical work, indicating that an adequate account of the electron correlation is essential for describing such systems. The structural distortion of a charged soliton (half-width calculated as L ) 13) is shown to be more extended than that of a neutral soliton (L ) 6); the geometrical distortion is even more extended in a polaron. In bipolarons, our calculations show that the coupling of the soliton-antisoliton pair might be longer ranged than expected. The phase transition from a bipolaron to a separated soliton-antisoliton pair occurs when n is close to 100.

I. Introduction Since it was discovered in the 1970s that electrical conductivity of trans-polyacetylene (PA) can be improved significantly through doping with strong electron acceptors or donors,1-3 trans-PA received sustained attention from both academic and industrial researchers. So far it has been widely accepted that the most important electric current carriers in trans-PA are solitons in the polymer chain, which can originate either from self-localized nonlinear excitations in small amount or from charge doping in abundance. As one of the most fundamental and inherent features of quasi-one-dimensional conducting polymer, self-localized nonlinear excitation can lead to various structural defects, such as neutral or charged solitons, polarons, or bipolarons, depending on the properties of the polymer chain as well as on the concentration of charge doping. While an undisturbed trans-PA can adopt either of the two bond length alternation patterns, which are energetically degenerate, opposite bond length alternation patterns can coexist in a disturbed transPA with soliton excitation serving as domain walls. At very low doping concentrations, the defects often appear as separated charged solitons. As the doping concentration goes up, more solitons will be formed, with the distance between a pair of solitons becoming shorter and shorter. When it reaches a stage where the coupling between neighboring solitons becomes strong enough, bipolarons (coupled soliton-antisoliton pairs) are formed. Unlike undisturbed trans-PA chains, which show regular single and double bond alternation due to the Peierls distortion, the above-mentioned defects bring about strong delocalization over several or even tens of (CH) unit cells, delimiting trans-PA regions that present different carboncarbon bond length alternation patterns. With such a great interest in existing conducting conjugated organic polymers and further structural modulations, it is very desirable to acquire a thorough understanding of the role of * Address correspondence to this author. E-mail: [email protected]. † E-mail: [email protected].

defects such as solitons and bipolarons in electrical conductivity. For this reason, in the past few decades, many theoretical studies were devoted to the geometrical and electronic characterization of pure and charged trans-PA oligomers.4-34 The relationship between the topological electronic defect and the geometrical distortion was first established by Su, Schrieffer, and Heeger (SSH).4,5 Their SSH model was based on the Hu¨ckel molecular orbital theory, in which π electrons are coupled to distortions in the polymer backbone through electron-phonon interaction. The SSH model has been found to be very useful in a qualitative interpretation of many physical properties of conjugated systems, but, as an empirical method without the inclusion of explicit electron-electron interactions, it fails to give reliable prediction of the molecular structure and the spin densities distributions in trans-PA. To improve upon this situation, some semiempirical and ab initio self-consistent-field (SCF) methods were employed to study the defects in neutral and charged trans-PA oligomers.7-19 Among the earlier work, the most impressive one is a detailed semiempirical MNDO (modified neglect of diatomic overlap) investigation by Boudreaux et al.7 on neutral and charged solitons, polarons, and bipolarons. Later, several studies emphasized the importance of an adequate account of electron correlations for a more accurate description of trans-PA.20-29 The importance of including electron correlation was also emphasized very recently by Fonseca et a.l30,31 and Champagne et al.32,33 in their studies of singly and doubly charged transPA chains using the second-order Møller-Plesset theory (MP2). It was revealed that neglect of electron correlations will lead to an underestimation of the defect width as well as an overestimation of the averaged bond length alternation (BLA, defined as the difference between consecutive C-C bond lengths) in both neutral and charged trans-PA oligomers. Unfortunately, due to the high computational demand, high level ab initio electron correlation methods such as MP2 cannot be routinely applied to large sized molecules. Meanwhile, density functional theories (DFT), which is a popular choice in many other problems and has a very favorable cost, have been shown to be inappropriate

10.1021/jp0602528 CCC: $33.50 © 2006 American Chemical Society Published on Web 11/30/2006

Nonlinear Elementary Excitations in trans-Polyacetylene

J. Phys. Chem. B, Vol. 110, No. 51, 2006 26489

for treating large conjugated molecules,35-38 most likely as a result of an overestimation of locality with current-generation exchange-correlation functionals. So it remains computationally challenging to take proper account of electron correlation in the study of elementary excitations in trans-PA systems.39 In our recent studies,40-42 we used the semiempirical PariserParr-Pople (PPP) model to study medium to large sized transPA oligomers. The π electron correlations were fully taken into account by employing the density matrix renormalization group (DMRG) method.43,44 With our PPP-DMRG method, reasonable geometries as well as energies for both ground and several lowest excited states were predicted. Moreover, spin density distribution around a neutral defect in trans-PA was also reasonably reproduced. In this work, to make a systematic investigation of various elementary excitations in trans-PA, we extend the application of PPP-DMRG to the study of neutral and charged trans-PA oligomers with up to 101 carbon atoms.

choice for numerical studies of such systems.50 Here we give a very brief description of the DMRG technique as well as how it is utilizied for geometry optimization. In DMRG, the interaction between different fragments is taken into account by working on a “superblock” AB, which is composed of one “system block” A and one “environment block” B. The importance of each basis function in the system block for a high-quality description of a specific state of the superblock is defined through the reduced density matrix (Fii′) ) F, which is defined as

Fii′ )

H)

Ψ)

1

∑∑σ βijaiσ+ajσ + 2∑ij γij(ni - 1)(nj - 1)

(1)

〈ij〉

+ where aiσ and aiσ are the creation and annihilation operators, respectively, which creates or annihilates an electron with spin + σ in the mutually orthogonal atomic π orbital pi; ni ) ∑σaiσ aiσ is the π-electron number operator; βij is the exchange or hopping integral; γij is the effective electron-electron Coulomb repulsion integral; and 〈ij〉 denotes summation restricted to the nearest neighbors. The expression for the resonance integral was taken from the paper of Schulten et al.46 as,

βij ) -2.43 + 3.21(rij - 1.397)

(2)

where rij is the distance between the nuclear sites i and j. γij is parametrized based on Ohno-formula47 through the on-site Hubbard interaction Uii, which is set to be 11.13 eV according to ref 48.

γij ) 14.397

[

]

28.794 + rij2 (Uii + Ujj)2

-1/2

(3)

B. DMRG Method. The difficulties of obtaining the exact solution of the multielectron models for large conjugated systems come from the exponential increase in the number of configurations with the system size. However, physical intuition suggests that not all of the enormous number of configurations contribute equally to the low-lying states, and therefore various methods of truncations are available for cutting down the computational cost while yielding reasonable results. Unlike the conventional post Hatree-Fock (HF) methods built on the SCF wave function, the DMRG method uses position-dependent matrix-product states as many-body wave functions directly to account for the electron correlation effect. It is essentially a recursive transformation of the basis for the full CI Hilbert space, followed by a selection of a limited number of relevant states, which are chosen to be eigenvectors of a particular density matrix F. This algorithm has achieved unprecedented precision in the description of one-dimensional quantum systems, being demonstrated to be able to yield results very close to the exact full CI results.49 It has therefore quickly become the method of

(4)

where ψij is the contribution of a direct product of basis |i〉 in the system block and |j〉 in the environment block to the specific state of interest. In general, this state is denoted as,

II. Computational Methodology A. Model Hamiltonian. The PPP model has been widely used in the context of conjugated polymers45 and it can be written as (energies in eV and distances in Å)

∑j ψijψi′j

∑ij ψij|i〉|j〉

(5)

The diagonalization of F leads to a set of eigenvalues ωR and eigenvectors uR. Following the definition of the density matrix, the states corresponding to larger eigenvalues of density matrix F are more probable configurations for the system block. Accordingly, the m largest eigenstates are retained, spanning a truncated space for the system block. All the operators (e.g., H) are transformed into this new representation through an orthogonal transformation, H′ ) OHO+, where each column in matrix O is a retained eigenvector of F. After the transformation, the Hilbert space for the “system block” A receives a very compact representation, and then at the next step we can increase the size of the system block by one site and repeat the calculation for a new “superblock” AB. A typical real space DMRG computation consists of two stages. The first stage employs the infinite system algorithm, during which, starting from a small fragment of a targeting system, the system is elongated by a few atoms at each iteration, until the superblock reaches the size of the targeting system. The second stage adopts the finite system algorithm, where the basis set for the system and environment blocks are optimized to further improve the accuracy, while the size of the superblock is kept constant. In our scheme, the geometry optimization is incorporated into each iteration of the DMRG steps. For the geometry of transPAs, a planar, all-trans configuration is assumed, and all bond angles are fixed to 120°. The conjugated carbon-carbon bond lengths are calculated according to the widely used bond orderbond length relationship51

rij ) 1.517 - 0.18Pij

(6)

with π bond orders Pij calculated as

〈 |∑

Pij ) Ψ

1

2

σ

|



+ + (aiσ ajσ + ajσ aiσ) Ψ

(7)

By using an initial guess of the molecular geometry, Pij are calculated with the DMRG method, which will be used to obtain a new geometry. In this way, the molecular geometry is optimized iteratively until it is converged. It should be pointed out that, due to the feature of the model Hamiltonian, the geometry optimization procedure used in our work does not account for the chain curvature that emerges in charged long oligoenes.

26490 J. Phys. Chem. B, Vol. 110, No. 51, 2006

Ma et al.

TABLE 1: Comparison of the Carbon-Carbon Bond Lengths (in Å) in C10H122+ Optimized with Various Theoretical Methodsa PPP-DMRG MP2b MP4 CASSCFc HF

1-1′

1-2

2-3

3-4

4-5

1.444 1.452 1.465 1.454 1.456

1.383 1.391 1.382 1.369 1.350

1.430 1.434 1.443 1.431 1.430

1.405 1.416 1.415 1.397 1.393

1.392 1.392 1.390 1.376 1.366

a The backbone carbon atoms are numbered sequentially from the central atom (prime indicates the symmetrically equivalent atom). b From ref 31. c From ref 34.

Figure 1. Calculated carbon-carbon bond length (in units of Å) for C10H122+, using different theoretical methods. MP2 data are from ref 31 and CASSCF data are from ref 34.

III. Results and Discussion A. Validation. In our previous studies, PPP-DMRG results were compared to those of other advanced ab initio electron

correlation methods as well as available experimental data for undoped trans-PA oligomers, and a good agreement was found for many different cases.40-42 Here, we shall demonstrate that PPP-DMRG is also suitable for charge-doped trans-PA oligomers. First we shall compute the optimized geometry of C10H122+ with PPP-DMRG, and compare the results against Hatree-Fock (HF), MP2, MP4, and complete active space self-consistent field (CASSCF) methods. All computed C-C bond lengths along the carbon chain are listed in Table 1 and displayed in Figure 1. Clearly, bond lengths from PPP-DMRG calculations fall within 0.02 Å of those from post-HF ab initio methods, which is a fairly good agreement. Very similar bond length alternation patterns along the chain are also obtained. In contrast, HF calculations underestimate several bond lengths, with discrepancies as large as 0.03-0.04 Å, leading to exaggerated bond length oscillation in the middle of the oligomer chain. Optimized bond lengths and averaged BLA in some other singly or doubly charged trans-PA oligomers from HF, MP2, and PPP-DMRG calculations are presented in Table 2. Just like C10H122+, PPPDMRG results agree well with MP2 ones, while the bond length alternations from HF calculations are quite overestimated for closed shell systems and underestimated for open shell systems. In agreement with previous studies,21-32 the above comparisons once again confirm that electron correlation plays an important role in the description of charged defects in trans-PA, which effectively reduces the magnitude of bond length alternation. What is particularly encouraging to us is that PPP-DMRG appears to work equally in charged trans-PA oligomers. Bond length alternation in the middle of the long oligomer chain could be taken as an approximation for the bond length alternation for the entire polymer chain (noted as ∆r∞), which is one of the important structure parameters for the study of defects in trans-PA chains. Structure defects on the chain usually reduce bond length alternation around the center of the defects,

TABLE 2: Carbon-Carbon Bond Lengths and Averaged BLA (in Å) of Some Charged trans-PA Oligomers Optimized with Different Theoretical Methodsa C11H13+ C13H15+ C15H17+ C12H14+ C14H16+ C16H18+ C12H142+ C14H162+ C16H182+

c

HF/6-31G MP2/6-31Gb PPP-DMRG HF/6-31G MP2/6-31Gb PPP-DMRG HF/6-31G MP2/6-31Gb PPP-DMRG HF/6-31G MP2/6-31G PPP-DMRG HF/6-31G MP2/6-31G PPP-DMRG HF/6-31G MP2/6-31G PPP-DMRG HF/6-31G MP2/6-31Gc PPP-DMRG HF/6-31G MP2/6-31Gc PPP-DMRG HF/6-31G MP2/6-31Gc PPP-DMRG

1-1′

1-2

2-3

3-4

4-5

5-6

1.390 1.406 1.408 1.391 1.408 1.407 1.390 1.406 1.408 1.344 1.389 1.382 1.451 1.443 1.440 1.343 1.389 1.383

1.383 1.404 1.403 1.396 1.411 1.410 1.384 1.405 1.404 1.392 1.403 1.406 1.392 1.404 1.409 1.391 1.403 1.406 1.448 1.444 1.439 1.346 1.390 1.384 1.447 1.439 1.437

1.411 1.421 1.418 1.372 1.399 1.398 1.406 1.416 1.415 1.398 1.408 1.415 1.393 1.393 1.403 1.395 1.402 1.412 1.358 1.396 1.389 1.437 1.435 1.432 1.351 1.393 1.387

1.359 1.389 1.389 1.420 1.426 1.424 1.365 1.395 1.394 1.392 1.378 1.394 1.401 1.411 1.420 1.396 1.387 1.399 1.413 1.422 1.420 1.368 1.400 1.394 1.426 1.427 1.426

1.442 1.447 1.439 1.353 1.386 1.385 1.427 1.430 1.428 1.418 1.442 1.436 1.394 1.373 1.390 1.404 1.415 1.424 1.404 1.421 1.412 1.400 1.413 1.412 1.376 1.404 1.399

1.335 1.368 1.366 1.446 1.450 1.442 1.349 1.383 1.383 1.373 1.345 1.369 1.420 1.447 1.440 1.395 1.370 1.387 1.357 1.386 1.385 1.413 1.427 1.418 1.389 1.408 1.407

6-7

7-8

1.333 1.366 1.365 1.449 1.453 1.444

1.332 1.365 1.363

1.373 1.343 1.367 1.421 1.450 1.443

1.373 1.342 1.365

1.351 1.381 1.380 1.420 1.431 1.422

1.346 1.378 1.377

BLA 0.068 0.047 0.040 0.069 0.045 0.040 0.070 0.045 0.040 0.017 0.040 0.028 0.015 0.042 0.030 0.014 0.040 0.031 0.061 0.033 0.034 0.062 0.034 0.033 0.063 0.033 0.033

a The carbon backbone is numbered sequentially from the central atom (prime indicates the symmetrically equivalent atom). b From ref 30. From ref 31.

1.358 1.454 1.358 1.454 1.371 1.358 1.454 1.371 1.448 1.359 1.454 1.371 1.448 1.374 1.359 1.454 1.371 1.448 1.374 1.446 1.359 1.454 1.371 1.448 1.374 1.446 1.375 1.359 1.454 1.371 1.448 1.374 1.446 1.376 1.445 1.359 1.454 1.372 1.447 1.375 1.445 1.376 1.444 1.377 1.359 1.454 1.372 1.447 1.375 1.445 1.377 1.444 1.377 1.443 1.359 1.454 1.372 1.447 1.376 1.444 1.377 1.443 1.378 1.443 1.378 1.359 1.453 1.373 1.446 1.377 1.443 1.378 1.442 1.379 1.442 1.379 1.441 a

1.439 1.393 1.427 1.398 1.423 1.400 1.421 1.401 1.420 1.401 1.419 1.402 1.419 1.402 1.419 1.402 1.419 1.402 1.404 1.418 1.406 1.416 1.407 1.415 1.407 1.414 1.407 1.414 1.407 1.413 1.408 1.413 1.408 1.413 1.408 1.413 C7H9 C9H11 C11H13 C13H25 C15H17 C17H19 C19H21 C21H23 C23H25 C25H27 C27H29 C29H31 C31H33 C33H35 C35H37 C37H39 C39H41 C41H43

The carbon backbone is numbered sequentially from the central atom.

1.360 1.453 1.374 1.445 1.378 1.442 1.379 1.441 1.380 1.440 1.381 1.440 1.381 1.360 1.452 1.375 1.444 1.379 1.441 1.381 1.439 1.382 1.439 1.382 1.438 1.383 1.438 1.361 1.451 1.376 1.442 1.381 1.439 1.383 1.437 1.384 1.436 1.385 1.436 1.385 1.436 1.385 1.363 1.450 1.378 1.440 1.383 1.436 1.386 1.434 1.387 1.433 1.388 1.433 1.388 1.432 1.388 1.432 1.365 1.448 1.381 1.437 1.387 1.433 1.389 1.431 1.391 1.430 1.392 1.429 1.392 1.429 1.392 1.428 1.392

19-20 18-19 17-18 16-17 15-16 14-15 13-14 11-12 10-11 9-10 8-9 7-8 6-7 5-6 4-5 3-4 2-3 1-2

where ∆rn is the BLA around the nth atom counted from the center of a soliton. In our calculations, ∆r∞ is found to be 0.069 Å with the PPP-DMRG method. L is the half-width of the soliton, whose value yields a best fitting of ∆rn against n. Optimized C-C bond lengths for small to medium sized trans-PA radicals CnHn+2 (n ) 7, 9, ..., 41) are listed in Table 3. Only half of the bonds are tabulated, because the soliton distribution is symmetrical with respect to the center of the chain. Since the soliton district serves as a domain wall to separate the chain into two halves, where the terminal carboncarbon bonds are always double bonds, the oligomers can be classified into two series: one is composed of 4n + 1 carbon atoms, and the other is composed of 4n + 3 carbon atoms. In the C4n+1H4n+3 series, the centermost carbon-carbon bond is slightly longer than its nearest neighbors; in the C4n+3H4n+5 series, the centermost one is shorter than its nearest neighbors. In both cases the bond length alternations will increase gradually when moving away from the middle of the chain, which can be well described with eq 8. The formation of a soliton leads to extra C-C bond delocalization around the center of this soliton. At first, the conjugated region will be extended gradually with increasing chain length. But once reaching a certain size, the width of the soliton converges to a finite value. The soliton width is predicted by fitting BLA vs n for C61H63 and C101H103, as shown in Figure 2. Due to the terminal effects, slight deviations from the hyperbolic curve can be noticed on several ending carbon atoms. The data points in these regions were excluded in the data fitting. Several integral values of L for C101H103 were tried, and when L equals 6 the best fitting was obtained. This optimal value for L was also supported by results for moderately sized radicals, such as C61H63. We noticed that different values for the half-width were suggested in previous studies: for example, MNDO SCF calculations yielded more localized soliton with L ) 3,7 while SSH4,5 and HF9 calculations led to slightly more delocalized solitons with L ) 7.

12-13

(8)

TABLE 3: Optimized Carbon-Carbon Bond Lengths (in Å) for Neutral trans-PA Oligomersa

(Ln)

∆rn ) (-1)n∆r∞ tanh

20-21

and the extent and region of such deviations is a reflection of the nature of the defect. In our recent paper,42 we applied the PPP-DMRG method to a series of oligoenes C2nH2n+2 (n ) 2, 3,. .., 50), and ∆r∞ was predicted as 0.069 Å, which is in reasonable accordance with the value (0.08 Å) measured by nutation NMR experiment.52 B. Neutral Solitons. A neutral soliton automatically appears on a conjugated segment of a trans-PA chain with an odd number of carbon atoms. In chemical terms, it is called a neutral free radical (charge Q ) 0 and spin S ) 1/2). In our previous work, we correctly predicted the spin distribution picture in neutral trans-PA with the PPP-DMRG method.41 In this work, we will focus on geometrical characters of neutral solitons in the trans-PA chain. Starting from one side of a soliton, the length of the double (single) bonds is known to gradually increase (decrease) such that, when arriving at the other side of the defect, the single and double bond characters are fully exchanged. Naturally, in the middle of the chain there is a critical point where two neighboring carbon-carbon bonds have exactly the same length. If we take this point as a phase boundary, then on one side of the boundary the BLAs are positive; on the other side they are negative. This geometrical distortion from the regular single and double bond alternation is generally written as a hyperbolic tangent function4,5

1.358

J. Phys. Chem. B, Vol. 110, No. 51, 2006 26491

1.368 1.444 1.386 1.433 1.391 1.429 1.394 1.426 1.395 1.425 1.396 1.425 1.397 1.424 1.397 1.424 1.397 1.424

Nonlinear Elementary Excitations in trans-Polyacetylene

Ma et al.

1.359 1.453 1.373 1.445 1.376 1.360 1.453 1.373 1.445 1.377 1.442

15-16

1.360 1.452 1.373 1.444 1.378 1.441 1.380

14-15

1.360 1.452 1.374 1.444 1.378 1.440 1.380 1.438

13-14

1.360 1.451 1.374 1.443 1.379 1.439 1.381 1.437 1.383

12-13

1.361 1.451 1.375 1.442 1.380 1.438 1.382 1.435 1.384 1.434

11-12

1.361 1.450 1.376 1.441 1.381 1.436 1.383 1.434 1.385 1.432 1.386

10-11

a

1.363 1.446 1.379 1.436 1.385 1.431 1.388 1.428 1.390 1.426 1.391 1.425 1.392 1.424 1.365 1.444 1.381 1.434 1.387 1.429 1.390 1.426 1.392 1.424 1.393 1.423 1.394 1.422 1.394 1.366 1.442 1.383 1.431 1.389 1.426 1.392 1.423 1.394 1.421 1.395 1.420 1.396 1.420 1.396 1.419 1.369 1.439 1.385 1.428 1.391 1.423 1.394 1.420 1.396 1.419 1.397 1.418 1.398 1.417 1.399 1.417 1.399 1.372 1.434 1.389 1.424 1.394 1.419 1.397 1.417 1.398 1.416 1.400 1.415 1.400 1.414 1.401 1.414 1.401 1.414 1.428 1.393 1.418 1.398 1.415 1.400 1.413 1.401 1.412 1.402 1.412 1.403 1.412 1.403 1.411 1.404 1.411 1.404 1.401 1.411 1.403 1.410 1.404 1.409 1.405 1.409 1.405 1.409 1.406 1.409 1.406 1.409 1.406 1.409 1.406 1.409 C7H9 C9H11+ C11H13+ C13H15+ C15H17+ C17H19+ C19H21+ C21H23+ C23H25+ C25H27+ C27H29+ C29H31+ C31H33+ C33H35+ C35H37+ C37H39+ C39H41+ C41H43+

The carbon backbone is numbered sequentially from the central atom.

1.362 1.448 1.378 1.438 1.383 1.433 1.386 1.430 1.388 1.428 1.389 1.427 1.390

1.362 1.449 1.377 1.439 1.382 1.435 1.385 1.432 1.386 1.430 1.388 1.429

9-10 8-9 7-8 6-7 5-6 4-5 3-4 2-3 1-2

C. Charged Solitons. In addition to neutral solitons, solitons in trans-PA can also be positively or negatively charged. They are different from neutral ones in that these defects are spinless (S ) 0). Charged solitons can be generated through reduction or oxidation of trans-PA, especially at low concentrations. Earlier studies revealed that geometrical character of negatively and positively charged solitons can be slightly different, because a positive charge on the chain tends to provoke the planar sp2 hybridization, while a negative charge tends to invoke the nonplanar sp3 hybridization.7 Unfortunately, as a parametrized semiempirical method for π conjugated systems, the PPP model cannot be used to distinguish between positively and negatively charged solitons,7 so our results for charged defects in transPA only apply to positively charged, planar species. Optimized carbon-carbon bond lengths for some small to medium sized singly charged trans-PA oligomers with an odd number of carbon atoms CnHn+2+ (n ) 7, 9, ..., 41) are listed in Table 4. BLA values at each carbon atom on the chain are displayed in Figure 3 for C61H63+ and C101H103+, where a hyperbolic tangent function fitting is also displayed and the halfwidth of the soliton is determined to be L ) 13. MNDO SCF results7 also suggested more extended charged solitons (the halfwidth increased from L ) 3 for neutral solitions to L ) 5 for charged solitons). Villar et al.’s RHF/6-31G calculation9 indicated that both neutral and charged solitons have similar half-width (L ) 7). In the more recent work of Champagne et al., different values of L ) 5 (HF/STO-3G), 8 (HF/6-31G), 7 (HF/6-31G*), 9 (MP2/STO-3G), 18 (MP2/6-31G), and 22 (B3LYP/STO-3G) were reported.32 Obviously, these predicted L values depend strongly on the basis set and on the way electron correlation is being taken into account, with the two factors balanced in a good ab initio model. Our semiempirical calculations yielded a value of L ) 13, lying between those from HF calculations and those from MP2 and B3LYP calculations. The underestimation of L by ab initio HF and MNDO SCF methods, with respect to those obtained from PPPDMRG, MP2, as well as B3LYP calculations, reflects that that electron correlation tends to reduce bond alternations in charged trans-PA. It was pointed out that neglecting the electron correlation in HF calculations led to less π electron conjugation, which in turn resulted in less relaxation of solitons, or in other words, an underestimation of the length of the geometrical distortion region.12 Since the available MP2 and B3LYP results were obtained at small basis sets, it is inappropriate for us to comment on the reliability of the values determined with these

+

Figure 2. Bond length alternation (in units of Å) pattern around the neutral defect site for C61H63 and C101H103. The solid curve represents a hyperbolic function given by ∆rn ) ∆r∞ tanh(n/L) with L ) 6.

TABLE 4: Optimized Carbon-Carbon Bond Lengths (in Å) for Singly Charged trans-PA Oligomers Having Odd Carbon Atomsa

16-17

1.359 1.453 1.372 1.446

17-18

1.359 1.453 1.372

18-19

1.359 1.453

19-20

1.359

20-21

26492 J. Phys. Chem. B, Vol. 110, No. 51, 2006

Nonlinear Elementary Excitations in trans-Polyacetylene

Figure 3. Bond length alternation (in units of Å) pattern around the charged defect site for C61H63+ and C101H103+. The solid curve is drawn with the hyperbolic function given by ∆rn ) ∆r∞ tanh(n/L) with L ) 13.

Figure 4. Net charge per CH unit for C61H63+ and C101H103+.

methods. We note earlier calculations37 indicated that DFT calculations tended to overestimate the spin distribution region of the neutral solitons in the trans-PA chain. For all charged solitons under consideration, an oscillation of net charge per CH unit is observed when approaching the defect site, with the largest charge on the central unit and then gradually decreasing charges toward the chain ends. This damping charge wave is shown for C61H63+ and C101H103+ in Figure 4. Similar phenomena have been reported,7,9,32,33 and can be easily understood based on chemical intuition. For instance, Boudreaux et al.7 reasoned that a positively charged carbon would induce a negative charge on the nearby carbon atom, resulting in polarization of the π bonds next to the charged site. When different chains lie parallel to each other, each with a set of sufficiently large electric dipoles, the interchain coupling can become very strong, which can then lead to interchain hopping, another important mode of charge transportation (besides moving along the chains). As a consequence of large soliton width, the charge is spread over a broad region while the amplitude of the charge oscillation is very small, which is consistent with Champagne et al.’s MP2 results,32,33 but contrary to Boudreaux et al.’s MNDO SCF7 and Villar et al.’s HF9 results. D. Polarons. Just like charged solitons in trans-PA, a polaron could also be generated by chemical reduction or oxidation. A

J. Phys. Chem. B, Vol. 110, No. 51, 2006 26493 polaron is different from a charged soliton in that it involves an unpaired electron and thus a nonzero spin S ) 1/2. A polaron is also different from a neutral soliton in that it spreads over an even number of distorted atoms, in other words, a polaron is symmetrical with respect to the bond in the middle of the chain. A polaron can be roughly described as a bound state of a neutral and a charged soliton. We list in Table 5 optimized carbon-carbon bond lengths for trans-PA cations CnHn+2+ (n ) 6, 8, ..., 40), and an illustration of BLA at each carbon atom position is presented in Figure 5 for C60H62+ and C100H102+. It could be seen that the bond length alternation decreases when approaching the center of the polaron; however, unlike in the case of a soliton, the bond length alternation pattern in a polaron is not reversed when going from one side to the othersBLAs are positive on both sides. It is also found that, when moving away from the middle toward either end of the chain, BLA increases slower than in the case of charged solitons, which implies a polaron is more extended along the chain than a charged soliton. The charge distribution patterns for C60H62+ and C100H102+ are shown in Figure 6. Interestingly, as noticed previously by Boudreaux et al.,7 while the central atom in a charged soliton possesses the largest charge, in a polaron the largest positive charge is found away from the center of the chain. On the basis of the MNDO SCF calculations, Boudreaux predicted that the largest positive charge located on the sixth atom away from the middle and the largest negative charges are even farther away. In our work, the largest positive charge is found at eight sites away from the center. The charge sign oscillation is also different from that of a charged soliton, instead of beginning from the center of a soliton, Boudreaux et al. predicted that it starts at the third site away from the center of the chain, while our calculation shows that it starts at the eleventh site. E. Bipolarons and Soliton-Antisoliton Pairs. A bipolaron is defined as a pair of charges (a dication in this work) associated with a strong local lattice distortion, and can also be regarded as a confined soliton-antisoliton pair. It can be generated by taking away one more electron from a singly charged trans-PA chain. Normally, the defect excitation in doubly charged structures ranges from a bipolaron in short oligomers to a separated soliton-antisoliton pair in long oligomers.7,31 From Table 6, which lists carbon-carbon bond lengths of bipolarons from C6H82+ to C40H422+, it could be seen that there exist two domain walls (centers for soliton and antisolitons) symmetrically located away from the center of the chain. BLA of each carbon atom in C60H622+ and C100H1022+, as well as two curves which simulate isolated charged soliton, are presented in Figure 7. From Table 6 and Figure 7, it can be seen that BLA in the middle of the chain increases very slowly, they are calculated as 0.057 Å, and 0.063 Å in C40H422+ and C60H622+, respectively. Both values are smaller than ∆r∞ (0.069 Å), indicating that the soliton-antisoliton pair are still not quite independent even in oligomers as long as C60H622+. In contrast, BLAs in the centermost region in C100H1022+ have reached ∆r∞ (0.069 Å) and become very stable in the central region of 20 carbon atoms, implying that the transition from bipolaron to completely separated soliton-antisoliton pair occurs at a chain length with no more than 100 carbon atoms. This is also supported by Figure 7, which shows that the C100H1022+ curve is a good match to the ∆rn ) ∆r∞ tanh(n/13) curve for charged solitons curve, whereas the C60H622+ curve deviates significantly. These results are quite different from MNDO SCF7 and HF9 ones. MNDO SCF calculations showed that the BLAs in the centermost region of the C40H422+ chain already reached

The carbon backbone is numbered sequentially from the central atom (a prime indicates symmetry equivalent atoms). a

1-2

2-3

3-4

4-5

5-6

6-7

7-8

8-9

9-10

10-11

11-12

12-13

13-14

14-15

15-16

16-17

17-18

18-19

19-20 1.411 1.409 1.408 1.408 1.407 1.408 1.407 1.408 1.406 1.409 1.406 1.409 1.406 1.410 1.405 1.410 1.405 1.411

1.417 1.405 1.411 1.406 1.409 1.406 1.409 1.406 1.409 1.406 1.409 1.405 1.410 1.405 1.410 1.405 1.411 1.405

1.386 1.425 1.399 1.415 1.403 1.412 1.404 1.411 1.404 1.411 1.404 1.411 1.404 1.411 1.404 1.411 1.404 1.412

1.377 1.431 1.394 1.420 1.399 1.416 1.401 1.414 1.402 1.413 1.403 1.413 1.403 1.413 1.403 1.413 1.403

1.372 1.436 1.390 1.424 1.396 1.419 1.399 1.416 1.400 1.415 1.401 1.415 1.401 1.414 1.401 1.414

1.369 1.440 1.387 1.427 1.393 1.422 1.396 1.419 1.398 1.418 1.399 1.417 1.399 1.417 1.400

1.367 1.443 1.384 1.431 1.391 1.425 1.394 1.422 1.396 1.420 1.397 1.419 1.398 1.419

1.365 1.445 1.382 1.433 1.388 1.428 1.392 1.425 1.394 1.423 1.395 1.422 1.396

1.364 1.447 1.380 1.435 1.387 1.430 1.390 1.427 1.392 1.425 1.393 1.424

1.363 1.448 1.379 1.437 1.384 1.432 1.388 1.429 1.390 1.427 1.391

1.362 1.449 1.377 1.439 1.383 1.434 1.386 1.431 1.388 1.430

1.361 1.450 1.376 1.440 1.382 1.436 1.384 1.433 1.386

1.361 1.451 1.375 1.442 1.380 1.437 1.383 1.435

1.360 1.451 1.375 1.443 1.379 1.439 1.382

1.360 1.452 1.374 1.444 1.379 1.440

1.360 1.452 1.374 1.444 1.378

1.360 1.453 1.373 1.445

1.359 1.453 1.373

1.359 1.453

1.359 C6H8 C8H10+ C10H12+ C12H14+ C14H16+ C16H18+ C18H20+ C20H22+ C22H24+ C24H26+ C26H28+ C28H30+ C30H32+ C32H34+ C34H36+ C36H38+ C38H40+ C40H42+

1-1′

26494 J. Phys. Chem. B, Vol. 110, No. 51, 2006

TABLE 5: Optimized Carbon-Carbon Bond Lengths (in Å) for Singly Charged trans-PA Oligomers Having Even Carbon Atomsa +

Ma et al.

Figure 5. Bond length alternation (in units of Å) pattern around the charged defect site for C60H62+ and C100H102+.

Figure 6. Net charge per CH unit for C60H62+ and C100H102+.

MNDO-determined ∆r∞ (0.106 Å), while HF calculations predicted that interaction of the soliton-antisoliton pair might start to disappear with an even smaller C22H242+, where BLA of the centermost atom was computed to be 0.107 Å as opposed to a HF-based ∆r∞ (0.112 Å). The disagreement between our conclusion and those from MNDO SCF and ab initio HF indicates that the coupling between the soliton-antisoliton pair might be stronger and more extended than what one expects based on SCF calculations. This viewpoint was supported by MP2 calculations of Oliveira et al. for doubly charged transPA oligomers with lengths up to 24 carbon atoms.31 The BLA of the centermost atom was also found to change slowly with the size of the oligomer. For C24H262+, it was computed to be 0.048 Å, which is much smaller than ∆r∞ (0.08 Å) calculated with MP2. Further understanding of the interaction between the solitonantisoliton pair could be reached by inspecting the charge distribution patterns for C60H622+ and C100H1022+ in Figure 8. Two positive charge peaks are located ca. 34 (for C60H622+) or 62 (for C100H1022+) carbon atoms apart, and the charge population on the central carbon atom decreases from about 0.01 (in C60H622+) to 0.001 (in C100H1022+). Different magnitudes of these two values imply again that the phase transition from bipolaron to separated soliton-antisoliton pair has been completed at the length of 100 carbon atoms.

1.363 1.446

18-19

1.364 1.444 1.382 1.430 1.365 1.443 1.383 1.429 1.392

1.364 1.445 1.381

17-18 16-17 15-16

Figure 7. Bond length alternation (in units of Å) pattern for C60H622+ and C100H1022+. Curve 1 is drawn with the hyperbolic function given by ∆rn ) ∆r∞ tanh((n - 17)/13) and Curve 2 is drawn with the function ∆rn ) ∆r∞ tanh((n - 31)/13).

1.366 1.442 1.384 1.427 1.393 1.419

14-15

1.366 1.441 1.386 1.426 1.394 1.417 1.400

13-14

a

1.369 1.437 1.389 1.422 1.398 1.413 1.404 1.408 1.408 1.370 1.435 1.391 1.419 1.400 1.411 1.406 1.406 1.410 1.402 1.372 1.433 1.393 1.417 1.402 1.409 1.408 1.404 1.413 1.400 1.416 1.374 1.430 1.395 1.414 1.405 1.406 1.410 1.401 1.415 1.398 1.419 1.395 1.377 1.427 1.398 1.411 1.408 1.403 1.413 1.399 1.418 1.396 1.421 1.394 1.424 1.380 1.423 1.402 1.407 1.411 1.401 1.416 1.397 1.420 1.394 1.423 1.392 1.426 1.390 1.385 1.418 1.407 1.403 1.415 1.397 1.420 1.394 1.423 1.392 1.426 1.390 1.428 1.388 1.430 1.392 1.412 1.412 1.399 1.420 1.394 1.424 1.391 1.427 1.389 1.429 1.388 1.431 1.387 1.432 1.385 1.401 1.405 1.420 1.394 1.426 1.390 1.428 1.389 1.430 1.387 1.432 1.386 1.433 1.385 1.434 1.384 1.435 1.418 1.396 1.430 1.389 1.432 1.387 1.433 1.386 1.433 1.385 1.434 1.385 1.435 1.384 1.436 1.383 1.436 1.382 1.386 1.442 1.383 1.439 1.384 1.437 1.384 1.436 1.384 1.436 1.384 1.436 1.383 1.437 1.383 1.437 1.382 1.438 1.455 1.381 1.444 1.382 1.440 1.383 1.438 1.383 1.437 1.383 1.437 1.383 1.437 1.383 1.437 1.382 1.438 1.381 C6H8 C8H102+ C10H122+ C12H142+ C14H162+ C16H182+ C18H202+ C20H222+ C22H242+ C24H262+ C26H282+ C28H302+ C30H322+ C32H342+ C34H362+ C36H382+ C38H402+ C40H422+

The carbon backbone is numbered sequentially from the central atom (a prime indicates symmetry equivalent atoms).

1.367 1.439 1.387 1.424 1.396 1.415 1.402 1.410

12-13 11-12 10-11 9-10 8-9 7-8 6-7 5-6 4-5 3-4 2-3 1-2 1-1′

2+

TABLE 6: Optimized Carbon-Carbon Bond Lengths (in Å) for Doubly Charged trans-PA Oligomersa

J. Phys. Chem. B, Vol. 110, No. 51, 2006 26495

1.363

19-20

Nonlinear Elementary Excitations in trans-Polyacetylene

Figure 8. Net charge per CH unit for C60H622+ and C100H1022+.

IV. Conclusion Structural distortions induced by neutral and charged solitons, polaron, and bipolaron in trans-PA oligomers with up to 101 carbon atoms were studied with the PPP-DMRG method. Within this method, the π electron correlations receive nearly full account by virtue of DMRG, which yielded wider distribution of all kinds of defects in trans-PA than previous SCF theoretical work. It is also found that structural distortion of the charged solitons (L ) 13) is much more extended than that of the neutral solitons (L ) 6). Meanwhile, in polaron, the geometrical distortion is more extended than in solitons and the largest charge population appears symmetrically on two sites which are several carbon atoms away from the center of the chain. In the case of doubly charged oligomers, our calculations covered the transition from bipolaron to separated soliton-antisoliton pair with increasing oligomer size. We found that this phase transition has been completed for oligomers with less than 100 carbon atoms. However, it should be mentioned that properties of charge doped trans-PA can vary greatly with counterions used in charge doping; as discussed in refs 13 ,19, and 32 a more accurate description of the phase transition should take into account the nature of the counterions. Acknowledgment. This work is supported by China NSF under Grant Nos. 20273030, 20433020, and 20573051. Portions of the computations were done on Origin 3800 and Dawning

26496 J. Phys. Chem. B, Vol. 110, No. 51, 2006 3000 computers at Nanjing University. C.L. is very grateful to Dr. Y. H. Shao for correcting and revising the manuscript. References and Notes (1) Chiang, C. K.; Fincher, C. R., Jr.; Park, Y. W.; Heeger, A. J.; Shirakawa, H.; Louis, E. J.; Gau, S. C.; MacDiarmid, A. G. Phys. ReV. Lett. 1977, 39, 1098. (2) Shirakawa, H.; Louis, E. J.; MacDiarmid, A. G.; Chiang, C. K.; Heeger, A. J. J. Chem. Soc., Chem. Commun. 1977, 16, 579. (3) Chiang, C. K.; Druy, M. A.; Gau, S. C.; Heeger, A. J.; Louis, E. J.; MacDiarmid, A. G.; Park, Y. W.; Shirakawa, H. J. Am. Chem. Soc. 1978, 100, 1013. (4) Su, W. P.; Schrieffer, J. R.; Heeger, A. J. Phys. ReV. Lett. 1979, 42, 1698. (5) Su, W. P.; Schrieffer, J. R.; Heeger, A. J. Phys. ReV. B 1980, 22, 2099. (6) Heeger, A. J.; Kivelson, S.; Schrieffer, J. R.; Su, W. P. ReV. Mod. Phys. 1988, 60, 781. (7) Boudreaux, D. S.; Chance, R. R.; Bre´das, J. L.; Silbey, R. Phys. ReV. B 1983, 28, 6927. (8) Fo¨rner, W.; Wang, C. L.; Martino, F.; Ladik, J. Phys. ReV. B 1988, 37, 4567. (9) Villar, H. O.; Dupuis, M.; Clementi, E. Phys. ReV. B 1988, 37, 2520. (10) Villar, H. O.; Dupuis, M.; Watts, J. D.; Hurst, G. J. B.; Clementi, E. J. Chem. Phys. 1988, 88, 1003. (11) Kirtman, B.; Toto, J. L.; Robins, K. A.; Hasan, M. J. Chem. Phys. 1995, 102, 5350. (12) Champagne, B.; Deumens, E.; O ¨ hrn, Y. J. Chem. Phys. 1997, 107, 5433. (13) Champagne, B.; Spassova, M.; Jadin, J.; Kirtman, B. J. Chem. Phys. 2002, 116, 3935. (14) An, Z.; Wong, K. Y. J. Chem. Phys. 2001, 114, 1010. (15) An, Z.; Wong, K. Y. J. Chem. Phys. 2003, 119, 1204. (16) de Melo, C. P.; Fonseca, T. L. Chem. Phys. Lett. 1996, 261, 28. (17) Larsson, S.; Rodrı´guez-Monge, L. Int. J. Quantum Chem. 1997, 63, 655. (18) Stafstro¨m, S.; Bre´das, J. L.; Lo¨gdlund, M.; Salaneck, W. R. J. Chem. Phys. 1993, 99, 7938. (19) Spassova, M.; Champagne, B.; Kirtman, B. Chem. Phys. Lett. 2005, 412, 217. (20) Yonemitsu, K.; Ono, Y.; Wada, Y. J. Phys. Soc. Jpn. 1988, 57, 3875. (21) Sim, F.; Salahub, D. R.; Chin, S.; Dupuis, M. J. Chem. Phys. 1991, 95, 4317. (22) Suhai, S. Int. J. Quantum Chem. 1992, 42, 193. (23) Villar, H. O.; Dupuis, M. Theor. Chim. Acta 1992, 83, 155.

Ma et al. (24) Bally, T.; Roth, K.; Tang, W.; Schrock, R. R.; Knoll, K.; Park, L. Y. J. Am. Chem. Soc. 1992, 114, 2440. (25) Rodrı´guez-Monge, L.; Larsson, S. J. Chem. Phys. 1995, 102, 7106. (26) Hirata, S.; Torii, H.; Tasumi, M. J. Chem. Phys. 1995, 103, 8964. (27) Fu¨lscher, M. P.; Matzinger, S.; Bally, T. Chem. Phys. Lett. 1995, 236, 167. (28) Guo, H.; Paldus, J. Int. J. Quantum Chem. 1997, 63, 345. (29) Perpe`te, E. A.; Champagne, B. J. Mol. Struct. (THEOCHEM) 1999, 487, 39. (30) Fonseca, T. L.; Castro, M. A.; Cunha, C.; Amaral, O. A. V. Synth. Met. 2001, 123, 11. (31) Oliveira, L. N.; Amaral, O. A. V.; Castro, M. A.; Fonseca, T. L. Chem. Phys. 2003, 289, 221. (32) Champagne, B.; Spassova, M. Phys. Chem. Chem. Phys. 2004, 6, 3167. (33) Monev, V.; Spassova, M.; Champagne, B. Int. J. Quantum Chem. 2005, 104, 354. (34) Kawshima, Y.; Nakayama, K.; Nakano, H.; Hirao, K. Chem. Phys. Lett. 1997, 267, 82. (35) Champagne, B.; Perpe`te, E. A.; van Gisbergen, S. J. A.; Baerends, E. J.; Snijders, J. G.; Soubra-Ghaoui, C.; Robins, K. A.; Kirtman, B. J. Chem. Phys. 1998, 109, 10489. (36) van Gisbergen, S. J. A.; Schipper, P. R. T.; Gritsenko, O. V.; Baerends, E. J.; Snijders, J. G.; Champagne, B.; Kirtman, B. Phys. ReV. Lett. 1999, 83, 694. (37) Bally, T.; Hrovat, D. A.; Borden, W. T. Phys. Chem. Chem. Phys. 2000, 2, 3363. (38) Cai, Z.-L.; Sendt, K.; Reimers, J. R. J. Chem. Phys. 2002, 117, 5543. (39) Heeger, A. J. J. Phys. Chem. B 2001, 105, 8475. (40) Ma, H.; Liu, C.; Jiang, Y. J. Chem. Phys. 2004, 120, 9316. (41) Ma, H.; Cai, F.; Liu, C.; Jiang, Y. J. Chem. Phys. 2005, 122, 104909. (42) Ma, H.; Liu, C.; Jiang, Y. J. Chem. Phys. 2005, 123, 84303. (43) White, S. R. Phys. ReV. Lett. 1992, 69, 2863. (44) White, S. R. Phys. ReV. B 1993, 48, 10345. (45) Linderberg, J.; O ¨ hrn, Y. J. Chem. Phys. 1968, 49, 716. (46) Schulten, K.; Ohmine, I.; Karplus, M. J. Chem. Phys. 1976, 64, 4422. (47) Ohno, K. Theor. Chim. Acta 1964, 2, 219. (48) Dewar, M. J. S.; de Llano, C. J. Am. Chem. Soc. 1969, 91, 789. (49) Yaron, D.; Moore, E. E.; Shuai, Z.; Bre´das, J. L. J. Chem. Phys. 1998, 108, 7451. (50) Schollwo¨ck, U. ReV. Mod. Phys. 2005, 77, 259 and references therein. (51) Coulson, C. A.; Golebiewski, A. Proc. Phys. Soc., London 1961, 78, 1310. (52) Yannoni, C. S.; Clarke, T. C. Phys. ReV. Lett. 1983, 51, 1191.