Atomistic Modeling of Quantum Dots at ... - ACS Publications

primary factors: first, the large computational cost associated with simulating QDs of experimentally relevant sizes. (typically ... Computational cos...
0 downloads 14 Views 2MB Size
Subscriber access provided by READING UNIV

Article

Atomistic Modeling of Quantum Dots at Experimentally Relevant Scales Using Charge Equilibration Nathan Weeks, and Kevin Tvrdy J. Phys. Chem. A, Just Accepted Manuscript • DOI: 10.1021/acs.jpca.7b09060 • Publication Date (Web): 07 Nov 2017 Downloaded from http://pubs.acs.org on November 8, 2017

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

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

Page 1 of 27

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

The Journal of Physical Chemistry

Atomistic Modeling of Quantum Dots at Experimentally Relevant Scales Using Charge Equilibration Nathan Weeks† and Kevin Tvrdy†* †

Department of Chemistry & Biochemistry, University of Colorado at Colorado Springs, Colorado Springs, Colorado *[email protected]

Abstract Quantum dots (QDs) have been successfully employed within a vast array of fundamental and applied studies spanning all subdisciplines of chemistry. However, ab initio models of QD behavior are inherently limited by computational cost due to the large number of atoms within QDs of experimentally relevant size. This work builds upon the method of charge equilibration (qEQ) to account for system interactions unique to QDs (QD-qEQ) and demonstrates accuracy through calculated per-QD energies and dipole moments that agree generally with ab initio calculations and experimental observation, respectively. By forgoing electronic structure information, QD-qEQ exhibits a distinct advantage in its exceptionally low computational cost, which affords consideration of over 35,000 unique spherical wurtzite CdSe structures with radii ≤12.5Å.

A comparison of QD-qEQ calculations with

experimental data relating to the phenomenon of CdSe magic size crystals (MSCs) affords statistical and structural insight into why MSCs are observed. Consideration of structures ≤12.5Å reveals QD sizes corresponding with local minima in QD energy, correlating closely with experimentally observed MSCs. The physical origin of observed energy minima is assigned to QD structures with surfaces exhibiting large fractions of highly coordinated atoms, a physical trait postulated to yield fewer reaction sites for stepwise growth, resulting in MSC stability. The low computational cost along with the per-atom and per-structure electrostatic data afforded by QD-qEQ makes this method an enticing approach to address dynamic QD behavior and has the potential to impact a broad range of fields concomitant to those in which QD inclusion has already proven useful.

ACS Paragon Plus Environment

1

The Journal of Physical Chemistry

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

Page 2 of 27

Introduction Semiconducting quantum dots (QDs) exhibit tunable size dependent optical and electronic properties that serve as both a model for quantized behavior at the nanoscale1,2 and as active components within industrial scale displays,3 biological imaging and sensing schemes,4 therapeutics,5 and next generation photovoltaics.6,7 Synthetic benchmarks for QDs include the first observation of size dependent behavior in colloidal nanocrystals,8 the development of methods affording samples with narrow size distribution,9 and the refinement of such methods to enable growth from less toxic precursors.10 Further, the introduction of low temperature growth methods has revealed the existence of so called “magic-sized” QD crystals (MSCs) that exhibit exceptional thermodynamic stability, ultra-narrow emission linewidths, and may lead to the bulk scale creation of single-structure nanocrystals.11–14 While there have been extensive experimental efforts to identify reaction conditions affording precise control over the size, shape, and surface composition of QDs, such studies are not typically theory-driven, instead relying on an expensive and time consuming trial-and-error approach.15 Given this, the mechanism of QD growth steps has been the focus of many theoretical investigations16–20 but remains a challenging target owing to two primary factors: first, the large computational cost associated with simulating QDs of experimentally relevant sizes (typically hundreds of atoms) using ab initio methods; and second, the atomic arrangement of bulk and surface atoms within QDs of experimentally relevant sizes varies so significantly that sufficient sampling of possible QD structures requires approximately 102-106 calculations—depending on the size range under investigation. This work presents a new approach toward modeling QDs of experimentally relevant size that utilizes a modification of charge equilibration methods (qEQ) first developed by Rappé and Goddard21 to address the unique environments encountered by QDs (QD-qEQ). While not quantum mechanical in nature—and thus foregoing electronic structure information—QD-qEQ affords per-atom charge and energy (allowing QD energy and dipole moment to be considered) and is sufficiently computationally inexpensive to allow consideration of an unprecedented number of unique QD structures relevant to experimental phenomena, such as MSC growth. Methods for modelling QD growth generally fall into one of two categories: ab initio or kinetic, each having unique benefits and limitations. Computational cost typically limits ab initio studies of QD behavior to either a handful of representative structures whose properties are extrapolated to describe an experimentally continuous growth process, or to consideration of pristine crystallographic surfaces composed of QD-relevant

ACS Paragon Plus Environment

2

Page 3 of 27

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

The Journal of Physical Chemistry

atoms.22,23 For example, high level density functional theory simulating adatom and admolecule adsorption on QD surfaces has been used to afford insight into conditions favoring rod vs. dot growth.24 However, such methods remain prohibitively costly when addressing large numbers of unique QD structures of experimentally relevant size. In contrast, kinetic methods addressing reactions involved in QD growth have accurately predicted QD size distributions throughout the synthetic process.25 However, exclusive use of rate equations to describe QD behavior does not afford atomistic and mechanistic insight into how and why certain QD surface sites may favor or limit growth under a given condition. The development of an atomistic, accurate, and computationally efficient method to model QD dynamics can offer currently unavailable insight into the relation between experimental growth conditions and resultant sample characteristics. Such a method has the potential to impact nearly every area where QD characteristics have proven theoretically interesting or experimentally advantageous. The remainder of this work begins with a brief introduction to qEQ and its modification for the treatment of QDs. This method is applied to stoichiometric CdSe wurtzite QDs of high symmetry, the results of which are compared with ab initio calculations and experimental measurements. To afford a statistically significant number of unique QD structures for analysis, a 2-atom growth algorithm is presented that employs a spherical restriction parameter to generate all possible unique wurtzite CdSe QD structures within a given size range. Finally, QD-qEQ is used to analyze a comprehensive set of CdSe QD structures, resulting in strong correlation between theoretical predictions and experimental MSC observation. This analysis provides atomistic insight into the origins of MSC existence using methodology not currently available within out-of-the-box qEQ software. Collectively, this work considers per-atom charge and energy data from the consideration of over 35,000 unique CdSe QDs and lays the foundation for a new method that could significantly impact how QD growth dynamics are addressed.

Results and Discussion Computational Approach Foundation of QD-qEQ. Introduced by Rappé and Goddard,21 qEQ methods take advantage of matrix solution algorithms to globally minimize the energy of a spatially fixed system of point charges exhibiting atom like behavior. The point charge characteristic of this approach is derived by consideration of interatom electrostatic interactions, while the atom like nature of qEQ originates from the assignment of a charge dependent ionization energy to each atom. Since its inception, qEQ has been modified to describe systems exhibiting long range

ACS Paragon Plus Environment

3

The Journal of Physical Chemistry

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

Page 4 of 27

periodicity,26 to reduce the number of ad hoc parameters necessary for solution,27 and to supplement experimental ionization state data with ab initio calculations to enable accurate inclusion of cationic nuclei.28 This work modifies qEQ to describe nanoscale crystals of finite size and spherical shape (QD-qEQ). A brief description of the central equations and assumptions associated with QD-qEQ is presented here, while a full and pedagogical derivation is provided in Sections S1-S3 of the Supporting Information (SI). The foundational assumption for any charge equilibration method is that a spatially fixed system of

N

atoms, each of which can take on an independent charge, exhibits system behavior based on two central conditions. The first condition takes total system energy to be the sum of various single atom properties and pairwise atomic interactions. Denoting individual nuclei within the system as i such that

v

1≤ i ≤ N ,

total system energy is

v

expressed as a function of atom charge ( qi ) and position ( ri ) within a spatially variant dielectric ( ε ): N N 1 N N v v v v v Esystem ( q1 , q2 ,..., qN ; r1 , r2 ,..., rN ; ε ) = ∑ Ei ( qi ) + ∑∑ qi q j J ij ( rij ) + ∑ qi ∆ i ( R, ε ) 2 i =1 j =1 i =1 i =1

(Eq. 1)

where the ionization energy term, first on the righthand-side (RHS), sums over the charge dependent internal energy possessed by each atom in the system. The electrostatic energy term (second on RHS) accounts for all pairwise electrostatic interactions, where rij is the interatomic distance and the factor of one-half is included to compensate for overcounting within a complete sum. Finally, the dielectric energy term (third on RHS) accounts for the placement of each charged atom at position R within a spatially variant dielectric. Note that the dielectric term distinguishes QD-qEQ from other charge equilibration methods, owing to the unique nature of QDs existing as spherically symmetric dielectrics imbedded within a medium of mismatched dielectric. The physical origins of each term contributing to QD-qEQ total system energy are depicted in Fig. 1A. The second central condition of qEQ is a requirement that all terms on the RHS of Eq. 1 are represented by a truncated Taylor expansion. When truncated at the quadratic term this condition ensures that differentiation with respect to a system of single atom charges,

qi , generates a system of independent linear equations. Together, these

two conditions afford the use of existing algorithms to efficiently solve for the per-atom partial charges that globally minimize total system energy.

ACS Paragon Plus Environment

4

Page 5 of 27

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

The Journal of Physical Chemistry

Functional forms of Energy Terms.

The first

contribution to E system considered by QD-qEQ is the ionization energy of each atom, which is described using a truncated Taylor expansion of the ionization energy term centered about ionization state Q i0 :

Ei ( qi ) Q0 ≈ Ei ( Qi0 ) + ( qi − Qi0 ) i

dEi ( qi ) dqi Q0 i

+

1 ( qi − Q 2

)

0 2 i

d Ei ( qi ) dqi2 Q0 2

.

i

(Eq. 2) Note that the domain of per-atom charge of any atom within the QD-qEQ system is continuous; however, experimental or theoretical data for single-atom energy are only available at integer values of charge, restricting

Q i0 to integers.

Further discussion

regarding the functional form of Eq. 2 is provided in Section S3.1 of the SI. For a given element, Eq. 2 is related to experimental or computationally calculated ionization energies through algebraic combinations of Eq. 2 where q i = Qi0 + 1 and q i = Qi0 − 1 , yielding threepoint first and second derivatives centered about Q i0 :

Fig. 1. Schematic (A) and plots (B-D) of the three energyrelated factors considered within QD-qEQ. (A) Schematic depictions of the charge and spatial relationships constituting ionization, electrostatic, and dielectric mismatch energetic contributions. (B) Calculated single-atom energetics as a function of single-atom charge for both Cd and Se (open circles) along with quadratic fits of calculated data for truncated Taylor expansions centered about integer charges between -2 and +2. The Cd (green) and Se (purple) curves centered about charges of +2 and -2 respectively are highlighted as these curves were utilized in this work for consideration of CdSe QDs by QD-qEQ. (C) Electrostatic energy vs. interatomic separation for the three types of interactions considered—Cd-Se, Cd-Cd, and Se-Se— determined using Eq. 2 with a system dielectric of ε=6.25, where qi=qj=1. Internuclear separations shorter than one CdSe wurtzite bond length are not encountered in this work. (D) Dielectric mismatch energy resulting from the placement of a single charge of q=1 at various radial positions within dielectric spheres (ε2=6.25) of radii R=10, R=15, R=20, and R=25Å embedded within an infinite dielectric of ε1=2.6, calculated using Eq. 6.

 Ei ( Qi0 + 1) − Ei ( Qi0 )  −  Ei ( Qi0 − 1) − Ei ( Qi0 )  dEi ( qi )    ≡ χ Qi0 = i , dqi Q0 2

(Eq. 3)

i

0 d 2 Ei ( qi ) =  Ei ( Qi0 + 1) − Ei ( Qi0 )  +  Ei ( Qi0 − 1) − Ei ( Qi0 )  ≡ ηiQi .(Eq. 4) 2 dqi Q0 i

ACS Paragon Plus Environment

5

The Journal of Physical Chemistry

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

Page 6 of 27

In the special case of Qi0 = 0 , which was employed in the original qEQ formulation,21 the first and second bracketed terms on the RHS of both Eqs. 3 and 4 are the ionization potential and the electron affinity of atom i , respectively, resulting in the assignment of Eq. 3 as the electronegativity, χ i0 , and Eq. 4 as the idempotential (i.e. 0

atomic hardness), η i0 . For the case of Qi0 ≠ 0 , as is used here and elsewhere,27,28 the notation of χ iQi

0

and η iQi

represent the ionic electronegativity and ionic idempotential, respectively. Given the focus of this work on the properties of CdSe QDs, implementation of Eq. 2 requires use of computational methods to estimate ionization energies for Cd and Se at charges relevant to the system. While adequate experimental data exist for elemental cationic ionization states (i.e. Cd2+), such data do not exist for anionic ionization states (i.e. Se2- and Cdn-).29 In order to achieve consistency in this method across all nuclei types considered, the ionization energies of negatively charged Cd species were calculated even though the experimental interpretation of such energies is of limited value. Utilizing the methods presented by Wells and coworkers in their development of formal charge equilibration (FC-qEQ),28 this work employs single-atom ionization energies for both Cd and Se calculated using the coupled-cluster with single, double and perturbative triple excitations [CCSD(T)] method and the augmented cc-pVQZ basis set within the quantum chemistry package Gaussian 09.30 Atomic spin states were assigned to minimize the energy associated with each ionization state. Single atom ionization energy vs. charge curves for Cd and Se obtained from insertion of Eqs. 3 and 4, at various values of Q i0 , into Eq. 2 are presented in Fig. 1B while numerical lists of calculated parameters are provided in the Section S5 of the SI. The second contribution to E system considered by QD-qEQ is the electrostatic interaction between all atom pairs. In the original formulation of qEQ, Rappé and Goddard assumed that the charge on each atom was spread over a Slater type atomic orbital and described electrostatic interactions and shielding effects collectively with two center coulomb integrals.21 The complexity of this approach led to reports focused on identifying alternative empirical approximations that both increase accuracy and reduce computational cost.31 This work utilizes the DasGupta-Huzinaga empirical two-center Coulomb integral approximation:32

 2K J ij ( rij ) = K  rij + 0 k r Q0 k r  ηiQi e i ij + η j j e j ij 

−1

  ,  

(Eq. 5)

ACS Paragon Plus Environment

6

Page 7 of 27

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

The Journal of Physical Chemistry

where K is a collection of physical constants that includes system dielectric and has units of energy·distance,

ki = k j = 0.4 is the Klondike parameter, and pairwise interactions between an atom and itself ( rij → 0 ) are taken to be J ii ( 0 ) ≡ 0 . Further discussion of the implementation of this functional form within QD-qEQ is provided in Section S3.2 of the SI. Partial atomic charges calculated using qEQ and Eq. 5 have been found to be in good agreement with charges assigned using either Mulliken population analysis or electrostatic potentials calculated at the HF/6-31G (d,p) level of theory.31,33

A plot of

Jij ( rij )

for each of three interatomic electrostatic interactions

relevant to this study (Cd-Se, Cd-Cd, and Se-Se) is shown in Fig. 1C. The third and final contribution to E system considered by QD-qEQ originates from the placement of system charges within a medium of dielectric mismatch. Specifically, QDs exist as spherical crystals surrounded by capping ligands and immersed within a solvent—with each layer exhibiting a unique dielectric. Following the formulation of Brus,8 this model assumes the case where system kinetic energy is a fraction of an eV, and thus employs high-frequency dielectrics due to atomic core electrons.

For example, the high-frequency dielectric

constant of CdSe is 6.25,34 while that of the common CdSe QD capping ligand and synthetic solvent Trioctyl Phosphine Oxide (TOPO) is

2.6.35

The dielectric mismatch between CdSe QDs and their surrounding

environment requires consideration of the energetics of placing charges within an environment of discontinuous dielectric. Surrounding the QDs with the dielectric environment of TOPO is an approximation of ligand effects on QD properties. The energy associated with imbedding a single charge at radial distance sphere of radius R and dielectric

∆i ( S , R, ε1, ε 2 ) = K where

ε1

ε2

S

from the center of a

is given by:8

( ε 2 − ε1 ) 2Rε1

F 2 , ε1ε+1ε 2 ;1 + ε1ε+1ε2 ; RS 2  , 2

2 1

is the dielectric of the material surrounding the sphere and

2

(Eq. 6)

F1 [ a , b ; c ; z ] is a Gauss hypergeometric

function. A plot of the dielectric mismatch contributions for CdSe spheres of various radii imbedded within TOPO is shown in Fig. 1D. Note that both a derivation of Eq. 6 and a justification for the exclusion of screened pairwise interactions are provided in the Section S3.3 of the SI.

ACS Paragon Plus Environment

7

The Journal of Physical Chemistry

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

Page 8 of 27

Simultaneous Solution and Single-Parameter Calculations. Given functional forms of each term on the RHS of Eq. 1, a system of linear equations must be generated to simultaneously solve for the per-atom charges that globally minimize system energy. Briefly, the system chemical potential with respect to the charge on atom

µk =

∂Esystem ∂qk

N 0 0 =  qk χ kQk + qk − Qk0 ηkQk  + ∑ q j J kj ( rkj )    j =1

(

{

)

}

k

is written as:

 N v   + ∑ 2qk ∆k ( R, ε )  ,   j =1 

(Eq. 7)

where the first, second, and third bracketed terms in Eq. 7 are derived from the ionization energy, electrostatic energy, and dielectric mismatch terms, respectively. At equilibrium, the system chemical potential with respect to

µ1 = µ 2 = L = µ N , a condition that yields N − 1 independent equalities. An additional

each charge is equal,

equality can be written restricting total system charge to N

∑q

i

qtotal such that:

= qtotal , (Eq. 8)

i =1

which, when considered alongside the equalities formulated through application of Eq. 7, yields a total of N independent equations and N unknown per-atom charges. This system of equations can be expressed in terms of the matrix formulation:

D⋅Q = C , where the N × N

(Eq. 9)

D matrix has elements:

(

)

(

 J − J + δ χ Qi0 + 2∆ − δ η QN0 + 2∆ jN ij i i jN N N Dij =  ij 1  where

)

for 1 ≤ i ≤ N − 1  for i = N 

(Eq. 10)

δ ij is the Kronecker delta, the N × 1 C matrix has elements:  − χ Qi + χ NQN + Qi0ηiQN − QN0 η NQN Ci =  i qtotal  0

and the N × 1

0

0

0

for 1 ≤ i ≤ N − 1 , for i = N 

Q matrix contains the unknown partial atomic charges.

(Eq. 11)

Equation 9 can be solved efficiently using

built-in functions within commonly available computational software packages, yielding values for partial atomic charges that globally minimize the system energy according to Eq. 1.

ACS Paragon Plus Environment

8

Page 9 of 27

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

The Journal of Physical Chemistry

A clear example of system energy minimization using QD-qEQ considers the one-, two-, and threedimensional CdSe structures illustrated in Fig. 2A. Given the symmetry of these structures, restricting the total charge of each system to

qtot = 0 forces each Cd atom to take on an equal and opposite charge as that adopted by

each Se atom such that

qCd = − qSe . This unique feature allows these hypothetical systems to be reduced from

energy minimization problems in 2-, 4-, and 8-dimensional charge-space to minimizations with respect to the single variable q . A plot of total system energy with respect to q , along with contributing components: total Cd atom ionization energy, total Se atom ionization energy, sum of all pairwise electrostatic interactions, and sum of all dielectric mismatch interactions, is shown for each system in Fig. 2B. These examples demonstrate that each contributing factor to E system is quadratic and that total energy in each case takes the form of a concave up parabola, resulting in each system exhibiting a single (and global) energy minimum. For ease of comparison with other charge equilibration algorithms, numerical specifications of the D- and C-matrices (from Eqs. 10 and 11) used to compute qmin for each system are provided in Section S4 of the SI.

Fig. 2. Geometric setup (A) and QD-qEQ based energy minimization (B) for the three special cases in one-, two-, and threedimensions in which a system of charged atoms restricted to qtot=0 imparts equal and opposite charge on each Cd (+q) and Se (-q) atom. Due to this symmetry, these three systems afford straightforward visualization of individual-component and totalsystem energy minimization along a single abscissa of absolute per-atom charge. Single-atom energy vs. charge curves for Cd and Se were obtained from Eq. 2 centered about QCd0 and QSe0 of +2 and -2, respectively. A full symbolic and numerical detail of these example QD-qEQ calculations, including for intermediate steps, is provided in Section S4 of the SI.

ACS Paragon Plus Environment

9

The Journal of Physical Chemistry

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

Application

of

QD-qEQ

to

Page 10 of 27

Highly

Symmetric Wurtzite CdSe Spheres This work is limited to consideration of CdSe nanocrystals with stoichiometric ratios of Cd and Se atoms affording straightforward comparison of per-structure energies because of atomic arrangement and QD radius as opposed to an

imbalance

consideration

of

of only

nuclei.

Additionally,

stoichiometric

CdSe

clusters enables the rational implementation of

qtotal = 0 in Eqs. 8 and 11. Such stoichiometric crystals can be generated through by selective inclusion of wurtzite CdSe atoms falling inside of a spherical cutoff located at one of four wurtzite symmetry points (denoted here as origins O1, O2, O3, and O4) illustrated in Fig. 3A. Note that, while O1 and O4 within Fig. 3A are both located at the midway point along a Cd-Se bond, they are geometrically unique because of the slight difference in experimentally determined bond lengths along the short and long axes of wurtzite CdSe.36

Unique crystals were identified by

calculating the number of atoms included within spheres centered at each origin in radial increments of 10-5 Å ranging from 1.0-12.5Å, as shown in Fig. 3B. For a given series, a vertical discontinuity in number of included atoms

Fig. 3.

Description of the identification and analysis of unique wurtzite CdSe crystals using QD-qEQ. Crystals “cut” from spheres within the CdSe wurtzite lattice centered at one of four symmetry points (A) affords the generation of 165 unique structures with stoichiometric numbers of Cd and Se atoms (B). Calculation of perQD energy (C) and dipole (D). reveal radial dependent results in general agreement with ab initio methods (LDA-DFT) and experimental observation, respectively. Single-QD per-atom analysis of atom-energy and atom-charge for relatively small (E) and intermediate (F) sized QDs reveals energy and charge stabilization afforded by larger numbers of atoms.

ACS Paragon Plus Environment

10

Page 11 of 27

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

The Journal of Physical Chemistry

represents a unique QD structure. This method resulted in the identification of 72, 33, 31, and 29 unique CdSe QDs with origins O1, O2, O3, and O4, respectively, totaling 165 unique spherical stoichiometric CdSe structures with radii ≤12.5Å. A description of the method used in this work to identify structural uniqueness of CdSe crystals is provided in Section S6 of the SI. To identify radially dependent trends in either charge or QD energy, QD-qEQ was applied to each of the 165 unique crystals, taking total QD energy to be the normalized sum of all pairwise electrostatic interactions within a given structure:

EQD =

1 N N ∑∑ qi q j Jij ( rij ) , (Eq. 12) N i =1 j =1

where the zero-point energy for this analysis is taken to be the charged atoms at infinite separation. The dipole moment of each QD was determined using charge-vector addition originating at the geometric center (GC, mean Cartesian coordinate position of all atoms) of each crystal and corrected for effects of electric depolarization caused by spherical system dielectric mismatch:37

µQD = L where

3ε1 ε1 + ε 2

N

v v ∑q ( r − r ) , i

i

GC

(Eq. 13)

i =1

v v ri and rGC are the positions of atom i and the system GC, respectively, and

L is a collection of physical

constants with units of Debye·C-1·m-1. The calculated energy and dipole moment of each unique crystal identified in Fig. 3B is shown in Figs. 3C and 3D, respectively, plotted against QD radius. Formulation of the plots in Figs. 3C and 3D requires a definition of nanocrystal radius that enables straightforward comparison between the simulated QDs investigated here and experimental measurements determined through statistical analysis of high-resolution TEM images. Assuming that QD edge identification via microscopy identifies atoms lacking full coordination (a fully-coordinated wurtzite site has four nearest neighbors), this work took the radius of a given simulated QD to be the average distance from the GC of the crystal to each surface atom (singly, doubly, or triply coordinated atoms). The radial dependence of per-CdSe pair energies calculated using QD-qEQ (Fig. 3C) agrees generally with studies conducted using ab initio methods. Ramprasad and coworkers used the local density approximation (LDA) within density functional theory (DFT) to model 56 unique wurtzite CdSe QDs formed using spherical cutoffs centered about O1, O2, and O3, reporting per-pair energies ranging from -0.65 to -1.0eV, agreeing generally with this

ACS Paragon Plus Environment

11

The Journal of Physical Chemistry

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

work.38

Page 12 of 27

Further, the general trend of DFT-calculated pair-energy decreasing with radius (a ~50% reduction in per-

pair energy with increasing radii over the range: ~3Å-7Å, and another ~10% reduction over the range: ~7Å-12Å)39 is captured by QD-qEQ. This general agreement demonstrates that QD-qEQ is a viable method for modelling aspects of QD behavior driven by per-atom charge and energy. Additionally, there is agreement between permeant QD dipole moments calculated from QD-qEQ and those measured experimentally.

Given that the per-QD dipole calculated within a given size range varies

significantly, individual QD data points in Fig. 3D were binned within radial groups of width 0.5Å and averaged (open black circles in Fig. 3D), affording a more straightforward comparison with experimental reports. GuyotSionnest and coworkers used dielectric dispersion to measure the size dependence of the permeant dipole moment of CdSe QDs, reporting a gradual increase in such from 41-98D over the size range of 13-28Å,40 a trend captured by QD-qEQ. Additionally, QD-qEQ qualitatively agrees with experimentally measured dipole moments, predicting an average dipole of 56D at 12Å radii QDs versus an experimentally determined value of 41D for 13Å radii QDs.40 The structural origin of this trend has been shown to be due to the deviation from ideal tetrahedral geometry along the z-axis (c wurtzite lattice vector direction) in experimental structures (experimental bond lengths were used in this study).36,37 A striking aspect of calculated QD-qEQ data is its general agreement with theoretical and experimental reports despite the exclusion of a complete treatment of surface ligand or solvent effects, which are approximated by the inclusion of the dielectric mismatch term. While this may be due to comparisons of QD-qEQ with QD behavior driven by the nature of crystal structure (and not surface ligands), an alternative interpretation is that the exclusive consideration of stoichiometric structures in this work generates an inherent chemical balance within each QD that is otherwise delegated to surface ligand interactions. Regardless, the treatment of non-stoichiometric QDs and ligand effects with QD-qEQ remain potential topics of future studies. In addition to providing accurate values of QD energy and charge, QD-qEQ also affords atomistic insight in the form of individual-atom energies and charges. Images of representative QDs with radii of 5.25Å and 10.25Å are shown in Fig. 3E, while the per-atom energy and per-atom charge of each structure are shown in Figs. 3F and 3G, respectively. Each plot shows both atom type (as symbol color) and atom coordination (as symbol shape) as a function of atomic radial position from the crystals’ geometric center. The atomistic analyses shown in Figs. 3F and 3G provide detailed information that can be used to explain experimental QD behavior. For example, while the per-

ACS Paragon Plus Environment

12

Page 13 of 27

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

The Journal of Physical Chemistry

atom energy exhibited by a 5.25Å QD is within the same general range as that of a 10.25Å QD, the spread of peratom energy is much greater in the smaller crystal, a phenomenon likely due to a fewer number of fully coordinated (core) atoms contributing to surface stability. This broad distribution of per-atom energy leads to the presence of both energetically unstable (highly reactive) and stable (less reactive) surface sites on smaller crystals, agreeing with a size dependent experimental reports of QD stability.41 Additionally, the 10.25Å QD has a spread of per-atom energies that increases at greater radial distances from the GC, with surface atoms exhibiting both greater and lower energies relative to their core counterparts. A detailed analysis of surface site energy distribution in QDs may afford insight into site-dependent reactivities relevant to catalysis and sensing schemes. Such an analysis is delegated to a future study focused on the application of QD-qEQ to dynamic chemical processes and will be compared with other computational and experimental results relating to such processes.42 Given its low computational cost, QD-qEQ can be used to analyze a far greater number of structures than the 165 presented in Fig. 3. The following subsection outlines a method that affords (given a definition of “spherical” crystals) generation of all possible atomic arrangements of spherical wurtzite CdSe crystals within a specified diameter range, enabling a complete analysis of structurally dependent QD behavior.

Complete-set Generation of Spherical Wurtzite Structures Despite limiting wurtzite crystals to spherical shapes, there still exists a significant number of unique structural combinations within experimentally relevant QD size regimes. To ensure that calculations performed on a set of unique QD structures is representative of the possible morphologies grown within a colloidal synthesis, methods beyond those presented in Fig. 3 are needed. This work introduces a 2-atom stepwise structure generation algorithm driven by two axioms: first, rejection of structures exhibiting non-spherical geometry; and second, rejection of structures exhibiting chained atoms (singly coordinated atoms bound to doubly coordinated atoms). Utilizing growth increments of one Cd atom and one Se atom per step, this algorithm exclusively generates stoichiometric QDs. To increase the number of unique structures considered, Cd and Se atoms added within a given growth-step are not restricted to locations at neighboring lattice sites, although consideration of such may provide insight into the mechanistic aspects of QD growth.43 The 2-atom growth algorithm is summarized in the following three steps:

ACS Paragon Plus Environment

13

The Journal of Physical Chemistry

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

1.

Each of

Page 14 of 27

Nseed number of unique stoichiometric wurtzite CdnSen seed structures is considered separately to

identify its respective adatoms (vacant wurtzite lattice positions nearest-neighbor to seed structure atoms), which are categorized into Cd-adatom sites (totaling N aa , Cd ) and Se-adatom sites (totaling N aa , Se ). 2.

Each seed-structure is grown to form all possible combinations of Cdn+1Sen+1 by generating growthstructures through combinatorial additions of one Cd atom and one Se atom within respective adatom vacancies, totaling N aa ,Cd ⋅ N aa , Se growth-structures generated from each seed-structure.

3.

All generated growth-structures of common Cdn+1Sen+1 are collectively screened for: first, spherical nature using the spherical restriction parameter

Θ s detailed below; second, repeated identical structures detailed

in Section S.6 of SI); and third, presence of chained atoms. All unique spheres lacking chained atoms, which may be utilized as By seeding growth with

Nseed structures within Step 1 of a subsequent growth iteration are retained.

N seed = 2 Cd1Se1 structures (one each with internuclear separation corresponding to the

short and long bond in wurtzite CdSe),36 structural-sets of unique, spherical, and stoichiometric wurtzite CdnSen can be generated to arbitrarily large n (corresponding to arbitrarily large QD radius). A central component of this algorithm is choice of spherical restriction parameter,

Θ s , which quantifies

the spherical symmetry of a crystalline structure. Specifically, the maximum distance from the structural GC to all atoms in the QD, rc ,max , is compared with the minimum distance of all adatoms, raa ,min , to yield a single structurespecific spherical parameter:

Θ = raa,min − rc ,max ,

(Eq. 14)

a visual representation of which is shown in Fig. 4A. Through the specification of a spherical restriction parameter

Θ s and rejection of structures with Θ < Θ s , spherical growth within Step 3 of the algorithm is enforced. Contrasting with the spherical cutoff method (Fig. 3), the 2-atom growth algorithm generates every possible wurtzite structural combination that lacks chained atoms and meets the spherical condition enforced by choice of

Θs ,

thereby creating a complete structure set within a given radius range. A histogram demonstrating the relative abundance of unique spherical CdSe structures at each radius (with 0.05Å resolution) is shown in Fig. 4B for

Θs = 0 , while a table listing the number of structures deemed spherical,

ACS Paragon Plus Environment

14

Page 15 of 27

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

The Journal of Physical Chemistry

unique, and lacking chained atoms for system sizes

2 ≤ N ≤ 400 is provided in Section S7.1 of

the SI. Using

Θs = 0 , a total of 8,051 unique structures

are identified over the range of which exhibit radii ≤12.5Å. within

the

Θs = 0 for CdSe

2 ≤ N ≤ 400

2 ≤ N ≤ 400 ,

7,881 of

Consideration of structures range

affords

complete

consideration of QDs with radii ≤12.5Å, as shown in radii vs.

N

plots within Section S7.2 of the SI. Implementation of spherical restriction parameters

Θs < 0 relaxes the spherical symmetry necessary for structural inclusion and affords the generation of structuralsets containing a larger number of unique crystals. Given experimental observations of slightly asymmetric QDs formed using colloidal synthetic methods,44 consideration of

Θs < 0 is more representative of experimental conditions than

Θs = 0 . This work uses the bond-length-difference

between the short and long bonds of wurtzite CdSe,

l∆ = 1.423 × 10 −3 Å, as a natural length within this system, and thus considers

Θ s at integer values of l∆ . Enforcing

spherical restriction parameters of

Θ s = −10 l∆ and

Θ s = −15 l∆ within Step 3 of the growth algorithm results

Fig. 4. Schematic (A) demonstrating the definition of “spherical” for crystals utilized in this study and the resultant radial histograms (B-D, bin size of 0.05Å) of complete spherical wurtzite CdSe structure sets generated with a spherical restriction parameter, Θs, set to various integer values of the bond-length-differential of wurtzite CdSe. The strictest value of Θs=0.0Å rejects crystals whose innermost adatom lies at a smaller radial distance from the geometric structure center than that of the outermost contained atom during the stepwise 2-atom growth process. In comparison, Θs