Faraday cage screening reveals intrinsic aspects of the van der Waals attraction.

General properties of the recently observed screening of the van der Waals (vdW) attraction between a silica substrate and silica tip by insertion of graphene are predicted using basic theory and first-principles calculations. Results are then focused on possible practical applications, as well as an understanding of the nature of vdW attraction, considering recent discoveries showing it competing against covalent and ionic bonding. The traditional view of the vdW attraction as arising from pairwise-additive London dispersion forces is considered using Grimme's "D3" method, comparing results to those from Tkatchenko's more general many-body dispersion (MBD) approach, all interpreted in terms of Dobson's general dispersion framework. Encompassing the experimental results, MBD screening of the vdW force between two silica bilayers is shown to scale up to medium separations as 1.25 d e/d, where d is the bilayer separation and d e is its equilibrium value, depicting antiscreening approaching and inside d e Means of unifying this correlation effect with those included in modern density functionals are urgently required.

Materials in zero, one, two, and three dimensions of relevance to conceived future fabrication and electronics technologies are often held together by the van der Waals (vdW) dispersion force (1). Methods of measuring and first-principles simulations of the free energies of formation of such systems are becoming available (2). Often the critical issues involve situations in which the forces holding systems together become nonadditive, that is, the interaction between two parts of a system is modulated by the presence of nearby matter, with dispersion and other aspects all contributing (3, 4). At short distances typical of chemical bonding, it is now being recognized that dispersion forces can sometimes compete with traditional chemical covalent and ionic bonding forces to control outcomes by a number of means (5–13). Related parallel work demonstrates how ionic forces can control typical scenarios associated with dispersion (14), as well as scenarios in which general solvent effects including dispersion control structure (15, 16). Alternatively, at long distances, the Casimir dispersion effect becomes critical (17–19), as well as other exotic phenomena associated with the wavelike nature of charge polarization in nanoscale objects (4). While answers to each of the issues raised can be formed in isolation, a generally useful understanding of the vdW force remains elusive. Indeed, how different computational methods perceive dispersion forces at long and short distances have been found to be uncorrelated, raising fundamental questions concerning the nature of the force (12). What happens to long-range phenomena at vdW separations and then as chemical bond distances are reached will form a key part of future understanding.

To initiate discussion, we consider the extremely nonadditive vdW interactions observed by Tsoi et al. (20) in systems involving a silica substrate, a silica atomic force microscope (AFM) tip, and an intervening conducting graphene layer. The remarkable result from this work is that a large interobject dispersion force was switched off by the insertion of graphene in between the objects. Here, the origin and basic properties of this effect are elucidated using first-principles computational methods applied to a model 2D system. Discussion is considered using the framework for understanding vdW phenomena developed recently by Dobson (21).

Studies of 2D materials are currently very prevalent, with first-principle simulations providing powerful tools to facilitate understanding (11, 22–25). While dispersion interactions are critical for determining the structure and properties of such systems, the most commonly used method applied for first-principles materials simulations, density functional theory (DFT) using a conventional generalized-gradient approximation (GGA), improperly treats its contribution. As a result, a wide range of empirical correction schemes are commonly added to GGA calculations so as to produce a realistic description of the critical interactions (24, 25). The vdW dispersion interactions described by these schemes typically involve sums over interatomic interactions, each described by the London force (21, 26–28), with only small corrections. Related variants include replacing the atomic sums by electron-density integrals. In either case, these approaches are intrinsically pairwise additive, meaning that adding more atoms to the system just systematically increases the net dispersion interaction.

We select two widely applied types of methods to investigate the Tsoi et al. (20) experiment and its wider implications. First, we use pairwise-additive approaches based on “D3”-type schemes of Grimme that have achieved wide-ranging success (2, 29–33), particularly when applied to understand chemical vdW structures and energetics at equilibrium separations. These methods are also being found useful for understanding dispersion contributions down to chemical bond-length scales (10, 12, 31, 34). Second, we apply the more general many-body dispersion (MBD) approaches of Tkatchenko and others (4, 24, 35–45) that have also been shown to be widely successful, especially when large objects are involved at both equilibrium and larger separations. Of general interest is how different methods treat the vdW force during small excursions from equilibrium (46). More significantly, the observation that the vdW force can be switched off severely challenges the traditional pairwise-additive point of view and hence provides an external reference frame for this discussion. Here, combining the results from the different computational methods reveals that the switching off of the vdW interaction is an effect somewhat analogous to the “Faraday cage” effect that is well known throughout physics and engineering. The model system used for these calculations is shown in Fig. 1 and contains two bilayer silicon sheets separated by vacuum from an interposed graphene layer.

Fig. 1.

Geometries used for determining the screening of the vdW interaction between two silica bilayers induced by inserting an intermediary graphene sheet. (A) Two silica bilayers, (B) just one silica bilayer, (C) full system with graphene inserted midway between the two silica bilayers, (D) graphene only, and (E) graphene plus just a single silica bilayer.

Dobson’s (21) framework for understanding the dispersion force conceptually links traditional ideas concerning pairwise additivity to “top-down” Lifshitz theory in which screening emerges as a natural consequence of the quantum fluctuations that facilitate the dispersion interaction (47, 48). In the traditional view, each interatomic interaction is specified by the London contribution (49), which scales as r−6, where r is the interatomic distance, with only small corrections. The origin of the basic interaction is the polarization of one atom by spontaneous quantum dipole fluctuations occurring inside another.

At the simplest level, all atoms of the same type could be treated as having the same vdW attractions (the gas-phase values) independent of chemical environment. Dobson classified corrections to this approach into three types (21). The first type, here called “Dobson A,” focuses on the quantal effect of the local chemical environment, insofar as it modifies the polarizability of each atom (e.g., via orbital compression) and hence modifies its vdW interaction with other atoms. Just bringing two chemically noninteracting atoms close together induces a Dobson-A effect as the polarizability of one atom will change owing to the presence of the other. Most modern computational methods treat such terms to high accuracy, sometimes through empirical parameterization and sometimes through explicit environmental modeling.

In the second type, here called “Dobson B,” the fluctuating electric fields that mediate the vdW interaction between a pair of atoms are disrupted (screened or antiscreened) by the sympathetic induced fields caused by polarization of other atoms. This produces long-ranged “many-body” (many-atom) vdW effects (21, 27).

A third type of correction, called “Dobson C,” involves long-ranged fluctuating charge transfer occurring on a length scale larger than the size of an atom. This becomes particularly significant for long-ranged interactions between low-dimensional metallic conductors (26, 50, 51), but, in undoped graphene at T = 0 K, it is less important (22). It leads on (21) to the Casimir effect (17–19), something that can become very important in large systems (52).

Hence, in the context of the Tsoi et al. (20) experiment, Dobson-B effects take on a central focus: (i) What are they? (ii) how reliably do modern computational methods treat them? (iii) can they account for the observed screening of the silica–silica dispersion attraction by graphene? (iv) how do they influence critical equilibrium properties of vdW heterostructures? (v) how do they influence the way that dispersion forces are perceived at shorter distances of the order of chemical bond lengths? (vi) how do they behave asymptotically at long range and at very short range in the unified-atom limit? and (vii) can they be manipulated to make new functional materials and devices?

Results

Nature of the Dobson-B Screening Effect and Its Relationship to a Faraday Cage.

Fig. 2 illustrates the significance of type B (21) effects in the situation that we will model here, namely, a monolayer graphene sheet midway between two well-separated silica sheets. The graphene sheet will turn out to reduce the dispersion interaction between the two SiO2 sheets in a process analogous to the operation of a Faraday cage. Faraday caging arises within classical electrostatics, preventing an external electric field from penetrating through an infinite conducting sheet or a continuous metallic surface or “cage.” Typical ramifications of this effect include the blocking of radio waves by a steel-framed bridge and the reduction of mobile-phone signals by typically a million-fold by placing a MRI machine inside a room with continuous metallic walls. Polarization of the atoms within the metal by the external field cooperate across the whole extent of the metal to prevent penetration of the electric and magnetic fields. The effect applies independent of the source of the electric field, for example, permanent charges, electromagnetic radiation, or, as considered here, the spontaneous quantum dipole fluctuations involved in vdW attractions. Hence, screening induced by the Dobson-B effect can be thought of as a generalization of the abovementioned traditional notions of Faraday caging. However, vdW screening involves electronic fluctuations happening on the attosecond–femtosecond timescales, whereas Faraday caging formally involves responses at zero frequency, with microwave frequencies being so low in comparison that they are considered “zero frequency” in this context. Our application of the expression “Faraday cage” herein is also a little loose as its meaning is usually associated with continuous metallic conductors, whereas here we discuss atomistic graphene sheets. Certainly, the shielding of the fluctuating fields by an undoped, cold graphene sheet should be less complete than the screening of static fields by a true metal sheet, but it is still substantial, as will be demonstrated later. With these qualifications, the Faraday cage concept captures the essential essence of an important aspect of Dobson B, providing simple understanding to the Tsoi et al. (20) experiment.

Fig. 2.

MBD interactions. Electric field lines (solid) produced by a short-lived spontaneous dipole on a top-layer atom induce polarization of other atoms: These atoms therefore produce reaction fields (dashed or dotted). Weaker fields and dipoles are here shown as fainter. See the text for a detailed analysis of the diagrams. A full many-body treatment was required for the Faraday cage reported here.

As depicted in Fig. 2, Left, pairwise dispersion theories that neglect Dobson-B effects give the direct SiO2–SiO2 dispersion energy via interaction of the original top dipole with the dashed reaction field generated by the induced bottom dipole (53). This process misses simultaneous interactions involving both silica bilayers and graphene, and so cannot describe Faraday caging. The simplest way to include such interactions is through the three-atom Axilrod–Teller interaction (also known as the Axilrod–Teller–Muto interaction) terms (middle diagram in Fig. 2). These arise because the solid field lines from the top dipole also induce dipoles in the graphene layer. A graphene dipole then produces the dotted field line, which tends to cancel the solid direct field line, at the position of the atom in the lower SiO2 layer. This reduces the dipole induced on the bottom atom, and thereby reduces the dashed reaction field (fainter dashed lines). This in turn reduces the vdW interaction between atoms in the top and bottom layers. It is genuinely a three-atom effect, additional to the pairwise interactions of all species. The middle atom has “screened” the pair interaction between top and bottom atoms, yielding a partial Faraday cage effect.

Another way to view the three-atom Axilrod–Teller effect is the following. The polarizability of atoms in the combined subsystem (graphene plus bottom SiO2 layer) differs from the sum of the individual polarizabilities owing to the Coulomb interaction that occurs between the bottom two layers. This change of polarizability then affects the vdW interaction between atoms in the top Si layer and those in the combined (graphene plus lower SiO2) system.

However, this Axilrod–Teller term is only the simplest level of the Dobson-B electrostatic screening phenomenon as unique contributions also appear when four atoms are simultaneously considered, etc. Many-body terms (Fig. 2, Right) lead to a more complete description of Faraday cage screening. Here, any atom can screen (or enhance, “antiscreen”) the vdW interaction between any others, in accordance with the global geometry of the system. In this many-body case, for simplicity we do not distinguish initiating dipoles and fields from induced ones, showing all field lines as thin solid lines in the right-hand diagram of Fig. 2.

Inclusion of Dobson-B Effects in Electronic Structure Calculations of Molecules and Materials.

Among ab initio computational approaches, the simplest that includes dispersion is second-order Møller–Plesset theory (54) (MP2). This includes dispersion up to Dobson A only. Higher-order methods like coupled-cluster singles and doubles theory (CCSD) (55) and the random-phase approximation (RPA) (21) include also elements of both Dobson B and Dobson C.

In DFT approaches, most treatments of dispersion, including double-hybrid functionals, are based intrinsically on the assumption of pairwise additivity of the dispersion energy (21, 25), thus excluding explicit Dobson-B effects. These include the D3 family [except D3(ABC), which has some Axilrod–Teller terms] (29, 30), the XDM method (56), and many older approaches (57–60). All practical computational methods involve many assumed equations and their parameters, with some methods (like D3) being fully empirical with many parameters, while others are semiempirical, specifying equations that automatically generate required properties. Some methods may go beyond pure pairwise additivity by various means; D3(ABC) explicitly includes some three-body contributions, while most other methods include an uncontrolled amount of Dobson-B contributions into the Dobson-A terms (25) by modifying the pairwise-additive energies based only on the immediate chemical environment. D3 appears as an exemplary method demonstrating great success when screening and antiscreening effects can be treated empirically (2, 29–33), while D3(ABC) appears as the first step in a perturbative improvement strategy.

Going beyond these levels of treatment, of the efficient, semiempirical DFT+dispersion approaches, only the MBD method of Tkatchenko et al. (35, 41) and its descendants (24, 42) explicitly include Dobson-B corrections to all orders. It does so in a way that depends on global geometries and so includes Faraday screening and other many-body effects. This approach is known to be very successful in treating screening and antiscreening interactions at short range, that is, separations at and just beyond vdW contact (4, 24, 35–45). It is inappropriate for small-gap systems at very long ranges (21), but the separation used in the Tsoi et al. (20) experiment is just 20 nm, an intermediate distance at which MBD would still be expected to apply.

Very recently, a complementary interpretation of the Tsoi et al. (20) experiment was presented using a specialist method embodying Dobson B. In contrast to MBD-based approaches, this new approach was developed only to treat only very long-range interactions (17) and is inappropriate at vdW contact distances and thereabouts. It uses a model dielectric for graphene that does not involve the necessary cutoff to the intraband “Dirac cone” transitions and also neglects important interband transitions. Consequently, it is applicable only at very large distances (61).

An alternative method that encompasses the effects included in both MBD (35, 41) and this new work (17), as well as Dobson C and other effects, is RPA applied using DFT orbitals (62) (“dRPA”). This includes, at a realistic level, all dispersive effects operating over all length scales (53, 63). However, at the moment it is not computationally feasible to apply it to the problem at hand; nevertheless, efficient RPA methods are emerging (64–68), presenting a variety of possible future options.

Hence, at the moment, only MBD and its descendants appear as feasible approaches for capturing screening effects in the complex chemical environment of the silica–graphene system from the recognized important regime of close contact up to the separations of primary interest herein. Demonstration that this approach can account for a large reduction of the vdW attraction between silica bilayers on insertion of graphene, whereas methods that do not fully include Dobson B cannot, would present a result indicative that Faraday cage screening is responsible for the effect seen in the Tsoi et al. (20) experiments, connecting this effect to problems of great interest in biochemistry, chemistry, and condensed-matter physics.

Evidence for Faraday Cage Screening by Graphene.

The bilayer silica (20) used in our computational model constitutes a 2D layered material that has similar properties to usual 3D silica. However, it provides a test system small enough to allow calculations to very high precision (10 μeV), as is required for the evaluation of the dispersion force between bilayers over a large range in separations. A variety of specialized techniques are required to converge calculations to this precision, as detailed in . To further enhance precision, the geometries of the silica bilayers and the graphene sheet were frozen at their individually optimized values. This is a poor approximation at and inside vdW contact, but is very adequate in the critical intermediate-range region. Similarly, we choose a commensurate lattice with lattice parameter of 5.18 Å for silica and graphene. This preserves key qualitative properties of both materials, such as band structures, but ensures the cell size is reasonable; key results are insensitive to the chosen lattice parameter.

To determine the dispersion interaction between two silica bilayers and the effect of intervening graphene on it, calculations are performed for five systems, as described in Fig. 1. Fig. 1 shows two silica bilayers separated by interfacial distance d; its energy is labeled E. To get the unscreened dispersion energy ΔE at this distance, the energy of the corresponding isolated single bilayers, E, obtained from the structure shown in Fig. 1 is subtracted, yielding the following:To get the screened dispersion interaction ΔE, we insert a graphene layer in the center between two silica bilayers, as shown in Fig. 1. The energy of this structure is labeled E, from which must be subtracted the energies of each component and the energies of each silica–graphene dispersion interaction. To do this, we evaluate the energy of an isolated graphene sheet E (Fig. 1) and the energy of the bilayer–graphene interaction at distances d/2, E (Fig. 1), yielding the following:Of particular interest is the ratio specifying the reduction of the intersilica–bilayer dispersion force induced by the insertion of graphene. In addition, the total binding energy between the three layers is given by the following:The D3 and MBD computational approaches both have multiple variants. For D3, we select the methods knows as “D3(BJ)” and “D3(BJ,ABC),” both involving Becke–Johnson damping, and where ABC indicates that some three-body Axilrod–Teller terms are included (29, 69). For MBD (41, 70), we employ the fractional ionic (FI) method (43), which incorporates a superior polarizability model for oxides and a more physical treatment of many-body screening effects compared with the original method. Full details of the computational methods used are provided in (71–75).

Fig. 3 shows the total binding energy ΔE as a function of the interlayer spacing d, revealing that the predicted equilibrium separations are d = 6.65 Å both for MBD-FI and D3(BJ,ABC) or else 6.50 Å for D3(BJ), all in the range reported by Włodarczyk et al. (76). Next, Fig. 3 shows the unscreened vdW dispersion interaction ΔE. All three methods considered give similar results for small d, with D3(BJ) predicting ∼30% larger values than the other methods at equilibrium and hence also ∼30% larger binding ΔE. While the magnitude of the additional binding strength predicted by D3(BJ) is significant, the overall similarity of the results reflects the observation that all of these methods can describe most systems in vdW contact (2, 24, 25, 77). Most significantly, the results from the very different approaches D3(BJ,ABC) and MBD-FI shown in Fig. 3 (for ΔE) and Fig. 3 (for its small component ΔE) are mostly in excellent agreement.

Fig. 3.

Screening at work. (A) The total interaction energy ΔE at interlayer distance d (Fig. 1). (B) The contributions to this arising from the unscreened silica–silica vdW energies ΔE obtained without an intervening graphene. (C) The screening ratio ΔE/ΔE obtained after insertion of graphene. (D) The total dispersion energy for the BGB system, with the Inset showing that for the H2 molecule on a log–log scale. The green line indicates the equilibrium separation of de = 6.65 Å as calculated using MBD-FI and D3(BJ,ABC); D3(BJ) gives 6.60 Å instead.

However, how the various approaches treat the small bilayer–bilayer contribution to the total energy upon insertion of graphene, a quantity analogous to that measured in the Tsoi et al. (20) experiments, is highlighted in Fig. 3, where the ratio is shown. A purely pairwise-additive method would yield = 1 at all geometries, and the results for both D3 variants are indeed very close to this value at vdW contact and beyond. Small deviations occur through explicit inclusion of three-body corrections and the environment dependence of the parameters used in the method and, for larger distances, through rounding errors. In striking contrast, the MBD-FI ratio shown in Fig. 3 are fitted to yield the following:reducing toward zero with a d−1 dependence. The screening is half the equilibrium value at twice the equilibrium bond length, but the value of 1.25 at equilibrium indicates that antiscreening (factor > 1) is occurring at this critical geometry. One feature relevant here is that perfect Faraday cage screening only happens for continuous metallic conductors, and at short distances close to interatomic spacings, electric fields can penetrate through materials. Furthermore, the electron correlation effects manifested as the dispersion interaction must smoothly connect with the strong electron correlation effects manifested at short range within the GGA approximation, making the atomic nature of matter an unlikely cause of the antiscreening. In any case, a key feature of interest is that, at equilibrium distances, the net effect of Dobson-B terms becomes small (here, a 25% correction, with also a kink in the curve making the derivative small at this critical geometry), allowing sensibly parameterized pairwise-additive models to depict realistic many important features (2, 31, 32, 78).

However, the screening increases rapidly as the layers separate, reducing the unscreened vdW force to half at about twice the equilibrium separation. In this way, Dobson-B effects can become important even at or close to equilibrium separations, and approaches like MBD have been shown to have significant advantages in a range of applications (4, 24, 35, 37–43).

At moderately long range, extrapolation the MBD-FI results using Eq. predicts, at the separation of d = 20 nm as used in the Tsoi et al. (20) experiments, that the screening ratio should be 4.4%, close to the observed value of 5%. The calculations thus provide an explanation for the observation.

Anticipated Properties at Long Range.

The MBD-FI calculations and the Tsoi et al. (20) experiments consider screening in the 0.2- to 3-nm range and at 20 nm, respectively. Anticipated properties asymptotically at long range are also worthy of consideration. As an example of Faraday cage screening, we note that a 2D continuous metallic layer will completely reflect long-wavelength static (ω = 0) electric fields, but even for a metal, dynamic (high-frequency) fields are not fully reflected. The dispersion interaction depends on the exchange of such dynamic (high ω) electric fields. An undoped, zero-temperature graphene sheet is a good electrical conductor at zero frequency, but graphene is not a metal in the usual sense, having a zero electronic density of states at the Fermi energy. Correspondingly, it is known (68) that such a single undoped graphene layer, at T = 0 K as modeled herein, in contrast to a 2D metal, reflects most but not all of any long-wavelength static field, and this refection coefficient will be reduced at the finite frequencies of interest for the dispersion interaction. It is therefore not obvious how complete the Faraday caging effect will become for an undoped zero-temperature graphene sheet. The fitted d−1 dependence (Eq. ) found at short to intermediate distances may therefore in reality be replaced asymptotically by a constant minimum value. This is consistent with recent long-range calculations (17).

How the MBD method is expected to perform at very long range is also of interest. Rather than manifesting this anticipated theoretical limit, the MBD approach actually treats graphene as an insulator. Hence, asymptotically, MBD is expected to produce no screening at all (if it could be carried out numerically at asymptotic separations), and, like D3, would yield = 1.

Insight into the vdW Force and Its Means of Computation.

Most chemical discussion of the vdW force considers it as an example of the effects of electron correlation operating at distances above two to three times chemical bond lengths. Table 1 summarizes understanding of its properties as well as more “chemical” properties associated with smaller length scales. The simplest approach for describing electronic structure, Hartree–Fock theory (HF), treats electrons as effective noninteracting particles, ignoring all electron correlation. It is qualitatively descriptive of covalent and ionic bonding but omits the subtle features often responsible for quantitative analysis of chemical reaction property differences and, accordingly, does not include any vdW effects.

Table 1.

Overview of the effects of electron correlation

Example method	Correlation type	Effect at chemical bond lengths	Effect at near vdW separations	Long-range Dobson-B effects like screening and antiscreening
HF	None	Qualitatively descriptive of covalent and ionic bonds but poor quantitative accuracy	Not included	Not included
HF+D3	None	Akin to that from GGAs	Mostly quantitative	Not included
MP2	Fully pairwise	Quantitatively useful	Descriptive	Not included
CCSD, RPA	Fully many-body	Quantitative	Quantitative	Quantitative
GGA	Local many-body	Quantitatively useful	Improperly described	Improperly described
GGA+D3	Local many-body + pairwise vdW	Quantitativity improved by D3	Mostly quantitative	Not included
GGA+MBD	Local many-body + many-body vdW	Quantitativity reduced by MBD	Quantitative	Adequate at short to intermediate range*
Ref. 17	Long-range only	Inappropriate	Inappropriate	Designed for long range

Except for low-dimensional conductors for which Dobson C becomes critical, as may also be the related Casimir effect.

Its ab initio improvement, MP2 theory, includes both short-range and long-range correlation in a pairwise fashion, greatly improving chemical reaction properties as well as describing vdW bonding, but Dobson-B effects like screening and antiscreening are not included. The excluded effects are well known to influence chemical bonding at short range, and here we see that they are also can be profound at long range. This shows that Dobson-B effects are influential over all length scales. At the next level of ab initio improvement to HF theory, methods like CCSD and RPA (approximately) include all effects, providing bonding descriptions to chemical accuracy as well as including screening and antiscreening (21). While these methods do not directly provide insight in terms of intuitive descriptors like “chemical bonding” as distinct from “vdW bonding,” they embody a smooth link between these concepts.

DFT approaches like GGAs with vdW corrections are in much more widespread use owing to their computational efficiency and their ready applicability to both periodic systems and systems with small or no band gap, but the great challenge for the present concerns how to link together their disparate descriptions of chemical bonding and vdW bonding. GGAs provide an exact description of electron correlation in a free-electron gas over all length scales and so in principle fully include dispersion. However, they do so in a way that yields very poor results in any atomic system, defining the core problem needing to be solved (Table 1). Of note is that they embody a fully many-body description of short-range electron correlation that typically provides a better description of chemical properties than does a purely pairwise method like MP2.

So the situation is that many-body effects are important at short range and, as demonstrated by the Faraday cage screening effect, can also be critical at intermediate to long range. However, at the short-range distances typical of vdW equilibria, many-body effects often appear to cancel, allowing the combination of GGAs with pairwise-additive schemes like D3 to become widely successful (2, 31, 32, 78) (Table 1). Nevertheless, Table 1 indicates that replacing pairwise-additive methods like D3 with many-body ones like MBD allows all Dobson-B effects to be included at short to intermediate range, often leading to improved descriptions at vdW contact (4, 24, 35, 37–45). One would hope that such improvements continued down to the very short length scales of chemical bonding and beyond.

In Fig. 3, properties of the perceived dispersion interaction at short to intermediate separations are expanded from inside vdW contact (d ∼ de = 6.65 Å) to the very short distances typical of covalent and ionic bonding (d ∼ 4 Å if alignment allows) down toward the unified atom limit at which nuclear fusion would occur (d → 0). The MBD-FI screening curve follows the 1/d law well until near close contact (∼9 Å), undergoing some wiggles down to chemical bonding distances (5 Å) and then going flat (Fig. 3). In contrast, the D3 and B3(ABC) screening is near unity until ∼5 Å and then behaves very erratically. This arises as the London C6 coefficients are coordination dependent, allowing implicit inclusion of Dobson-B effects. However, the implementation of this concept in D3 and D3(ABC) is purely empirical. It is intended to describe typical chemical bonding scenarios, but simply compressing the BGB system together as is done in Fig. 3 does not correspond to any known chemical process and hence is not something that these methods were designed to do, leading to the erratic results. Our focus here is not in these troubles with the more empirical approaches, but rather the differences between them and MBD-FI in the region below vdW contact heading toward covalent bonding for which the D3 and D3(ABC) descriptions may be thought to be realistic: MBD-FI predicts strong antiscreening, an effect that flows from its dramatically improved description at intermediate distances exemplified by Faraday cage screening and other effects (46), while D3 and D3(ABC) do not.

A design criterion for D3 was that, using a predefined algorithm, it can be automatically applied to any density functional, as it has been done to at least 115 different density functionals (30, 31, 33). How the vdW correction links to the underlying description of electron correlation embodied in each functional has been a focus of much of the empiricism of the method. Significantly, it can also be applied to HF, to which it adds a description of the short-range electron correlation that can be competitive with descriptions provided by GGA-type density functionals such as Perdew–Burke–Ernzerhof (PBE) (30, 69, 79, 80). The central idea here is one of compatibility: The expectation is that adding D3 to any density functional largely just improves that functional’s description of covalent bonding.

In comparison, the development of MBD-based approaches are just in their infancy. MBD has been applied now to eight density functionals (45). Like D3, MBD assumes an empirical damping function specifying the changeover between the covalent and dispersive bonding regimes. However, for MBD, its form and parameter are chosen considering only data pertinent to dispersive interactions (45).

To bring focus to this discussion, we note that recently it was found that vdW forces can control the internal chemical structure of CuInP2S6 and CuBiP2Se6 monolayers (12), that specific Cu–Cu vdW interactions outcompete ionic forces to control antiferroelectric/ferroelectric arrangements of stacked monolayers (13), and that the nature of the vdW forces controlling stacking is uncorrelated with the nature of the vdW forces controlling internal chemical structure (12). Furthermore, D3-based approaches give results for the internal chemical structure that closely relate to experiment, being much improved on compared with the underlying density functional used alone, whereas MBD performs very poorly in this regard, degrading predictions.

In Fig. 3, the total dispersion energy (energy for method including dispersion less energy using only the underlying density functional) is shown for the BGB trilayer system. At short to intermediate distances in which Faraday cage screening is important, D3 and D3(ABC) overestimate the vdW attraction compared with MBD-FI. However, this situation reverses at very short distances comparable to covalent bond lengths with d < 5 Å. At short to intermediate distances, the D3 and D3(ABC) total dispersion energies scale as expected (21) as d−4, with MBD-FI scaling similarly as d−3.8. However, at extremely short distances with d < 4 Å, the MBD-FI scaling reducing to d−0.94, but the D3 and D3(ABC) interactions are being even more greatly dampened. This dampening of D3 dispersion energies compared with the continuing increase of MBD-type energies is a general effect, highlighted in the Inset in the figure, which shows the total dispersion energy for the H2 molecule from intermediate distances through vdW contact and then covalent contact to the unified atom limit. The vdW force must approach a constant value in the unified-atom limit (81), as exemplified only by D3.

Discussion

Our results show that the dramatic switching off of the vdW force demonstrated by the experiments of Tsoi et al. (20) is explained by Faraday cage-type screening effects involving the many-body response of all atoms in graphene to spontaneous quantum fluctuations in neighboring matter (Dobson-B effects). Modern computational methods can quantitatively reproduce these effects, provided that they include all responses to the fluctuations, as does MBD. However, subtleties, including correct asymptotic limits, needed to establish full quantitative agreement may require metallic (Dobson-C) models (68), this becoming critical whenever fully metallic low-dimensional systems are involved (21). Extension to include very long-range effects (17, 21) may also be required.

Methods like D3(BJ,ABC) that only partially include these responses are found to be inadequate for a general description of vdW screening in layered systems. These inadequacies are profound when screened silica bilayers are separated at distances of order twice the vdW separation and scale up with d−1 at larger (but still intermediate) distances. However, the screening effect cancels as the equilibrium position is reached, allowing simplistic pairwise-additive methods to often yield good results for this important class of problems (2, 29–33, 56–60).

We see here that the Tsoi et al. (20) experiment reveals fundamental physics critical to the understanding of systems at equilibrium as neglected screening (and other Dobson-B) effects appear as strong functions of separation around equilibrium structures. To better understand the dispersion force from asymptotic distances down to contact, there is great need for this experiment to be repeated so as to determine accurately the screening ratio as a function of distance d.

In terms of the conceptual understanding of the vdW force and its links to chemical forces, the results provide an understanding as to how Dobson-B terms behave at very short length scales below vdW equilibrium separations and their link with covalent bonding. Poor results found recently from MBD calculations for the competition between covalent and vdW forces controlling CuInP2S6 and CuBiP2Se6 monolayers (12, 13), poor results demonstrated for the dispersion force perceived in the H2 molecule in the unified-atom limit, may be related. These results therefore provide insights into how one might improve vdW models at atomic bond-length scales via simultaneous improvements to Dobson-A and -B effects, tailoring atom-scaled vdW forces. Clearly, modification of the design and parameterization of MBD approaches (45) to include these effects is highly desired.

More fundamentally, the results presented illustrate the competitive (and hard-to-define) role played by Dobson-A and -B terms at atomic scales; for example, the effect of the horizontal field lines in the right diagram of Fig. 2 can by mimicked by local changes to the response function. Can methods be developed that seamlessly include Faraday cage screening and covalent bonding? To make progress on this issue, an analytical understanding of how the RPA perceives Dobson-B effects at short range is urgently required. We feel that a numerical RPA correlation energy calculation at all separations, if feasible for the small energy differences investigated here, would be close to the current state of the art and would also very illuminating for these problems. Simplified versions of RPA (64–67) may help here.

In terms of nanotechnological applications, the consequences of Faraday cage screening may become significant in future devices made from 2D materials. Already, multiple layered systems are being made, allowing screening to occur (1). However, the screening effect acts to modify long-range interactions that are themselves small, so its influence may not be competitive with other effects in asserting control of self-assembly and device structure. This is in contrast to possible profound screening effects that conducting layers could have on functional properties housed inside layers that they separate, acting like a “ground plane” in electronics circuit boards.

To make the effects of vdW screening competitive, the objects being screened from each other need to themselves be large. The silica substrate and AFM tip used in the Tsoi et al. (20) experiment are examples of this effect. In this way, a conducting intermediary layer could dramatically reduce the force between two nanostructures of considerable thickness. If some external signal (e.g., transmitted current, applied electric or magnetic field, or an optical signal) could switch the conductance on or off, then the friction force between the nanoscale materials would also be modulated, blocking or facilitating slippage or separation.

The vdW force is not traditionally viewed as being controllable, as ionic forces, hydrogen bonding, and covalent bonding are. We now know that it can compete with these traditional chemical motifs to control structure and function (11). Developing the means of controlling and manipulating it therefore opens up new ways for designing functional materials and devices. The screening reported here shows that such manipulation is not only possible but can be modeled using existing theories, thereby enabling computer-led design of such materials.

Materials and Methods

We place all structures in a hexagonal periodic box with dimensions of 5.17 × 5.17 Å in the bilayer plane and 50 Å long, as indicated in Fig. 1. Initial density functional calculations are then performed using full 3D periodicity, as the long box length results in negligible PBE interactions between silica bilayer sheets. For the subsequent vdW dispersion corrections, however, only 2D periodic imaging is performed as these interactions remain significant on the 50-Å scale.

All calculations are performed by VASP (71, 72) using 24 × 24 × 1 k-points. The PBE density functional (73) is used along with PAW pseudopotentials describing core electrons (74). A plane-wave basis set is used set to an energy cutoff of 500 eV. The shortest-wavelength plane wave used in the direction perpendicular to the bilayers is 0.1 Å, and all calculated structures have layers translated by integral multiples of it. This ensures that the relationship between the grid points used in the electronic structure calculations and the precise nuclear coordinates remains invariant, seriously increasing the precision to which the calculations are performed.

In calculations of this type, the graphene lattice must be commensurate with that used for the silica bilayers. The two lattices are of the same type, but the PBE-calculated lattice constant is 5.18 Å for silica bilayers and 4.92 Å for graphene. In the calculations, a compromise value of 5.17 Å is used. The calculated band structures for the silica bilayer and for graphene at this value closely mimic those at the equilibrium geometries, manifesting the silica bilayer as a high bandgap insulator, as well as the graphene Dirac cones.

Dispersion correction options available in VASP allow for variations of the D3 and MBD methods. In all cases, we use PBE as the base functional. For D3, we select the methods known as “D3(BJ)” (VASP option “IVDW = 12”) and “D3(BJ,ABC)” (VASP option “IVDW = 12” and “noabc = .false.”) involving Becke–Johnson damping (both) and three-body Axilrod–Teller corrections (the latter only) (29, 69). For MBD (41, 70), we employ the FI method (43), which incorporates a superior polarizability model for oxides and uses eigenvalue rescaling to avoid the polarization catastrophe that can occur in highly polarizable systems (VASP option “IVDW = 263”) (42, 75, 82). Also, to enhance precision for the MBD calculations, we use large summation grid (VASP option “VDW_MBD_SIZE = 8 8 1”), which can fully capture the long-range screening.