Simulation of Electronic Transport in Silicon Nanocrystal Solids

Apr 30, 2012 - IMEP-LAHC, MINATEC Campus, 3 parvis Louis Néel, 38016 Grenoble Cedex 1, France. §. INL, INSA Lyon, 7 avenue Jean Capelle, 69621 ...
1 downloads 0 Views 3MB Size
Article pubs.acs.org/JPCC

Simulation of Electronic Transport in Silicon Nanocrystal Solids Hadrien Lepage,*,†,§ Anne Kaminski-Cachopo,‡ Alain Poncet,§ and Gilles le Carval† †

CEA, LETI, MINATEC Campus, 17 rue des Martyrs, 38054 Grenoble Cedex 9, France IMEP-LAHC, MINATEC Campus, 3 parvis Louis Néel, 38016 Grenoble Cedex 1, France § INL, INSA Lyon, 7 avenue Jean Capelle, 69621 Villeurbanne, France ‡

S Supporting Information *

ABSTRACT: Despite the emergence of quantum dot optoelectronics, so far very few attempts have been made to estimate theoretically the carrier mobility in nanocrystal solids. This paper presents a method to simulate hopping transport in 3D nanocrystal solids which is applied to silicon nanocrystals embedded in silicon oxide. Electronic states and charging energies are calculated separately for each nanocrystal. Electron transfer is calculated in the frame of Marcus theory with weak electronic couplings given by Bardeen’s formula. Additionally, the long-range nature of electrostatic interactions is explicitly taken into account. An accelerated kinetic Monte Carlo algorithm is employed to enable the simulation of an arbitrary degree of disorder. Electron mobility is found to be rather weakly impacted by size and positional disorder at room temperature, which is consistent with recent experimental results performed on PbSe nanocrystal solids. Electron mobility is also found to increase with decreasing nanocrystal size despite the increase of hops. Temperature dependent measurements, which have led to various interpretations in the literature, are indeed found to be tricky, as thermally activated mechanisms are numerous: inelastic hopping between nanocrystals, polaronic effects, or Coulombic repulsion, for instance. This suggests that simulation is a powerful tool to understand transport measurements and move toward the application of NC solids in optoelectronic devices. charge filling, and packing geometry. Its major contributions, compared with previous studies,9,10 are the ability to compute absolute mobility and to treat an arbitrary degree of disorder. The NC solid is modeled as a three-dimensional cubic lattice with periodic boundary conditions along the three space dimensions. The electronic structure and the electrostatic properties are calculated independently for each NC. The coupling between electronic states on adjacent NCs is calculated in the weak coupling limit, as the electron−phonon scattering rate, responsible for phase memory loss, is much higher than the transfer rate between two NCs. In other words, the thermal energy kT is much higher than the electronic coupling hΓ.11 The transfer rates Γ are then implemented in an accelerated superbasin kinetic Monte Carlo (ASKMC)12 code to calculate the carrier mobility. The range of calculations is performed for NC sizes giving significant as well as reasonable quantum confinement and for oxide thicknesses and degrees of disorder small enough to maximize the mobility while being experimentally achievable.

1. INTRODUCTION For the past decade, semiconductor quantum dots embedded in an insulating matrix have been investigated for potential optoelectronic applications1 such as solar cells2 or light emitting diodes. Bandgap widening due to quantum confinement (QC) may indeed lead to tunable materials. While the QC effect has been extensively studied, the electrical transport mechanisms through three-dimensional networks of semiconductor quantum dots are still under investigation.3 The dispersion of electrical results has so far prevented reliable characterization of the respective impact of disorder, quantum confinement, Coulomb blockade, and temperature on the mobility in such materials. In particular, percolation paths may shunt the major part of the current when measuring the resistivity.4 Therefore, simulation is mandatory to isolate the impact of each parameter and to elucidate the electrical characterization. For instance, the electrical behavior of various semiconductor nanocrystal solids has often been ascribed to disorder in the frame of variable range hopping (VRH) 5 or nearest-neighbor hopping (NNH),6,7 but recent experiments show that carrier mobility is independent from size disorder.8 This paper presents a comprehensive study of the carrier mobility in silicon nanocrystal solids in the weak coupling limit with the help of a kinetic Monte Carlo (KMC) code. It is dedicated to incoherent transport between nanocrystals (NC) where the carriers diffuse and drift by hopping. It discusses the mobility dependence with temperature, NC size, size disorder (diagonal disorder), position disorder (off-diagonal disorder), © 2012 American Chemical Society

2. COMPUTATIONAL METHODS The electronic wave function of a NC is described by the product of a function which is dependent on the crystal periodicity and an envelope function which is dependent on the Received: February 21, 2012 Revised: April 30, 2012 Published: April 30, 2012 10873

dx.doi.org/10.1021/jp301713v | J. Phys. Chem. C 2012, 116, 10873−10880

The Journal of Physical Chemistry C

Article

size of the NCs.13 For the problem of electron transfer between NCs, the relevant information is contained within the envelope function. The Schrödinger equation of the envelope wave function Ψ for an electron in a spherical NC of radius R is, in the effective mass approximation,14 ⎛ ℏ2 ⎞ −1 ⎜ − ∇⃗[m] ∇⃗ + U (r )⎟Ψ( r ⃗) = E Ψ( r ⃗) ⎝ 2 ⎠

(1)

where U(r) is the finite barrier potential at the interface between the silicon and the insulating matrix

⎧0 r ℏωphonon).22 The tunnelling Hamiltonian Hij, treated as a perturbation (nonadiabatic electron transfer), is given by the Bardeen formula23

ΔEij0 = Ej − Ei

(16)

The second one is related to the applied electric field Fext along z. F ΔEab = qFext(zb − za)

(17)

And the third one is the change in the electrostatic interactions C ΔEab = (Σb + nbEbe − e) − (Σa + (na − 1)Eae − e) + qΔVab

(18)

It is itself the sum of two contributions, one arising from the on-site charging energies calculated before and the other arising from the Coulomb interactions between carriers localized on different nanocrystals. The potential difference ΔVab to evaluate for a carrier hopping from the NC a to the NC b is 10875

dx.doi.org/10.1021/jp301713v | J. Phys. Chem. C 2012, 116, 10873−10880

The Journal of Physical Chemistry C N

ΔVab = Vb − Va − V ba =

Article N

The master equation associated with the Markovian dynamics described previously is usually solved using a rejection free kinetic Monte Carlo algorithm. At each iteration, a first uniformly distributed random number is drawn on the [0, Γtot=ΣΓab] interval to select the event to be executed, a second random number is drawn to increment the time t by the elapsed time Δt, and the transition probabilities per unit time are updated. However, a major limitation of this standard KMC appears when the time scale between rare and frequent events becomes large. Due to the exponential dependence of the tunnelling probabilities Γab with the oxide thickness between adjacent NCs, this situation occurs as soon as positional and/or size disorder is introduced. In such a situation, KMC becomes unable to solve the problem because the algorithm gets stuck in the short time scale events for a very large number of iterations so that the time t progress very slowly. In other words, a carrier hops many times among a group of NCs which are efficiently coupled to one another. To circumvent this problem, Voter and Chaterjee have introduced the so-called accelerated superbasin kinetic Monte Carlo (ASKMC)12 method. This method provides acceleration by gradually lowering the rate constants observed many times during the simulation. As a consequence, fast processes are slowed down and slow processes are executed more often. Inside a SB of efficiently coupled NCs where the electron hops many times, a local equilibrium is achieved before the electron escapes from it; the ASKMC method is designed to maintain this equilibrium while lowering the transition rates belonging to the SB. In this way, the error remains reasonable12 (see Figure 13). The flowchart of the algorithm is given in Scheme 1. Once an event is observed more than Nf times, a search algorithm is performed to determine whether or not the selected process Γab is part of a superbasin (SB) and whether or not the processes Γa′b′ belonging to the SB were all observed more than Nf times too. The processes Γa′b′ belong to the selected process’ SB if they are connected and if they satisfy the SB criterion

∑ niVbi − ∑ niVai − Vba i=1

i=1

⎛ ⎞ ⎜ q ⎟ 1 ⎟ = ∑ ni⎜ ∑ ⃗ πε ε 4 ⎜ ⎟⎟ | + | r R ⃗ out 0 ⃗ ib i=1 ⎜ R | rib⃗ + R⃗| ≠ 0 ⎝ ⎠ ⎛ ⎞ N ⎜ q 1 ⎟⎟ − ∑ ni⎜ ∑ | ria⃗ + R⃗| ⎟⎟ i=1 ⎜ R⃗ ⎜ 4πεoutε0 ⃗ | + | ≠ r R 0 ⃗ ⎝ ⎠ ia q 1 − ∑ 4πεoutε0 R⃗ | rab⃗ + R⃗| N

⎡ ⎤ ⎢ ⎥ N ⎢ ⎞⎥ ⎛n − δ q n i ⎟⎥ = ∑⎢ ∑ ⎜ i ia − ⃗ | + | | + r R r R ⃗ | ⎠⎥ ⃗ ⃗ ⎝ ⃗ ib ia i = 1 ⎢ 4πεoutε0 R ⎢ ⎥ | rib⃗ + R⃗| ≠ 0 ⎢⎣ ⎥⎦ | ria⃗ + R⃗| ≠ 0 q(na − nb − 1) + 4πεoutε0| rab⃗ | (19)

where Vb is the potential created in the NC b by the charges out of the NC b. Vij is the potential created in the NC j by one charge in the NC i, rib = (xb − xi, yb − yi, zb − zi), and N is the total number of NCs contained in the simulation box. Because of the 3D periodic boundary conditions, the summation is over the replica of the simulation box denoted R = (mxX, myY, mzZ) with integers mx,y,z. The summations to compute Vb and Va are divergent due to the non-neutral nature of the system, but ΔVab converges as the alternate series in square brackets is conditionally convergent. However, this series converges very slowly due to the long-range nature of Coulomb interactions. To circumvent this problem which makes the computation too demanding, the summation for the calculation of potentials Vij is replaced by the Ewald sum24 ⎛ erfc(β| rij⃗ + R⃗|) q ⎜ 4π V ji = ∑ + ⎜ 4πεoutε0 ⎝ R⃗ Ω | rij⃗ + R⃗| ⎞ ⎛ |k |⃗ 2 ⎞ exp⎜ − 2 ⎟ cos(k ⃗· rij⃗ )⎟⎟ ⎝ 4β ⎠ ⎠

∑ k ⃗ ≠ 0⃗

γπa ′Γa ′ b ′ > πa Γab

(21)

where γ is the superbasin tolerance chosen by the user and πa is the probability of residing in the NC a at equilibrium πa =

1 |k |⃗ 2

⎛ ⟨E ⟩ ⎞ 1 exp⎜ − a ⎟ ZSB ⎝ kBT ⎠

(22)

πa follows the classical Maxwell−Boltzmann statistic because the NCs are weakly coupled. ZSB is the partition function for the SB (20)

ZSB =

where β is the Ewald parameter, k = (2πkx/X, 2πky/Y, 2πkz/Z) with integers kx,y,z is the reciprocal space vector, and Ω = XYZ is the simulation box volume. The restriction of real-space and reciprocal-space sums within two images (mx,y,z∈[−2, 2] and kx,y,z∈[−2, 2]) along with β = 2/Ω1/3 gives a good accuracy. The potentials Vij given in eq 20 are evaluated and stored once before the KMC iterations begin. Then, after every hopping event, the potential Va on each site a is easily updated, the stored contributions Via of the carrier at its initial position i are removed, and the stored contributions Vfa at its final position f are added. In this way, the computational effort is low, while the procedure is equivalent to solve the Poisson equation after every KMC event.25 This procedure actually benefits from the discrete nature of the nanocrystal lattice.

⎛ ⟨E ⟩ ⎞ exp⎜ − a ⎟ ⎝ kBT ⎠ a ∈ SB



(23)

And ⟨Ea⟩ is the average energy of one electron in the NC a ⟨Ea⟩ =

∑i gifi (1)Ei ∑i gi

(24)

If the processes Γa′b′, identified by the search algorithm as belonging to the selected process’ SB, are observed more than Nf times, all the transition probabilities belonging to the SB are lowered by the parameter α chosen by the user and the number of times they were observed is reset. Otherwise, the SB criterion is not checked and the transition probabilities are not lowered. Also, the search algorithm must ensure that the SB is neither the all simulation box nor a percolation path connecting 10876

dx.doi.org/10.1021/jp301713v | J. Phys. Chem. C 2012, 116, 10873−10880

The Journal of Physical Chemistry C

Article

3. RESULTS AND DISCUSSION The NCs were located on a three-dimensional 3 × 3 × 3 cubic lattice with periodic boundary conditions (see insets of Figure 11) in order to consider the specific 3D percolation effect. Only the nearest neighbors are coupled, the coordination number being equal to 6 in a cubic lattice. The electric field is applied along the Z direction, and the mobility is given by

Scheme 1. Flowchart of the KMC Algorithm (White Blocks Only) and the ASKMC Algorithm (Both White and Green Blocks)a

μe = q × net carrier transition through the simulation box × Z t × total number of carriers in the simulation box × Fext (25)

The mobility decreases exponentially with increasing oxide thickness between adjacent NCs as expected. The quantitative mobility calculated in Figure 6 is several orders of magnitude

Figure 6. Electron mobility as a function of the oxide thickness between adjacent NCs.

below the one calculated by Jiang and Green;26 this difference is ascribed to the band-like transport hypothesis employed in their paper which may be valid for kT < hΓ,11 that is at low temperature only for the material considered. Contrary to experimental results of Liu et al.8 and Kang et al.,27 the results in Figure 7 demonstrate that electron mobility

a

The superscript m refers to the number of times the transition probability has been lowered.

Figure 7. Electron mobility as a function of the NC radius.

increases as the electronic coupling increases with decreasing NC diameter (see Figure 5). Despite the different electronic structures of Si and PbSe, this general feature should be observed for both. As a consequence, the increase in total hops needed for carriers to flow through NC solid with decreasing NC diameter cannot solely account for their observation. Figure 8 shows two resonant peaks in the electric field dependence. The first peak appears when the driving electric field counterbalances the reorganization energy to be paid at

the two edges of the simulation box. Otherwise, it would corrupt the true dynamics. This algorithm enables the simulation of an arbitrary degree of disorder. The ASKMC calculations were performed with an integration time of t = 1 s along with α = 1.5, γ = 2, and δ = 0.1 which defines Nf (see Scheme 1). 10877

dx.doi.org/10.1021/jp301713v | J. Phys. Chem. C 2012, 116, 10873−10880

The Journal of Physical Chemistry C

Article

Figure 8. Electron mobility as a function of the applied electric field for different temperatures.

Figure 9. Electron mobility as a function of the applied electric field, for different degrees of size (energetic) disorder and with fixed oxide thickness between adjacent NCs.

each charge transfer |qFext(zb − za)| = |λab|. The second peak of mobility appears when the electron mainly flows from a LUMO state to a LUMO+1 state |qFext(zb − za)| = |λab + (E(LUMO+1) − ELUMO)|. The position of these peaks is then a function of the packing geometry (lattice geometry, NC size, and distance between adjacent NCs) as well as a function of the materials (NC and matrix). The relative amplitude of these peaks is linked to a better penetration of the wave functions in the oxide with the level of excitation (see Figure 2). At ambient temperature, this effect becomes less pronounced as electronic levels are more broadened (see Figure 4). To assess the influence of size disorder, the diameter of each NC is randomly drawn from a log-normal distribution f (d ) =

⎛ (ln d − μ)2 ⎞ exp⎜ − ⎟ 2σ 2 ⎝ ⎠ d 2πσ 2 1

(26)

with a relative size dispersion D and a mean radius R as input parameters D=

exp(σ 2) − 1

(27)

R=

⎛ 1 σ2 ⎞ exp⎜μ + ⎟ 2 2⎠ ⎝

(28)

Figure 10. NC size distribution (red lines) compared with NC size distribution contributing to carrier flow (black bars).

Two situations are addressed: either the oxide thickness is fixed or the oxide thickness is let free. The former option is a numerical artifice to separate energetic (diagonal) and coupling (off-diagonal) disorder. To achieve a good statistical convergence whatever the disorder strength, the calculation of the mobility is repeated until standard deviation of μe mean value of μe ×

number of draws

≤ 0.1 (29)

Figure 9 shows that resonance effects are logically smoothed by disorder and that energetic disorder has a narrow influence on electron mobility. Indeed, at room temperature, the carriers flow through most of the NCs, leading to a weak reduction of the carrier mobility (see Figures 10c and 11). These results are consistent with the recent observation of Liu et al.8 on PbSe NCs, since the sensitivity to disorder is assumed to be comparable. However, at low electric field and low temperature (see Figures 10a and 11), the smallest NCs are actually avoided, leading to a dramatic reduction of carrier mobility as the carriers experience slower hops (higher activation energies between NCs) and take winding pathways. To compare with experimental observations in the literature,5−7,28−31 the mobility dependence on the temperature

Figure 11. Arrhenius plot of the electron mobility as a function of the temperature for different degrees of size (energetic) disorder with fixed oxide thickness between adjacent NCs. Lines are exponential fits.

at low electric field is plotted in an Arrhenius plot (Figure 11). First of all, the polaronic effects induce a thermally activated behavior in the absence of any disorder. The addition of energetic disorder reinforces this thermally activated behavior, especially with fixed oxide thickness, but it does not lead to sublinear behavior typical of Mott-VRH (log μe ∼ −1/T1/2), 10878

dx.doi.org/10.1021/jp301713v | J. Phys. Chem. C 2012, 116, 10873−10880

The Journal of Physical Chemistry C

Article

Efros-Shklovskii-VRH (log μe ∼ −1/T1/4), or percolation-NNH (log μe ∼ −1/T1/2). Note that the VRH mechanism is unlikely in Si/SiO2 NC solids due to the high tunnel barrier in SiO2, but it is theoretically feasible in other binary systems.32 Indeed, in order to test this hypothesis, second nearest-neighbor hopping probabilities were computed without any impact on the results. If the distance between NCs is let free, the carriers flow through the biggest NCs (see Figure 10d), which are the closest, leading to efficient percolation paths. As a consequence, the mobility is higher because tunnelling barriers are thinner within these percolation paths (Figure 12). Additionally, the second peak appears at lower electric field, since the energy difference between LUMO and LUMO+1 states is narrower among the largest NCs (see Figure 1).

Figure 14. Electron mobility as a function of the applied electric field, for different levels of charge filling.

handled by basic arguments. However, when charge filling reaches one charge per NC, the electron mobility in Figure 15

Figure 12. Electron mobility as a function of the applied electric field with and without fixed oxide thickness between adjacent NCs. Figure 15. Electron mobility as a function of the applied electric field, for different temperatures with one electron per NC.

To assess the effect of positional disorder, each NC coordinate is shifted from its ideal position with respect to the lattice geometry. This shift follows a Gaussian distribution described by its full-width at half-maximum (fwhm). Figure 13

demonstrates a clear Coulomb blockade effect similar to single electron transistors (SET).33 Note that charge filling induces strong activated behavior at low electric field. As a consequence, conductivity measurement as a function of temperature is not self-sufficient to distinguish between VRH, NNH, polaronic, or Coulomb blockade effects.

4. CONCLUSION This paper presents a general method to numerically describe hopping transport through NC solids. It provides insight onto the factors impacting the mobility to prevent from misinterpretation. In particular, for silicon nanocrystals embedded in SiO2, diagonal and/or off-diagonal disorder induces rather weak reduction of mobility at room temperature as observed by Liu et al. for PbSe NCs.8 On the other hand, the mobility increases with decreasing NC diameter. Indeed, basic calculation shows that electronic coupling increases with decreasing NC diameter as long as the oxide material is not altered; it is also confirmed in recent studies.34 Furthermore, the results suggest that experimental conductivity as a function of temperature is not sufficient to make the difference between disorder, polaronic, or even Coulomb blockade effects, especially if we assume that a real device brings interfacial defects35,36 and injection/ejection phenomena into play. Finally, this paper provides quantitative results about the electron mobility in silicon NCs embedded in SiO2 which could be useful inputs to estimate the potential applications.

Figure 13. Electron mobility as a function of the applied electric field, for different degrees of positional disorder.

shows the results for different degrees of positional disorder without size disorder. As soon as the fwhm is higher than 0.1 nm, the ASKMC algorithm becomes necessary to keep the calculation time reasonable. It is worth noting that the comparison between ASKMC and usual KMC calculation in Figure 13 demonstrates the validity of the accelerated algorithm. 3D percolation paths actually enable the carrier to flow without a dramatic drop of mobility, and resonant peaks are smoothed like in the case of energetic disorder. Finally, the effect of charge filling is assessed. As can be seen in Figure 14, the effect of charge filling cannot be on the whole 10879

dx.doi.org/10.1021/jp301713v | J. Phys. Chem. C 2012, 116, 10873−10880

The Journal of Physical Chemistry C



Article

(25) Meng, L.; Wang, D.; Li, Q.; Yi, Y.; Bredas, J.-L.; Shuai, Z. J. Chem. Phys. 2011, 134, 124102. (26) Jiang, C.-W.; Green, M. A. J. Appl. Phys. 2006, 99, 114902. (27) Kang, M. S.; Sahu, A.; Norris, D. J.; Frisbie, C. D. Nano Lett. 2010, 10, 3727−3732. (28) Rafiq, M. A.; Durrani, Z. A. K.; Mizuta, H.; Hassan, M. M.; Oda, S. J. Appl. Phys. 2008, 104, 123710. (29) Porter, V. J.; Mentzel, T.; Charpentier, S.; Kastner, M. A.; Bawendi, M. G. Phys. Rev. B 2006, 73, 155303. (30) Mentzel, T. S.; Porter, V. J.; Geyer, S.; MacLean, K.; Bawendi, M. G.; Kastner, M. A. Phys. Rev. B 2008, 77, 075316. (31) Romero, H. E.; Drndic, M. Phys. Rev. Lett. 2005, 95, 156801. (32) Liu, H.; Pourret, A.; Guyot-Sionnest, P. ACS Nano 2010, 4, 5211−5216. (33) See, J.; Dollfus, P.; Galdin, S. IEEE Trans. Electron Devices 2006, 53, 1268−1273. (34) Chu, I.-H.; Radulaski, M.; Vukmirovic, N.; Cheng, H.-P.; Wang, L.-W. J. Phys. Chem. C 2011, 115, 21409−21415. (35) Zimina, A.; Eisebitt, S.; Eberhardt, W.; Heitmann, J.; Zacharias, M. Appl. Phys. Lett. 2006, 88. (36) Godefroo, S.; Hayne, M.; Jivanescu, M.; Stesmans, A.; Zacharias, M.; Lebedev, O. I.; Van Tendeloo, G.; Moshchalkov, V. V. Nat. Nanotechnol. 2008, 3, 174−178.

ASSOCIATED CONTENT

S Supporting Information *

Figures showing Coulomb repulsion charging energy as a function of the oxide thickness between adjacent NCs in a cubic lattice for different NC radius (jp301713v_si_001.pdf), self-energy as a function of the oxide thickness between adjacent NCs in a cubic lattice for different NC radius (jp301713v_si_002.pdf), and electronic coupling |H| of two LUMO states on identical NCs separated by SiO2 for different NC radius (SI correc.tif). This material is available free of charge via the Internet at http://pubs.acs.org.



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS The research was supported by the CEA LETI. We thank F. Triozon, P. Blaise, and F. de Crécy for many fruitful discussions.



REFERENCES

(1) Talapin, D. V.; Lee, J.-S.; Kovalenko, M. V.; Shevchenko, E. V. Chem. Rev. 2010, 110, 389−458. (2) Kramer, I. J.; Sargent, E. H. ACS Nano 2011, 5, 8506−8514. (3) Balberg, I. J. Appl. Phys. 2011, 110, 061301. (4) Strelniker, Y. M.; Havlin, S.; Berkovits, R.; Frydman, A. Phys. Rev. E 2005, 72, 016121. (5) Yu, D.; Wang, C.; Wehrenberg, B. L.; Guyot-Sionnest, P. Phys. Rev. Lett. 2004, 92, 216802. (6) Fujii, M.; Mamezaki, O.; Hayashi, S.; Yamamoto, K. J. Appl. Phys. 1998, 83, 1507−1512. (7) Zhou, X.; Usami, K.; Rafiq, M. A.; Tsuchiya, Y.; Mizuta, H.; Oda, S. J. Appl. Phys. 2008, 104, 024518. (8) Liu, Y.; Gibbs, M.; Puthussery, J.; Gaik, S.; Ihly, R.; Hillhouse, H. W.; Law, M. Nano Lett. 2010, 10, 1960−1969. (9) van de Lagemaat, J. Phys. Rev. B 2005, 72, 235319. (10) Chandler, R. E.; Houtepen, A. J.; Nelson, J.; Vanmaekelbergh, D. Phys. Rev. B 2007, 75, 085325. (11) Vanmaekelbergh, D.; Liljeroth, P. Chem. Soc. Rev. 2005, 34, 299−312. (12) Chatterjee, A.; Voter, A. F. J. Chem. Phys. 2010, 132, 194101. (13) Tisdale, W. A.; Zhu, X.-Y. Proc. Natl. Acad. Sci. U.S.A. 2011, 108, 965−970. (14) Anchala; Purohit, S. P.; Mathur, K. C. Appl. Phys. Lett. 2011, 98, 043106. (15) Lannoo, M.; Delerue, C.; Allan, G. Phys. Rev. Lett. 1995, 74, 3415−3418. (16) Franceschetti, A.; Williamson, A.; Zunger, A. J. Phys. Chem. B 2000, 104, 3398−3401. (17) Wang, L.-W.; Zunger, A. Phys. Rev. Lett. 1994, 73, 1039−1042. (18) Merrill, W. M.; Diaz, R. E.; LoRe, M. M.; Squires, M. C.; Alexopoulos, N. G. IEEE Trans. Antennas Propag. 1999, 47, 142−148. (19) Barbara, P. F.; Meyer, T. J.; Ratner, M. A. J. Phys. Chem. 1996, 100, 13148−13168. (20) Brus, L. Phys. Rev. B 1996, 53, 4649−4656. (21) Coropceanu, V.; Cornil, J.; da Silva Filho, D. A.; Olivier, Y.; Silbey, R.; Brédas, J.-L. Chem. Rev. 2007, 107, 926−952. (22) Jdira, L.; Overgaag, K.; Stiufiuc, R.; Grandidier, B.; Delerue, C.; Speller, S.; Vanmaekelbergh, D. Phys. Rev. B 2008, 77, 205308. (23) Valentin, A.; Galdin-Retailleau, S.; Dollfus, P. J. Appl. Phys. 2009, 106, 044501. (24) Casalegno, M.; Raos, G.; Po, R. J. Chem. Phys. 2010, 132, 094705. 10880

dx.doi.org/10.1021/jp301713v | J. Phys. Chem. C 2012, 116, 10873−10880