Optimal artificial boundary condition for random elliptic media
Jianfeng Lu Felix Otto




[2010]35B27; 65N99

[
    Jeffrey Hatley
    ?; revised ?; accepted ?. - To be entered by editorial office
=================================================================




We are given a uniformly elliptic coefficient field 
that we regard as a realization of a stationary and finite-range (say, range unity) ensemble of coefficient fields.
Given a (deterministic) right-hand-side supported in a ball of size ℓ≫ 1 and of vanishing average, 
we are interested in an algorithm to compute the (gradient of the) solution near the origin, 
just using the knowledge of the (given realization of the) coefficient field in some large box of size L≫ℓ. 
More precisely, we are interested in the most seamless (artificial) boundary condition
on the boundary of the computational domain of size L. 



Motivated by the recently introduced multipole expansion in random media, we propose an algorithm. 
We rigorously establish an error estimate (on the level of the gradient) in terms of L≫ℓ≫ 1, 
using recent results in quantitative stochastic homogenization. More precisely,
our error estimate has an a priori and an a posteriori aspect: With a priori overwhelming probability,
the (random) prefactor can be bounded by a constant that is computable without much further effort, 
on the basis of the given realization in the box of size L.



We also rigorously establish that the order of the error estimate in both L and ℓ is optimal,

where in this paper we focus on the case of d=2.
This amounts to a lower bound on the variance of the quantity of interest when conditioned on
the coefficients inside the computational domain, 
and relies on the deterministic insight that a sensitivity analysis wrt a defect 
commutes with (stochastic) homogenization.
Finally, we carry out numerical experiments that show that this optimal convergence rate already sets in at only moderately large L,
and that more naive boundary conditions perform worse both in terms of rate and prefactor.








§ INTRODUCTION AND MAIN RESULTS


Let the dimension d≥ 2 and the ellipticity ratio λ>0 be fixed.
We will be considering symmetric tensor fields a on d-dimensional space ℝ^d
that are uniformly elliptic:


    |ξ|^2≥ξ· a(x)ξ≥λ|ξ|^2  ,


where wlog we've set the upper bound to unity. Symmetry is notationally convenient at a few places, but by no
means essential; while we use scalar notation and language, all results hold for systems.
For some localized rhs g, say near the origin, 
we are interested in the decaying (ie Lax-Milgram) whole-space solution u of


    -∇· a∇ u=∇· g.


More precisely, we are interested in -∇ u near the origin, say,
at the origin: -∇ u(0). In the language of electrostatics, we
are interested in the electric field generated by the neutral and
localized charge distribution ∇· g.  We pose the question
to which precision -∇ u(0) can be inferred without solving a
PDE in whole-space. Let us denote by ℓ the (linear) size of the
support of g.  In case of constant coefficients a_h, the
explicit fundamental function G_h allows to reduce the determination
of -∇ u(0) to the evaluation of an integral over B_ℓ, the
centered ball of radius ℓ.  



In our case of variable coefficients, we ask the
question of whether one can do better than solving a boundary value
problem with homogeneous boundary data, say the Dirichlet problem


    -∇· a∇ u_0=∇· g  Q_L,  
    u_0=0  ∂ Q_L


on the centered cube Q_L:=(-L,L)^d for some large scale L≥ℓ (where we take cubes instead
of balls for computational convenience). Under the assumption
that ℓ is the only scale of g in the sense that there exists a function ĝ such that


    g(x)=ĝ(x/ℓ)   ĝ   B_1  ∇̂ĝ ,


one expects and we experimentally show in Section <ref>
that the approximation (<ref>) is no better than what generically holds in the constant-coefficient case, namely
∇(u_0-u)(0)=O((ℓ/L)^d).
More precisely, we ask the question whether one can do better without knowing the coefficients 
outside Q_2L, which is hopeless when a has no further structure. In this paper, we thus consider the
case when a comes from a stationary finite-range ensemble ⟨·⟩ of uniformly elliptic coefficient fields;
wlog we assume the range to be unity. We recall that stationary means that a and its shifted version a(·+z) for
any shift vector z∈ℝ^d have the same distribution under ⟨·⟩;
unit range means that for two subset D,D'⊂ℝ^d with dist(D,D')>1,
the restrictions a_|D and a_|D' of the coefficient field a are independent under ⟨·⟩.






 §.§ Approximation algorithm and error bound


Loosely speaking, our first main result Theorem <ref> states that there is an algorithm, the outcome of which is u^(L),
that



  * only involves knowing the realization a restricted to Q_2L,

  * next to the solution of a Dirichlet problem on Q_L only requires
the solution of 2 (in d=2) respectively 9 (in d=3) further Dirichlet problems on Q_2L, 


  * improves upon (<ref>) by (almost) a factor L^-d/2, 
albeit with a random prefactor,

  * with overwhelming probability in L≫ 1, this prefactor is
  dominated by a constant that can be computed at the cost of further
  2 (for d=2) respectively 27 (for d=3) (constant-coefficient) Dirichlet problems on
  Q_2L.



Note that Theorem <ref> is a mixture of an a priori result, namely the
probabilistic estimate on when the approach is successful at all, and
an a posteriori result, namely the domination of the prefactor by the computable quantity r_*^(L).
Loosely speaking, the second main result, Theorem <ref>, states that in terms of
scaling (both in L and ℓ), there is no better algorithm since
for a relevant class of ensembles ⟨·⟩, the square
root of the variance of ∇ u(0) conditioned on a_|Q_2L is
of order (ℓ/L)^dL^-d/2.  The argument shows
that the factor L^-d/2 is of CLT-type and loosely speaking
arises as the inverse of the square root of the volume of the
neighboring “annulus” Q_4L-Q_2L.



This paper only discusses the algorithm that gives the (near) optimal result in case of d=2;
the optimal algorithm for the case of d=3 (and higher d in general)
would require a refinement, namely the second-order corrector, but no new concepts. 
More precisely, the theory of dipoles developed in <cit.>
would have to be replaced by its systematic generalization to multipoles (in particular
quadrupoles) developed in <cit.>, relying on second-order correctors. 
It is the exponent β 
on 1/L in (<ref>), which can be taken arbitrarily close to 1 (=d/2 for d=2), 
that provides the near-optimal CLT improvement over the homogeneous boundary value problem (<ref>).



Since we do not assume (Hölder-) continuity of the realization a, we do not have access
to the point evaluation (at the origin) of the gradient. Hence in all statements, pointwise control
is replaced by L^2-control over an order-one ball.



Correctors ϕ_i. Not surprisingly, our algorithm makes use of the correctors ϕ_i, which
for every coordinate direction i=1,⋯,d provide a-harmonic coordinates x_i+ϕ_i by satisfying


    -∇· a(e_i+∇ϕ_i)=0,


where e_i denotes the unit vector in direction i.
According to the classical qualitative theory of stochastic homogenization (by “qualitative” theory 
we mean the one only relying on ergodicity and stationarity of the ensemble ⟨·⟩), 
for almost every realization a, a corrector ϕ_i of sublinear growth, that is,


    lim_r↑∞1/r(_B_r(ϕ_i-_B_rϕ_i)^2)^1/2=0


can be constructed. Here and in the sequel, _B_r denotes the average over the ball B_r (of radius r centered at the origin).
Moreover, again almost surely, the homogenized coefficients a_h may be inferred from


    a_he_i=lim_L↑∞_B_Lq_i     q_i:=a(e_i+∇ϕ_i).


Hence a naive guess would be to replace the approximation u_0 defined through (<ref>) by
the approximation u_I defined through


    -∇· a∇ u_I=∇· g  Q_L,  
    u_I=ũ_h  ∂ Q_L,


where ũ_h is the decaying solution of the homogenized problem in the whole space, that is,


    -∇· a_h∇ũ_h=∇· g.


This can be seen to generically yield no improved scaling of the error
in L (ie at fixed ℓ); it only improves the scaling of the error
when ℓ≫ 1.  Incidentally, it would not fall in the class of the
algorithms we consider, since inferring the homogenized coefficient
a_h requires solving the whole-space problem (<ref>).  In the
context of multiscale method, this is the approach taken in
<cit.> (with additional steps to approximate a_h).



In both periodic and random homogenization, it is known that the so-called two-scale expansion


    (1+ϕ_i∂_i)ũ_h,


where we use Einstein's convention of summation over repeated indices,
provides a better approximation to u than ũ_h itself; in particular,
this approximation is necessary to get closeness of the gradients (in the regime ℓ≫ 1).
Hence a second attempt would be to replace the approximation u_I defined through (<ref>) by
the approximation u_II defined through


    -∇· a∇ u_II=∇· g  Q_L,  
    u_II=(1+ϕ_i∂_i)ũ_h  ∂ Q_L,


As our numerical experiments in Section <ref> show, this generically yields no improved scaling of the error in L.



Dipoles. The problem with all three approaches, (<ref>), (<ref>), and (<ref>),
is that as soon as ℓ≪ L,
the far field of ũ_h generically has the wrong dipole behavior. This phenomenon
was observed in <cit.>, where the right-hand side g in (<ref>)
was replaced by g_i(e_i+∂_iϕ) in order for the gradient of 
the two scale-expansion ∇(1+ϕ_i∂_i)ũ_h to be O( (ℓL)^d (1L)^β)-close to ∇ u. 

This is the right strategy for concentrated g, ie for ℓ∼ 1. In order
to also treat a more spread rhs, ie ℓ≫ 1, we hold on to ũ_h but correct it through
an a_h-dipole δ u_h coming from the first moments of g_i∂_iϕ. In formulas,
we pass to 


    u_h:=ũ_h+(∫∇ϕ_i· g)∂_iG_h,


where G_h is the fundamental solution for -∇· a_h∇.
Hence our more educated ansatz is to replace the approximation u_II defined through (<ref>) by
the approximation u_III defined through


    -∇· a∇ u_III=∇· g  Q_L,  
    u_III=(1+ϕ_i∂_i)u_h  ∂ Q_L.


Corollary <ref> shows that this approximation indeed reduces the (generic) error of (<ref>)
by a factor with the desired L-scaling 1/L^β. Corollary <ref> relies on Lemma <ref>,
which is a minor modification of <cit.>.



Flux correctors σ_i. As can be seen from Lemma <ref>, the prefactor comes in form of r_*^β for some 
length scale r_*. This length scale has the interpretation that for larger scales r≥ r_*,
the quantified sublinear growth of the correctors ϕ_i sets in, cf (<ref>). 
The sublinear growth is quantified through the exponent β
(with β close to 0 meaning almost linear growth, and β close to 1 meaning
almost no growth). However, for quantitative results like Lemma <ref>, which closely follows
<cit.>, itself inspired by <cit.>, it is not sufficient
to monitor just ϕ_i. In fact, the harmonic vector field e_i+∇ϕ_i (the electric field
in the language of electrostatics) is not just a closed 1-form; but through the flux 
q_i:=a(e_i+∇ϕ_i) (the electric current in the language of electrostatics) it provides a 
closed (d-1)-form. Hence there is not just the 0-form (a scalar potential) x_i+ϕ_i, 
or rather its correction ϕ_i,
but there naturally is also a (d-2)-form (a vector potential in the 3-d language, or a stream function
in the 2-d language), or rather its correction σ_i, which we can write as a skew-symmetric
tensor field σ_i={σ_ijk}_j,k=1,⋯,d. In view of (<ref>),
this correction should satisfy


    q_i=a_he_i+∇·σ_i,


where (∇·σ_i)_j:=∂_kσ_ijk. Note that by skew symmetry of σ_i we have ∇·∇·σ_i=0
so that (<ref>) contains the familiar (<ref>), as it implies ∇· q_i = 0.
Clearly, (<ref>) determines σ_i only up to a (d-3)-form, ie the freedom of the choice of a gauge.
A particularly simple choice of gauge is


    -Δσ_ijk=∂_jq_ik-∂_kq_ij.


This skew-symmetric field is not uncommon in periodic homogenization; in 
qualitative stochastic homogenization <cit.> 
it has been shown to almost surely exist with sublinear growth:


    lim_r↑∞1/r(_B_r|σ_i-_B_rσ_i|^2)^1/2=0.
 



Radius r_*. Loosely speaking, starting from the length scales r
at which the lhs expression (<ref>) and (<ref>) drop below a threshold only depending on d
and λ, the operator -∇· a∇ inherits the
regularity theory of -∇· a_h∇, both for Schauder
theory on the C^1,α level (then the threshold depends in
addition on 0<α<1) <cit.>, and for the
Calderon-Zygmund theory on the Ḣ^1,p level (then the threshold
depends in addition on 1<p<∞) <cit.>. This type of
theory had been developed by Avellaneda & Lin <cit.> for the
periodic case; Armstrong & Smart <cit.> were first to extent this
to the random case.  As mentioned, r_* quantifies (<ref>) and
(<ref>); it is defined as the length scale starting from which
β-sublinear growth of the scalar and vector correctors kick in,
that is,


    1/r(_B_r|(ϕ,σ)-_B_r(ϕ,σ)|^2)^1/2≤(r_*/r)^β   r≥ r_*,


where (ϕ,σ) is the collection of all components {ϕ_i,σ_ijk}_i,j,k.
While the form (<ref>) is more natural, it is only seemingly weaker than


    1/r(_B_r|(ϕ,σ)-_B_r_*(ϕ,σ)|^2)^1/2≤
    C(β)(r_*/r)^β   r≥ r_*,


as we shall show at the beginning of the proof of Lemma <ref>. 
In <cit.> it is shown that r_* satisfies (optimal) stretched exponential bounds
even under weak correlation decay of ⟨·⟩.



Clearly, this notion of r_* singles out the origin (it can and will be defined for other bounds
r_*(y), making it a stationary random field). The origin plays a special role in our analysis in 
two ways: It is where we want to monitor the error (on the level of the gradient) and where
the rhs g is concentrated (on scale ℓ). In view of the above-mentioned
C^1,α-regularity theory that kicks in (only) from scales r_* onwards, it is not
surprising that we can localize the error (on the level of the gradients) only to scales
R≥ r_*, see Proposition <ref> (which in Theorem <ref> is expressed in terms of the proxy r_*^(L)).
For the same reason, we need the condition ℓ≥ r_* in Proposition <ref>. 



Algorithm. The “algorithm” (<ref>) is not admissible, since it involves solving 
the d whole-space problems (<ref>).
The natural idea is to replace these whole-space problems by Dirichlet problems on Q_2L;
for reasons that will become clearer later, we do the same for the d×d(d-1)/2
whole-space problems (<ref>) (while having constant coefficients they feature an extended supported rhs). 

For any coordinate direction i=1,⋯,d, let the function _i 
and the skew-symmetric tensor field _i={_ijk}_j,k=1,⋯,d 
be determined through


    -∇· a(e_i+∇_i)=0  Q_2L,     _i=0  ∂ Q_2L,
    
    -Δ_ijk=∂_j_ik-∂_k_ij  Q_2L,     _ijk=0  ∂ Q_2L,


where we have set for abbreviation _i:=a(e_i+∇_i). While
an easy calculation shows that (<ref>) & (<ref>) imply
Δ(q_i-∇·σ_i)=0, this does in general not yield
_i-_Q_2L_i =∇·_i. The latter would be
automatic in case of periodic boundary conditions, in which case we
would replace the homogenized coefficient a_h, cf (<ref>),
by _h e_i =_Q_2L_i. In our (more ambitious) case of
Dirichlet boundary conditions, we pick a mask ω̂ of an
averaging function with


    ω(x)=1/L^dω̂(x/L)   ω̂   Q_1, ∫ω̂=1, and ∇̂ω̂ ,


and set


    _he_i:=∫ω_i;


it is a consequence of (<ref>) in Lemma <ref> that _h is elliptic.
We now make the corresponding changes on the level of ũ_h and u_h:
We substitute ũ_h defined in (<ref>) by the decaying solution _h of


    -∇·_h∇_h=∇· g


and u_h defined in (<ref>) by _h defined through


    _h:=_h+(∫∇_i· g)∂_i_h,


where _h denotes the fundamental solution of -∇·_h∇.



As mentioned above, Theorem <ref> is a mixture of an a priori and
an a posteriori result: Through the scale L_0, which only depends on the dimension d<∞, 
the ellipticity ratio λ>0, the sublinear growth exponent β<1, 
and the stretched exponential exponent s<2(1-β),
Theorem <ref> provides an a priori estimate on the
probability that the random pick of a realization a is so bad that
the algorithm fails. Through (<ref>), which characterizes the
scale _* in a computable fashion, it provides an a posteriori
estimate on the constant (_*)^β in the
1/L^β-improvement of the error estimate (<ref>).  We
think of (<ref>) as an a posteriori estimate, since it
relies on an auxiliary computation based on the given realization a.




Let ⟨·⟩ be a stationary ensemble of uniformly elliptic coefficient fields a,
cf (<ref>), that is of unit range. 
Then for any exponents β∈(0,1) and s∈(0,2(1-β)), there exists a scale L_0=L_0(d,λ,β,s)
so that for any scale L≥ L_0, with probability 1-exp(-(L/L_0)^s) 
a realization a has the following property: 



Let {_i}_i=1,⋯,d, _i={_ijk}_i,j,k=1,⋯,d,
and _h be defined through (<ref>), (<ref>), and (<ref>).
Suppose that (,) is quantitatively sublinear 
in the sense of that there exists a scale 1≤_*≤ L with


    1/r(_Q_r|(,)|^2)^1/2≤(_*/r)^β   2L≥ r≥_*.




Let g be of the form (<ref>) for some scale ℓ∈[_*,L] and mask ĝ,
and let u be the decaying solution of (<ref>).
Consider the solution  of


    -∇· a∇=∇· g  Q_L,
      =(1+_i∂_i)_h  ∂ Q_L,


where _h is defined through  (<ref>) and (<ref>).
Then we have 


    (_B_R|∇(-u)|^2)^1/2≤ C(ℓ/L)^d(_*/L)^β   L≥ R≥_*


with a constant C that only depends on d, λ, β, and the Hölder norm
of ∇̂ĝ as well as the supremum norm of ∇̂ω̂.




Theorem <ref> has three ingredients, the deterministic a priori
error estimate provided by Proposition <ref>, the stochastic
ingredient Lemma <ref>, and the deterministic Lemma <ref>, that
allows to pass from an a priori to an a posteriori error estimate. We
say that Proposition <ref> provides an a priori and deterministic error estimate since it is formulated in terms of
r_* characterizing the sublinear growth of the augmented corrector
(ϕ,σ), cf (<ref>).  In addition, it starts from a given
uniformly elliptic coefficient field a, which might but does not have
to be a realization under ⟨·⟩.  The only assumption
on the Dirichlet proxy (,) is that it is well-behaved on
the large scale 2L, cf (<ref>), but not necessarily on smaller
scales as in (<ref>).


Consider a given λ-uniformly elliptic coefficient field a on ℝ^d,
cf (<ref>). Suppose that there exists a tensor a_h and, for i=1,⋯,d, a scalar field 
ϕ_i and a skew-symmetric tensor field σ_i such that (<ref>) holds.
Suppose that for given β∈(0,1) there exists a radius r_* such that (<ref>) holds.



For given L≥ r_*, let _i, _i, and _h be defined through
(<ref>), (<ref>), and (<ref>). We assume that


    1/2L(_Q_2L|(,)|^2)^1/2≤(r_*/2L)^β.




Given ĝ and ℓ∈[r_*,L] let u and  be defined as in Theorem <ref>.
Then we have 


    (_B_R|∇(-u)|^2)^1/2≤ C(ℓ/L)^d(r_*/L)^β   L≥ R≥ r_*,


where the constant C is of the same type as in Theorem <ref>.





The following lemma is the key ingredient for Proposition <ref>,
it shows that indeed the multipole has to be corrected in the sense of
(<ref>). Its proof essentially follows <cit.>.






Let a, a_h, ϕ_i, σ_i, β, and r_* be as in Proposition <ref>.
Let g be of the form (<ref>) for some ℓ≥ r_* and ĝ,
and let u be the decaying solution of (<ref>).
Let u_h be defined through  (<ref>) and (<ref>).
Then we have


    (_B_R^c|∇(u-(1+ϕ_i∂_i)u_h)|^2)^1/2≤ C(ℓ/R)^d(r_*/R)^β   R≥ r_*,


where ϕ_i is normalized through


    _B_r_*ϕ_i=0.


Here C is a constant of the same type as in Theorem <ref>. Furthermore
_B_R^c is the abbreviation of 1/R^d∫_B_R^c, where B_R^c
is the complement of B_R.





Equipped with Lemma <ref>, we may assess the effect of a computational domain Q_L
endowed with the Dirichlet conditions given by (1+ϕ_i∂_i)u_h, cf (<ref>).


Let a, a_h, ϕ_i, σ_i, β, and r_* be as in Proposition <ref>.
Let g be of the form (<ref>) for some ℓ≥ r_* and ĝ,
and let u be the decaying solution of (<ref>).



For L≥ℓ let u_III be the solution of the Dirichlet problem (<ref>).
Then we have


    (_B_R|∇(u_III-u)|^2)^1/2≤ C(ℓ/L)^d(r_*/L)^β   L≥ R≥ r_*,


where C is of the type as in Theorem <ref>.




In order to pass from Corollary <ref> to Proposition <ref>, that is, from u_III to u^(L),
we need to replace (ϕ,σ) by the computable (,)
in order to pass from u_III to .


Let a, a_h, ϕ_i, and σ_i be as in Proposition <ref>.



Suppose that for given β∈(0,1) there exists a radius r_* such that (<ref>) holds in the weaker
form of


    1/2L(_Q_2L|(ϕ,σ)-_Q_2L(ϕ,σ)|^2)^1/2≤(r_*/2L)^β,


for given L≥ r_*. Let _i, _i, and _h be defined through
(<ref>), (<ref>), and (<ref>). Suppose that (<ref>) holds.
Then we have 


    (_Q_3/2L|∇(-ϕ)|^2)^1/2   ≤ C(r_*/L)^β,
    (_Q_L|∇(-σ)|^2)^1/2   ≤ C(r_*/L)^β,
      |_h-a_h|   ≤ C(r_*/L)^β,


where C is of the same type as in Theorem <ref>.




As mentioned, the following Lemma <ref> allows to pass from the deterministic a priori estimate
of Proposition <ref> to the deterministic a posteriori estimate of Theorem <ref>.
It shows that sublinear growth of (ϕ,σ) on scales R≥ L with (some) pre-factor (_*)^β, cf (<ref>),
and sublinear growth of the Dirichlet proxy (,) on scales 2L≥ r≥_*, cf (<ref>),
implies sublinear growth of (ϕ,σ) on all scales r≥ r_*^(L), that is, r_*≲_*.


Let a, a_h, ϕ_i, and σ_i be as in Proposition <ref>.



For given L≥ r_*, let _h, _i, and _i be defined through
(<ref>), (<ref>), and (<ref>). We assume that for a given β∈(0,1)
there exists a scale _*≤ L such that we have on the one side


    1/R(_B_R|(ϕ,σ)-_B_R(ϕ,σ)|^2)^1/2≤(_*/R)^β   R≥ L


and on the other side


    1/r(_B_r|(,)|^2)^1/2≤(_*/r)^β    2L≥ r≥_*.


Then we have


    1/r(_B_r|(ϕ,σ)-_B_r(ϕ,σ)|^2)^1/2≤ C(_*/r)^β   r≥_*


with the constant C being of the same type as in Theorem <ref>.





The only stochastic ingredient, which we “take from the shelf”,
for Theorem <ref> is part ii) of the following lemma;
part i) is needed in the argument for the lower bound, Theorem <ref> below. Lemma <ref> provides
stochastic bounds on the (augmented) corrector (ϕ,σ), and essentially amounts to
saying that it is bounded with overwhelming probability in d>2 and almost so in d=2
(the form of the statement of Lemma <ref> is marginally weakened by not distinguishing
the cases d=2 and d>2).
By now, there are several approaches to such a result: The (historically) first result of this type is
<cit.>,
and is based on functional inequalities (at first in case of a discrete medium;
see <cit.> for an extension to the continuum case)
and thus (indirectly) relies on an underlying product structure of the ensemble ⟨·⟩
(eg a Gaussian field or Poisson point process like in Definition <ref> below).
Here, we follow a more recent approach that is based on a finite, say unit, range assumption.


There are two possible references for this second approach: <cit.> based
on the variational approach of <cit.>, and <cit.>, based
on a semi-group approach. For reasons of familiarity, we opt for the second reference,
which also has the advantage of treating σ next to ϕ.


Let ⟨·⟩ be a stationary ensemble of uniformly elliptic coefficient fields, cf (<ref>),
that has unit range of dependence. Then for every i=1,⋯,d, there exist a (random) scalar field ϕ_i
and a (random) skew symmetric tensor field σ_i such the gradient fields ∇ϕ_i and ∇σ_i
are stationary, of finite second moments, and of vanishing expectation, and such that (<ref>) (and thus (<ref>))
and (<ref>) hold. Moreover,



  i) For every exponent 0<β<1 there exists a (random) radius r_* such that (<ref>) holds
and which satisfies


    ⟨exp(r_*^s)⟩≤ C


for any exponent 0<s<2 with C=C(d,λ,β,s).


  ii) For every exponents 0<β<1 and 0<s<2(1-β), there exists
a scale L_0=L_0(d,λ,β,s) such that for every L≥ L_0 the statement


    1/R(_B_R|(ϕ,σ)-_B_R(ϕ,σ)|^2)^1/2≤(1/R)^β   R≥ L


fails with probability ≤exp(-(L/L_0)^s).









 §.§ The lower bound


Our second main result states that the fluctuations of ∇ u(0), when conditioned on the
coefficient field inside the ball B_L are at least of the order of (ℓ/L)^d(1/L)^d/2. 
The first factor is a (deterministic) consequence of fact that ∇ u(0), in view of the rhs supported in B_ℓ,
is not very sensitive in the coefficients in B_L^c. 
The second factor scales as the inverse of the square root of the volume 
of the annulus B_2L-B_L, and thus has a CLT-flavor to it. More precisely, instead of ∇ u(0) we monitor
smooth averages of ∇ u on a sufficiently large but order-one scale R near the origin. 
Similarly to (<ref>) we fix a (universal) averaging function on B_R through


    ω(x)=1/R^dω̂(x/R)    ω̂=0  B_1^c,  ∫ω̂=1,  ∇̂ω̂


and consider ∫ω∇ u.
We establish this lower bound on the fluctuations under convenient assumptions on the ensemble: 


Let ⟨·⟩ denote the distribution of the Poisson point process X in ℝ^d of
unit intensity. We assume that there exists a measurable map A from the space of point configurations into
the space of λ-uniformly elliptic coefficient fields; in other words, we
consider the ensemble of such coefficient fields given by


    a(x)=A(x,·)


and thus make the following assumptions on A:
For all points x∈ℝ^d, point configurations X,X'⊂ℝ^d, shift vectors z∈ℝ^d
we impose



  * shift-invariance, that is,


    A(z+x,z+X)=A(x,X),



  * locality, that is,


    A(x,X)=A(x,X')   X∩ B_1(x)=X'∩ B_1(x),



  * monotonicity, that is,


    A(x,∅)<a_h  .




In particular, ⟨·⟩ is stationary and of range unity, so that by Lemma <ref> i), for every β<1
and s<2, there exists r_* with (<ref>) and such that


    ⟨exp(r_*^s)⟩≤ C(d,λ,β,s).





















Let the ensemble be as in Definition <ref>.
Consider the solution u of (<ref>) with rhs of the form (<ref>) for a given ĝ.
Then there exists a radius R such for all scales L, ℓ with C≤ℓ≤1/CL we have


    ⟨|∫ω∇ u-⟨∫ω∇ u|B_L⟩|^2⟩^1/2≥1/C(ℓ/L)^d(1/L)^d/2|∫ĝ|,


where ω is defined as in (<ref>). Here the radius R and the constant C depend on the ensemble,
on the sup norm of ∇̂ω̂,
and on the Hölder norm of ∇̂ĝ.




Like Theorem <ref>, Theorem <ref> relies on a purely deterministic result, namely Proposition <ref> which is
of independent interest. Proposition <ref> monitors the effect of a “defect” in the medium a; 
more precisely, we consider the medium a' given by


    a'=a_0   B_R(y)     a'=a   B_R(y),


where a_0 is some other λ-uniformly elliptic coefficient field, cf (<ref>),
and y is some point and R some radius. We are interested in the effect on the solution
u of our whole-space problem (<ref>) with localized rhs g.
Hence we compare u to u' given by the decaying solution of


    -∇· a'∇ u'=∇· g.


Here, we think of B_R(y) as being far from the origin where g is localized, cf (<ref>).
Proposition <ref> states that to leading order, the effect of the inclusion is captured by its
effect on the level of the homogenized coefficients. 
More precisely, we consider the coefficient field a_h' given by


    a'_h=a_0   B_R(y)     a'_h=a_h   B_R(y)


and the decaying solution u_h' of


    -∇· a'_h∇ u'_h=-∇· a_h∇ u_h,


where u_h is given by (<ref>).
Proposition <ref> states that indeed u'-u≈ u_h'-u_h.
As for the other results, Proposition <ref> does so on the level
of the gradient (and thus involves the two-scale expansion) and in a localized way.
It is the positivity of the exponent α in (<ref>) that ensures that
indeed to leading order, ∇(u'-u) behaves like ∂_i(u_h'-u_h)(e_i+∇ϕ_i), 
since, as a classical argument shows, ∇(u'_h-u_h)(0) generically scales as
(ℓ/|y|)^d(R/|y|)^d. Loosely speaking, Proposition <ref>
states that sensitivity analysis (the dependence of a solution on the coefficient field) and
homogenization commute.


Consider a given λ-uniformly elliptic coefficient field a on ℝ^d,
cf (<ref>). Suppose that there exists a tensor a_h and, for i=1,⋯,d, a scalar field
ϕ_i and a skew-symmetric tensor field σ_i such that (<ref>) holds.
Suppose that for given β∈(0,1) there exists a radius r_*(0) such that (<ref>) holds
and a radius r_*(y) such that (<ref>) holds with the origin replaced by some point y.



We are given another λ-uniformly elliptic coefficient field a_0, cf (<ref>),
and consider a' and a'_h given through (<ref>) and (<ref>).
Consider u, u', u_h, and u_h' defined through (<ref>), (<ref>), (<ref>), and (<ref>).



Under the assumption that ℓ≥ r_*(0) and provided the radius R satisfies 1/2|y|≥ R≥ r_*(0),r_*(y)
we have 


    (_B_R|∇((u'-u)   -(1+ϕ_i∂_i)(u_h'-u_h))|^2)^1/2
       ≤ C(ℓ/|y|)^d(R/|y|)^d((R/|y|)^β+(r_*(y)/R)^α)


with an exponent α>0 that only depends on d, λ and β
and a constant C<∞ that in addition depends on the Hölder norm of ĝ.




Proposition <ref> relies on two ingredients: The first ingredient
is Lemma <ref> which to leading order characterizes u'-u in terms
of the homogenized solution u_h and a tensor δ a_ij,
defined in (<ref>), that captures how the corrector ϕ is
affected by the defect. Since we are interested in u'-u near the
origin, the Green's function G_h of the homogenized operator is also
involved. Since we characterize u'-u on the level of gradients, the
estimate involves the correction ∇ϕ_k.  The second
ingredient is Lemma <ref> that establishes a version of Proposition
<ref> with the general solution u of an equation involving
-∇· a∇ replaced by a specific one, namely the
corrector ϕ_i. With help of this lemma, we establish Corollary
<ref> that characterizes the tensor δ a_ij in terms of
its counterpart δ a_h.



How is the solution of the corrector equation (<ref>) affected by replacing a by a', cf (<ref>)?
We denote by ϕ'_i the solution of


    -∇· a'(e_i+∇ϕ'_i)=0


that behaves as ϕ_i at infinity (in the sense that ∇(ϕ_i'-ϕ_i) is square integrable).
Lemma <ref> compares ϕ_i' to ϕ_ih', the decaying solution of


    -∇· a'_h(e_i+∇ϕ_ih')=0.


More precisely, Lemma <ref> compares the two solutions on the level of the potentials and of flux averages.
Similarly to (<ref>) we fix a (universal) averaging function on B_2R(y) through 




    (x)=1/R^dω̂(x-y/R)  
       ω̂=1/2|B_1|  B_1,  ω̂=0  B_2,  ∫ω̂=1,  ∇̂ω̂ .




In terms of arguments, Lemma <ref> takes inspiration from <cit.>.


We make the same assumptions on a, a_h, ϕ_i, σ_i, r_*(y) as in Proposition <ref>.
Let a_0, a', and a_h' be as in Proposition <ref>. For i=1,⋯,d we consider the solutions
ϕ_i' and ϕ_ih' of (<ref>) and (<ref>). Then for R≥ r_*(y) 
and , cf (<ref>), we have


    1/R(_B_2R(y)(ϕ_i'-ϕ_ih'-_B_2R(y)(ϕ_i'-ϕ_ih'))^2)^1/2≤ C(r_*(y)/R)^α,
    |∫ a'(e_i+∇ϕ_i')-∫ a_h'(e_i+∇ϕ_ih')|
    ≤ C(r_*(y)/R)^α,


where α and C only depend on d, λ and β.




As we shall see in Lemma <ref>, it is the tensor 


    δ a_ij:=∫(e_j+∇ϕ_j')·(a'-a)(e_i+∇ϕ_i)


that to leading order governs the sensitivity of the solution operator of -∇· a∇ under changing
the coefficients from a to a'.  It is an easy consequence of Lemma <ref> that to leading order, we may
replace the medium described by a (in both a itself and in its perturbation a', and the corresponding
correctors) by its homogenized version a_h in the expression (<ref>).


Under the assumptions of Lemma <ref> we have for R≥ r_*(y)


    |1/R^d∫(e_j+∇ϕ_j')·(a'-a)(e_i+∇ϕ_i)
       -1/R^d∫(e_j+∇ϕ_jh')·(a'_h-a_h)e_i|≤ C(r_*(y)/R)^α,


where α and C are as Lemma <ref>.




The following lemma heavily relies on <cit.>, indirectly through
expanding on Lemma <ref>, but also more directly for localization.


We make the same assumptions on a, a_h, ϕ_i, σ_i, r_*, r_*(y) as in Proposition <ref>.
Let a_0, a', a_h', u, and u' be as in Proposition <ref>. Then for R≥ r_*, r_*(y) we have


    (_B_R|∇(u'-u)   
    -∂_iu_h(y)δ a_ij∂_j∂_kG_h(-y)(e_k+∇ϕ_k)|^2)^1/2
       ≤ C(ℓ/|y|)^d(R/|y|)^d+β,


where C depends on d, λ, β and the Hölder norm of ĝ.




Equipped with Lemma <ref> and Corollary <ref>, we now may show that ∇ u(0), or rather
a smooth average ∫ω∇ u of ∇ u near the origin, substantially reacts to a change
in the medium at B_R(y). This reaction is characterized in terms of the tensor 


    δ a_hij:=∫(e_j+∇ϕ_jh')·(a_h'-a_h)e_i


appearing in (<ref>).
In fact, Theorem <ref> is not inferred from Proposition <ref>,
but rather directly from Corollary <ref> and Lemma <ref>,
which we combine for that purpose to


We make the same assumptions on a, a_h, ϕ_i, σ_i, r_*, r_*(y) as in Proposition <ref>.
Let a_0, a', a_h', u, and u' be as in Proposition <ref>. Then provided R≤ℓ≤|y| we have


    |∫ω∇(u'-u)
    -∫ g_m∂_m∂_iG_h(y)δ a_hij∂_j∂_kG_h(-y)e_k|≤ C
    (ℓ/|y|)^d(R/|y|)^d
       ×((R/|y|)^β+(r_*/R)^β+(r_*(y)/R)^α
    +ℓ/|y|+(r_*/ℓ)^β)(1+r_*+r_*(y)/R)^3/2d+β.


where ω is defined in (<ref>).
Here C depends on d, λ, β, and the L^2-norm of ∇̂ĝ.


The first rhs summand in (<ref>) comes from Lemma <ref>, the second one
from eliminating ∇ϕ_k in (<ref>), the third one comes from
Corollary <ref>, the fourth comes from the dipole expansion of u_h,
and the fifth from simplifying the dipole moment. The last factor on the rhs of (<ref>) arises because
we avoid any smallness condition on max{r_*,r_*(y)} in terms of R or ℓ.



Under our assumptions on the ensemble, cf Definition <ref>, and
equipped with the deterministic result of Corollary <ref>, we obtain
a lower bound on the variance of the expectation of
∫ω∇ u conditioned on the Poisson process restricted to
B_R(y).  Note that we first average over the Poisson process in the
complement B_R^c(y), and then consider the variance wrt the Poisson process on B_R(y)
 — as opposed to the opposite order, 

which would amount to a weaker result, in
particular too weak for the purpose of Theorem <ref>.  This lower
bounds holds provided the order-one radius R, which also governs the
average through ω, cf (<ref>), is sufficiently large.  The lower bound is
optimal in terms of the scaling in the ratios of the length scales 1
(the range of dependence of the random coefficient field a), ℓ (the
scale of the source g), and |y| (the distance between the source
and the site at which the variance is probed).


Under the assumptions of Theorem <ref> there exists an R 
so that provided C≤ℓ≤1/C|y| we have 


    ⟨|⟨∫ω∇ u|B_R(y)⟩-⟨∫ω∇ u⟩|^2⟩^1/2   ≥1/C(ℓ/|y|)^d(1/|y|)^d|∫ĝ|,


where R and C are as in Theorem <ref>.



The following auxiliary lemma will be used in the proof of Lemma <ref>
and states that the random (scalar) field r_*=r_*(a,x) is not very sensitive to local changes in the coefficient field a.
Here we think of r_* as being defined as the minimal radius with the property (<ref>);
we denote by r_*' the radius for the medium a', and by r_*(y), r_*'(y)
the corresponding radii with the origin replaced by the point y (in the neighborhood of which
a' differs from a). Here, and in analogy to ϕ_i', we think of σ_ijk' as the solution
of -Δσ_ijk'=∂_jq'_ik-∂_kq'_ij, where q_i':=a'(e_i+∇ϕ_i'),
that behaves at ∞ as σ_ijk (ie ∫|∇(σ_i'-σ_i)|^2<∞). We note that
by construction we have Δ(q_i'-∇·σ_i')=0, so that by the decay of q_i'-q_i
and ∇σ_i'-∇σ_i, the relation (<ref>) is preserved, that is,
q_i'=a_he_i+∇·σ'_i.


We have 


    r_*'(y)+r_*'≤ C(d,β)(r_*(y)+r_*+R).









§ NUMERICAL RESULTS


In our numerical tests, we consider a random ensemble according to
Definition <ref>. Let X be the Poisson point process on
ℝ^2 with unit intensity, the coefficient field a is given by 

    a(x) = max_ξ_i ∈ X {φ (x - ξ_i) } + 1/2,

where φ is a bump function supported in B_1/2

    φ(x̂) = exp(- 1/1 - 4 x̂^2).


The elliptic PDE is 

    - ∇· a ∇ u = f

where the rhs is given by 

    f(x) = x_2 exp(- 5/5 - x^2)

and thus f is compactly supported inside the ball B_√(5)
with average 0. It is clear that we can rewrite f = ∇· g with g supported in B_√(5).

To numerically approximate the solution of the Dirichlet problems
(<ref>), (<ref>), and (<ref>), we use a standard
second-order centered finite-difference scheme with mesh size
δ x = 0.1. 

In Figure <ref>, we plot the numerically obtained _h as
a function of L for 4 independent realizations of the coefficient
field a. It shows that the Dirichlet approximation _h converges as L increases, validating (<ref>) in Lemma <ref>.


We now consider our algorithm to approximate ∇ u(0) based on
the a |_Q_2L. To validate the algorithm, we check the
numerical convergence rate and compare it with two approximations that
with a slower convergence rate. Recall that our algorithms consists of
the steps:


  * Solve the equation (<ref>) for approximate correctors _i;

  * calculate the homogenized coefficient _h as in (<ref>);

  * obtain _h from the homogenized equation (<ref>) and the dipole correction (<ref>); 

  * solve the equation (<ref>) for . 

In comparison, we consider two other algorithms: 1) Solving the
equation with homogeneous Dirichlet boundary condition, ie
(<ref>), referred as “Dirichlet algorithm”; 2) Dropping the
dipole correction, ie, instead of (<ref>), one solves (<ref>)
using approximate homogenized coefficient and correctors

    - ∇· a ∇_II  = ∇· g  in Q_L,   _II = ( 1 + _i ∂_i) _h  on ∂ Q_L.

Thus in the boundary condition, _h is replaced by _h so the
dipole correction is dropped. This will be referred as the “no dipole
algorithm”. In Figure <ref>, we compare the numerical
convergence for the three algorithms with two different realization of
the random media. The difference
∇ u^(2L)(0) - ∇ u^(L)(0) is plotted for various L
for the three algorithms. As the plot is in loglog scale, a straight
line shows algebraic convergence with the slope indicating the
convergence rate. We observe from the numerical results that the
proposed algorithm achieves almost O(L^-3) convergence, which is
consistent with the theoretical rate
O((ℓL)^2 (1L)^β) with β close to
1. The wiggled line is due to the randomness of the coefficients
added in the annulus Q_4L - Q_2L when L is doubled.  We also
observe that the other two algorithms have the slower rate of
O(L^-2), the main reason being that they do not capture the
correct dipole in the far field. It can be also seen that by using the
information from homogenization, the “no dipole algorithm” indeed
reduces the error compared with the simple minded Dirichlet
approximation.



Finally, we study the error committed by the proposed algorithm with
the uncertainty of ∇ u(0) due to the unknown coefficient
outside Q_2L. In the second row of Table <ref>, we
compare the approximation ∇(0) with the reference solution
obtained by a calculation with a large domain (L = 768 is
chosen). In the remaining rows, we show how ∇ u(0) varies 
when we re-sample the random medium outside Q_2L (where again, we take ∇ u^(768)(0) as a proxy for ∇ u(0)). These thus show the sensitivity of ∇ u(0) with 
respect to the change of the media. We observe from the numerical
results that while the error made by our algorithm is larger than the
typical fluctuations  of ∇ u(0), the decay rate matches and the error is
also comparable to the sensitivity for various L.












§ PROOFS



We start with the remark that (<ref>) actually implies (<ref>).







Indeed, what separates (<ref>) from (<ref>) is 


    |_B_r(ϕ,σ)-_B_r'(ϕ,σ)|
    ≲ r(r_*/r)^β   r≥ r'≥ r_*.


Since β<1, we may apply dyadic decomposition to reduce this to


    |_B_r(ϕ,σ)-_B_r'(ϕ,σ)|
    ≲ r(r_*/r)^β   2r'≥ r≥ r'≥ r_*,


which by the triangle inequality is a consequence of (<ref>). In view of (<ref>),
we may change (ϕ,σ) by an additive constant, and thus without affecting (<ref>)
to upgrade (<ref>) to


    1/r(_B_r|(ϕ,σ)|^2)^1/2≲(r_*/r)^β   r≥ r_*.


This is in line with (<ref>).



For later reference we also note that


    (_B_R|e_i+∇ϕ_i|^2)^1/2≲ 1   R≥ r_*,


which is a consequence of Caccioppoli's estimate on the a-harmonic function x_i+ϕ_i, cf (<ref>),
and (<ref>) in form of (_B_R(x_i+ϕ_i)^2)^1/2
≲ R(r_*/R)^β.



We now follow the steps of the proof of <cit.> and like there may wlog assume that r_*=1.
In analogy to Step 2 in the proof of <cit.> we claim


    (∫|∇ u|^2)^1/2≲ℓ^d/2,


and


    |x|^d|∇ u_h(x)|+|x|^d+1|∇^2u_h(x)|≲ℓ^d   x≠0.


We furthermore claim that the u and u_h have not only vanishing “constant invariant”


    ∫∇η_r· a∇ u=∫∇η_r· a_h∇ u_h=0   r≥ℓ,


but also identical “linear invariants”, that is,


    ∫∇η_r·((x_i+ϕ_i)a∇ u-ua(e_i+∇ϕ_i))
       =∫∇η_r·(x_ia_h∇ u_h-u_ha_he_i)   r≥ℓ,


where η_r(x)=η(x/r), and η is a cut-off function for B_1 in B_2.
Indeed, (<ref>) follows from the energy estimate for (<ref>) and the form (<ref>)
of the rhs. By the triangle inequality, we split (<ref>) into the corresponding estimate
for ũ_h and for δ u_h:=ξ_i∂_iG_h with ξ_i:=∫∇ϕ_i· g.
We first turn to ũ_h defined through (<ref>) which in view of (<ref>) is of
the form ∇ũ_h(x)=∇̂û_h(x/ℓ) with -∇̂· a_h∇̂û_h
=∇̂·ĝ, so that the desired estimate follows from


    |∇̂û_h|≲1/|x̂|^d+1    
    |∇̂^2 û_h|≲1/|x̂|^d+1+1,


which in turn are a consequence of standard Schauder theory based on the Hölder continuity of ∇̂·ĝ
for the near-field, and the Green's function representation for the far-field. 
We now turn to δ u_h and note that the estimates


    |∇δ u_h|≲ℓ^d/|x|^d    
    |∇^2 δ u_h|≲ℓ^d/|x|^d+1


follow from the homogeneity of ∇^2 G_h of degree -d and the fact that |ξ_i|≲∫_B_ℓ|∇ϕ_i|, 
Jensen's inequality together with (<ref>).
We now turn to (<ref>). The identity for u follows directly from (<ref>) and the fact that g is
supported in B_ℓ, cf (<ref>). For the same reason, the contribution to the constant invariant
of u_h coming from ũ_h vanishes; the one coming from δ u_h, which is constant in r>0 by a_h-harmonicity
(thus the name “invariant”),
must vanish by its homogeneity of order 1-d.
We finally turn to (<ref>). By (<ref>) and (<ref>) we learn from an integration by parts and the support condition
on g that for r≥ℓ,


    ∫∇η_r·((x_i+ϕ_i)a∇ u-ua(e_i+∇ϕ_i))
    =∫(e_i+∇ϕ_i)· g.


Likewise, we obtain from (<ref>) that


    ∫∇η_r·(x_ia_h∇ũ_h-ũ_ha_he_i)=∫ e_i· g.


Finally, by definition of the Green's function G_h we have


    ∫∇η_r·(x_ia_h∇∂_j G_h-∂_j G_ha_he_i)=δ_ij.


By u_h=ũ_h+ξ_j∂_jG_h and the definition of ξ, the last two identities combine to


    ∫∇η_r·((x_ia_h∇ u_h-u_ha_he_i)=∫(e_i+∇ϕ_i)· g.


This and (<ref>) yield the identity (<ref>).



Following Step 3 in the proof of <cit.> (with η=0, the additional multiplicative factor ℓ^d,
and β playing the role of 1-α)
we consider the error in the two-scale expansion w:=u-(1+ϕ_i∂_i)u_h we claim that


    (∫_B_1^c|∇ w|^2)^1/2≲ℓ^d


and that


    -∇· a∇ w=∇· h     h:=(ϕ_ia-σ_i)∇∂_i u_h,


where 


    (1/R^d∫_B_r^c|h|^2)^1/2≲ℓ^d(1/r)^d+β   r≥ 1.


We first turn to (<ref>) and note
∇ w=∇ u-∂_iu_h(e_i+∇ϕ_i)-ϕ_i∇∂_iu_h. Hence for the ∇ u-contribution
the desired estimate follows from (<ref>) (and ℓ≥ 1). According to (<ref>), the contribution
from ∂_iu_h(e_i+∇ϕ_i) is
estimated by ℓ^d(∫_B_1^c|x|^-2d|e_i+∇ϕ_i|^2)^1/2; dividing the integral into
dyadic annuli and using (<ref>), we see that it is ≲ 1. Still according to (<ref>),
the contribution from ϕ_i∇∂_iu_h is estimated by 
ℓ^d(∫_B_1^c|x|^-2(d+1)ϕ_i^2)^1/2; dividing also this integral into
dyadic annuli and using (<ref>), we see that it is ≲ 1.

As -∇· a∇ u=∇· g=∇· a_h∇ u_h (in view of (<ref>), (<ref>),
and the a_h-harmonicity of G_h and thus δ u_h away from the origin), 
the identity (<ref>) follows by a straightforward calculation
as in Step 3 of <cit.>. It is the merit of σ_i that the
rhs can be written in divergence-form ∇· h. 

We finally turn to the estimate (<ref>); by (<ref>) and (<ref>) we have
|h|≲ℓ^d|x|^-(d+1)|(ϕ,σ)|. Hence the estimate follows once more from (<ref>) 
and division into dyadic annuli.



Following Step 4 and Step 5 of <cit.>
(still with α replaced by 1-β and the additional factor of ℓ^d)
we obtain


    (1/R^d∫_B_r^c|∇ w|^2)^1/2≲ℓ^d(1/r)^d+β   r≥ 2,


which (in view of r_*=1) turns into (<ref>). Note that the outcome of Step 5 of <cit.>
is worse by a log R; which however can be easily avoided, cf <cit.>
with k=1 and r=r_*=1.



Let w:=u-(1+ϕ_i∂_i)ũ_h denote the error in the two-scale convergence.
A crucial but elementary calculation based on the equation (<ref>), (<ref>),
and (<ref>), see for instance Step 2 of the proof of <cit.>, shows that
that the related residuum can be written in divergence form:


    -∇· a∇ w=∇· h     h:=(ϕ_ia-σ_i)∇∂_iũ_h.


This is the main merit of σ.
We note that the rhs h has good decay:


    (_B_2r-B_r|h|^2)^1/2≲(r_*/r)^β(ℓ/r)^d
       r≥ r_*.


Indeed, by (<ref>) and (<ref>), this follows from 


    sup_B_r^c|∂_iũ_h|+rsup_B_r^c|∇∂_iũ_h|
    ≲(ℓ/r)^d,


which itself is a consequence of the form (<ref>) and standard regularity theory for the
constant-coefficient equation (<ref>).



We now argue that we may construct a solution w̃ of the equation (<ref>) that inherits the good decay of h, that is,


    (_B_R^c|∇w̃|^2)^1/2≲(r_*/R)^β(ℓ/R)^d
       R≥ r_*,


which however satisfies the equation only in an exterior domain, that is,


    -∇· a∇w̃=∇· h   B_r_*^c.


Indeed, to this purpose we resort to a dyadic decomposition: For every dyadic multiple r of r_* let
w_r denote the Lax-Milgram solution of


    -∇· a∇ w_r=∇·( I(B_2r-B_r)h),


where I(B_2r-B_r) denotes the characteristic function of the annulus B_2r-B_r.
By the energy estimate and the bound (<ref>) on h we have


    (1/r^d∫|∇ w_r|^2)^1/2≲(r_*/r)^β(ℓ/r)^d.


Fix an exponent α∈(β,1).
By <cit.> there exists an a-harmonic function δ w_r in B_r_*^c such that on the one hand,
it captures the large-scale behavior of w_r, which is a-harmonic in B_2r^c,  in the sense of


    (_B_R^c|∇(w_r-δ w_r)|^2)^1/2≲(r/R)^d+α(1/r^d∫|∇ w_r|^2)^1/2   R≥ r.


On the other hand, it behaves like a dipole in the sense that we have the inverse estimate 


    (_B_R^c|∇δ w_r|^2)^1/2≲(r/R)^d(_B_r^c|∇δ w_r|^2)^1/2   r_*≤ R≤ r.


From this and (<ref>) for R=r we obtain


    (_B_R^c|∇δ w_r|^2)^1/2≲(r/R)^d(1/r^d∫|∇ w_r|^2)^1/2   r_*≤ R≤ r,


so that by the triangle inequality in L^2


    (_B_R^c|∇(w_r-δ w_r)|^2)^1/2≲(r/R)^d(1/r^d∫|∇ w_r|^2)^1/2   r_*≤ R≤ r.


Inserting (<ref>) into (<ref>) and (<ref>) we obtain


    (_B_R^c|∇(w_r-δ w_r)|^2)^1/2   ≲
    (ℓ/R)^d(r_*/R)^β(r/R)^α-β   R≥ r,
    (_B_R^c|∇(w_r-δ w_r)|^2)^1/2   ≲
    (ℓ/R)^d(r_*/R)^β(R/r)^β   r_*≤ R≤ r.


Thanks to both β<α and β>0, which leads to the convergence of the geometric series in r,
we now see that w̃:=∑_r(w_r-δ w_r) satisfies not only (<ref>), cf (<ref>)
and the a-harmonicity of δ w_r in B_r_*^c, but also (<ref>).



We now argue that there exists an a_h-dipole δ u_h, that is, a function of the form 


    δ u_h=ξ_i∂_iG_h   ξ,


such that 


    (_B_R^c|∇(w-w̃-(1+ϕ_i∂_i)δ u_h)|^2)^1/2≲ (ℓ/R)^d(r_*/R)^β   R≥ r_*.


Indeed, by (<ref>) and (<ref>), w-w̃ is a-harmonic in B_r_*^c. Hence by <cit.>, it remains
to show that it has (at least) dipole-like decay, that is,


    (_B_R^c|∇(w-w̃)|^2)^1/2≲ (ℓ/R)^d   R≥ r_*.


In view of (<ref>), by the triangle inequality in L^2, it is enough to consider w. Again
by the triangle inequality on ∇ w =∇ u -∂_iũ_h(e_i+∇ϕ_i)
-ϕ_i∇∂_iũ_h, and in view of (<ref>), this reduces to


    (_B_R^c|∇ u|^2)^1/2≲ (ℓ/R)^d/2≤ (ℓ/R)^d,
    (_B_R^c(1/|x|^d|e_i+∇ϕ_i|)^2)^1/2≲1/R^d,
    
    (_B_R^c(1/|x|^d+1ϕ_i)^2)^1/2≲1/R^d
    (r_*/R)^β≤1/R^d.


Estimate (<ref>) follows from the energy estimate for (<ref>) and the form
(<ref>) of the rhs:


    (_B_R^c|∇ u|^2)^1/2≤(1/R^d∫|∇ u|^2)^1/2≲(1/R^d∫|g|^2)^1/2≲ (ℓ/R)^d/2.


Estimate (<ref>) follows by decomposition into dyadic annuli from (<ref>). Likewise, estimate
(<ref>) follows by dyadic decomposition from


    (_B_R|e_i+∇ϕ_i|^2)^1/2≲ 1   R≥ r_*,


which is a consequence of Caccioppoli's estimate on the a-harmonic function x_i+ϕ_i, cf (<ref>),
and (<ref>) in form of (_B_R(x_i+ϕ_i)^2)^1/2 
≲ R(r_*/R)^β.



We now are in the position to conclude: The combination of (<ref>) and (<ref>) yields (<ref>).
The claimed representation (<ref>) reduces to the statement that ξ_i=∫∇ϕ_i· g, cf (<ref>),
to which we turn now.
We note that by the support condition on g, in the exterior domain B_ℓ^c,
u_h is a_h-harmonic and u is a-harmonic. Since for any coordinate direction j,
the functions x_j and x_j+ϕ_j are a_h-harmonic and a-harmonic, respectively, the vector fields 
x_ja_h∇ u_h-u_ha_he_j and (x_j+ϕ_j)a∇ u-ua(e_j+∇ϕ_j) are divergence-free
in B_ℓ^c. Fixing a smooth function η̂(x̂) that is 1 on {|x̂|≤ 1} and vanishes
on {|x̂|≥ 2} and setting η_r(x):=η̂(x/r) we thus have that


    ∫∇η_r·(x_ja_h∇ u_h-u_ha_he_j),
    ∫∇η_r·((x_j+ϕ_j)a∇ u-ua(e_j+∇ϕ_j))


do not depend on r≥ℓ. According to <cit.>, these invariants are well-behaved under
the two-scale expansion:


    ∫∇η_r·(x_ja_h∇ u_h-u_ha_he_j)
       =lim_R↑∞∫∇η_R·((x_j+ϕ_j)a∇ (1+ϕ_i∂_i)u_h-
    ((1+ϕ_i∂_i)u_h)a(e_j+∇ϕ_j)).


We note that (<ref>) implies that


    ∫∇η_R   ·((x_j+ϕ_j)a∇ (1+ϕ_i∂_i)u_h-
    ((1+ϕ_i∂_i)u_h)a(e_j+∇ϕ_j))
    
    -∫∇η_R   ·((x_j+ϕ_j)a∇ u-ua(e_j+∇ϕ_j))


converges to zero as R↑∞: For the first summand, this follows directly from
(<ref>) and the (L^2-averaged) linear growth of x_j+ϕ_j on scales R≥ r_*, cf
(<ref>). For the second summand, we note that the integrals on both sides do not change under modifying u_h and u by an additive constant, 
so that we may post-process (<ref>) by Poincaré's inequality on the dyadic annulus {R≤|x|≤ 2R} to obtain smallness of
u-(1+ϕ_i∂_i)u_h-c, which we combine with (<ref>). 
Hence the two invariants coincide:


    ∫∇η_r·(x_ja_h∇ u_h-u_ha_he_j)=∫∇η_r·((x_j+ϕ_j)a∇ u-ua(e_j+∇ϕ_j)),


which we rewrite as


    ∫∇η_r·(x_ja_h∇δ u_h-δ u_ha_he_j)
       =∫∇η_r·((x_j+ϕ_j)a∇ u-ua(e_j+∇ϕ_j))
    -∫∇η_r·(x_ja_h∇ũ_h-ũ_ha_he_j).


Using the equations (<ref>),  (<ref>), 
and (<ref>) in the distributional form of -∇· a_h∇δ u_h
=ξ_i∂_iG_h, the latter identity turns into the desired


    ξ_j=∫(e_j+∇ϕ_j)· g-∫ e_j· g.










 

Let us introduce a tool we will often use, namely the C^0,1-estimate for an a-harmonic function w
in B_R, the crucial role of which was recognized in <cit.>. As is obvious from its characterization (<ref>), r_* 
dominates, up to a multiplicative constant C=C(d,λ,β), the “minimal radius” introduced in <cit.>
and from which on the “mean-value property”, as an estimate, holds for |∇ w|^2, see <cit.> for
the proof. We record this


    _B_r|∇ w|^2≤ C(d,λ,β)_B_R|∇ w|^2   r_*≤ r≤ R.




We note that the combination of (<ref>), (<ref>) and (<ref>) may be rephrased in terms of w:=u_III-u 
and w_D:=u-(1+ϕ_i∂_i)u_h as


    -∇· a∇ w=0  Q_L     w=w_D  ∂ Q_L.


Once we show that this implies


    (∫_Q_L|∇ w|^2)^1/2≲(∫_Q_2L-Q_L|∇ w_D|^2)^1/2,


we see that (<ref>) follows from (<ref>) for R=L, and a subsequent application of the C^0,1-estimate
for w to get from Q_L to B_R for R≥ r_*. 



We now turn to the argument that (<ref>) implies (<ref>) under
the mere assumption of uniform ellipticity, cf (<ref>). Hence by rescaling, we may wlog
assume L=1. Since (<ref>) and (<ref>) are oblivious to additive
constant we may wlog assume ∫_Q_2-Q_1w_D=0, which allows to construct an extension
w̅_D of w_D on Q_1 such that


    (∫_Q_2|∇w̅_D|^2)^1/2≲(∫_Q_2-Q_1|∇ w_D|^2)^1/2.


This extension allows us to reformulate (<ref>) as


    -∇· a∇(w-w̅_D)=∇· a∇w̅_D  Q_1,
       w-w̅_D=0  ∂ Q_1.


Hence (<ref>) follows by the energy estimate on the latter, followed by a triangle
inequality in L^2 yielding (∫_Q_1|∇ w|^2)^1/2
≲(∫_Q_1|∇w̅_D|^2)^1/2, into which we insert
(<ref>).





We first turn to (<ref>) and, for notational simplicity, drop the index i=1,⋯,d.
By (<ref>), which is as mentioned a consequence of (<ref>), and (<ref>) 
we have -∇· a∇(-ϕ) =0 in Q_2L,
so that by Caccioppoli's estimate and the triangle inequality in L^2 we have


    (_Q_3/2L   |∇(-ϕ)|^2)^1/2≲1/L(_Q_2L((-ϕ)-_Q_2Lϕ)^2)^1/2
       ≤1/L(_Q_2L()^2)^1/2
    +1/L(_Q_2L(ϕ-_Q_2Lϕ)^2)^1/2.


Hence (<ref>) follows from the ϕ-parts in (<ref>) and (<ref>) .



For (<ref>), we argue in a similar way: According to (<ref>), which yields (<ref>), and (<ref>)
we have -Δ(_jk-σ_jk) =∂_j(_k-q_k)-∂_k(_j-q_j) in Q_2L, 
where we recall q=a(e+∇ϕ), =a(e+∇) and thus -q =a∇(-ϕ).
By Caccioppoli's estimate and the upper bound on a provided by (<ref>) we therefore obtain


    (_Q_L|∇(-σ)|^2)^1/2
       ≲(_Q_3/2L|∇(-ϕ)|^2)^1/2
    +1/L(_Q_2L||^2)^1/2
    +1/L(_Q_2L|σ-_Q_2Lσ|^2)^1/2.


Hence (<ref>) follows from (<ref>) and the σ-parts of both (<ref>) and (<ref>).



We now turn to (<ref>) and recall our averaging function ω, cf (<ref>).


We first claim that also for the whole-space case,
the effective coefficient can approximately be recovered from averaging the flux wrt to ω:


    |a_he-∫ω a(e+∇ϕ)|≲(r_*/L)^β.



The argument for this is based on (<ref>), which yields the identity
a_he-∫ω a(e+∇ϕ) =∫(σ-_B_Lσ) ∇ω, and thus by (<ref>) the estimate 
|a_he -∫ω a(e+∇ϕ)| ≲1/L_B_L|σ-_B_Lσ|. Hence we obtain
(<ref>) from (<ref>) and Jensen's inequality. 



We now may conclude: By the upper bound provided in (<ref>), we obtain from
(<ref>) and the definition (<ref>) by the triangle inequality


    |(_h-a_h)e|≲_Q_L|∇(-ϕ)|+(r_*/L)^β,


so that (<ref>) follows from (<ref>) and Jensen's inequality.





We first compare the two solutions _h and ũ_h of (<ref>) and (<ref>), respectively,
and claim that


    |∇(_h-ũ_h)|+L|∇^2(_h-ũ_h)|
    ≲ (r_*/L)^β(ℓ/L)^d   Q_L^c.


Indeed, the difference w:=_h-ũ_h satisfies


    -∇·_h∇ w=∇· (_h-a_h)∇ũ_h.


The form (<ref>) of the rhs g in (<ref>) transmits to the solution ũ_h:
∇ũ_h(x) =∇̂ũ_h(x/ℓ), and then also to w:
∇ w(x) =∇̂ŵ(x/ℓ). Recall that ĝ is compactly supported
in B_1 with Hölder continuous derivatives; in particular, ĝ is in the class of vector fields 
decaying as O(|x|^-d) with derivatives decaying as O(|x|^-d-1) and differences of the derivatives
decaying as O(|x|^-d-1-α|x-y|^α), in line with α-Hölder continuity for some α∈(0,1).
This class is preserved under the constant-coefficient Helmholtz projection (recall that _h satisfies
(<ref>)). Hence we obtain in particular


    |∇̂ŵ|≲|_h-a_h|1/(|x̂|+1)^d,     
    |∇̂^2ŵ|≲|_h-a_h|1/(|x̂|+1)^d+1.


Translating back to the microscopic variables and using (<ref>) in Lemma <ref> we have


    |∇ w|≲ (r_*/L)^β(ℓ/|x|+ℓ)^d,     ℓ|∇^2w|≲(r_*/L)^β(ℓ/|x|+ℓ)^d+1,


from which we extract (<ref>).



We now compare _h and u_h defined in (<ref>) and (<ref>) and claim that


    |∇(_h-u_h)|+L|∇^2(_h-u_h)|
    ≲ (r_*/L)^β(ℓ/L)^d   Q_L^c.


For later purpose, we also record


    |∇_h|+L|∇^2 _h|
    ≲(ℓ/L)^d   Q_L^c.


In order to pass from (<ref>) to (<ref>), and from (<ref>) and (<ref>) to (<ref>), it remains to control
the dipole contributions:


    |∇(_i∂_i_h-ξ_i∂_iG_h)|   +L|∇^2(_i∂_i_h-ξ_i∂_iG_h)|
    ≲ (r_*/L)^β(ℓ/L)^d
      
    |∇ξ_i∂_i_h|   +L|∇^2ξ_i∂_i_h|
    ≲(ℓ/L)^d   Q_L^c,


where we have set for abbreviation _i :=∫∇_i· g and (as above) ξ_i =∫∇ϕ_i· g.
Because of the obvious estimates on the constant-coefficient (and thus homogeneous) fundamental solution


    |∇∂_i_h|+L|∇^2∂_i_h|≲1/L^d   Q_L^c


and


    |∇∂_i(_h-G_h)|+L|∇^2∂_i(_h-G_h)|   ≲|_h-a_h|1/L^d
       ≲ (r_*/L)^β1/L^d   Q_L^c,


it suffices to show


    |ξ|≲ℓ^d    |-ξ|≲ (r_*/L)^βℓ^d.


By definition of ,ξ, the form (<ref>) of g, and Jensen's inequality the latter two estimates follow from


    (_B_ℓ|∇ϕ|^2)^1/2≲ 1
        (_B_ℓ|∇(-ϕ)|^2)^1/2≲ (r_*/L)^β.


The first estimate in (<ref>) follows from (<ref>) since ℓ≥ r_*. The second estimate 
in (<ref>) follows from (<ref>) in Lemma <ref> and the C^0,1-estimate for an a-harmonic function
w in Q_L, cf (<ref>),




which requires L≥ℓ≥ r_*.



We finally compare  defined through (<ref>) with u_III defined through (<ref>) and claim that


    (_B_R|∇(-u_III)|^2)^1/2≲ (ℓ/L)^d(r_*/L)^β   L≥ R≥ r_*.


This allows us to pass from (<ref>) in Corollary <ref> to this proposition's statement (<ref>).
We note that because w:=-u_III satisfies


    -∇· a∇ w=0   Q_L,   w=w_D   ∂ Q_L,


where w_D:=(1+_i∂_i)_h-(1+ϕ_i∂_i)u_h. Hence we have by (a slight adaptation of)
the argument at the beginning of the proof of Corollary <ref> that


    (_Q_L|∇ w|^2)^1/2≲(_Q_3/2L-Q_L|∇ w_D|^2)^1/2.


We combine this with the C^0,1-estimate, cf (<ref>), for a-harmonic functions,






so that it remains to show


    (_Q_3/2L-Q_L|∇ w_D|^2)^1/2≲ (ℓ/L)^d(r_*/L)^β.


We appeal to the triangle inequality in L^2 to split ∇ w_D into the four contributions


    ∇ w_D   =∂_i(_h-u_h)(e_i+∇ϕ_i)+ϕ_i∇∂_i(_h-u_h)
       +∂_i_h∇(_i-ϕ_i)+(_i-ϕ_i)∇∂_i_h.


For the first contribution, (<ref>) follows from the first part of (<ref>) and (<ref>). 
For the second contribution, (<ref>) follows from the second part of (<ref>) and (<ref>), recall
the normalization of ϕ_i through _B_r_*ϕ_i=0, cf (<ref>).
For the third contribution, we appeal to the first part of (<ref>) and (<ref>).



For the fourth contribution we appeal once more to (<ref>), and now also to (<ref>),








which we combine to (_Q_3/2L(_i-ϕ_i)^2)^1/2 ≲ L(r_*/L)^β. 
In connection with the second part of (<ref>), we see that also that this last contribution is controlled as stated in (<ref>).







As mentioned in the introduction, (<ref>) and (<ref>) do not ensure the analogue of (<ref>),
that is, _i :=a(e_i+∇_i) =_he_i+∇·_i.
However, we will need the C^1,α-theory for a from <cit.> on scales L/2≥ r≥_*, 
which relies on this identity next to the controlled sublinearity (<ref>). Hence
for arbitrary but fixed i=1,⋯,d,
we need to modify _i to (a still skew symmetric) σ̃_i with


    _i=_h e_i+∇·σ̃_i   Q_L,


while retaining the σ-part of (<ref>) in (the seemingly weaker) form of


    1/r(_Q_r|σ̃_i-_Q_rσ̃_i|^2)^1/2≲ (_*/r)^β   L≥ r≥_*.


To this purpose, we first argue that the “defect” _i:=_i-_h e_i-∇·_i, 
which as mentioned in the introduction is component-wise harmonic on Q_2L by (<ref>) and (<ref>), satisfies


    (_Q_L|_i|^2)^1/2≲(_*/L)^β.


Indeed, this follows immediately from the identity (<ref>) in conjunction with the three estimates
(<ref>), (<ref>), and (<ref>) in Lemma <ref> (applied with r_* replaced by _*).



We now turn to (<ref>) and (<ref>). In view of the definition of the defect _i, Poincaré's inequality, 
and (<ref>), it is enough to construct δσ̃_i such that


    _i=∇·δσ̃_i   Q_L


and


    (_Q_r|∇δσ̃_i|^2)^1/2≲ (_*/L)^β   L≥ r≥_*.


To this purpose, we extend the restriction of _i to Q_L periodically (without changing the notation)
and let δ_i be the Q_L-periodic solution of


    -Δδ_ijk=∂_j_ik-∂_k_ij    _Q_Lδ_ijk=0.


By a remark in the introduction, periodic boundary conditions ensure


    _i-_Q_L_i=∇·δ_i.


We claim that (<ref>) holds for δ_i:


    (_Q_r|∇δ_i|^2)^1/2≲ (_*/L)^β   L≥ r≥_*.


Indeed, since by the energy estimate for the periodic problem (<ref>), that is,


    (_Q_L|∇δ_i|^2)^1/2≲(_Q_L|_i|^2)^1/2,


for (<ref>) it suffices to establish the stronger estimate


    sup_Q_L/4|∇δ_i|≲(_Q_L|_i|^2)^1/2
    +(_Q_L|∇δ_i|^2)^1/2


and to appeal to (<ref>). Estimate (<ref>) is a consequence of Sobolev's estimate


    sup_Q_L/4|∇δ_i|≲(_Q_L/4|∇^n+1δ_i|^2)^1/2
    +(_Q_L/4|∇δ_i|^2)^1/2,


where n is an integer larger than d/2, of a localized energy estimate for 
the constant-coefficient equation (<ref>) in form of


    (_Q_L/4|∇^n+1δ_i|^2)^1/2
    +(_Q_L/4|∇δ_i|^2)^1/2≲(_Q_L/2|∇^n_i|^2)^1/2
    +(_Q_L/2|∇δ_i|^2)^1/2,


and of a localized energy estimate for -Δ_i=0


    (_Q_L/2|∇^n_i|^2)^1/2
    ≲(_Q_L|_i|^2)^1/2.


In order to pass from (<ref>) to (<ref>), 
it is enough to (explicitly) construct an affine δσ̃_ijk
with _Q_L_i =∇·δσ̃_ijk and such that |∇δσ̃_ijk|
≲|_Q_L_i|.
 


We now note that we control the gradient of the proxy (,) down
to scales _*:


    (_B_r|∇(,)|^2)^1/2≲ 1   L/4≥ r≥_*.


Indeed, since by (<ref>), x_i+_i is a-harmonic in B_L, we have by Caccioppoli's estimate


    (_B_L/2|e_i+∇_i|^2)^1/2   ≲1/L(_B_L(x_i+_i-_B_L(x_i+_i))^2)^1/2
       ≲ 1+1/L(_B_L(_i-_B_L_i)^2)^1/2,


so that by (<ref>) for r=L, we get (<ref>) for r=L/2. The remaining range
of L/2≥ r≥_* follows since there, 
according to our hypotheses (<ref>) for ϕ^(L) and to (<ref>) and (<ref>) for σ̃, 
the medium a is well-behaved and thus admits the C^0,1-estimate, cf (<ref>), 
which applied to x_i+_i yields


    (_B_r|e_i+∇_i|^2)^1/2≲(_B_L/2|e_i+∇_i|^2)^1/2   L/2≥ r≥_*.


This establishes the -contribution to (<ref>), which in particular implies for the flux
_i=a(e_i+∇_i)


    (_B_r|_i|^2)^1/2≲ 1   L/2≥ r≥_*.


By the equation for _i with rhs given by the curl of _i, cf (<ref>), we obtain from
Caccioppoli's estimate


    (_B_r|∇_ijk|^2)^1/2≲1/r(_B_2r(_ijk)^2)^1/2
    +(_B_2r|_i|^2)^1/2,


so that we obtain the -part of (<ref>) from the -part of (<ref>) and
from (<ref>).




Next we argue that, again on the level of the gradients, the difference between proxy and
true corrector are small on scale L:


    (_B_L/2|∇(-ϕ)|^2)^1/2
    +(_B_L/4|∇(-σ)|^2)^1/2≲(_*/L)^β.


Indeed, by (<ref>) and (<ref>), -ϕ is (component-wise) a-harmonic in B_L.
Hence by the argument at the beginning of the proof of Lemma <ref> we obtain the ϕ-part
of (<ref>) from (<ref>) for R=L and (<ref>) for r=L. For the σ-part we note
that by (<ref>) and (<ref>), _i-σ_i solves
a (tensorial) Poisson equation with a rhs given by the curl of a∇(_i-ϕ_i). Hence by Caccioppoli
and the triangle inequality in L^2


    (_B_L/4|∇(-σ)|^2)^1/2
       ≲(_B_L||^2)^1/2
    +(_B_L|σ-_B_Lσ|^2)^1/2
    +(_B_L/2|∇(-ϕ)|^2)^1/2,


so that the σ-part of (<ref>) follows from the ϕ-part of it,
together with the σ-parts of both (<ref>) and (<ref>) for
R=r=L.



We now come to the central piece, namely that the differences between proxy and true corrector are
small down to scales _*:


    (_B_r|∇(-ϕ,-σ)|^2)^1/2≲(_*/L)^β   L/4≥ r≥_*.


As noticed above, _i-ϕ_i is a-harmonic in B_L so that passing from r=L/2,
cf (<ref>) and (<ref>) in Lemma <ref> with r_* replaced by _*, 
to L/2≥ r≥_* follows from the C^0,1-estimate already used
earlier, (<ref>). This settles the ϕ-part of (<ref>); in order to deal with the
σ-part, we need the full C^1,α-theory for a-harmonic functions from <cit.>
(for an α∈(0,1) fixed, say α=1/2), which holds for radii L≥ r≥_* since, as already
remarked above, the medium is well behaved there in the sense that there exist a tensor a_h^(L), scalar fields  ϕ_i^(L),
and skew symmetric vector fields σ̃_i, related by (<ref>) and satisfying the estimates
(<ref>) & (<ref>).


Applied to the a-harmonic function ϕ_i-_i in B_L this yields a vector ξ (which should carry an index i) 
such that 


    (_B_r|∇(ϕ_i-_i-ξ_j(x_j+_j))|^2)^1/2   ≲(r/L)^α(_B_L/2|∇(ϕ_i-_i)|^2)^1/2
     L/2≥ r≥_*   |ξ|   ≲(_B_L/2|∇(ϕ_i-_i)|^2)^1/2.


Recalling the definition of the fluxes q_i=a(e_i+∇ϕ_i) and _i=a(e_i+∇_i)
and inserting (<ref>), we obtain


    (_B_r|q_i-_i-ξ_j_j|^2)^1/2≲(r/L)^α(_*/L)^β   L/2≥ r≥_*,
      |ξ|≲(_*/L)^β.


By (<ref>) and (<ref>) we have that w:=σ_i-_i-ξ_j_j solves
a Poisson equation with the curl of q_i-_i-ξ_j_j as rhs. Hence the first part of (<ref>)
translates into


    -Δ w=∇· h    (_B_r|h|^2)^1/2≲(r/L)^α(_*/L)^β   L/2≥ r≥_*.


In the next paragraph, we shall argue that thanks to α>0, this implies


    (_B_r|∇ w|^2)^1/2≲(_*/L)^β+(_B_L/4|∇ w|^2)^1/2   L/4≥ r≥_*.


By definition of w and the triangle inequality, this yields


    (_B_r|∇(σ_i-_i)|^2)^1/2≲(_*/L)^β+(_B_L/4|∇(σ_i-_i)|^2)^1/2
       +|ξ_j|((_B_r|∇_j|^2)^1/2
    +(_B_L/4|∇_j|^2)^1/2)
       L/4≥ r≥_*.


Inserting the estimate on |ξ| from (<ref>), the localized estimates on
 from (<ref>), and the large-scale estimate on σ- from (<ref>),
we get the localized estimate on σ- stated in (<ref>).



It remains to argue that (<ref>) implies (<ref>). Like in the proof of Lemma <ref>,
we resort to a construction via a decomposition into dyadic annuli. For any dyadic multiple r
of _* with r≤L/4 we consider the Lax-Milgram solution w_r of


    -Δ w_r=∇·(I(B_2r-B_r)h);


the solution of the Poisson equation with rhs ∇·(I(B__*)h) is denoted by w__*/2
for notational consistency.
By the energy estimate for the Poisson equation and (<ref>) we have for all dyadic 
_*/2≤ r≤L/4


    (1/r^d∫|∇ w_r|^2)^1/2≲(r/L)^α(_*/L)^β.


From this and the (standard) mean-value property (note that unless r=_*/2, w_r is harmonic in B_r) 
we obtain for every radius R≥_* (the lower bound arises because of w__*/2)


    (_B_R|∇ w_r|^2)^1/2≲(r/L)^α(_*/L)^β.


Because of α>0 we obtain for w̃:=∑__*/2≤ r≤L/4w_r that


    (_B_R|∇w̃|^2)^1/2≲(_*/L)^β   R≥_*.


Since by construction, w-w̃ is harmonic in B_L/4, we obtain by the (standard) mean-value
property


    (_B_R|∇(w-w̃)|^2)^1/2≲(_B_L/4|∇(w-w̃)|^2)^1/2   R≤L/4.


By the triangle inequality in L^2, the two last inequalities yield (<ref>).



Equipped with the localized estimates on ∇(ϕ-,σ-) from (<ref>),
it is easy to conclude: Because of (<ref>), it is sufficient to establish (<ref>) for
L/4≥ r≥_*. In fact, we shall show that (<ref>) entails (<ref>) in this range.
This follows instantly from Poincaré's inequality in form of


    1/r(_B_r|(ϕ,σ)   -_B_r(ϕ,σ)|^2)^1/2
    -1/r(_B_r|(,)|^2)^1/2
       ≲(_B_r|∇(-ϕ,-σ)|^2)^1/2,


into which we plug (<ref>).








Since as a finite-range ensemble, ⟨·⟩ is in particular ergodic, <cit.>
applies, and yields the existence of ϕ_i and σ_i with the stated properties.
We set for abbreviation F_R:=(_B_R|(ϕ,σ)-_B_R(ϕ,σ)|^2)^1/2.



For given 1≤ s'<2, we start by extracting the stochastic bound


    ⟨exp(F_R/Clog^1/2R)^s'⟩≤ C   R≥ 2


from <cit.>.
Here C=C(d,λ,s') denotes a generic constant the value of which may change from line to line.
We focus on the case of d=2 (the result is stronger for d>2),
in which case the statement of <cit.> takes the form of


    log⟨exp((1/Clog(2+|x|)∫ G_1(·-x)|(ϕ,σ)
    -∫ G_1(ϕ,σ)|^2)^s'/2)⟩≤ 1


for all x∈ℝ^d,
where the Gaussian G_1(x):=(2π)^-d/2exp(-|x|^2/2) plays the role of a spatial averaging function.
We rewrite this 


    ⟨exp((1/Clog(2+|x|)∫ G_1(·-x)|(ϕ,σ)
    -∫ G_1(ϕ,σ)|^2+1)^s'/2)⟩≲ 1,


where ≲ stands for ≤ C(d,λ,s'). Since [0,∞)∈ F↦exp((F+1)^s'/2) is convex
for s'≥ 1, we obtain for arbitrary R≥ 2 from averaging over x∈ B_R


    ⟨exp((1/Clog R_B_R|(ϕ,σ)
    -∫ G_1(ϕ,σ)|^2+1)^s'/2)⟩≲ 1,


where now ∫ G_1(ϕ,σ) may be replaced by _B_R(ϕ,σ), so that this turns into (<ref>).



We now turn to part i) of the lemma, ie (<ref>), and fix β<1 and let the random radius
r_*∈[1,+∞] be minimal with (<ref>).
Because of 1/r'(_B_r'·)^1/2
≤(r/r')^d/2+11/r(_B_r·)^1/2,
r_* is dominated by the smallest dyadic radius, which for simplicity we call
again r_*, with the property that


    1/RF_R≤ 2^-d/2-1(r_*/R)^β   R>r_*.


Hence we have for any (deterministic) threshold r≥ 2


    ⟨ I(r_*≥ r)⟩≤∑_R>r ⟨ I(F_R≥ 2^-d/2-1 r^βR^1-β)⟩,


and obtain from (<ref>) for s<s'<2 (say s'=s+2/2) by Chebyshev's inequality


    ⟨ I(F_R≥ 2^-d/2-1 r^βR^1-β)⟩≤exp(-(r^β R^1-β/Clog^1/2r)^s'+1),


where C=C(d,λ,s,β) denotes a generic constant the value of which may change line by line.
We now claim that β<1 implies that the first term in the dyadic sum dominates, ie,


    ∑_R>r exp(-(r^β R^1-β/Clog^1/2r)^s')
    ≲exp(-(r/Clog^1/2r)^s'),


where ≲ stands short for ≤ up to a generic multiplicative constant C=C(d,λ,s,β).
Indeed, setting for abbreviation C(r):=(r/Clog^1/2r)^s',
this amounts to showing ∑_R>rexp(-C(r)(R/r)^s'(1-β)) ≲exp(-C(r)),
which because of C(r)≥ 1 (for r≫ 1, which means that r is larger than some constant only depending
on d,λ,s,β) reduces to the elementary ∑_n=1^∞exp(-((2^n)^s'(1-β)-1))≲ 1,
where we used β<1 and s'>0.
The combination of (<ref>), (<ref>), and (<ref>) yields


    ⟨ I(r_*≥ r)⟩≲exp(-(r/Clog^1/2r)^s')   r≫ 1,


which implies (<ref>) thanks to s<s'.



We finally turn to part ii) of the lemma and fix β<1 and 0<s<2(1-β). We have to show
the existence of L_0 such that


    ⟨ I(∃ R≥ L   F_R>R^1-β)⟩≤exp(-(L/L_0)^s)   L≥ L_0;


by the same argument as for part i) this reduces to


    ∑_R≥ L ⟨ I(F_R>2^-d/2-1R^1-β)⟩
       ≤exp(-(L/L_0)^s)   L≥ L_0.


From (<ref>) for s/1-β<s'<2 (say s':=1/2(s/1-β+2))
we obtain for each summand with R≥ 2


    ⟨ I(F_R≥ 2^-d/2-1 R^1-β)⟩≤exp(-(R^1-β/Clog^1/2R)^s'+1),


and as for (<ref>) we find for the sum


    ∑_R>L exp(-(R^1-β/Clog^1/2R)^s')
    ≲exp(-(L^1-β/Clog^1/2L)^s'),


Since s<s'(1-β) the last two statements imply (<ref>) for some L_0≲ 1.









According to Lemma <ref> ii) and with probability 1-exp(-(L/L_0)^s),
the hypothesis (<ref>) of Lemma <ref> is satisfied with r_*^(L)=1; the second hypothesis
(<ref>) is satisfied by assumption (<ref>) of the theorem. Hence by (<ref>) in Lemma <ref>, (<ref>) holds with
r_*≲max{r_*^(L),1}=r_*^(L), so that we may apply Proposition <ref> with r_*^(L)
playing the role of r_*. Hence (<ref>) turns into the desired (<ref>).











For notational simplicity, we drop the index i.
We start by collecting some estimates on ϕ_h'. Rewriting (<ref>) as
-∇· a_h'∇ϕ_h' =∇·(a_h'-a_h)e and noting that the rhs (a_h'-a_h)e
is bounded and supported in B_R(y), cf (<ref>), we obtain from the energy estimate 


    (1/R^d∫|∇ϕ_h'|^2)^1/2≲ 1.


It is convenient to introduce


    v'_h:=x· e+ϕ'_h    ∇· a_h'∇ v_h'(<ref>)=0,
 

and to reformulate (<ref>) as


    (_B_3R(y)|∇ v_h'|^2)^1/2≲ 1,


which in view of (<ref>), by Meyer's estimate, upgrades to


    (_B_2R(y)|∇ v_h'|^p)^1/p≲ 1   p=p(d,λ)>2.


Rewriting (<ref>) once more, this time as -∇· a_h∇ϕ_h' =∇·(a_h-a_h')∇ v_h',
and noting that a_h is a constant coefficient and that (a_h-a_h')∇ v_h' is supported in B_R(y), cf (<ref>),
we see that ϕ_h' decays like the gradient of the fundamental solution G_h. Hence, from (<ref>)
we obtain the estimate


    Rsup_x∈ B_2R^c(y)(|x-y|/R)^d+1|∇^2 v_h'|(<ref>)=Rsup_x∈ B_2R^c(y)(|x-y|/R)^d+1|∇^2ϕ_h'|≲ 1.


We may even get closer to the boundary of B_R^c(y) at the expense of a bad constant: For any boundary layer width 0<ρ≤ R we
obtain from (<ref>) and the fact that v_h' is constant-coefficient harmonic in B_R^c(y), cf (<ref>),


    sup_B_2R(y)-B_R+ρ(y)(ρ|∇^2 v_h'|+|∇ v_h'|)≲ (R/ρ)^d/2.




We now turn to the comparison of ϕ' and ϕ_h', at first in the strong topology on the level of gradients.
We carry this out in terms of the harmonic functions


    v':=e· x+ϕ'    ∇· a'∇ v'(<ref>)=0


and v'_h, by monitoring the error in the two-scale expansion, that is,



    w:=v'-(1+ζϕ_i∂_i)v_h'.


Here ζ denotes a cut-off function with


    ζ=0  B_R+ρ,  ζ=1  B_R+2ρ,  |∇ζ|≲ρ^-1


for a boundary layer width 0<ρ≤ R to be optimized at the end of the proof.
Based on the equations (<ref>), (<ref>), and on (<ref>), we obtain the following formula


    -∇· a'∇ w=∇·(h_far+h_near)  
       h_far:=(ϕ_ia-σ_i)∇(ζ∂_iv_h')
         h_near:=(1-ζ)(a'-a_h')∇ v_h';


which is the same (elementary) calculation as in Step 2 of the proof of <cit.>
and a slight variation of (<ref>).
For the proof of this lemma, we normalize (ϕ,σ) by _B_r_*(y)(y)(ϕ,σ) =0
so that by (<ref>) (with the origin replaced by y) we have


    1/r(_B_r(y)|(ϕ,σ)|^2)^1/2≲(r_*(y)/r)^β   r≥ r_*(y).




The near-field rhs h_near is supported in the thin annulus B_R+2ρ-B_R of thickness 2ρ, 
as a consequence of (<ref>), (<ref>) and (<ref>). With help of Meyer's estimate (<ref>),
we may capitalize on this by Hölder's inequality:


    (1/R^d∫|h_near|^2)^1/2≲ (ρ/R)^1/2-1/p.


In the complement of B_2R(y), in view of (<ref>) the far-field term assumes the simpler form 
h_far =(ϕ_ia-σ_i)∇∂_iv_h' and thus is easily estimated:


    (1/R^d∫_B_2R^c(y)|h_far|^2)^1/2
       (<ref>)≲1/R(1/R^d∫_B_2R^c(y)((R/|x-y|)^d+1|(ϕ,σ)|)^2)^1/2(<ref>)≲(r_*(y)/R)^β;


the last estimate can be seen by a decomposition into dyadic annuli; here we also use our assumption that R≥ r_*(y).
In B_2R(y) by (<ref>) the far-field term is supported in B_R+ρ^c(y) and estimated by
|h_far| ≲ |(ϕ,σ)|(|∇^2 v_h'|+1/ρ|∇ v_h'|). Hence we obtain from (<ref>)
and (<ref>) 


    (1/R^d∫_B_2R(y)|h_far|^2)^1/2≲(R/ρ)^1+d/2(r_*(y)/R)^β.


Using (<ref>), (<ref>), and (<ref>) in the energy estimate for (<ref>) we obtain


    (1/R^d∫|∇ w|^2)^1/2≲ (ρ/R)^1/2-1/p+(R/ρ)^1+d/2(r_*(y)/R)^β.




We are interested in the difference of the potentials, cf (<ref>), and in the difference of flux averages, 
cf (<ref>), and thus need to post-process (<ref>).
We first turn to the difference of the potentials and note that by definitions
(<ref>), (<ref>), and (<ref>) of v_h', v' and w we have
ϕ'-ϕ_h' =w+ζϕ_i∂_iv_h',
so that by Poincaré's inequality and the support condition on ζ, cf (<ref>),


    1/R(_B_2R(y)(ϕ'-ϕ_h'-_B_2R(y)(ϕ'-ϕ_h'))^2)^1/2
       ≲(1/R^d∫|∇ w|^2)^1/2
    +1/R(_B_2R(y)|ϕ|^2)^1/2sup_B_2R(y)-B_R+ρ(y)|∇ v_h'|,


so that from (<ref>) for the first rhs term and (<ref>) & (<ref>) for the second rhs term
we obtain


    1/R(_B_2R(y)(ϕ'-ϕ_h'   -_B_R(y)(ϕ'-ϕ_h'))^2)^1/2
       ≲ (ρ/R)^1/2-1/p+(R/ρ)^1+d/2(r_*(y)/R)^β.




We now turn to the differences of flux averages.
The first post-processing step is based on the formula


    a'∇ v'-a_h'∇ v_h'+∇·(ζ∂_iv_h'σ_i)=a'∇ w+h_far+h_near,


which in fact is the basis for the formula in (<ref>).
Hence from (<ref>), (<ref>), (<ref>), and (<ref>) we obtain


    (1/R^d∫|a'∇ v'-a_h'∇ v_h'   +∇·(ζ∂_iv_h'σ_i)|^2)^1/2
    ≲    (ρ/R)^1/2-1/p+(R/ρ)^1+d/2(r_*(y)/R)^β.


The second post-processing step consists in noting that the additional term
∇·(ζ∂_iv_h'σ_i) in the integrand has small average: From integration by parts,


    ∫∇·(ζ∂_iv_h'σ_i)=-∫∇·(ζ∂_iv_h'σ_i),


and since  is supported in B_2R(y), cf (<ref>), this rhs is bounded by


    sup|∇|   (∫_B_2R(y)|∇ v_h'|^2)^1/2(∫_B_2R(y)|σ|^2)^1/2
       (<ref>),(<ref>)≲1/R(_B_2R(y)|σ|^2)^1/2(<ref>)≲(r_*(y)/R)^β.


Combining this with (<ref>) yields 


    |∫(a'∇ v'-a_h'∇ v_h')|
    ≲ (ρ/R)^1/2-1/p+(R/ρ)^1+d/2(r_*(y)/R)^β.




The last task is to upgrade (<ref>) and (<ref>) by 
optimizing in the boundary layer thickness ρ. Indeed, choosing ρ/R =(r_*(y)/R)^β/1/2-1/p+1+d/2 
we obtain (<ref>) and (<ref>) with the potentially small but positive exponent
α=β1/2-1/p/1/2-1/p+1+d/2. 







We start by noting that


    (_B_2R(y)|e_j+∇ϕ_j'|^2)^1/2≲ 1.


Indeed, by (<ref>) in form of


    (_B_2R(y)|e_j+∇ϕ_j|^2)^1/2≲ 1


and the triangle inequality in L^2, it suffices to show


    (_B_2R(y)|∇(ϕ_j'-ϕ_j)|^2)^1/2≲ 1.


Since -∇· a'∇(ϕ_j'-ϕ_j)=∇·(a'-a)(e_j+∇ϕ_j), cf (<ref>) and (<ref>),
and since a'-a is supported in B_R(y), cf (<ref>), we have by the energy estimate


    (1/R^d∫|∇(ϕ_j'-ϕ_j)|^2)^1/2≲(_B_2R(y)|e_j+∇ϕ_j|^2)^1/2,


so that (<ref>) is a consequence of (<ref>).



Since a'-a and a_h'-a_h are supported in B_R(y), cf (<ref>) and (<ref>),
we may smuggle in the averaging function , cf (<ref>), into the
lhs of the desired (<ref>) so that it is enough to show


    |∫(e_j+∇ϕ_j')·(a'-a)(e_i+∇ϕ_i)
       -∫(e_j+∇ϕ_jh')·(a'_h-a_h)e_i|≲(r_*(y)/R)^α.


By integration by parts, and using (<ref>) and (<ref>), we obtain for the first integral in (<ref>)


    ∫(e_j+∇ϕ_j')·(a'-a)(e_i+∇ϕ_i)
       =∫(e_j+∇ϕ_j')· a'e_i
    -∫ e_j· a(e_i+∇ϕ_i)
       -∫(ϕ_i-_B_2R(y)ϕ_i)(e_j+∇ϕ_j')· a'∇
       +∫(ϕ_j'-_B_2R(y)ϕ_j')∇· a(e_i+∇ϕ_i).


On the first rhs term in (<ref>) we apply (<ref>) of Lemma <ref>:


    |∫(e_j+∇ϕ_j')· a'e_i-∫(e_j+∇ϕ_jh')· a'_he_i|≲(r_*(y)/R)^α.


Using (<ref>), the second rhs term in (<ref>) (without the minus sign) can be rewritten as


    ∫ e_j· a_h e_i-∫ e_j·(σ_i-_B_2R(y)σ_i)∇,


where the second contribution is estimated as follows


    |∫ e_j·(σ_i-_B_2R(y)σ_i)∇|
    (<ref>)≲1/R(_B_2R(y)|σ_i-_B_2R(y)σ_i|^2)^1/2(<ref>)≲(r_*(y)/R)^β,


so that we obtain


    |∫ e_j· a(e_i+∇ϕ_i)-∫ e_j· a_h e_i|≲(r_*(y)/R)^β.


The third rhs term in (<ref>) is estimated as follows


    |∫(ϕ_i-_B_2R(y)ϕ_i)(e_j+∇ϕ_j')· a'∇|
       (<ref>)≲1/R(_B_2R(y)(ϕ_i-_B_2R(y)ϕ_i)^2)^1/2(_B_2R(y)|e_j+∇ϕ_j'|^2)^1/2   (<ref>),(<ref>)≲(r_*(y)/R)^β.


We now turn to the last rhs term in (<ref>). We first apply (<ref>) in Lemma <ref>
to the effect of


    |∫(ϕ_j'-_B_2R(y)ϕ_j')∇· a(e_i+∇ϕ_i)
    -∫(ϕ_jh'-_B_2R(y)ϕ_jh')∇· a(e_i+∇ϕ_i)|
       (<ref>)≲1/R(_B_2R(y)(ϕ_j'-ϕ_jh'-_B_2R(y)(ϕ_j'-ϕ_jh'))^2)^1/2(_B_2R(y)|e_i+∇ϕ_i|^2)^1/2
       (<ref>),(<ref>)≲(r_*(y)/R)^α.


We then note that by (<ref>), the skew symmetry of σ_i and two integration by parts


    ∫(ϕ_jh'-_B_2R(y)ϕ_jh')∇· a(e_i+∇ϕ_i)
       =-∫∇ϕ_jh'· a_he_i-∫∇·(σ_i-_B_2R(y)σ_i)∇ϕ_jh',


so that


    |∫(ϕ_jh'-_B_2R(y)ϕ'_jh)∇· a(e_i+∇ϕ_i)
    +∫∇ϕ_jh'· a_he_i|
       (<ref>)≲1/R(_B_2R(y)|σ_i-_B_2R(y)σ_i|^2)^1/2(_B_2R(y)|∇ϕ_jh'|^2)^1/2(<ref>),(<ref>)≲(r_*(y)/R)^β.


The combination of (<ref>) and (<ref>) yields for the last term in (<ref>):



    |∫(ϕ_j'-_B_2R(y)ϕ_j')
    ∇· a(e_i+∇ϕ_i)+∫∇ϕ_jh'· a_he_i|
    ≲(r_*(y)/R)^α.


Inserting the four estimates (<ref>), (<ref>), (<ref>) and (<ref>) into (<ref>),
we obtain (<ref>).







Starting point is Lemma <ref>, more precisely (<ref>) for L:=1/2|y| and (<ref>) in form of


    (_B_L^c|∇(u-(1+(ϕ_i-_B_r_*ϕ_i)∂_i)u_h)|^2)^1/2≲ (ℓ/|y|)^d (r_*/|y|)^β.


The first post-processing step is to replace _B_r_*ϕ_i by _B_r_*(y)ϕ_i
in (<ref>). Indeed, by (<ref>) (for r=L=1/2|y| and r'=r_*) we have


    1/L|_B_L   (ϕ,σ)   -_B_r_*   (ϕ,σ)|
    +1/L|_B_L(y)(ϕ,σ) -_B_r_*(y)(y)(ϕ,σ)|
       ≲ (max{r_*,r_*(y)}/|y|)^β.


As a direct consequence of (<ref>) (for r=3L) and B_L∪ B_L(y)⊂ B_3L we have 


    1/L|_B_L(ϕ,σ)-_B_L(y)(ϕ,σ)|
    ≲(r_*/|y|)^β,


so that (<ref>) implies


    1/L|_B_r_*(ϕ,σ)-_B_r_*(y)(y)(ϕ,σ)|
    ≲ (max{r_*,r_*(y)}/|y|)^β.


We combine this with (<ref>) in form of


    sup_B_L^c|∇ u_h|+sup_B_L^cL|∇^2u_h|≲(ℓ/|y|)^d


to the desired


    (_B_L^c|∇(u-(1+(ϕ_i-_B_r_*(y)(y)ϕ_i)∂_i)u_h)|^2)^1/2≲ (ℓ/|y|)^d (max{r_*,r_*(y)}/|y|)^β.


Hence wlog we assume _B_r_*(y)(y)ϕ_i=0 so that on the one hand, the above simplifies to


    (_B_L^c|∇(u-(1+ϕ_i∂_i)u_h)|^2)^1/2≲ (ℓ/|y|)^d (max{r_*,r_*(y)}/|y|)^β


and (<ref>), with the origin replaced by y, assumes the form of


    1/r(_B_r(y)|(ϕ,σ)|^2)^1/2≲(r_*(y)/r)^β   r≥ r_*(y).




Following Step 6 in the proof of <cit.>, we now localize (<ref>) around y, 
making use of (<ref>).
To this purpose, we appeal once more to the formula (<ref>) for the error in the two-scale convergence


    -∇· a∇(u-(1+ϕ_i∂_i)u_h)=∇· h    
    h:=(ϕ_ia-σ_i)∇∂_iu_h.


The combination of (<ref>) and (<ref>) yields for the rhs


    (_B_r(y)|h|^2)^1/2≲(ℓ/|y|)^dr/L(r_*(y)/r)^β   L≥ r≥ r_*(y).


We now argue, starting from the large-scale anchoring of (<ref>) in form of


    (_B_L(y)|∇(u-(1+ϕ_i∂_i)u_h)|^2)^1/2≲ (ℓ/|y|)^d (max{r_*,r_*(y)}/|y|)^β,


that (<ref>) allows for the desired localization to our scale of interest R≥ r_*(y):


    (_B_R(y)|∇(u-(1+ϕ_i∂_i)u_h)|^2)^1/2≲(ℓ/|y|)^d(max{r_*,r_*(y)}/|y|)^β.


To this purpose, for every radius r=R,2R,⋯ with r≤ L  
we consider the Lax-Milgram solution of


    -∇· a∇ u_r=∇·(I(B_r(y)-B_r/2(y))h)


with the understanding that for r=R, the rhs is given by I(B_R(y))h. From the
energy estimate and (<ref>) we obtain 


    (1/r^d∫|∇ u_r|^2)^1/2≲(ℓ/|y|)^dr/L(r_*(y)/r)^β.


Since for r≥ 2R, u_r is a-harmonic in B_r/2(y), we may apply the C^0,1-estimate, cf (<ref>),
to localize the above to


    (_B_R(y)|∇ u_r|^2)^1/2≲(ℓ/|y|)^dr/L(r_*(y)/r)^β.


Since the exponent 1-β of r is positive, ũ:=∑_L≥ r≥ Ru_r satisfies


    (_B_R(y)|∇ũ|^2)^1/2≲(ℓ/|y|)^d(r_*(y)/|y|)^β.


Likewise, we obtain directly from (<ref>)


    (1/L^d∫|∇ũ|^2)^1/2≲(ℓ/|y|)^d(r_*(y)/|y|)^β.


Since by construction, u-(1+ϕ_i∂_i)u_h-ũ 
is harmonic in B_L/2(y), we may apply the C^0,1-estimate
to the effect of


    (_B_R(y)|∇(u-(1+ϕ_i∂_i)u_h-ũ)|^2)^1/2
       ≲(_B_L(y)|∇(u-(1+ϕ_i∂_i)u_h-ũ)|^2)^1/2.


By the triangle inequality in L^2, we obtain (<ref>) from 
(<ref>), (<ref>), (<ref>), and (<ref>).



We now post-process (<ref>), making use of r_*,r_*(y)≤ R, to


    (_B_R(y)|∇ u-∂_i u_h(y)(e_i+∇ϕ_i)|^2)^1/2≲(ℓ/|y|)^d(R/|y|)^β


and also note for later purpose


    (_B_R(y)|∇ u|^2)^1/2,
    (_B_R(y)|∂_i u_h(y)(e_i+∇ϕ_i)|^2)^1/2≲(ℓ/|y|)^d.


In deriving (<ref>) from (<ref>), 
we first replace ∇(1+ϕ_i∂_i)u_h by ∂_iu_h(e_i+∇ϕ_h),
which means that by the triangle inequality in L^2, we have to estimate ϕ_i∇∂_iu_h.
We have


    (_B_R(y)|ϕ_i∇∂_iu_h|^2)^1/2≤sup_B_R(y)|∇^2u_h|(_B_R(y)|ϕ|^2)^1/2;


by B_R(y)⊂ B_L^c and (<ref>), the first factor is estimated by 1/|y|(ℓ/|y|)^d.
By (<ref>), the second factor is estimated by R(r_*(y)/R)^β.
Hence because of R≤|y|, this contribution is contained in the rhs of (<ref>).
We now replace ∂_iu_h(e_i+∇ϕ_h) by ∂_iu_h(y)(e_i+∇ϕ_h). By
the triangle inequality in L^2, we are lead to estimating 


    (_B_R(y)|(∂_iu_h-∂_iu_h(y))(e_i+∇ϕ_i)|^2)^1/2
       ≤ Rsup_B_R(y)|∇∂_iu_h|(_B_R(y)|e_i+∇ϕ_i|^2)^1/2.


As above, the first factor is estimated by R/|y|(ℓ/|y|)^d. According to
(<ref>) (with y playing the role of the origin and using R≥ r_*(y)), the
second factor is estimated by 1. Because of β≤ 1, 
also this contribution is contained in the rhs of (<ref>). The same argument leads to the
second estimate in (<ref>). The first estimate in (<ref>) follows from the second
one and (<ref>) via the triangle inequality in L^2.



We now turn to the sensitivity estimate and consider


    w:=u'-u,    -∇· a'∇ w=∇· h'    
    h':=(a'-a)∇ u.


We want to apply the localized homogenization error estimate of <cit.>, 
for the medium given by a'. This means comparing w to the solution of 


    -∇· a_h∇ w_h=∇· (h_i'(e_i+∂_iϕ')).


Note that by assumption, and by Lemma <ref> we have r_*'(y),r_*'≲ R.
Since h' is supported in B_R(y), by <cit.> 
(for the medium a', R playing the role of r_* there, the roles of y and the origin exchanged,
and β replacing 1-α there)  we have


    (_B_R|∇ w-∂_kw_h(e_k+∇ϕ_k')|^2)^1/2≲(R/|y|)^d+β(_B_R(y)|h'|^2)^1/2.


(As explained in the proof of Lemma <ref>, the logarithm in <cit.>
can be avoided.)
Among other ingredients, this estimate is based on the following estimate of w_h,
cf Step 2 of the proof of <cit.>,


    sup_B_L^c(y)(|∇ w_h|+L|∇^2 w_h|)
    ≲(R/L)^d(_B_R(y)|h'|^2)^1/2.


In view of the definition of h', cf (<ref>), and (<ref>), these two estimates turn into


    (_B_R|∇ w-∂_kw_h(e_k+∇ϕ_k')|^2)^1/2≲(ℓ/|y|)^d(R/|y|)^d+β,
    sup_B_L^c(y)(|∇ w_h|+L|∇^2 w_h|)
    ≲(ℓ/|y|)^d(R/L)^d.


Like for the passage from (<ref>) to (<ref>), (<ref>) may be post-processed 
with help of (<ref>) to


    (_B_R|∇ w-∂_kw_h(0)(e_k+∇ϕ'_k)|^2)^1/2≲(ℓ/|y|)^d(R/|y|)^d+β.


In this argument, we just have to replace the medium a by the medium a' and y by the origin.

















We continue with post-processing and argue that we may replace ϕ_k' by ϕ_k in (<ref>):


    (_B_R|∇ w-∂_kw_h(0)(e_k+∇ϕ_k)|^2)^1/2≲(ℓ/|y|)^d(R/|y|)^d+β.


In fact, we claim that the error term is of (substantially) higher order:


    |∂_kw_h(0)|(_B_R|∇(ϕ_k'-ϕ_k)|^2)^1/2≲(ℓ/|y|)^d(R/|y|)^2d,


which by (<ref>) reduces to


    (_B_R|∇(ϕ_k'-ϕ_k)|^2)^1/2≲(R/|y|)^d.


The latter can be seen noting that 


    -∇· a∇(ϕ_k'-ϕ_k)=∇·(a'-a)(e_k+∇ϕ_k').


From (<ref>) we learn at first that ϕ_k'-ϕ_k is a-harmonic in B_R^c(y) so
that by the C^0,1-estimate, cf (<ref>), and since R≥ r_*


    (_B_R|∇(ϕ_k'-ϕ_k)|^2)^1/2≲(_B_L|∇(ϕ_k'-ϕ_k)|^2)^1/2,


where we recall our abbreviation L=1/2|y|. Moreover, since the rhs of (<ref>) is in
divergence form, we have by the dualized C^0,1-estimate (see Step 5 in the proof of 
<cit.> for such a duality argument) and since R≥ r_*(y)


    (1/L^d∫_B_L^c(y)|∇(ϕ_k'-ϕ_k)|^2)^1/2≲(R/L)^d(1/R^d∫_B_R^c(y)|∇(ϕ_k'-ϕ_k)|^2)^1/2.


Finally, by the energy estimate for (<ref>) we have


    (1/R^d∫|∇(ϕ_k'-ϕ_k)|^2)^1/2≲(_B_R(y)|e_k+∇ϕ_k'|^2)^1/2≲ 1,


where in the very last estimate, we've used (<ref>) (with (a,0) replaced by (a',y)),
which we may since by Lemma <ref> we have r_*'(y)≲ R. Since by definition of L,
B_L⊂ B_L^c(y) the last three estimates combine to (<ref>).



In the remainder of the proof we argue that we may replace w_h by a more
explicit expression. More precisely, in order to pass from (<ref>) to (<ref>) we have to
replace ∂_kw_h(0) by


    ∂_j∂_kG_h(-y)∂_iu_h(y)∫(e_j+∇ϕ_j')·(a'-a)(e_i+∇ϕ_i),


and then appeal to the definition (<ref>) of δ a_ij. The basis for this is the representation


    ∂_kw_h(0)=∫∂_j∂_kG_h(-x)((e_j+∇ϕ_j')·(a'-a)∇ u)(x)dx,


which is a consequence of the definition (<ref>) of w_h, yielding the representation


    ∂_kw_h(0)   =∫∇∂_kG_h(-x)·(h_i'(e_i+∂_iϕ')))(x)dx
       =∫∂_j∂_kG_h(-x)(h'·(e_j+∇ϕ_j')))(x)dx,


into which we insert the definition (<ref>) of h'. We split the passage from (<ref>) to
(<ref>) into two steps. In view of (_B_R|e_k+∇ϕ_k|^2)^1/2≲ 1,
which is a consequence of R≥ r_*, cf (<ref>), it is enough to show


    R^dsup_x∈ B_R(y)|∂_j∂_kG_h(-x)-∂_j∂_kG_h(-y)|
    _B_R|e_j+∇ϕ_j'||∇ u|   ≲(ℓ/|y|)^d(R/|y|)^d+1,
    
    R^d|∂_j∂_kG_h(-y)|
    _B_R|e_j+∇ϕ_j'||∇ u-∂_iu_h(y)(e_i+∇ϕ_i)|
       ≲(ℓ/|y|)^d(R/|y|)^d+β,


which by the obvious decay properties of G_h reduces to


    _B_R|e_j+∇ϕ_j'||∇ u|   ≲(ℓ/|y|)^d,
    _B_R|e_j+∇ϕ_j'||∇ u-∂_iu_h(y)(e_i+∇ϕ_i)|
       ≲(ℓ/|y|)^d(R/|y|)^β.


By the Cauchy-Schwarz inequality, because of (_B_R|e_j+∇ϕ_j'|^2)^1/2≲ 1,
which is a consequence of R≳ r_*', this reduces to (<ref>) and (<ref>).







We may apply Lemma <ref> to the heterogeneous medium a replaced by the homogeneous a_h;
in which case we may choose (ϕ,σ)≡ 0 and thus in particular r_*=r_*(y)=0, whereas in view of
(<ref>) and (<ref>), ϕ' is being replaced by ϕ_h' and thus δ a by δ a_h, cf (<ref>). 
An inspection of the proof of Lemma <ref> shows that me may replace w=u'-u by w_h:=u_h'-u_h,
since the crucial property of that difference was, when translated to the homogeneous medium,


    w_h:=u'_h-u_h,    -∇· a'_h∇ w_h=∇· h'_h    
    h'_h:=(a'_h-a_h)∇ u_h,


cf (<ref>) for w, which for w_h follows immediately from the 
definition (<ref>) of u_h'. Hence (<ref>) assumes the form


    (_B_R|∇(u'_h-u_h)
    -∂_iu_h(y)δ a_hij∂_j∂_kG_h(-y)e_k|^2)^1/2≲(ℓ/|y|)^d(R/|y|)^d+β.


Into this estimate, we insert (<ref>) in form of 1/R^d|δ a_ij-δ a_hij|
≲(r_*(y)/R)^α, combined with (<ref>) below to


    (_B_R|∂_k(u'_h-u_h)   
    -∂_iu_h(y)δ a_ij∂_j∂_kG_h(-y)|^2)^1/2
       ≲(ℓ/|y|)^d(R/|y|)^d((R/|y|)^β+(r_*(y)/R)^α).


In combination with (<ref>) in form of (_B_R|e_k+∇ϕ_k|^2)^1/2≲ 1,
this assumes the form


    (_B_R|∂_k(u'_h-u_h)(e_k+∇ϕ_k)   
    -∂_iu_h(y)δ a_ij∂_j∂_kG_h(-y)(e_k+∇ϕ_k)|^2)^1/2
       ≲(ℓ/|y|)^d(R/|y|)^d((R/|y|^β+(r_*(y)/R)^α).


By the triangle inequality in L^2 with (<ref>) in its original form, we obtain (<ref>).






Clearly, it suffices to establish (<ref>) at the origin (which here plays the role of the general point), that is,


    r_*'≲ r_*+r_*(y)+R.


According to the definition of r_* as the minimal radius with the property (<ref>), 
by the triangle inequality in L^2, it is enough to show


    1/r(_B_r|(ϕ'-ϕ,σ'-σ)   -_B_r(ϕ'-ϕ,σ'-σ)|^2)^1/2
       ≲(r_*(y)+R/r)^β   r≥ r_*(y)+R,


which by Poincaré's inequality follows from


    (_B_r|∇(ϕ'-ϕ,σ'-σ)|^2)^1/2≲(r_*(y)+R/r)^β   r≥ r_*(y)+R.


In view of β<1≤d/2 and of (_B_rf^2)^1/2 
≲(R/r)^d/2(1/R^d∫ f^2)^1/2
this follow from


    (1/R^d∫|∇(ϕ'-ϕ,σ'-σ)|^2)^1/2≲ (1+r_*(y)/R)^d/2,


which we shall establish now.



In order to tackle (<ref>), we fix a coordinate direction i=1,⋯,d and suppress the index i in (ϕ_i,σ_i)
and (ϕ_i',σ_i'). Since -∇· a'∇(ϕ'-ϕ) =∇·(a'-a)(e+∇ϕ), we have
by the energy inequality


    (1/R^d∫|∇(ϕ'-ϕ)|^2)^1/2≲(_B_R(y)|e+∇ϕ|^2)^1/2.


We now note that Caccioppoli's estimate (<ref>) implies


    (_B_R(y)|e+∇ϕ|^2)^1/2≲ (1+r_*(y)/R)^d/2;


This is immediate in case of R≥ r_*(y) and follows from 
(_B_R(y)f^2)^1/2 ≤(R/r_*(y))^d/2 (_B_r_*(y)f^2)^1/2
in the other case. The combination of (<ref>) and (<ref>) yields (<ref>)
for the ϕ-part.



We now turn to the σ-part of (<ref>) and note that -Δ(σ_jk'-σ_jk) =∂_j(q'-q)_k
-∂_k(q'-q)_j where q'-q =a'∇(ϕ'-ϕ) +(a'-a)(e+∇ϕ). Hence we obtain by
the energy estimate


    (1/R^d∫|∇(σ'-σ)|^2)^1/2≲(1/R^d∫|∇(ϕ'-ϕ)|^2)^1/2+(_B_R(y)|e+∇ϕ|^2)^1/2,


so that it remains to appeal to (<ref>) and (<ref>) once more.









Motivated by (<ref>), let us introduce for abbreviation


    δ a_hij:=∫(e_j+∇ϕ_jh')·(a_h'-a_h)e_i


so that the tensor δ a_h={δ a_hij}_ij=1,⋯,d satisfies


    ξ̃·δ a_hξ=∫(ξ+∇ϕ_hξ')·(a_h'-a_h)ξ̃    ϕ_hξ':=ξ_jϕ_jh'.


We first claim that the non-degeneracy condition (<ref>) translates into


    -ξ·δ a_hξ≥λ_0|B_R||ξ|^2   ξ∈ℝ^d.


Indeed, since by definition (<ref>) of ϕ_hξ' and (<ref>) we have
-∇· a_h'(ξ+∇ϕ_hξ')=0, which we rewrite as
-∇· a_h'∇ϕ_hξ'=∇·(a_h'-a_h)ξ, we obtain from (<ref>) the representation


    ξ̃·δ a_hξ=∫ξ·(a_h'-a_h)ξ̃-∫∇ϕ_hξ'· a_h'∇ϕ_hξ̃'


and thus the inequality for the quadratic part


    ξ·δ a_hξ≤∫_B_R(y)ξ·(a_0-a_h)ξ.


Hence (<ref>) indeed translates into (<ref>).


We start by post-processing Lemma <ref> by getting rid of the constraint that R≥ r_*,r_*(y).
Indeed, in case of R≤max{r_*,r_*(y)} we apply Lemma <ref> with max{r_*,r_*(y)} playing
the role of R so that (<ref>) turns into 


    (_B_max{r_*,r_*(y)}|∇(u'-u)
       -∂_iu_h(y)δ a_ij∂_j∂_kG_h(-y)(e_k+∇ϕ_k)|^2)^1/2
       ≲ (ℓ/|y|)^d(max{r_*,r_*(y)}/|y|)^d+β,


which because of _B_R ≤(max{r_*,r_*(y)}/R)^d_B_max{r_*,r_*(y)} yields


    (_B_R|∇(u'-u)
       -∂_iu_h(y)δ a_ij∂_j∂_kG_h(-y)(e_k+∇ϕ_k)|^2)^1/2
       ≲ (ℓ/|y|)^d
    (R/|y|)^d+β
    (max{r_*,r_*(y)}/R)^d+β+d/2.


Hence we obtain in either case 


    (_B_R|∇(u'-u)
       -∂_iu_h(y)δ a_ij∂_j∂_kG_h(-y)(e_k+∇ϕ_k)|^2)^1/2
       ≲ (ℓ/|y|)^d
    (R/|y|)^d+β(1+r_*+r_*(y)/R)^d+β+d/2.




We now pass from a strong-L^2-estimate to an estimate of averages, which allows us to get rid
of the corrector in (<ref>):


    |∫ω(∇(u'-u)
    -∂_iu_h(y)δ a_ij∂_j∂_kG_h(-y)e_k)|
       ≲(ℓ/|y|)^d(R/|y|)^d((R/|y|)^β+(r_*/R)^β)
    (1+r_*+r_*(y)/R)^3/2d+β.


Indeed, this follows from the obvious estimate |∫ω g| ≲(_B_R|g|^2)^1/2
applied to the field g :=∇(u'-u)-∂_iu_h(y)δ a_ij∂_j∂_kG_h(-y)(e_k+∇ϕ_k),
combined with an argument that the contribution from ∫ω∇ϕ_k is negligible. For the latter we first
note that


    R^d|∇^2 G_h(-y)|≲ (R/|y|)^d,  
    |∇ u_h(y)|(<ref>)≲(ℓ/|y|)^d


and that


    1/R^d|δ a_ij|≲ (1+r_*(y)/R)^d/2.


In case of R ≥ max{r_*(y),r_*'(y)}, (<ref>) follows immediately from the
definition (<ref>), the fact that the integral is restricted to B_R(y) by (<ref>),
and the Caccioppoli estimate (<ref>) with y playing the role of the origin and
also applied to the perturbed medium a'. In case of R ≤ max{r_*(y),r_*'(y)}
we argue as in the previous paragraph, that is, we apply the above estimate with 
max{r_*(y),r_*'(y)} playing the role of R, which is legitimate in view of the only
constraint (<ref>) on R. This yields (<ref>) with a rhs given by
(1+max{r_*(y),r_*'(y)}/R)^d/2, so that it remains to appeal to Lemma <ref>.
We now may turn to ∫ω∇ϕ_k itself: By integration by parts we obtain


    |∫ω∇ϕ_k|=|∫(ϕ_k-_B_Rϕ_k)∇ω|
    (<ref>)≲1/R(_B_R(ϕ_k-_B_Rϕ_k)^2)^1/2.


In case of R≥ r_*, the rhs is controlled by (r_*/R)^β according to (<ref>).
In the other case, we argue as before to obtain an estimate by (r_*/R)^d/2+1. 
Hence we obtain in either case


    1/R(_B_R(ϕ_k-_B_Rϕ_k)^2)^1/2≲(r_*/R)^β(1+r_*/R)^d/2+1-β.


Inserting (<ref>) into (<ref>) and combining with (<ref>) and (<ref>) yields


    |∂_iu_h(y)δ a_ij∂_j∂_kG_h(-y)∫ω∇ϕ_k|
       ≲
    (ℓ/|y|)^d(R/|y|)^d(1+r_*(y)/R)^d/2(r_*/R)^β(1+r_*/R)^d/2+1-β,


which completes the argument for (<ref>) since (1+r_*(y)/R)^d/2(1+r_*/R)^d/2+1-β
≲ (1+r_*+r_*(y)/R)^3/2d+β.



We now bring Corollary <ref> into play, which in our abbreviations (<ref>) and (<ref>) reads
1/R^d|δ a_ij-δ a_hij|≲(r_*(y)/R)^α in case of R≥ r_*(y).
In the other case, we use the triangle inequality 1/R^d|δ a_ij-δ a_hij| ≤1/R^d|δ a_ij| +1/R^d|δ a_hij| recall the argument for (<ref>) that gave
1/R^d|δ a_ij| ≲(max{r_*(y),r_*'(y)}/R)^d/2
≲(1+r_*(y)/R)^d/2, while we easily get
1/R^d|δ a_hij|≲ 1 from (<ref>). Hence in either case, we have


    1/R^d|δ a_ij-δ a_hij|≲(r_*(y)/R)^α(1+r_*(y)/R)^d/2-α.


This allows us to upgrade (<ref>) to


    |∫ω∇(u'-u)-∂_iu_h(y)δ a_hij∂_j∂_kG_h(-y)e_k|
       ≲(ℓ/|y|)^d(R/|y|)^d((R/|y|)^β+(r_*/R)^β+(r_*(y)/R)^α)
    (1+r_*+r_*(y)/R)^3/2d+β,


where we absorbed both (r_*/R)^β and (r_*(y)/R)^α into
(max{r_*,r_*(y)}/R)^α.



Finally, we will substitute u_h by (∫ g)·∇ G_h in (<ref>).














Starting point is the following representation of u_h which follows from its definition in (<ref>) and (<ref>):


    u_h(x)   =∫∇ G_h(x-x')· g(x')dx'+(∫∇ϕ_m· g)∂_mG_h(x)  
    ∂_i u_h(y)   =∫∂_i∂_mG_h(y-x') g_m(x')dx'+(∫∇ϕ_m· g)∂_i∂_mG_h(y).


Hence we have


    |∂_i u_h(y)-(∫(e_m+∇ϕ_m)· g)∂_i∂_mG_h(y)|
       (<ref>)≤ℓ^dsup_x'∈ B_ℓ|∂_i∂_mG_h(y-x')-∂_i∂_mG_h(y)|
    ≲(ℓ/|y|)^d+1.


In order to argue that the contribution from ∫∇ϕ_m· g is negligible we argue as for
(<ref>) and (<ref>) to obtain


    |∫∇ϕ_m· g|
    
    ≲ℓ^d(r_*/ℓ)^β(1+r_*/R)^d/2+1-β.


Hence with help of |∂_i∂_mG_h(y)|≲1/|y|^d we may upgrade (<ref>) to


    |∂_i u_h(y)-(∫ g_m)∂_i∂_mG_h(y)|
    ≲(ℓ/|y|)^d(ℓ/|y|+(r_*/ℓ)^β(1+r_*/R)^d/2+1-β).


Since |δ a_hij|≲ R^d, see above, and |∂_j∂_kG_h(-y)|≲1/|y|^d
this allows to pass from (<ref>) to (<ref>).





In view of (<ref>), in order to conclude, it suffices to establish that for some non-zero ξ∈ℝ^d


    |(∫ g_m)   ∂_m∂_iG(y)δ a_hij∂_j∂_kG_h(-y)ξ_k|
    ≳|∫ĝ||ξ|(ℓ/|y|)^d(R/|y|)^d
         |y/|y|-∫ĝ/|∫ĝ||≪ 1.


By homogeneity and rotational invariance, it suffices to show this for ∫ĝ=e_1, in which case we
show the above for ξ=e_1 so that by ∫ g=ℓ^d∫ĝ this reduces to


    |ℓ^d∂_1∂_iG_h(y)δ a_hij∂_j∂_1G_h(-y)|
       ≳(ℓ/|y|)^d(R/|y|)^d
         |y/|y|-e_1|≪ 1.


Applying (<ref>) with ξ_i =∂_1∂_iG_h(y) =∂_1∂_iG_h(-y) we see
that it is enough to establish


    |∇∂_1 G_h(y)|^2≳(1/|y|)^2d    |y/|y|-e_1|≪ 1,


which by the (-d)-homogeneity of ∇^2G_h reduces to


    |∇∂_1 G_h(ŷ)|
    ≳ 1    |ŷ-e_1|≪ 1,


and by the continuity of ∇^2G_h away from the origin to the obvious statement of


    0≠∇∂_1G_h(e_1)=(1-d)∇ G_h(e_1),


where the last identity follows from the (1-d)-homogeneity of ∇ G_h.





In view of the independence properties of the Poisson point process, 
it follows from elementary probability theory that for any random variable F
(ie a function F=F(X) of the point configuration X) we have


    ⟨(⟨ F|B_R-1(y)⟩-⟨ F⟩)^2⟩
    =1/2⟨⟨(⟨ F'⟩_out-⟨ F⟩_out)^2
    ⟩_in⟩_in',


where ⟨·⟩_out:=⟨·|B_R-1(y)⟩,
where ⟨·⟩_in and ⟨·⟩_in' denote two independent
copies of the Poisson process restricted to B_R-1(y), 
and where F':=F(X') with X' denoting the realization of
the Poisson point process that arises from concatenating X_|B_R-1^c(y) and X'_|B_R-1(y).
We apply this to F=∫ω∇ u and note that F'=∫ω∇ u',
provided a':=A(·,X'), cf Definition <ref>. By definition of X and the locality assumption (<ref>) on A, 
this is consistent with the second condition in (<ref>); 
the first condition in (<ref>) comes for free by setting a_0:=a'_|B_R(y).
Hence in order to establish (<ref>) (with R replaced by R-1), we have to show for R≫ 1 that


    ⟨⟨(⟨∫ω∇ u'⟩_out-⟨∫ω∇ u⟩_out)^2
    ⟩_in⟩_in'
    ≳((ℓ/|y|)^d(1/|y|)^d)^2,


where by homogeneity and rotational invariance, we have assumed wlog


    ∫ĝ=e_1.




Since by the above remark, we have that a'=a outside of B_R(y), we may apply Corollary <ref>.
More precisely, we apply ⟨·⟩_out to (<ref>) and obtain by Jensen's inequality
in probability


    |⟨∫ω∇ u'⟩_out-⟨∫ω∇ u⟩_out
       -(∫ g_m)∂_m∂_iG_h(y)⟨δ a_hij⟩_out∂_j∂_kG_h(-y)e_k|≲⟨ err⟩_out,


where we have set for abbreviation


    err:=(ℓ/|y|)^d(R/|y|)^d
       ×((R/|y|)^β+(r_*/R)^β+(r_*(y)/R)^α
    +ℓ/|y|+(r_*/ℓ)^β)(1+r_*+r_*(y)/R)^3/2d+β.


Thanks to the normalization (<ref>), which by (<ref>) turns into
∫ g=ℓ^de_1, we obtain from just considering the first component, and using the symmetry and homogeneity
of ∇^2G_h


    |(∫ g_m)∂_m∂_iG_h(y)⟨δ a_hij⟩_out∂_j∂_kG_h(-y)e_k|
       ≥ (ℓ/|y|)^d(R/|y|)^d
    |∂_i∂_1G_h(ŷ)⟨1/R^dδ a_hij⟩_out∂_j∂_1G_h(ŷ)|.


Hence under the proviso


    -ξ·⟨1/R^dδ a_hij⟩_outξ≥λ_0|ξ|^2   ξ∈ℝ^d


for some λ_0>0 to be chosen later, we obtain


    |(∫ g_m)∂_m∂_iG_h(y)⟨δ a_hij⟩_out∂_j∂_kG_h(-y)e_k|
    ≥λ_0(ℓ/|y|)^d(R/|y|)^d|∇^2 G_h(ŷ)e_1|^2,


which by the invertibility of the Hessian matrix ∇^2G_h(ŷ) and its continuity on {|ŷ|=1} implies


    |(∫ g_m)∂_m∂_iG_h(y)⟨δ a_hij⟩_out∂_j∂_kG_h(-y)e_k|
    ≳λ_0(ℓ/|y|)^d(R/|y|)^d.


Therefore, always under the proviso (<ref>), we obtain the following lower bound from (<ref>)


    λ_0(ℓ/|y|)^d(R/|y|)^d
    ≲|⟨∫ω∇ u'⟩_out-⟨∫ω∇ u⟩_out|
    +⟨ err⟩_out.




We note that by definition, cf (<ref>), (<ref>) and (<ref>), 
δ a_hij only depends on a_0_|B_R(y)=a'_|B_R(y)=A(·,X)_|B_R(y), which by the locality of A,
cf (<ref>), implies that δ a_hij and ⟨δ a_hij⟩_out are
in particular independent of X_|B_R-1(y). Hence the proviso (<ref>) is not affected by
(squaring and) applying ⟨·⟩_in to (<ref>): Under the proviso (<ref>), we have


    (λ_0(ℓ/|y|)^d(R/|y|)^d)^2
    ≲⟨|⟨∫ω∇ u'⟩_out-⟨∫ω∇ u⟩_out|^2⟩_in
    +⟨ err^2⟩,


where we've used Jensen's inequality in probability on the last term. In particular, we may multiply this estimate with
the characteristic function I(-1/R^d⟨δ a_h⟩_out≥λ_0 id) of the event
(<ref>) and then apply ⟨·⟩_in', which yields


    prob(λ_0(ℓ/|y|)^d(R/|y|)^d)^2
       ≲⟨⟨|⟨∫ω∇ u'⟩_out
    -⟨∫ω∇ u⟩_out|^2⟩_in⟩_in'
    + prob⟨ err^2⟩,


where we have set for abbreviation


    prob:=⟨ I(-1/R^d⟨δ a_h⟩_out≥λ_0 id)⟩_in'.


Since r_* has finite moments, cf (<ref>), (and thus also r_*(y) in view of the shift-covariance
(<ref>), which translates into shift covariance of ∇(ϕ,σ) and thus of r_*, 
in conjunction with the stationarity of ⟨·⟩) we have


    ⟨ err^2⟩≲((ℓ/|y|)^d(R/|y|)^d((R/|y|)^β+(1/R)^β+(1/R)^α
    +ℓ/|y|+(1/ℓ)^β))^2.


Hence provided


    |y|≫λ_0^-1/βR,  
    R≫λ_0^-1/α,  
    |y|≫λ_0^-1ℓ,  ℓ≫λ_0^-1/β,


we may absorb the second rhs term to obtain


    prob(λ_0(ℓ/|y|)^d(R/|y|)^d)^2
    ≲⟨⟨|⟨∫ω∇ u'⟩_out
    -⟨∫ω∇ u⟩_out|^2⟩_in⟩_in'.




In order to derive (<ref>) for a fixed, but sufficiently large R, from (<ref>), 
it remains to argue that there exists a λ_0>0 only depending on the ensemble such that


    prob=⟨ I(-1/R^d⟨δ a_h⟩_out≥λ_0 id)⟩_in'≥exp(-CR^d)


for some constant C=C(d), where it only matters that the rhs of (<ref>) only depends on R
and is positive for every finite R.
To this purpose, we first claim that


    1/(R-2)^d∫_B_R-2(y)(a_h-a')≥2λ_0 id  ⟹  
    -1/R^dδ a_h≥λ_0 id.


Here comes the argument: by definition (<ref>) we have


    ξ̃·δ a_hξ=∫(ξ+∇ϕ_hξ')·(a_h'-a_h)ξ̃    ϕ_hξ':=ξ_jϕ_jh'.


Indeed, since by definition (<ref>) of ϕ_hξ' and (<ref>) we have
-∇· a_h'(ξ+∇ϕ_hξ')=0, which we rewrite as
-∇· a_h'∇ϕ_hξ'=∇·(a_h'-a_h)ξ, we obtain from (<ref>) and (<ref>) the representation


    ξ̃·δ a_hξ=∫_B_R(y)ξ·(a_0-a_h)ξ̃-∫∇ϕ_hξ'· a_h'∇ϕ_hξ̃'.


Hence for the quadratic part we have the inequality


    ξ·δ a_hξ≤∫_B_R(y)ξ·(a_0-a_h)ξ,


which in view of a_0=a'_|B_R(y) implies


    1/R^d∫_B_R(y)(a_h-a')≥λ_0 id  ⟹  
    -1/R^dδ a_h≥λ_0 id.


This yields (<ref>) since 


    |1/(R-2)^d∫_B_R-2(y)(a_h-a')-1/R^d∫_B_R(y)(a_h-a')|(<ref>)≲1/R


and since we are in the regime (<ref>), which in view of α≤ 1 includes 1/R≪λ_0.



Based on (<ref>) we now argue that


    ⟨ I(-1/R^d⟨δ a_h⟩_out≥λ_0 id)⟩_in'
    ≥⟨ I(1/(R-2)^d∫_B_R-2(y)(a_h-a)≥2λ_0 id)⟩.


This follows immediately from (<ref>) and the observation that thanks to the locality assumption
(<ref>), the event
1/(R-2)^d∫_B_R-2(y)(a_h-a')≥2λ_0 id only depends on X' via its restriction to 
B_R-1(y).
Hence by the independence property of the Poisson point process, on the one hand, this event implies
-1/R^d⟨δ a_h⟩_out≥λ_0 id, 
and on the other hand, the conditional probability ⟨·⟩_in' of this event agrees with the unconditional one
⟨·⟩.



In view of (<ref>), in order to establish (<ref>), it is enough to argue that
there exists a λ_0>0 such that


    ⟨ I(1/(R-2)^d∫_B_R-2(y)(a_h-a)≥2λ_0 id)⟩≥exp(-CR^d).


Here comes the argument for (<ref>):
Letting λ_0>0 be such that a_h-A(0,∅)≥2λ_0/|B_1| id, 
cf the monotonicity assumption (<ref>), we have by the locality assumption (<ref>) and shift-invariance assumption
(<ref>) that for any realization X of the Poisson point process:


    X∩ B_R-1(y)=∅   ⟹    a_h-a≥2λ_0/|B_1| id  B_R-2(y)
       ⟹   1/(R-2)^d∫_B_R-2(y)(a_h-a)≥2λ_0 id.


Hence (<ref>) follows from the defining property of the Poisson point process:


    ⟨ I(X∩ B_R-1(y)=∅)⟩=exp(-|B_R-1(y)|).








As a consequence of the independence property of the ensemble ⟨·⟩ of the Poisson point process,
we have for any random variable, in particular F:=∫ω∇ u, that


    ⟨|F-⟨ F|B_L⟩|^2⟩≥∑_i=1^N⟨|⟨ F|B_R(y_n)⟩-⟨ F⟩|^2⟩,


provided we have for the sets


    B_L, B_R(y_1), ⋯ ,B_R(y_N)  .


Hence Theorem <ref> follows immediately from Lemma <ref>, since under the assumptions of the theorem
there exists a family {y_n}_n=1,⋯,N of points such that (<ref>) holds while


    |y_n|≤ 2L     N≳ L^d,


where R is the order-one radius given by the lemma.
More precisely, we use (<ref>) with y_n playing the role of y and which by the
first property in (<ref>) assumes the form of


    ⟨|⟨ F|B_R(y_n)⟩-⟨ F⟩|^2⟩≳((ℓ/L)^d(1/L)^d|∫ĝ|)^2.


By the second property in (<ref>) we obtain for the sum


    ∑_n=1^N⟨|⟨ F|B_R(y_n)⟩-⟨ F⟩|^2⟩≳1/L^d((ℓ/L)^d|∫ĝ|)^2,


so that (<ref>) follows via (<ref>).




Acknowledgments.
The work of JL is supported in part by the National Science Foundation under award DMS-1454939. FO learned from Jim Nolen the general strategy (<ref>) for a lower bound on the variance.


Here comes the main stochastic ingredient: For our random variable F:=∫ω∇ u we claim that


    ⟨(F-⟨ F|B_L⟩)^2⟩≥1/2∑_n=1^N⟨⟨(F'-F)^2|B_R^c(y_n)⟩'⟩.


Here F' stands short for F(a') where a' is distributed according to the conditional
expectation ⟨·|B_R^c(y_n)⟩', in particular a' coincides with a outside of B_R(y_n).
Inequality (<ref>) is a consequence of the monotonicity


    ⟨(F-⟨ F|B_L⟩)^2⟩≥⟨(F-⟨ F|(⋃_n=1^NB_R(y_n))^c⟩)^2⟩,


the sub-additivity


    ⟨(F-⟨ F|(⋃_n=1^NB_R(y_n))^c⟩)^2⟩≥∑_n=1^N⟨(F-⟨ F|B_R^c(y_n)⟩)^2⟩


of conditional variances and the formula


    ⟨(F-⟨ F|B_R^c(y_n)⟩)^2⟩
    =1/2⟨(F'-F)^2|B_R^c(y_n)⟩'⟩.


While  (<ref>) relies on the elementary formula


    ⟨(F-⟨ F⟩)^2⟩
    =1/2⟨⟨(F'-F)^2⟩'⟩


applied to the conditional measure ⟨·|B_R^c(y_n)⟩ (and then we apply ⟨·⟩
to the outcome), both (<ref>) and (<ref>) start from the decomposition of the variance


    ⟨(F-⟨ F⟩)^2⟩=⟨(F-⟨ F|D^c⟩)^2⟩+⟨(⟨ F|D^c⟩-⟨ F⟩)^2⟩.


The monotonicity follows from applying (<ref>) with ⟨·⟩ replaced by ⟨·|B_L⟩
and D=⋃_n=1^NB_R(y_n)⊂ B_L^c, throwing away the second rhs term, and then applying
⟨·⟩ to the outcome.  Missing argument for (<ref>)



plain


