Metal oxide nanoparticles have wide range of commercial and technological applications ranging from biomedical engineering to chemical catalysis. [1][2][3][4][5][6][7] Oxide nanoparticles have a number of favorable properties including a high surface-area to volume ratio, tunable optical properties, and a range of surface reactivities. For example, oxide nanoparticles can be used as inert support materials in catalysis, 8 or provide reactive surfaces that interact with supported catalysts, 9 or act directly as catalysts. 8,10 In particular, oxide nanoparticles are commonly used for photocatalysis since it requires high surface areas, bandgaps that are aligned with the redox potentials of products and reactants, and surface chemistry that facilitates chemisorption and reaction of intermediates. 11,12 Both the band structure and chemical reactivity of nanoparticles are controlled by their atomic structure and size, as well as the surrounding environment such as solvents, capping agents and supports. However, atomic-scale simulations of nanoparticles are challenging due to the number of atoms and electrons present, especially if their environment is considered. Therefore, it is a natural starting point to study the structure and reactivity of isolated nanoparticles to establish an understanding of the key factors governing band structure and reactivity to help design nanoparticles with specific optical and catalytic properties. 13 In particular, TiO 2 is of paramount technological importance. The applications of TiO 2 include hydrogen production, photo-voltaic cells, degradation of harmful organic compounds in the environment, medical applications such as bone implants, and supports for other catalytic materials. 1,[14][15][16][17][18][19] For this reason, TiO 2 is widely studied in the scientific literature and is often referred to as a "model oxide". [20][21][22] A large number of ab-initio studies focus on models of bulk polymorphs, [23][24][25][26] extended surfaces [27][28][29] and small nanoparticles but relatively few studies have been done on investigation of finite size e↵ects arising in small TiO 2 nanoparticles, [30][31][32] especially for the rutile polymorph of TiO 2 . These finite size e↵ects can be classified into geometric e↵ects which arise due to unique atomic configurations in the nanoparticle, electronic e↵ects, which arise due to features in the electronic structure of the nanoparticle system, and quantum e↵ects, which arise due to quantization of energy levels in small particles. 33 Since rutile TiO 2 has a large number of applications in the field of heterogeneous catalysis and photocatalysis, [34][35][36] it represents a natural starting point for evaluating the finite size e↵ects on surface and electronic properties. In this work, we study the finite size e↵ects in surface energy which is an important property that dictates the surface structure and stability. The electronic property in consideration is the band structure, which is crucial in applications of photocatalysis because nanoparticles of TiO 2 are often used, but the interplay between nanoparticle structure and photon absorption is not well understood.
Previous works in the literature have investigated the intrinsic particle size e↵ect with increasing size of metallic nanoparticles using adsorption energy to characterize surface catalytic properties. Kleis et al. 37 performed a DFT study using the revised Perdew-Burke-Ernzerhof (RPBE) functional 38 on gold metal nanoparticles ranging from 13 to 1,415 atoms to show how surface properties varied with system size at two local geometries that resemble surfaces of ( 111) and (211) slabs. They show that surface properties converge to the slab limit at a characteristic length of 27 Å (560 atoms). The generality of findings for other transition metals were confirmed with a similar study on freestanding cuboctahedral platinum metal nanoparticles which show analogous convergence with size but at a smaller characteristic length of 16 Å (147 atoms). 33 The main limitation of the prior work described above is that they focus on metallic systems and are limited to semilocal GGA exchange-correlation functionals. One exception is work by Lamiel-Garcia et al., 30 that studied finite size e↵ects in anatase nanoparticles at the hybrid level of theory with localized basis sets. However, the convergence of surface properties toward the slab limit was not investigated, and the particles were constructed to minimize the presence of geometric defects which may play a significant role in catalytic systems.
Jinnouchi and Asahi 39 have also shown that metal alloy nanoparticles can contain heterogeneous atomic configurations such as atomic-scale defects that cannot be explained by the single-crystal surface models, and that these defects dominate the catalytic activity of nanoparticles. They propose a machine-learning scheme which adopts a descriptor-based approach to map slab-based surface models to the nanoparticle surfaces and defects. The machine-learning model is trained on single crystal slabs with various defects and compositions. The model is then able to accurately predict energetics of complex nanoparticles, which indicates that the geometric defects dominate the nanoparticle behavior. However, it is unclear if a similar strategy will work for (mixed) metal oxide particles, particularly for photocatalysts where the electronic band structure is important.
Evaluating finite size e↵ects on band structure properties is complicated by the need for hybrid-level functionals. Exchange-correlation functionals based on the local-density approximation (LDA) and GGA significantly underestimate the bandgap because of a derivative discontinuity in the exchange-correlation potential. 40 However, the hybrid functionals that are implemented within the generalized Kohn-Sham (KS) formalism instead of standard KS formalism incorporate part of the discontinuity which leads to a bandgap that is in good agreement with experimental values. 41 Despite the fact that planewave codes are the most widely used method for solving the KS equations, it becomes impractical to perform planewave hybrid calculations on systems with > 100 atoms due to the early onset of the cubic scaling bottleneck and the large associated prefactor. In particular, the limited scalability on parallel computing platforms restricts the time to solution that can be achieved.
Another disadvantage of planewave codes is the nature of the Fourier basis that restricts the method to periodic boundary conditions, hence making it necessary to add artificial periodicity using vacuum and dipole corrections, which limits accuracy in the study of systems with dipole moments as well as charged systems. [42][43][44] The use of a finite di↵erence basis set can overcome these challenges by enabling ideal parallelization and application of non-periodic boundary conditions.
In this work, we investigate finite size e↵ects in rutile TiO 2 , which is the most stable polymorph of titania. 45 We utilize a collection of model systems including bulk, surface slabs, and nanoparticles to deconvolute e↵ects arising from various types of defects. We use these model systems to elucidate finite size e↵ects on the surface energy and band structure, and we explore the influence of the exchange correlation functional on these e↵ects by performing calculations with both semilocal (PBE) and hybrid (PBE0) functionals. The study utilizes the new finite-di↵erence Simulation Package for Ab-initio Real-space Calculations (SPARC) 42 code to enable hybrid-level calculations on nanoparticle systems. We also utilize a numerical convergence technique to assess the length scale of electronic interactions, and apply a geometric fingerprinting scheme based on machine-learning 46 to deconvolute geometrical and electronic finite size e↵ects 33 in surface energy. The results indicate that the surface energy converges quickly for both functionals, approaching the infinite particle limit at particle sizes of ⇠ 10 Å. However, we find that the band structure is highly sensitive to functional choice, particle size, and particle shape, indicating that hybrid-level calculations are required to assess the band structure of nanoparticles well beyond the sizes investigated here.
To study the finite size e↵ects in surface properties and band structure of nanoparticles, we construct the "TiO 2 model space" which is the collection of all model systems using the equilibrium structural parameters. This is illustrated in Fig. 1. For the bulk rutile system, we find equilibrium lattice constants of a = 4.636 Å and c = 2.967 Å. These lattice constants are in good agreement with the literature values and are summarized in Table 1 in the SI.
The slab models include a mix of symmetric and asymmetric slabs, and the facet energy of asymmetric slabs is calculated using the average between the two di↵erent facets. The cuboidal nanoparticles are non-periodic slabs, with vacuum in all directions instead of periodic boundary conditions. These cuboidal isolated nanoparticles are stoichiometric, and the faces of these cuboidal nanoparticles have an atomic structure consistent with the low-index facets of TiO 2 slabs. The faces of cuboidal particles will also have defect sites arising due to edges, corners and sub-edge atoms at edges (sub-edges) that are not present in the facets of slab models. These nanoparticles are a convenient way of studying the trends of surface properties and comparing with the extended slab models since their geometries are most similar to the slab models, making it easier to isolate electronic finite size e↵ects. However, the unphysical nature of the cuboidal geometries may cause artifacts when studying systemlevel properties such as the bandgap and DOS. Hence, non-stoichiometric nanoparticles that are closer to realistic systems are also created using the Wul↵-construction algorithm. 47 To facilitate comparison between the particle models, we restrict the Wul↵-constructed particles to only contain the same low-index facets that appear on the cuboidal particles (100, 010, and 001 facets).
Finite size e↵ects can be categorized as geometric e↵ects, electronic e↵ects, and quantum e↵ects. 33 In the case of energetic properties such as surface or adsorption energies, the geometric finite size e↵ects are expected to be local due to the nearsightedness of electrons. 48 Therefore, it is possible to deconvolute the geometric and electronic finite size e↵ects by partitioning the energy to specific types of geometric defects. Any deviation from this partitioning can be assumed to arise from electronic or quantum e↵ects that are implicit in the electronic structure of the particle. Here, we group quantum e↵ects with electronic structure e↵ects since they are both inherent to the behavior of electrons in the system.
To quantify the finite size e↵ects on the nanoparticle surface energy, we utilize the total surface energy of the nanoparticle obtained from DFT as the ground truth:
where E nanoparticle is the total energy of nanoparticle from DFT, E bulk is the bulk energy extracted from the linear interpolation model for slab surface energies (Eq. 5) and n surface Ti ,total
is the total number of TiO 2 units belonging to low-index facets and surface defects. In the absence of any finite size e↵ects, the nanoparticle surface energy could be approximated by a linear combination of semi-infinite slab energies:
where E facet,i is the facet energy and n facet,i is the number of surface-like TiO 2 units resembling each facet type. The semi-infinite slab energies are obtained by the linear interpolation method as described in the methods section and are referred to as "facet energies" to distinguish them from the energy of other types of surface defects. The converged slab energies are summarized in Table 1 and are shown in Fig. 2.
To deconvolute the influence of geometric finite size e↵ects, we utilize a regression model (Eq. 3) to obtain the energy contribution of edge, corner, and sub-edge geometric defects present within the system. The total surface energy of the particle based on the regression model can be obtained by adding the energy contributions of all surface defects that constitute of atoms resembling facet types and defects such as corners, edges and sub-edge defects and normalizing to the number of surface Ti atoms:
where the term P defect types j E defect,j n defect,j sums over the contributions of corner, edge, and sub-edge Ti atoms and accounts for the geometric finite size e↵ects. The results from the regression model are summarized in Table 1.
The trends in surface energy convergence for all three models are plotted in Fig. 3.
Fig. 3(a) clearly shows that the linear slab model poorly approximates the surface energy, which indicates that finite size e↵ects play a significant role at all particle sizes investigated. However, the regression model provides a much more accurate approximation of the nanoparticle surface energy, indicating that geometric finite size e↵ects are dominant. The residual between the regression model and the ground truth energy can be interpreted as the contribution due to electronic finite size e↵ects towards the surface energy. The electronic finite size e↵ects are quantified on a log scale in Fig. 3(b). At very small particle sizes (<4 Å) the electronic finite size e↵ects are very significant (⇠0.03 Ha/atom), which exceeds the typical exchange-correlation error (⇠0.01 Ha/atom). At particle sizes from ⇠4-8 Å (⇠30-200 atoms) the electronic finite size e↵ects are below the typical exchange-correlation error, but still exceed the numerical error (⇠0.001 Ha/atom). Finally, beyond ⇠10 Å the electronic finite size e↵ects are reliably at or below the threshold of numerical accuracy, suggesting that they can safely be neglected. This trend is consistent between both PBE and PBE0 functionals. This indicates that the surface properties of TiO 2 particles larger than ⇠ 10 Å can be studied using a combination of slab models, as long as slab models that capture the relevant geometric defects are included.
To validate the regression model, we apply it to particles not included in the regression analysis. These particles are generated by extending four of the cuboidal nanoparticles (24, 72, 108, and 162 atoms) in the 100, 010, and 001-directions and computing their surface energies with PBE. These particles will have the same types of geometric defects (edges, corners, and sub-edge atoms), but the proportions will be di↵erent from the original cuboidal particles. The electronic finite size e↵ects (the absolute di↵erence between E surface,DFT and E surface,regression ) as a function of particle size are shown in the SI (Fig. 1). Similar to Fig. 3, the residual decreases to ⇠1e-3 Ha/atom beyond ⇠ 10 Å, which confirms that the regression model can be generalized to systems that were not included in the training procedure. This also provides further support for the conclusion that the surface properties of nanoparticles with a minimum dimension greater than 10 Å can be accurately modeled using appropriate semi-infinite slabs.
The results of this work indicate that electronic finite size e↵ects have a negligible impact on the surface energetic properties of TiO 2 particles larger than ⇠10 Å. However, this conclusion is material dependent, with the characteristic length scale of finite size e↵ects depending on the localization of electrons in the material. TiO 2 is a large bandgap semiconductor, which suggests that electronic finite size e↵ects will decay quickly compared to more delocalized systems such as metals. Indeed, this is consistent with the fact that adsorption energies on metal nanoparticles tend to converge around ⇠2-3 nm. 33,37 Yet, a quantitative assessment of this length scale is extremely computationally demanding with the methods presented here and in prior work, since it requires the calculation of many large particles with DFT.
An alternative approach is to directly characterize the length scale that a↵ects the electronic convergence in the bulk solid. This is achieved by the "nearsightedness analysis" discussed in the methods section, which is based on the convergence of the density error as a function of the length scale of the interaction region around each point in the system. The results are shown in Fig. 4 for the bulk crystal and smallest nanoparticle. The box plots for the bulk crystal indicate that at a length scale of 10.6 Å the normalized electron density error decays to under 0.01 in all cases and below 0.001 in most cases. To correlate this to energy convergence, the normalized electron density error is compared to the total energy error at a variety of SCF tolerances for bulk TiO 2 and a linear fit is used to provide an upper bound based on the systems studied in this work (see SI Fig. 2). The results indicate that a normalized electron density error of 0.01 -0.001 corresponds to an SCF energy error of 10 4 10 6 Ha/atom. To ensure that these findings hold for the nanoparticle systems we also perform the nearsightedness analysis for the smallest nanoparticle, Ti 8 O 16 . The results indicate that the interactions decay as fast or faster in the nanoparticle system, with the error in SCF energy decaying to under 0.001 Ha/atom at a length scale of 8.0 Å in most cases.
The results of the nearsightedness analysis for both bulk and nanoparticle systems are consistent with the analysis of convergence of surface energy for nanoparticles, where we observe that the electronic finite size e↵ects disappear at a length scale >10.0 Å. The nearsightedness analysis and the corresponding energy error calibration can be performed using only the bulk system, and can be applied to any material. This suggests that the nearsightedness analysis is a far more computationally e cient route for evaluating the characteristic length scale of electronic finite size e↵ects for a given material.
The band structure of TiO 2 plays an important role in photocatalysis, especially the bandgap and defect states within the gap. To assess the finite size e↵ects on the band structure of the particles, we analyze the bandgap and DOS of stoichiometric TiO 2 nanoparticles and compare the results with bulk and slab models. The experimental bandgap for bulk rutile TiO 2 is ⇠ 3.0-3.1 eV, 45,49 while the bulk gap predicted by PBE and PBE0 are 1.89 and 4.19 eV, respectively. This is consistent with the well-known underestimation of band gaps by PBE, 40,50,51 and prior reports of over-estimation of TiO 2 band gaps with hybrid functionals. 52,53 Nonetheless, the hybrid results are considered to be more reliable since they accurately incorporate exchange interactions that lead to electron localization.
The bandgaps for all nanoparticles are shown in Fig. 5(a) as a function of particle size.
The plot indicates that bandgaps for all nanoparticles are much smaller than the bulk for both PBE and PBE0. Moreover, the bandgap values decrease, rather than increase as a function of particle size, indicating that the particle sizes investigated are far from the bulk limit. The PBE results suggest that the particles are metallic, indicating a qualitative failure of PBE. For this reason, PBE is omitted from subsequent analyses.
(a) (b) To understand the reason for smaller bandgaps in the nanoparticles, we analyze the full DOS of two of the cuboidal particles constituting of 162 and 288 atoms. We compare the DOS of these particles with that of the bulk crystal and slab models and these are plotted in Fig. 5(b). The DOS is normalized with respect to the total number of electrons in each of the particles. The vertical dotted lines highlight the bandgap of the bulk crystal. In the DOS of cuboidal particles, we observe that the bandgaps are much lower compared to the bulk bandgap and defect states within the bandgap are forming a continuum above the HOMO of the bulk crystal. This continuum leads to a very small gap of ⇠0.2 eV.
We compare the DOS of the cuboidal particles with the DOS of slab models. We observe that the (100)-slab model has a bandgap very close to the bulk and there are no defect states within the gap, whereas the DOS of the ( 010) and (001)-slabs reveal the defect states present within the gap resulting in a lower bandgap. The (010)-facet which is highly unstable also has a continuum of states within the gap, which suggests that the slab models provide some insight into the defect states within the gap for these cuboidal particles. However, these cuboidal particle models are unrealistic due to the presence of high energy surfaces in equal proportions. Therefore, we hypothesize that the nature of DOS of cuboidal particles is an artifact of the unphysical nature of the particles.
To test our hypothesis, we plot the DOS of Wul↵-constructed TiO2 particles which are more realistic, since the area occupied by each low-index facet is inversely proportional to the surface energy of the low-index facet. These particles have larger gaps with discrete defect states within the gaps. However, these gaps do not seem to converge to the bulk limit as the particle size increases and the location of the defect states cannot be predicted from the slab DOS plots. Since the minimum length scales of the Wul↵ particle with 136 and 282 atoms are 9.27 Å and 11.08 Å, this leads us to conclude that even at ⇠ 10 Å, the band structure is highly sensitive to the shape and size of the nanoparticles. Therefore, hybrid calculations with specific nanoparticle morphologies are required to understand their band structure.
We also visualize the orbitals around highest occupied molecular orbital (HOMO) and lowest unoccupied molecular orbital (LUMO) levels of the cuboidal and Wul↵-constructed particles which are included in the SI (Fig. 3). The orbitals are localized around ( 010) and (001) low-index facets of the cuboidal particles, whereas they appear to be much more delocalized throughout the Wul↵ particles. This behavior of orbitals supports the conclusion that presence of high energy facets results in defect states within the gap localized around those facets in cuboidal particles, whereas in the case of Wul↵ particles, the defects are delocalized throughout the nanoparticle, and will depend on the details of the particle size and shape.
In this work, we study the convergence of surface energy and band structure as a function of nanoparticle size for rutile TiO 2 . We utilize the PBE GGA functional and PBE0 hybrid functional, and evaluate cuboidal and Wul↵-constructed particles with a maximum size of 384 atoms (⇠13 Å). The results indicate that geometric finite size e↵ects play a significant role in the particle surface energy, but that they can be accounted for with a simple linear regression algorithm. In contrast, electronic finite size e↵ects on the surface energy decay rapidly, and can be neglected for particles with a minimum dimension of >10 Å. These findings are consistent for both the PBE and PBE0 functionals. Analysis of the nearsightedness of the electronic interactions provides further confirmation of the findings and suggests that the characteristic length scale of electronic interactions for a given material can be computed directly from the bulk electronic structure. This analysis may negate the need for expensive nanoparticle simulations to evaluate the length scale at which finite size e↵ects no longer impact energetic properties.
In the case of band structure, the findings indicate that the nanoparticle systems studied are well below the size at which the bandgap or band structure converges to the bulk limit.
The cuboidal particles are found to have a qualitatively di↵erent band structure from bulk TiO 2 , and the band structure does not change significantly with particle size. This can be explained by the band structure of the slab models, since the high-energy (010) surface slab exhibits a similar band structure to the cuboidal nanoparticles. The Wul↵ construction is used to generate more realistic particles without highly unstable surfaces, and their band structures are more similar to the bulk, exhibiting a large gap with more discrete defects within the gap. However, the size of the gap and the location of the defects varied considerably with particle size, and the bandgap decreased with increasing particle size. These results suggest that the bandgap and band structure of TiO 2 nanoparticles are highly sensitive to both the particle size and morphology at the length scales investigated here (<15.0 Å). This indicates that hybrid-level simulations of specific TiO 2 nanoparticle morphologies are required to elucidate their band structures.
The findings have implications for the electronic structure theory of nanoparticle systems and nanoparticle catalysis. The findings demonstrate that the finite-di↵erence SPARC DFT code is capable of performing hybrid-level DFT simulations for large (>350 atoms) systems, and that a nearsightedness analysis can be used to rapidly assess the characteristic length scale of electronic interactions. For catalysis, the findings show that the finite size e↵ects of surface properties (e.g. surface energy and adsorption energies) are dominated by geometric defects, and can be simulated using appropriate semi-infinite slab models. On the other hand, in the case of photocatalysis or other applications where the details of the band structure are important, explicit nanoparticle models are likely required to describe very small nanoparticles below ⇠ 20 Å. Further work is necessary to evaluate the band structure of more realistic particle morphologies and elucidate the role of solvents and adsorbates.
However, the emergence of highly parallelized hybrid DFT codes like SPARC, along with the increasing prevalence of petascale computing resources, represents promising progress toward evaluating these complex phenomena.
We perform Kohn-Sham DFT calculations using the state-of-the-art "Simulation Package for Ab-Initio Real-space Calculations" (SPARC) software 42,54,55 -a real-space DFT code that has comparable accuracy to established planewave codes, while requiring walltimes that are more than an order of magnitude lower. In all calculations, we neglect spin and choose optimized norm-conserving Vanderbilt (ONCV) pseudopotentials 56 from the SG15 57 collection. In addition, we employ the Perdew-Burke-Ernzerhof (PBE) 58 and PBE0 59,60 exchange-correlation functionals for performing GGA and hybrid level calculations, respectively.
In all simulations, we choose the twelfth-order finite-di↵erence approximation. For bulk calculations, where periodic boundary conditions are prescribed in all three coordinate directions, we employ a mesh-size of 0. For the nanoparticle calculations, where Dirichlet boundary conditions are employed in all three coordinate directions, we employ a mesh-size of 0.3 bohr and vacuum of 8 bohr in each coordinate direction. These and other parameters have been chosen to provide an accuracy of 0.001 Ha/atom in the energy. Note that we perform structural relaxation for only the bulk system, while the atoms in the slab systems and nanoparticles are held fixed, i.e., only the electronic ground state is computed for the given atomic positions. All atomic positions and structures are defined using the Atomic Simulations Environment (ASE) package. 62
We adopt the following decomposition for the nanoparticle energy:
where E bulk is the energy per TiO 2 unit for the crystal, n bulk is the number of bulk-like TiO 2 units in the nanoparticle, E facet,i is the average facet energy per TiO 2 unit for the i th slab type, n facet,i is the number of i th slab-like surface TiO 2 units in the nanoparticle, E defect,j is the defect energy per TiO 2 unit of the j th defect type, and n defect,j is the number of j th defect-like TiO 2 units in the nanoparticle. We have separated surface-like defects from other defects, since the energy of surface defects can be obtained from slab calculations. In the current work, we consider 3 slab types, i.e., (100)-symmetric, (010)-asymmetric, and (001)-symmetric, and 3 defect types, i.e., edge, corner, and sub-edge. The values for E defect,j are determined via least-squares regression, with the values for the remaining quantities computed using the methodology outlined below.
The energies E nanoparticle and E bulk are immediately available from DFT calculations for the nanoparticle and bulk, respectively. To calculate E surface,i , we adopt the extrapolation scheme of Fiorentini and Methfessel. 63 In particular, for each of the slab types, we first calculate the energy of the slab as a function of the number of layers N, which we denote by E slab,i (N). Next, we determine the average surface energy by fitting the data to the relation:
where Ẽbulk,i is the extrapolated bulk energy. In the current work, we have used N = 12 layers, which results in E surface,i values converged to within 0.0002 Ha/atom. Note that for the (010)-asymmetric slab, though the computed E surface,i is the average over the two di↵erent surfaces on either side of the slab, it can be used in our formulation since all nanoparticles considered in this work have both surfaces of the (010)-asymmetric slab, if present at all.
To be able to systematically determine n bulk , n slab,i , and n defect,i , we develop a machine learning scheme that determines the classification of each TiO 2 unit in the nanoparticles based on the local atomic environment. Specifically, we use the atomic descriptors proposed by Behler and Parrinello, 46 details of which can be found in the SI. We generate the atomic descriptors for fingerprinting the model space by using the Atomistic Machine-learning Package (AMP) 64 and select the hyperparameters of the featurization scheme such that there is maximum separation between di↵erent geometrical configurations in the model space. We apply dimensionality reduction using kernel principal component analysis (kPCA) 65 on the scaled descriptors to generate linearly independent features and use MeanShift clustering, 66 a density-based clustering algorithm to form clusters of atomic configurations for di↵erent types of titanium atoms in nanoparticles. The dimensionally reduced features facilitate visualization in a lower-dimensional space by using the interactive visualization tool Elec-troLens. 67 The scikit-learn software package 68 is used for all machine-learning models. We train the clustering algorithm with the model space and use ElectroLens to assign clusters to class labels for the categories of interest.
We perform the nearsightedness analysis for electronic interactions using the real-space Spectral Quadrature (SQ) method, 69-71 a technique developed for performing large-scale linear-scaling Kohn-Sham DFT calculations. In particular, for the electronic ground state computed using standard diagonalization-based schemes in SPARC, 42,54,55 we determine the convergence in electron density with size of interaction region -quadrature order to be large enough to make associated errors significantly smaller than those considered in this workanalogous to previous such results obtained for the energy and atomic forces in bulk aluminum at various temperatures. 72 Specifically, for a given interaction length scale, we express the electron density at any point in space as a bilinear form in terms of the Hamiltonian, which is then approximated by a Gauss quadrature rule that remains spatially localized to the interaction region by exploiting the locality of electronic interactions in real-space, 73 i.e, the exponential decay of the density matrix in real-space for insulators as well as metals at finite temperature. Indeed, in the limit of infinite size for the interaction region, the ground state electron density computed using diagonalization is exactly recovered.