Locality of Interatomic Interactions in Self-Consistent Tight Binding Models

Abstract

A key starting assumption in many classical interatomic potential models for materials is a site energy decomposition of the potential energy surface into contributions that only depend on a small neighbourhood. Under a natural stability condition, we construct such a spatial decomposition for self-consistent tight binding models, extending recent results for linear tight binding models to the nonlinear setting.

1 Introduction

Electronic structure models are widely used to calculate many optical, magnetic and mechanical properties of materials and molecules (Finnis 2003; Martin 2004; Kittel 2004; Szabo and Ostlund 1996). Self-consistent (nonlinear) tight binding models are simple electronic structure models that are interesting in their own right but also provide convenient prototypes for the much more complicated density functional theories. A paradigm example is the density functional tight binding (DFTB) method (Koskinen and Mäkinen 2009; Seifert and Joswig 2012; Elstner and Seifert 2014).

In contrast to previous works on the linear tight binding model (Chen et al. 2018; Chen and Ortner 2016; Chen et al. 2019), the present work allows for the “on-site” Hamiltonian matrix entries to depend on the local electron density through a potential term. This results in a nonlinear model since the electron density itself depends on the Hamiltonian via a self-consistency condition. The present work therefore represents a key stepping stone between linear tight binding and more accurate nonlinear models such as Kohn–Sham density functional theory.

The self-consistency introduces the interesting issue of stability of the electronic structure problem; the main additional technical difficulty when compared to the linear model. Therefore, under a suitable stability condition (E and Lu 2010), we show that the potential energy surface in this model can be decomposed into exponentially localised site contributions, thus justifying many classical interatomic potential (IP) and multi-scale methods.

Despite the relative simplicity of tight binding models, a naive implementation demands $O(N^3)$ computational cost, where N is the number of particles in the system. It may therefore be advantageous to instead implement an IP model. In this case, the relevant parameters can be fitted to high accuracy by machine learning techniques together with theoretical data resulting from a high-fidelity model (Bartók et al. 2018, 2010; Behler and Parrinello 2007; Shapeev 2016). In most IP models for materials, a necessary starting assumption is that the potential energy surface can be decomposed into localised site contributions. That is, for atomic positions $y = \{y_n\}$, the total energy $E = E(y)$ may be written

$$\begin{aligned} E(y) = \sum _\ell E_\ell (y) \quad \text {with} \quad \left| \frac{\partial ^j E_\ell (y)}{\partial y_{n_1}\dots \partial y_{n_j}}\right| \lesssim e^{-\eta \sum _{l=1}^j |y_\ell - y_{n_l}|}, \end{aligned}$$

(1.1)

for some $\eta > 0$.

Typically, classical IP models are short-ranged and thus only justified if the exponent $\eta$ in (1.1), which measures the interatomic interaction range, is not too small. Therefore, in the present paper, not only are we interested in obtaining a site energy decomposition as in (1.1), but we also wish to describe the exponent.

Also, in the context of QM/MM multi-scale methods, a minimal starting assumption on the QM model is the strong locality of forces (Csányi et al. 2005). This condition is implied by the strong energy locality condition (1.1) that we consider in the present paper. Moreover, the exponent $\eta$ from (1.1) gives a guide for the size of the (computationally more expensive) QM region that must be imposed (Chen and Ortner 2016, 2017; Csányi et al. 2005).

In addition to the partial justification for IP and multi-scale models, our results also allow the thermodynamic and zero Fermi-temperature limit results of Ortner and Thomas (2020), Chen et al. (2018) to be extended to the nonlinear setting. We sketch the main ideas in the concluding remarks in Sect. 5.

1.1 List of Assumptions

The analysis requires a number of assumptions on the tight binding model, the atomistic configuration and associated electronic densities.

In general, we require the following assumptions:

(L):

p. 4 Uniform non-interpenetration with non-interpenetration constant $\mathfrak {m}>0$

(TB):

p. 5 Nonlinear tight binding Hamiltonian $\mathcal {H}(u;\rho )$

(SC):

p. 7 Self-consistent electronic density $\rho$

(STAB):

p. 8 Stability of the configuration with stability operator $\mathscr {L}(u;\rho )$

In addition to the above, we require the following assumptions for the particular case of point defects in insulating multi-lattice materials:

(REF):

p. 9 Multi-lattice reference configuration $\varLambda ^{\mathrm {ref}}$ with reference density $\rho ^{\mathrm {ref}}$

(GAP):

p. 10 Insulating material (band gap in the system) with gap $\mathsf {g}^{\mathrm {ref}}>0$

$\mathbf (STAB ^\mathbf{ref }{} \mathbf )$:

p. 11 Stability of the reference configuration with reference stability operator $\mathscr {L}^{\mathrm {ref}}$

(P):

p. 11 Point defect reference configuration $\varLambda$

(FF):

p. 11 Far-field assumption on the electronic density

1.2 Summary of Results

The results of this paper are divided into two sections: we discuss general locality estimates in Sect. 2 and improve these results for the specific case of point defects in insulating multi-lattice materials in Sect. 3.

1.2.1 General Locality Estimates

The previous works (Chen et al. 2018; Chen and Ortner 2016) provide partial justification for (1.1) in the setting of linear tight binding models at finite Fermi-temperature. However, in general, we can only expect

$$\begin{aligned} \eta \sim \beta ^{-1} \quad \text {as} \quad \beta \rightarrow \infty \end{aligned}$$

where $\beta$ is the inverse Fermi-temperature, meaning that for low-temperature regimes the practical value of (1.1) is limited. However, in the case of insulating multi-lattice materials (where there is a spectral gap in the system), the locality estimates are improved (Chen et al. 2019), and extended to the zero Fermi-temperature case. In this setting, the exponent $\eta$ is linear in the spectral gap.

In Sect. 2, we simultaneously extend both Chen and Ortner (2016); Chen et al. (2018) (in the case of finite Fermi-temperature) and Chen et al. (2019) (for insulators at finite or zero Fermi-temperature) to the nonlinear setting.

1.2.2 Point Defects in Insulators

Simulating local defects in materials remains an issue of great interest in the solid state physics and materials science communities (Pisani 1994; Stoneham 2001). See Cancès and Le Bris (2013) for a mathematical review of some works related to the modelling of point defects in materials science.

When considering point defects in the material, “pollution” of the spectrum may enter band gap which a priori affects the exponent $\eta$ in (1.1). However, by approximating the defect as a low rank perturbation, we show in Sect. 3 that the locality results only weakly depend on the defect and the estimates resemble the defect-free estimates for sites away from the defect core. That is, the exponents in the locality estimates depend only on a local environment of the particular atomic site. This extends results of Chen et al. (2019) to the nonlinear setting.

1.3 Notation

For an operator T, the discrete spectrum (isolated eigenvalues of finite multiplicity) and essential spectrum are denoted by $\sigma _{\mathrm {disc}}(T)$ and $\sigma _{\mathrm {ess}}(T) {:=}\sigma (T) \setminus \sigma _{\mathrm {disc}}(T)$, respectively.

For sequences $\psi$, the $j^\mathrm {th}$ entry is written $[\psi ]_j$ and the $\ell ^2$ norm is $\Vert \psi \Vert _{\ell ^2}$. For bounded $\mathbb {C}$-valued functions, we denote by $\Vert \cdot \Vert _\infty$ the supremum norm. For matrices (or operators with matrix entries only nonzero on a finite sub-matrix) M, we denote by $\Vert M\Vert _{\mathrm {max}}{:=}\max _{ij} |M_{ij}|$ the maximum-norm. For operators T on $\ell ^2$, $\Vert \cdot \Vert _{\mathrm{F}}$ denotes the Frobenius or Hilbert–Schmidt norm and $\Vert \cdot \Vert _{\ell ^2 \rightarrow \ell ^2}$ the operator norm.

On $\mathbb R^n$ or $\mathbb {C}$, we will denote the Euclidean norm by $|\cdot |$ and the open balls of radius $\delta$ about a and 0 by $B_\delta (a)$ and $B_\delta {:=}B_\delta (0)$, respectively. For a subset A of $\mathbb R^n$ or $\mathbb {C}$, we write $B_\delta (A) {:=}\{ x :\mathrm {dist}(x,A) < \delta \}$ where $\mathrm {dist}(b,A) {:=}\inf _{a\in A}|a - b|$. The Hausdorff distance between two sets A and B is denoted $\mathrm {dist}(A,B) {:=}\max \{ \sup _{a \in A} \mathrm {dist}(a, B), \sup _{b \in B} \mathrm {dist}(b, A) \}$.

We write $b + A = \{ b + a :a \in A \}$, $A - b {:=}\{a - b :a \in A\}$ and $r A {:=}\{r a :a \in A\}$. If A is finite then $\#A$ denotes the cardinality of A. For an index set A, we denote by $\delta _{ij}$ the Kronecker delta function for $i,j \in A$. The set of non-negative real numbers will be denoted by $\mathbb R_+$.

The symbol C will denote a generic positive constant that may change from one line to the next. In calculations, C will always be independent of Fermi-temperature. The dependencies of C will be clear from context or stated explicitly. When convenient to do so we write $f \lesssim g$ to mean $f \hbox {\,\,\char 054\,\,}C g$ for some generic positive constant as above.

2 Results: General Locality Estimates

2.1 Tight Binding Model

For a locally finite reference configuration $\varLambda \subset \mathbb R^d$ and displacement $u :\varLambda \rightarrow \mathbb R^d$, we write $\varvec{r}_{\ell k}(u) {:=}\ell + u(\ell ) - k - u(k)$ and $r_{\ell k}(u) {:=}|\varvec{r}_{\ell k}(u)|$. We consider displacements u satisfying the following uniform non-interpenetration condition:

(L):

There exists a non-interpenetration constant $\mathfrak {m}>0$ such that $r_{\ell k}(u) \hbox {\,\,\char 062\,\,}\mathfrak {m}|\ell - k|$ for all $\ell , k \in \varLambda$.

We consider $N_{\mathrm{b}}$ atomic orbitals per atom, indexed by $1\hbox {\,\,\char 054\,\,}a,b \hbox {\,\,\char 054\,\,}N_{\mathrm{b}}$. For given displacements u, and corresponding electronic densities, $\rho :\varLambda \rightarrow \mathbb R_+$, we define the two-centre tight binding Hamiltonian as follows:

(TB):

For $\ell , k \in \varLambda$ and $1\hbox {\,\,\char 054\,\,}a,b \hbox {\,\,\char 054\,\,}N_{\mathrm{b}}$, we suppose that the entries of the Hamiltonian take the form

$$\begin{aligned} \mathcal {H}(u;\rho )_{\ell k}^{ab} {:=}h^{ab}_{\ell k}\left( \varvec{r}_{\ell k}(u)\right) + v(\rho (\ell )) \delta _{\ell k} \delta _{ab} \end{aligned}$$

(2.1)

where $h^{ab}_{\ell k}:\mathbb R^d \rightarrow \mathbb R$ are $\nu$ times continuously differentiable for some $\nu \hbox {\,\,\char 062\,\,}1$ and v is a bounded smooth function on $(0,\infty )$ with bounded derivatives.

Further, we assume that there exist $h_0, {\gamma _0}>0$ such that, for each $1 \hbox {\,\,\char 054\,\,}j \hbox {\,\,\char 054\,\,}\nu$,

$$\begin{aligned} \left| h^{ab}_{\ell k}(\xi )\right| \hbox {\,\,\char 054\,\,}h_0\,e^{-{\gamma _0}|\xi |} \quad \text {and} \quad \left| \partial ^\alpha h_{\ell k}^{ab}(\xi )\right| \hbox {\,\,\char 054\,\,}h_0\, e^{-{\gamma _0}|\xi |} \quad \forall \xi \in \mathbb R^d \end{aligned}$$

(2.2)

for all multi-indices $\alpha \in \mathbb N^d$ with $|\alpha |_1 = j$.

Finally, we suppose that $h^{ab}_{\ell k}(\xi ) = h_{k\ell }^{ba}(-\xi )$ for all $\xi \in \mathbb R^d$, $1\hbox {\,\,\char 054\,\,}a,b\hbox {\,\,\char 054\,\,}N_{\mathrm{b}}$ and $\ell , k \in \varLambda$.

Remark 1

(Tight binding Hamiltonian)

(i) :

Despite being able to define the tight binding Hamiltonian for any fixed electronic density $\rho$, we will be concerned with the case that $\rho$ satisfies a self-consistency condition $\rho = F(u;\rho )$ (see (SC) on page 7 for the definition of F). This equation introduces the nonlinearity to the model; the main departure from the previous works (Chen and Ortner 2016; Chen et al. 2019).

(ii) :

The constants $h_0$ and $\gamma _0$ in (2.2) are independent of the atomic sites.

(iii) :

The symmetry assumption in (TB) means that the Hamiltonian, $\mathcal {H}(u;\rho )$, is symmetric and thus the spectrum is real. Moreover, $\sigma (\mathcal {H}(u;\rho )) \subset [\underline{\sigma },\overline{\sigma }]$ for some $\underline{\sigma },\overline{\sigma }$ depending on $\mathfrak {m}, d, h_0, \gamma _0, \Vert v\Vert _\infty$ and are independent of the size of the system and the displacement u satisfying (L) with the constant $\mathfrak {m}$. In fact, by generalising the proof of (Chen and Ortner 2016, Lemma 4) to the setting we consider here, we obtain explicit bounds on the width of the spectrum: $\sigma (\mathcal {H}(u;\rho )) \subset \Vert v\Vert _\infty + C_d h_0 (\gamma _0 \mathfrak {m})^{-d} [-1,1]$.

(iv) :

The pointwise bound on $|h^{ab}_{\ell k}|$ in (2.2) is more general than many tight binding models which impose a finite cut-off radius. The assumption for $|\alpha | = 1$, states that there are no long range interactions in the model and so, in particular, we assume that Coulomb interactions have been screened, which is typical in practical tight binding models (Cohen et al. 1994; Mehl and Papaconstantopoulos 1996; Papaconstantopoulos et al. 1997).

Following (Chen and Ortner 2016; Chen et al. 2018; Ortner and Thomas 2020; Chen et al. 2019), we consider finite energy displacements. For $u :\varLambda \rightarrow \mathbb R^d$, $\ell \in \varLambda$ and $\sigma \in \varLambda - \ell$, we define the finite difference $D_\sigma u(\ell ) {:=}u(\ell + \sigma ) - u(\ell )$ and the full finite difference stencil $Du(\ell ) {:=}(D_\sigma u(\ell ))_{\sigma \in \varLambda - \ell }$. We then let $\dot{\mathscr {W}}^{1,2}(\varLambda )$ be the set of finite energy displacements: for ${\Upsilon }> 0$,

$$\begin{aligned} \dot{\mathscr {W}}^{1,2}(\varLambda ) {:=}\bigg \{ u :\varLambda \rightarrow \mathbb R^d :\Vert Du\Vert _{\ell ^2_{\Upsilon }}^2 {:=}\sum _{\ell \in \varLambda } \sum _{\sigma \in \varLambda - \ell } e^{-2{\Upsilon }|\sigma |} |D_\sigma u(\ell )|^2 < \infty \bigg \} \end{aligned}$$

which is well defined since all the semi-norms $\Vert D\cdot \Vert _{\ell ^2_{\Upsilon }}$ are all equivalent for ${\Upsilon }>0$ (Chen et al. 2019). For the remainder of this paper, we fix an exponent ${\Upsilon }>0$. In this section, we require perturbations to have finite energy (see Lemma 2), whereas when considering the improved estimates in the case of point defects in Sect. 3, we also require the configurations to be given by finite energy displacements.

2.2 Spatial Decomposition of Quantities of Interest

Let $u :\varLambda \rightarrow \mathbb R^d$ satisfy (L) and suppose $\rho$ is an associated electronic density.

2.2.1 Local Quantities of Interest

We suppose $\mathfrak {o}:\mathbb R\rightarrow \mathbb R$ is a function that extends analytically to an open neighbourhood of $\sigma (\mathcal {H}(u;\rho ))$ in $\mathbb {C}$ and define the corresponding local quantities, $O_\ell$ (for $\ell \in \varLambda$), by

$$\begin{aligned} O_\ell (u;\rho ) {:=}- \frac{1}{2\pi i} \sum _a \oint _{\mathscr {C}_\mathfrak {o}} \mathfrak {o}(z) \left[ ( \mathcal {H}(u;\rho ) - z )^{-1} \right] _{\ell \ell }^{aa} \mathrm {d}z \end{aligned}$$

(2.3)

where $\mathscr {C}_\mathfrak {o}$ is a simple, closed, positively oriented contour contained within the region of holomorphicity of $\mathfrak {o}$, encircling $\sigma (\mathcal {H}(u;\rho ))$ and such that

$$\begin{aligned} \mathsf {d}_\mathfrak {o} {:=}\min _{z \in \mathscr {C}_\mathfrak {o}} \mathrm {dist}\big ( z, \sigma (\mathcal {H}(u;\rho )) \big ) > 0. \end{aligned}$$

(2.4)

For finite systems, we may diagonalise the Hamiltonian, $\mathcal {H}(u;\rho ) = \sum _s \lambda _s \vert \psi _s\rangle \langle \psi _s\vert$ (where $(\lambda _s,\psi _s)$ are normalised eigenpairs), and note that many quantities of interest, including the Helmholtz free energy, grand potential and the particle number functional (Chen et al. 2018; Chen and Ortner 2016), may be written as a sum of the local contributions (2.3) for some appropriate choice of $\mathfrak {o}$:

(2.5)

This decomposition is well known and follows from a spatial decomposition of the total density of states (Ercolessi 2005; Chen and Ortner 2016; Finnis 2003; Chen et al. 2018).

Our motivation comes from viewing (2.5) as the total energy of the system and (2.3) as a spatial decomposition of this energy. The aim of the present paper is to show the exponential localisation of (2.3) with respect to perturbing atomic positions, which, in the case of site energies, justifies IP and multi-scale models as discussed in the introduction.

2.2.2 Self-consistency

The electronic structure of the system is obtained by assigning electrons to the eigenstates of lowest energy, according to the Fermi-Dirac occupation distribution and subject to Pauli’s exclusion principle. We now make this precise for (possibly infinite) atomic systems.

After fixing an inverse Fermi-temperature $\beta >0$, we define $f(z - \mu ) = (1 + e^{\beta (z - \mu )})^{-1}$ to be the Fermi-Dirac occupation distribution where $\mu$ is a fixed chemical potential. We let ${\mathscr {C}}_f$ be a simple closed contour encircling $\sigma (\mathcal {H}(u;\rho ))$ and avoiding the singularities of $f(\,\cdot - \mu )$. That is, avoiding $\mu + i \pi \beta ^{-1} \left( 2\mathbb Z + 1\right)$ and such that $\mathsf {d}_f {:=}\min _{z \in \mathscr {C}_f} \mathrm {dist} \left( z, \sigma (\mathcal {H}(u;\rho )) \right) > 0$. In general, we may choose $\mathscr {C}_f$ so that $\mathsf {d}_{f} \hbox {\,\,\char 062\,\,}\tfrac{\pi }{2\beta }$. However, for insulators, as we shall see in Sect. 3, this constant may be chosen to be linear in the spectral gap at $\mu$.

For the case of zero Fermi-temperature, we take the pointwise $\beta \rightarrow \infty$ limit and define $f(z-\mu ) {:=}\chi _{(-\infty ,\mu )}(z) + \frac{1}{2}\chi _{\{\mu \}}(z)$. For insulating systems, there is a spectral gap at $\mu$ (that is, with $\mathsf {g} {:=}\inf [\mu ,+\infty ) \cap \sigma (\mathcal {H}(u;\rho )) - \sup (-\infty , \mu ] \cap \sigma (\mathcal {H}(u;\rho )) > 0$) and so $f(\,\cdot \, - \mu )$ is analytic in a neighbourhood of $\sigma (\mathcal {H}(u;\rho ))$. In this case, we let $\mathscr {C}_f$ be a simple closed contour encircling $\sigma (\mathcal {H}(u;\rho )) \cap (-\infty , \mu )$ and avoiding $\sigma (\mathcal {H}(u;\rho )) \cap [\mu , \infty )$ with $\mathsf {d}_f = \min _{z \in \mathscr {C}_f}\mathrm {dist} \left( z, \sigma (\mathcal {H}(u;\rho )) \right) \hbox {\,\,\char 062\,\,}\frac{1}{2}\mathsf {g}$.

For fixed $\beta \in (0,\infty ]$ as above, we can therefore define

$$\begin{aligned} F_\ell (u;\rho ) {:=}O_\ell (u;\rho ) \quad \text {with} \quad \mathfrak {o} = f(\,\cdot - \mu ) \quad \text {and} \quad \mathscr {C}_\mathfrak {o} = \mathscr {C}_f \end{aligned}$$

(2.6)

as in (2.3). To simplify notation in the following, we will write $F(u;\rho ) = (F_\ell (u;\rho ))_{\ell \in \varLambda }$.

For a given displacement, we consider corresponding self-consistent electronic densities, giving rise to the nonlinearity of the problem:

(SC):

For a displacement $u:\varLambda \rightarrow \mathbb R^d$ satisfying (L), we say that $\rho$ is an associated self-consistent electronic density if $\rho = F(u;\rho )$.

Remark 2

(Self-consistency) For a finite system, the self-consistency equation (SC) takes the following form:

$$\begin{aligned} \rho (\ell ) = \sum _{s,a} f(\lambda _s - \mu ) [\psi _s]_{\ell a}^2 \end{aligned}$$

(2.7)

where $\mathcal {H}(u;\rho ) =\sum _s \lambda _s \vert \psi _s\rangle \langle \psi _s\vert$ for normalised eigenpairs $(\lambda _s, \psi _s)$. Written in this form, we can see that $f(\lambda _s - \mu )$ represents the occupation number for the electronic state s.

2.2.3 Stability

We wish to show that, for fixed $u:\varLambda \rightarrow \mathbb R^d$ and an associated self-consistent electronic density $\rho$, the quantities $O_\ell (u;\rho )$ are exponentially localised. As shown in Chen et al. (2019) for the linear model (that is, neglecting the $v(\rho (\ell )) \delta _{\ell k}$ term in (2.1)), the exponent in these locality results are linear in $\mathsf {d}_\mathfrak {o}$ from (2.4). For the more complicated, nonlinear model that we consider here, the locality also depends on the stability of the model, which we will now discuss.

Supposing that $(u,\rho )$ satisfies a natural stability condition (defined in (STAB) on page 8, below), it is possible to rewrite the local quantities of interest as a function of the displacements. That is, for $\widetilde{u}$ in a neighbourhood of u, there exists a locally unique $\widetilde{\rho } = \widetilde{\rho }(\widetilde{u})$ in a neighbourhood of $\rho$ such that $(\widetilde{u},\widetilde{\rho }(\widetilde{u}))$ satisfies (SC) and we can therefore write

$$\begin{aligned} O_\ell ^\mathrm {sc}(\widetilde{u}) {:=}O_\ell \big (\widetilde{u};\widetilde{\rho }(\widetilde{u})\big ). \end{aligned}$$

(2.8)

Moreover, the mapping $\widetilde{u} \mapsto O^\mathrm {sc}_\ell (\widetilde{u})$ is $\nu$ times continuously differentiable in a neighbourhood of u. See Lemma 2 for the rigorous statement.

We can now consider the derivatives of the local quantities of interest with respect to the perturbation of atomic positions. Using (2.3), it is sufficient to consider the derivatives of the resolvent operators. Since the linear contribution has been studied in Chen et al. (2019), we are only concerned with the additional nonlinear part, which involves derivatives of the electronic density. Due to the self-consistency, we obtain the formula

$$\begin{aligned} \frac{ \partial \rho (\ell ) }{ \partial [u(m)]_i } = \left[ (I - \mathscr {L}(u;\rho ))^{-1} \phi ^{(m)} \right] _\ell \end{aligned}$$

(2.9)

where the stability operator, $\mathscr {L}(u;\rho )$, is the Jacobian of $F(u;\rho )$ with respect to $\rho$ and $\phi ^{(m)} \in \ell ^2(\varLambda )$ satisfies $|\phi _\ell ^{(m)}| \lesssim e^{-\eta r_{\ell m}(u)}$ for some $\eta > 0$. Therefore, the following stability condition, which we take from E and Lu (2010, (2012), is the minimal starting assumption required for the analysis:

(STAB):

We say $(u,\rho )$ is stable if $I - \mathscr {L}(u;\rho )$ is invertible as an operator on $\ell ^2(\varLambda )$ where $\mathscr {L}(u;\rho )$ is the Jacobian of $F(u;\rho )$ with respect to $\rho$. For a stable configuration $(u,\rho )$, we write

$$\begin{aligned} \mathsf {d}_\mathscr {L} {:=}\mathrm {dist}(1, \mathscr {L}(u;\rho ) ) > 0. \end{aligned}$$

(2.10)

Remark 3

(i) :

Equivalently, (2.10) states that $\Vert \left( I - \mathscr {L}(u;\rho ) \right) ^{-1} \Vert _{\ell ^2\rightarrow \ell ^2} \hbox {\,\,\char 054\,\,}\mathsf {d}_\mathscr {L}^{-1}$.

(ii) :

A simple calculation reveals that $\mathscr {L}(u;\rho ):\ell ^2(\varLambda ) \rightarrow \ell ^2(\varLambda )$ has matrix entries

$$\begin{aligned} \mathscr {L}(u;\rho )_{\ell k} = \frac{1}{2\pi i} \oint _{\mathscr {C}_f} f(z-\mu ) \sum _{a,b = 1}^{N_{\mathrm{b}}} \left( \left[ \left( \mathcal {H}(u;\rho ) - z \right) ^{-1} \right] ^{ab}_{\ell k} \right) ^2 \mathrm {d}z \, v^\prime (\rho (k)). \end{aligned}$$

(2.11)

2.2.4 Locality

We are now in a position to state general locality estimates for $\beta \in (0,\infty )$ and for insulators in the case $\beta = \infty$:

Theorem 1

(General Locality Estimates) Suppose $(u,\rho )$ satisfies (L), (SC), (STAB) and let $\mathsf {d}_{\mathfrak {o}}$, $\mathsf {d}_{f}$, $\mathsf {d}_\mathscr {L}>0$ be the constants from (2.4) and (2.10). Then, for $1 \hbox {\,\,\char 054\,\,}j \hbox {\,\,\char 054\,\,}\nu$, $\ell \in \varLambda$, $\varvec{m} = (m_1,\dots ,m_j) \in \varLambda ^j$ and $1\hbox {\,\,\char 054\,\,}i_1,\dots ,i_j \hbox {\,\,\char 054\,\,}d$, there exists $C_{j}>0$ such that

$$\begin{aligned} \left| \frac{\partial ^j O^{\mathrm {sc}}_{\ell }(u)}{\partial [u(m_1)]_{i_1}\dots \partial [u(m_j)]_{i_j}} \right| \hbox {\,\,\char 054\,\,}C_{j} e^{-\eta \sum _{l=1}^j r_{\ell m_l}(u) } \end{aligned}$$

where $\eta {:=}c \min \left\{ 1, \mathsf {d}_{\mathfrak {o}}, \mathsf {d}_{f}, c_f \mathsf {d}_{\mathscr {L}}\right\}$, $c_f {:=}\mathsf {d}_{f}^2 \min \{ 1, \mathsf {d}_{f}^{d+1} \}$ and c is a positive constant depending on $\gamma _0$, $h_0$, $\mathfrak {m}$, $N_{\mathrm{b}}$, d, j, $\Vert Du\Vert _{\ell ^2_{\Upsilon }}$, $\Vert v^\prime \Vert _\infty$ and on the length of $\mathscr {C}_f$.

Proof

The proof of this result follows the analogous proof in the linear case (Chen et al. 2019), together with bounds on the nonlinear contribution (2.9). Full details are presented in Sect. 6. $\square$

3 Results: Point Defects in Insulating Materials

Now we consider the specific example of point defect reference configurations. In this case, we show “improved” locality estimates in which the pre-factors and exponents behave like the corresponding reference constants.

3.1 Point Defect Configurations

We suppose that $\varLambda ^{\mathrm {ref}} \subset \mathbb R^d$ is a multi-lattice and $\rho ^{\mathrm {ref}}$ is a corresponding self-consistent electronic density:

(REF):

We suppose that there exist a non-singular matrix $\mathsf {A} \in \mathbb R^{d\times d}$ and a unit cell $\varGamma \subset \mathbb R^d$ such that $\varGamma$ is finite, contains the origin and

$$\begin{aligned} \varLambda ^{\mathrm {ref}} {:=}\bigcup _{\gamma \in \mathbb Z^d} \left( \varGamma + \mathsf {A}\gamma \right) . \end{aligned}$$

Moreover, we require the Hamiltonian to satisfy the following translational invariance property: $h^{ab}_{\ell + \mathsf {A}\gamma _1, k + \mathsf {A}\gamma _2}(\xi ) = h^{ab}_{\ell k}(\xi )$ for all $\gamma _1, \gamma _2 \in \mathbb Z^d$, $\xi \in \mathbb R^d$ and $1\hbox {\,\,\char 054\,\,}a,b \hbox {\,\,\char 054\,\,}N_{\mathrm{b}}$. Further, we suppose $\rho ^{\mathrm {ref}}$ is a translational invariant, self-consistent electronic density: $(x^{\mathrm {ref}}, \rho ^{\mathrm {ref}})$ satisfies (SC) where $x^{\mathrm {ref}}:\varLambda ^{\mathrm {ref}}\rightarrow \varLambda ^{\mathrm {ref}}$ is the identity configuration on $\varLambda ^{\mathrm {ref}}$ and $\rho ^{\mathrm {ref}}(\ell + \mathsf {A}\gamma ) = \rho ^{\mathrm {ref}}(\ell )$ for all $\ell \in \varGamma$ and $\gamma \in \mathbb Z^d$.

It is possible that there are many electronic densities $\rho ^{\mathrm {ref}}$ satisfying the conditions of (REF). However, this is an issue that we will not concern ourselves with here, and, for the remainder of this paper, we will simply fix any electronic density $\rho ^{\mathrm {ref}}$ satisfying (REF).

The translational invariance property of the Hamiltonian states that all $\ell + \mathsf {A}\gamma$ are of the same atomic species.

To simplify notation, when we consider the reference configuration ($\varLambda = \varLambda ^{\mathrm {ref}}$), we will write,

$$\begin{aligned} \mathcal {H}^{\mathrm {ref}} {:=}\mathcal {H}(\varvec{0};\rho ^{\mathrm {ref}}) \quad \text {and} \quad \mathscr {L}^{\mathrm {ref}} {:=}\mathscr {L}(\varvec{0};\rho ^{\mathrm {ref}}). \end{aligned}$$

By exploiting the translational invariance of the reference configuration, and by use of the Bloch transform (Kittel 2004), we may conclude that $\sigma (\mathcal {H}^{\mathrm {ref}})$ and $\sigma (\mathscr {L}^{\mathrm {ref}})$ can be written as a union of finitely many spectral bands:

$$\begin{aligned} \sigma (\mathcal {H}^{\mathrm {ref}}) = \bigcup _n \lambda _n(\varGamma ^\star ) \qquad \text {and} \qquad \sigma (\mathscr {L}^{\mathrm {ref}}) = \bigcup _n \varepsilon _n(\varGamma ^\star ) \end{aligned}$$

(3.1)

where $\varGamma ^\star \subset \mathbb R^d$ is a compact, connected set and $\lambda _n$ and $\varepsilon _n$ are continuous functions. Full details follow the calculations of E and Lu (2010) and are given in Appendix 1 for completeness.

In the following, we consider insulating materials and thus assume that there is a band gap in the reference Hamiltonian:

(GAP):

We assume $\mu \not \in \sigma (\mathcal {H}^{\mathrm {ref}})$ and define

$$\begin{aligned} \mathsf {g}^{\mathrm {ref}} {:=}\inf \left( \sigma (\mathcal {H}^{\mathrm {ref}}) \cap (\mu ,\infty ) \right) - \sup \left( \sigma (\mathcal {H}^{\mathrm {ref}}) \cap (-\infty ,\mu ) \right) > 0. \end{aligned}$$

Therefore, in the finite temperature case ($\beta < \infty$), we may consider a contour $\mathscr {C}_f$ as in (2.6) with $\mathsf {d}^{\mathrm {ref}}_f {:=}\min _{z \in \mathscr {C}_f} \mathrm {dist}(z,\sigma (\mathcal {H}^{\mathrm {ref}})) \hbox {\,\,\char 062\,\,}\frac{1}{2}\mathrm {dist}\left( \mu , \sigma (\mathcal {H}^{\mathrm {ref}})\right)$. In particular, we have $\mathsf {d}^{\mathrm {ref}}_f \not \rightarrow 0$ in the zero temperature limit. On the other hand, for zero temperature ($\beta = \infty$), we may choose $\mathscr {C}_f$ as in (2.6) with $\mathsf {d}^{\mathrm {ref}}_f = \frac{1}{2}\mathsf {g}^{\mathrm {ref}}$.

Further, we suppose that the reference configuration is stable:

$\mathbf (STAB ^\mathbf{ref }{} \mathbf )$:

$\exists \,\mathsf {d}_{\mathscr {L}}^{{\mathrm {ref}}}>0$ such that $I - \mathscr {L}^{\mathrm {ref}}$ is invertible with $\Vert (I - \mathscr {L}^{\mathrm {ref}})^{-1}\Vert _{ \ell ^2 \rightarrow \ell ^2 } \hbox {\,\,\char 054\,\,}(\mathsf {d}_{\mathscr {L}}^{{\mathrm {ref}}})^{-1}$.

Now, given a reference, $\varLambda ^{\mathrm {ref}}$, as above, we consider a point defect configuration, $\varLambda$, satisfying:

(P):

There exists $R_\mathrm {def}>0$ such that $\varLambda \setminus B_{R_\mathrm {def}} = \varLambda ^{\mathrm {ref}} \setminus B_{R_\mathrm {def}}$ and $\varLambda \cap B_{R_\mathrm {def}}$ is finite.

We will consider electronic densities $\rho : \varLambda \rightarrow \mathbb R_+$ satisfying the following mild technical assumption on the far-field behaviour:

(FF):

$\lim \limits _{|\ell | \rightarrow \infty }|v(\rho (\ell )) - v(\rho ^{\mathrm {ref}}(\ell ))| = 0$ where $\rho ^{\mathrm {ref}}$ is the fixed electronic density satisfying (REF).

This assumption is explained in more detail in Remark 4, below.

We now restrict the class of admissible configurations by considering finite energy displacements, $u \in \dot{\mathscr {W}}^{1,2}(\varLambda )$, and show that, for such displacements, the spectra can be described in terms of the reference spectra.

Lemma 1

(Perturbation of the Spectrum) Let $u \in \dot{\mathscr {W}}^{1,2}(\varLambda )$ satisfy $\mathbf (L)$ and $\rho$ be an associated electronic density satisfying (SC) and (FF).

Then, for all $\delta > 0$, there exists $S_\delta$ such that $\#\left( \sigma ( \mathcal {H}(u;\rho ) ) \setminus B_\delta ( \sigma (\mathcal {H}^{\mathrm {ref}}) ) \right) + \#\left( \sigma ( \mathscr {L}(u;\rho ) ) \setminus B_\delta ( \sigma (\mathscr {L}^{\mathrm {ref}}) ) \right) \hbox {\,\,\char 054\,\,}S_\delta$.

Proof

We may apply Lemmas 6 and 7 together with (Chen et al. 2019, Proof of Lemma 3) to conclude. $\square$

Supposing that (GAP) is satisfied, Lemma 1 states that there are at most finitely many isolated eigenvalues lying inside the band gap and bounded away from the band edges.

3.2 Improved Locality Estimates

Just as in the linear case (Chen et al. 2019), a Combes–Thomas resolvent estimate (Combes and Thomas 1973) applied to the spectral projections corresponding to the finitely many eigenvalues bounded away from the edges of the bands, together with a finite rank update formula (i.e. the Woodbury identity), allows us to approximate $(\mathcal {H}(u;\rho ) - z)^{-1}$ and $(I - \mathscr {L}(u;\rho ))^{-1}$ by finite rank updates of the reference resolvents $(\mathcal {H}^{\mathrm {ref}} - z)^{-1}$ and $(I - \mathscr {L}^{\mathrm {ref}})^{-1}$, respectively. Therefore, by applying the locality estimates of Theorem 1 on the reference spectrum, we obtain the following improved estimates for point defect configurations:

Theorem 2

(Improved Locality Estimates) Suppose that the reference configuration $\varLambda ^{\mathrm {ref}}, \rho ^{\mathrm {ref}}$ satisfies (REF), (GAP) and $\mathbf (STAB ^\mathbf{ref }{} \mathbf )$. Moreover, we fix $\varLambda$ satisfying (P), $u\in \dot{\mathscr {W}}^{1,2}(\varLambda )$ and $\rho$ satisfying (L), (SC), (STAB) and (FF). Then,

• (i) for $1 \hbox {\,\,\char 054\,\,}j \hbox {\,\,\char 054\,\,}\nu$, $\ell \in \varLambda$, $\varvec{m} = (m_1,\dots ,m_j) \in \varLambda ^j$ and $1\hbox {\,\,\char 054\,\,}i_1,\dots ,i_j \hbox {\,\,\char 054\,\,}d$, there exists positive constants $C_j = C_j(\ell ,\varvec{m})$, $\eta _j = \eta _j(\ell ,\varvec{m})$ such that

  $$\begin{aligned} \left| \frac{\partial ^j O^\mathrm {sc}_{\ell }(u)}{\partial [u(m_1)]_{i_1} \dots \partial [u(m_j)]_{i_j}} \right| \hbox {\,\,\char 054\,\,}C_{j} e^{-\eta _j \sum _{l=1}^j r_{\ell m_l}(u) } \end{aligned}$$

  where $\eta _j {:=}c_j \min \left\{ 1, \mathsf {d}^{\mathrm {ref}}_{\mathfrak {o}}, \mathsf {d}^{\mathrm {ref}}_{f}, c_{f}(\ell ,\varvec{m}) \mathsf {d}^{\mathrm {ref}}_{\mathscr {L}} \right\}$, $c_j>0$ depends on j, $\gamma _0$, $h_0$, $\mathfrak {m}$, d, $\Vert v^\prime \Vert _\infty$ and on the lengths of $\mathscr {C}_\mathfrak {o}$, $\mathscr {C}_f$ and $c_{f}$ is a constant depending only on $\ell$, $\varvec{m}$, ${\mathsf {d}^{\mathrm {ref}}_{f}}$ and ${\mathsf {d}_{f}}$.

• (ii) $C_j(\ell ,\varvec{m})$ is uniformly bounded and $c_{f}(\ell ,\varvec{m})$ is uniformly bounded away from zero independently of $\ell$ and $\varvec{m}$. Let $C^{\mathrm {ref}}_j{:=}C_j(\ell ,\varvec{m})$ and $c^{\mathrm {ref}}_{f}{:=}c_{f}(\ell ,\varvec{m})$ when $\varLambda = \varLambda ^{\mathrm {ref}}$, $u = 0$ and $\rho = \rho ^{\mathrm {ref}}$. If $\ell ,m_1, \dots , m_j \in B_R(\xi )$ for some $R > 0$, then $C_j(\ell ,\varvec{m}) \rightarrow C_j^{\mathrm {ref}}$ and $c_{f}(\ell ,\varvec{m}) \rightarrow c_{f}^{\mathrm {ref}}$ as $|\xi | \rightarrow \infty$, with exponential rates.

Remark 4

Here, we briefly give examples of when the assumption (FF) is satisfied and we show that it is equivalent to $\Vert \rho - \rho ^{\mathrm {ref}}\Vert _{\ell ^2} < \infty$.

(i):

Derivative of v sufficiently small. If $\Vert v^\prime \Vert _\infty$ is sufficiently small, we may treat the nonlinear tight binding model as a nonlinear perturbation of the corresponding linear model and would thus expect the main assumption $|v(\rho (\ell )) -v(\rho ^{\mathrm {ref}}(\ell ))| \rightarrow 0$ to be satisfied. In fact, we now show that if $\Vert v^\prime \Vert _\infty$ is sufficiently small, then the stronger condition $\Vert \rho - \rho ^{\mathrm {ref}}\Vert _{\ell ^2} < \infty$ is satisfied. Using the Combes–Thomas estimates (Lemma 3) for the resolvents we can conclude that there exists $\eta > 0$ such that

$$\begin{aligned} \left| \rho (\ell ) - \rho ^{\mathrm {ref}}(\ell ) \right|&\hbox {\,\,\char 054\,\,}C \Big | \sum _a \oint _{\mathscr {C}_f} f(z-\mu ) \big [ \mathscr {R}_z(u;\rho ) - \mathscr {R}_z^{\mathrm {ref}} \big ]^{aa}_{\ell \ell } \mathrm {d}z \Big | \nonumber \\&\hbox {\,\,\char 054\,\,}C \sup _{z \in \mathscr {C}_f} \sum _a \Big | \big [ \mathscr {R}_z(u;\rho ) \big ( \mathcal {H}^{\mathrm {ref}} - \mathcal {H}(u;\rho ) \big ) \mathscr {R}_z^{\mathrm {ref}} \big ]^{aa}_{\ell \ell } \Big | \nonumber \\&\hbox {\,\,\char 054\,\,}C \sum _{\ell _1,\ell _2} e^{-\eta \left( |\ell - \ell _1| + |\ell _1 - \ell _2| + |\ell _2 - \ell | \right) } |D_{\ell _2 - \ell _1} u(\ell _1)| \nonumber \\&\quad + C \sum _{\ell _1} e^{-2\eta |\ell - \ell _1|}|v(\rho (\ell _1)) - v(\rho ^{\mathrm {ref}}(\ell _1))| \end{aligned}$$

(3.2)

where $\mathscr {R}_z(u;\rho ) {:=}(\mathcal {H}(u;\rho ) - z)^{-1}$ and $\mathscr {R}^{\mathrm {ref}}_z {:=}(\mathcal {H}^{\mathrm {ref}} - z)^{-1}$. Here, we have abused notation as the operators $\mathscr {R}_z(u;\rho )$ and $\mathscr {R}_z^{\mathrm {ref}}$ are defined on different spatial domains. This issue is resolved in (Chen et al. 2019, Sect. 4.3) and also briefly explained on page 22. After squaring (3.2) and summing over $\ell \in \varLambda$, we obtain

$$\begin{aligned} \Vert \rho - \rho ^{\mathrm {ref}}\Vert _{\ell ^2}^2&\hbox {\,\,\char 054\,\,}C_1+ {C}_2\Vert v^\prime \Vert ^2_\infty \Vert \rho - \rho ^{\mathrm {ref}} \Vert _{\ell ^2}^2. \end{aligned}$$

Therefore, if $C_2\Vert v^\prime \Vert ^2_\infty < 1$, then $\Vert \rho - \rho ^{\mathrm {ref}}\Vert ^2_{\ell ^2} \hbox {\,\,\char 054\,\,}C_1(1 - C_2\Vert v^\prime \Vert _\infty ^2)^{-1}$.

(ii):

Stability of the electronic structure. Another approach involves integrating along a path between $\rho$ and $\rho ^{\mathrm {ref}}$. In order to compare u and reference configuration, we must assume that (FF) is satisfied for $\rho (\varvec{0})$, a self-consistent electronic density associated with the identity configuration on $\varLambda$. By the translational invariance of the Hamiltonian (i.e. for all $c\in \mathbb R^d$, $\mathcal {H}(u;\rho ) = \mathcal {H}(u + c;\rho )$ where $(u+c)(\ell ) = u(\ell ) + c$), we obtain translational invariance of the quantities of interest (as in Chen et al. 2018). In particular, the quantities of interest may be written as functions of the finite difference stencil $Du(\ell )$ for some $\ell \in \varLambda$. Therefore, the electronic density solving $\rho = F(u;\rho )$ can also be written as a function of $Du(\ell )$. Now since $\left| \frac{\partial \rho (\ell )}{\partial u(m)} \right| \lesssim e^{-\eta r_{\ell m}}$ (see (6.10)), we formally obtain

$$\begin{aligned} |v(\rho (\ell )) - v(\rho ^{\mathrm {ref}}(\ell ))|&\hbox {\,\,\char 054\,\,}\Vert v^\prime \Vert _\infty \left| \int _0^1 \sum _{\sigma \in \varLambda - \ell } \frac{\partial \rho _t(\ell )}{\partial D_\sigma u(\ell )} \cdot D_\sigma u(\ell ) \mathrm {d}t \right| \nonumber \\&\quad +|v(\rho (\varvec{0})(\ell )) - v(\rho ^{\mathrm {ref}}(\ell ))| \nonumber \\&\lesssim \sum _{\sigma \in \varLambda - \ell } e^{-\eta |\sigma |} |D_\sigma u(\ell )| + |v(\rho (\varvec{0})(\ell )) - v(\rho ^{\mathrm {ref}}(\ell ))| \end{aligned}$$

(3.3)

where $\rho _t {:=}\rho (tDu(\ell ))$. Therefore, by taking $\ell \rightarrow \infty$ we may conclude. However, in (3.3), we have assumed that along the linear path between u and $\varvec{0}$, the electronic density is a well-defined differentiable function of the displacement. Generalising the argument above, we only need to assume that there exists a sequence of displacements such that we can integrate along a piecewise linear path between u and $\varvec{0}$. That is, along the path, we need unique self-consistent electronic densities, the uniform non-interpenetration assumption to be satisfied with a uniform constant and, in the case of zero Fermi-temperature, the spectrum of the Hamiltonian must avoid the chemical potential.

(iii):

An equivalent assumption. We claim that (FF) is equivalent to the (a priori stronger) condition that $\Vert \rho - \rho ^{\mathrm {ref}}\Vert _{\ell ^2} < \infty$. Indeed, by assuming that $|v(\rho (\ell )) - v(\rho ^{\mathrm {ref}}(\ell ))| \rightarrow 0$, the diagonal operator defined by $D_{\ell \ell } {:=}v(\rho (\ell )) - v(\rho ^{\mathrm {ref}}(\ell ))$ is compact and so Lemma 6 (given below) allows us to approximate the Hamiltonian $\mathcal {H}(u;\rho )$ with a finite rank update of $\mathcal {H}^{\mathrm {ref}}$. We can therefore use (3.2) to obtain the following stronger bound: for all $\delta > 0$, there exists a Hilbert–Schmidt operator $P^\delta$ such that $\Vert P^\delta \Vert _{\mathrm{F}} \hbox {\,\,\char 054\,\,}\delta$ and $|\rho (\ell ) - \rho ^{\mathrm {ref}}(\ell )| \hbox {\,\,\char 054\,\,}C_\delta e^{-\eta |\ell |} + P^\delta _{\ell \ell }$. This is an argument similar to (Chen et al. 2019, (4.18)–(4.20)). We have simply written $|v(\rho (\ell )) -v(\rho ^{\mathrm {ref}}(\ell ))| \rightarrow 0$ as an assumption in Lemma 1 and Theorem 2 to simplify the presentation and avoid the technical issues detailed above. We briefly remark here that this is the minimal assumption needed for our analysis to hold. Indeed, if (FF) is not satisfied, then the operator $\mathcal {H}(u;\rho ) - \mathcal {H}^{\mathrm {ref}}$ is not compact and thus the compact perturbation results which we rely on in the proofs cannot be applied.

4 Applications: Geometry Optimisation

As well as justifying a number of interatomic potential and multi-scale methods (Chen and Ortner 2016, 2017; Csányi et al. 2005), we may use the locality results of this paper to formulate limiting variational problems for infinite systems. That is, for a fixed configuration $u_0$ with associated self-consistent electronic density $\rho _0$ such that $(u_0, \rho _0)$ is stable, we can renormalise the total energy and define

$$\begin{aligned} \mathcal G^\beta (u) {:=}\sum _{\ell \in \varLambda } \left( G^\beta _\ell (u) - G^\beta _\ell (u_0) \right) \end{aligned}$$

(4.1)

where $G^\beta _\ell$ is given by (2.8) with $\mathfrak {g}^\beta (z) = \frac{2}{\beta } \log ( 1 - f_\beta (z-\mu ) )$ for finite Fermi-temperature and $\mathfrak {g}^\infty (z) = 2(z - \mu ) \chi _{(-\infty ,\mu )}(z)$ in the case of zero Fermi-temperature. By the stability of the configuration $(u_0,\rho _0)$, it follows from the locality results of this paper together with Chen et al. (2019) that (4.1) is well defined in a $\Vert D\cdot \Vert _{\ell ^2_{\Upsilon }(\varLambda )}$-neighbourhood of $u_0$.

We can then consider the following geometry relaxation problems

$$\begin{aligned} \overline{u} \in \arg \min \left\{ \mathcal G^\beta (u) :u \in B_\delta (u_0;\Vert D\cdot \Vert _{\ell ^2_{\Upsilon }}) \text { satisfies }{} \mathbf (L) \right\} \end{aligned}$$

(4.2)

where “$\arg \min$” denotes the set of local minimisers. We emphasise here that, in order to define these problems, we require the differentiability of the site energies and so can only define these problems locally around stable configurations. We may follow the proofs of Ortner and Thomas (2020), to extend the results to the case of nonlinear tight binding models. For example, we may show the following:

Theorem 3

Suppose that $\mu$ is fixed such that (GAP) is satisfied and that $\overline{u}$ solves (4.2) for $\beta = \infty$ such that $\langle \delta ^2\mathcal {G}^\infty (\overline{u})v, v\rangle \hbox {\,\,\char 062\,\,}c_0\Vert Dv \Vert ^2_{\ell ^2_{\Upsilon }}$ for all $v \in \dot{\mathscr {W}}^{1,2}(\varLambda )$ and some $c_0>0$. Then, there exists $\overline{u}_\beta$ solving (4.2) with $\beta <\infty$ such that $\Vert D(\overline{u}_\beta - \overline{u})\Vert _{\ell ^2_{\Upsilon }} \lesssim e^{-c\beta }$.

Sketch of the Proof

Here, to distinguish between the finite and zero Fermi-temperature cases, we will write $F^{\beta }(u;\rho )$, $\mathscr {L}^\beta (u;\rho )$ and $F^\infty (u;\rho )$, $\mathscr {L}^\infty (u;\rho )$, respectively.

Firstly, we note that there exists a locally unique electronic density $\overline{\rho } = F^\infty (\overline{u};\overline{\rho })$. By stability, $I - \mathscr {L}^\infty (\overline{u};\overline{\rho })$ is invertible and so it follows from zero Fermi-temperature limit results (see Ortner and Thomas 2020, Lemma 5.9) that $I - \mathscr {L}^\beta (\overline{u};\overline{\rho }) = \left( I - \mathscr {L}^\infty (\overline{u};\overline{\rho }) \right) - (\mathscr {L}^\beta (\overline{u};\overline{\rho }) - \mathscr {L}^\infty (\overline{u};\overline{\rho }))$ is also invertible for all sufficiently large $\beta$. This means that for u in a neighbourhood of $\overline{u}$, there exists a locally unique electronic density ${\rho }_\beta$ satisfying ${\rho }_\beta = F^\beta ({u};{\rho }_\beta )$. Therefore, for u in a neighbourhood of $\overline{u}$, we may write $G^\beta _\ell ({u}) {:=}G^\beta _\ell ({u};{\rho }_\beta )$. This means that, for $\beta$ sufficiently large, we may apply the inverse function theorem on $\delta \mathcal G^\beta$ about $\overline{u}$ as in (Ortner and Thomas 2020, Theorem 2.3). $\square$

In Theorem 3, we restrict ourselves to the grand-canonical ensemble where there is a fixed chemical potential. By following the proofs of Ortner and Thomas (2020), one can also show analogous results for the canonical ensemble where the Fermi-level arises as a Lagrange multiplier for the particle number constraint.

5 Conclusions

We have extended the locality results of Chen et al. (2019) to nonlinear tight binding models. More specifically, the results of this paper are twofold: (i) we have written analytic quantities of interest (which includes the total energy of the system) as the sum of exponentially localised site contributions. Moreover, (ii) under a mild assumption on the electronic density, we have shown that point defects in the material only weakly affect the locality estimates. That is, away from the defect, where the local atomic environment resembles that of the corresponding defect-free configuration, the locality estimates resemble that of the defect-free case.

We have formulated variational problems for (possibly infinite) atomic systems and extended the zero Fermi-temperature limit results of Ortner and Thomas (2020) to the nonlinear setting. We believe that the thermodynamic limit results of Ortner and Thomas (2020) can also be extended to nonlinear tight binding models. The only additional technical detail is to show that the limiting configuration gives rise to stable configurations defined along the sequence of finite domain approximations. This means the choice of boundary condition and the number of electrons imposed plays a key role in the analysis. While we do not see any problem in extending the results of Ortner and Thomas (2020) for a supercell model, it is much less clear how and when the boundary effects may inhibit the stability of the electronic structure when considering clamped boundary conditions, for example.

The results of this paper represent a first natural stepping stone between the linear tight binding results of (Chen and Ortner 2016; Chen et al. 2018, 2019) towards more accurate electronic structure models, such as Kohn-Sham density functional theory. The strong energy locality estimates of this paper and the previous works (Chen and Ortner 2016; Chen et al. 2019) require an implicit screening assumption to justify the short-ranged Hamiltonian. Therefore, the key question of incorporating long-range Coulomb interactions into our analysis remains an important direction for future study.

6 Proofs: General Locality Estimates

In order to simplify the notation in the following, we will write $\mathcal {H}(u;\rho ) = \mathcal {H}^\mathrm {L}(u) + \mathcal {H}^\mathrm {NL}(\rho )$ where $\mathcal {H}^\mathrm {L}(u)_{\ell k}^{ab} {:=}h_{\ell k}^{ab}( \varvec{r}_{\ell k}(u) )$ and $\mathcal {H}^\mathrm {NL}(\rho )_{\ell k}^{ab} = v(\rho (\ell )) \delta _{\ell k} \delta _{ab}$. Further, we denote the resolvent operator by $\mathscr {R}_z(u;\rho ) {:=}(\mathcal {H}(u;\rho ) - z)^{-1}$.

6.1 Preliminaries

Firstly, we prove that we may write the local quantities of interest as a function of the displacement, (2.8):

Lemma 2

Suppose that $(u,\rho )$ satisfies (L), (SC) and (STAB). Then, there exist $\delta _u, \delta _\rho > 0$ such that for all $\widetilde{u} :\varLambda \rightarrow \mathbb R^d$ with $\Vert D(\widetilde{u} - u)\Vert _{\ell ^2_{\Upsilon }(\varLambda )} < \delta _u$, there exists a unique electronic density $\widetilde{\rho } = \widetilde{\rho }(\widetilde{u})$ satisfying $\Vert \widetilde{\rho } - \rho \Vert _{\ell ^2(\varLambda )} < \delta _\rho$ and $\widetilde{\rho }= F(\widetilde{u};\widetilde{\rho })$.

Further, the mapping $\widetilde{u} \mapsto \widetilde{\rho }$ is smooth.

Sketch of the Proof

We apply the implicit function theorem on $T(\widetilde{u};\widetilde{\rho }) {:=}\widetilde{\rho }- F(\widetilde{u};\widetilde{\rho })$, a smooth map in a neighbourhood of $(u,\rho )$. Since $\rho$ satisfies (SC), we have $T(u;\rho ) = 0$. By (STAB), we have: for each $\widetilde{u}$ in a neighbourhood of u, there exists a locally unique $\widetilde{\rho }= \widetilde{\rho }(\widetilde{u})$ with $\widetilde{\rho }= F(\widetilde{u};\widetilde{\rho })$.

The fact that F is indeed a smooth map in a neighbourhood of $(u,\rho )$ follows from the fact that small perturbations in u and $\rho$ lead to small perturbations in $\sigma (\mathcal {H}(u;\rho ))$Kato (1995). This means that the fixed contour $\mathscr {C}_f$, which depends on $(u,\rho )$, can be used in the definition of F in a neighbourhood of $(u,\rho )$.

For full details in a slightly different setting, see (E and Lu 2010, Theorem 5.3). $\square$

We now state a Combes–Thomas type estimate Combes and Thomas (1973) for the resolvent:

Lemma 3

(Combes–Thomas Resolvent Estimates) Given $u \in \dot{\mathscr {W}}^{1,2}(\varLambda )$ satisfying (L), suppose that T(u) is an operator on $\ell ^2(\varLambda \times \{1,\dots ,N_{\mathrm{b}}\})$ given by

$$\begin{aligned}{}[T(u)w](\ell ;a) {:=}\sum _{k \in \varLambda }\sum _{1\hbox {\,\char 054\,}b\hbox {\,\char 054\,}N_{\mathrm{b}}} T(u)_{\ell k}^{ab}w(k;b) \quad \text {where} \quad |T(u)_{\ell k}^{ab}|\hbox {\,\,\char 054\,\,}c_T e^{-\gamma _T r_{\ell k}(u)} \end{aligned}$$

for some $c_T,\gamma _T > 0$. Then, if $z \in \mathbb {C}$ with $\mathfrak {d} {:=}\mathrm {dist}(z,\sigma (T(u))) > 0$, we have,

$$\begin{aligned} \left| \left[ (T(u) - z)^{-1} \right] _{\ell k}^{ab} \right| \hbox {\,\,\char 054\,\,}2 \mathfrak {d}^{-1} e^{-\gamma _\mathrm {CT}r_{\ell k}(u)} \end{aligned}$$

where $\gamma _\mathrm {CT}{:=}c\gamma _T \min \{1, c_T^{-1}\gamma _T^d \mathfrak {d}\}$ where $c>0$ depends only on $\Vert Du\Vert _{\ell ^2_{\Upsilon }}$, $\mathfrak {m}$ and d.

Proof

The proof is analogous to (Chen and Ortner 2016, Lemma 6) where the claimed $\mathfrak {d}$-dependence can be obtained by following the same proof and calculating the pre-factor in (Chen et al. 2019, (4.4)).

That is, it can be shown that

$$\begin{aligned} \sup _{\ell \in \varLambda } \sum _{k \in \varLambda } c_T e^{-\gamma _T r_{\ell k}(u)} (e^{\gamma _\mathrm {CT} r_{\ell k}(u)} - 1) \hbox {\,\,\char 054\,\,}C \frac{c_T}{\gamma _T^{d+1}} \gamma _\mathrm {CT}, \qquad \forall \,\gamma _\mathrm {CT} \hbox {\,\,\char 054\,\,}\frac{1}{2}\gamma _T \end{aligned}$$

(6.1)

for some $C>0$ depending only on $\Vert Du\Vert _{\ell ^2_{\Upsilon }}$, $\mathfrak {m}$ and d. The proof follows by choosing $\gamma _\mathrm {CT} > 0$ sufficiently small such that the right hand side of (6.1) is less than $\frac{1}{2}\mathfrak {d}$. $\square$

Remark 5

More careful analysis reveals that the above proof gives

$$\begin{aligned} \gamma _\mathrm {CT}{:=}\frac{1}{2}\gamma _T \min \left\{ 1, \left( \frac{\mathfrak {m}}{2} \right) ^{d+1} \frac{\gamma _T^{d} \mathfrak {d}}{d! \Vert Du\Vert _{\ell ^\infty } c_T} \right\} \quad \text {where} \quad \Vert Du\Vert _{\ell ^\infty } {:=}\sup _{\ell \in \varLambda } \sup _{\rho \in \varLambda - \ell } \frac{|D_{\rho }u(\ell )|}{|\rho |}. \end{aligned}$$

Here, we have used the fact that $\Vert D\cdot \Vert _{\ell ^\infty }$ defines a semi-norm that is equivalent to $\Vert D\cdot \Vert _{\ell ^2_{\Upsilon }}$Chen et al. (2019).

By applying Lemma 3 to $\mathcal {H}(u;\rho )$, we obtain locality estimates for the resolvents $\mathscr {R}_z(u;\rho )$: for $z \in \mathbb {C}$ with $\mathrm {dist} \left( z,\sigma ( \mathcal {H}(u;\rho ) ) \right) \hbox {\,\,\char 062\,\,}\mathfrak {d} > 0$ we have

$$\begin{aligned} \left| \mathscr {R}_z(u;\rho )_{\ell k}^{ab} \right| \hbox {\,\,\char 054\,\,}2 \mathfrak {d}^{-1} e^{-\gamma _{\mathrm {r}}(\mathfrak {d}) r_{\ell k}(u) } \end{aligned}$$

(6.2)

where $\gamma _{\mathrm {r}}(\mathfrak {d}) {:=}c \min \{1,\mathfrak {d}\}$ and c is a positive constant that depends only on $h_0$, $\gamma _0$, $\Vert Du\Vert _{\ell ^2_{\Upsilon }}$, $\mathfrak {m}$ and d. We will apply (6.2) for both $z \in \mathscr {C}_f$ (with $\mathfrak {d} = \mathsf {d}_f$) and $z \in \mathscr {C}_\mathfrak {o}$ (with $\mathfrak {d} = \mathsf {d}_\mathfrak {o}$).

Therefore, by (2.11), we have

$$\begin{aligned} \left| \mathscr {L}(u;\rho )_{\ell k}^{ab} \right| \hbox {\,\,\char 054\,\,}C \mathsf {d}_{f}^{-2} e^{-2\gamma _{\mathrm {r}}(\mathsf {d}_f) r_{\ell k}(u) }. \end{aligned}$$

(6.3)

Therefore, applying Lemma 3 again with T replaced with $\mathscr {L}(u;\rho )$ (and with $z = 1$), we obtain,

$$\begin{aligned} \left| \left[ (I - \mathscr {L}(u;\rho ))^{-1} \right] _{\ell k}^{ab} \right| \hbox {\,\,\char 054\,\,}2 \mathsf {d}_{\mathscr {L}}^{-1} e^{-\gamma _{\mathrm {s}} r_{\ell k}(u) } \end{aligned}$$

(6.4)

where $\mathsf {d}_\mathscr {L}$ is the constant from (2.10) and $\gamma _{\mathrm {s}} {:=}c_1 \gamma _{\mathrm {r}}(\mathsf {d}_f) \min \{1,\mathsf {d}_{f}^2\gamma _{\mathrm {r}}^d \mathsf {d}_{\mathscr {L}}\}$. By expanding $\gamma _{\mathrm {r}}$ in terms of $\mathsf {d}_f$, we obtain $\gamma _{\mathrm {s}} = c_2 \min \left\{ 1, \mathsf {d}_{f}, \mathsf {d}_{f}^2\mathsf {d}_{\mathscr {L}}, \mathsf {d}_{f}^{d+3}\mathsf {d}_{\mathscr {L}} \right\}$ for some $c_2 >0$ depending only on $h_0, \gamma _0, \Vert Du\Vert _{\ell ^2_{\Upsilon }}, \mathfrak {m}, d$, the length of $\mathscr {C}_f$, $\Vert v^\prime \Vert _\infty$ and $N_{\mathrm{b}}$.

6.2 Proof of Theorem 1: General Locality Estimates

We are now in a position to prove the locality estimates. Since we may write $O^\mathrm {sc}_\ell$ as an integral of the resolvent operator, derivatives of $O^\mathrm {sc}_\ell$ can be written as derivatives of the resolvent operators.

We start with the case $j = 1$: for $z \in \mathbb {C}$ with $\mathrm {dist}(z,\sigma (\mathcal {H}(u;\rho ))) \hbox {\,\,\char 062\,\,}\mathfrak {d}>0$, we have

$$\begin{aligned} \frac{ \partial [\mathscr {R}_z(u;\rho )]^{aa}_{\ell \ell } }{ \partial [u(m)]_i }&= - \left[ \mathscr {R}_z(u;\rho ) \frac{\partial \left[ \mathcal {H}(u;\rho ) \right] }{\partial [u(m)]_i} \mathscr {R}_z(u;\rho ) \right] ^{aa}_{\ell \ell } \nonumber \\&= - \left[ \mathscr {R}_z(u;\rho ) \frac{\partial \left[ \mathcal {H}^\mathrm {L}(u;\rho )\right] }{\partial [u(m)]_i} \mathscr {R}_z(u;\rho ) \right] ^{aa}_{\ell \ell }\nonumber \\&\quad - \sum _{k\in \varLambda } \sum _{b=1}^{N_{\mathrm{b}}} \left( [ \mathscr {R}_z(u;\rho ) ]^{ab}_{\ell k} \right) ^2 v^\prime (\rho (k)) \frac{\partial \rho (k)}{\partial [u(m)]_i}. \end{aligned}$$

(6.5)

Here, we have used the fact that,

$$\begin{aligned} \frac{\partial \left[ \mathcal {H}(u;\rho ) \right] _{\ell k}^{ab}}{\partial [u(m)]_i}&= \frac{\partial \mathcal {H}^\mathrm {L}(u)_{\ell k}^{ab}}{\partial [u(m)]_i} + v^\prime (\rho (\ell )) \frac{\partial \rho (\ell )}{\partial [u(m)]_i} \delta _{\ell k}\delta _{ab}. \end{aligned}$$

The first contribution in (6.5) can be treated by applying (6.2) as in (Chen and Ortner 2016; Chen et al. 2019):

$$\begin{aligned} \left| \left[ \mathscr {R}_z(u;\rho ) \frac{\partial [\mathcal {H}^\mathrm {L}(u)]}{\partial [u(m)]_i} \mathscr {R}_z(u;\rho ) \right] _{\ell \ell }^{aa} \right|&\hbox {\,\,\char 054\,\,}C\mathfrak {d}^{-2} e^{-\min \{ \gamma _{\mathrm {r}}(\mathfrak {d}), \gamma _0 \} r_{\ell m}(u)}. \end{aligned}$$

(6.6)

Now we move on to consider the nonlinear contribution in (6.5). By taking derivatives in the self-consistency equation for $\rho$ (that is, $\rho = F(u;\rho )$ from (SC)), we obtain the following identity,

$$\begin{aligned} \frac{\partial \rho (\ell )}{\partial [u(m)]_i}&= - \frac{1}{2\pi i} \sum _{a} \oint _{\mathscr {C}_f} f(z-\mu ) \frac{ \partial [\mathscr {R}_z(u;\rho )]_{\ell \ell } }{ \partial [u(m)]_i } \mathrm {d}z\\&= \frac{1}{2\pi i} \sum _a \oint _{\mathscr {C}_f} f(z-\mu ) \left[ \mathscr {R}_z(u;\rho ) \frac{\partial \left[ \mathcal {H}^\mathrm {L}(u)\right] }{\partial [u(m)]_i} \mathscr {R}_z(u;\rho ) \right] ^{aa}_{\ell \ell } \mathrm {d}z \\&\quad + \left[ \mathscr {L}(u;\rho ) \frac{\partial \rho }{\partial [u(m)]_i} \right] _{\ell } \end{aligned}$$

where $\mathscr {L}(u;\rho )$ is the stability operator given in (2.11). That is,

$$\begin{aligned} \frac{\partial \rho (\ell )}{\partial [u(m)]_i} = \left[ (I - \mathscr {L}(u;\rho ))^{-1} \phi ^{(m)} \right] _\ell \end{aligned}$$

(6.7)

where $\phi ^{(m)} \in \ell ^2(\varLambda )$ is given by

$$\begin{aligned} \phi ^{(m)}_\ell {:=}\frac{1}{2\pi i} \sum _{a = 1}^{N_{\mathrm{b}}} \oint _{\mathscr {C}_f} f(z-\mu ) \left[ \mathscr {R}_z(u;\rho ) \frac{\partial \left[ \mathcal {H}^\mathrm {L}(u;\rho ) \right] }{\partial [u(m)]_i} \mathscr {R}_z(u;\rho ) \right] ^{aa}_{\ell \ell } \mathrm {d}z. \end{aligned}$$

(6.8)

Applying (6.6) and using the fact that f is uniformly bounded, we have

$$\begin{aligned} \big | \phi ^{(m)}_\ell \big | \hbox {\,\,\char 054\,\,}C \mathsf {d}_{f}^{-2} e^{-\min \{ \gamma _{\mathrm {r}}, \gamma _0 \} r_{\ell m}(u)}. \end{aligned}$$

(6.9)

Combining (6.9) with the resolvent estimate for $(I - \mathscr {L}(u;\rho ))^{-1}$ from (6.4), we obtain

$$\begin{aligned} \left| \frac{\partial \rho (\ell )}{\partial [u(m)]_i} \right|&\hbox {\,\,\char 054\,\,}C \mathsf {d}_{f}^{-2} \mathsf {d}_{\mathscr {L}}^{-1} \sum _{k \in \varLambda } e^{-\gamma _{\mathrm {s}} r_{\ell k}(u)} e^{-\min \{ \gamma _{\mathrm {r}}, \gamma _0 \} r_{km}(u)} \nonumber \\&\hbox {\,\,\char 054\,\,}C \mathsf {d}_{f}^{-2} \mathsf {d}_{\mathscr {L}}^{-1} e^{-\frac{1}{2}\min \{ \gamma _{\mathrm {s}}, \gamma _{\mathrm {r}}, \gamma _0 \} r_{\ell m}(u) }. \end{aligned}$$

(6.10)

Therefore, we may bound the second term in (6.5): for $z \in \mathbb {C}$ with $\mathrm {dist}(z,\sigma (\mathcal {H}(u;\rho ))) \hbox {\,\,\char 062\,\,}\mathfrak {d}$, we have

$$\begin{aligned}&\sum _{k \in \varLambda } \sum _{1\hbox {\,\char 054\,}b \hbox {\,\char 054\,}N_{\mathrm{b}}} \left( \mathscr {R}_z(u;\rho )_{\ell k}^{ab} \right) ^2 v^\prime (\rho (k)) \frac{\partial \rho (k)}{\partial [u(m)]_i} \nonumber \\&\quad \hbox {\,\,\char 054\,\,}C \mathfrak {d}^{-2}\mathsf {d}_{f}^{-2} \mathsf {d}_{\mathscr {L}}^{-1} \sum _{k \in \varLambda } |v^\prime (\rho (k))| e^{-2\gamma _{\mathrm {r}}(\mathfrak {d}) r_{\ell k}(u)} e^{-\frac{1}{2}\min \{ \gamma _{\mathrm {s}}, \gamma _{\mathrm {r}}(\mathsf {d}_f), \gamma _0\} r_{k m}(u)} \nonumber \\&\quad \hbox {\,\,\char 054\,\,}C \mathfrak {d}^{-2}\mathsf {d}_{f}^{-2} \mathsf {d}_{\mathscr {L}}^{-1} \Vert v^\prime \Vert _{\infty } e^{-\frac{1}{4}\min \{ \gamma _{\mathrm {s}}, 4\gamma _{\mathrm {r}}(\mathfrak {d}), \gamma _{\mathrm {r}}(\mathsf {d}_f), \gamma _0 \} r_{\ell m}(u)}. \end{aligned}$$

(6.11)

Combining (6.6) and (6.11) with $\mathfrak {d} = \mathsf {d}_\mathfrak {o}$ and using the fact that $\mathfrak {o}$ is uniformly bounded along the contour $\mathscr {C}_\mathfrak {o}$, we can conclude the proof for $j=1$.

Higher derivatives can be treated by taking derivatives of (6.5). The first contribution in (6.5) is what arises in the linear case and so derivatives of this term can be treated in the same way as in Chen et al. (2019). We sketch the argument for $j=2$ for the second contribution in (6.5). We fix $k \in \varLambda$ and $b\in \{1,\dots ,N_{\mathrm{b}}\}$ and note

$$\begin{aligned}&\frac{\partial }{\partial u(n)} \left\{ \left( [\mathscr {R}_z(u;\rho )]^{ab}_{\ell k} \right) ^2 v^\prime (\rho (k)) \frac{\partial \rho (k)}{\partial u(m)} \right\} \nonumber \\&\quad = 2 [\mathscr {R}_z(u;\rho )]^{ab}_{\ell k} \frac{\partial [\mathscr {R}_z(u;\rho )]^{ab}_{\ell k}}{\partial u(n)} v^\prime (\rho (k)) \frac{\partial \rho (k)}{\partial u(m)} \nonumber \\&\qquad + \left( [\mathscr {R}_z(u;\rho )]^{ab}_{\ell k} \right) ^2 v^{\prime \prime }(\rho (k)) \frac{\partial \rho (k)}{\partial u(n)} \frac{\partial \rho (k)}{\partial u(m)}\nonumber \\&\qquad + \left( [\mathscr {R}_z(u;\rho )]^{ab}_{\ell k} \right) ^2 v^{\prime }(\rho (k)) \frac{\partial ^2 \rho (k)}{\partial u(n) \partial u(m)}. \end{aligned}$$

(6.12)

After summing over $k \in \varLambda$, the first two contributions in (6.12) may be bounded above by a constant multiple of $e^{-\eta (r_{\ell m}(u) + r_{\ell n}(u))}$ for some $\eta > 0$ depending only on the exponents in (6.2), (6.6), and (6.11). The final contribution in (6.12) involves the second derivative of the electronic density which may be bounded above as follows: using (6.7), we have

$$\begin{aligned} \bigg | \frac{\partial ^2 \rho (k)}{\partial u(n) \partial u(m)} \bigg |&\hbox {\,\,\char 054\,\,}\bigg | \bigg [ \frac{\partial (I - \mathscr {L}(u;\rho ))^{-1}}{\partial u(n)} \phi ^{(m)} \bigg ]_k \bigg | + \bigg | \bigg [ (I - \mathscr {L}(u;\rho ))^{-1} \frac{\partial \phi ^{(m)}}{\partial u(n)} \bigg ]_k \bigg | \nonumber \\&\hbox {\,\,\char 054\,\,}C e^{-\eta (r_{kn}(u) + r_{km}(u)) } \end{aligned}$$

(6.13)

where $\eta > 0$ depends only on the exponents in (6.2), (6.4), (6.10) and in the locality estimates of the first contribution in (6.5). The estimate in (6.13) is easy to prove but is lengthy and very similar to the calculations above and so is omitted. Using (6.13) and summing over $k \in \varLambda$ in (6.12) we can conclude.

7 Proofs: Improved Locality Estimates

Before we prove Theorem 2, we require an improved Combes–Thomas type estimate for the resolvent operators; see Lemma 5. In the following section, we discuss this result and explain how we can use it despite the fact that the reference and defect Hamiltonians are defined on different spatial domains. Then, in Sect. 7.2, we show that the operators $\mathcal {H}(u;\rho )$ and $\mathscr {L}(u;\rho )$ satisfy the conditions of Lemma 5 and thus prove Theorem 2.

7.1 Preliminaries

We show an improved resolvent estimate for operators on $\ell ^2(\varLambda \times \{1,\dots ,N_{\mathrm{b}}\})$ that can be decomposed into a reference operator and two perturbations that are small in the sense of rank and Frobenius norm (Lemma 5), respectively. First, we need a basic identity for the inverse of an updated operator (Hager 1989):

Lemma 4

(Woodbury) Suppose that A and P are operators on a Banach space such that A and $A + P$ are invertible. Then, $I + PA^{-1}$ and $I + A^{-1}P$ are invertible and

$$\begin{aligned} (A + P)^{-1}&= A^{-1} - A^{-1} (I + PA^{-1})^{-1} P A^{-1} \\&= A^{-1} - A^{-1}P (I + A^{-1}P)^{-1} A^{-1}. \end{aligned}$$

Proof

Firstly, $I + PA^{-1} = (A + P)A^{-1}$ is invertible with inverse $A(A+P)^{-1}$. Therefore, we have $(A+P)^{-1} - A^{-1} = (A+P)^{-1}[A - (A + P)] A^{-1} = - [(I+PA^{-1})A]^{-1}PA^{-1}$.

The second formulation can be shown similarly.

Using this Woodbury identity, we may prove the following “improved” Combes–Thomas estimate:

Lemma 5

(Improved Combes–Thomas Resolvent Estimate) Suppose $\delta , R > 0$, $T^{\mathrm {ref}}, T^\mathrm {FR}, T^\delta$ are operators on $\ell ^2(\varLambda \times \{1,\dots ,N_{\mathrm{b}}\})$ and define $T {:=}T^{\mathrm {ref}} + T^\mathrm {FR} + T^\delta$. Further, suppose that:

• $\left| [ T^{\mathrm {ref}} ]_{\ell k}^{ab} \right| \hbox {\,\,\char 054\,\,}c_T e^{-\gamma _T |\ell - k|}$ for some $c_T, \gamma _T > 0$,

• $[T^\mathrm {FR}]_{\ell k}^{ab} = 0$ if $\ell \not \in \varLambda \cap B_{R}$ or $k \not \in \varLambda \cap B_{R}$,

• $\Vert T^\delta \Vert _{\mathrm{F}} \hbox {\,\,\char 054\,\,}\delta$ and,

• $z \in \mathbb {C}$ with $\mathfrak {d} {:=}\mathrm {dist}(z,\sigma (T)) > 0$ and $\mathfrak {d}^{\mathrm {ref}} {:=}\mathrm {dist}(z,\sigma (T^{\mathrm {ref}})) - \delta > 0$.

Then, there exists a constant C, depending on $\delta$, R, $c_T$, $\gamma _T$, $\Vert T^\mathrm {FR}\Vert _{\mathrm {max}}$, $\mathfrak {d}$, $\mathfrak {d}^{\mathrm {ref}}$ and d, such that

$$\begin{aligned} \left| \left[ (T - z)^{-1} \right] _{\ell k}^{ab} \right| \hbox {\,\,\char 054\,\,}C_{\ell k} e^{-\gamma _\mathrm {CT}(\mathfrak {d}^{\mathrm {ref}}) |\ell - k|} \quad \text {where}\\ C_{\ell k} {:=}2 (\mathfrak {d}^{\mathrm {ref}})^{-1} + C(1+|z|)^{2} e^{-\gamma _\mathrm {CT}(\mathfrak {d}^{\mathrm {ref}}) \left( |\ell | + |k| - |\ell - k|\right) }, \end{aligned}$$

and $\gamma _\mathrm {CT}(\mathfrak {d}^{\mathrm {ref}}) {:=}c \gamma _T \min \{ 1, c_T^{-1}(\gamma _T)^d \mathfrak {d}^{\mathrm {ref}} \}$ is the constant from Lemma 3.

Proof

This proof closely follows the ideas of (Chen et al. 2019, Sect. 4.4) but for more general operators T. We sketch the argument for completeness.

After defining $\mathscr {R}_z {:=}(T - z)^{-1}$ and $\mathscr {R}_z^\delta {:=}(T^{\mathrm {ref}} + T^\delta - z)^{-1}$, we apply Lemma 4 to obtain:

$$\begin{aligned} \mathscr {R}_z&= \mathscr {R}^\delta _z - \mathscr {R}^\delta _z (I + T^\mathrm {FR}\mathscr {R}^\delta _z)^{-1} T^\mathrm {FR} \mathscr {R}^\delta _z = \mathscr {R}^\delta _z - \mathscr {R}^\delta _z T^{\mathrm {FR}} (I + \mathscr {R}^\delta _zT^\mathrm {FR})^{-1} \mathscr {R}^\delta _z. \end{aligned}$$

(7.1)

We will consider the two terms in (7.1) separately.

Firstly, since $\mathrm {dist}\left( \sigma (T^{\mathrm {ref}}), \sigma (T^{\mathrm {ref}} + T^\delta )\right) \hbox {\,\,\char 054\,\,}\Vert T^\delta \Vert _{\mathrm {F}} \hbox {\,\,\char 054\,\,}\delta$ Kato (1995), we may apply Lemma 3 directly to conclude that $| [\mathscr {R}_z^\delta ]_{\ell k}^{ab} | \hbox {\,\,\char 054\,\,}2 (\mathfrak {d}^{{\mathrm {ref}}})^{-1} e^{-\gamma _\mathrm {CT}(\mathfrak {d}^{\mathrm {ref}}) |\ell - k| }$.

Next, we note that $(I + T^\mathrm {FR}\mathscr {R}^\delta _z)^{-1}T^\mathrm {FR}$ is a finite rank operator with $\big [ (I + T^\mathrm {FR}\mathscr {R}^\delta _z)^{-1} T^\mathrm {FR} \big ]_{\ell k}^{ab} = \big [ T^\mathrm {FR} (I + \mathscr {R}^\delta _zT^\mathrm {FR})^{-1} \big ]_{\ell k}^{ab} = 0$ for all $(\ell ,k) \not \in (\varLambda \cap B_R)^2$. Therefore,

$$\begin{aligned} \left| [ \mathscr {R}^\delta _z (I + T^\mathrm {FR}\mathscr {R}^\delta _z)^{-1} T^\mathrm {FR}\mathscr {R}^\delta _z ]_{\ell k}^{ab} \right|&\hbox {\,\,\char 054\,\,}\Vert (I + T^\mathrm {FR}\mathscr {R}^\delta _z)^{-1} T^\mathrm {FR} \Vert _{\mathrm {max}}\nonumber \\&\quad \sum _{\ell _1,\ell _2 \in \varLambda \cap B_R} e^{-\gamma _\mathrm {CT}(\mathfrak {d}^{\mathrm {ref}}) \left( |\ell -\ell _1| + |\ell _2 - k| \right) }\nonumber \\&\hbox {\,\,\char 054\,\,}C \Vert (I + T^\mathrm {FR}\mathscr {R}^\delta _z)^{-1} T^\mathrm {FR} \Vert _{\mathrm {max}} e^{-\gamma _\mathrm {CT}(\mathfrak {d}^{\mathrm {ref}}) \left( |\ell | + |k| \right) }. \end{aligned}$$

(7.2)

We obtain the claimed z-dependence after showing that $\Vert (I + T^\mathrm {FR}\mathscr {R}^\delta _z)^{-1} T^\mathrm {FR} \Vert _{\mathrm {max}} \lesssim (1+|z|)^2$. By (7.1), we have

$$\begin{aligned} (I + T^\mathrm {FR}\mathscr {R}^\delta _z)^{-1} T^\mathrm {FR}&= (T^{\mathrm {ref}} + T^\delta - z) (\mathscr {R}_z - \mathscr {R}_z^\delta ) (T^{\mathrm {ref}} + T^\delta - z) \nonumber \\&= -(T^{\mathrm {ref}} + T^\delta - z) \mathscr {R}_z T^\mathrm {FR} \mathscr {R}_z^\delta (T^{\mathrm {ref}} + T^\delta - z). \end{aligned}$$

(7.3)

Using the fact $[T^{\mathrm {FR}}]_{\ell _1 \ell _2}^{ab} = 0$ unless $\ell _1, \ell _2 \in \varLambda \cap B_R$, we have

$$\begin{aligned} |[ \mathscr {R}_z T^\mathrm {FR} \mathscr {R}_z^\delta ]_{\ell k}^{ab}| \hbox {\,\,\char 054\,\,}C \sum _{\ell _1,\ell _2 \in \varLambda \cap B_R} e^{-\gamma _\mathrm {CT}(\mathfrak {d}^{\mathrm {ref}})|\ell - \ell _1|} e^{-\gamma _\mathrm {CT}(\mathfrak {d}) |\ell _2 - k|} \hbox {\,\,\char 054\,\,}C e^{-\eta (|\ell | + |k| - 2R)} \end{aligned}$$

(7.4)

where $\eta {:=}\frac{1}{2} \min \{ \gamma _\mathrm {CT}(\mathfrak {d}), \gamma _\mathrm {CT}(\mathfrak {d}^{\mathrm {ref}}) \}$. Therefore, by combining (7.3) and (7.4), we obtain

$$\begin{aligned} \left| \big [ (I+T^\mathrm {FR}\mathscr {R}_z^\delta )^{-1} T^\mathrm {FR} \big ]_{\ell k}^{ab} \right|&\hbox {\,\,\char 054\,\,}C \bigg ( \sum _{\ell _1 \in \varLambda } \left( c_T + \delta + |z| \right) e^{-\eta (|\ell _1| - R)} \bigg )^2 \hbox {\,\,\char 054\,\,}C(1 + |z|)^2 \end{aligned}$$

and can thus conclude. $\square$

We will now show that $\mathcal {H}(u;\rho )$ and $\mathscr {L}(u;\rho )$ can be written as in the statement of Lemma 5 so that we may apply these improved resolvent estimates.

Since $\mathcal {H}(u;\rho )$ and $\mathcal {H}^{\mathrm {ref}}$ are defined on different spatial domains, we cannot directly compare the Hamiltonian with the corresponding reference operator. In order to alleviate this issue, we follow the arguments of Chen et al. (2019). Firstly, we shift the operator by a constant multiple of the identity cI and replace the contour and chemical potential by $\mathscr {C}_\mathfrak {o} + c$ and $\mu + c$, respectively, so that 0 is not encircled by $\mathscr {C}_\mathfrak {o} + c$. By changing variables in the integration, we can conclude that this shift does not affect the quantities defined by (2.3). We then add zero rows and columns so that the operators are defined on the same spatial domain $\varLambda \cup \varLambda ^{\mathrm {ref}}$. For example, for $\ell ,k \in \varLambda \cup \varLambda ^{\mathrm {ref}}$, if $\ell \in \varLambda ^{\mathrm {ref}} \setminus \varLambda$ or $k \in \varLambda ^{\mathrm {ref}} \setminus \varLambda$, we redefine $\widetilde{\mathcal {H}}(u;\rho )_{\ell k}^{ab} {:=}0$. This only affects the spectrum by adding zero as an eigenvalue of finite multiplicity and so, because 0 is not encircled by the contour $\mathscr {C}_\mathfrak {o}$, the value of (2.3) is unchanged. For full details, see Chen et al. (2019).

By replacing ${\mathcal {H}}(u;\rho )$ by $\widetilde{\mathcal {H}}(u;\rho )$ in (2.3) we obtain $\widetilde{O}_\ell (u;\rho ) {:=}O_\ell (u;\rho ) \chi _\varLambda (\ell )$. In particular, if we write $\widetilde{F}(u;\rho )$ as a function of electronic densities defined on $\varLambda \cup \varLambda ^{\mathrm {ref}}$, we find that the Jacobian of $\widetilde{F}(u;\rho )$ with respect to $\rho$, which we denote by $\widetilde{\mathscr {L}}(u;\rho )$, is obtained from ${\mathscr {L}}(u;\rho )$ by inserting finitely many additional zero rows and columns. Therefore, $I - \widetilde{\mathscr {L}}(u;\rho )$ is invertible with $\Vert \big ( I-\widetilde{\mathscr {L}}(u;\rho ) \big )^{-1} \Vert _{\ell ^2\rightarrow \ell ^2} \hbox {\,\,\char 054\,\,}\max \{ 1, \mathsf {d}_{\mathscr {L}}^{-1} \}$.

For the remainder of this paper, we consider the redefined quantities $\widetilde{\mathcal {H}}(u;\rho ), \widetilde{\mathcal {H}}^{\mathrm {ref}}$ and $\widetilde{\mathscr {L}}(u;\rho ), \widetilde{\mathscr {L}}^{\mathrm {ref}}$ and drop the tilde in the notation.

7.2 Proof of Theorem 2: Improved Locality Estimates

We now show that we can apply Lemma 5 to the Hamiltonian and the stability operators and thus conclude the proof of Theorem 2.

Lemma 6

(Perturbation of the Hamiltonian) Let $u \in \dot{\mathscr {W}}^{1,2}(\varLambda )$ satisfy $\mathbf (L)$ and $\rho$ be an associated electronic density satisfying (SC) and (FF).

Then, for all $\delta > 0$, there exist operators $\mathcal {H}^\delta = \mathcal {H}^\delta (u;\rho )$ and $\mathcal {H}^\mathrm {FR} = \mathcal {H}^\mathrm {FR}(u;\rho )$ such that

$$\begin{aligned} \mathcal {H}(u;\rho ) = \mathcal {H}^{\mathrm {ref}} + \mathcal {H}^\mathrm {FR} + \mathcal {H}^\delta \end{aligned}$$

(7.5)

where $\Vert \mathcal {H}^\delta \Vert _{\mathrm {F}} \hbox {\,\,\char 054\,\,}\delta$ and there exists an $R>0$ such that $[\mathcal {H}^\mathrm {FR}]_{\ell k}^{ab} = 0$ for all $(\ell ,k) \not \in (\varLambda \cap B_R)^2$.

Proof

Applying (Chen et al. 2019, Lemma 9), we may conclude that (7.5) holds for the linear Hamiltonian $\mathcal {H}^\mathrm {L}(u)$:

$$\begin{aligned} \mathcal {H}(u;\rho ) = \mathcal {H}^{\mathrm {ref}} + P(u) + Q(u) + D(\rho ) \end{aligned}$$

where $\Vert P(u)\Vert _{\mathrm {F}} \hbox {\,\,\char 054\,\,}\delta$, there exists an $R>0$ such that $Q(u)_{\ell k}^{ab} = 0$ for all $(\ell ,k) \not \in (\varLambda \cap B_R)^2$ and $D(\rho )$ is a diagonal operator with $D_{\ell \ell }^{aa} = v(\rho (\ell )) - v(\rho ^{\mathrm {ref}}(\ell ))$ for all $\ell \in \varLambda \setminus B_{R_{\mathrm {ref}}}$. After defining,

$$\begin{aligned} D^R(\rho )_{\ell \ell }^{aa} {:=}{\left\{ \begin{array}{ll} D_{\ell \ell }^{aa}(\rho ) &{}\text {if } \ell \in \varLambda \cap B_R \\ 0 &{}\text {otherwise,} \end{array}\right. } \quad \text {we have} \\ \lim _{R\rightarrow \infty } \Vert D(\rho ) - D^R(\rho )\Vert _{\ell ^2 \rightarrow \ell ^2} = \limsup _{|\ell | \rightarrow \infty } |v(\rho (\ell )) - v(\rho ^{\mathrm {ref}}(\ell ))| = 0. \end{aligned}$$

That is, $D(\rho )$ may be approximated with appropriate finite rank operators.

Remark 6

We remark here that if (FF) is not satisfied then $D(\rho )$ from the Proof of Lemma 6 is not compact and thus $\mathcal {H}(u;\rho ) - \mathcal {H}^{\mathrm {ref}}$ is also not compact. This means that, as noted at the end of Remark 4, the main techniques used in the Proof of Theorem 2 cannot be applied.

We now use Lemma 6 to show an analogous result for the stability operator.

Lemma 7

(Perturbation of the Stability Operator) Let $u \in \dot{\mathscr {W}}^{1,2}(\varLambda )$ satisfy $\mathbf (L)$ and $\rho$ be an associated electronic density satisfying (SC) and (FF). Then, for all $\delta > 0$, there exist operators $\mathscr {L}^\delta$ and $\mathscr {L}^\mathrm {FR}$ such that

$$\begin{aligned} \mathscr {L}(u;\rho ) = \mathscr {L}^{\mathrm {ref}} + \mathscr {L}^\mathrm {FR} + \mathscr {L}^\delta . \end{aligned}$$

where $\Vert \mathscr {L}^\delta \Vert _{\mathrm {F}} \hbox {\,\,\char 054\,\,}\delta$ and there exists an $R>0$ such that $[\mathscr {L}^\mathrm {FR}]_{\ell k} = 0$ for all $(\ell ,k) \not \in (\varLambda \cap B_R)^2$.

Proof

Using the notation from Lemma 6, we may apply Lemma 4 and obtain: for $z \in \mathscr {C}_f$,

$$\begin{aligned} \mathscr {R}_z(u;\rho )&= \big ( \mathcal {H}^{\mathrm {ref}} + \mathcal {H}^\mathrm {FR} - z \big )^{-1} + \big ( \mathscr {R}_z(u;\rho ) - \big ( \mathcal {H}^{\mathrm {ref}} + \mathcal {H}^\mathrm {FR} - z \big )^{-1} \big ) \nonumber \\&= \mathscr {R}^{\mathrm {ref}}_z - \mathscr {R}^{\mathrm {ref}}_z \big ( I + \mathcal {H}^\mathrm {FR} \mathscr {R}^{\mathrm {ref}}_z \big )^{-1} \mathcal {H}^\mathrm {FR} \mathscr {R}^{\mathrm {ref}}_z \nonumber \\&\quad - \mathscr {R}_z(u;\rho ) \mathcal {H}^\delta \big ( \mathcal {H}^{\mathrm {ref}} + \mathcal {H}^\mathrm {FR} - z \big )^{-1}\nonumber \\&{=:}\mathscr {R}^{\mathrm {ref}}_z + P_z + Q_z. \end{aligned}$$

(7.6)

Therefore, using (2.11), we have

$$\begin{aligned} \mathscr {L}(u;\rho )_{\ell k} - [\mathscr {L}^{\mathrm {ref}}]_{\ell k}= & {} \frac{1}{2\pi i} \oint _{\mathscr {C}_f} f(z-\mu ) \sum _{ab} \big [ P_z + Q_z \big ]_{\ell k}^{ab} \nonumber \\&\quad \big [ 2 \mathscr {R}_z^{\mathrm {ref}} + P_z + Q_z \big ]_{\ell k}^{ab} \mathrm {d}z\, v^\prime (\rho (k)) \end{aligned}$$

(7.7)

and so, using the fact that $\Vert 2 \mathscr {R}_z^{\mathrm {ref}} + P_z + Q_z \Vert _{\mathrm {max}} < \infty$ and arguing as in (7.2), we have

$$\begin{aligned} \Vert \mathscr {L}(u;\rho ) - \mathscr {L}^{\mathrm {ref}}\Vert _{\mathrm {F}}^2&\hbox {\,\,\char 054\,\,}C\sup _{z \in \mathscr {C}_f} \Vert P_z + Q_z\Vert _{\mathrm {F}}^2 \hbox {\,\,\char 054\,\,}C \sup _{z \in \mathscr {C}_f} (\Vert P_z\Vert _{\mathrm {F}} + \Vert Q_z\Vert _{\mathrm {F}})^2 \nonumber \\&\hbox {\,\,\char 054\,\,}C\sum _{\ell ,k \in \varLambda } e^{-\gamma _1 (|\ell | + |k|)}\nonumber \\&\quad + C \sum _{\ell , \ell _1, \ell _2, k \in \varLambda } \sum _{a_1,a_2} e^{-\gamma _2 r_{\ell \ell _1}} \big | [\mathcal {H}^\delta ]^{a_1a_2}_{\ell _1 \ell _2} \big |^2 e^{-\gamma _3 r_{\ell _2 k}} \nonumber \\&\hbox {\,\,\char 054\,\,}C\bigg ( \sum _{\ell \in \varLambda } e^{-\gamma _1 |\ell |} \bigg )^2 + C \Vert \mathcal {H}^\delta \Vert ^2_{\mathrm {F}} \hbox {\,\,\char 054\,\,}C \end{aligned}$$

(7.8)

where $\gamma _l {:=}2\gamma _{\mathrm {r}}(\mathfrak {d}_l)$ and $\gamma _{\mathrm {r}}$ is the constant (6.2) with appropriate choices of $\mathfrak {d}_l$ for $l=1,2,3$. Here, we have implicitly assumed that $\delta$ is sufficiently small such that for all $z \in \mathscr {C}_f$, we have $\mathrm {dist} (z, \sigma ( \mathcal {H}^{\mathrm {ref}} + \mathcal {H}^\mathrm {FR} ) ) \hbox {\,\,\char 062\,\,}\mathfrak {d}_3 > 0$. Therefore, by applying (7.8), for sufficiently large R, we can define

$$\begin{aligned} {[}\mathscr {L}^\mathrm {FR} ]_{\ell k} {:=}{\left\{ \begin{array}{ll} [\mathscr {L}(u;\rho ) - \mathscr {L}^{\mathrm {ref}}]_{\ell k} &{} \quad \text {if }\ell , k \in \varLambda \cap B_R \\ 0 &{} \quad \text {otherwise} \end{array}\right. } \end{aligned}$$

and conclude $\Vert \mathscr {L}(u;\rho ) - \mathscr {L}^{\mathrm {ref}} - \mathscr {L}^\mathrm {FR}\Vert _{\mathrm {F}} \hbox {\,\,\char 054\,\,}\delta$.

For fixed $(u,\rho )$ and $\delta > 0$ sufficiently small, we fix operators $\mathcal {H}^\mathrm {FR}, \mathcal {H}^\delta$ and $\mathscr {L}^\mathrm {FR}, \mathscr {L}^\delta$ as in Lemmas 6 and 7 and apply Lemma 5 to obtain: for $z \in \mathbb {C}$ with $\mathrm {dist}(z,\sigma (\mathcal {H}(u;\rho ))) \hbox {\,\,\char 062\,\,}\mathfrak {d} > 0$ and $\mathrm {dist}(z,\sigma (\mathcal {H}^{\mathrm {ref}})) - \delta \hbox {\,\,\char 062\,\,}\mathfrak {d}^{\mathrm {ref}} > 0$ we have

$$\begin{aligned} \big | \mathscr {R}_z(u;\rho )_{\ell k}^{ab} \big |&\hbox {\,\,\char 054\,\,}C_{\ell k} e^{-\gamma _{\mathrm {r}}(\mathfrak {d}^{\mathrm {ref}}) |\ell - k| } \quad \text {where} \nonumber \\ {C}_{\ell k}&{:=}2(\mathfrak {d}^{\mathrm {ref}})^{-1} + C (1 + |z|)^2 e^{-\gamma _{\mathrm {r}}(\mathfrak {d}^{\mathrm {ref}}) \left( |\ell | + |k| - |\ell - k|\right) } \end{aligned}$$

(7.9)

and $\gamma _{\mathrm {r}}(\mathfrak {d}^{\mathrm {ref}}) = c \min \{ 1, \mathfrak {d}^{\mathrm {ref}} \}$ is the constant from (6.2). We will apply (7.9) with $z \in \mathscr {C}_f$ (with $\mathfrak {d} = \mathsf {d}_f$ and $\mathfrak {d}^{\mathrm {ref}} = \mathsf {d}_f^{\mathrm {ref}} {:=}\min _{z \in \mathscr {C}_f} \mathrm {dist}( z,\sigma (\mathcal {H}^{\mathrm {ref}}) ) - \delta$) and $z \in \mathscr {C}_\mathfrak {o}$ (with $\mathfrak {d} = \mathsf {d}_\mathfrak {o}$ and $\mathfrak {d}^{\mathrm {ref}} = \mathsf {d}_\mathfrak {o}^{\mathrm {ref}} {:=}\min _{z \in \mathscr {C}_\mathfrak {o}} \mathrm {dist}(z, \sigma (\mathcal {H}^{\mathrm {ref}})) - \delta$).

By (2.11), we have $\big | \mathscr {L}(u;\rho )_{\ell k}\big | \lesssim {C}_{\ell k}^2 e^{-2\gamma _{\mathrm {r}}(\mathsf {d}_f^{\mathrm {ref}}) |\ell - k| }$ and so, applying Lemma 5 again but now with T replaced with $\mathscr {L}(u;\rho )$ (and with $z = 1$), we obtain,

$$\begin{aligned} \left| \left[ (I - \mathscr {L}(u;\rho ))^{-1} \right] _{\ell k}^{ab} \right| \hbox {\,\,\char 054\,\,}\widetilde{C}_{\ell k} e^{-\gamma ^{\mathrm {ref}}_{\mathrm {s}} r_{\ell k} } \quad \text {where} \nonumber \\ \widetilde{C}_{\ell k} {:=}2 (\mathsf {d}_{\mathscr {L}}^{\mathrm {ref}})^{-1} + C e^{-\gamma ^{\mathrm {ref}}_{\mathrm {s}} \left( |\ell | + |k| - |\ell - k|\right) } \end{aligned}$$

(7.10)

and $\gamma ^{\mathrm {ref}}_{\mathrm {s}} {:=}c_1 \gamma _{\mathrm {r}}(\mathsf {d}^{\mathrm {ref}}_f) \min \{1,C_{\ell k}^{-2}\gamma _{\mathrm {r}}(\mathsf {d}^{\mathrm {ref}}_f)^d \mathsf {d}^{\mathrm {ref}}_{\mathscr {L}}\}$ with $\mathsf {d}^{\mathrm {ref}}_{\mathscr {L}} {:=}\mathrm {dist}(1,\sigma (\mathscr {L}^{\mathrm {ref}})) - \delta$. By expanding $\gamma _{\mathrm {r}}(\mathsf {d}_f^{{\mathrm {ref}}})$ in terms of $\mathsf {d}_f^{\mathrm {ref}}$, we obtain $\gamma _{\mathrm {s}}^{\mathrm {ref}} = c_2 \min \big \{ 1, \mathsf {d}^{\mathrm {ref}}_f, C_{\ell k}^{-2} \min \{ 1, (\mathsf {d}_f^{\mathrm {ref}})^{d+1} \} \mathsf {d}^{{\mathrm {ref}}}_{\mathscr {L}} \big \}$.

Proof of Theorem 2

The arguments in the Proof of Theorem 1 can be applied with the resolvent estimates of (6.2) and (6.4) replaced with the corresponding improved estimates (7.9) and (7.10). This means that the exponents $\gamma _{\mathrm {r}}(\mathfrak {d})$ and $\gamma _{\mathrm {s}}$ can be replaced with the improved exponents $\gamma _{\mathrm {r}}(\mathfrak {d}^{\mathrm {ref}})$ and $\gamma ^{\mathrm {ref}}_{\mathrm {s}}$, respectively, and the pre-factors can be replaced with constants that depend on the atomic sites. These constants decay exponentially to the constants in the defect-free case as the subsystem moves away from the defect core together. This can be seen by noting that $C_{\ell k} \rightarrow 2(\mathfrak {d}^{\mathrm {ref}})^{-1}$ (where $C_{\ell k}$ is the constant from (7.9)) and $\widetilde{C}_{\ell k} \rightarrow 2(\mathsf {d}_{\mathscr {L}}^{\mathrm {ref}})^{-1}$ (where $\widetilde{C}_{\ell k}$ is the constant from (7.10)) as $|\ell | + |k| - |\ell -k| \rightarrow \infty$ with exponential rates. See (Chen et al. 2019, (4.21)–(4.23)) for the analogous argument in the linear case that can be readily adapted to the setting we consider here. $\square$

Appendix A: Band Structure of \(\mathcal {H}^{\mathrm {ref}}\) the Reference Hamiltonian and \(\mathscr {L}^{\mathrm {ref}}\) the Reference Stability Operator

Recall that the unit cell $\varGamma \subset \varLambda ^{\mathrm {ref}}$ is finite and satisfies $\varLambda ^{\mathrm {ref}} = \bigcup _{\gamma \in \mathbb Z^d} (\varGamma + \mathsf {A}\gamma )$ and $\varGamma + \mathsf {A}\gamma$ pairwise disjoint for each $\gamma \in \mathbb Z^d$. Suppose $\varGamma ^\star \subset \mathbb R^d$ is a bounded connected domain containing the origin and such that $\mathbb R^d = \bigcup _{\eta \in \mathbb Z^d} (\varGamma ^\star + 2\pi \mathsf {A}^{-\mathrm {T}}\eta )$ and the $\varGamma ^\star + 2\pi \mathsf {A}^{-\mathrm {T}}\eta$ are disjoint. Therefore, for each $\xi \in \mathbb R^d$, there exist unique $\xi _0 \in \varGamma ^\star$ and $\eta \in \mathbb Z^d$ such that $\xi = \xi _0 + 2\pi \mathsf {A}^{-\mathrm {T}}\eta$, and since $\mathsf {A}\gamma \cdot \mathsf {A}^{-\mathrm {T}}\eta = \gamma \cdot \eta$, we have $e^{-i\mathsf {A}\gamma \cdot \xi } = e^{-i\mathsf {A}\gamma \cdot \xi _0}$ for $\gamma \in \mathbb Z^d$.

Let us define the unitary operator $U :\ell ^2(\varLambda ^{\mathrm {ref}} \times \{1,\dots ,N_{\mathrm{b}}\} ) \rightarrow L^2\left( \varGamma ^\star ;\ell ^2(\varGamma \times \{1,\dots ,N_{\mathrm{b}}\}) \right)$ by

$$\begin{aligned} (U\psi )_\xi (\ell ;a) = \sum _{\gamma \in \mathbb Z^d} \psi (\ell + \mathsf {A}\gamma ;a) e^{-i(\ell + \mathsf {A}\gamma )\cdot \xi }. \end{aligned}$$

Here, $L^2\left( \varGamma ^\star ;\ell ^2(\varGamma \times \{1,\dots ,N_{\mathrm{b}}\})\right)$ is a Hilbert space with inner product

$$\begin{aligned} \langle \varPsi ,\varPhi \rangle _{L^2(\varGamma ^\star ;\ell ^2)} {:=}\frac{1}{|\varGamma |} \int _{\varGamma ^\star } \langle \varPsi _\xi ,\varPhi _\xi \rangle _\ell \mathrm {d}\xi = \frac{1}{|\varGamma |} \sum _{\ell \in \varGamma } \sum _{1\hbox {\,\,\char 054\,\,}a \hbox {\,\,\char 054\,\,}N_{\mathrm{b}}} \int _{\varGamma ^\star } \varPsi _\xi (\ell ;a)\overline{\varPhi _\xi (\ell ;a)} \mathrm {d}\xi . \end{aligned}$$

A simple calculation reveals that $( U\mathcal {H}^{\mathrm {ref}} \psi )_\xi = \mathcal {H}^{\mathrm {ref}}_\xi ( U \psi )_\xi$ where

$$\begin{aligned} {[} \mathcal {H}^{\mathrm {ref}}_\xi ]_{\ell k}^{ab} = \sum _{\gamma \in \mathbb Z^d} h_{\ell k}^{ab}(\ell - k + \mathsf {A}\gamma ) e^{-i(\ell - k + \mathsf {A}\gamma )\cdot \xi } + \delta _{\ell k}\delta _{ab} \sum _{\gamma \in \mathbb Z^d} v(\rho ^{\mathrm {ref}}(\ell )) e^{-i\mathsf {A}\gamma \cdot \xi } . \end{aligned}$$

Letting $\lambda _n(\xi )$ be the ordered eigenvalues of $\mathcal {H}^{\mathrm {ref}}_\xi$ for $\xi \in \varGamma ^\star$ and $n = 1,\dots ,N_{\mathrm{b}}\cdot \#\varGamma$, we can use the fact that U is unitary to conclude that,

$$\begin{aligned} \sigma (\mathcal {H}^{\mathrm {ref}}) = \bigcup _n \bigcup _{\xi \in \varGamma ^\star } \lambda _n(\xi ). \end{aligned}$$

Since $\lambda _n :\varGamma ^\star \rightarrow \mathbb R$ are continuous, we may conclude that $\sigma (\mathcal {H}^{\mathrm {ref}})$ is composed of finitely many spectral bands.

Similarly,

$$\begin{aligned} \sigma (\mathscr {L}^{\mathrm {ref}})= & {} \bigcup _{\xi \in \varGamma ^\star } \sigma (\mathscr {L}^{\mathrm {ref}}_\xi ), \quad \text {where} \\ \big [\mathscr {L}^{\mathrm {ref}}_\xi \big ]_{\ell k}= & {} \sum _{\gamma \in \mathbb Z^d} \frac{1}{2\pi i} \oint _{\mathscr {C}_f} f(z-\mu ) \sum _{a,b = 1}^{N_{\mathrm{b}}} \left( [\mathscr {R}_z^{\mathrm {ref}}]_{\ell +\mathsf {A}\gamma ,k}^{ab} \right) ^2 \mathrm {d}z\, v^\prime (\rho (k)) e^{-i(\ell - k + \mathsf {A}\gamma )\cdot \xi }. \end{aligned}$$