Ab Initio Calculation of the Dispersion Interaction between a

Oct 22, 2009 - Present address: Institut für Physik, Technische Universität Ilmenau, Postbox 10 05 65, D-98684 Ilmenau, Germany. Cite this:J. Phys. Ch...
0 downloads 0 Views 2MB Size
Download PDF
J. Phys. Chem. C 2009, 113, 19897–19904

19897

Ab Initio Calculation of the Dispersion Interaction between a Polyaromatic Molecule and a Noble Metal Substrate: PTCDA on Ag(110) Afshin Abbasi†,‡,§ and Reinhard Scholz*,†,| Walter Schottky Institut und Physik Department, Technische UniVersita¨t Mu¨nchen, D-85748 Garching, Germany, and Institut fu¨r Physik, Technische UniVersita¨t Chemnitz, D-09107 Chemnitz, Germany ReceiVed: March 17, 2009; ReVised Manuscript ReceiVed: September 23, 2009

The geometry of PTCDA adsorbed on a (110)-oriented silver crystal is optimized on a finite metal cluster used as a rigid substrate. In a comparison between second-order Møller-Plesset perturbation theory (MP2) and density functional theory (DFT), we find pronounced differences; MP2 gives a nearly flat adsorbate, in agreement with the geometries deduced from X-ray standing wave studies on similar substrates, whereas the lack of dispersion interactions in DFT results in a strongly bent geometry. In MP2, the van-der-Waals attraction counterbalances the overlap repulsion, resulting in an average height of the PTCDA adsorbate above the topmost substrate layer of 2.68 Å for C atoms, comparing favorably with recent measurements based on X-ray standing wave techniques. The carboxylic oxygens interact more strongly with the substrate than the anhydride groups, so that we predict a height difference of 0.13 Å between both types of oxygen atoms. Introduction From a technological point of view, the adsorption geometry and the resulting energetics at the molecule-metal interface define the contact properties in devices like organic field effect transistors and light-emitting diodes.1 Therefore, various electronic and spectroscopic methods have been applied to detailed investigations of such interfaces.2 The aromatic compound PTCDA (3,4,9,10-perylene tetracarboxylic dianhydride) analyzed in the present work is among the most intensively studied molecular semiconductors, and its tendency to form well-ordered films has helped to elucidate various details of its electronic and spectroscopic properties. On (111)-oriented noble metals, PTCDA forms an adsorbate layer with two molecules per unit cell, in a geometry resembling the bulk phases.2–6 The height of a PTCDA molecule above the topmost monolayer of the Ag(111) substrate has been determined with X-ray standing wave (XSW) techniques,5–7 and it was argued that DFT calculations reproduce a similar geometry.5 The latter finding resulted from overcompressed basis orbitals reducing the overlap repulsion, and a revised calculation with more extended basis orbitals gave only a shallow potential minimum at a much too large distance between adsorbate and substrate.8,9 These well-known shortcomings of DFT when applied to systems involving dispersion interactions can be overcome by semiempirical corrections to DFT or by quantum chemical techniques including the main contributions to the vander-Waals interaction at an ab initio level. More reactive surfaces like Ag(110) induce a much stronger interaction with organic adsorbates, so that both PTCDA and the similar compound NTCDA (1,4,5,8-naphthalene tetracar* To whom correspondence should be addressed. E-mail: reinhard. [email protected]; [email protected]. † Technische Universita¨t Mu¨nchen. ‡ Technische Universita¨t Chemnitz. § Present address: Max-Planck-Institut fu¨r Eisenforschung GmbH, MaxPlanck-Strasse 1, D-40237 Du¨sseldorf, Germany. | Present address: Institut fu¨r Physik, Technische Universita¨t Ilmenau, Postbox 10 05 65, D-98684 Ilmenau, Germany.

boxylic dianhydride) prefer an adsorbate unit cell containing only a single molecule.10,11 For PTCDA, the adsorption site was determined by manipulation of Ag adatoms,12 whereas the most favorable position of chemisorbed NTCDA can be deduced from DFT calculations.13 In the bulk phases and in the surface unit cell of PTCDA on Ag(111), the area per molecule is only 119 Å2, so that rather short intermolecular C-H · · · O hydrogen bridges give a substantial contribution to the cohesive energy.14 For PTCDA chemisorbed to Ag(110), the area per molecule increases to a much larger value of 141 Å2. Placing a molecule in its optimized monomer geometry into each adsorbate unit cell, it turns out that the distance between oxygen and hydrogen atoms on different molecules increases to about 3.2 Å. At these large distances, hydrogen bridges give no significant intermolecular interaction, reducing the latter to the Coulomb potential between electric multipoles. Therefore, due to the strongly reduced intermolecular interactions in the adsorbate layer, PTCDA on Ag(110) is an ideally suited model system for quantum chemical studies with advanced methods including dispersion interactions; first, an investigation of the interaction between a single molecule and the substrate is sufficient, second, the adsorption site is known from previous experimental studies,12 and third, XSW gives a rather small distance of only 2.59 ( 0.01 Å between the carbon core and the topmost substrate layer, indicating a quite strong chemisorption.15 The adsorption of small molecules on the surface of transition metals has been analyzed in great detail with density functional theory (DFT), resulting in a quantitative understanding of the chemical bonding and of subsequent reactions catalyzed by the metal.16,17 The success of these investigations depends crucially on the large binding energy of the adsorbates studied, so that the absence of dispersion interactions in standard DFT does not induce significant shortcomings. For aromatic molecules adsorbed to transition metals, DFT can still reproduce the large chemisorption energies resulting from the interaction of the adsorbate with a partly filled d band.18–21

10.1021/jp902370b CCC: $40.75  2009 American Chemical Society Published on Web 10/22/2009

19898

J. Phys. Chem. C, Vol. 113, No. 46, 2009

In sharp contrast to this behavior, the binding potential of aromatic molecules on noble metals is rather shallow, so that the attractive van-der-Waals interaction between adsorbate and substrate may give a major contribution to the adsorption energy. For pentacene on gold, it was found that DFT with the gradientcorrected PBE functional understimates the adsorption energy by a factor of 4, and plane waves and localized basis sets give equivalent results once a calculation with the latter accounts for the counterpoise (CP) correction of the basis set superposition error (BSSE) according to the Boys-Bernardi scheme.22–24 In reference cases where dispersion interactions dominate, it was found that 25 common density functionals fail completely to predict both reasonable geometries and interaction energies,25 just like Hartree-Fock. Therefore, it was suggested that Hartree-Fock or DFT augmented by dispersion interactions included in the form of a sum over pair potentials may give a significant improvement with respect to standard DFT.26,27 Even though it cannot be denied that these semiempirical schemes reduce systematic deviations of DFT, they become questionable even for dimers of rare gas atoms where a pair potential for the dispersion interaction should be adequate.27 Nevertheless, in other cases like rare gas atoms physisorbed to noble metals, DFT augmented by a sum over a van-der-Waals pair potential can give quantitatively meaningful results.28 For polyaromatic molecules, it is evident that the calculated dispersion interaction has to be based on molecular orbitals, as opposed to a sum over pair potentials whose validity cannot be guaranteed. Second-order Møller-Plesset perturbation theory (MP2) contains the dispersion interaction arising from the repulsion between the Hartree-Fock ground state and double excited states.29 For rather small systems such as benzene dimers, the validity of MP2 can be checked with coupled cluster schemes like CCSD(T), revealing that MP2 overestimates the dispersion interaction by a small amount, whereas triple excitations included at the CCSD(T) level give a repulsive contribution to the dispersion interaction.30 Due to the different scaling with the size N of the system, that is, N3 for DFT, N5 for MP2, and N7 for CCSD(T), for large systems, the most advanced CCSD(T) scheme would require a prohibitive amount of computation time. In theses cases, MP2 emerges as the ideal compromise for an ab initio study of dispersion interactions, as demonstrated for the adsorption of rare gas atoms on copper31 and for benzene on noble metals.32,33 The most time-consuming step of a MP2 calculation consists of the evaluation of the correlation energy defined at the level of second-order perturbation theory, scaling as the fifth power of the system size. Within the turbomole 5.734 suite of programs used throughout the present work, this time-limiting step is sped up by expressing products of basis orbitals and products of molecular orbitals by a set of auxiliary atom-centered basis functions.35,36 This step resembles the insertion of a completeness relation, so that it is called resolution of identity (RI). Even though the required auxiliary basis set is much larger than the variational basis set used for the electronic orbitals, this RIMP2 scheme is about one order of magnitude faster than the respective MP2 calculation.36 Test calculations for DNA pairs have revealed that the difference between the total energies obtained with the MP2 and RI-MP2 computational schemes is very small.37 These investigations have demonstrated that the systematic deviations introduced by the RI scheme are orders of magnitude smaller than the residual errors introduced by the finite size of the variational basis set for the electronic orbitals. Hence, for systems which are too large to treat with MP2, it can be expected that RI-MP2 still yields reliable results.

Abbasi and Scholz

Figure 1. Intermolecular potential along the stacking direction of R-PTCDA, obtained from different ab initio approaches, with (b, 9) and without (O, 0) counterpoise correction;24 Hartree-Fock/SV(P) (magenta b and O, annoted HF), BLYP/SV(P) (green 9 and 0, B3LYP/ SV(P) (blue 9 and 0), and MP2/SV(P) (red b and O), and MP2/TZVP (orange b and O), as annoted. The stacking distance a ) 3.72 Å is indicated by a vertical line.

Methods We apply RI-MP2 calculations to the interaction of PTCDA with a Ag(110) substrate, and we use DFT calculations without any semiempirical correction for an analysis of the adsorbate geometry resulting in the absence of dispersion interactions. Before investigating chemisorption of this large polyaromatic molecule, we test the suitability of different ab initio methods on two simpler problems dominated by dispersion interactions, the intermolecular potential along the stacking direction in R-PTCDA and physisorption of benzene on (111)-oriented noble metals. Then, MP2 and DFT are applied to the optimization of a PTCDA molecule interacting with a (110)-oriented silver substrate, revealing substantial differences between the optimized geometries obtained with these computational schemes. A subsequent analysis of the resulting charge distribution in the adsorbate can be interpreted in more detail in terms of specific orbitals highlighting selected interaction mechanisms. The main achievements of the present work are summarized in the conclusion. Stacked PTCDA Pair In order to assess the suitability of different ab initio approaches for our purposes, we start with an analysis of the interaction between stacked PTCDA molecules in a geometry compatible with the crystalline R-phase.38 The substantial molecular quadrupole moments influence the geometric arrangement of adjacent molecules, so that the tilted stack geometry reflects a compromise between van-der-Waals attraction, electrostatic contributions, and overlap repulsion. Figure 1 shows the intermolecular potential along the stacking direction of PTCDA, calculated with the turbomole 5.734 software package and a split-valence basis including polarization functions, SV(P).36,39 From the comparison between Hartree-Fock, the gradient-corrected density functional BLYP, the hybrid functional B3LYP,40,41 and MP2,42 we find clear evidence that only in the latter case are the long-range correlations responsible for the van-der-Waals interaction included in an appropriate way, as discussed earlier for rare gas dimers and pairs of benzene molecules.30,43,44 For a pair of stacked benzene molecules, calculations without CP correction have resulted in intermolecular potentials resembling the respective curves in Figure 1, revealing a failure of Hartree-Fock and common density functionals for quantifying the van-der-Waals attraction, as opposed to MP2 including dispersion interactions.25

PTCDA on Ag(110) After applying the CP correction, a more flexible triple-ζ valence-polarized basis (TZVP)45 does not significantly modify the shape of the MP2 potential energy surface in Figure 1, so that a SV(P) basis turns out to be an adequate compromise for our purposes. For these intermediate-size basis sets, we find that the BSSE increases more rapidly than the net interaction after the counterpoise correction. This situation resembles earlier studies of benzene and naphthalene dimers with similar basis sets46,47 and transition-metal carbonyls treated with different contraction schemes for the orbital basis.48 Therefore, the quantitative validity of the counterpoise correction has been questioned,48 but independently of the basis size, it is generally accepted that the counterpoise correction is a reliable way of significantly reducing the BSSE error, corroborated by an analytic proof relying on full configuration interaction.49 In the MP2 calculation of the PTCDA dimer, the SV(P) basis is already large enough to produce a small overbinding, so that the minimum of the intermolecular potential occurs at a stacking distance of 3.66 Å, or about 2% below the lattice constant of a ) 3.72 Å in R-PTCDA. In the larger TZVP basis, the MP2 overbinding increases further, resulting in a potential minimum at 3.60 Å (compare Figure 1). Benzene on (111)-Oriented Noble Metals The adsorption of benzene on transition metals and noble metals is an interesting model system for understanding their suitability as catalysts for hydrogenation of aromatic compounds. In previous MP2 studies of benzene on Cu, Ag, and Au, it was found that the interaction is dominated by dispersion interactions, whereas charge donation or back-donation does not play a siginificant role.32,33 In the following, we calculate the potential energy surface between a planar benzene molecule and a rigid (111)-oriented metal cluster composed of 32 atoms, with 19 in the topmost layer, 12 in the second, and 1 in the third, corresponding exactly to the model geometries analyzed with MP2 before.32,33 In this seminal work, the 7 metal atoms directly underneath of the adsorbate and the 25 atoms in the surrounding part of the metal cluster were handled with different effective core potentials. This approach reduced the number of electrons per atom entering the Hartree-Fock scheme to 11 in the relevant region of the topmost layer and 1 in the embedding part of the cluster. The SV(P) basis set used in the following treats 19 electrons of each metal atom explicitly, whereas the electrons in the core (10 for Cu; 28 for Ag) are replaced by pseudopotentials.50,51 In a comparison between all-electron calculations and pseudopotential calculations for Ag and Ag+, it was proven that the energetic changes for the electronic configurations of interest (d9s2, d10s1, d9s1, d10) do not vary by more than 0.05 eV, indicating that the use of pseudopotentials has no significant influence on our results.50 For C6H6 adsorbed on Cu(111), Figure 2 shows the resulting MP2 potential energy surfaces together with a decomposition into the Hartree-Fock contribution and the long-range correlation energy. The CP correction of the BSSE shifts the potential minimum from 3.08 to 3.73 Å, where the latter value is close to a distance of 3.59 Å reported earlier.32 Our BSSE-corrected binding energy of -0.30 eV corresponds reasonably well to the previous value of -0.35 eV.32 This agreement is quite remarkable because the MP2 methods applied are not entirely equivalent. We treat a larger number of electrons per metal atom explicitly, but our orbital basis is less flexible. At the rather large adsorbate distance realized for benzene, we do not expect that the explicit inclusion of the filled sp shell

J. Phys. Chem. C, Vol. 113, No. 46, 2009 19899

Figure 2. Potential energy surfaces for benzene adsorbed on a 32 atom cluster of Cu(111), with (b) and without (O) counterpoise correction;24 Hartree-Fock/SV(P) (magenta), MP2/SV(P) (red), and the correlation energy included at the MP2 level (blue), as annoted. The positions of the potential minima with CP correction (3.73 Å) and without CP correction (3.08 Å) are indicated as vertical dashed lines.

below the d shell has a significant impact on the results. On the other hand, excitations from filled d bands to sp states above the Fermi energy dominate the optical response of silver.52,53 Hence, the explicit inclusion of d-symmetric basis states on substrate atoms at a larger distance from the adsorbate changes the resulting orbitals of the entire cluster and the lowest dipoleactive transitions, modifying in turn the contributions to the dispersion interaction included at the MP2 level. For the much lower height realized for adsorbed PTCDA, the increased overlap between the orbitals in the adsorbate and substrate will give a larger Pauli repulsion. As a result, the orthogonalization of the single-particle orbitals and the subsequent antisymmetrization of the many-electron wave function will push a larger amount of electronic charge density arising from the tails of the metal orbitals back toward the substrate. Therefore, in order to have a realistic shape of substrate orbitals compatible with the point group of the system composed of substrate and adsorbate, we think it is adequate to include d states on all metal atoms. Optimized Geometries of PTCDA Chemisorbed to Ag(110) Substrate Clusters. In the previous section, we have shown that the variational basis set and the effective core potential underlying our MP2 calculations are adequate for a quantitative investigation of weak dispersion interactions between an aromatic molecule and a noble metal substrate. For PTCDA, the 12 Ag atoms in the topmost layer dominate the chemisorption behavior, whereas the second Ag layer allows a redistribution of excess charges induced by the adsorbate, including favorable image charges compensating the charge distribution in the adsorbate. Due to the steep scaling of MP2 with system size, it is unavoidable that the substrate cluster has to be a compromise between computational requirements and a fair representation of the surroundings of the Ag atoms in the topmost layer. In the following, we use the two model clusters shown in Figure 3, Ag22 with 12 atoms in the first layer and 10 in the second and Ag32, where the number of atoms in the second layer is increased to 20. For the silver atoms interacting with the carboxylic groups of PTCDA, the Ag22 cluster realizes three nearest neighbors and the Ag32 cluster five, compared to seven for an extended (110) surface. When comparing the midgap energy (EHOMO + ELUMO)/2 of each silver cluster with the respective value for PTCDA, we find the molecular energy at a lower value. This misalignment

19900

J. Phys. Chem. C, Vol. 113, No. 46, 2009

Abbasi and Scholz

Figure 3. Optimized geometries of PTCDA on a (110)-oriented Ag22 cluster (left) and on a (110)-oriented Ag32 cluster (right). The geometries have been optimized with MP2/SV(P) or B3LYP/SV(P), as annoted. In each case, the distance from the substrate has been shifted rigidly to the minimum of the counterpoise-corrected potential surface (compare Table 1 and Figure 4 for the MP2/SV(P) calculation of PTCDA on Ag32).

TABLE 1: Calculated Height (Å) of Different Atoms above the Topmost ML of the Ag Substrate, with the Distance Shift Arising from the Counterpoise Correction of the BSSE Included in All Entriesa Substrate Ag32

Ag22 method

B3LYP

MP2

B3LYP

MP2

average C deviation C (calc. - exp.15) perylene core C central ring C end groups carboxylic O anhydride O average H BSSE distance shift E (eV), including BSSE E (eV), CP-corrected

3.15 0.56 3.23 3.33 2.73 2.49 2.57 3.30 0.04 -0.99 -0.27

2.81 0.22 2.83 2.84 2.72 2.63 2.68 2.98 0.20 -5.43 -2.12

3.05 0.46 3.13 3.24 2.66 2.37 2.53 3.15 0.04 -2.27 -1.41

2.68 0.09 2.69 2.62 2.64 2.50 2.63 2.86 0.24 -7.94 -3.74

a The last two lines give the binding energies, without and with counterpoise correction of the BSSE. All calculations have been performed using a SV(P) orbital basis, and the calculated average height of the carbon atoms is compared to the measured value of 2.59 ( 0.01 Å.15

in the region of the frontier orbitals favors transfer of electronic charge from the metal into the molecule, involving especially the former LUMO of the adsorbate; see the section on electronic orbitals for further details. Nevertheless, due to compensating mechanisms which will be discussed below, the net charge of chemisorbed PTCDA remains rather small. Relaxed Adsorbate Geometries. Figure 3 visualizes the relaxed geometry of an adsorbed PTCDA molecule obtained at the MP2/SV(P) and B3LYP/SV(P) levels.35,42,54 In all cases, the geometry of the metal cluster was kept fixed at a Ag-Ag bond length of 2.88 Å as in the bulk crystal. The heights of different atoms with respect to the topmost layer of silver are given in Table 1, including a rigid shift toward the potential minimum after applying the counterpoise correction (compare Figure 4 for the MP2 calculation of PTCDA on Ag32). The distance dependence of the relatively large BSSE found for the correlation part of the MP2 calculation modifies the slope of the BSSECP-corrected surface, so that the minimum of the potential surface shifts by 0.24 Å to larger distances. Due to the increased curvature of the MP2 potential surface with respect to the behavior of the PTCDA dimer in Figure 1 and that of physisorbed benzene in Figure 2, the distance shift arising from

Figure 4. Scan of the MP2/SV(P) adsorption energy of PTCDA on a (110)-oriented Ag32 cluster along the distance between the molecule and substrate, where the reference distance corresponds to the geometry optimized without counterpoise correction; b: with BSSE-CP correction;24 O: without BSSE-CP correction; lines without symbols: size of BSSE correction for the different energetic contributions (magenta: Hartree-Fock; blue: correlation at the MP2 level; red: MP2). The minimum of the counterpoise-corrected MP2 potential is indicated by a vertical line.

the BSSE-CP correction remains somewhat smaller. The BSSE of the Hartree-Fock potential depends only weakly on the distance, so that the shape of the potential surface after the BSSE-CP correction is hardly modified. Similarly, the shortrange correlations included at the B3LYP level show only a small BSSE depending weakly on the distance between the adsorbate and substrate, so that the distance shift of the minimum of the BSSE-CP-corrected potential surface remains as small as 0.04 Å (compare Table 1). B3LYP: Contributions to Chemisorption. Similar to the interaction between stack neighbors shown in Figure 1, the counterpoise correction of the BSSE produces a substantial reduction of the binding energy (compare Table 1). The B3LYP adsorption energies remain below the cohesive energy in R-PTCDA,14 falling instead into the range of previous DFT calculations for similar systems.9,13,55 If the small adsorption energy found in B3LYP were realistic, it should be possible to thermally desorb both multilayers and the last monolayer in a similar temperature range, but instead, it was found that the PTCDA monolayer cannot be thermally desorbed. Therefore, an exact monolayer coverage can be achieved by deposition of a multilayer and subsequent thermal desorption of all molecules which are not directly interacting with the metal substrate.4 On the Ag32 cluster, B3LYP/SV(P) gives a deviation of 0.46 Å between the calculated average height of the carbon atoms with respect to the measured value of 2.59 ( 0.01 Å,15 in the same range as the deviation observed for DFT calculations of PTCDA on Ag(111).5,8,9 As DFT does not include dispersion interactions, only the following four ingredients contribute to the binding between the adsorbate and substrate: (i) charge donation via directional bonds between oxygen atoms and the substrate, (ii) backdonation through a double occupancy of the former LUMO of PTCDA, (iii) a positive image charge in the region below the negatively charged molecule giving an attractive Coulomb contribution, and (iv) exchange repulsion between filled orbitals of the molecule and substrate. Due to the lack of a long-range nonlocal correlation term required for the van-der-Waals attraction, the overlap repulsion increases the distance between the perylene core and substrate, as in earlier DFT studies of similar systems.13,55 Nevertheless, the O-Ag bonds are so strong that they keep the height of the carbon core significantly below the sum of the van-der-Waals radii of silver and carbon, rvdW(C)

PTCDA on Ag(110) + rvdW(Ag) ) 1.70 Å + 1.72 Å ) 3.42 Å.56 This indicates a substantial strength of the chemical interaction, as opposed to physisorption, where the position of the adsorbate would essentially reproduce a value expected from the sum of the vander-Waals radii, a behavior observed for PTCDA on gold57 or for benzene on noble metals (compare Figure 2). MP2: Influence of Dispersion Interaction. In the MP2/ SV(P) calculation, the above phenomena are still included at the Hartree-Fock level, but in addition, the leading term of the dispersion interaction is taken into account. Therefore, exchange repulsion and van-der-Waals attraction compensate for each other to a large extent, resulting in nearly flat geometries and large binding energies (compare Figure 3 and Table 1). In fact, for the larger substrate cluster, the calculated average height of 2.68 Å of the carbon atoms is in excellent agreement with the measured position of 2.59 ( 0.01 Å.15 The residual deviation between calculated and measured height can be related to unavoidable shortcomings of our geometry optimization; obviously, the Ag32 cluster is still rather small, so that the precise geometry of the adsorbate may be influenced by boundary effects. Moreover, with the software used, the geometry could not be optimized at the counterpoise-corrected level. Instead, we have performed a counterpoise-corrected scan along the most important configuration coordinate, defined by the distance between the adsorbate and substrate, as visualized in Figure 4 for MP2 on the larger silver cluster. The MP2 geometry shows a bending of the C-H bonds away from the substrate, indicating a substantial repulsion between hydrogen and silver. As the distance between the H atoms and the topmost substrate plane is only 2.86 Å (compare Table 1), slightly below the sum of the van-der-Waals radii rvdW(H) + rvdW(Ag) ) 1.20 Å + 1.72 Å ) 2.92 Å, such a behavior has to be expected. In the B3LYP geometry, on the other hand, the H atoms remain at such a large distance that the geometric distortion with respect to the adjacent carbon atoms remains much smaller. Adsorption Energies. In B3LYP, the adsorption energy on the Ag32 substrate cluster increases by 1.14 eV with respect to Ag22, and in MP2, this difference is as large as 1.62 eV. The geometry of the strong O-Ag bonds is very similar on both substrate models, defining in turn a similar contribution to the adsorption energy. However, the larger substrate cluster allows a more flexible distribution of favorable image charges, with a large impact on the adsorption energy. In MP2, the chemisorption energies for both substrate models are much larger than the interaction between stack neighbors in R-PTCDA depicted in Figure 1, rationalizing the observation that the first ML of PTCDA on Ag(110) cannot be desorbed thermally without destroying the molecule. When comparing the two substrate clusters, the Hartree-Fock energy changes only from + 0.10 to -0.22 eV, but the MP2 correlation energy grows from -2.22 to -3.52 eV. The convergence of the dispersion interaction can be analyzed by a calculation of a sum over pair potentials scaling as the long-range part of the van-der-Waals interaction proportional to 1/r6. On each substrate cluster, we estimate the dispersion interaction for the reference position in the center of the central aromatic ring of PTCDA. On the Ag22 cluster, this gives 83% of the sum over 1/r6 pair potentials with respect to a silver substrate of infinite size. For the Ag32 cluster, this fraction increases to 87.5%, and due to the significant decrease of the height of the central ring on the larger substrate cluster, the lattice sum over 1/r6 is larger by a factor 1.31 when compared to the smaller substrate cluster. The substantially larger increase

J. Phys. Chem. C, Vol. 113, No. 46, 2009 19901

Figure 5. Partial charges in a free PTCDA molecule (black dashed), in the B3LYP geometry of the adsorbate visualized in Figure 3 (blue) and in the RIMP2 geometry (red). In each case, the partial charges are defined as Mulliken charges obtained at the B3LYP/SV(P) level, and the labeling of the atoms is according to the scheme displayed in the inset.

of the MP2 correlation by a factor of 1.59 indicates that the lower energetic cost for single excitations of the larger silver cluster plays a significant role, underlining the limitations of a simple sum over pair potentials. The large change of the adsorption energy between different substrate clusters demonstrates that a direct comparison may be problematic. However, in terms of the calculated height above the substrate, in MP2, the deviation between the relaxed geometries on both substrate models remains rather small, 0.13 Å for the average height of the carbon atoms and the carboxylic oxygens and 0.05 Å for the anhydride oxygens. The strong interaction of the four carboxylic oxygens with the substrate outweighs the attraction between the anhydride groups and the Ag surface, resulting in a height difference of 0.13 Å between the two kinds of O atoms on Ag32. To some extent, this difference in height can be rationalized from the Coulomb interaction between each oxygen atom and the silver atom underneath. Both kinds of oxygen atoms carry a similar negative net charge, but the silver atom below the anhydride oxygen has a particularly large negative charge, so that the respective O-Ag distance is increased by Coulomb repulsion. Charge Balance The above optimized geometries raise the question how the charge distribution in a chemisorbed molecule differs from the reference of a free molecule. In order to guarantee a comparable interpretation of different reference geometries, we have to choose the same computational scheme for all of them, for example, Hartree-Fock or B3LYP. As the latter hybrid functional is known to be quite reliable for various properties of organic systems, we shall apply it in the following. In Figure 5, we show the distribution of Mulliken charges obtained at the B3LYP/SV(P) level. In both adsorbate geometries, all atoms in the functional OdC-O-CdO end groups (atoms 1-5) increase their net charge, and the fact that no substantial difference between both adsorbate geometries occurs demonstrates that the Ag-O bonds formed in both cases are very similar. The central region of the perylene core (atoms 13-19) shows a small increase of electronic charges corresponding to a decrease of the positive charges close to the long axis, with a very similar behavior of both adsorbate geometries. In B3LYP, the charges on the hydrogen atoms (atoms 13 and 19) remain essentially unchanged with respect to the free molecule, but in the MP2 geometry, the energetic cost for the repulsion between H atoms and the substrate is reduced by a

19902

J. Phys. Chem. C, Vol. 113, No. 46, 2009

Abbasi and Scholz

Figure 6. Selected occupied orbitals in the MP2 model geometry, calculated with B3LYP/SV(P), ordered according to their energy. For each orbital, the representation in the point group C2V is indicated.

decrease of the electronic charges, reflected in an increased positive charge. In the periphery of the perylene core, the B3LYP geometry of the adsorbate follows the charge distribution of the free molecule, but in MP2, this region is subject to a charge redistribution around the long axis, without major impact on the total charge of the entire group. In MP2, the negative charge on the silver atom below the carbon atoms 9 and 16 is larger than that in the B3LYP calculation, indicating an increased influence of C-Ag interactions via hybrid orbitals arising from molecular orbitals mixing with substrate states. Electronic Orbitals In the MP2 geometry, the distance between the adsorbate and substrate is significantly below the sum of van-der-Waals radii, a clear indication of chemisorption. Therefore, in addition to the discussion of the global charge distribution in the previous section, it is interesting to establish the main interaction mechanisms via an analysis of the electronic orbitals involved. A few key examples are visualized in Figure 6, and their population analysis is summarized in Table 2. The most striking feature is the fact that the LUMO of the free molecule is now occupied by two electrons, resulting in STM pictures of the adsorbate reproducing the LUMO pattern.4 From this modified occupation of a previously empty molecular orbital, one would expect a large negative charge in the adsorbate. However, the atomic charges visualized in Figure 5

sum to a much smaller net charge of only about q ) -0.41 e on the entire molecule, indicating several compensating mechanisms transferring electronic charge into the substrate. Even though the visualization of the former LUMO in Figure 6 seems to indicate a localization exclusively on the adsorbate, the population analysis in Table 2 reveals a significant contribution of 16% from the metal substrate. Similar to the LUMO of the free molecule, its HOMO shows some hybridization with substrate orbitals, with an even smaller contribution of only 9% from the silver cluster. Nevertheless, the moderate hybridization of both orbitals reduces the negative net charge of the adsorbate by about half of an electronic charge. Due to the smaller point group of the adsorbate, that is, C2V instead of D2h for the free molecule, pairs of representations in D2h can be superimposed. Among the electronic states visualized in Figure 6, the b2 orbital at -11.44 eV has mainly σ character in the core region, hybridizing with π states on the end groups. Such a hybridization would be forbidden in D2h, but in the smaller point group of the adsorbate, it allows transfer of σ-symmetric electronic charge density from the core region into mainly π-symmetric contributions on the OdC-O-CdO functional groups. In the adsorbate, the slightly reduced electronic charge on the hydrogen atoms shown in Figure 5 can be understood quantitatively from two such σ-π hybrids transferring electronic charge from the core region toward the end groups. Besides the occupation of the former LUMO of PTCDA, the main interaction mechanism consists of the formation of bonding hybrid orbitals between the oxygen atoms and the substrate. In the free molecule, π orbitals on O atoms occur only in the form of conjugated states extending over a major part of the core region, but in the adsorbate, this conjugation between the core region and functional groups is modified, so that π orbitals localized mainly on O may occur. These states can then form bonding hybrids with the underlying substrate atoms, involving metal orbitals which would not contribute in the same way to occupied states of the substrate cluster alone. Such bonding orbitals based on previously occupied oxygen π states and on contributions from both previously empty and previously occupied orbitals in the substrate donate electronic charge from the adsorbate into the metal, so that they contribute to a reduction of the negative net charge in the adsorbate. As an example of such a bonding O-Ag hybrid, Figure 6 visualizes a b1 orbital at -9.63 eV. Of course, most of the π states still extend over the entire molecule, and in selected cases, the part in the core region of perylene can hybridize with orbitals in silver, as shown for a b2 state at -10.28 eV. Among the mechanisms defining the net charge of the molecule, this type of hybridization is the less efficient, mainly because most of these states occur as hybrid pairs, differing only in the signs of the wave function on the molecule and substrate. Therefore, when summing over both states in such a pair, a bonding and an antibonding one, only a small amount of electronic charge can be transferred between the adsorbate and substrate.

TABLE 2: Population Analysis of the Orbitals Calculated with B3LYP/SV(P) in the MP2 Geometry (compare their visualization in Figure 6) energy eV b2 LUMO(PTCDA) a2 HOMO(PTCDA) b1 π(O)-Ag hybrid b2 π(core)-Ag hybrid b2 σ-π hybrid

-4.97 -6.23 -9.63 -10.28 -11.44

C σ(s, p)

C π(p)

Cd

O σ(s, p)

O π(p)

0.02 0.02

0.01

0.12 0.15 0.11

0.02 0.51

0.69 0.72 0.01 0.10 0.17

0.03

0.15

Hσ

Ag

0.09

0.16 0.09 0.87 0.87 0.04

PTCDA on Ag(110)

J. Phys. Chem. C, Vol. 113, No. 46, 2009 19903

TABLE 3: Kohn-Sham Energies of the Orbitals in the MP2 Geometry Visualized in Figure 6 (in eV), of Similar Electronic States Occuring in the Bent B3LYP Model Geometry, And of the Respective States in a Free Molecule in D2h Symmetrya substrate -

Ag32 model geometry

MP2

B3LYP

B3LYP

b2 LUMO(PTCDA) a2 HOMO(PTCDA) b1 π(O)-Ag hybrid b2 π(core)-Ag hybrid b2 σ-π hybrid

-4.97 -6.23 -9.63 -10.28 -11.44

-4.65 -5.88 -9.62 -10.20 -11.18

-4.01 -6.57 -10.72 -11.77

a In all model geometries, the electronic states have been calculated with B3LYP/SV(P).

A comparison between the MP2 and B3LYP model geometries raises the questions which of the above interaction mechanisms can occur both in the flat and in the bent adsorbate geometries, and which are severely modified. In the bent B3LYP geometry, electronic states resembling the orbitals visualized in Figure 6 occur at similar energies (compare Table 3), but states localized mainly on the adsorbate have a significantly lower binding energy with respect to the flat MP2 geometry, a difference related to the larger net charge of q ) -0.49 e in the bent adsorbate. Even though the orbitals in both geometries have similar patterns, there may be significant differences. For the σ-π hybrid shown at the bottom of Figure 6, in the B3LYP geometry, the degree of hybridization in a similar orbital is three times smaller, so that this state is less efficient in transferring a small amount of electronic charge density from the H atoms and the perylene core toward the end groups. As a further σ-π hybrid shows a similar reduced hybridization, these minor differences between the orbitals in both model geometries have a major impact on the net charge of the H atoms shown in Figure 5. In the bent adsorbate geometry, this type of hybridization is so small that the electron density around the H atoms resembles the free molecule. Concerning π states hybridizing with substrate states via their core region, it is clear that the much larger distance in the bent geometry strongly reduces this specific hybridization mechanism, contributing in turn to a smaller transfer of electronic charge toward the substrate, so that the negative charge on the adsorbate remains somewhat larger. Due to the negative net charge of PTCDA in both adsorbate geometries, orbitals localized mainly on the molecule have a reduced binding energy with respect to a free molecule (compare Table 3), an effect which is more pronounced for the B3LYP geometry with a larger negative charge in the adsorbate. For the former LUMO, this difference between the two adsorbate geometries survives, but in both cases, the orbital energy moves below the Fermi energy of the substrate, that is, to a much larger binding energy with respect to the free molecule. In the MP2 model geometry, the redistribution of charge in the peryphery of the perylene core (atoms 8, 9, and 10 in Figure 5) seems to be related to a larger hybridization of the core region of molecular π states with the substrate. However, this region of the adsorbate is involved in so many σ and π states that a clear assignment of this modified charge distribution to a small number of electronic states could not be achieved, as opposed to the orbitals visualized in Figure 6, where a clear correspondence between orbitals, charge distribution, and interaction mechanisms can easily be established.

Conclusion In conclusion, the chemisorption of PTCDA on Ag(110) occurs via charge donation from the oxygen atoms into the metal substrate, back-donation from the metal into the former LUMO of the free molecule, and Coulomb attraction between the small negative net charge of the adsorbate and the induced positive image charge in the substrate. The occupation of the former LUMO with two electrons does not favor a specific geometric configuration with respect to the topmost substrate atoms, in sharp contrast to the formation of bonding orbitals between oxygen atoms and silver, playing a decisive role for the definition of the adsorption site. The compensation between van-der-Waals attraction and overlap repulsion is responsible for a nearly flat adsorption geometry. From previous MP2 studies of the interaction of rare gas atoms or benzene with noble metals, it is evident that this computational scheme is likely to provide quantitative insight into the influence of dispersion interactions on the geometry of the adsorbate. In the present work, a similar MP2-based scheme was applied for the first time to a much larger polyaromatic molecule, resulting in a deviation between calculated and measured heights below 0.1 Å. MP2 calculations like the one presented can serve as a benchmark for the test of improved density functionals accounting explicitly for nonlocal correlations required for the vander-Waals interaction.58,59 Once such tests for small substrate clusters would be convincing, the improved DFT schemes would eventually allow one to circumvent the unavoidable size limits arising from the N5 scaling of MP2, so that in the future, density functionals including nonlocal correlation might help to elucidate the impact of dispersion interactions on the chemisorption of large molecules on extended noble metal substrates. References and Notes (1) Hill, I. G.; Rajagopal, A.; Kahn, A.; Hu, Y. Appl. Phys. Lett. 1998, 73, 662–664. (2) Umbach, E.; Sokolowski, M.; Fink, R. Appl. Phys. A 1996, 63, 565–576. (3) Schmitz-Hu¨bsch, T.; Fritz, T.; Sellam, F.; Staub, R.; Leo, K. Phys. ReV. B 1997, 55, 7972–7976. (4) Glo¨ckler, K.; Seidel, C.; Soukopp, A.; Sokolowski, M.; Umbach, E.; Bo¨hringer, M; Berndt, R.; Schneider, W. D. Surf. Sci. 1998, 405, 1–20. (5) Hauschild, A.; Karki, K.; Cowie, B. C. C.; Rohlfing, M.; Tautz, F. S.; Sokolowski, M. Phys. ReV. Lett. 2005, 94, 036106. (6) Gerlach, A.; Sellner, S.; Schreiber, F.; Koch, N.; Zegenhagen, J. Phys. ReV. B 2007, 75, 045401. (7) Henze, S.; Bauer, O.; Lee, T.-L.; Sokolowski, M.; Tautz, F. Surf. Sci. 2007, 601, 1566–1573. (8) Rurali, R.; Lorente, N.; Ordejon, P. Phys. ReV. Lett. 2005, 95, 209601. (9) Hauschild, A.; Karki, K.; Cowie, B. C. C.; Rohlfing, M.; Tautz, F. S.; Sokolowski, M. Phys. ReV. Lett. 2005, 95, 209602. (10) Fink, R.; Gador, D.; Stahl, U.; Zou, Y.; Umbach, E. Phys. ReV. B 1999, 60, 2818–2826. (11) Seidel, C.; Awater, C.; Liu, X. D.; Ellerbrake, R.; Fuchs, H. Surf. Sci. 1997, 371, 123–130. (12) Bo¨hringer, M.; Schneider, W. D.; Glo¨ckler, K.; Umbach, E.; Berndt, R. Surf. Sci. 1998, 419, L95–L99. (13) Alkauskas, A.; Baratoff, A.; Bruder, C. Phys. ReV. B 2006, 73, 165408. (14) Scholz, R.; Kobitski, A.; Zahn, D.; Schreiber, M. Phys. ReV. B 2005, 72, 245208. (15) Bauer, O.; Hauschild, A.; Soubatch, S.; Henze, S. K. M.; Temirov, R.; Lassise, A.; Tautz, F. S.; Sokolowski, M. The adsorption of PTCDA on coin metal surfaces: a correlation between the adsorption height and the chemisorptive nature of the adsorbate-substrate bonding. DPG (VI)43/1, Berlin, Germany; 2008; Poster O.18.31. (16) Hammer, B.; Morikawa, Y.; Nørskov, J. K. Phys. ReV. Lett. 1996, 76, 2141–2144. (17) Alavi, A.; Hu, P.; Deutsch, T.; Silvestrelli, P. L.; Hutter, J. Phys. ReV. Lett. 1998, 80, 3650–3653.

19904

J. Phys. Chem. C, Vol. 113, No. 46, 2009

(18) Campbell, J. M.; Seimanides, S.; Campbell, C. T. J. Phys. Chem. 1989, 93, 815–826. (19) Xu, C.; Tsai, Y. L.; Koel, B. E. J. Phys. Chem. 1994, 98, 585– 593. (20) Saeys, M.; Reyniers, M. F.; Marin, G. B.; Neurock, M. J. Phys. Chem. B 2002, 106, 7489–7498. (21) Morin, C.; Simon, D.; Sautet, P. J. Phys. Chem. B 2004, 108, 5653– 5665. (22) Lee, K.; Yu, J.; Morikawa, Y. Phys. ReV. B 2007, 75, 045402. (23) Perdew, J. P.; Burke, K.; Ernzerhof, M. Phys. ReV. Lett. 1996, 77, 3865–3868. (24) Boys, S. F.; Bernardi, F. Mol. Phys. 1970, 19, 553–566. (25) Johnson, E. R.; Wolkow, R. A.; DiLabio, G. A. Chem. Phys. Lett. 2004, 394, 334–338. (26) Douketis, C.; Scoles, G.; Marchetti, S.; Zen, M.; Thakkar, A. J. J. Chem. Phys. 1982, 76, 3057–3063. (27) Wu, X.; Vargas, M. C.; Nayak, S.; Lotrich, V.; Scoles, G. J. Chem. Phys. 2001, 115, 8748–8757. (28) Lazic´, P.; Crljen, v.; Brako, R.; Gumhalter, B. Phys. ReV. B 2005, 245407. (29) Szabo, A.; Ostlund, N. S. Modern Quantum Chemistry: Introduction to AdVanced Electronic Structure Theory; Courier Dover Publications: Mineola, NY, 1996. (30) Sinnokrot, M. O.; Sherrill, C. D. J. Phys. Chem. A 2006, 110, 10656–10668. (31) Bagus, P. S.; Staemmler, V.; Wo¨ll, C. Phys. ReV. Lett. 2002, 89, 096104. (32) Bagus, P. S.; Hermann, K.; Wo¨ll, C. J. Chem. Phys. 2005, 123, 184109. (33) Caputo, R.; Prascher, B. P.; Staemmler, V.; Bagus, P. S.; Wo¨ll, C. J. Phys. Chem. A 2007, 111, 12778–12784. (34) Ahlrichs, R.; Ba¨r, M.; Ha¨ser, M.; Horn, H.; Ko¨lmel, C. Chem. Phys. Lett. 1989, 162, 165–169. (35) Weigend, F.; Ha¨ser, M. Theor. Chem. Acc. 1997, 97, 331–340. (36) Weigend, F.; Ha¨ser, M.; Patzelt, H.; Ahlrichs, R. Chem. Phys. Lett. 1998, 294, 143–152. (37) Jure cˇka, P.; Nachtigall, P.; Hobza, P. Phys. Chem. Chem. Phys. 2001, 3, 4578. (38) Lovinger, A. J.; Forrest, S. R.; Kaplan, M. L.; Schmidt, P. H.; Venkatesan, T. J. Appl. Phys. 1984, 55, 476.

Abbasi and Scholz (39) Scha¨fer, A.; Horn, H.; Ahlrichs, R. J. Chem. Phys. 1992, 97, 2571– 2577. (40) Becke, A. D. Phys. ReV. A 1988, 38, 3098–3100. (41) Lee, C.; Yang, W.; Parr, R. G. Phys. ReV. B 1988, 37, 785–789. (42) obtained with RI-MP2 scheme35,36 in TURBOMOLE 5.7, a development of University of Karlsruhe and Forschungszentrum Karlsruhe GmbH, 1989-2007, TURBOMOLE GmbH, since 2007; available from http://www.turbomole.com. (43) Kristyan, S.; Pulay, P. Chem. Phys. Lett. 1994, 229, 175–180. (44) Tsuzuki, S.; Lu¨thi, H. P. J. Chem. Phys. 2001, 114, 3949–3957. (45) Scha¨fer, A.; Huber, C.; Ahlrichs, R. J. Chem. Phys. 1994, 100, 5829–5835. (46) Sinnokrot, M.; Valeev, E.; Sherrill, C. J. Am. Chem. Soc. 2002, 124, 10887–10893. (47) Gonzalez, C.; Lim, E. C. J. Phys. Chem. A 2003, 107, 10105– 10110. (48) Persson, B. J.; Taylor, P. R. Theor. Chem. Acc. 2003, 110, 211– 217. (49) van Duijneveldt, F. B.; van Duijneveldt-van de Rijdt, J. G. C. M.; van Lenthe, J. H. Chem. ReV. 1994, 94, 1873–1885. (50) Andrae, D.; Ha¨ussermann, U.; Dolg, M.; Stoll, H.; Preuss, H. Theor. Chem. Acc. 1990, 123–141. (51) Eichkorn, K.; Weigend, F.; Treutler, O.; Ahlrichs, R. Theor. Chem. Acc. 1997, 97, 119–124. (52) Fehrenbach, G. M.; Bross, H. Phys. Status Solidi B 1996, 193, 231– 238. (53) Johnson, P. B.; Christy, R. W. Phys. ReV. B 1972, 6, 4370–4379. (54) The RI-MP2 optimization of PTCDA on Ag32 takes 60 days of CPU time on an AMD Opteron processor with 2.0 GHz and 8 GB of RAM. (55) Du, S. X.; Gao, H. J.; Seidel, C.; Tsetseris, L.; Ji, W.; Kopf, H.; Chi, L. F.; Fuchs, H.; Pennycook, S. J.; Pantelides, S. T. Phys. ReV. Lett. 2006, 97, 156105. (56) Bondi, A. J. Phys. Chem. 1964, 68, 441–451. (57) Duhm, S.; Gerlach, A.; Salzmann, I.; Bro¨ker, B.; Johnson, R.; Schreiber, F.; Koch, N. Org. Electron. 2008, 9, 111–118. (58) Andersson, Y.; Langreth, D. C.; Lundqvist, B. I. Phys. ReV. Lett. 1996, 76, 102–105. (59) Langreth, D. C.; Dion, M.; Rydberg, H.; Schroder, E.; Hyldgaard, P.; Lundqvist, B. I. Int. J. Quantum Chem. 2005, 101, 599.

JP902370B