Nanostructures can be bound together at equilibrium by the van der Waals (vdW) effect, a small but ubiquitous many-body attraction that presents challenges to density functional theory. How does the binding energy depend upon the size or number of atoms in one of a pair of identical nanostructures? To answer this question, we treat each nanostructure as a whole object, not as a collection of atoms. Our calculations start from an accurate static dipole polarizability for each considered nanostructure, and an accurate equilibrium center-to-center distance for the pair (the latter from experiment or from the vdW-DF-cx functional). We consider the competition in each term −C2k/d2k (k = 3, 4, 5) of the long-range vdW series for the interaction energy, between the size dependence of the vdW coefficient C2k and that of the 2kth power of the center-to-center distance d. The damping of these vdW terms can be negligible, but in any case, it does not affect the size dependence for a given term in the absence of non-vdW binding. To our surprise, the vdW energy can be size-independent for quasi-spherical nanoclusters bound to one another by vdW interaction, even with strong nonadditivity of the vdW coefficient, as demonstrated for fullerenes. We also show that, for low-dimensional systems, the vdW interaction yields the strongest size-dependence, in stark contrast to that of fullerenes. We illustrate this with parallel planar polycyclic aromatic hydrocarbons. The size dependences of other morphologies or bonding types lie between, as shown by sodium clusters.

Conventional Kohn-Sham density functional theory (DFT) has reached a high level of sophistication and achieved practical success, due to the good balance between efficiency and achievable accuracy. In recent years, many reliable semilocal density functionals have been proposed1–9 and some of them have been widely used in electronic structure calculations. However, these conventionally constructed DFT methods often produce large errors for molecular complexes and solids,10 interface problems,11 and ionic solids.12 A fundamental reason is that while conventional DFT methods can give an accurate description of the short-range part, the long-range vdW interaction is missing in these methods.

The long-range vdW interaction is an important nonlocal correlation due to instantaneous electric charge fluctuations. It affects many properties of molecular complexes and solids,13–21 layered materials,10,22,23 and ionic solids.12,24 It has been shown that performances of conventional DFT methods can be substantially improved with a vdW correction. A number of accurate vdW corrections have been developed. Most of them are based on atom pairwise effective models,25–29 while a few of them26,30–35 go beyond atom pairwise models. These models are very accurate for small or mid-size intermolecular interaction and have been widely used in electronic structure calculations, but errors may grow with system size or number of atoms in a system and can seriously affect the performances of vdW-corrected methods, as system size approaches the nanoscale. For example, it has been shown that while the dispersion-corrected atom pairwise model PBE+D2 is accurate for the binding of a pair of small molecules, the errors become large for fullerenes.21,36 It is known that an important source of errors is the nonadditivity of the vdW coefficients,37 due to electron delocalization and many-body (e.g., three-body) interactions.38,39

Additivity means that the multipole polarizability of a nanostructure scales linearly with N, the number of atoms in it, and that each vdW coefficient C2k between identical nanostructures scales as N2, the behavior predicted by atom pairwise interactions. As we will see, this expectation is rarely exact even for the dipole polarizability and for C6, and never for higher-order contributions. A proper treatment of non-additivity requires treating each nanostructure as a whole object, and not as a collection of atoms.

In solids, a particle at I not only interacts with another particle at J, but also with pairs JK of other particles, etc. The energy of vdW interaction of a particle at I with all others can be written as a many-body expansion39–43 

EvdW(I)=JEvdW(2)(IJ)+JKEvdW(3)(IJK)+.
(1)

where EvdW(2) is the two-body contribution, while EvdW(3) accounts for the three-body contribution.

A model for the two-body vdW interaction energy is usually developed from the asymptotic expansion of the vdW interaction at large separation, which can be written as

EvdW(2)=(C6/d6)fd,6(d/dvdW)(C8/d8)fd,8(d/dvdW)(C10/d10)fd,10(d/dvdW),
(2)

where d is the distance between the centers of two interacting density fragments and fd,2k is the damping function,44–48 with dvdW being the sum of vdW radii. In the development of vdW corrections, there are two important tasks. One is the calculation of vdW coefficients and the other is the design of a proper damping function, according to the short-range interaction from a semilocal DFT. The former involves important many-body effects and has received much attention. It was shown21 that in nanostructured materials such as fullerenes, fd,2k(d/dvdW ∼ 1) ≈ 1. The sum of vdW radii can define the minimum separation between the centers of density fragments for nanostructured materials without formation of a covalent bond. The center-to-center distance can be written as d = dvdW + Δ, with Δ being determined by the nature of the interaction. It can be positive, zero (i.e., direct contact), or negative. For nanostructures bound at equilibrium by only the vdW interaction, d/dvdW is of order unity, and nearly independent of system size (as we will show), so even when damping is important for the interaction energy, it is not important for the size-dependence of each term in Eq. (2).

The two-body vdW coefficients can be calculated from the dynamic multipole polarizability via the Casimir-Polder formula,49C2kAB=[(2k2)!/(2π)]l1=1k2[1/(2l1)!(2l2)!]0duαl1A(iu)αl2B(iu), where l2 = kl1 − 1. The required dynamic multipole polarizability can be modeled with the spherical-shell model. Since the electron density in nanostructures is nearly independent of system size or number of atoms in a system, we may simplify the spherical-shell model with the single-frequency approximation (SFA),50,51 in which we assume that (i) only valence electrons in the outermost subshell of an atom in a molecule are polarizable and (ii) the density is uniform inside the effective or vdW radius Rl and zero otherwise. This model is particularly useful for nanostructures or larger systems, in which the electron density is slowly varying.50–52 The only required input is the static polarizability which can be obtained from accurate ab initio methods.

Within the SFA, the model dynamic multipole polarizability takes the simple expression

αlSFA(iu)=Rl2l+1ωl2ωl2+u21ρl1βlρl,
(3)

where Rl is the effective outer radius of the shell defined below [Eq. (4)], βl=ωl2ω̃l2/[(ωl2+u2)(ω̃l2+u2)] describes the coupling of the sphere and cavity plasmon oscillations, and ρl=(1tl/Rl)2l+1 describes the shape of the shell, with tl being the shell thickness.50,53 In the static limit, the model dynamic multipole polarizability reduces to the true static polarizability, i.e., αlSFA(0)=αl0. For fullerenes, we set tl = 3.4 bohrs. ωl=ωpl/(2l+1) is the average sphere plasmon frequency, ω̃l=ωp(l+1)/(2l+1) is the cavity plasmon frequency, and ωp=4πn¯ is the average plasmon frequency of the extended electron gas, with n¯=N/Vl and Vl being the effective vdW volume and N being the total number of valence electrons in the outermost subshell (N=2 for carbon atom, while for B and N atoms, N=1, 3).

For a classical conducting shell of uniform density, the static multipole polarizability is related to the static dipole polarizability by αl(0)=[α1(0)](2l+1)/3. The vdW radius is defined by

Rl=[αl(0)]1/(2l+1),
(4)

with R1=R2=R3=[α1(0)]1/3. The preceding formula, which predicts the static higher-order polarizabilities from the static dipole polarizability, is valid for slowly varying densities, but only approximately true for strongly inhomogeneous densities such as atoms and ions.12 It can be shown that the average electron density n¯ within the shell is nearly a constant for nanostructures such as fullerenes. Therefore, we can write

EvdW(2)=α1A(0)α1B(0)f11(n¯A,n¯B)([α1A(0)]1/3+[α1B(0)]1/3+Δ)6α1A(0)α2B(0)f12(n¯A,n¯B)+P21([α1A(0)]1/3+[α1B(0)]1/3+Δ)8α1A(0)α3B(0)f13(n¯A,n¯B)+P22+P31([α1A(0)]1/3+[α1B(0)]1/3+Δ)10,
(5)

where fl1l2(n¯A,n¯B) represents the integral over the imaginary frequency iu, whose explicit expression can be extracted from Ref. 36. P21, P22, and P31 represent the terms containing α2A(0)α1B(0), α2A(0)α2B(0), and α3A(0)α1B(0), respectively. Note that Eq. (5) is valid for any two nanostructures, no matter whether they are identical or not. Consider a sequence of systems in which system size or number of atoms in a system increases from the initial (i) to final (f) size. The static multipole polarizability also changes from αli(0) to αlf(0)=Nf1+δl[αli(0)/Ni] for NfNi. Here δl is the nonadditivity of the static multipole polarizability measuring the deviation of the conventional value αlconv,f(0)=Nf[αli(0)/Ni] from the true value αlf(0). For example, δ1 = 0.2 for the dipole polarizability of fullerenes54,55 and −0.084 for sodium clusters evaluated from the dipole polarizability of Ref. 56. As seen below, even when δ1 = 0, δl > 0 for l ≥ 2.

In the preceding paragraph, δl is a measure of nonadditivity of the static multipole polarizability, as explained above. It is given by54δl=[(2l+1)(1+δ1)3]/3(1/3)[(2l+1)3](lnNi/lnNf), where the second term is a correction to the first term. For Ni=1 or NiNf, the second term vanishes. Then, we obtain δ2 = (2 + 5δ1)/3 and δ3 = (4 + 7δ1)/3. Substituting the expression for the nonadditivity of the static multipole polarizability into Eq. (5) and considering A = B leads to a dramatically simplified size-dependence of the vdW interaction energy between nanostructures,

EvdW(2)=[α1i(0)]2f1(n¯)(2[α1i(0)]1/3+Δ/R)6α1i(0)α2i(0)f2(n¯)(2[α1i(0)]1/3+Δ/R)8α1i(0)α3i(0)f3(n¯)(2[α1i(0)]1/3+Δ/R)10,
(6)

where R=[Nf(1+δ1)/Ni]1/3. This is the main result for the two-body vdW interaction. We can see from Eq. (6) that the size-dependence of the two-body vdW interaction depends upon the difference Δ between the true center-to-center distance and the sum of the vdW radii. When Δ is close to 0, the vdW attraction reduces to the value it takes initially for d=2[α1i(0)]1/3 and we are almost unable to observe the size dependence of the vdW interaction. When Δ is large and positive, the vdW attraction is rather weak. When Δ is negative, a stronger (e.g., covalent) bond may form between nanostructures. But in both cases, the size-dependence of the vdW interaction should be observed. We will demonstrate these observations with fullerenes,54 sodium clusters, and polycyclic aromatic hydrocarbons (PAHs) as follows.

In solids, a molecule or ion core not only interacts with another molecule or ion core but also simultaneously interacts with other two or more species. Three-body vdW interaction is much weaker12,38,39,42 than two-body vdW interaction. Here we only consider the most important lowest-order dipole-dipole-dipole interaction. The asymptotic form of the three-body interaction energy is given by EvdW(3)(IJK)=C9[3cos(θI)cos(θJ)cos(θK)+1]/(dIJdIKdJK)3, where dIJ = |IJ| are the sides of a triangle and θI, θJ, θK are the internal angles of the triangle formed by dIJdIKdJK. C9 is the three-body vdW coefficient, which can be calculated with the dipole polarizability from

C9=3π0duα1I(iu)α1J(iu)α1K(iu).
(7)

Substituting the dipole polarizability of Eq. (3) into the expression for the three-body vdW coefficient [Eq. (7)] and performing the integration over the imaginary frequency, we can obtain an expression for the three-body coefficient. Next we consider identical nanostructures. Making use of an analysis similar to that for the two-body interaction, we can easily obtain the final expression for the size dependence of three-body vdW interaction energy

EvdW(3)=[α1i(0)]3f(n¯,θI,θJ,θK)(2[α1i(0)]1/3+ΔIJ/R)×(2[α1i(0)]1/3+ΔIK/R)(2[α1i(0)]1/3+ΔJK/R)3,

which is similar to the two-body vdW interaction, with R being defined below Eq. (6). If ΔIJ, ΔIK, and ΔJK are all small, the three-body interaction is also size independent. Otherwise, it is size-dependent. This is very similar to the two-body interaction, but its effect on the total vdW interaction is small.

Summarizing, our size-dependence of the vdW interaction energy for identical nanostructures is formulated in terms of the hollow-sphere model within the SFA of Eq. (3). This model is valid for both spherical and non-spherical nanostructures of slowly varying densities because non-sphericity50,51 can enter the model via the input static polarizability αl(0) = (αl,xx + αl,yy + αl,zz)/3 obtained from accurate ab initio calculations. But Eq. (6) is only instructive for quasi-spherical nanostructures at equilibrium center-to-center distance between them. It is least instructive for finite low-dimensional parallel nanostructures, for which the center-to-center distance at equilibrium is nearly independent of system size. In this case, the size-dependence of the vdW interaction energy mainly arises from the size dependence of vdW coefficients. As such, we can employ the asymptotic formula of Eq. (2) directly to study its size dependence, as exemplified here with polycyclic aromatic hydrocarbons. For nanostructures of infinite length such as nanotubes, the situation is more complicated and will not be discussed.

Finally, it is worth pointing out that the vdW interaction energy can also display a size-dependence in the case of nonzero Δ, even when the vdW coefficients are additive. This suggests that these phenomena can be also described with an atom pairwise model.20,25,26

Fullerene is an important material with many remarkable properties such as great chemical stability and high sublimation or cohesive energy, leading to a variety of applications.57,58 Its properties resemble those of graphene in the large-size limit, such as zero-energy gap and binding energy.54 In fullerene solids, the vdW interaction is dominantly important, while covalent interaction between fullerene molecules is negligibly small.21 The shape of fullerenes is quasispherical and the electron density on the surface of fullerenes is nearly uniform. Therefore, they are ideal model systems for the study of the vdW interaction energy.

Table I for fullerene solids shows that both the leading-order and higher-order vdW interactions, −C2k/d2k, are nearly size independent. This is because Δ values are all small, compared to the sum of vdW radii dvdW=2[α1(0)]1/3. Table I also shows the comparison of expensive TDHF (time-dependent Hartree-Fock) calculations and model polarizability-based calculations for C6. From Table I, we can see the good agreement of our model C6 with TDHF values. The size independence of fullerene pair interactions may not be valid for other near neighbor (NN) pair interactions because for other NN pair interactions, Δ increases while dvdW is a constant. But the influence of other NN pair interactions is small, leading to the near size-independence of fullerene pair interaction. Table I also shows the very slow convergence of the vdW series for the interaction energy between two fullerenes at equilibrium, as anticipated in Refs. 21 and 59.

TABLE I.

Terms of the van der Waals interaction energy for a pair of identical fullerenes in a fullerene solid21 with the nearest neighbor center-to-center separation dcc obtained from fcc experimental lattice constants60,61 and for a pair of identical sodium clusters. The dipole polarizabilities of Na4 and Na8 are from Ref. 56, while that of Na19 is obtained as α1(0) = (R + δ)3, with R = N1/3rs and δ = 1.5 bohrs.62dcc is calculated as the distance between the centers of mass of two sodium clusters by putting them side by side63 [AA stacking for (Na4)2] with Quantum ESPRESSO64 using the vdW-DF-cx functional.31 All quantities are in atomic units. The reference values for C6 are from Ref. 55 for fullerenes and Ref. 56 for sodium clusters, except for Na19, which is taken as an average of Na18 and Na20. The calculated vdW coefficients50,52 are obtained from the hollow-sphere model within SFA of Eq. (3). A minus sign “−” in front of all the vdW interactions has been suppressed.

α1C6ref/103C6/103C8/105C10/108C9/105dccdvdWdccdvdWΔ103C6/d6103C8/d8103C10/d10(C9/d9)103
C60–C60 536.6 100.1 98.91 356.9 105.9 396.0 18.9 16.3 1.16 2.6 2.2 2.4 2.0 0.13 
C70–C70 659.1 141.6 144.7 601.8 205.7 711.1 20.1 17.4 1.16 2.7 2.2 2.5 2.2 0.13 
C78–C78 748.3 178.2 184.2 836.1 311.9 1027 20.8 18.2 1.14 2.6 2.2 2.5 2.2 0.14 
C84–C84 806.1 207.7 213.3 1019 400.2 1281 21.5 18.6 1.16 2.9 2.2 2.5 2.2 0.13 
Na4–Na4 511.5 16.80 17.28 51.40 13.75 57.89 13.0 16.0 0.81 −2.7 3.6 6.3 10 0.55 
Na8–Na8 883.9 52.48 55.68 251.3 98.47 342.5 15.7 19.2 0.82 −3.3 3.7 6.8 11 0.59 
Na19–Na19 1804 241.9 250.5 1941 1249 3389 16.5 24.3 0.68 −7.8 12.4 35.3 84 3.74 
α1C6ref/103C6/103C8/105C10/108C9/105dccdvdWdccdvdWΔ103C6/d6103C8/d8103C10/d10(C9/d9)103
C60–C60 536.6 100.1 98.91 356.9 105.9 396.0 18.9 16.3 1.16 2.6 2.2 2.4 2.0 0.13 
C70–C70 659.1 141.6 144.7 601.8 205.7 711.1 20.1 17.4 1.16 2.7 2.2 2.5 2.2 0.13 
C78–C78 748.3 178.2 184.2 836.1 311.9 1027 20.8 18.2 1.14 2.6 2.2 2.5 2.2 0.14 
C84–C84 806.1 207.7 213.3 1019 400.2 1281 21.5 18.6 1.16 2.9 2.2 2.5 2.2 0.13 
Na4–Na4 511.5 16.80 17.28 51.40 13.75 57.89 13.0 16.0 0.81 −2.7 3.6 6.3 10 0.55 
Na8–Na8 883.9 52.48 55.68 251.3 98.47 342.5 15.7 19.2 0.82 −3.3 3.7 6.8 11 0.59 
Na19–Na19 1804 241.9 250.5 1941 1249 3389 16.5 24.3 0.68 −7.8 12.4 35.3 84 3.74 

The same analysis will also apply to the three-body vdW interaction energy. To demonstrate the size independence of the three-body interaction energy, we first evaluate the three-body vdW coefficient in each fullerene solid with Eq. (7). Then we evaluate the three-body vdW interaction. The results are also displayed in Table I. From Table I, we can observe that, similar to the two-body vdW interaction, the three-body interaction is also size-independent, confirming our prediction.

Sodium cluster is a simple-metal cluster. It was shown that a sodium cluster can form a giant atom.65,66 Therefore, the pair interaction between sodium clusters can form a covalent bond. This leads to an intermolecular distance between centers of two sodium clusters significantly shorter than the sum of the vdW radii of sodium clusters, dvdW. To confirm this, we have calculated the distance between the centers of mass of two identical sodium clusters (Na4)2 (D2h), (Na8)2 (Td), and (Na19)2 (D5h) by putting them side by side, with Quantum ESPRESSO using the vdW-DF-cx functional,31 which has proven accurate in nanostructures.54 In this calculation, the energy cutoff is 30 hartree and only Γ points are included in the k-mesh. The results are displayed in Table I. From Table I, we can see that cancellation of the size dependence between vdW coefficients and the vdW radii is incomplete. We have also repeated the calculation of the distance between the centers of mass of sodium cluster pair (Na8)2 by putting them head to head and head to tail (see the supplementary material). The results are shorter, compared to that for side to side, suggesting that the cancellation of the size dependence between vdW coefficients and the vdW radii is incomplete in all the cases. In other words, the vdW interaction between sodium clusters is indeed size dependent. Table I also shows the comparison of the vdW coefficients obtained from the ab initio method and the hollow-sphere model within the SFA. From Table I, we see that the model calculation is in good agreement with more expensive ab initio values.

Finally, we have studied the size dependence of the vdW interaction energy between two identical PAHs with AA stacking. For such a geometry, a recent calculation with the vdW-DF-cx has shown54 that the plane-to-plane distance is only shrinking a little with system size, due to the fact that the vdW force also increases with system size, but the vdW coefficients per carbon atom increase rapidly with system size as we go from C6H6, C10H8, C14H10, to C18H12, owing to the nonadditivity. This leads to a rapid increase of the vdW interaction with system size. From Table II, we can observe that, even when the molecules are highly non-spherical, the model vdW coefficients [Eq. (3) with tl = Rl] are still accurate, compared to the more expensive time-dependent DFT-B3LYP calculations, suggesting the reliability of our results. In our calculations, the plane-to-plane distance of 2,3-benzanthracene54 was used for C18H12–C18H12.

TABLE II.

Terms of the van der Waals interaction energy between two identical polycyclic aromatic hydrocarbons (PAHs) of AA stacking with the nearest neighbor center-to-center separation d = dcc54 obtained from the vdW-corrected vdW-DF-cx nonlocal functional.31 All quantities are in atomic units (hartrees for energy, bohrs for distance). The vdW coefficients are evaluated from the hollow-sphere model within the SFA of Eq. (3) with tl = Rl, where the average valence electron density n¯=N/v.51 The reference values of C6 are from TDDFT (time-dependent DFT) calculations.67,68 The input static dipole polarizabilities are taken from TDHF (time-dependent Hartree-Fock) calculations.68 The higher-order multipole polarizabilities are estimated from the conventional formula of Eq. (4). We have used a density n¯=0.0468. A minus sign “−” in front of all the vdW interactions has been suppressed.

α1C6ref/102C6/102C8/104C10/106dccdccdcc(C6H6C6H6)(C6/d6)102(C8/d8)102(C10/d10)102
C6H6–C6H6 68.23 17.73 17.93 15.65 11.35 7.69 1.00 0.87 1.28 1.57 
C10H8–C10H8 117.3 48.67 50.42 63.16 65.71 7.47 0.98 2.90 6.51 12.1 
C14H10–C14H10 176.6 100.3 108.5 178.5 244.0 7.37 0.97 6.77 20.5 51.6 
C18H12–C18H12 244.8 175.1 199.0 407.1 691.8 7.30 0.96 13.1 50.5 161.0 
α1C6ref/102C6/102C8/104C10/106dccdccdcc(C6H6C6H6)(C6/d6)102(C8/d8)102(C10/d10)102
C6H6–C6H6 68.23 17.73 17.93 15.65 11.35 7.69 1.00 0.87 1.28 1.57 
C10H8–C10H8 117.3 48.67 50.42 63.16 65.71 7.47 0.98 2.90 6.51 12.1 
C14H10–C14H10 176.6 100.3 108.5 178.5 244.0 7.37 0.97 6.77 20.5 51.6 
C18H12–C18H12 244.8 175.1 199.0 407.1 691.8 7.30 0.96 13.1 50.5 161.0 

To see the damping effect on the vdW interaction energy, we calculate the ratio of the center-to-center distance over the sum of the vdW radii dccdvdW and the center-to-center distance of PAHs over the benzene-to-benzene distance dccdcc(C6H6C6H6). The former are all displayed in Table I, while the latter are displayed in Table II. From Table I, we see that all the ratios are nearly independent of system size, except for Na19–Na19. This means that the damping function should not change the size-dependence of the vdW interaction energy. The reason for this is the cancellation of the size-dependences of the equilibrium center-to-center distance and vdW radii. From Table II, we can also see that the ratios are nearly the same for PAHs because the plane-to-plane distance of PAHs should not be size-dependent.

Figure 1 shows the comparison of the dipole-dipole interaction energy term −C6/d6 for identical fullerene pairs, sodium cluster pairs, and PAH pairs. From Fig. 1, we observe that the size dependence of fullerene pairs is nearly a constant, while that of PAH pairs yields the strongest size dependence. The size-dependence of the vdW interaction between sodium clusters is between these two extreme cases. From Tables I and II, we also see similar size dependences for higher-order interactions.

FIG. 1.

Comparison of the size dependences for the absolute value of the vdW interaction energy term −C6/d6 between identical fullerenes in fcc solids (upper left panel), sodium clusters (upper right panel), and PAHs (lower panel).

FIG. 1.

Comparison of the size dependences for the absolute value of the vdW interaction energy term −C6/d6 between identical fullerenes in fcc solids (upper left panel), sodium clusters (upper right panel), and PAHs (lower panel).

Close modal

In conclusion, we employ the model polarizability and the experimental center-to-center distance (for fullerenes) or the distance from the vdW-DF-cx functional (for sodium clusters and PAHs) to study the size dependence of the vdW energy at the equilibrium distance. The former offers a good description of vdW coefficients via the Casimir-Polder formula, while the latter is accurate in the prediction of the center-to-center equilibrium distance, but inaccurate in the vdW coefficients for some nanostructures such as fullerene.36,54 The dependence of the vdW interaction energy on system size or number of atoms is a common feature for nanostructures. It arises from the competition between the size-dependences of the vdW coefficients and the center-to-center distance. In this work, starting from the asymptotic long-range vdW interaction, we have derived an expression for the size-dependence of the vdW interaction, which is valid for identical nanostructures. We have studied the size-dependence of the vdW interaction for fullerenes, sodium clusters, and PAHs. Our calculations show that for two identical nearest neighbor fullerenes in a fullerene solid, the vdW interaction is size-independent. This is unexpected, given that the vdW coefficients of fullerenes have very strong nonadditivity50 or non-linear effects, due to the electron delocalization. However, for low-dimensional nanostructures, the vdW interaction shows the strongest size-dependence. We illustrate this with planar PAHs. For sodium clusters, the size-dependence of the vdW interaction is between those of fullerenes and PAHs.

See supplementary material for the details of different molecular geometries of sodium cluster.

We thank Andrey Solov’yov for kindly sending us Ref. 65 with the atomic coordinates of sodium clusters. J.T. was supported by the NSF under Grant No. CHE 1640584. J.T. was also supported on Temple start-up from J.P.P. J.P.P. and C.S. were supported by the NSF under Grant No. DMR-1607868. H.T. acknowledges support from the DOE Office of Science, Basic Energy Sciences (BES), under Grant No. DE-SC0018194. Computational support was provided by HPC of Temple University and NERSC.

1.
J. P.
Perdew
,
K.
Burke
, and
M.
Ernzerhof
,
Phys. Rev. Lett.
77
,
3865
(
1996
).
2.
J.
Tao
,
J. P.
Perdew
,
V. N.
Staroverov
, and
G. E.
Scuseria
,
Phys. Rev. Lett.
91
,
146401
(
2003
).
3.
Y.
Zhao
and
D. G.
Truhlar
,
J. Chem. Phys.
125
,
194101
(
2006
).
4.
J.
Sun
,
A.
Ruzsinszky
, and
J. P.
Perdew
,
Phys. Rev. Lett.
115
,
036402
(
2015
).
5.
J.
Sun
,
R. C.
Remsing
,
Y.
Zhang
,
Z.
Sun
,
A.
Ruzsinszky
,
H.
Peng
,
Z.
Yang
,
A.
Paul
,
U.
Waghmare
,
X.
Wu
,
M. L.
Klein
, and
J. P.
Perdew
,
Nat. Chem.
8
,
831
(
2016
).
6.
J.
Tao
and
Y.
Mo
,
Phys. Rev. Lett.
117
,
073001
(
2016
).
7.
G.
Tian
,
Y.
Mo
, and
J.
Tao
,
J. Chem. Phys.
146
,
234102
(
2017
).
8.
M.
Ernzerhof
and
G. E.
Scuseria
,
J. Chem. Phys.
110
,
5029
(
1999
).
9.
A. D.
Becke
,
J. Chem. Phys.
104
,
1040
(
1996
).
10.
H.
Peng
,
Z.-H.
Yang
,
J. P.
Perdew
, and
J.
Sun
,
Phys. Rev. X
6
,
041005
(
2016
).
11.
J.
Tao
and
A. M.
Rappe
,
Phys. Rev. Lett.
112
,
106101
(
2014
).
12.
J.
Tao
,
F.
Zheng
,
J.
Gebhardt
,
J. P.
Perdew
, and
A. M.
Rappe
,
Phys. Rev. Mater.
1
,
020802(R)
(
2017
).
13.
D. A.
Egger
and
L.
Kronik
,
J. Phys. Chem. Lett.
5
,
2728
(
2014
).
14.
W. A.
Al-Saidi
,
V. K.
Voora
, and
K. D.
Jordan
,
J. Chem. Theory Comput.
8
,
1503
(
2012
).
15.
O.
Hod
,
J. Chem. Theory Comput.
8
,
1360
(
2012
).
16.
A.
Otero-de-la-Roza
and
E. R.
Johnson
,
J. Chem. Phys.
137
,
054103
(
2012
).
17.
A.
Tkatchenko
,
R. A.
DiStasio
,
R.
Car
, and
M.
Scheffler
,
Phys. Rev. Lett.
108
,
236402
(
2012
).
18.
T.
Risthaus
and
S.
Grimme
,
J. Chem. Theory Comput.
9
,
1580
(
2013
).
19.
J. G.
Brandenburg
and
S.
Grimme
,
Top. Curr. Chem.
345
,
1
24
(
2014
).
20.
A.
Otero-de-la-Roza
and
E. R.
Johnson
,
J. Chem. Theory Comput.
11
,
4033
(
2015
).
21.
J.
Tao
,
J.
Yang
, and
A. M.
Rappe
,
J. Chem. Phys.
142
,
164302
(
2015
).
22.
B.
Sachs
,
T. O.
Wehling
,
M. I.
Katsnelson
, and
A. I.
Lichtenstein
,
Phys. Rev. B
84
,
195414
(
2011
).
23.
T.
Bučko
,
S.
Lebègue
,
J.
Hafner
, and
J. G.
Ángyán
,
Phys. Rev. B
87
,
064110
(
2013
).
24.
T.
Gould
,
S.
Lebègue
,
J. G.
Ángyán
, and
T.
Bućko
,
J. Chem. Theory Comput.
12
,
5920
(
2016
).
25.
A.
Tkatchenko
and
M.
Scheffler
,
Phys. Rev. Lett.
102
,
073005
(
2009
).
26.
S.
Grimme
,
J.
Antony
,
S.
Ehrlich
, and
H.
Krieg
,
J. Chem. Phys.
132
,
154104
(
2010
).
27.
A. D.
Becke
and
E. R.
Johnson
,
J. Chem. Phys.
127
,
154108
(
2007
).
28.
S. N.
Steinmann
and
C.
Corminboeuf
,
J. Chem. Phys.
134
,
044117
(
2011
).
29.
S. N.
Steinmann
and
C.
Corminboeuf
,
J. Chem. Theory Comput.
7
,
3567
(
2011
).
30.
M.
Dion
,
H.
Rydberg
,
E.
Schröder
,
D. C.
Langreth
, and
B. I.
Lundqvist
,
Phys. Rev. Lett.
92
,
246401
(
2004
).
31.
K.
Berland
and
P.
Hyldgaard
,
Phys. Rev. B
89
,
035412
(
2014
).
32.
K.
Berland
,
V. R.
Cooper
,
K.
Lee
,
E.
Schröder
,
T.
Thonhauser
,
P.
Hyldgaard
, and
B. I.
Lundqvist
,
Rep. Prog. Phys.
78
,
066501
(
2015
).
33.
O. A.
Vydrov
and
T.
Van Voorhis
,
Phys. Rev. Lett.
103
,
063004
(
2009
).
34.
J.
Granatier
,
P.
Lazar
,
M.
Otyepka
, and
P.
Hobza
,
J. Chem. Theory Comput.
7
,
3743
(
2011
).
35.
J.
Klimes
,
D. R.
Bowler
, and
A.
Michaelides
,
Phys. Rev. B
83
,
195131
(
2011
).
36.
A.
Ruzsinszky
,
J. P.
Perdew
,
J.
Tao
,
G. I.
Csonka
, and
J. M.
Pitarke
,
Phys. Rev. Lett.
109
,
233203
(
2012
).
37.
J. F.
Dobson
,
A.
White
, and
A.
Rubio
,
Phys. Rev. Lett.
96
,
073201
(
2006
).
38.
M. L.
Klein
and
R. J.
Munn
,
J. Chem. Phys.
47
,
1035
(
1967
).
39.
O. A.
von Lilienfeld
and
A.
Tkatchenko
,
J. Chem. Phys.
132
,
234109
(
2010
).
40.
R. J.
Bell
and
A. E.
Kingston
,
Proc. Phys. Soc.
88
,
901
(
1966
).
41.
42.
M. B.
Doran
and
I. J.
Zucker
,
J. Phys. C: Solid State Phys.
4
,
307
(
1971
).
44.
J.-D.
Chai
and
M.
Head-Gordon
,
Phys. Chem. Chem. Phys.
10
,
6615
(
2000
).
45.
K. T.
Tang
and
J. P.
Toennies
,
J. Chem. Phys.
80
,
3726
(
1984
).
46.
Q.
Wu
and
W.
Yang
,
J. Chem. Phys.
116
,
515
(
2002
).
47.
E. R.
Johnson
and
A. D.
Becke
,
J. Chem. Phys.
124
,
174104
(
2006
).
48.
S.
Grimme
,
S.
Ehrlich
, and
L.
Goerigk
,
J. Comput. Chem.
32
,
1456
(
2011
).
49.
S. H.
Patil
and
K. T.
Tang
,
J. Chem. Phys.
106
,
2298
(
1997
).
50.
J.
Tao
and
J. P.
Perdew
,
J. Chem. Phys.
141
,
141101
(
2014
).
51.
J.
Tao
and
A. M.
Rappe
,
J. Chem. Phys.
144
,
031102
(
2016
).
52.
J.
Tao
,
Y.
Mo
,
G.
Tian
, and
A.
Ruzsinszky
,
Phys. Rev. B
94
,
085126
(
2016
).
53.
G. K.
Gueorguiev
,
J. M.
Pacheco
, and
D.
Tománek
,
Phys. Rev. Lett.
92
,
215501
(
2004
).
54.
J.
Tao
,
Y.
Jiao
,
Y.
Mo
,
Z.-H.
Yang
,
J.-X.
Zhu
,
P.
Hyldgaard
, and
J. P.
Perdew
, “
First-principles study of the binding energy between nanostructures and its scaling with system size
,”
Phys. Rev. B
(submitted).
55.
J.
Kauczor
,
P.
Norman
, and
W. A.
Saidi
,
J. Chem. Phys.
138
,
114107
(
2013
).
56.
A.
Jiemchooroj
,
P.
Norman
, and
B. E.
Sernelius
,
J. Chem. Phys.
125
,
124306
(
2006
).
57.
R.
Loutfy
and
E.
Wexler
, “
Ablative and flame-retardant properties of fullerenes
,” in
Perspectives of Fullerene Nanotechnology
, edited by
E.
Ōsawa
(
Kluwer Academic Publishers
,
New York
,
2002
), Part VII, pp.
275
280
.
58.
Nanoclusters, Volume 1: A Bridge across Disciplines
, edited by
P.
Jena
and
A. W.
Castleman
, Jr.
(
Elsevier
,
Amsterdam
,
2010
).
59.
J. P.
Perdew
,
A.
Ruzsinszky
,
J.
Sun
,
S.
Glindmeyer
, and
G. I.
Csonka
,
Phys. Rev. A
86
,
062714
(
2012
).
60.
M. N.
Magomedov
,
High Temp.
43
,
379
(
2005
).
61.
R. D.
Bendale
and
M. C.
Zerner
,
J. Phys. Chem.
99
,
13830
(
1995
).
62.
63.
I. A.
Solov’yov
,
A. V.
Solov’yov
, and
W.
Greiner
,
Phys. Rev. A
65
,
053203
(
2002
).
64.
P.
Giannozzi
 et al,
J. Phys.: Condens. Matter
21
,
395502
(
2009
).
65.
S.
Saito
and
S.
Ohnishi
,
Phys. Rev. Lett.
59
,
190
(
1987
).
66.
Y.
Ishii
,
S.
Saito
, and
S.
Ohnishi
,
Z. Phys. D: At., Mol. Clusters
7
,
289
(
1987
).
67.
M. A. L.
Marques
,
A.
Castro
,
G.
Malloci
,
G.
Mulas
, and
S.
Botti
,
J. Chem. Phys.
127
,
014107
(
2007
).
68.
A.
Jiemchooroj
,
P.
Norman
, and
B. E.
Sernelius
,
J. Chem. Phys.
123
,
124312
(
2005
).

Supplementary Material