Microscopic Characterization of the Interface between Aromatic

Apr 2, 2008 - John Kestell , Joshua Walker , Yun Bai , J. Anibal Boscoboinik , Michael Garvey , and ... Jennifer E. Laaser , Wei Xiong , and Martin T...
0 downloads 0 Views 391KB Size
Download PDF
J. Phys. Chem. C 2008, 112, 6413-6421

6413

Microscopic Characterization of the Interface between Aromatic Isocyanides and Au(111): A First-Principles Investigation Yan Li,*,† Deyu Lu,† Sally A. Swanson,‡ J. Campbell Scott,‡ and Giulia Galli† Department of Chemistry, UniVersity of California, DaVis, California 95616, and IBM Research DiVision, Almaden Research Center, San Jose, California 95120 ReceiVed: NoVember 21, 2007

Molecular assemblies comprising aromatic rings have received widespread attention as possible components of molecular electronic devices. An essential prerequisite to understand their stability and transport properties is the microscopic characterization of the interface formed with metallic leads. We present a comprehensive study, based on density functional theory (DFT) of the interface between the Au(111) surface and 1,4-phenylenediisocyanide, the simplest representative of aromatic isocyanides. We find that low coverage is favored at low temperature; this prediction is accessible to experimental validation as the metal work function shift upon adsorption has opposite signs, depending on the coverage. The computed binding energy between the isocyanide groups and the surface (∼0.5 eV) is smaller than that obtained for benzenethiols, but the presence of surface defects may considerably increase binding, by about 0.8 eV. The computed NtC stretching frequencies on atop sites show a blue shift, compared with the gas phase, consistent with experiments. In addition, our calculations show that energy differences between geometries with straight and those with tilted molecules on the surface are small compared with room temperature, which may explain the disordered patterns recently inferred from spectroscopy measurements. This suggests that longer isocyanide chains may be better suited to produce ordered self-assembled monolayers (SAMs) on Au(111). Finally, we discuss the electronic properties at the organic/metal interface by including self-energy corrections through many-body perturbation theory and surface polarization effects. Our results indicate that electronic structure calculations beyond DFT are required to make a correct, qualitative assessment of energy level alignment between SAMs and the metallic leads.

1. Introduction The search for molecular assemblies with promising transport properties for electronic devices has been an active field of research for the past two decades.1 The performance of such devices is strongly influenced by the electronic interactions, charge exchange and energy level alignment at the molecule/ electrode interface. Therefore, a microscopic characterization of this interface is an essential prerequisite to understanding and predicting molecular transport properties and to discerning promising candidates for molecular assemblies. One of the most extensively studied systems is composed of thiols on gold surfaces;2 this assembly is particularly attractive because of the strong chemical bond formed between the thiol radical (-S) and the surface gold atoms, with binding strengths ranging from 1 to 2 eV. However, the stability of this system may be hampered by oxidation, and this has motivated the search for alternative molecular assemblies. Among those, isocyanides have received widespread attention lately: the presence of a triple bond in the isocyanide group (-NC) may effectively connect pπ orbitals residing on the aromatic moiety and the dπ orbitals pointing out of the gold surface, thus providing a good network of delocalized electrons and thereby lowering the barrier for charge injection from the metal into the organic semiconductor.3 Several experimental efforts have been concentrated on determining the vibrational properties of isocyanide self* Corresponding author. E-mail: [email protected]. † University of California. ‡ Almaden Research Center.

assembledmonolayers (SAMs) usinginfrared(IR) spectroscopy.4-7 The shifts in the NtC stretching frequency upon SAM adsorption on metal surfaces may provide useful information, including the charge exchange mechanism between the organic molecule and the metal (e.g., whether donation-dominated or backdonation-dominated), binding site preference, and adsorption geometries. Temperature-dependent surface-enhanced Raman spectroscopy (SERS) and near-edge X-ray absorption fine structure spectroscopy (NEXAFS) have also been used to probe the binding and structure of aromatic isocyanide SAMs on Ag and Au surfaces;8 these studies suggested low surface coverage and orientational disorder at room temperature, as indicated, for example, by the negligible change in NEXAFS spectra at two different incident angles (90° and 20°). However, a definite assessment of the coverage and geometry of these SAMs is not yet available. Electrical conductivity measurements through aromatic isocyanide SAMs have been reported by various experimental groups.9-14 Nevertheless, no consensus has been reached in the literature about whether isocyanide groups serve as better molecular links to the gold electrodes and produce a higher conductance, compared with their thiol-terminated counterparts. For example, Chu et al. studied the conductance of alkanes, biphenyl, and diphenyl acetylene bridges with thiol and isocyanide terminals and found that the conductivity was enhanced by an order of magnitude when the thiol/gold linkage is replaced with an isocyanide/gold linkage,14 while the opposite trend was reported by other groups.13 Only few theoretical investigations of transport properties of isocyanide molecules on metals have so far accompanied

10.1021/jp7111044 CCC: $40.75 © 2008 American Chemical Society Published on Web 04/02/2008

6414 J. Phys. Chem. C, Vol. 112, No. 16, 2008 experimental studies, and these have been based on density functional theory (DFT) in the local density (LDA) or generalized gradient corrected (GGA) approximations, and on nonequilibrium Green functions calculations.15,16 One of the studies15 reported that the gold Fermi level (EF) falls between the highest occupied molecular orbital (HOMO) and the lowest unoccupied molecular orbital (LUMO) of phenylenediisocyanide, and is in near resonance with the delocalized LUMO state, leading to a large zero-bias conductance (45.8 µS), much higher than that of benzenedithiol (4.8 µS) calculated using the same method.17 On the other hand, another DFT study18 of phenylenediisocyanide molecule sandwiched between metal electrodes reported an energy difference between the LUMO and the Fermi level, E(LUMO) - EF, of 1.08 and 1.14 eV for Pt and Pd electrodes, respectively, and suggested that a bias of about 1 eV would be sufficient to significantly enhance the density of states (DOS) near EF. However, we note that approaches using DFT within LDA or GGA usually overestimate the conductance because of underestimated energy gaps between occupied and empty states in the molecule and possibly because of an incorrect alignment of energy levels at the molecule/metal interface.19,20 The binding of isocyanides to a gold surface has been studied within DFT, using both finite cluster models8,21,22 and periodic slabs;22,23 yet, the results of these calculations are not in agreement with each other. By studying the connection of phenylisocyanide molecule to one, two, and three gold atoms and using a B3PW91/LANL2DZ level of theory, Seminario et al. concluded that binding is only possible to a single Au atom.21 However, using a fixed 14-atom single-layer cluster, a LACVP** basis set for Au, and 6-31G** basis sets for other atoms, Kim et al. found that the same molecule binds more strongly at fcchollow site (1.0 eV) than at atop site (0.4 eV).8 These results disagree with those of Zhou et al. who studied the same system using a two-layer Au10 cluster at a B3lYP/6-31G**/LANL2DZ level of theory; these authors reported that adsorption is possible only on the atop site.22 The adsorption energy, however, changed dramatically, from 0.03 to 1.13 eV when the Au cluster was allowed to relax. In addition, when using a periodic three-layer slab model, Zhou et al. found that adsorption became possible at all sites, with adsorption energies of 0.87, 0.40, 0.37 eV at the atop, bridge, and hcp site, respectively.22 The minimal-size clusters and the type of basis sets used in these studies highlight the need for more accurate and systematic calculations, in order to resolve controversies about existing theoretical results. Finally, we mention a recent study carried out by Gilman et al., who studied the adsorption of CNH and CNCH3 on Au(111). A four-layer periodic slab model was employed, and the calculations were performed using DFT within GGA and a plane-wave basis set. The molecules were found to only adsorb on the atop site and at low coverage, because of strong dipoledipole repulsion. However, it is not straightforward to generalize these results to aromatic isocyanides, as the type and strength of intermolecular interactions within the SAM are expected to differ from those found in CNH and CNCH3. The existing theoretical and experimental controversies and the lack of consensus in the literature on the basic properties of isocyanide/ Au(111) interfaces calls for accurate and comprehensive theoretical investigations capable of providing a unified understanding of existing measurements. Here, we present a first-principles theoretical study of isocyanide on gold, aimed at a microscopic characterization of the interface and at providing a comprehensive understanding of its structural, vibrational, and electronic properties. Our results

Li et al. provide important information regarding the strength and weaknesses of these compounds as possible candidates for molecular electronic devices. In particular, we have studied in detail the adsorption of 1,4-phenylenediisocyanide (PDI) molecules on the Au(111) surface. These are the simplest yet representative molecules of the family of aromatic isocyanides. The effects of aromaticity, number of benzene rings, and different binding groups have also been investigated by comparing results for various adsorbate systems, including hydroisocyanic acid (CNH), 4,4′-biphenyldiisocyanide (BPDI), and 1,4benzenedithiol (BDT). We find that the isocyanide group binds to the gold surface with a moderate strength (the binding energy is ∼0.5 eV), and the most favorable structure for the SAM is at low coverage. Our computed vibrational properties are fully consistent with available experimental data. Although our calculations indicate disordered, low coverage pattern for PDI SAMs, there are possibilities of improving the packing density and inducing ordering through chemical modifications, which we briefly discuss in the paper. Finally, our description of the electronic states, using many-body perturbation theory in the GW approximation24 (G, Green’s function; W, screened Coulomb interaction), highlights the need for an accurate description of energy level alignments between organics and the metal, beyond the one provided by DFT in local or semilocal density approximations. The rest of the paper is organized as follows. In section 2, we describe the method and geometrical models used in our calculations. In section 3, we discuss our results for structural, vibrational, and electronic properties of PDI SAMs adsorbed on Au(111). Section 4 concludes the paper. 2. Computational Methods Plane-Wave Pseudopotential Approach. Structural, electronic, and vibrational properties have been obtained by performing plane-wave pseudopotential calculations within DFT in the Perdew-Burke-Ernzerhof (PBE) generalized-gradient approximation,25 using the PWSCF package.26 Norm-conserving pseudopotentials have been used for Au, H, C, N, and S, all generated with the PBE exchange and correlation functionals and including scalar relativistic effects for Au. For gold, we have generated pseudopotentials using the FHI98PP package27 according to the Troullier-Martin procedure,28 with cutoff radii 2.35 au (1.24 Å) for s and p orbitals and 1.35 au (0.71 Å) for d orbitals; the reference energy to pseudize the scatting p states was set at 10 eV, which was found to give an excellent transferability in the energy range of interest. Pseudopotentials for other elements were generated and extensively tested in previous studies.29 Wave functions were expanded in plane waves with a kinetic energy cutoff Ecut ) 50 Ry, which was found to yield well-converged results for most of the properties studied here. All calculations are spin-unpolarized, except those for insulating systems with an odd number of electrons, for example, the AuPDI complex or isolated thiol radicals. Slab Model of Au(111) Surface. Before carrying out calculations for isocyanides absorbed on Au(111), we computed properties of the bulk fcc gold lattice, using a 10 × 10 × 10 Monkhorst-Park (MP) k-point grid30 with Methfessel-Paxton smearing occupation (Gaussian spreading of 0.01 Ry) of the Bloch states. The equilibrium lattice constant, a0, and bulk modulus, B, were calculated to be 4.16 Å and 1439 kbar (experimental values are a0 ∼ 4.06 Å and B ∼ 1720 kbar), respectively, in good agreement with all-electron FPLAPW results;31 our findings are consistent with the general tendency of gradient-corrected DFT to slightly overestimate the lattice

Interface between Aromatic Isocyanides and Au(111) constants and underestimate the bulk moduli of solids. Nevertheless, GGA-PBE is expected to accurately account for the properties of systems in which homogeneous and inhomogeneous electron densities coexist, compared with LDA.32 The Au(111) surface was modeled with a repeated slab of four atomic layers with calculated equilibrium lattice constant of 4.16 Å, separated by nine layers of vacuum (∼22 Å). For SAMs at full coverage (Θ ) 1), defined as one adsorbate molecule per three surface gold atoms, the (x3 × x3)R30° unit cell was used, and the Brillouin zone was sampled with a 6 × 6 × 1 MP grid of k-points. For the low coverage case (Θ ) 1/3), the (3 × 3)R30° unit cell was used and sampled with a 3 × 3 × 1 MP grid. During structure optimization, the bottom Au layer was kept fixed while the top three Au layers and the adsorbed molecules were fully relaxed until the minimum force was smaller than 0.026 eV/Å. The calculated surface energy for the clean Au(111) is 767 erg/cm2, in agreement with previous DFT calculations.33 This value is smaller than the reported experimental values, in the range of 1045-1410 erg/cm2.34 The underestimate of the surface energy is a well-known shortcoming of gradient-corrected DFT approaches.35,36 Because of the periodicity of the supercell in the z direction imposed by the slab geometry, an artificial electric field is induced by the two-dimensional dipole layer formed at the molecule/gold interface.37 We have verified that, for all systems studied here, the vacuum space is large enough to render dipole corrections to the adsorption geometry and energetics negligible. In the following, dipole corrections were only applied after geometry relaxation, to yield quantitative results for the dipole moment and work function shift. Phonon Calculation. Vibrational frequencies of gas-phase molecules and the adsorbate systems were calculated using density functional perturbation theory (DFPT).38 For gas-phase molecules, the full vibrational spectra were computed and compared with the partial vibrational spectra where only the C and N atoms in the isocyanide groups were displaced from their equilibrium positions. The difference in the NtC vibrational frequency obtained in the two ways is marginal (about a few cm-1 or 0.1%). For computational efficiency, we calculated only the partial phonon spectra for isocyanides on Au(111), while fixing the equilibrium positions of the surface atoms and of those atoms not belonging to the isocyanide groups. Since gold has quite large atomic mass compared with that of first row atoms, the resulting stretching frequencies are expected to have small deviations from those obtained from the full phonon calculation. GW Self-Energy Corrections. It is well-known that DFT (LDA or GGA) calculations tend to underestimate the electron removal energy and to overestimate the electron addition energy of small molecules, which is approximated by the energy levels of HOMO and LUMO, respectively. An accurate way of correcting these energies is to compute self-energy corrections to the Kohn-Sham eigenvalues through many-body perturbation theory, using the GW approximation.24 When adsorbed on a metal surface, the HOMO and LUMO energies are further modified by polarization from the metal, which is not fully captured at the LDA or GGA level, because of the incorrect asymptotic tail of the surface exchange-correlation potential.39 We computed the energy alignment of molecular levels relative to the gold Fermi level by employing GW calculations for the free-standing SAMs, and we added surface polarization effects by using a classical image potential model. The GW calculations were performed using a parallel version of the ABINIT code40 modified by the authors; computational details can be found in Supporting Information.

J. Phys. Chem. C, Vol. 112, No. 16, 2008 6415 TABLE 1: Calculated Energetic and Structural Parameters for Molecular SAM Structures Adsorbed on the Au(111) Surfacea coverage

Emol-Au (eV)

Emol-mol (eV)

Ead (eV)

|EvdW| (eV)

dAu-R (Å)

1 1 (Herringbone) 1/3

PDI/Au(111), atop 0.44 -0.38 0.06 0.46 -0.22 0.24 0.59 -0.05 0.54

0.29 0.24 0.01

2.04 2.04 2.03

1 1 (Herringbone) 1/3

BDT/Au(111), fcc 1.13 -0.27 0.86 1.12 -0.09 1.03 1.15 -0.03 1.12

0.27 0.25 0.01

2.55 2.55 2.45

1 1/3

NH3/Au(111), atop 0.18 -0.05 0.13 0.37 -0.02 0.35

0.01 0.0

2.57 2.41

a Contributions to the adsorption energy (Ead) include interactions between the molecules and the Au(111) surface (Emol-Au), and interactions within the SAM (Emol-mol), both calculated within density functional theory. Empirical vdW interactions (EvdW) are listed for comparison. dAu-R is the shortest gold-adsorbate distance (R ) C, S, and N for PDI, BDT radical and NH3, respectively). The molecular axes were set normal to the surface in the initial configurations, and the optimized structures may not reflect the most stable geometries, especially at low coverage.

3. Results and Discussion 3.1. Adsorption Geometries and Energetics. The adsorption energy per molecule for PDI on Au(111) is defined as

Ead ) -(EPDI/Au(111) - EAu(111) - EPDI)

(1)

where EPDI/Au(111), EAu(111), and EPDI are the total energies of a PDI SAM adsorbed on Au(111), of a clean Au(111) surface and of a gas-phase PDI molecule, respectively. According to this definition, an exothermic adsorption corresponds to a positive value of Ead. Since experimental evidence from IR measurements seems to indicate that PDI molecules bind to Au(111) surface with a terminal (η1) geometry,4-7 we first studied the case of PDI molecules vertically placed on the atop site. The aromatic rings of the PDI molecules were arranged to be parallel to the next-nearest-neighbor (NNN) direction so as to minimize steric repulsion, which, at high coverage, was estimated to amount to as much as 5 eV per molecule if the rings are parallel to the nearest-neighbor (NN) direction. The system was then fully relaxed, and the calculated bond lengths and adsorption energies are listed in Table 1. At full coverage (Θ ) 1), it was found that a herringbone arrangement of the PDI SAM with a (2x3 × x3) unit cell is energetically more favorable than the (x3 × x3) structure (see Figure 1). By decomposing the adsorption energy Ead into contributions from molecule-surface interactions (Emol-Au) and molecule-molecule interactions (Emol-mol), one finds that, although the first term is always positive, the significant intermolecular repulsion calculated at the DFT level almost cancels the molecule-surface attraction for the (x3 × x3) structure. Most of the repulsion comes from quadrupole-quadrupole interactions as the PDI molecule has a rather large quadrupole moment along the molecular axis, |Qzz| ≈ 23 Debye‚Å, contributed by the isocyanide groups at both ends. This is much larger than |Qzz| of the benzene molecule, which was calculated to be less than 3 Debye‚Å. As a result, Emol-mol for a benzene (x3 × x3) SAM is only of the order of 10 meV. By rearranging the orientations of the PDI molecules into a herringbone structure, the intermolecular repulsion is effectively

6416 J. Phys. Chem. C, Vol. 112, No. 16, 2008

Li et al. TABLE 2: Density Functional Theory Results for Isolated PDI and AuPDI Molecules, and PDI Adsorbed on Au(111) Surface at Low Coverage (Θ ) 1/3), Including the Binding Energy (Emol-Au), the Shortest Gold-Adsorbate Distance (dAu-C), Bond Lengths of the Isocyanide Groups (dNtC), and the NtC Stretching Frequency (νCtN) νC≡N (cm-1)

dC≡N (Å) system PDI AuPDI

Emol-Au dAu-C (eV) (Å) unbound bound 0.99

2.01

0.59 0.71 0.71

2.03 2.15 2.21

1.20 1.20

1.19

unbound 2040 2036

PDI/Au(111), Θ ) 1/3 1.20 1.19 2039 1.20 1.21 2026 1.20 1.22 2033 2120-2122

bound 2102

Figure 1. Top view of the optimized atomic configurations of (x3 × x3) structure (left) and (2x3 × x3) or herringbone structure (right) at full coverage.

atop bridge fcc expt4-7

reduced, while the binding strength between the SAM and the gold surface remains almost the same. In order to assess the accuracy of computed adsorption energies, one needs to take into account that van der Waals (vdW) forces may contribute to the interaction energy between aromatic molecules and these are not treated accurately within the GGA approximation adopted here. In particular, at high coverage, this inaccurate treatment may result in an underestimate of the adsorption energy and surface stability. We used an empirical dispersion model41 to estimate long-range dispersion corrections to the total energy (details in Supporting Information). By adding this contribution, denoted as EvdW, to the total adsorption energy Ead determined within DFT, the intermolecular interaction is much reduced at high coverage, but the herringbone geometry remains more stable than the (x3 × x3) structure (see Table 1). For comparison, we also considered adsorption geometry and energetics of BDT SAMs at the same intermolecular spacing, which were found to be in excellent agreement with the results of previous DFT calculations on benzenethiol.42 Overall, the binding of thiol molecules to gold was found to be stronger than that of the isocyanide molecules. Nevertheless, the calculated binding energy of PDI is more than sufficient to form a stable monolayer at room temperature. In our calculations, the most favorable geometry was found at low coverage with a (3 × 3) unit cell, where the strength of molecule-molecule interaction is reduced to the order of a few tens of millielectronvolts, while the molecules are drawn closer to the surface, resulting in a strong binding to the gold atoms. From Table 1, one finds that the binding strength of the PDI molecule with Au(111) surface is intermediate between that of the BDT radical (with the hydrogen atom desorbed from one end and forming a strong covalent bond to gold through the thiol group), and that of the NH3 molecule, which is mainly bound by dispersion forces.43 We note that the presence of surface defects may considerably increase the binding strength. For instance, the adsorption energy of the PDI molecule on top of an Au adatom was calculated to be 1.3 eV, about 0.8 eV higher than that at atop site on the flat surface. Table 2 lists the DFT results for PDI adsorbed on different sites of the Au(111) surface. Binding at bridge and fcc sites was found to be energetically more favorable than binding at atop site, contrary to the conclusions of experiments reported in refs 4-7. Similar discrepancies were found in the literature for CO adsorption on noble and transition metals (this is sometimes referred to as the CO/Pt(111) puzzle44), with most DFT calculations predicting the wrong adsorption sites. This is likely due to the exchangecorrelation functionals used at the LDA and GGA levels, which tend to overestimate π* backdonation contributions and therefore

favor adsorption sites with higher coordination numbers. One solution may be the introduction of an orbital-dependent term U to the GGA energy functional, that is, the use of the socalled GGA+ molecular U approach.45,46 An alternative one might be to employ hybrid functionals. For example, a recent study of CO was performed using hybrid DFT implemented using periodic boundary conditions and plane-wave basis,47 which was shown to capture the key aspects of the CO/Pt(111) interaction, including the identification of the adsorption site and the CO frequency shift upon adsorption. By replacing PDI molecules with BPDI (which contains two aromatic rings), one finds almost identical geometrical parameters (e.g., dAu-C and Emol-Au), showing that the number of aromatic rings does not have a significant influence on the binding strength between the isocyanide group and the gold surface. The vdW contribution, however, becomes much stronger at high coverage with multiple aromatic rings (|EvdW| ∼ 0.5 eV for BPDI in the herringbone structure), making it favorable to form closely packed SAMs. Dhirani et al. investigated the structures of monolayers of oligo(phenylethynyl) benzenethiols on Au(111) and revealed that the degree of order in these systems increases with the number of aromatic rings.48 Similar tendency may also be expected for aromatic isocyanides, but so far no experimental evidence for well-ordered BPDI or TPDI SAMs has been reported. Although definite, robust conclusions on the preferred coverage cannot be drawn based only on energetics, we will see below that a comparison of our calculation with measured vibrational spectra provide compelling evidence in favor of low coverage being a more stable geometric configuration for PDI SAMs. 3.2. Vibrational Frequencies of NtC Stretching Mode. Several recent experimental studies4-7 were reported on the adsorption of aromatic isocyanide SAMs on noble or transition metals. The frequency of the metal-coordinated isocyanide stretching mode was measured by infrared spectroscopy, and it was concluded that the isocyanide molecules bind to the Au(111) surface with an end-on or terminal (η1) geometry, indicated by the appearance of an additional vibrational peak about 50-60 cm-1 higher than that of the bulk sample. Although the relative binding strength at different adsorption sites was not correctly predicted at the GGA level (Table 2), it was recently demonstrated that the vibrational frequencies can be predicted fairly accurately using GGA, as they are not controlled directly by the π* backdonation strength.46 Since the isocyanide group is isoelectric with CO, the bonding configurations studied here are similar to those of CO. Therefore, we expect a DFT/ GGA description of the vibrational frequency to be accurate. As shown in Table 2, the NtC stretching frequency of the bound isocyanide group decreases with increasing coordination

2099 1923 1847 2170-2181

Interface between Aromatic Isocyanides and Au(111) number, while that of the unbound isocyanide group is almost independent of the adsorption sites. These results are in good agreement with physical intuition and experimental findings,4-7 supporting the conclusion that the appearance of a higher vibrational frequency peak corresponds to the NtC stretching mode of the metal-coordinated isocyanide group, bound to the Au(111) surface with an end-on coordination, that is, on the atop site. We also computed νNC for the high coverage cases, which show a red-shift of 50-100 cm-1 compared with the gas phase. The decrease of νNC is most likely caused by the strong intermolecular interactions at close distance, which weaken the stretching of the NtC bond, and by the decrease of σ donation to metal because of competing effects from neighboring molecules. Since no red-shift was observed in the measured IR spectra, this finding supports our previous conclusion that low coverage is favored over high coverage. We note that a different tendency was observed for CO adsorbed on metals, where the CO stretching frequency actually increases with increasing CO coverage; this was attributed to dynamic dipole-dipole coupling49 or weakening of adsorbate-substrate bond.50 Although CO is isoelectric to the isocyanide group, the dominant mechanisms for both the intermolecular interaction and the charge exchange at the adsorbate/metal interface are different in the two cases, and therefore it is not entirely surprising that the shifts in frequency as a function of coverage are not the same for the two groups. From Table 2, one also sees that the DFT determined vibrational frequencies are systematically lower than the experimental values. To test the convergence of νNC with respect to the kinetic energy cutoff Ecut (see Computational Methods), we used the isolated AuPDI molecule as our model system. When Ecut is increased from 50 Ry up to 100 Ry, the bond lengths of both isocyanide groups decrease only by ∼0.01 Å, and Emol-Au differs only by less than 10 meV. However, νNC of the unbound (bound) isocyanide group increases from 2036 cm-1 (2102 cm-1) to 2106 cm-1 (2165 cm-1) and becomes much closer to the values found in IR spectra of PDI on Au(111). This tendency is consistent with the fact that convergence of vibrational frequencies usually requires a higher Ecut than that of geometrical properties. All of the phonon calculation results presented here are not scaled to account for cutoff convergence, although in the literature a scaling factor is sometimes used to fit the theoretical results to the experimentally observed vibrational frequencies.7 3.3. Tilting at Low Coverage. In the calculations discussed so far, PDI molecules were initially placed vertically on the Au(111) surface along the NNN direction. Because of mirror symmetry, the relaxed geometry remains vertical. It is interesting to see how large are the energy differences between configurations with the molecules tilted away from the normal to the surface. We therefore repeated geometrical optimization by starting with the lower NtC bond of the PDI molecules tilted at 30° and 60° while keeping the molecular geometry unperturbed. The initial and final tilted angle, denoted as θ0 and θ respectively, are given by the angle between the lower NtC bond and the surface normal in the corresponding geometries; both in-plane and out-of-plane (relative to the aromatic ring) configurations were studied. All calculations were carried out at low coverage (Θ ) 1/3). The results are shown in Table 3, and the optimized structures are plotted in Figure 2. The bond length of the bound isocyanide group is very similar after geometry relaxation for the five different tilted configurations, while the Au-C bond length varies in the range 2.032.09 Å. The adsorption energy differs by less than 50 meV,

J. Phys. Chem. C, Vol. 112, No. 16, 2008 6417 TABLE 3: Calculated Structural, Energetic, and Vibrational Parameters for Tilted PDI Molecules on the Au(111) Surfacea νN≡C (cm-1) θ0 (°) θ (°) ∆Ead (eV) dAu-C (Å) unbound bound in-plane out-of-plane

0 30 60 30 60

0 16 32 21 41

[≡0] 0.04 0.02 ∼0 -0.04

2.03 2.04 2.08 2.05 2.09

2039 2034 2043 2034 2041

2099 2067 2038 2077 2027

a θ0 and θ are the initial and final tilt angles between the bound NtC bond and the surface normal. The variation of dNtC is negligible and therefore not listed.

which implies that, at room temperature, these configurations are energetically indistinguishable and may coexist at low coverage. This may make it difficult to form well-organized structures at low coverage as the molecules can take many different orientations without suffering serious intermolecular repulsion. This observation is in agreement with the disordered patterns recently suggested by NEXAFS spectroscopy.8 The height of these molecules, measured from top of the clean Au(111) surface, varies from 5.8 to 10.0 Å, which may also make it difficult to form molecular arrays with uniform height for use in transport devices. However, with proper chemical modification or functionalization, it may still be possible to force these molecules to arrange vertically on the surface. For example, by increasing the number of aromatic rings, one can effectively enhance the π-π attraction, and this may drive the adsorbed molecules into a well-ordered geometry, in order to maximize the energy gain. Steric hindrance, for example by m-methyl substitution, could also favor vertical alignment and promote ordering within the layer. New classes of aryl isocyanides may also be worth exploring. For example, isocyanides and diisocyanides based on the nonbenzenoid azulenic framework have been synthesized and investigated recently, which were shown to be remarkably air- and thermally stable and may provide the possibility of forming SAMs of controlled orientation and packing.51,52 One interesting phenomenon is that νNC(bound) decreases as the tilt angle increases and may even come close to that of the unbound isocyanide group. These findings are consistent with the existence of two peaks found in experimental vibrational spectra, one centered at the free C-N vibrational frequency and the other at higher frequencies, and we speculate that the lowfrequency peak may include contributions from bound isocyanide groups with large tilt angles. Experimentally, it was found that the PDI vibrational spectrum shows poor resolution between the two peaks whereas the BPDI spectrum goes down to the baseline between the peaks.6 The latter is broadened asymmetrically, sharp on the high-energy side and broadened toward low energy. This may be indicative of higher order for BPDI (because of the stronger aromatic interactions), although with some remaining tendency to tilt. In addition, the free and bound isocyanide bands are better resolved in tetramethyl-PDI than in the parent molecule, again suggesting a higher degree of order in the sterically constrained molecule. It would be interesting to see whether the relative intensities of these two peaks vary as a function of temperature, thus indicating a tendency toward ordering when the sample is cooled down. 3.4. Charge Transfer and Energy Level Alignment. The charge transfer between molecules and metal surfaces is sometimes used to quantify interfacial interactions. Changes of atomic orbital populations, for example, evaluated from Lowdin charges, may give a qualitative idea of this transfer. The loss

6418 J. Phys. Chem. C, Vol. 112, No. 16, 2008

Li et al.

Figure 2. Optimized structures of PDI molecules on Au(111) with initial tilt angle of θ0 ) 0°, 30° (in-plane), 60° (in-plane), 30° (out-of-plane), and 60° (out-of-plane).

Figure 3. Plane-integrated electron density difference ∆F (a) and planeaveraged electrostatic potential difference (b) before and after the PDI molecule adsorbs on Au(111) surface. z ) 0 corresponds to the top layer of Au(111) surface. Also shown is the isosurface of charge difference for the (3 × 3) structure (red indicates electron gain and blue indicates electron loss).

of electron populations in the σ orbitals and the gain in the π orbitals of the PDI molecule found in our calculations are consistent with the Blyholder model.53 On the other hand, the overall charge transfer between the whole PDI molecule and the Au(111) surface is less than 0.1 e, which is within the error of the atomic projection technique. Instead of examining projected atomic charges, which are only accurate from a qualitative point of view, it may be more useful to look at the spatial distribution of the induced charge density, shown in Figure 3. Clearly, the induced charge density decays from the molecule-surface interface in an oscillating fashion. Because of the different strength of surface interaction at different coverages, the induced axial dipole moments turn out to be of opposite signs: µ ) -0.37 Debye for the full coverage herringbone structure and µ ) 0.85 Debye for Θ ) 1/3, respectively. Consequently, the work function of the gold surface is increased by 0.63 eV for the former and decreased by 0.47 eV for the latter case (relative to the vacuum level on the SAM side). This difference in work function is accessible to experimental verification, for example, through photoelectron spec-

Figure 4. Projected density of states for PDI adsorbed on Au(111) at Θ ) 1/3, onto the surface Au atom before and after binding to the PDI molecule (a) and onto the bound (b) and unbound (c) isocyanide groups. The Fermi level is located at 0 eV. Representative hybridized orbitals are marked at the corresponding PDOS peaks (see text) and plotted in Figure 5.

troscopy or scanning Kelvin probe, and may be used to probe the predictions about coverage provided by our calculations. In the following, the bonding nature of the isocyanide group to Au(111) is discussed in detail by analyzing the projected density of states (PDOS). Plotted in Figure 4a is the PDOS of the surface Au atom binding to the PDI molecule, where one finds many extra features compared with that of the free surface. In a similar fashion, the splitting and shift of PDOS peaks and delocalization of electron densities near the Fermi energy for the bound isocyanide group in Figure 4b are signatures of chemical bonding with the Au atom. We first discuss σ interactions. The σ state of the isocyanide group (containing 2s, 2pz orbitals of C and N atoms), originally at -3 eV below the Fermi level, interacts strongly with the gold 5dz2 and 6s

Interface between Aromatic Isocyanides and Au(111)

J. Phys. Chem. C, Vol. 112, No. 16, 2008 6419

Figure 5. Isosurfaces of square moduli for the hybridized d˜ πy, π˜ y, d˜ σ, and π˜ /y orbitals of PDI adsorbed on Au(111), as indicated in Figure 4b. Red and blue colors correspond to opposite phases of the orbital wave functions. Note that d˜ πy orbital contains mixed contributions from πy and π/y orbitals of the bound isocyanide group in such a ways that the wave function is vanishing at the bottom C atom.

orbitals. As a result, a new peak appears at about -8.5 eV, corresponding to a localized σ˜ state, indicating an adsorbate state of σ bonding with the surface metal atoms. At the same time, resonance states with antibonding character are formed with energies close to EF. These states are denoted as d˜ σ, since most contributions come from the metal d band,54 as indicated by the isosurface plot in Figure 5. Since the square moduli of these states are confined to the interface and do not interact with orbitals in the benzene ring, they are not expected to contribute to the electron transport along the molecule. We now turn to the description of π interactions. These are more complex, as they contain contributions from both in-plane (πx) and out-of-plane (πy) orbitals of the isocyanide group as well as from the benzene ring. In addition to forming bonding and antibonding states with gold atoms individually, the π and π* states of the isocyanide group also hybridize when binding to the gold atom. For example, the d˜ πy state contains both πy and π/y state contributions, which are mixed in such a way that a nonbonding state is formed between C and Au atoms, with an almost zero amplitude of the wave function on the C atom (see Figure 5). Because of the involvement of π orbitals on the benzene ring, the charge density is effectively delocalized over the whole π network of the molecule. A recent study of CNH adsorption on Au(111) surface yielded results similar to ours for the electronic structure,23 although that of PDI molecules are more complex because of the hybridization of isocyanide orbitals with those of the benzene ring. In the field of molecular transport, the barriers for electron and hole transport are usually quantified by |EF - EH| and |EF - EL|, which measure the energy difference between HOMO and LUMO of the adsorbed molecules and the metal Fermi level. The peaks of πy and π/y states on the PDOS plot in Figure 4 thus correspond to a 0.7 eV barrier for electron conduction and a 2.9 eV barrier for hole conduction. However, the Kohn-Sham energy eigenvalues of DFT (LDA or GGA) calculations may severely underestimate the electron removal energy while the opposite trend is expected for the electron addition energy, yielding a HOMO-LUMO energy gap smaller than experiment. This is mainly due to the overestimate of screening within the molecules, which is much weaker compared with that of the free electron gas, on which the derivation of most exchangecorrelation functionals is based. For example, our DFT calculaDFT DFT tion predicted EHOMO ) -6.9 eV and ELUMO ) -3.2 eV for a

gas-phase PDI molecule. A more accurate estimate from total energy differences of the neutral molecule and corresponding cation/anion yield EHOMO ) -9.5 eV and ELUMO ) -0.9 eV instead. (This was obtained using Gaussian03,55 spin-restricted PBE with aug-cc-pvTZ basis set.) Errors introduced by the GGA approximation are expected to be smaller for molecules adsorbed on a metal surface because of screening from the metal, which stabilizes the added electron or hole and reduces the energy gap of the isolated molecule. However, to properly account for the screening effects from both the substrate and the neighboring molecules in the SAM, one needs an accurate determination of energy alignment beyond LDA/GGA. A more accurate way to compute molecular energy levels is through the use of many-body perturbation theory and in particular the GW approximation56 (see computational details in Supporting Information). We employed this approximation and split the problem into two parts: self-energy corrections for a free-standing PDI SAM were first computed through a GW calculation; the surface polarization effect was then included through a classical image potential model. Both terms were added to the GGA-determined eigen energies for the adsorbate system, with two assumptions: (1) the energy shift due to an induced surface dipole is already well described at the GGA level; and (2) for the electronic states of interests, the wave functions are close enough to those of the free-standing SAM so that surface polarization can be treated as an external perturbation. Since the GGA calculation already yielded a rather large lower bound for |EF - EH| (∼3 eV), here, we will focus on the unoccupied π/y states. For a π/y state, an overlap of ∼90% is found with LUMO of the SAM. The energies of the adsorbate levels can then be approximately written as:

Ead(j) ≈ Ead,GGA(j) + [ESAM,GW(j) - ESAM,GGA(j)] + ∆Psurface(j) (2) where the terms on the right-hand side are the Kohn-Sham (GGA) energy for an adsorbate orbital j, the self-energy correction for the corresponding orbital in the SAM and the energy correction from the surface polarization, respectively. The ∆ symbol on the right-hand side of eq 2 is due to the fact that Ead,GGA(j) already contains contributions from the metal surface potential, and one only has to correct for the asymptotic tail of that potential. The HOMO-LUMO gap was found to

6420 J. Phys. Chem. C, Vol. 112, No. 16, 2008 increase dramatically from 3.8 to 8.3 eV. The corrections of occupied and unoccupied energy levels are asymmetric, with an upward shift of 3.1 eV for LUMO and a downward shift of 1.4 eV for HOMO. For the high-density PDI SAM, the HOMO-LUMO gap was found to be about 1 eV smaller, indicating a weak but finite screening within the SAM when molecules are closely packed. The screening from the surface polarization is estimated to reduce the LUMO energy level by about 0.5 eV. We note that a larger surface polarization energy (1-1.5 eV) is found for LUMO of benzene physisorbed on the graphite surface.56 In our case, the aromatic ring is separated further away from the surface because of the presence of the isocyanide group; therefore, the metal polarization effect plays a less important role in reducing the HOMO-LUMO gap. In calculating the GW energy corrections, we employed the random phase approximation (RPA) for the dielectric matrix, and the frequency dependence was described using a plasmonpole model (PPM). PPM works well for semiconductors and insulators but is not fully justified for molecular systems. To test the reliability of this model for molecular systems, we carried out a GW calculation for a gas-phase benzene molecule since experimental data and results of other calculations are available for comparison. The resulting HOMO and LUMO levels are in good agreement with the GW calculation by Niehaus et al. using Gaussian-type orbitals,57 but they are about 1 eV higher compared with experimental values and the results reported by Tiago et al.58 In ref 58, the frequency dependence was explicitly described by the Casida formalism;59 that is, no PPM was used, and contributions beyond RPA were added through an exchange-correlation kernel. Taking into account possible errors coming from the use of PPM in our GW calculations, overall, one expects about 1.5-2 eV upward shift for the unoccupied adsorbate states compared with the computed GGA value. This implies that a quite large bias would be required to make PDI on gold a good molecular conductor. There is still a debate in the literature regarding over the dominating electronic transport mechanism of the Au-PDIAu system. The temperature dependence of conductance for such a system in a nanopore structure was measured by Chen et al.,9 which suggests that non-Ohmic thermionic emission is the dominant process with a thermionic barrier of 0.35 eV. Other groups reported results consistent with a tunneling mechanism,10,13,14 but different estimates of the energy offsets near EF were predicted by different experiments. Therefore, an accurate prediction of energy level alignments at molecule/metal interface is extremely important in helping to design and interpret transport measurements. Calculations are underway to carry out full (metal plus adsorbate) GW calculations and thus obtain a more precise value for the correction to the energy levels obtained within DFT. 4. Conclusion In summary, by using first-principles calculations, we have studied the energetic, electronic, and vibrational properties of aromatic isocyanide SAMs adsorbed on the Au(111) surface. Our computed total energies, together with an analysis of vibrational spectra, indicate that the most stable geometry corresponds to low coverage, where the adsorption energy to a perfect surface is about 0.5 eV. This value may be considerably increased (up to 1.3 eV) by the presence of defects, for example, adatoms. As the SAM coverage increases, the adsorption energy is reduced because of increased intermolecular repulsion. Our predicted, stable coverage is amenable to experimental verification, as the shift of the metal work function upon SAM

Li et al. adsorption has opposite signs, depending on the coverage (Θ ) 1 or 1/3). Overall, the binding of isocyanide to gold is found to be weaker than that of thiols. In order to further validate our proposed structural model, we have also computed the vibrational frequencies of the isocyanide NtC stretching mode on different adsorption sites, and we have found trends as a function of the geometry in good agreement with experiment. In particular, the computed coverage dependence of the NtC stretching mode strongly supports the low coverage as being the most stable geometric configuration of isocyanide on Au(111). An important part of our work has been devoted to studying energy alignment at the molecule/metal interface beyond the common DFT/LDA or GGA used in the literature. By using many-body perturbation theory under the so-called GW approximation and a classical image potential model, we have estimated corrections to the Kohn-Sham (GGA) single particle energies of the adsorbate system and thus obtained a more realistic description of the energy levels at the molecule/ metal interface. The LUMO states are found to be shifted upward from EF by about 1.5-2 eV, yielding a qualitatively different picture of electron transport barriers, as compared with that extracted from single particle GGA energy levels.15,18 In particular, the estimated barrier appears to be rather high for building efficient devices. In addition, our calculations indicate an unexpected tendency of the SAM to disorder through tilting of the metal-carbon bond. This may hinder the formation of ordered monolayer films and thus complicate the experimental verification of improved charge injection. Overall our results highlight the need for further functionalization and/or SAM modifications, in order to build a promising interface for molecular transport, when using isocyanide on gold. Work is in progress to investigate modifications that may lead to an improved stability as well as to a more favorable energy level alignment at the interface. Acknowledgment. This work was funded by NSF/CPIMA Grant DMR-0213618 and by DOE/BES Grant DE-FG0206ER46262. Some of the calculations were performed at the NERSC and the SDSC facilities. We thank G. Y. Liu, A. Riposan, Z. Bao, and I. Dabo for useful discussions. Supporting Information Available: Estimate of van der Waals interactions between molecules, computational details of quasiparticle energy calculations, and Cartesian coordinates of optimized geometries. This material is available free of charge via the Internet at http://pubs.acs.org. References and Notes (1) Wang, W.; Lee, T.; Reed, M. A. Rep. Progr. Phys. 2005, 68, 523544. (2) Reed, M. A.; Zhou, C.; Muller, C. J.; Burgin, T. P.; Tour, J. M. Science 1997, 278, 252-254. (3) Chen, J.; Wang, W.; Klemic, J.; Reed, M. A.; Axelrod, B. W.; Kaschak, D. M.; Rawlett, A. M.; Price, D. W.; Dirk, S. M.; Tour, J. M.; Grubisha, D. S.; Bennett, D. W. Ann. N. Y. Acad. Sci. 2002, 960, 69-99. (4) Henderson, J. I.; Feng, S.; Bein, T.; Kubiak, C. P. Langmuir 2000, 16, 6183-6187. (5) Murphy, K. L.; Tysoe, W. T.; Bennett, D. W. Langmuir 2004, 20, 1732-1738. (6) Swanson, S. A.; McClain, R.; Lovejoy, K. S.; Alamdari, N. B.; Hamilton, J. S.; Scott, J. C. Langmuir 2005, 21, 5034-5039. (7) Stapleton, J. J.; Daniel, T. A.; Uppili, S.; Cabarcos, O. M.; Naciri, J.; Shashidhar, R.; Allara, D. L. Langmuir 2005, 21, 11061-11070. (8) Kim, S.; Ihm, K.; Kang, T.-H.; Hwang, S.; Joo, S.-W. Surf. Interface Anal. 2005, 37, 294-299. (9) Chen, J.; Calvet, L. C.; Reed, M. A.; Carr, D. W.; Grubisha, D. S.; Bennett, D. W. Chem. Phys. Lett. 1999, 313, 741-748. (10) Hong, S.; Reifenberger, R.; Tian, W.; Datta, S.; Henderson, J. I.; Kubiak, C. P. Superlattices Microst. 2000, 28, 289-303.

Interface between Aromatic Isocyanides and Au(111) (11) Dupraz, C. J. F.; Beierlein, U.; Kotthaus, J. P. Chem. Phys. Chem. 2003, 4, 1247-1252. (12) Lee, J.; Lientschnig, G.; Wiertz, F.; Struijk, M.; Janssen, R.; Egberink, R.; Reinhoudt, D.; Hadley, P.; Dekker, C. Nano Lett. 2003, 3, 113-117. (13) Kim, B.; Beebe, J.; Jun, Y.; Zhu, X.-Y.; Frisbie, C. J. Am. Chem. Soc. 2006, 128, 4970-4971. (14) Chu, C.; Ayres, J.; Stefanescu, D.; Walker, B.; Gorman, C.; Parsons, G. J. Phys. Chem. C 2007, 111, 8080-8085. (15) Xue, Y.; Ratner, M. A. Phys. ReV. B 2004, 69, 085403. (16) Bai, P.; Li, E.; Neerja Collier, P. IEEE Trans. Nanotechnol. 2005, 4, 422-429. (17) Xue, Y.; Ratner, M. A. Phys. ReV. B 2003, 68, 115406. (18) Morari, C.; Rignanese, G.-M.; Melinte, S. Phys. ReV. B 2007, 76, 115428. (19) Toher, C.; Filippetti, A.; Sanvito, S.; Burke, K. Phys. ReV. Lett. 2005, 95, 146402. (20) Toher, C.; Sanvito, S. Phys. ReV. Lett. 2007, 99, 056801. (21) Seminario, J.; Zacarias, A.; Tour, J. J. Am. Chem. Soc. 1999, 121, 411-416. (22) Zhou, J.-H.; Shi, L.-W.; Zhang, T.; Chen, M.-B. Chin. J. Chem. 2007, 25, 1223-1228. (23) Gilman, Y.; Allen, P. B.; Hybertsen, M. S., arXiv:cond-mat/ 0411320v2. (24) Hedin, L.; Lundquist, S. Solid State Physics. In Vol. 23; Seitz, F., Turnbull, D., Ehrenreich, H., Eds.; Academic Press: New York, 1969. (25) Perdew, J. P.; Burke, K.; Ernzerhof, M. Phys. ReV. Lett. 1996, 77, 3865-3868. (26) Baroni, S.; Corso, A. D.; de Gironcoli, S.; Giannozzi, P.; Cavazzoni, C.; Ballabio, G.; Scandolo, S.; Chiarotti, G.; Focher, P.; Pasquarello, A.; Laasonen, K.; Trave, A.; Car, R.; Marzari, N.; Kokalj, A. PWSCF and PHONON: plane-wave pseudo-potential codes, http://www.pwscf.org. (27) Fuchs, M.; Scheffler, M. Comput. Phys. Commun. 1999, 119, 6798. (28) Troullier, N.; Martins, J. L. Phys. ReV. B 1992, 46, 1754-1765. (29) Kanai, Y.; Cicero, G.; Selloni, A.; Car, R.; Galli, G. J. Phys. Chem. B 2005, 109, 13656-13662. (30) Monkhorst, H. J.; Pack, J. D. Phys. ReV. B 1976, 13, 5188-5192. (31) Khein, A.; Singh, D. J.; Umrigar, C. J. Phys. ReV. B 1995, 51, 4105-4109. (32) Gro¨nbeck, H.; Curioni, A.; Andreoni, W. J. Am. Chem. Soc. 2000, 122, 3839-3842. (33) Vargas, M. C.; Giannozzi, P.; Selloni, A.; Scoles, G. J. Phys. Chem. B 2001, 105, 9509-9513. (34) Cosandey, F.; Madey, T. E. Surf. ReV. Lett. 2001, 8, 73-93. (35) Yan, Z.; Perdew, J. P.; Kurth, S.; Fiolhais, C.; Almeida, L. Phys. ReV. B 2000, 61, 2595-2598.

J. Phys. Chem. C, Vol. 112, No. 16, 2008 6421 (36) Carling, K.; Wahnstro¨m, G.; Mattsson, T. R.; Mattsson, A. E.; Sandberg, N.; Grimvall, G. Phys. ReV. Lett. 2000, 85, 3862-3865. (37) Bengtsson, L. Phys. ReV. B 1999, 59, 12301-12304. (38) Baroni, S.; de Gironcoli, S.; Dal Corso, A.; Giannozzi, P. ReV. Mod. Phys. 2001, 73, 515-562. (39) White, I. D.; Godby, R. W.; Rieger, M. M.; Needs, R. J. Phys. ReV. Lett. 1998, 80, 4265-4268. (40) The ABINIT code is a common project of the Universite´ Catholique de Louvain, Corning Incorporated, and other contributors (URL http:// www.abinit.org). (41) Grimme, S. J. Comput. Chem. 2004, 25, 1463-1473. (42) Nara, J.; Higai, S.; Morikawa, Y.; Ohno, T. J. Chem. Phys. 2004, 120, 6705-6711. (43) Bilic´, A.; Reimers, J. R.; Hush, N. S.; Hafner, J. J. Chem. Phys. 2002, 116, 8981-8987. (44) Feibelman, P. J.; Hammer, B.; Nørskov, J. K.; Wagner, F.; Scheffler, M.; Stumpf, R.; Watwe, R.; Dumesic, J. J. Phys. Chem. B 2000, 105, 40184025. (45) Kresse, G.; Gil, A.; Sautet, P. Phys. ReV. B 2003, 68, 073401. (46) Dabo, I.; Wieckowski, A.; Marzari, N. J. Am. Chem. Soc. 2007, 129, 11045-11052. (47) Wang, Y.; de Gironcoli, S.; Hush, N.; Reimers, J. J. Am. Chem. Soc. 2007, 129, 10402-10407. (48) Dhirani, A.-A.; Zehner, R.; Hsung, R.; Guyot-Sionnest, P.; Sita, L. J. Am. Chem. Soc. 1996, 118, 3319-3320. (49) Eichler, A. Surf. Sci. 2003, 526, 332-340. (50) Ko¨hler, L.; Kresse, G. Phys. ReV. B 2004, 70, 165405. (51) Holovics, T.; Robinson, R.; Weintrob, E.; Toriyama, M.; Lushington, G.; Barybin, M. J. Am. Chem. Soc. 2006, 128, 2300-2309. (52) DuBose, D.; Robinson, R.; Holovics, T.; Moody, D.; Weintrob, E.; Berrie, C.; Barybin, M. Langmuir 2006, 22, 4599-4606. (53) Blyholder, G. J. Phys. Chem. 1964, 68, 2772-2777. (54) Fo¨hlisch, A.; Nyberg, M.; Bennich, P.; Triguero, L.; Hasselstro¨m, J.; Karis, O.; Pettersson, L. G. M.; Nilsson, A. J. Chem. Phys. 2000, 112, 1946-1958. (55) Frisch, M. J. et al. Gaussian 03, Revision C.02; Gaussian, Inc.: Wallingford, CT, 2004. (56) Neaton, J. B.; Hybertsen, M. S.; Louie, S. G. Phys. ReV. Lett. 2006, 97, 216405. (57) Niehaus, T. A.; Rohlfing, M.; Sala, F. D.; Carlo, A. D.; Frauenheim, T. Phys. ReV. A 2005, 71, 022508. (58) Tiago, M. L.; Chelikowsky, J. R. Phys. ReV. B 2006, 73, 205334. (59) Casida, M. E. Recent AdVances in Density Functional Methods, Part I; Chong, D., Ed.; World Scientific: Singapore, 1995.