Green's function for elliptic systems: moment bounds

We study estimates of the Green's function in $\mathbb{R}^d$ with $d \ge 2$, for the linear second order elliptic equation in divergence form with variable uniformly elliptic coefficients. In the case $d \ge 3$, we obtain estimates on the Green's function, its gradient, and the second mixed derivatives which scale optimally in space, in terms of the"minimal radius"$r_*$ introduced in [Gloria, Neukamm, and Otto: A regularity theory for random elliptic operators; ArXiv e-prints (2014)]. As an application, our result implies optimal stochastic Gaussian bounds in the realm of homogenization of equations with random coefficient fields with finite range of dependence. In two dimensions, where in general the Green's function does not exist, we construct its gradient and show the corresponding estimates on the gradient and mixed second derivatives. Since we do not use any scalar methods in the argument, the result holds in the case of uniformly elliptic systems as well.


Introduction
This paper is a contribution to recently very active area of quantitative stochastic homogenization of second order uniformly elliptic operators, the main goal of which is to quantify how close is the large scale behavior of the heterogeneous operator (−∇ · A(x)∇) −1 to the behavior of the constant-coefficient solution operator (−∇ · A hom ∇) −1 . Here A(x) stands for the nonconstant (random) coefficient field defined on R d and A hom is called the matrix of homogenized coefficients.
As originally realized in their seminal papers by Papanicolaou and Varadhan [19] and, independently, by Kozlov [16], the central object in the homogenization of elliptic operators with random coefficients is the corrector φ ξ , defined for each direction ξ ∈ R d as a solution of the following elliptic problem − ∇ x · (A(x)∇ x (x · ξ + φ ξ (A, x))) = 0 in the whole space R d . The function φ is called corrector since it corrects the linear function x · ξ, which is clearly solution to the constant-coefficient equation, to be a solution of the equation with heterogeneous coefficients. Since φ serves as a correction of a linear function, it should naturally be smaller, i.e., sublinear. Assuming the distribution of random coefficient fields A is stationary (meaning the joint distribution of A at any two points in R d is the same) and ergodic (meaning any shift-invariant random variable is almost surely constant, a property encoding decorrelation of coefficient fields over large scales), they showed that correctors are almost surely sublinear and can be used to define the homogenized coefficient A hom e i := A(e i + ∇φ e i ) .
Since the problem is linear, it clearly suffices to study the d correctors φ i := φ e i for i = 1, . . . , d. Second, borrowing notation from the statistical physics, · stands for the ensemble average (expected value) with respect to a probability distribution on the space of coefficient fields A.
Here and also later, we will often drop the argument A in random quantities like the corrector as well as the argument x related to the spatial dependence of quantities like the coefficient field A or the corrector φ.
Both mentioned works [16,19] were purely qualitative in the sense that they showed the sublinearity of the corrector in the limit of large scales without any rate. Assuming the correlation of the coefficient fields decays with a specific rate (either encoded by some functional inequality like the Spectral Gap estimate or the Logarithmic Sobolev Inequality, or by some mixing conditions or even assuming finite range of dependence), one goal of quantitative theory is to quantify the sublinearity (smallness) of the corrector and consequences thereof.
Though the present result is purely deterministic in the sense that it translates the fact that the energy of any A-harmonic function satisfies "mean value property" from some scale on (a fact that follows from the sublinearity of the corrector, see [13]) into estimates on Green's function and its derivatives, we will first mention some recent results related to sublinearity of the corrector without going too much into details.
Here and in what follows B R stands for a ball of radius R centered at the origin and ffl denotes the average integral. In [13] they showed that for small enough δ = δ(d, λ) the sublinearity of the corrector implies the mean value property, meaning that for R ≥ r ≥ r * and for any A-harmonic function u on B R (i.e., a solution of −∇ · A∇u = 0 in B R ) one has Moreover, assuming that the ensemble on the coefficient fields satisfies a coarsened version of the Logarithmic Sobolev Inequality, they showed that the minimal radius r * has stretched exponential moments where 0 ≤ β < 1 appearing in the exponent is related to the coarsening rate in the Logarithmic Sobolev Inequality. Observe that in the case β = 0, i.e., the case when we consider the Logarithmic Sobolev Inequality without coarsening, this is the optimal Gaussian bound.
Recently, reviving the parabolic approach used in the discrete setting [14], which has a benefit of conveniently disintegrating contributions to the corrector from different scales, Gloria and Otto [15] obtained a similar results assuming the coefficient fields have finite range of dependence (it is a known fact that the assumption of the finite range of dependence does not imply any of the functional inequalities we have mentioned above). As a by-product, assuming finite range of dependence Gloria and Otto get the estimates for the minimal radius r * with optimal stochastic integrability of the form Finally, on the other side of the spectrum, Fischer and Otto [11] combined Meyer's estimate together with sensitivity analysis to show that for strongly correlated coefficient fields (more precisely, they consider coefficient fields which are 1-Lipschitz images of a stationary Gaussian field with correlations bounded by |x| −β , where 0 < β ≪ 1 is coming from Meyer's estimates) it holds For the sake of completeness, without discussing any details, let us also mention the work of Armstrong and Smart [3] which predates the previously mentioned works of Otto and coauthors (see also subsequent works [1,2] for more general results) and which contains estimates with the optimal stochastic integrability for some quantity related but different from r * .
In the present paper we will obtain deterministic estimates for the Green's function based on the minimal radii r * at different points. More precisely, fixing two points x 0 , y 0 ∈ R d , we take as the input the coefficient field A and the corresponding minimal radii r * (x 0 ), r * (y 0 ) (here and in what follows r * (x 0 ) stays for the minimal radius r * of the shifted coefficient field A(·− x 0 ) ), and produce estimates on the Green's function G and its derivatives ∇ x G, ∇ y G, ∇ x ∇ y G, averaged over small scale around the points x 0 and y 0 . This averaging is necessary since we do not assume any smoothness of the coefficient fields.
Our only goal in this paper is to obtain bounds, and not to show existence (or other properties) of the Green's function. In fact, a well known counterexample of De Giorgi [8] shows that there are uniformly elliptic coefficient fields for which the Green's function does not exist. Nevertheless, as recently shown in [7] by Conlon, Otto, and the second author, this is not a generic behavior. More precisely, in [7] they show that for any uniformly elliptic coefficient field A the Green's function G = G(A; x, y) exists at almost every point y ∈ R d , provided the dimension d ≥ 3. Therefore, in the case d ≥ 3, we will assume that the Green's function G(A; ·, y) ∈ L 1 loc (R d ) exists, at least in the almost everywhere sense (i.e., for a.e. y ∈ R d ), and focus solely on the estimates. Since in R 2 the Green's function does not have to exist, but its "gradient" can possibly exists, using a reduction from 3D (where the Green's function exists) in Section 4 we construct and estimate ∇G.
There are several works studying estimates on the Green's function in the context of uniformly elliptic equations with random coefficients. Using De Giorgi-Nash-Moser approach for a parabolic equation (which is naturally restricted to the scalar case), Delmotte and Deuschel [9] obtained annealed estimates on the first and second gradient of the Green's function, in L 2 and L 1 in probability respectively, under mere assumption of stationarity of the ensemble (see also [17] for a different approach). Using different methods, Conlon, Otto, and the second author [7] recently obtained similar estimates, together with other properties of the Green's function, but without the restriction to the scalar case. For a single equation and in the discrete case, assuming that the spatial correlation of the coefficient fields decays sufficiently fast to the effect that the Logarithmic Sobolev Inequality is satisfied, Marahrens and Otto [18] upgraded the Delmotte-Deuschel bounds to any stochastic moments. Recently, for a single equation this work was extended by Gloria and Marahrens [12] into the continuum setting.
Before we state the main result, let us mention other works relating the smallness of the corrector and the properties of solutions to the heterogeneous equation. Together with Otto [6], we compare the finite energy solution u of − ∇ · A∇u = ∇ · g, with g ∈ L 2 (R d ; R d ) being supported in a unit ball around the origin, with twice corrected solution u hom of the homogenized equation − ∇ · A hom ∇u hom = ∇ ·g.
Here by twice corrected we mean that first the right-hand side g from the heterogeneous equation is replaced byg = g(Id + ∇φ) in the constant-coefficient equation, and second, we compare u with (Id + φ i ∂ i )v at the level of gradients. Using duality argument together with a compactness lemma, this gives an estimate of the difference between ∇ x ∇ y G(x, y) and ∂ i ∂ j G hom (e i + ∇φ i (x)) ⊗ (e j + ∇φ j (y)) (averaged over small balls both in x and y). In order to get such estimates, it is not enough to assume that the corrector is at most linear with small slope (as in (1)), but rather we need to assume that for some β ∈ (0, 1), it grows in the L 2 -sense at most like |x| 1−β : Hence, in comparison with the present work, in [6] we get a stronger statement (since we estimate the difference between the heterogeneous Green's function and corrected constantcoefficient Green's function while in the present paper we only control the heterogeneous Green's function alone), at the expense of stronger assumption on the smallness of the corrector and more involved proof. More precisely, here we show that the second mixed derivative of the Green's function ∇ x ∇ y G behaves like C|x − y| −d (clearly this estimate is sharp in scaling since it agrees with the behavior of the constant-coefficient Green's function), while in [6] we show that the homogenization error, i.e., the difference between ∇ x ∇ y G and twice corrected mixed second derivative of the constant-coefficient Green's function, is estimated by C|x − y| −(d+β)that means we gain a factor of |x − y| −β , where the exponent β ∈ (0, 1) is the one appearing in (2). Last, let us mention the work of Otto and the authors [5], where we push farther the results of [6] using higher order correctors. The second and higher order correctors were introduced into the stochastic homogenization setup by Fischer and Otto [10], in order to extend the C 1,α regularity estimates on large scales [13] to C 2,α estimates and C k,α estimates respectively. In [5], under the assumption of smallness of the corrector, we obtain two results about A-harmonic function in exterior domains: first, for any integer k we construct a finite dimensional space of functions such that the distance between any A-harmonic function in the exterior domain and this space is bounded by C|x − y| −(d+k) (this statement can be seen as an analogue of Liouville statements for finite energy solutions in the exterior domain). Second, assuming smallness of both the first and second order augmented correctors (i.e., also including the second order vector potential, which we had to introduce), compared with [6] we improve by 1 the exponent in the estimate between the solution of the heterogeneous equation in the exterior domain and of some corrected solution of the constant-coefficient equation.
The paper is organized as follows: In the next section we will state our assumptions together with the main result, Theorem 1, and its corollaries, Corollary 1, Corollary 2, and Corollary 3. In Section 3 we prove Theorem 1 and in Section 4 we give the argument for Corollary 1, which is the only corollary which does not immediately follow from the theorem.
Notation. Throughout the article, we denote by C a positive generic constant which is allowed to depend on the dimension d and the ellipticity contrast λ, and which may be different from line to line of the same estimate. By we will mean ≤ C. Finally, the integrals without specified domain of integration are meant as integrals over the whole space R d .

The main result
We fix a coefficient field A ∈ L ∞ (R d ; R d×d ), which we assume to be uniformly elliptic in the sense thatˆR where λ ∈ (0, 1) is fixed throughout the paper. Then we have the following result: Let G = G(A; x, y) be the Green's function defined through assuming it exists for a.e. y ∈ R d . Then we havê where r ′ * (y) = r * (A t , y) denotes the minimal radius for the adjoint coefficient field A t at a point y.
Though the Green's function does not have to exist in 2D, with the help of the Green's function in 3D we can at least define and estimate "its gradient & second mixed derivatives": Let A be a uniformly elliptic coefficient field on R 2 in the sense of (3), such that for its extension into R 3 of the form there exists two pointsX,Ȳ ∈ R 3 so that the minimal radii forĀ t andĀ at those points are finite, respectively (i.e., r * (Ā t ,X) < ∞, r * (Ā,Ȳ ) < ∞). Then for a.e. y ∈ R 2 there exists a function on R 2 , which we denote ∇G(A; ·, y), so that it satisfies in a weak sense − ∇ x · A∇G(A; ·, y) = δ(· − y).
Moreover, given x 0 , y 0 ∈ R 2 with |x 0 − y 0 | ≥ 10, we have estimates for ∇G as well as for ∇ y ∇G:ˆB Remark 1. Assuming that the coefficient field A in the statement of Corollary 1 is chosen at random with respect to a stationary and ergodic ensemble, by the standard ergodic argument (see, e.g., [13]), applied in 3D to the ensemble obtained as a push-forward of the 2D ensemble through (9), the assumption on the finiteness of the minimal radii is almost surely satisfied (say withX =Ȳ = (0, 0, 0)).

Remark 2.
It is clear from the proof of Theorem 1 that all the above estimates, i.e. (5)- (11), are true also if the domains of integration B 1 (x 0 ) and B 1 (y 0 ) are replaced by larger balls with the corresponding radii r * . Moreover, the radii of the balls could be even larger than the minimal radii (as long as these new radii are not larger than one third of a distance between centers of those balls), in which case we need to replace the minimal radii on the right-hand sides of those estimates with the actual radii of the balls.
Remark 3. The appearance of different minimal radii in (10) and (11) (in (10) the minimal radii are related to the equation in 2D, while in (11) they are the minimal radii for the equation in 3D) is not a typo. The reason is that while (10) is proved directly in 2D, the proof of (11) passes through 3D -hence the need to consider the minimal radii in 3D. In view of the relation r * (A, x) ≤ r * (Ā,x), which easily follows from the fact that any A-harmonic function in B R ⊂ R 2 can be trivially extended to anĀ-harmonic function inB R ⊂ R 3 , the estimate (11) seems to be less optimal.
For notational convenience we state the result for a single equation. Since in the proof of Theorem 1 we do not use any scalar methods (like for example De Giorgi-Nash-Moser iteration), the result holds also in the case of elliptic systems -for that one just considers that u has values in some finite-dimensional Hilbert space. Naturally, in that case all the constants will depend on the dimension of this Hilbert space.
Using the Gaussian bounds on r * for the case of coefficient fields with finite range of dependence, which were obtained recently in [15], Theorem 1 implies the following bounds: Corollary 2. Suppose · is an ensemble of λ-uniformly elliptic coefficient fields which is stationary and of unity range of dependence, and let d ≥ 2. Then there exist C(d, λ) such that for every two points x 0 , y 0 ∈ R d , |x 0 − y 0 | ≥ 10, and every ǫ > 0 we have In the case of coefficient fields with stronger correlations we can use the result from [13]: Corollary 3. Suppose d ≥ 2, and that the ensemble · is stationary and satisfies a logarithmic Sobolev inequality of the following type: There exists a partition {D} of R d not too coarse in the sense that for some 0 ≤ β < 1 it satisfies Moreover, let us assume that there is 0 < ρ ≤ 1 such that for all random variables F where the carré-du-champ of the Malliavin derivative is defined as Then there exists a constant 0 < C < ∞, depending only on d, λ, ρ, β, such that

Proof of Theorem 1
The proof is inspired by a duality argument of Avellaneda and Lin [4, Theorem 13], which they used to obtain Green's function estimates in the periodic homogenization. After stating and proving two auxiliary lemmas, we first prove the estimate on the second mixed derivative (5). Then, (6) will follow from (5) using Poincaré inequality and one additional estimate. Next we observe that (7) can be obtained from (6) by replacing the role of x and y, which can be done by considering the adjoint A t instead of A. Finally, (8) will follow from (7) in a similar way as (6) follows from (5). We thus start with the following two auxiliary lemmas. The first one is very standard: Lemma 1 (Caccioppoli inequality). Let ρ > 0, δ > 0, and let u be a solution of a uniformly elliptic equation −∇ · A∇u = 0 in B (1+δ)ρ . Then for any c ∈ R.
Lemma 2. Let R 0 ≥ r * (0), and let u be an A-harmonic function in B R 0 . Then we have Proof. Throughout the proof we write r * instead of r * (0). We assume that 2r * < R 0 , since otherwise (13) is trivial. For r ∈ [r * , R 0 ] we denote u r := ffl Br u. We have Hence, to prove (4) it is enough to show To prove it, we use the following estimate which in fact holds for any function u ∈ W 1,2 (B 2r ). We first argue how to obtain (14) the proof thanks to estimate (15): Let n ≥ 0 be the largest integer that satisfies 2 n r * ≤ R 0 /2; using (15) multiple times we get |u r * − u 2 n r * | ≤ n−1 k=0 u 2 k r * − u 2 k+1 r * Using Jensen's inequality and the fact that R 0 ≤ 2 n+2 r * we get Combination of the two previous estimates then gives (14). It remains to prove (15). Using Jensen's and Poincaré's inequalities we get

Proof of (5)
We denote R 0 := |x 0 − y 0 |/3. We split the proof of (5) into 4 steps. In the first step we show thatˆB for any ρ ∈ [1, 2] and any functional F ρ on L 2 (B ρ (x 0 )) which satisfies for any v ∈ W 1,2 (B ρ (x 0 )). In the second step, using Neumann eigenfunctions on a ball, we define a family of functionals F k satisfying (17), which will play the role of Fourier coefficients. Using these F k we then estimate´|v| 2 with a sum of N terms |F k (v)| 2 plus a residuum in the form 1 λ N´| ∇v| 2 . Here, λ N denotes the N -th Neumann eigenvalue of Laplacian on a ball. Using (17) together with the second step we get an estimate on´|∇ x ∇ y G| 2 in terms of a (good) term and a small prefactor times´|∇ x ∇ y G| 2 , integrated over a slightly larger ball. In the last step, we apply iteratively this estimate.
Let f ∈ L 2 (B R 0 (y 0 ); R d ), and let u be the finite energy solution of Then on the one hand, the Green's function representation formula yields If r * (x 0 ) ≤ R 0 (w. l. o. g. we assume r * (x 0 ) ≥ ρ), we use this in (4) to get If r * (x 0 ) ≥ R 0 , we simply have Since F ρ is linear, using (19) we have where the dot means that F ρ acts on the first variable. The previous relations then give Using definition of the norm L 2 (B R 0 (y 0 )) by duality we get Let r ′ * (y 0 ) play the same role as r * (x 0 ) but for the adjoint equation. In the case r ′ * (y 0 ) ≤ R 0 , since y → F ρ (∇ x G(·, y)) solves the adjoint equation −∇ · A t ∇F ρ (∇ x G(·, y)) = 0 in B R 0 (y 0 ), an analogue of (4) implieŝ If r ′ * (y 0 ) ≥ R 0 , we get the same conclusion for free: Step 2. Let ρ ∈ [1, 2) and δ > 0 be fixed such that (1 + δ)ρ ≤ 2. For n ≥ 0, let F n denote the functional on L 2 (B (1+δ)ρ ), defined as an inner product of the n-th eigenfunction of Neumann Laplacian with v, and let λ n be the associated eigenvalue. By writing any v ∈ W 1,2 (B (1+δ)ρ ) in terms of Neumann eigenfunctions (which can be done since they form an orthonormal basis in L 2 ) we getˆB Step 3. Combination of Step 1 and Step 2 (applied to ∇ y G(·, y)) and use of (21) yieldŝ Step 4. For a given sequence δ k > 0 such that ρΠ ∞ k=1 (1 + δ k ) ≤ 2 we consider the following iteration procedure. Let ρ 0 := 1, and for k ≥ 1 set ρ k := (1 + δ k )ρ k−1 . We denote For any N k ≥ 1, estimate (23) in Step 3 yields where the values of δ k and N k are at our disposal. We choose δ k := (2k) −2 and N k := αk 2d 2 d . Since Π ∞ k=1 (1 + 1 4k 2 ) ∼ 1.46 ≤ 2, for this choice of δ k for all k ≥ 1 we have ρ k ∈ [1,2]. Using lower bound on the Neumann eigenvalues for the ball in the form λ k ≥ Ck 2 d (in the case of a cube one can use trigonometric functions to explicitly write down the formula for eigenfunctions and eigenvalues; for the ball, one uses the monotonicity of the eigenvalues with respect to the domain, which follows from the variational formulation of the eigenvalues), we can find large enough α such that the prefactor in front of M k+1 above satisfies

For this choice (24) turns into
Iterating this we get Assuming we have sup k M k < ∞, we send K → ∞ to get Since the sum on the right-hand side is summable, we get that M 1 1. It remains to justify the assumption sup k M k < ∞. For any Λ ≥ 1, let χ Λ (y) be the characteristic function of the set y ∈ B 1 (y 0 ) :´B 2 (x 0 ) |∇ x ∇ y G(x, y)| 2 ≤ Λ . Using the previous arguments, applied to ∇ x ∇ y G(x, y)χ Λ (y), we get that where the right-hand side does not depend on Λ. Now we send Λ → ∞, and get by the Monotone Convergence Theorem. This completes the proof of (5).

Proof of (6)
We first observe that using Poincaré's inequality we can control the difference between ∇ y G and its averages over B 1 (x 0 ) by the L 2 -norm of ∇ x ∇ y G, which we already control by (5). Hence, to obtain (6) it is enough to estimate averages ffl B 1 (x 0 ) ∇ y G(x, y) dx. Such estimate will follow from an analogue of (16) applied to one particular functional F . Since in this setting we need to work with´|u| 2 and not with previously used´|∇u| 2 , we will need to use Lemma 2.
Step 1. By Poincaré inequality in the x-variable we havê By the triangle inequality we havê and so (6) follows from (25) provided we shoŵ (26) Step 2. Proof of (26). Similarly as for (5), consider arbitrary f ∈ L 2 (B R 0 (y 0 ); R d ) and the finite energy solution u of − ∇ · A∇u = −∇ · f in R d , which satisfies the energy estimatê Let F be a linear functional on If otherwise r * (x 0 ) > R 0 , we do not need anymore to appeal to Lemma 2 and may directly bound and proceed as in the previous inequality. As before, we use linearity of F and write F (∇ y G(·, y))f (y) dy .
Since f ∈ L 2 (B R 0 (y 0 ); R d ) was arbitrary, combination of the two previous estimates yieldŝ As before, it remains to argue that by going from´B R 0 (y 0 ) to´B 1 (y 0 ) we gain a factor R −d 0 . We define v(y) := F (G(·, y)), and observe that −∇A t ∇v = 0 in B R 0 (y 0 ), where A t denotes the adjoint coefficient field. Then by definition of v estimate (28) implieŝ

Proof of (8)
Similarly to the proof of (6), we use Poincaré's inequality (Step 1) to show that (8) follows from (7) provided we control averages of G (Step 2).
Step 1. By Poincaré's inequality in the x-variable we havê Then by the triangle inequality we havê G(x, y) dx 2 dy, and so (8) follows provided we shoŵ (30) Step 2. Proof of (30). Similarly as for (6), consider arbitrary f ∈ L 2 (B R 0 (y 0 )), but this time u being a finite energy solution of − ∇ · A∇u = f in R d . In order to get the energy estimate, we test the equation with u to obtain: We point out that compared to the proof of (5) or (6), we got additional R 2 0 due to the right-hand side of the equation being f and not ∇ · f .
Let F be a linear functional on If otherwise r * (x 0 ) > R 0 , then we directly bound and proceed analogously to the other case. Using the Green's function representation formula we have u(x) =´B R 0 (y 0 ) G(x, y)f (y) dy, and thus the linearity of F yields F (G(·, y))f (y) dy .
Since f ∈ L 2 (B R 0 (y 0 )) was arbitrary, we may combine the two previous estimates and concludê As before, it remains to argue that by going from´B R 0 (y 0 ) to´B 1 (y 0 ) we gain a factor R −d 0 . We define v(y) := F (G(·, y)), and observe that −∇A t ∇v = 0 in B R 0 (y 0 ). Now we use Lemma 2 with v to get (33) (33) is exactly (30).

Proof of Corollary 1
We provide a generalization of (6)- (7) in the two-dimensional case. When d = 2, the Green's function for the whole space R 2 does not have to exist; nevertheless, we may give a definition for ∇G via the Green's function on R 3 . To this purpose we introduce the following notation: If x ∈ R 3 , we writex = (x, x 3 ) ∈ R 2 × R and denote byB r ⊂ R 3 and B r ⊂ R 2 the balls of radius r and centered at the origin. For a given bounded and uniformly elliptic coefficient field A in R 2 , recall that its trivial extensionĀ to R 3 was defined in (9)  It will become clear below that the argument for the representation formula for ∇G through ∇ xḠ calls for the notion of pointwise existence inȳ ∈ R 3 of the Green's functionḠ(Ā; ·,ȳ). As mentioned in Section 1, in the case of systems we may only rely on a definition of the Green's function for almost every singularity pointȳ. Therefore, differently from the previous sections, we need to bear in mind this weaker notion of existence ofḠ.
Step 1. We argue that for almost every y ∈ R 2 the function ∇G (since G does not exist, ∇G should be understood as a symbol for a function and not as a gradient of some function G), defined through i.e., in a weak sense it solves −∇ x · A∇G(A; ·, y) = δ(· − y).
We now choose a sequence {ζ n } n∈N of test functionsζ n = η n ζ, with ζ = ζ(x) ∈ C ∞ 0 (R 2 ) and η n = η n (x 3 ) smooth cut-off function for {|x 3 | < n} in {|x 3 | < n + 1}: From the previous identity and definition (9) it followŝρ We now want to send n → +∞ in the previous identity : By our assumptions onρ andζ n , if we show thatˆs then by the Dominated Convergence Theorem we may conclude thatρ and thus (35) by the arbitrariness of the test functionρ and the separability of C ∞ 0 (R 2 ). To argue inequality (36) we proceed as follows: We define a finite radius M such that M ≥ max(r * (Ā t ,X), r * (A,Ȳ )) and supp(ρ) ⊂B M (Ȳ ), supp(ζ) ⊂ B M/2 (X), and observe that inequality (36) is implied bŷB SinceĀ is translational invariant in x 3 , the minimal radius r * (Ā t , ·) is independent of x 3 . Then, by the definition of M and Remark 2 we havêB We now cover the cylinder B M/2 (X) × R with countably many balls of radius M centered at the points (X, ±M n) ∈ R 3 . By translational invariance we can w. l. o. g. assume that Y 3 = 0. We thus bound the integral in (37) bŷB We claim that ∇xḠ(Ā; ·, ·) ∈ L 1 loc (R 3 ×R 3 ), and so the first integral on the r. h. s. of the previous identity is finite.
Here we only sketch the idea why∇xḠ ∈ L 1 loc (R 3 × R 3 ); for the proof with all the details we refer to the proof of [7,Theorem 1]. To show that∇xḠ ∈ L 1 loc it suffices to show thatB R (0)´B R (0) |∇xḠ| < ∞. In order to do that we observe that for given two distinct points x,ȳ ∈ R 3 , the proof of Theorem 1 (without the use of r * to go to smaller scales; see also Remark 2) implies in 3D where r = |x −ȳ|/3, which by Hölder's inequality turns intôB Using a simple covering argument, the above estimate holds also in the case when the balls are replaced by cubes. SinceB R (0) ×B R (0) can be written as a null-set plus a countable union of pairs of open cubesQ rn (x n ) ×Q rn (ȳ n ), each with size r n := |x n −ȳ n |/3 and such that each pair of points (x,ȳ) ∈B R (0) ×B R (0) belongs to at most one such pair of cubes, we concludêB where we used that for (x,ȳ) ∈Q rn (x n ) ×Q rn (ȳ n ) we have |x −ȳ| ∼ r n . Going back to the second term on the right-hand side of (39), an application of Hölder's inequality in both variablesx andȳ yields for the the sum over n We have established (36).
We now appeal to Lebesgue's theorem and conclude (38).
Concerning (10), there are two possible ways how to proceed. For the first we observe that (35) implies for every test function φ ∈ C ∞ c (R 2 ) where f ∈ L 2 (R 2 ; R 2 ) and u is a solution of −∇ · A∇u = −∇ · f . Therefore we have that ∇u(x) =ˆ∇ y ∇G(x, y) · f (y) dy, and the proof of (5) applies verbatim. A different way would be to mimic the argument for (11), i.e., to define ∇ y ∇G as an integral of the second mixed derivative of the Green's function in three dimension. Unfortunatelly, this way we would obtain the estimate where the minimal radii in 2D appearing on the right-hand side of (10) would need to be replaced (with possibly larger) minimal radii for 3D.