Scott Armstrong
email: armstrong@ceremade.dauphine.fr

Tuomo Kuusi
email: tuomo.kuusi@aalto.fi

Jean-Christophe Mourrat
email: jean-christophe.mourrat@ens-lyon.fr

MESOSCOPIC HIGHER REGULARITY AND SUBADDITIVITY IN ELLIPTIC HOMOGENIZATION

Keywords: . 2010 Mathematics Subject Classification. 35B27, 35B45. Key words and phrases. stochastic homogenization, higher regularity, subadditivity, sublinear corrector estimate

We introduce a new method for obtaining quantitative results in stochastic homogenization for linear elliptic equations in divergence form. Unlike previous works on the topic, our method does not use concentration inequalities (such as Poincaré or logarithmic Sobolev inequalities in the probability space) and relies instead on a higher (C k , k ≥ 1) regularity theory for solutions of the heterogeneous equation, which is valid on length scales larger than a certain specified mesoscopic scale. This regularity theory, which is of independent interest, allows us to, in effect, localize the dependence of the solutions on the coefficients and thereby accelerate the rate of convergence of the expected energy of the cell problem by a bootstrap argument. The fluctuations of the energy are then tightly controlled using subadditivity. The convergence of the energy gives control of the scaling of the spatial averages of gradients and fluxes (that is, it quantifies the weak convergence of these quantities) which yields, by a new "multiscale" Poincaré inequality, quantitative estimates on the sublinearity of the corrector.

Introduction

In this paper, we introduce a new method for obtaining quantitative results in elliptic homogenization. It is based on a new regularity theory for higher derivatives of solutions which is valid on mesoscopic scales. The regularity estimates essentially localize the dependence of the solutions on the coefficients, allowing for the mixing assumptions to dictate the rate of homogenization. This idea is formalized by a novel bootstrap argument which uses the regularity estimates to accelerate the convergence of the natural subadditive and superadditive quantities associated to the variational formulation of the equation. As an application, we give explicit estimates on the sublinearity of the correctors and the weak convergence of their rescaled gradients.

1.1. Motivation and informal summary of results. We consider the linear elliptic equation (1.1) -∇ ⋅ (a(x)∇u) = 0 in bounded open subsets of R d . The coefficient a(⋅) is a random field valued in the set of d-by-d symmetric matrices with eigenvalues belonging to the interval [1, Λ] for a fixed ellipticity constant Λ ≥ 1. The law of a(⋅) is given by a probability measure P which is assumed to satisfy a quantitative mixing condition.

We are interested in obtaining quantitative information on the statistical properties of solutions of (1.1) on very large scales, a concern which lies within the realm of homogenization. The qualitative theory of homogenization for such equations was developed by [START_REF] Papanicolaou | Boundary value problems with rapidly oscillating random coefficients[END_REF][START_REF] Kozlov | Averaging of differential operators with almost periodic rapidly oscillating coefficients[END_REF][START_REF] Yurinskiȋ | Averaging of symmetric diffusion in a random medium[END_REF][START_REF] Jikov | Homogenization of differential operators and integral functionals[END_REF] and the first quantitative results are due Yurinskii [START_REF] Yurinskiȋ | Homogenization error estimates for random elliptic operators[END_REF]. Developing a quantitative theory of stochastic homogenization for linear elliptic equations has received a lot of attention since the groundbreaking work of Gloria and Otto [START_REF] Gloria | An optimal variance estimate in stochastic homogenization of discrete elliptic equations[END_REF][START_REF] Gloria | An optimal error estimate in stochastic homogenization of discrete elliptic equations[END_REF][START_REF] Gloria | Quantitative results on the corrector equation in stochastic homogenization[END_REF] and Gloria, Neukamm and Otto [START_REF] Gloria | An optimal quantitative two-scale expansion in stochastic homogenization of discrete elliptic equations[END_REF][START_REF] Gloria | Quantification of ergodicity in stochastic homogenization: optimal bounds via spectral gap on Glauber dynamics[END_REF], who proved an array of optimal estimates under certain very strong mixing assumptions on the coefficients (in particular, that P satisfies some form of a spectral gap inequality). At the core of their results are the use of sensitivity estimates and the spectral gap assumption to derive moment bounds on the gradient of the corrector and the Green's functions, which culminated in the work of Marahrens and Otto [START_REF] Marahrens | Annealed estimates on the Green's function[END_REF]. These gradient estimates, following Naddaf and Spencer [START_REF] Naddaf | Estimates on the variance of some homogenization problems[END_REF], give control of the sensitivity of the correctors themselves to changes in the coefficient field, and thereby yield optimal estimates after another application of the spectral gap inequality.

A regularity theory for stochastic homogenization was recently introduced by Armstrong and Smart [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF]. In particular, [4, Theorem 1.2] provides a gradient bound, of exactly the sort required for the quantitative theory of homogenization, which applies to arbitrary solutions of (1.1) and with much stronger (and essentially optimal) stochastic integrability. The techniques of [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF] were further developed by Armstrong and Mourrat [START_REF] Armstrong | Lipschitz regularity for elliptic equations with random coefficients[END_REF], for general (nonlinear) divergenceform equations and systems, and, in a similar spirit by Gloria, Neukamm and Otto [START_REF] Gloria | A regularity theory for random elliptic operators[END_REF] for linear equations and systems. These works were inspired by the celebrated papers of Avellaneda and Lin [START_REF] Avellaneda | Compactness methods in the theory of homogenization[END_REF][START_REF] Avellaneda | L p bounds on singular integrals in homogenization[END_REF], who proved in particular uniform Lipschitz estimates for equations with periodic coefficients. As in [START_REF] Avellaneda | Compactness methods in the theory of homogenization[END_REF][START_REF] Avellaneda | L p bounds on singular integrals in homogenization[END_REF], the philosophy is to show that solutions of the heterogeneous equation (1.1) inherit the regularity of the limiting constant-coefficient equation due to homogenization. The method of [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF] is however different from these previous works, which were based on compactness arguments and used strong bounds on the correctors. Instead, the idea is to mimic more closely the proof of the classical Schauder estimates, using quantitative methods to replace compactness: one shows that, at each scale, a solution of (1.1) may be approximated by the solution of the homogenized equation and then iterates the resulting estimate over dyadic scales. If the error in the approximation is small enough (an algebraic or Dini rate of homogenization is needed), then the argument yields a uniform Lipschitz estimate. This argument is more robust than the compactness argument of [START_REF] Avellaneda | Compactness methods in the theory of homogenization[END_REF][START_REF] Avellaneda | L p bounds on singular integrals in homogenization[END_REF], since it separates the approximation step from the iteration step. Indeed, in addition to its applicability in the stochastic setting, it has yielded new results even in the periodic and almost periodic settings (cf. [START_REF] Armstrong | Lipschitz estimates in almost-periodic homogenization[END_REF][START_REF] Shen | Boundary estimates in elliptic homogenization[END_REF]).

Another idea from [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF][START_REF] Armstrong | Lipschitz regularity for elliptic equations with random coefficients[END_REF] is that subadditive arguments are the way to prove quantitative homogenization results that are optimal in stochastic integrability (under essentially any mixing assumption). This leads to a regularity theory which is optimal in terms stochastic integrability. Subadditive arguments are so effective for this purpose because they reduce questions of stochastic integrability for random variables with complicated dependence on the random environment (such as solutions of the PDE) to estimates of a finite sum of bounded random variables which satisfy the same mixing properties as the coefficients. Under assumptions on the coefficient field which provide a spectral gap inequality, the regularity theory (specifically the Lipschitz estimate) is strong enough to quickly recover the optimal quantitative estimates on the corrector proved by Gloria and Otto. Thus an important consequence of these new regularity estimates on the quantitative theory was to separate the gradient estimates, which can now be proved under very general mixing assumptions, from the rest of the quantitative theory, which until now requires spectral gap-type assumptions.

In the present paper, we propose to take this program one step further by deriving quantitative bounds for the correctors (with explicit exponents) without using the spectral gap or other Poincaré-type concentration inequalities. This is the first step in obtaining a quantitative theory of stochastic homogenization applicable to general coefficient fields. Instead of relying on concentration of measure, we show by a bootstrap argument that the regularity theory itself can be used to improve the rate of homogenization. In addition to providing a new point of view in the quantitative theory of stochastic homogenization and allowing for more general mixing assumptions, our arguments yield estimates which are much stronger in stochastic integrability (under strong mixing assumptions, we get exponential moments rather than just pth moments). This is because, as in [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF][START_REF] Armstrong | Lipschitz regularity for elliptic equations with random coefficients[END_REF], our method allows us to use subadditivity to control the stochastic fluctuations. So far, the estimates we can obtain (see Theorem 1.2, below) are unfortunately suboptimal in their scalings compared to what can be proved under spectral gap assumptions. This is due to the presence of boundary layers encountered in the analysis which eventually force the bootstrap to halt before desired. We hope to address this issue in the near future.

In addition to the uniform Lipschitz estimate, our approach requires higher C k,1 estimates for k ≥ 1, which were left essentially implicit in [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF][START_REF] Armstrong | Lipschitz regularity for elliptic equations with random coefficients[END_REF]. Of course, unlike the Lipschitz estimate, such higher derivative estimates for k ≥ 1 cannot hold uniformly on the unit scale, but they are valid on mesoscopic scales: the assertion, which is stated precisely in Theorem 2.1 below, is roughly that a solution of (1.1) on the ball B R , with R ≫ 1, can be well-appoximated by a kth degree polynomial on any mesoscopic ball B r (x) ⊆ B R 2 with radius r ≥ R 1-ε , for a specified exponent ε > 0 which depends in particular on k. By "well-approximated", we mean that the quality of the approximation of the solution by the polynomial is, up to a constant, as good as one would have for a harmonic function with the same oscillation. This result is presented below in Theorem 2.1. Fischer and Otto [START_REF] Fischer | A higher-order large-scale regularity theory for random elliptic operators[END_REF] have very recently developed a higher regularity theory along somewhat different lines (see also Remark 2.3).

While we state and prove our results for scalar equations with symmetric coefficients under the assumption of finite range of dependence, we emphasize that these choices are not imposed on us by any limitation of our method. Indeed, the arguments work verbatim for systems with just changes to the notation, the techniques of [START_REF] Armstrong | Lipschitz regularity for elliptic equations with random coefficients[END_REF] can be used to get rid of the symmetry assumption, and it is precisely under weaker mixing conditions that our estimates would be optimal, since the weak mixing stops the bootstrap before it sees boundary layers. Our choice to make these additional assumptions rather reflects a desire to maintain the readability of the paper by presenting the ideas in their simplest setting. 1.2. Assumptions. We work in the Euclidean space R d in dimension d ≥ 2 and with a fixed ellipticity constant Λ ≥ 1. We consider the space of coefficient fields a(⋅) valued in the symmetric d-by-d matrices satisfying, for all

ξ ∈ R d , (1.2) ξ 2 ≤ ξ ⋅ a(x)ξ ≤ Λ ξ 2 .
We define Ω to be the set of all such coefficient fields:

Ω ∶= a(⋅) ∶ a ∶ R d → R d×d is Lebesgue measurable, satisfies (1.
2) and a t = a .

We endow Ω with the translation group {T y } y∈R d , which acts on Ω via

(T y a)(x) ∶= a(x + y),
and the family {F (U)} of σ-algebras on Ω, with F (U) defined for each Borel subset U ⊆ R d by F (U) ∶= σ-algebra on Ω generated by the family of maps

a ↦ U q ⋅ a(x)p ϕ(x) dx, p, q ∈ R d , ϕ ∈ C ∞ c (R d ).
Roughly, F (U) contains the information about the coefficients restricted to U.

We denote the largest of these σ-algebras by F ∶= F (R d ). The translation group may be naturally extended to F itself by defining

T y A ∶= {T y a ∶ a ∈ A} , A ∈ F
and to any random element X by setting (T z X)(a) ∶= X(T z a).

Throughout the paper, we consider a probability measure P on (Ω, F ) which is assumed to satisfy the following two conditions: (P1) P is stationary with respect to Z d -translations: for every z ∈ Z d and A ∈ F ,

P [A] = P [T z A] .
(P2) P has a unit range of dependence: for every pair of Borel subsets U, V ⊆ R d with dist(U, V ) ≥ 1, F (U) and F (V ) are P-independent.

The expectation of an F -measurable random variable X with respect to P is denoted by E [X].

1.3. Notation. We continue with some notation used throughout the paper. For a measurable set E ⊆ R d , we denote the Lebesgue measure of E by E unless E is a finite set, in which case E denotes the cardinality of E. For a bounded Lipschitz domain

U ⊆ R d with U < ∞ and p ∈ [1, ∞), we denote the normalized L p (U) norm of a function f ∈ L p (U) by (1.3) f L p (U ) ∶= ⨏ U f (x) p dx 1 p . It is also convenient to denote f L ∞ (U ) ∶= f L ∞ (U ) . For a vector-valued F ∈ L p (U; R d ), we write F L p (U ) ∶= F L p (U ) . The average of a function f ∈ L 1 (U) on U is denoted by (f ) U ∶= ⨏ U f (x) dx. For U ⊆ R d with U < ∞ and vector fields F, G ∈ L 2 (U; R d ), we denote (1.4) ⟨F, G⟩ U ∶= ⨏ U F (x) ⋅ G(x) dx.
Since we also wish to quantify the convergence of various functions in the weak L p topology, it is natural to work with W -1,p norms. For this purpose we introduce the normalized

W -1,p norm of F ∈ L p (U; R d ) by F W -1,p (U ) ∶= sup ⟨F, η⟩ U ∶ η ∈ W 1,p ′ (U; R d ), (η) U = 0, ∇η(x) L p ′ (U ) = 1 .
As usual, p ′ denotes the Hölder conjugate of an exponent p ∈ [1, ∞]. We also use the shorthand notation

(1.5) F H -1 (U ) ∶= F W -1,2 (U ) .
Note that W -1,p ′ and H -1 are (somewhat unconventionally) used to denote the duals of W 1,p R and H 1 R, rather than W 1,p 0 and H 1 0 . We denote cubes of side length R > 0 by

◻ R = ◻ R (0) ∶= - 1 2 R, 1 2 R d , ◻ R (x) ∶= x + ◻ R .
The family of cubes of side length at least one is

(1.6) C ∶= ◻ R (x) ∶ x ∈ R d , R ≥ 1 .
The set of (real-valued) polynomials on R d with degree at most k ∈ N is denoted by P k . The (random) vector space of solutions of (

1.1) in U ⊆ R d is (1.7) A(U) ∶= u ∈ H 1 loc (U) ∶ ∀v ∈ H 1 0 (U), U ∇v(x) ⋅ a(x)∇u(x) dx = 0 .
Recall that, for each p ∈ R d , the corrector Φ(⋅, p) is defined for P-almost every a ∈ Ω as the unique function (up to a constant) in H 1 loc (R d ) with a Z dstationary, mean-zero gradient and which satisfies the equation

-∇ ⋅ (a(x) (p + ∇Φ(⋅, p))) = 0 in R d .
See [START_REF] Jikov | Homogenization of differential operators and integral functionals[END_REF] for details.

1.4. Subadditive quantities and main results. We now recall some objects from [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF][START_REF] Armstrong | Lipschitz regularity for elliptic equations with random coefficients[END_REF] which play a central role in the paper. Given a bounded Lipschitz domain U ⊆ R d and p, q ∈ R d , we denote the linear function of slope p by ℓ p (x) ∶= p ⋅ x and define

ν(U, p) ∶= inf v∈ℓp+H 1 0 (U ) ⨏ U 1 2 ∇v(x) ⋅ a(x)∇v(x) dx,
which is the natural subadditive quantity representing the energy of the solution of the "cell problem," introduced (in a more general form) by Dal Maso and Modica [START_REF] Maso | Nonlinear stochastic homogenization[END_REF][START_REF] Maso | Nonlinear stochastic homogenization and ergodic theory[END_REF] in their proof of qualitative homogenization of convex integral functionals. We also define

µ(U, q) = inf u∈H 1 (U ) ⨏ U 1 2 ∇u(x) ⋅ a(x)∇u(x) -q ⋅ ∇u(x) dx,
which is the natural superadditive quantity introduced in [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF] (also in a more general form), and is dual to ν(U, p).

The quantity ν( ⋅ , p) is subadditive in the sense that if U, U 1 , . . . , U k are bounded domains satisfying

U 1 , . . . , U k are pairwise disjoint and U ∖ (U 1 ∪ ⋯ ∪ U k ) = 0, then (1.8) ν(U, p) ≤ k i=1 U i U ν(U i , p).
This is immediate from the fact that a candidate minimizer for ν(U, p) can be obtained by assembling the minimizers of the quantities ν(U i , p), as these agree on the boundaries. (We remark that we are using the term "subadditive" in a nonconventional way, due to our normalization, as it would be more standard to say that U ↦ U ν(U, p) is subadditive.) Conversely, the quantity µ( ⋅ , q) is superadditive (that is, -µ( ⋅ , q) is subadditive), since the minimizer for µ(U, q) provides with a candidate minimizer for each µ(U i , q) by restriction. Note that the minimizer for µ(U, q) will be the solution of a Neumann problem for (1.1) in U, while of course the minimizer for ν(U, p) is the solution of a Dirichlet problem. Thus we may think of ν(U, p) as "the subadditive quantity for gradients" and of µ(U, q) as "the superadditive quantity for fluxes."

The homogenized matrix a can be defined as the symmetric matrix satisfying

(1.9) 1 2 p ⋅ ap = lim R→∞ E[ν(◻ R , p)].
We view the quantities µ and ν as central to the quantitative theory of homogenization. The main step in the quantitative arguments in [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF][START_REF] Armstrong | Lipschitz regularity for elliptic equations with random coefficients[END_REF] is to show roughly that there exists α(d, Λ) > 0 such that, for every ◻ ∈ C, ν(◻, p) -1 2 p ⋅ ap ≤ C p 2 ◻ -α with overwhelming probability, and a similar estimate for µ, µ(◻, q) + 1 2 q ⋅ a -1 q ≤ C q 2 ◻ -α with overwhelming probability.

The first main result of this paper improves the rate of convergence from the unspecified and tiny α(d, Λ) > 0 in [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF] to any exponent α < 1 d .

Theorem 1.1 (Rate of convergence of subadditive quantities). Let α < 1 d .

There exists C(d, Λ, α) < ∞ such that for every p, q ∈ R d , ◻ ∈ C and λ ∈ R,

(1.10) ⎧ ⎪ ⎪ ⎨ ⎪ ⎪ ⎩ log E exp λ p -2 ◻ α ν(◻, p) -1 2 p ⋅ ap ≤ C(1 + λ 2 ),
and

log E exp λ q -2 ◻ α µ(◻, q) + 1 2 q ⋅ a -1 q ≤ C(1 + λ 2 ).
Testing the definition of µ(U, q) with the minimizer of ν(U, p) yields (1.11) J(U, p, q) ∶= ν(U, p)µ(U, q)p ⋅ q ≥ 0.

Moreover, as we will see, J(U, p, q) can be expressed in a variational form, and the optimizer is precisely the difference of the minimizers of µ(U, q) and ν(U, p).

In other words, the difference of the subadditive energies is an energy of the differences of the minimizers. The proof of Theorem 1.1 is based on a bootstrap argument, using the higher regularity estimates, applied to this difference, to accelerate the convergence of J(◻, p, ap) to zero. See Section 4.1 for a sketch of the bootstrap argument.

The estimate (1.10) is optimal in the sense that it is false for any α > 1 d . To see this, we note that minimizers for µ satisfy (oscillating) Neumann conditions while those of ν satisfy Dirichlet conditions. Therefore, the gradient of the difference of the minimizers of µ(◻, ap) and ν(◻, p) will in general be O(1) in a boundary layer of at least unit thickness. The proportion of volume of such a boundary layer relative to the whole cube ◻ is of order ◻ -1 d . Therefore J(◻, p, ap) should be at least c ◻ -1 d in general.

It is the presence of this boundary layer that so far limits our quantitative results to consequences of (1.10). On the other hand, one expects the minimizers of ν(U, p) and µ(U, ap) to be much closer in the interior of U than near the boundary, and so if the boundary layer could be neglected, there is hope to prove much more precise results.

Theorem 1.1 can be used to prove an array of quantitative estimates in homogenization with exponents which, under the strongest mixing assumptions, will typically differ from the optimal one by a square root. For example, one can show that the L 2 error in homogenization is O(ε 1 2-), where the microscopic length scale is ε, or that the H 1 norm of the two-scale expansion is O(ε 1 4-). The application we present here is an estimate of the sublinearity of the corrector and of the H -1 norm of its gradient (recall that the H -1 norm measures weak convergence in L 2 ). It roughly states that

Φ(⋅, p) -(Φ(⋅, p)) B R L 2 (B R ) + ∇Φ(⋅, p) H -1 (B R ) ≲ R 1 2 + ,
with very strong stochastic integrability. Note that estimates on the sublinearity of the corrector are intimately connected to estimates for the error in homogenization.

Theorem 1.2 (Sublinear growth of the corrector). Let β ∈ 0, 1 2 . There exists a constant C(d, Λ, β) < ∞ such that, for every p ∈ R d , R ≥ 1 and λ ∈ R,

log E exp λ p -2 R -2+2β Φ(⋅, p) -(Φ(⋅, p)) B R 2 L 2 (B R ) + ∇Φ(⋅, p) 2 H -1 (B R ) ≤ C(1 + λ 2 ).
Theorem 1.2 is new, even though it is almost certainly suboptimal and does not compare favorably with the results of Gloria and Otto [START_REF] Gloria | Quantitative results on the corrector equation in stochastic homogenization[END_REF] under spectral gap assumptions, who proved that correctors are almost bounded (see also Gloria, Neukamm and Otto [START_REF] Gloria | A regularity theory for random elliptic operators[END_REF] for more results under spectral gap assumptions). Indeed, it seems to be an open question (and there is some doubt) whether finite range of dependence implies a spectral gap-type inequality in d > 1. The best previous result for finite range of dependence was for β(d, Λ) > 0 very small (cf. [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF]). Moreover, as mentioned above, our methods are applicable under essentially any mixing condition and we believe yield essentially optimal estimates under weaker mixing conditions (such as mixing conditions so slow that spatial averages of the coefficients on scale R converge slower than R -1 to their mean). This will be explained in more details in future papers.

The novelty of Theorem 1.2 is not however in its statement, but in its proof.

1.5. Outline of the paper. Our results rely crucially on C k,1 estimates for solutions of (1.1) on mesoscopic scales. This was essentially proved in [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF], where the higher regularity theory in stochastic homogenization was introduced; the precise statements we require are proved in Section 2. In Section 3, we make the fundamental observation that the quantity J in (1.11) can be expressed as a "modulated energy" of the difference of the minimizers of µ and ν.

We break the proof of Theorem 1.1 into two main steps. First, in Section 4, we use an induction argument on the exponent α to show that, for every α < 1 d ,

p ∈ R d and cube ◻ ∈ C, (1.12) E[ν(◻, p)] - 1 2 p ⋅ ap ≤ C p 2 ◻ -α ,
and a similar estimate with µ in place of ν. This statement is a priori much weaker than the conclusion of Theorem 1.1, but we show in Section 5 that the subadditivity of -µ and ν can be used to upgrade the stochastic integrability from L 1 to exponential moments without sacrificing any of the exponent.

Finally, in Section 6, we present a functional inequality, which appears to be new and which we call multiscale Poincaré inequality. It has a wavelet flavor and converts control over the spatial averages of the gradient of a function in triadic subcubes into an estimate on the oscillation of the function itself. This, together with the regularity theory, reduces Theorem 1.2 to Theorem 1.1.

Higher regularity on mesoscopic scales

A cornerstone of the methods in this paper is a quenched C k,1 estimate for solutions of (1.1) on mesoscopic scales. This was essentially proved in [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF], where the higher regularity theory in stochastic homogenization was introduced.

The rough statement, which is given precisely in the following theorem, asserts that a solution of (1.1) can be approximated, just as well as a harmonic function, by kth degree harmonic polynomials -on length scales larger than a fixed mesoscopic scale (or the microscopic scale, if k = 0). We let A k denote the set of polynomials of degree at most k which are a-harmonic. Theorem 2.1 (C k,1 regularity on mesoscopic scales). Fix s ∈ (0, d) and k ∈ N. There exist an exponent δ(s, d, Λ) > 0, a constant C(s, k, d, Λ) < ∞ and an F -measurable random variable X ∶ Ω → [1, ∞) satisfying the estimate

(2.1) E [exp (X s )] < ∞ and such that, for every R ≥ 2X , v ∈ A(B R ) and r ∈ X ∨ R k k+δ , 1 2 R , (2.2) inf w∈A k v -w L 2 (Br) ≤ C r R k+1 v L 2 (B R ) .
Theorem 2.1 is essentially proved in [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF]. That paper only stated the (uniform) C 0,1 estimate and left (mesoscopic) higher regularity statements for k ≥ 1 implicit in the proof. To prove the proposition, we just need to expound the argument from [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF].

A second ingredient in the proof of Theorem 2.1, which also plays a key role in the bootstrap argument, is the following quenched error estimate for the Dirichlet problems proved in [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF]. It gives a quenched, deterministic estimate for the error in homogenization for the Dirichlet problem on length scales larger than a certain random scale of characteristic size O(1).

Proposition 2.2 ([4, Theorem 1.1]). Fix s ∈ (0, d), ε > 0, and let U ⊆ B 1 be a Lipschitz domain. There exists an exponent δ(d, Λ, s, ε) > 0, a constant C(s, d, Λ, U) < ∞ and an F -measurable random variable R ∶ Ω → [1, ∞)
, which depends only on (d, Λ, s, ε) and satisfies the estimate

(2.3) E [exp (R s )] < ∞, such that, for every r ≥ R, f ∈ W 1,2+ε (rU) and solutions u, u ∈ f + H 1 0 (rU) of -∇ ⋅ (a∇u) = 0 and -∇ ⋅ (a∇u) = 0 in rU,
we have the estimate

(2.4) 1 r u -u L 2 (rU ) ≤ r -δ ∇f L 2+ε (rU ) .
Remark 2.3. The mesoscopic regularity estimate in Theorem 2.1 implies Liouville theorems of all orders. These state roughly that, for each k ∈ N, the subspace of A(R d ) consisting of functions which grow at most like o ( x k+1 ) has the same dimension as A k . The latter functions correspond precisely to the kth order correctors. Note that, in view of the results of [START_REF] Armstrong | Lipschitz regularity for elliptic equations with random coefficients[END_REF], Theorem 2.1 and therefore these Liouville theorems hold in much greater generality than we present here (e.g., for nonlinear equations and with weaker mixing conditions).

A similar result was recently proved by Fischer and Otto [START_REF] Fischer | A higher-order large-scale regularity theory for random elliptic operators[END_REF], who developed a similar higher regularity theory in the general stationary ergodic setting. Similar to [START_REF] Gloria | Quantification of ergodicity in stochastic homogenization: optimal bounds via spectral gap on Glauber dynamics[END_REF], they proceed a bit differently: rather than measure the distance of an element of A(B R ) to A k , they measure its distance to the kth order correctors. We believe that while both approaches are of interest, the one presented here is more faithful to what a C k,1 estimate should be in this context (we would like, for example, to measure the C k,1 seminorm of a jth order corrector, for j ≤ k, and not get zero, otherwise we gain no information on the correctors themselves). See below Remark 4.5 in [START_REF] Armstrong | Lipschitz estimates in almost-periodic homogenization[END_REF] for some similar comments.

That a quantitative estimate for the error in homogenization for the Dirichlet problem implies higher regularity estimates via a Campanato iteration was an idea introduced in [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF]. Here we formulate a version in the following lemma, which is a variation of [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF]Lemma 5.1].

Lemma 2.4. Fix R ≥ 2, α > 0, p ∈ [1, ∞] and u ∈ L p (B R ). For each k ∈ N and s ∈ (0, R], denote D k (s) ∶= inf w∈A k u -w L p (Bs)
Assume that h ∈ 1, 1 2 R and have the property that, for every r ∈ h, 1 2 R , there exists an a-harmonic

function v ∈ C ∞ (B r ) such that (2.5) u -v L p (Br) ≤ r -α D 0 (2r).
Then, for each k ∈ N, there exists a constant

C(d, Λ, k, α) < ∞ such that, for every r ∈ h, 1 2 R , (2.6) D k (r) ≤ C r R k+1 D k (R) + Cr -α r R D 0 (R).
Proof. Fix k ∈ N. Throughout, we denote by C and c positive constants depending only on (k, α, Λ) which may vary in each occurrence. In the first few steps, we make some preliminary observations and introduce the notation we need in the main part of the argument, which begins in Step 3.

Step 1. We first observe that the triangle inequality and the hypothesis (2.5) imply that, for every r ∈ h, 1 2 R and s ∈ 0, 1 2 r , (2.7)

D k (s) ≤ C s r k+1 D k (r) + C r s d p r -α D 0 (r).
Selecting a harmonic function v ∈ C ∞ (B r ) to satisfy (2.5), we find that

D k (s) = inf w∈A k u -w L p (Bs) ≤ inf w∈A k v -w L p (Bs) + u -v L p (Bs) ≤ C s r k+1 inf w∈A k v -w L p (Br) + Cr -α D 0 (r) r s d p .
In the last line, we used the fact that any a-harmonic function v satisfies, for every 0

< s ≤ 1 2 r, inf w∈A k v -w L ∞ (Bs) ≤ C s r k+1 inf w∈A k v -w L p (Br) .
Next we use the triangle inequality and (2.5) a second time to get inf

w∈A k v -w L p (Br ) ≤ D k (r) + Cr -α D 0 (r).
Substituting into the inequality above, we get

D k (s) ≤ C s r k+1 (D k (r) + r -α D 0 (r)) + C r s d p r -α D 0 (r) ≤ C s r k+1 D k (r) + C r s d p r -α D 0 (r).
Step 2. We set up the rest of the argument. Fix θ = θ(k) ≥ c so small that Cθ = 1 4 , where C is the constant in (2.7), so that the latter inequality implies

Dk (θr) ≤ 1 2 Dk (r) + Cr -k-α D 0 (r),
where here and in what follows we set Dk (s) ∶= s -k D k (s). An iteration of the previous inequality gives

(2.8) Dk (θ m r) ≤ 2 -m Dk (r) + C m-1 j=0 2 j-m θ j r -k-α D 0 θ j r provided θ m-1 r ≥ h.
To shorten the notation, we denote r j ∶= θ j R, B j ∶= B r j and by w k,j the best kth degree polynomial approximation of u in B j , that is,

w k,j satisfies u -w k,j L p (B j ) = inf w∈P k u -w L p (B j )
Step 3. We now complete the proof of the proposition in the case k = 0. In fact, this has been already proved in [4, Lemma 5.1], but we give the argument for the sake of completeness.

First, since we may add constants to both u and v, we may assume without loss of generality that w 0,0 = 0. It follows that

w 1,0 L p (B R ) ≤ 2D 0 (R),
and hence we easily get (2.9)

w 1,0 L ∞ (B R ) + R ∇w 1,0 L ∞ (B R ) ≤ CD 0 (R).
Using (2.8) we obtain

1 r m+1 w 1,m+1 -w 1,m L p (B m+1 ) ≤ D1 (r m+1 ) + θ -1-d p D1 (r m ) ≤ 2 -m C D1 (R) + C m j=0 2 j-m r -1-α j D 0 (r j ) .
We thus deduce that

∇w 1,m+1 -∇w 1,m ≤ 2 -m C D1 (R) + C m j=0 2 j-m r -α j D 0 (r j ) r j ,
and it follows by summation that

∇w 1,n -∇w 1,0 ≤ C D1 (R) + C n j=0 r -α j D 0 (r j ) r j
provided that r n ≥ h. Combining this with (2.8) and (2.9) leads to

D 0 (r n ) r n ≤ D1 (r n ) + ∇w 1,n ≤ C D 0 (R) R + C n-1 j=0 r -α j D 0 (r j ) r j .
Taking supremum then gives sup m∈{0,...,n}

D 0 (r m ) r m ≤ C D 0 (R) R + C n j=0 r -α j sup m∈{0,...,n} D 0 (r m ) r m provided that r n ≥ h. Letting now n * = n * (d, α) be the largest integer such that C ∑ n * j=0 r -α j ≤ 1
2 and r n * ≥ h we obtain after straightforward manipulations that (2.10) sup

h≤r≤R D 0 (r) r ≤ C D 0 (R) R .
Step 4. We complete the argument in the case of general k ∈ N. Following the reasoning of the previous step, letting wk+1,j stand for the k th order polynomial part of w k+1,j , we get, by (2.8) and (2.10),

Dk (r j ) ≤ 1 r k j u -wk+1,j L p (B j ) (2.11) ≤ r j Dk+1 (r j ) + 1 r k j w k+1,j -wk+1,j L p (B j ) ≤ r j 2 -j Dk+1 (R) + Cr -k-α j D 0 (R) R + ∇ k+1 w k+1,j .
Thus we are left to estimate ∇ k+1 w k+1,j . Since

1 r k+1 j+1 w k+1,j+1 -w k+1,j+1 L p (B j+1 ) ≤ C Dk+1 (r j+1 ) + Dk+1 (r j ) ≤ C2 -j Dk+1 (R) + Cr -k-α j D 0 (R) R ,
we deduce that (2.12)

∇ k+1 w k+1,j+1 -∇ k+1 w k+1,j ≤ C2 -j Dk+1 (R) + Cr -k-α j D 0 (R) R .
Indeed, for any n ∈ N and polynomial φ ∈ P n , we have that

(2.13) ⨏ Br φ(x) p dx 1 p = ⨏ B 1 φ(rx) p dx 1 p ≥ c sup x∈B 1 φ(rx) ≥ cr n ∇ n φ for a constant c = c(d, n) > 0.
This is due to scaling and Lemma 2.6, which is stated and proved below. Now (2.12) implies after summation that

∇ k+1 w k+1,j -∇ k+1 w k+1,0 ≤ C Dk+1 (R) + Cr -k-α j D 0 (R) R .
Furthemore, since

w k,0 -w k+1,0 L p (B R ) ≤ D k+1 (R) + D k (R) ≤ 2D k (R), we obtain ∇ k+1 w k+1,0 ≤ CR -k-1 D k (R),
and consequently

∇ k+1 w k+1,j ≤ ∇ k+1 w k+1,j -∇ k+1 w k+1,0 + ∇ k+1 w k+1,0 ≤ CR -k-1 D k (R) + Cr -k-α j D 0 (R) R .
Combining this finally with (2.11) allows us to conclude with

D k (r j ) ≤ C r j R k+1 D k (R) + Cr -α j r j R D 0 (R) R ,
from which the statement (2.6) can be easily deduced.

Remark 2.5. Notice that since D k (R) ≤ D 0 (R), the estimate (2.6) implies

D k (r) ≤ C r R k+1 D 0 (R) + Cr -δ r R D 0 (R).
This is a C k,1 estimate on scales r for which the second term on the right side is smaller than the first term, that is, for r satisfying

r ≳ R k k+δ .
In the proof of Lemma 2.4, we used the following fact, which is a simple consequence of the equivalence of norms in finite-dimensional vector spaces.

Lemma 2.6. Fix k ∈ N. There exists a constant C(d, k) < ∞ such that, for every m ≤ k and w ∈ P k , ∇ m w L ∞ (B 1 ) ≤ C w L 1 (B 1 ) .
Proof of Theorem 2.1. The proof is essentially the same as -indeed, even considerably simpler than (due to homogeneity) -the argument given in the general nonlinear case considered in [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF][START_REF] Armstrong | Lipschitz regularity for elliptic equations with random coefficients[END_REF]. If we take U = B 1 and ε(d, Λ) > 0 to be the exponent in the interior Meyers estimate (cf. [2, Proposition B.6] for instance), then the solvability of the Dirichlet problem and the conclusion of Proposition 2.2, with R given there, ensures that the hypothesis of Lemma 2.4 is satisfied for h = R. The lemma then yields the result for X = R.

Properties of the modulated energy J

We begin this section with the simple but key observation that the difference between ν and µ can be expressed as a "modulated energy" of the difference of the minimizers. We define J (⋅, U, p, q) for functions w ∈ H 1 (U) by

J (w, U, p, q) ∶= ⨏ U - 1 2 ∇w(x) ⋅ a(x)∇w(x) -p ⋅ a(x)∇w(x) + q ⋅ ∇w(x) dx.
We denote the maximum of J (⋅, U, p, q) among solutions of the PDE by

J(U, p, q) ∶= max w∈A(U ) J (w, U, p, q).
To motivate the definition of J(U, p, q), we show that it is actually a familiar object: it can be decomposed into the quantities studied in [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF], which are the focus of the analysis in that paper.

Lemma 3.1. For every bounded Lipschitz domain U ⊆ R d and p, q ∈ R d , (3.1) J(U, p, q) = ν(U, p)µ(U, q)p ⋅ q.

Proof. Fix p, q ∈ R d , u ∈ A(U) and let v ∈ ℓ p + H 1 0 (U) be the minimizer in the definition of ν(U, p). Note that v ∈ A(U) and compute

ν(U, p) -⨏ U 1 2 (∇u(x) ⋅ a(x)∇u(x) -q ⋅ ∇u(x)) dx -p ⋅ q = ⨏ U 1 2 ∇v(x) ⋅ a(x)∇v(x) - 1 2 ∇u(x) ⋅ a(x)∇u(x) + q ⋅ ∇u(x) dx -p ⋅ q = ⨏ U - 1 2 (∇u(x) -∇v(x)) ⋅ a(x) (∇u(x) -∇v(x)) -p ⋅ a(x) (∇u(x) -∇v(x)) + q ⋅ (∇u(x) -∇v(x)) dx = J (u -v, U, p, q).
Here we used integration by parts twice, both taking advantage of the affine boundary condition for v to get

p = ⨏ U ∇v(x) dx
and another, which also uses that uv ∈ A(U) to get

⨏ U p ⋅ a(x) (∇u(x) -∇v(x)) dx = ⨏ U ∇v(x) ⋅ a(x) (∇u(x) -∇v(x)) dx.
We deduce that

J(U, p, q) = max w∈A(U ) J (w, U, p, q) = max u∈A(U ) J (u -v, U, p, q) = max u∈A(U ) ν(U, p) -⨏ U 1 2 (∇u(x) ⋅ a(x)∇u(x) -q ⋅ ∇u(x)) dx -p ⋅ q = ν(U, p) -µ(U, q) -p ⋅ q.
Notice that the proof of Lemma 3.1 gave more than its statement: namely, the maximizer in the definition of J(U, p, q) is precisely the difference of the minimizers of µ(U, q) and ν(U, p).

In the rest of this section, we present some basic properties of J which are needed in the bootstrap argument in the next section, and we fix a bounded Lipschitz domain U ⊆ R d throughout. We denote the unique (up to additive constants) maximizer in the definition of J(U, p, q) by u(⋅, U, p, q) ∶= maximizer of J (⋅, U, p, q) among functions in A(U).

Notice that u(⋅, U, 0, q) is the minimizer for µ(U, q) and -u(⋅, U, p, 0) is the minimizer for ν(U, p).

We next record the first and second variations of the optimization problem implicit in the definition of J(U, p, q). Henceforth, we make use of the notation (1.4). Lemma 3.2. For every p, q ∈ R d and v ∈ A(U),

(3.2)
⟨a∇u(⋅, U, p, q), ∇v⟩ U = ⟨-ap + q, ∇v⟩ U and

(3.3) J(U, p, q) -J (u(⋅, U, p, q) + v, U, p, q) = 1 2 ⟨∇v, a∇v⟩ U .
Proof. For convenience, for t ≥ 0, set u t ∶= u(⋅, U, p, q) + tv. Compute

0 ≤ J (u 0 , U, p, q) -J (u t , U, p, q) = t 2 ⨏ U 1 2 ∇v(x) ⋅ a(x)∇v(x) dx + t ⨏ U (a(x)∇u 0 (x) + a(x)p -q) ⋅ ∇v(x) dx.
Dividing by t and sending t → 0 gives

⨏ U (a(x)∇u 0 (x) + a(x)p -q) ⋅ ∇v(x) dx ≥ 0.
Repeating the argument with -v in place of v yields (3.2). Returning to the previous identity and taking t = 1 gives (3.3).

For reference, we observe that by taking v = u(⋅, U, p, q) in (3.2) we obtain the identities

J(U, p, q) = 1 2 ⟨a∇u(⋅, U, p, q), ∇u(⋅, U, p, q)⟩ U (3.4) = 1 2 ⟨-ap + q, ∇u(⋅, U, p, q)⟩ U .
It is easy to see from (3.2) that (p, q) ↦ u(⋅, U, p, q) is a linear map from R d ×R d into A(U); that is, for every p 1 , p 2 , q 1 , q 2 ∈ R d and s, t ∈ R,

(3.5) u(⋅, U, tp 1 + sp 2 , tq 1 + sq 2 ) = tu(⋅, U, p 1 , q 1 ) + su(⋅, U, p 2 , q 2 )
Likewise, J is quadratic: for every p, q ∈ R d and t > 0, (3.6) J(U, tp, tq) = t 2 J(U, p, q).

We next show that J (⋅, U, p, q) responds quadratically to perturbations near its maximum.

Lemma 3.3. For every p, q ∈ R d and v, w ∈ A(U),

(3.7) 1 4 ∇v -∇w 2 L 2 (U ) ≤ 2J(U, p, q) -J (w, U, p, q) -J (v, U, p, q) and (3.8) 2J (v, U, p, q) -J (w, U, p, q) -J(U, p, q) ≤ Λ 4 ∇v -∇w 2 L 2 (U ) . Proof. For any v 1 , v 2 ∈ A(U), 2J v 1 + v 2 2 , U, p, q -J (v 1 , U, p, q) -J (v 2 , U, p, q) = 1 4 ⟨a(∇v 1 -∇v 2 ), ∇v 1 -∇v 2 ⟩ U .
Then (3.7) follows by choosing v 1 = v and v 2 = w and using the maximality of J(U, p, q), and (3.8) follows similarly by choosing v 1 = w and v 2 = 2vw.

Lemma 3.4. For every p, q, p ′ , q ′ ∈ R d , (3.9) J(U, p + p ′ , q + q ′ ) -J(U, p, q) -J(U, p ′ , q ′ ) = ⟨-ap ′ + q ′ , ∇u(⋅, U, p, q)⟩ U .

Proof. By the first variation (3.2), ⟨a∇u(⋅, U, p, q), ∇u(⋅, U, p + p ′ , q + q ′ )⟩ U = ⟨-a(p + p ′ ) + q + q ′ , ∇u(⋅, U, p, q)⟩ U = 2J(U, p, q) + ⟨-ap ′ + q ′ , ∇u(⋅, U, p, q)⟩ U and, by (3.5), ⟨a∇u(⋅, U, p, q), ∇u(⋅, U, p

+ p ′ , q + q ′ )⟩ U = ⟨-ap + q, ∇u(⋅, U, p + p ′ , q + q ′ )⟩ U = 2J(U, p + p ′ , q + q ′ ) -⟨-ap ′ + q ′ , ∇u(⋅, U, p + p ′ , q + q ′ )⟩ U = 2J(U, p + p ′ , q + q ′ ) -2J(U, p ′ , q ′ ) -⟨-ap ′ + q ′ , ∇u(⋅, U, p, q)⟩ U .
Combining the two displays above gives (3.9).

We denote by ∇J(U, p, q) the gradient of J(U, ⋅) at (p, q). By (3.6), this is a linear mapping from R d × R d → R which, in view of the previous lemma, can be expressed by

(3.10) ∇J(U, p, q)(p ′ , q ′ ) = ⟨-ap ′ + q ′ , ∇u(⋅, U, p, q)⟩ U .
This relation between ∇J(U, p, q) and the spatial averages of the gradient and flux of u(⋅, U, p, q) will play an important role in the proof of Proposition 4.1.

For simplicity, the maps p ′ ↦ ∇J(U, p, q)(p ′ , 0) and q ′ ↦ ∇J(U, p, q)(0, q ′ ) are sometimes denoted by ∇ p J(U, p, q) and ∇ q J(U, p, q), respectively. Lemma 3.4 gives the identity

(3.11) 1 2 J(U, p 1 , q 1 ) + 1 2 J(U, p 2 , q 2 ) -J U, p 1 + p 2 2 , q 1 + q 2 2 = 1 4 J(U, p 1 -p 2 , q 1 -q 2 ) .
This readily implies the following upper convexity estimate for J in (p, q) and the lower convexity estimates for J in the variables p and q separately. Lemma 3.5. For every p, p 1 , p 2 , q, q 1 , q 2 ∈ R d ,

(3.12) 1 2 J(U, p 1 , q 1 ) + 1 2 J(U, p 2 , q 2 ) -J U, 1 2 (p 1 + p 2 ), 1 2 (q 1 + q 2 ) ≤ Λ p 1 -p 2 2 + q 1 -q 2 2 , (3.13) 1 2 J(U, p 1 , q) + 1 2 J(U, p 2 , q) -J U, 1 2 (p 1 + p 2 ), q ≥ 1 2 p 1 -p 2 2 and (3.14) 1 2 J(U, p, q 1 ) + 1 2 J(U, p, q 2 ) -J U, p, 1 2 (q 1 + q 2 ) ≥ 1 2Λ q 1 -q 2 2 .
Proof. We apply (3.11) together with (3.4) to obtain an upper bound for J(U, p 1 -p 2 , q 1 -q 2 ), and rearrange to obtain (3.12). We obtain (3.13) and (3.14) in a similar way, using the bounds

J(U, p, 0) = ν(U, p) ≥ 1 2 p 2 and J(U, 0, q) = -µ(U, q) ≥ 1 2Λ q 2 .
By subadditivity, stationarity and (1.9), for every ◻ ∈ C, we have

(3.15) 1 2 p ⋅ ap ≤ E[ν(◻, p)].
By (1.11), this implies

(3.16) 1 2 q ⋅ a -1 q ≤ -E[µ(◻, q)],
so that

(3.17) E[J(◻, p, q)] ≥ 1 2 p ⋅ ap + 1 2 q ⋅ a -1 q -p ⋅ q.
Due to Lemma 3.5, there exists a unique, deterministic matrix

Q = Q(U) such that, for every p ∈ R d , p = ∇ q E [J(U, 0, Qp)] = E [∇ q J(U, 0, Qp)] .
Indeed, for each fixed p ∈ R d , we can consider the minimum of the uniformly convex, quadratic function

q ↦ E [J(U, 0, q)] -p ⋅ q.
This defines a linear map from p to the minimum point q(p). The matrix Q is defined by q = Qp. It is easy to check that Q is symmetric, and the upper and lower uniform convexity of J(U, 0, ⋅) ensures that Q is positive and in particular invertible. In fact, we have (3.18)

I d ≤ Q ≤ ΛI d .
Notice that, for every ◻ ∈ C, we have

(3.19) Q(◻) ≤ a.
Indeed, by (3.17),

1 2 p ⋅ Q(◻)p = 1 2 Q(◻)p ⋅ E [∇ q J(U, 0, Q(◻)p)] = E [J(U, 0, Q(◻)p)] ≥ 1 2 Q(◻)p ⋅ a -1 Q(◻)p. Putting p ∶= Q(◻) -1 q gives Q(◻) -1 ≥ a -1
, which is equivalent to the claim. For future reference, we notice that since u(⋅, U, p, Qp) = u(⋅, U, p, 0) + u(⋅, U, 0, Qp) and -u(⋅, U, p, 0) is the minimizer of ν(U, p), we have that -u(⋅, U, p, 0) ∈ ℓ -p + H 1 0 (U) and ⨏ U ∇u(x, U, p, 0) dx = -p and thus by (3.10) that

∇ q E [J(U, p, Qp)] = ∇ q E [J(U, 0, Qp)] + ∇ q E [J(U, p, 0)] = p -p = 0.
Similarly, we denote by P = P (U) the deterministic, symmetric matrix defined via the relation q = ∇ p E [J(U, P q, 0)] . As above, we have the estimates

1 Λ I d ≤ P (U) ≤ I d ,
for every ◻ ∈ C,

P (◻) ≤ a -1 and (3.20)
∇ q E [J(U, P q, q)] = 0.

We next give an estimate for the difference between Q(◻) and the homogenized coefficients a in terms of the expected size of J(◻, p, Q(◻)p). Lemma 3.6. There exists C(d, Λ) < ∞ such that, for any ◻ ∈ C,

a -Q(◻) + a -1 -P (◻) ≤ C sup p∈R d p -2 E [J(◻, p, Q(◻)p)] .
Proof. We prove only the estimate for a -Q(◻) . The argument for the bound on a -1 -P (◻) is similar. We drop the dependence on ◻ and write Q = Q(◻) and J(p, q) = J(◻, p, q). Set η ∶= sup

p∈R d p -2 E [J(p, Qp)] .
By (3.17), we have, for every

p ∈ R d , 0 ≤ 1 2 p ⋅ ap + 1 2 Qp ⋅ a -1 Qp -p ⋅ Qp ≤ E [J(p, Qp)] ≤ η p 2 .
It follows by uniform convexity that, for every p ∈ R d ,

Qp -ap 2 ≤ Cη p 2
and thus, by upper uniform convexity of q ↦ E [J(p, q)] and the fact that this map achieves its minimum at q = Qp, we deduce

E [J(p, ap)] ≤ Cη p 2 .
Using (3.1), (3.15) and (3.16), we get

p ⋅ ap ≤ E [J(p, 0)] + E [J(0, ap)] ≤ p ⋅ ap + Cη p 2 .
In view of (3.15) and (3.16), this implies that, for every p, q ∈ R d ,

E [J(p, 0)] - 1 2 p ⋅ ap ≤ Cη p 2 and E [J(0, q)] - 1 2 q ⋅ a -1 q ≤ Cη q 2 .
From this, the definition of Q and the fact that J is quadratic, we get

1 2 p ⋅ Qp = 1 2 Qp ⋅ E [∇ q J(0, Qp)] = E [J(0, Qp)] ≤ 1 2 Qp ⋅ a -1 Qp + Cη p 2 . Thus Qa -1 Q ≥ Q -CηI d , that is, Q(a -1 -Q -1 )Q ≥ -CηI d . By (3.18), this implies Q ≥ a -CηI d .
In view of (3.19), the proof is now complete.

The bootstrap argument

This section is devoted to the proof of the following estimate for J.

E [J (◻, p, ap)] ≤ C p 2 ◻ -α .
An immediate consequence of Proposition 4.1 is Corollary 4.2. For every α ∈ 0, 1 d , there exist C(d, Λ, α) < ∞ such that, for every p, q ∈ R d and ◻ ∈ C,

(4.2) ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ E [ν(◻, p)] - 1 2 p ⋅ ap ≤ C p 2 ◻ -α ,
and

E [µ(◻, q)] + 1 2 q ⋅ a -1 q ≤ C q 2 ◻ -α .
Proof of Corollary 4.2. By (3.1) and Proposition 4.1,

(4.3) E [ν(◻, p)] -E [µ(◻, ap)] -p ⋅ ap ≤ C p 2 ◻ -α .
Combining this with (3.15) and (3.16) yields the result.

We begin the proof of Proposition 4.1 with some reductions. First, observe that it suffices to demonstrate (4.1) for cubes of the form ◻ R , R ≥ 1, since we can then apply the result to translations of the law P (recall that we only assume Z d -stationarity). Second, in view of Lemma 3.6, in order to prove (4.1), it suffices to show that, for each α < 1 d , there exists C(d, Λ, α) < ∞ such that, for every R ≥ 1 and p ∈ R d , (4.4)

E [J(◻ R , p, Q(◻ R )p)] ≤ C p 2 R -dα .
The proof of (4.4) is by induction. For each α ∈ (0, 1) and K ≥ 1, we let S(α, K) be the assertion that, for every R ≥ 1 and p ∈ R d ,

E [J(◻ R , p, Q(◻ R )p)] ≤ K p 2 R -dα .
In view of Lemma 3.1, the base case of our bootstrap was proved in [START_REF] Armstrong | Quantitative stochastic homogenization of convex integral functionals[END_REF]. S(α 0 , K 0 ) holds.

The previous proposition was proved in [4, Section 3] by showing that minimizers of µ(◻ n , q) for large n are expected to be flat-that is, close to a deterministic affine function-and this allows comparison to ν(◻ n , p) for an appropriate p (the slope of the affine function). Note that the result was proved with a in place of Q(◻ R ), but since Q(◻ R ) is the minimum of the map q ↦ E [J(◻ R , p, q)], we can make this replacement.

In view of Proposition 4.3, it therefore suffices to show that there exist ε(d, Λ) > 0 and C(d, Λ) < ∞ such that, for every α ∈ α 0 , 1 d and K ≥ 1, (4.6)

S(α, K) ⇒ S α + ε(1 -dα), CK 3 2
.

Indeed, an iteration of (4.6), starting from S(α 0 , K 0 ), yields (4.4). The rest of this section is focused on the proof of (4.6). Throughout, we fix α ∈ α 0 , 1 d and K ≥ 1, and assume that S(α, K) holds. We fix R ≥ 1 and p ∈ R d . By homogeneity, we may assume that p ≤ 1. We set (4.7)

q ∶= Q(◻ R )p,
where we recall that Q(◻ R ) is defined before Lemma 3.6. Throughout,

v ∶= u (⋅, ◻ R , p, q)
denotes the maximizer of J(◻ R , p, q). Notice that Lemma 3.3 gives

(4.8) ∇v 2 L 2 (◻ R ) ≤ 4J (◻ R , p, q) .
In particular, the induction hypothesis yields (4.9)

E ∇v 2 L 2 (◻ R ) ≤ 4KR -dα .
4.1. Rough sketch of the argument for (4.6). Before giving the complete details, we present an informal summary of the proof of (4.6). Choose a mesoscopic scale r ∈ (1, R) and partition the cube ◻ R into smaller cubes of the form ◻ r (y), y ∈ ◻ R . Observe that, using (3.2) in the last step,

J(◻ R , p, q) = R r -d y J (v, ◻ r (y), p, q) ≤ R r -d y ⟨-ap + q, ∇v⟩ ◻r(y) = R r -d y ⟨a∇ũ y , ∇v⟩ ◻r(y) ,
where ũy ∶= u(◻ r (y), p, q). If the induction hypothesis S(α, K) holds, then, using also (4.9) and the Lipschitz estimate, we have

(4.10) E[ ∇ũ y 2 L 2 (◻r(y)) ] ≤ Cr -dα and E[ ∇v 2 L 2 (◻r(y)) ] ≤ CR -dα .
No improvement on the estimate of E[J(◻ R , p, q)] can be obtained from these observations alone, as we have not yet used the mixing properties of the coefficients.

If v is sufficiently close (in say L 2 ) to an affine function of slope p y in the mesoscopic cube ◻ r (y), then, since ũy is a solution of (1.1), we may integrate by parts to get ⟨a∇ũ y , ∇v⟩ ◻r(y) = ⟨a∇ũ y , p y ⟩ ◻r(y) + a small error.

Recall that the dual vector q was chosen in (4.7) so that

E ⨏ ◻ R ∇v(x) dx = 0, that is, y E [p y ] = 0.
Hence, the sum The main insight is that the higher regularity theory guarantees that the slopes p y exist, that the approximation to v has a sufficiently small error, and that p y changes very slowly in y. Meanwhile, the spatial averages of the fluxes ⨏ ◻r(y) a∇ũ y are essentially independent as y varies. If we denote by s a mesoscale with r ≪ s ≪ R on which the vectors p y are approximately constant, then we can expect to bound (4.11) by (4.12) R -dα 2 r -dα 2 s r

-d 2 .
Indeed, by (4.10), R -dα 2 is the expected size of each p y and r -dα 2 is the expected size of the spatially averaged flux in a mesocube. The factor of (s r) -d 2 comes from the CLT scaling, as (s r) d is the number of smaller mesocubes of size r in each larger mesocube of size s.

The expression in (4.12) can be made smaller than R -dα by choosing the mesoscales r and s appropriately, provided that α < 1. This suggests that the correct estimate for E [J(◻ R , p, q)] should be (up to possible logarithmic factors) R -d . The argument we present below does not perform so well, and saturates at R -1 due to a boundary layer we have neglected in this sketch. 4.2. Declaration of parameters and mesoscales. We now proceed with the rigorous argument for (4.6). We take mesoscopic scales r, s, l ≥ 1 such that

r = CR 1-(m+2)ε , s = CR 1-2ε , l = CR 1-ε ,
for fixed parameters ε ∈ 0, 1 4m and m ≥ max{d, 2} to be selected below. Here the constants C in (4.14) are very close to 1 and are chosen so that the constraints (4.13) are satisfied.

The largest mesoscale l denotes the thickness of a boundary layer we remove from ◻ R in the first part of the argument; it also serves as a reference scale for the regularity estimates we will apply at the smaller scale s, where we apply Theorem 2.1 to obtain polynomial approximations to v. Finally, r denotes the size of the smallest mesoscale cells in which we compare v to local solutions (denoted by v z below) chosen to match the polynomial approximations made on the larger scale. Each of these mesoscales, even the smallest, will be chosen to be very close to the macroscale R.

The degree of the polynomials used in the mesoscopic approximation will be k ∈ N, also to be chosen. We let the exponent δ = δ(d, Λ) > 0 be the minimum of the exponent given in the statement of Theorem 2.1 for the choice s = 1 and the one given in the statement of Proposition 2.2 for the choices s = 1, U = ◻ 1 √ d and ε = 1. We may assume δ ≤ 1. We let X denote the maximum of the random variables X and R from Theorem 2.1 and Proposition 2.2, respectively, with the same choices of parameters.

In order to apply Theorem 2.1 from scale l to scale s, we must have that (4.15) s ≥ l k k+δ . For this, it suffices to impose the restriction

(4.16) ε ≤ δ 2(k + 1)
.

In the course of the proof, we will also find it convenient to fix the parameters m, k, ε as follows:

(4.17)

m ∶= 2(1 + α) 1 -α ∨ 2d(1 + α) 2 -dα ∨ d and (4.18) k ∶= ⌈ 10d(1 + θ)m θ ⌉ and ε ∶= α 2(m + 3) ∧ δ 2(k + 1) ∧ δθ 20d(1 + θ)m .
Here the positive parameter θ(d, Λ) is the exponent related to the local and global Meyers' estimates, which in particular give us that, for every ψ ∈ W 1,2(1+θ) (◻ r (y)), w ∈ A(◻ 2r (y)) and w ∈ (ψ

+ H 1 0 (◻ r (y))) ∩ A(◻ r (y)), (4.19) ∇w L 2(1+θ) (◻r(y)) ≤ C ∇w L 2 (◻ 2r (y)) and (4.20) ∇ w L 2(1+θ) (◻r(y)) ≤ C ∇ψ L 2(1+θ) (◻r(y)) .
Notice that, since α ∈ α 0 , 1 d , we have ε ≥ c(d, Λ) > 0. In fact, each of k, m, and ε are bounded above and below by constants depending only on (d, Λ). This implies that each of the constants C and c in the estimates below will depend only on (d, Λ) instead of (d, Λ, k, m, ε).

The improvement in the exponent α will be the result of a CLT scaling arising due to the fact that the v z 's inside each s sized cell are independent of each other. This is seen in the proof in Step 3, below. In the final step of the proof, we will define each of the mesoscales and other parameters in such a way that the error terms encountered in the argument will be at most CKR -dα-ε(1-dα) .

4.3.

Removal of boundary layer. The improvement in the exponent α is based on the application of regularity estimates to v. Since we do not have a boundary condition for v, these estimates are inapplicable near the boundary of the macroscopic cube ◻ R . Therefore we must remove a boundary layer, which is accomplished very simply by subadditivity. It is the error we make here that forces our bootstrap argument to halt before dα = 1.

As R l is an integer, the cubes {◻ l (y)} y∈lZ d ∩◻ R form a partition of ◻ R (up to a set of Lebesgue measure zero). We let ◻ ○ R denote the cube obtained from ◻ R after removing a mesoscopic boundary layer of thickness l:

◻ ○ R ∶= ◻ R ∖ ⋃ ◻ l (y) ∶ y ∈ lZ d , ∂◻ l ∩ ∂◻ R ≠ ∅ . Notice that ◻ R ∖ ◻ ○ R ◻ R ≤ C R l -1
.

Then by subadditivity, stationarity, Lemma 3.6 and the induction hypothesis, we obtain

E [J(◻ R , p, q)] (4.21) ≤ ◻ ○ R ◻ R E [J(◻ ○ R , p, q)] + ◻ R ∖ ◻ ○ R ◻ R E [J(◻ R ∖ ◻ ○ R , p, q)] ≤ ◻ ○ R ◻ R E [J(◻ ○ R , p, q)] + CK R l -1 l -dα .
In the second inequality in the display above, we used the fact that ◻ R ∖ ◻ ○ R is a union of cubes of the form y + ◻ l ; the induction hypothesis and Lemma 3.6 imply that (4.22)

E [ Q(◻ R ) -Q(y + ◻ l ) ] ≤ CKl -dα ,
so that another use of the induction hypothesis gives the desired estimate.

Throughout the rest of the argument we denote

Z ∶= rZ d ∩ ◻ ○ R .

4.4.

Local mesoscopic approximations of v. The next step in the argument is to introduce local approximations of v in each mesoscopic cube ◻ 3r (y), with y ∈ Z. These are also solutions of (1.1), that is, members of A(◻ 3r (y)). The advantage will be that the approximations, which we denote by v y , serve to localize the dependence on the environment and allow us to exploit the independence assumption.

We first introduce polynomial approximations of v in the larger mesoscopic cubes of the form ◻ 3s (z), z ∈ sZ d ∩ ◻ ○ R . With k ∈ N to be defined below, we select w z ∈ A k such that

(4.23) v -w z L 2 (◻ 3s (z)) = inf w∈A k v -w L 2 (◻ 3s (z)) .
That is, w z is the best approximation to v in L 2 (◻ 3s (z)) among a-harmonic polynomials of degree at most k.

We next present some basic estimates concerning the expected size of w z and the quality of the approximation to u. This is where we use the higher regularity estimates in Theorem 2.1.

(d, Λ, k) < ∞ such that, for each z ∈ sZ d ∩ ◻ ○ R , (4.24) E v -w z 2 L 2 (◻ 3s (z)) ≤ Cl 2 s l 2k+2 E ∇v 2 L 2 (◻ l (z)) + CR -2d-1 , (4.25) E ∇w z 2 L ∞ (◻ 3s (z)) ≤ CE ∇v 2 L 2 (◻ l (z)) + CR -2d-1
, and, for every m ∈ {2, . . . , k},

(4.26) E s 2(m-2) ∇ m w z 2 L ∞ (◻ 3s (z)) ≤ Cl -2 E ∇v 2 L 2 (◻ l (z)) + Cs -4 R -2d-1 .
Proof. We break the argument into three steps.

Step 1. The proof of (4.25). Since w z is an a-harmonic polynomial, we have

∇w z L ∞ (◻ 3s (z)) ≤ C s w z -(v) ◻ 3s (z) L 2 (◻ 3s (z)) (4.27) ≤ C s v -(v) ◻ 3s (z) L 2 (◻ 3s (z)) ≤ C ∇v L 2 (◻ 3s (z)) ,
where, in the above display, we used Lemma 2.6 in the first line, the optimality of w z in (4.23) tested against the constant function (v) ◻ 3s (z) and the triangle inequality to get the second line, and finally the Poincaré inequality to get the last line. Squaring and taking expectations, we get

E ∇w z 2 L ∞ (◻ 3s (z)) ≤ CE ∇v 2 L 2 (◻ 3s 1 {TzX ≤Cs} + CE ∇v 2 L 2 (◻ 3s (z)) 1 {TzX >Cs} ≤ CE ∇v 2 L 2 (◻ l (z)) + CE ∇v 2 L 2 (◻ 3s (z)) 1 {TzX >Cs} .
To bound the second term on the right side, we use (4.8) and the (deterministic) estimate J(◻ R , p, q) ≤ C to very crudely bound ∇v

2 L 2 (◻ 3s (z)) ≤ Cs -d R d ≤ CR d
and then combine this with strong integrability of X , which give the following estimate:

(4.28) E 1 {TzX >s} ∇v 2 L 2 (◻ 3s (z)) ≤ CR d+2 P [X > s] ≤ CR d+2 exp (-cs) ≤ CR -2d-1 .
In the last line we used s ≥ R 1 2 from (4.14). This completes the proof of (4.25).

Step 2. The proof of (4.24). According to Theorem 2.1 (recall that the assumption (4.16) implies (4.15)), we have

v -w z L 2 (◻ 3s (z)) 1 {TzX ≤Cs} ≤ C s l k+1 l ∇v L 2 (◻ l (z)) 1 {TzX ≤Cs} .
Recall that {T y } y∈R d is the translation group acting on the probability space Ω.

In the event that T z X is too large compared to s, we proceed differently. Applying (4.23) and the Poincaré inequality, we get

v -w z L 2 (◻ 3s (z)) 1 {TzX >Cs} ≤ v -(v) ◻ l (z) L 2 (◻ 3s (z)) 1 {TzX >Cs} ≤ C l s d 2 v -(v) ◻ l (z) L 2 (◻ l (z)) 1 {TzX >Cs} ≤ Cl l s d 2 ∇v L 2 (◻ l (z)) 1 {TzX >Cs} .
Combining the previous two displays yields

v -w z L 2 (◻ 3s (z)) ≤ C s l k+1 l ∇v L 2 (◻ l (z)) 1 {TzX ≤Cs} + Cl l s d 2 ∇v L 2 (◻ l (z)) 1 {TzX >Cs} .
Squaring and taking expectations, we obtain

(4.29) E v -w z 2 L ∞ (◻ 3s (z)) ≤ C s l 2k+2 l 2 E ∇v 2 L (◻ l (z)) + Cl 2 l s d E 1 {TzX >s} ∇v 2 L 2 (◻ l (z)) .
To estimate the second term on the right side, we proceed in the same way as the end of Step 1, above, to get

l 2 l s d E 1 {T [y] X >s} ∇v 2 L 2 (◻ l (z)) ≤ CR -2d-1 .
This completes the proof of (4.24).

Step 3. The proof of (4.26). We take ℓ z ∈ A 1 to be the best affine approximation to

v: v -ℓ z L 2 (◻ 3s (z) = inf ℓ∈A 1 v -ℓ L 2 (◻ 3s (z) .
Using the triangle inequality and (4.24) twice, once with k as above and once with k = 1, we get

E w z -ℓ z 2 L 2 (◻ 3s (z)) ≤ Cl 2 s l 4 E ∇v 2 L 2 (◻ l (z)) + CR -2d-1
. Since w zℓ z is an a-harmonic polynomial, this yields, by Lemma 2.6,

E sup m∈{2,...,k} s 2(m-2) ∇ m w z 2 L ∞ (◻ 3s (z)) ≤ CE s -4 w z -ℓ z 2 L ∞ (◻ 3s (z)) ≤ Cl -2 E ∇v 2 L 2 (◻ l (z)) + Cs -4 R -2d-1 , which is (4.26).
We next introduce the local mesoscopic approximations to v, which are based on w z . For each y ∈ Z, we denote by [y] the unique element of sZ d ∩ ◻ ○ R such that y ∈ ◻ s ([y]). For each y ∈ Z, we let v y denote the solution of the Dirichlet problem in ◻ r (y) with boundary condition w

[y] . That is, v y is the unique element of A(◻ r (y)) ∩ w [y] + H 1 0 (◻ r (y))
. Note that, for every φ ∈ A(◻ r (y)), (4.30) ⟨a ∇v y -∇w [y] , ∇φ⟩ ◻r(y) = 0.

We next give the estimate for the expected difference between the gradients of v and v y , using the previous lemma. Proof. Set η y ∶= vv y and let

(4.32) t ∶= r -r 1-ε , ε ∶= 6dmε(1 + θ) θ ,
where θ(d, Λ) is as in Meyers' estimates (4.19) and (4.20). Observe that the choice of ε in (4.18) guarantees that ε ≤ 1 2 . We first split each of the summands into two pieces as follows:

(4.33) ◻ r (y) ∇η y (x) 2 dx = ◻ t (y) ∇η y (x) 2 dx + y+◻r∖◻t ∇η y (x) 2 dx .
The first term we will estimate with the aid of the Caccioppoli estimate and Lemma 4.4, and the second one using the Meyers' estimates. Indeed, we have by the Caccioppoli estimate that 1 ◻ r ◻t(y)

∇η y (x) 2 dx ≤ C (r -t) 2 ⨏ ◻r(y) η y (x) 2 dx .
To estimate the term on the right, we have by the triangle inequality that

⨏ ◻ r (y) η y (x) 2 dx ≤ 2 v -w [y] 2 L 2 (◻r(y)) + 2 v y -w [y] 2 L 2 (◻r(y)) .
Now, for the first term we have

Z -1 y∈Z v -w [y] 2 L 2 (◻r(y)) = Z -1 z∈sZ d ∩◻ ○ R s r d v -w z 2 L 2 (◻s(z)) ≤ C s R d z∈sZ d ∩◻ ○ R v -w z 2 L 2 (◻s(z)) ,
and thus we obtain by (4.24) that Z

-1 y∈Z v -w [y] 2 L 2 (◻r(y)) ≤ C s R d z∈sZ d ∩◻ ○ R l 2 s l 2k+2 E ∇v 2 L 2 (◻ l (z)) + CR -2d-1 ≤ C l 2 s l 2k+2 E ∇v 2 L 2 (◻ R ) + R -2d-1 .
On the other hand, using Proposition 2.2 and (4.25) we get

E ⨏ ◻r(y) v y (x) -w [y] (x) 2 dx ≤ Cr 2-2δ E ∇w [y] 2 L ∞ (◻r(y)) ≤ Cr 2-2δ E ∇v 2 L 2 (◻ l (y)) + CR -2d-1 . Summing this over Z implies Z -1 y∈Z E ⨏ ◻r(y) v y (x) -w [y] (x) 2 dx ≤ Cr 2-2δ E ∇v 2 L 2 (◻ R ) + R -d-1 .
Connecting above estimates and using the fact, from (4.18), we have δ ≥ ε, we arrive at

(4.34) Z -1 y∈Z E ⨏ ◻ t (y) ∇η y (x) 2 dx ≤ C ⎛ ⎝ l r -t 2 s l 2k+2 + r 1-δ r -t 2 ⎞ ⎠ E ∇v 2 L 2 (◻ R ) + CR -d-1 .
By the choice of k and t in (4.18) and (4.32), respectively, we have

l r -t 2 s l 2k+2 = R 2(mε+ε-kε) = R -2ε k-6dm(1+θ) θ -m ≤ R -4dmε ,
and by the condition for ε in (4.18),

r 1-δ r -t 2 = R -2(1-ε(m+2)) δ-6dm(1+θ)ε θ ≤ R -4dmε .
Therefore the induction hypothesis yields by way of (4.9) that

(4.35) Z -1 y∈Z E ⨏ ◻ t (y) ∇η y (x) 2 dx ≤ CKR -dα-4dmε .
To treat the second term in (4.33), Meyers' estimates (4. [START_REF] Marahrens | Annealed estimates on the Green's function[END_REF]) and (4.20) provide us, via Hölder's inequality,

1 ◻ r y+◻r∖◻t ∇η(x) 2 dx ≤ ◻ r ∖ ◻ t ◻ r θ 1+θ ⨏ ◻r(y) ∇η(x) 2+2θ dx 1 1+θ ≤ C r -t r θ 1+θ ⨏ ◻ 3r (y) ∇v(x) + ∇w [y] (x) 2 dx .
Following Step 1 in the proof of Lemma 4.4 and applying (4.25) once more we obtain

E 1 ◻ r y+◻r∖◻t ∇η(x) 2 dx ≤ C r -t r θ 1+θ E ∇v 2 L 2 (◻ l (y)) + CR -2d-1 .
By the choice of t in (4.32) and ε in (4.18) we have

r -t r θ 1+θ = R -(1-ε(m+2)) εθ 1+θ ≤ R -4dmε .
Therefore the summation yields, as before,

(4.36) Z -1 y∈Z E 1 ◻ r y+◻r∖◻t ∇η(x) 2 dx ≤ CKR -dα-4dmε .
Combining (4.33) with (4.35) and (4.36) yields the statement of the lemma.

Next we show, using the previous lemma and the uniform convexity of J in p and q separately, that the expected difference between J (v, ◻ r (y), p, q) and J (v y , ◻ r (y), p, q) is small, after averaging over all y ∈ Z.

Lemma 4.6. Assume (4.16) and (4.18). There exists C(d, Λ, k) < ∞ such that

(4.37) Z -1 y∈Z E J (v, ◻ r (y), p, q) -J (v y , ◻ r (y), p, q) ≤ CKR -d(α+ε) .
Proof. For convenience, let ũy ∶= u(⋅, ◻ r (y), p, q) denote the minimizer in the definition of J(◻ r (y), p, q). Using (3.3) twice and summing, we find that (4.38) J (v, ◻ r (y), p, q) -J (v y , ◻ r (y), p, q) = 1 2 ⟨∇v -∇v y , a (2∇ũ y -∇v -∇v y )⟩ ◻ r (y) .

Thus (4.39) J (v, ◻ r (y), p, q) -J (v y , ◻ r (y), p, q)

≤ C ⨏ ◻r(y) ∇v(x) -∇v y (x) ( ∇ũ y (x) + ∇v(x) + ∇v y (x) ) dx.

We will estimate the term on the right. By (4.22) (with l replaced by r),

E ∇ũ y 2 L 2 (◻r(y)) ≤ CE [J(◻ r (y), p, q)] ≤ CE [J(◻ r (y), p, Q(◻ r )p)] + CKr -dα ≤ CKr -dα .
On the other hand, by the induction hypothesis and (4.25), we have

Z -1 y∈Z E ∇v y 2 L 2 (◻r(y)) ≤ C Z -1 y∈Z E ∇w [y] 2 L 2 (◻r(y))
≤ CKR -dα . We now use Hölder's inequality, Lemma 4.5 and the display above to get

Z -1 y∈Z E ⨏ ◻r(y) ∇v(x) -∇v y (x) ( ∇ũ y (x) + ∇v(x) + ∇v y (x) ) dx ≤ C E Z -1 y∈Z ∇v -∇v y 2 L 2 (◻r(y)) 1 2 CK 1 2 r -dα 2 ≤ CKR -dα 2-2dmε r -dα 2 .
Now we return to (4.39) and use the previous estimate to obtain

Z -1 y∈Z E J (v, ◻ r (y), p, q) -J (v y , ◻ r (y), p, q) ≤ CKr -dα 2 R -dα 2 R -2dmε .
By the choice of r in (4.14), we have r -dα 2 = CR -dα 2 R dα(m+2)ε 2 , and therefore we get the desired result using α ≤ 1 and m ≥ max{d, 2}.

4.5. Improving the exponent using independence. As in the proof of Lemma 4.6, we set ũy ∶= u(⋅, ◻ r (y), p, q). Applying (3.2) and (4.30), we find that, for every y ∈ Z, J (v y , ◻ r (y), p, q) ≤ ⟨-ap + q, ∇v y ⟩ ◻r(y)

(4.40) = ⟨a∇ũ y , ∇v y ⟩ ◻ r (y) = ⟨a∇ũ y , ∇w [y] ⟩ ◻r(y) .
This together with (4.21) and Lemma 4.6 imply that

(4.41) E [J (◻ R , p, q)] ≤ C R r -d E y∈Z ⟨a∇ũ y , ∇w [y] ⟩ ◻r(y) + CK R l -1 l -dα + CKR -d(α+ε) .
The sum on the right side of (4.41) allows us to see a CLT-type scaling because the terms are essentially independent for each y inside a single larger mesoscopic cube ◻ s (z). This is the mechanism which improves the exponent α. The precise statement we need is formalized in the following lemma. We remark that this is the only point in the proof that we use the choice q = Q(◻ R )p.

Lemma 4.7. Assume (4.16) and (4.18) hold. There exists

C(d, Λ, k) < ∞ such that, for every z ∈ sZ d ∩ ◻ ○ R , (4.42)
E ⎡ ⎢ ⎢ ⎢ ⎢ ⎣ sup w∈P k ⎧ ⎪ ⎪ ⎨ ⎪ ⎪ ⎩ ∇w L ∞ (◻s(z)) + l ∇ 2 w L ∞ (◻s(z)) -2 × ⎛ ⎝ y∈rZ d ∩◻s(z) ⟨a∇ũ y , ∇w⟩ ◻ r (y) ⎞ ⎠ 2 ⎫ ⎪ ⎪ ⎬ ⎪ ⎪ ⎭ ⎤ ⎥ ⎥ ⎥ ⎥ ⎦ ≤ CK s r 2d r -dα s r -d + l r -2
+ Kr -dα .

Proof. We fix z ∈ sZ d ∩ ◻ ○ R throughout the proof. We first prove a bound for each particular w ∈ P k , and then put the supremum over w inside the expectation, using linearity and that P k is a finite dimensional vector space. We begin with an estimate of the expectation of each term in the sum.

Step 1. We show that, for every w ∈ P k and y ∈ rZ d ∩ ◻ s (z),

(4.43) E ⟨a∇ũ y , ∇w⟩ ◻r(y) ≤ CKr -dα ∇w L ∞ (◻s(z)) + CK 1 2 r 1-dα 2 ∇ 2 w L ∞ (◻s(z)) .
We first observe that ⟨a∇ũ y , ∇w⟩ ◻r(y) -∇w(y) ⋅ ⨏ ◻ r (y) a(x)∇ũ y (x) dx

≤ C sup x∈◻r(y) ∇w(x) -∇w(y) ⨏ ◻ r (x) ∇ũ y (x) dx.
Thus, by (4.22) and the induction hypothesis,

E ⟨a∇ũ y , ∇w⟩ ◻ r (y) -∇w(y) ⋅ ⨏ ◻r(y) a(x)∇ũ y (x) dx (4.44) ≤ Cr ∇ 2 w L ∞ (◻r(y)) E ⨏ ◻r(x) ∇ũ y (x) 2 dx 1 2 ≤ CK 1 2 r 1-dα 2 ∇ 2 w L ∞ (◻r(y)) .
On the other hand, using Lemma 3.6 and the induction hypothesis again, we get, for every ◻ ∈ C, a -1 -P (◻) ≤ CK ◻ -α , from which we deduce, using also (4.22), that

p -P (◻ r (y))q ≤ CKr -dα . Therefore E ⨏ ◻r(y) a(x)∇ũ y (x) dx = E [∇ p J (◻ r (y), p, q)] ≤ E [∇ p J (◻ r (y), P (◻ r (y))q, q)] + CKr -dα = CKr -dα ,
as E [∇ p J (◻ r (y), P (◻ r (y))q, q)] = 0 by (3.20), and hence

(4.45) E ∇w(y) ⋅ ⨏ ◻ r (y) a(x)∇ũ y (x) dx ≤ CK ∇w(y) r -dα .
We now obtain (4.43) by combining (4.44) and (4.45) and the triangle inequality.

Step 2. We show using independence that, for every w ∈ P k , (

⎢ ⎢ ⎢ ⎢ ⎣ y∈rZ d ∩◻s(z) ⟨a∇ũ y , ∇w⟩ ◻r(y) ⎤ ⎥ ⎥ ⎥ ⎥ ⎦ ≤ CKr -dα s r d ∇w 2 L ∞ (◻s(z)) . 4.46) var ⎡
We expand the variance by writing

var ⎡ ⎢ ⎢ ⎢ ⎢ ⎣ y∈rZ d ∩◻s(z) ⟨a∇ũ y , ∇w⟩ ◻r(y) ⎤ ⎥ ⎥ ⎥ ⎥ ⎦ = y,y ′ ∈rZ d ∩◻s(z)
cov ⟨a∇ũ y , ∇w⟩ ◻r(y) , ⟨a∇ũ y ′ , ∇w⟩ ◻r(y ′ ) .

Using independence, the fact that ⟨a∇ũ y , ∇w⟩ ◻r(y) is F (◻ r (y))-measurable and each cube ◻ r+1 (y) has nonempty intersection with at most C cubes of the form ◻ r+1 (y ′ ) with y ′ ∈ rZ d ∩ ◻ s (z), we obtain

y,y ′ ∈rZ d ∩◻s(z) cov ⟨a∇ũ y , ∇w⟩ ◻ r (y) , ⟨a∇ũ y ′ , ∇w⟩ ◻ r (y ′ ) ≤ C y∈rZ d ∩◻s(z)
var ⟨a∇ũ y , ∇w⟩ ◻r(y) .

Finally, we observe that the induction hypothesis gives, for each y ∈ rZ d ∩◻ s (z),

var ⟨a∇ũ y , ∇w⟩ ◻ r (y) ≤ CE ∇ũ y 2 L 2 (◻r(y)) ∇w 2 L ∞ (◻s(z)) ≤ CKr -dα ∇w 2 L ∞ (◻s(z)) .
This completes the proof of (4.46).

Step 3. We complete the proof. Observe that combining the results of the first two steps gives, for every

w ∈ P k , (4.47) E ⎡ ⎢ ⎢ ⎢ ⎢ ⎣ ⎛ ⎝ y∈rZ d ∩◻s(z) ⟨a∇ũ y , ∇w⟩ ◻r(y) ⎞ ⎠ 2 ⎤ ⎥ ⎥ ⎥ ⎥ ⎦ ≤ CK s r 2d r -dα s r -d + l r -2 + Kr -dα × ∇w L ∞ (◻s(z)) + l ∇ 2 w L ∞ (◻s(z)) 2 .
To conclude, it remains to smuggle the supremum over w ∈ P k inside the expectation. As P k is a finite dimensional vector space with dimension depending only on (k, d), there exists an integer N(k, d) ∈ N and w 1 , . . . , w N ∈ P k such that {w 1 , . . . , w N } is a basis for P k , each w j satisfies

∇w j L ∞ (◻s(z)) + l ∇ 2 w j L ∞ (◻s(z)) = 1,
and for any w ∈ P k expressed in the form

w = c 1 w 1 + ⋯ + c N w N ,
we have, for some C(k, d) ≥ 1,

c 1 + ⋯ + c N ≤ C ∇w L ∞ (◻s(z)) + l ∇ 2 w L ∞ (◻s(z)) .
It follows that sup

w∈P k ∇w L ∞ (◻s(z)) + l ∇ 2 w L ∞ (◻s(z)) -1 y∈rZ d ∩◻s(z) ⟨a∇ũ y , ∇w⟩ ◻ r (y) ≤ C N j=1 y∈rZ d ∩◻s(z)
⟨a∇ũ y , ∇w j ⟩ ◻r(y) .

Squaring and taking the expectation of the previous inequality, using that N ≤ C and then applying (4.47), we get

E ⎡ ⎢ ⎢ ⎢ ⎢ ⎣ sup w∈P k ∇w L ∞ (◻s(z)) + l ∇ 2 w L ∞ (◻s(z)) -2 ⎛ ⎝ y∈rZ d ∩◻s(z) ⟨a∇ũ y , ∇w⟩ ◻ r (y) ⎞ ⎠ 2 ⎤ ⎥ ⎥ ⎥ ⎥ ⎦ ≤ CE ⎡ ⎢ ⎢ ⎢ ⎢ ⎣ ⎛ ⎝ N j=1 y∈rZ d ∩◻s(z) ⟨a∇ũ y , ∇w j ⟩ ◻r(y) ⎞ ⎠ 2 ⎤ ⎥ ⎥ ⎥ ⎥ ⎦ ≤ C N j=1 E ⎡ ⎢ ⎢ ⎢ ⎢ ⎣ ⎛ ⎝ y∈rZ d ∩◻s(z) ⟨a∇ũ y , ∇w j ⟩ ◻r(y) ⎞ ⎠ 2 ⎤ ⎥ ⎥ ⎥ ⎥ ⎦ ≤ CK s r 2d r -dα s r -d + l r -2 + Kr -dα .
This completes the proof.

Combining (4.40) and Lemma 4.7, we obtain y∈Z

E [J (v y , ◻ r (y), p, q)] ≤ CK 1 2 s r d r -dα 2 ⎛ ⎝ s r -d 2 + l r -1 + K 1 2 r -dα 2 ⎞ ⎠ × z∈sZ d ∩◻ ○ R E ∇w z 2 L ∞ (◻s(z)) + l 2 ∇ 2 w z 2 L ∞ (◻s(z)) 1 2
.

By Lemma 4.4 and (4.9), the sum in the line above is bounded by

C R s d K 1 2 R -dα 2 , so that y∈Z E [J (u y , ◻ r (y), p, q)] ≤ CK R r d R -dα 2 r -dα 2 ⎛ ⎝ s r -d 2 + l r -1 + K 1 2 r -dα 2 ⎞ ⎠ .
Combining this with (4.41) and using K ≥ 1, we obtain

E [J (◻ R , p, q)] ≤ CK 3 2 ⎛ ⎝ R -dα 2 r -dα 2 ⎛ ⎝ s r -d 2 + l r -1 ⎞ ⎠ + R -dα 2 r -dα + R l -1 l -dα + R -d(α+ε) ⎞ ⎠ .
Using the definitions of the mesoscales, 1) .

⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ R -dα 2 r -dα 2 s r -d 2 = CR -dα R -dε m 2 (1-α)-α R -dα 2 r -dα 2 l r -1 = CR -dα R -ε(m+1-dα(m+2) 2) R -dα 2 r -dα = CR -dα R dε(m+2)-dα 2 R l -1 l -dα = CR -dα R ε(dα-
It is easy to check that that the choices made in (4.16), (4.17) and (4.18) guarantee that the first three terms above are bounded by R -d(α+ε) . In view of these choices of the parameters, it is clear that the third term (the error in removing the boundary layer) is the limiting one, and we deduce that

E [J (◻ R , p, q)] ≤ CK R -dα R ε(dα-1) + R -d(α+ε) ≤ CKR -dα R ε(dα-1) .
Thus S(α + ε(1dα) d, CK) holds, for ε ≥ c(d, λ) and C(d, Λ) < ∞ which, after a redefinition of ε, completes the proof of the claim (4.6) and therefore the proof of Proposition 4.1.

Improvement of stochastic integrability by subadditivity

The aim of this section is to complete the proof of Theorem 1.1, by strengthening the stochastic integrability of Proposition 4.1 (or Corollary 4.2) from L 1 to exponential moments. The rough argument is as follows. We decompose ◻ R into subcubes (◻ r (y)). By subadditivity, we can bound the upper fluctuations of ν(◻ R , p) by the upper fluctuations of the average over y of ν(◻ r (y), p), up to an error controlled by the difference E[ν(◻ r , p)] -E[ν(◻ R , p)], which is small by Corollary 4.2. By independence, the average over y of ν(◻ r (y), p) is unlikely to be large. The same argument applied to -µ gives a control of the lower fluctuations of µ. By duality, we can then control upper and lower fluctuations of both quantities.

We make this idea precise in the following general statement. (Recall our slightly non-standard definition of subadditivity in (1.8), and that C is the set of cubes ◻ such that ◻ ≥ 1.) Theorem 5.1. Let δ < 1, R 0 < ∞, and let μ(⋅) ≤ ν(⋅) be respectively super-and subadditive quantities, such that for every R ≥ R 0 and x ∈ R d , μ(◻ R (x)) and ν(◻ R (x)) are F (◻ R+R δ (x))-measurable.

Let α > 0 and β ∈ (0, α) ∩ (0, 1 2]. Assume that there exists c < ∞ such that for every cube ◻ ∈ C,

(5.1) E[ν(◻)] -E[μ(◻)] ≤ c ◻ α .
Assume furthermore that there exists c ′ < ∞ such that for every cube ◻ satisfying 1 ≤ ◻ ≤ 2 d and every λ ∈ R,

(5.2) E[exp(λμ(◻))] ∨ E[exp(λν(◻))] ≤ c ′ (1 + λ 2 ).
Then there exists C = C(d, δ, R 0 , α, β, c, c ′ ) < ∞ such that for every cube ◻ ∈ C and every λ ∈ R,

(5.3) log E[exp(λ ◻ β (ν(◻) -E[ν(◻)]))] ≤ Cλ 2 ,
and the same estimate holds with ν replaced by μ.

Remark 5.2. The conclusion of Theorem 5.1 can be strengthened in several directions. First, one can replace λ 2 by λ 2 ∧ λ 1 1-β in the right side of (5.3). Moreover, we believe that (5.3) holds with α = β for α < 1 2 (and with a logarithmic correction when α = 1 2 ), and in fact, the proof given below can be adapted to show this stronger result when α is sufficiently small. On the other hand, even if μ and ν were additive, one could not improve the conclusion to β > 1 2 , by the central limit theorem.

The proof of Theorem 5.1 makes use of the following lemma.

Lemma 5.3 (parabolicity of log-Laplace).

There exists λ 0 ∈ (0, 1] such that the following holds. Let a ≥ 0 and X be a random variable such that E[exp ( X )] < ∞ and E[X] = 0. If the inequality (5.4) log E[exp (λX)] ≤ aλ 2 holds for every λ such that λ ∈ [λ 0 , 1], then it holds for every λ ∈ [-1, 1].

Proof. Since E[exp( X )] < ∞, the function Ψ ∶= λ ↦ log E[exp(λX)
] is infinitely differentiable on (-1, 1). Its value and first derivative at 0 vanish, while its second derivative at λ < 1 is bounded by E[X 2 exp(λX)]. We choose λ 0 > 0 sufficiently small that for every λ ∈ [-λ 0 , λ 0 ], ∀x ∈ R, x 2 exp (λx) ≤ exp(x) + exp(-x).

(By symmetry and monotonicity, it suffices to check the inequality for λ = λ 0 .) In particular, for every λ ∈ [-λ 0 , λ 0 ], we have

Ψ ′′ (λ) ≤ E[exp(X)] + E[exp(-X)] ≤ 2a.
The result then follows by integration.

Proof of Theorem 5.1. For any positive integer m and C ≥ 0, we denote by A m (C) the assertion that for every cube ◻ satisfying 1 ≤ ◻ ≤ 3 dm and every

λ ∈ R, (5.5
)

log E[exp(λ ◻ β (ν(◻) -E[ν(◻)]))] ≤ Cλ 2 ,
and that the same estimate holds with ν replaced by μ.

Step 1. We show that for any given m 0 , there exists C < ∞ such that A m 0 (C) holds. Every cube ◻ satisfying 1 ≤ ◻ ≤ 3 dm 0 can be decomposed into a finite number of subcubes of side length between 1 and 2. By subadditivity, the hypothesis and Hölder's inequality,

(5.6) log E[exp(λ(ν(◻) -E[ν(◻)]))] ≲ 1 + λ 2 ,
uniformly over cubes ◻ satisfying 1 ≤ ◻ ≤ 3 dm 0 and λ ∈ R. The conclusion for ν then follows by Lemma 5.3, and the reasoning for μ is identical.

Step 2. We show that there exists ε > 0 such that for every m sufficiently large and

C ≥ 1, A m-1 (C) ⇒ A m ((1 + 3 -εm )C).
Since ∏ m (1 + 3 -εm ) < ∞, this will complete the proof of the proposition.

Assuming A m-1 (C), we give ourselves a cube ◻ of side length R such that 3 m-1 < R ≤ 3 m . For notational convenience, we assume that ◻ = (0, R) d . We define a partition of ◻ into 3 d subcubes of side length L, each subcube being surrounded by a layer of smaller cubes of side length ℓ = 3 γm for some γ ∈ (δ, 1).

(We choose a triadic decomposition of ◻ for coherence with the rest of the paper, but a dyadic one would be fine too.) In order for the partition to be well-defined, we ask that L = (R -2ℓ) 3 be an integer multiple of ℓ. This requirement can easily be taken care of since ℓ ≪ R, so we will neglect it for clarity. We let (z i ) 1≤i≤3 d be the centers of the 3 d subcubes of side length L such that (0, L)

∪ (ℓ + L, ℓ + 2L) ∪ (2ℓ + 2L, R) d = 3 d ⋃ i=1 ◻ L (z i ),
and (z ′ j ) 1≤j≤N be the centers of the disjoint subcubes of side length ℓ such that

◻ ∖ 3 d ⋃ i=1 ◻ L (z i ) = N ⋃ j=1 ◻ ℓ (z ′ j ) up to a set of null measure, where N = (R d -(3L) d ) ℓ d . By subadditivity, ν(◻) ≤ 3 d i=1 ◻ L ◻ ν(◻ L (z i )) + N j=1 ◻ ℓ ◻ ν(◻ ℓ (z ′ j )).
Let r, s ∈ (1, ∞) be such that 1 r + 1 s = 1. By Hölder's inequality, for every

λ ≥ 0, log E exp λ ◻ β (ν(◻) -E[ν(◻)]) ≤ 1 r log E ⎡ ⎢ ⎢ ⎢ ⎢ ⎣ exp ⎛ ⎝ rλ ◻ β ◻ L ◻ 3 d i=1 (ν(◻ L (z i )) -E[ν(◻)]) ⎞ ⎠ ⎤ ⎥ ⎥ ⎥ ⎥ ⎦ (5.7) + 1 s log E exp sλ ◻ β ◻ ℓ ◻ N j=1 (ν(◻ ℓ (z ′ j )) -E[ν(◻)]) . (5.8)
We decompose the term in (5.7) into (5.9)

1 r log E ⎡ ⎢ ⎢ ⎢ ⎢ ⎣ exp ⎛ ⎝ rλ ◻ L ◻ 1-β 3 d i=1 (ν(◻ L (z i )) -E[ν(◻ L (z i ))]) ⎞ ⎠ ⎤ ⎥ ⎥ ⎥ ⎥ ⎦ + λ ◻ L ◻ 1-β 3 d i=1 (E[ν(◻ L (z i ))] -E[ν(◻)]),
and likewise, the term in (5.8) into

(5.10)

1 s log E exp sλ ◻ ℓ ◻ 1-β N j=1 (ν(◻ ℓ (z ′ j )) -E[ν(◻ ℓ (z ′ j ))]) + λ ◻ ℓ ◻ 1-β N j=1 (E[ν(◻ ℓ (z ′ j ))] -E[ν(◻)]).
It follows from sub-/superadditivity and (5.1) that

E[ν(◻)] ≥ lim ◻ →∞ E[ν( ◻)] = lim ◻ →∞ E[μ( ◻)] ≥ E[μ(◻ L (z i ))].
Hence, the second term in (5.9) can be estimated using (5.1):

E[ν(◻ L (z i ))] -E[ν(◻)] ≤ c ◻ L α ,
and similarly,

E[ν(◻ ℓ (z ′ j ))] -E[ν(◻)] ≤ c ◻ ℓ α .
By independence, the first term in (5.9) is equal to 1 r

3 d i=1 log E exp rλ ◻ L ◻ 1-β (ν(◻ L (z i )) -E[ν(◻ L (z i ))]) ,
which by the induction hypothesis is bounded by

(5.11) 3 d C r λr ◻ L 1-β ◻ 1-β 2 ≤ Crλ 2 ,
since β ≤ 1 2 and 3 d ◻ L ≤ ◻ . In order to estimate the first term in (5.10), we split the set {z

′ j , j ′ ≤ N} into 3 d subsets {z ′ j , j ∈ Z 1 }, . . . , {z ′ j , j ∈ Z 3 d } in such a way that if j 1 ∈ Z k 1 and j 2 ∈ Z k 2 with k 1 ≠ k 2 ,
then the cubes ◻ ℓ (z j 1 ) and ◻ ℓ (z j 2 ) are at distance at least ℓ from one another. By Hölder's inequality and independence, we get that the first term in (5.10) is bounded by

1 3 d s N j=1 log E exp 3 d sλ ◻ ℓ ◻ 1-β (ν(◻ ℓ (z ′ j )) -E[ν(◻ ℓ (z ′ j ))]) ,
which by the induction hypothesis is bounded by

(5.12) CN 3 d s 3 d sλ ◻ ℓ 1-β ◻ 1-β 2 ≤ C3 d s N ◻ ℓ ◻ λ 2 .
To sum up, we have shown that

log E exp λ ◻ β (ν(◻) -E[ν(◻)]) ≤ C r + 3 d s N ◻ ℓ ◻ λ 2 + λc ◻ β ◻ ℓ α
(recall that 3 d ◻ L + N ◻ ℓ = ◻ ), where r, s ∈ (1, ∞) such that 1 r + 1 s = 1 are arbitrary. Recall that ℓ = 3 γm with γ ∈ (δ, 1). Since β < α, we can choose γ sufficiently close to 1 that

◻ β ◻ ℓ α = 3 -mε
for some ε > 0. Moreover, there exists a constant

c d such that N ≤ c d 3 m(1-γ)(d-1) , hence log E exp λ ◻ β (ν(◻) -E[ν(◻)]) ≤ C r + sc d 3 d 3 -m(1-γ) λ 2 + λc 3 -mε .
We can now choose s = 3 m(1-γ) 2 and thus obtain that for every λ ≥ 0,

(5.13) log E exp λ ◻ β (ν(◻) -E[ν(◻)]) ≤ C (1 -3 -m(1-γ) 2 ) -1 + c d 3 d-m(1-γ) 2 λ 2 + λc 3 -mε .
The same reasoning applied to -μ shows that for every λ ≤ 0,

(5.14) log E exp λ ◻ β (μ(◻) -E[μ(◻)]) ≤ C (1 -3 -m(1-γ) 2 ) -1 + c d 3 d-m(1-γ) 2 λ 2 -λc 3 -mε .
Moreover, by (5.1), we can replace E[μ(◻)] by E[ν(◻)] in (5.14) provided we replace c by 2c in the right side of (5.14). Since ν ≥ μ, we deduce that for every

λ ∈ R, log E exp λ ◻ β (ν(◻) -E[ν(◻)]) ≤ C (1 -3 -m(1-γ) 2 ) -1 + c d 3 d-m(1-γ) 2 λ 2 + 2 λ c 3 -mε .
The result for ν follows by Lemma 5.3 and similar reasoning applies to μ.

We now complete the proof of Theorem 1.1 by combining Proposition 4.1 and Theorem 5.1.

Proof of Theorem 1.1. Let α < α ′ < 1 d, p ∈ R d be a unit vector, and q = ap.

By Proposition 4.1, there exists c < ∞ such that for every ◻ ∈ C,

(5.15) E[ν(◻, p)] -E[µ(◻, q)] -p ⋅ q ≤ c ◻ α ′ .
We apply Theorem 5.1 with ν = ν( ⋅ , p) and μ = µ( ⋅ , q) + p ⋅ q. The measurability assumption on ν and μ and the property that ν ≤ μ clearly hold. Assumption (5.2) is also satisfied, since ν and μ are bounded. We thus obtain the existence of a constant C < ∞ such that for every cube

◻ ∈ C and λ ∈ R, log E[exp(λ ◻ α (ν(◻, p) -E[ν(◻, p)]))] ≤ Cλ 2 .
Moreover, the constant C does not depend on the unit vector p (since the same is true of the constant c in (5.15)). By Corollary 4.2, we obtain, for every λ ∈ R,

log E exp λ ◻ α ν(◻, p) - 1 2 p ⋅ ap ≤ C(1 + λ 2 ).
The extension to arbitrary p ∈ R d then follows by homogeneity. The reasoning for µ is identical.

Sublinear growth of the correctors

In this section we prove Theorem 1.2. The proof naturally divides into two main steps: first we introduce a functional inequality, which appears to be new and which we term the multiscale Poincaré inequality. It converts control of spatial averages of gradients into control over the function itself. Then we show, using Theorem 1.1, that estimates of spatial averages of the gradient of the correctors can be reduced to the convergence of the subadditive energy quantities.

In this section it is convenient to work with triadic cubes, so we change the notation from the rest of the paper: for every m ∈ N, ◻ m ∶= -

1 2 3 m , 1 2 3 m d .
6.1. Multiscale Poincaré inequality. Here we present an inequality which gives an estimate of the H -1 norm of ∇u in terms of spatial averages of ∇u in cubes. This can be seen as a generalization of the usual Poincaré inequality giving the bound, for every u ∈ H 1 (◻ m ), (6.1)

⨏ ◻m u(x) -(u) ◻ m 2 dx ≤ C(d)3 2m ⨏ ◻m ∇u(x) 2 dx.
The sharpness of the scaling of the constant C3 2m in the Poincaré inequality is, of course, realized by considering an affine function. In the following proposition, we show that this scaling can be improved for functions with gradients having small spatial averages relative to their absolute size: roughly, if the gradient is canceling itself out, then the function has smaller oscillation. Proposition 6.1 (Multiscale Poincaré inequality). Fix m ∈ N and, for each

n ∈ N, n ≤ m, define Z n ∶= 3 n Z d ∩ ◻ m . There exists a constant C(d) < ∞ such that, for every u ∈ H 1 (◻ m ), (6.2) u -(u) ◻m L 2 (◻m) + ∇u H -1 (◻m) ≤ C ∇u L 2 (◻m) + C m-1 n=0 3 n Z n -1 y∈Zn (∇u) y+◻n 2 1 2
.

Proof. We first prove the estimate for ∇u H -1 (◻m) and then deduce the estimate for u -(u) ◻ m L 2 (◻m) as a simple consequence. Without loss of generality, we may suppose (u) ◻m = 0.

Step 1. The estimate for ∇u

H -1 (◻m) . Fix η ∈ H 1 (◻ m ; R d ) with ⨏ ◻m ∇η(x) 2 dx = 1.
We denote by w ∈ H 2 (◻ m ) the unique (up to additive constants) solution of the Neumann problem

-∆w = -∇ ⋅ η + b in ◻ m , ∂ ν w = 0 on ∂◻ m ,
where b ∶= ⨏ ◻m ∇ ⋅ η(x) dx is chosen to ensure solvability. Then according to [START_REF] Grisvard | Elliptic problems in nonsmooth domains[END_REF][START_REF] Adolfsson | L p -integrability of the second order derivatives for the Neumann problem in convex domains[END_REF], we have w ∈ H 2 (◻ m ) and

(6.3) ⨏ ◻m ∇ 2 w(x) 2 dx ≤ C ⨏ ◻m ∇η(x) 2 dx = C.
Testing the equation for w with u and using (u) ◻m = 0 yields (6.4)

⨏ ◻ m ∇u(x) ⋅ η(x) dx = ⨏ ◻ m ∇u(x) ⋅ ∇w(x) dx.
For every n ∈ {1, . . . , m} and z ∈ Z n , we have

z+◻n ∇u(x) ⋅ ∇w(x) -(∇w) z+◻n dx = y∈Z n-1 ∩(z+◻n) y+◻ n-1 ∇u(x) ⋅ ∇w(x) -(∇w) y+◻ n-1 dx + ◻ n-1 y∈Z n-1 ∩(z+◻n) (∇w) z+◻n -(∇w) y+◻ n-1 ⋅ (∇u) y+◻ n-1 .
By the Poincaré inequality,

y∈Z n-1 ∩(z+◻n) (∇w) z+◻n -(∇w) y+◻ n-1 2 ≤ C3 2n ⨏ z+◻n ∇ 2 w(x) 2 dx.
Therefore, after summing over z ∈ Z n and using Hölder's inequality, we get

z∈Zn z+◻n ∇u(x) ⋅ ∇w(x) -(∇w) z+◻n dx ≤ y∈Z n-1 y+◻ n-1 ∇u(x) ⋅ ∇w(x) -(∇w) y+◻ n-1 dx + C3 n(1+d 2) ◻ m ∇ 2 w(x) 2 dx 1 2 y∈Z n-1 (∇u) y+◻ n-1 2 1 2
.

Iterating this and using (6.3),

◻ m ∇u(x) ⋅ ∇w(x) dx ≤ z∈Z 0 z+◻ 0 ∇u(x) ⋅ ∇w(x) -(∇w) z+◻ 0 dx + C ◻ m 1 2 m-1 n=0 3 n(1+d 2) y∈Zn (∇u) y+◻n 2 1 2
.

By the Poincaré inequality and (6.3),

z∈Z 0 z+◻ 0 ∇w(x) -(∇w) z+◻ 0 2 dx ≤ C ◻ m ∇η(x) 2 dx = C ◻ m .
Thus using Hölder's and Young's inequalities, we obtain

◻m ∇u(x) ⋅ ∇w(x) dx ≤ C ◻ m 1 2 ◻ m ∇u(x) 2 dx 1 2 + C ◻ m 1 2 m-1 n=0 3 n(1+d 2) y∈Zn (∇u) y+◻n 2 1 2
.

Using (6.4) and rearranging the inequality yields the desired conclusion after taking the supremum over all such η.

Step 2. The estimate for u L 2 (◻m) . By the representation theorem, there exists φ ∈ L 2 (◻ m ) such that ⨏ ◻m φ(x) 2 dx = 1 and (6.5)

⨏ ◻ m u(x) 2 dx 1 2 = ⨏ ◻ m u(x)φ(x) dx = ⨏ ◻ m u(x) φ(x) -(φ) ◻m dx.
Denote by w ∈ H 2 (◻ m ) the unique (up to additive constants) solution of the Neumann problem

-∆w = φ -(φ) ◻ m in ◻ m , ∂ ν w = 0 on ∂◻ m .
We have that w ∈ H 2 (◻ m ) and The conclusion now follows from the previous step and (6.5).

6.2. Proof of Theorem 1.2. Throughout this subsection, we let v(x, U, p) ∶= -u(x, U, p, 0)p ⋅ x.

That is, v(⋅, U, p) ∈ H 1 0 (U) is the minimizer in the definition of ν(U, p) with the plane x ↦ p ⋅ x subtracted.

The multiscale Poincaré inequality motivates us to prove the sublinearity of the correctors by studying the spatial averages of its gradient in mesoscopic cubes. This is accomplished by a very simple energy comparison argument combined with the Lipschitz estimate, which reduce the needed estimates to the convergence of the subadditive quantities. The statements we need are given in the following two lemmas.

In the first lemma, we compare ∇v(⋅, ◻ m , p) in the large cube ◻ m to the gradient of the function obtained by gluing together the functions v(⋅, z+◻ m , p) on the mesoscopic subgrid {z + ◻ n ∶ z ∈ Z n } and bound the difference in terms of the cell problem energies in these cubes. This completes the proof.

We now use the Lipschitz estimate to upgrade the previous estimate. This will give us uniform local control as we blow up the macroscopic cube ◻ m to the whole space and thereby yield information on the stationary correctors. Note that the first inequality in the display above is obvious if l = k, and is a consequence of Theorem 2.1 for l ≥ k + 1. Taking the square root, summing over l = n, . . . , m -1 and using the triangle inequality yields (6.7).

We next combine the previous two inequalities with the multiscale Poincaré inequality to obtain estimates on the corrector itself in a large macroscopic cube terms of ν in triadic subcubes. .

For n ∈ N with n ≤ m, we use the fact that v(⋅, z + ◻ n , p) ∈ H 1 0 (z + ◻ n ) and use integration by parts, the Hölder inequality and Lemma 6. .

By the previous display, (6.9), the triangle inequality and the fact that the L 2 norm is stronger than the H -1 norm, we get ∇Φ(⋅, p) H -1 (◻m) ≤ ∇v(⋅, ◻ m , p) H -1 (◻m) + ∇Φ(⋅, p) -∇v(⋅, ◻ m , p) H -1 (◻m) ≤ ∇v(⋅, ◻ m , p) H -1 (◻m) + C ∇Φ(⋅, p) -∇v(⋅, ◻ m , p) . This is the desired estimate for ∇Φ(⋅, p) H -1 (◻m) . To complete the proof of (6. .

We are now ready to complete the proof of Theorem 1.2.

Proof of Theorem 1.2. By standard comparisons, it suffices to establish the result with B R replaced by ◻ m , m ∈ N. We may also assume p ≤ 1 by homogeneity. Fix β ∈ 0, 1 2 and define, for each n ∈ N,

s n ∶= ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ 3 -(m-n) -ν(◻ m , p) + Z n -1 z∈Zn ν(z + ◻ n , p) 1 2 if n ≤ m -1, ⎛ ⎝ -ν(◻ n+1 , p) + 3 -d z∈3 n Z d ∩◻ n+1 ν(z + ◻ n , p) ⎞ ⎠ 1 2
if n ≥ m.

By Lemma 6.4, the theorem is proved if we can show that

log E ⎡ ⎢ ⎢ ⎢ ⎢ ⎣ exp ⎛ ⎝ λ3 2mβ ∞ n=0 s n 2 ⎞ ⎠ ⎤ ⎥ ⎥ ⎥ ⎥ ⎦ ≤ C(1 + λ 2 ).
We let β ′ ∈ (β, 1 2), and apply Jensen's inequality with respect to the measure ∑ ∞ n=0 3 -n(β ′ -β) δ n and the convex function x ↦ exp(x 2 ) to get

exp ⎛ ⎝ λ3 2mβ ∞ n=0 s n 2 ⎞ ⎠ ≤ C ∞ n=0 3 -n(β ′ -β) exp λ3 2mβ+2n(β ′ -β) s 2 n .
We analyse first the sum over n ≥ m. By Theorem 1.1 with α = 2β ′ d,

∞ n=m 3 -n(β ′ -β) E exp λ3 2mβ+2n(β ′ -β) s 2 n ≤ ∞ n=m 3 -n(β ′ -β) exp C 1 + 3 -2β(n-m) λ 2 ≤ exp C(1 + λ 2 ) .
3 -n(β ′ -β) exp C 1 + 3 -2(1-β)(m-n) λ 2 ≤ exp C 1 + λ 2 .
This completes the proof.

Proposition 4 . 1 .

 41 For every α ∈ 0, 1 d , there exist C(d, Λ, α) < ∞ such that, for every p ∈ R d and ◻ ∈ C, (4.1)

Proposition 4 . 3 ([ 4 ,

 434 Theorem 3.1]). There exists α 0 (d, Λ) > 0 and K 0 (d, Λ) < ∞ such that(4.5)

( 4

 4 y , p y ⟩ ◻r(y) ] = R r -d y E p y ⋅ ⨏ ◻y(r)a∇ũ y (x) dx can be seen as a measure of the correlations between the (discrete) random fields ⨏ ◻ r (y) a∇ũ y and (p y ).

Lemma 4 . 4 .

 44 Assume that (4.16) holds. Then there exists C

Lemma 4 . 5 . 1 y∈ZE ∇v -∇v y 2 L 2

 45122 Assume that (4.16) and (4.18) hold. Then there exists a constant C(d, Λ, k) < ∞ such that (4.31) Z -(◻r(y)) ≤ CKR -dα-4dmε .

Lemma 4 .5 further gives Z - 1 y∈ZE ∇v -∇v y 2 L 2 1 y∈Z E ∇ũ y 2 L 2 2 L 2 2 L 2

 41221222222 (◻r(y)) ≤ CKR -dα .Combining above three displays and using the triangle inequality yield Z -(◻r(y)) + ∇v y (◻r(y)) + ∇v (◻r(y)) ≤ CKr -dα .

⨏ ◻m ∇ 2

 2 w(x) 2 dx ≤ C ⨏ ◻m φ(x) -(φ) ◻m 2 dx ≤ C.Testing the equation for w with u yields that⨏ ◻m u(x) φ(x) -(φ) ◻m dx = ⨏ ◻m ∇u(x) ⋅ ∇w(x) dx.

Lemma 6 . 2 . 2 2 L 2 2 ⟨p

 622222 For every p ∈ R d and m, n ∈ N with m ≥ n,(6.6) Z n -1 z∈Zn ∇v(⋅, ◻ m , p) -∇v(⋅, z + ◻ n , p) 2 L 2 (z+◻n) ≤ -ν(◻ m , p) + Z n -1 z∈Zn ν(z + ◻ n , p), where Z n ∶= 3 n Z d ∩ ◻ m .Proof. The lemma is a simple consequence of the uniform convexity of the energy functional and a proof can be found in [4, Lemma 2.1]. For the reader's convenience, we also provide the argument here. Define V ∈ H 1 0 (◻ m ) to be the function obtained by gluing together the functions v(⋅, z + ◻ n , p) for z ∈3 n Z d ∩ ◻ m . In other words, V ≡ v(⋅, z + ◻ n , p) in z + ◻ n , for each z ∈ 3 n Z d ∩ ◻ m . Since x ↦ p ⋅ x + v(x, ◻ m , p) is a solution of (1.1), we have ⟨a (p + ∇v(⋅, ◻ m , p)) , ∇V ⟩ ◻m = 0.Using this and a direct computation, we find that1 ∇v(⋅, ◻ m , p) -∇V , ◻ m , p) -∇V, a (∇v(⋅, ◻ m , p) -∇V )⟩ ◻m = 1 2 ⟨p + ∇V, a (p + ∇V )⟩ ◻m -1 + ∇v(x, ◻ m ,p), a (p + ∇v(x, ◻ m , p))⟩ ◻m = Z n -1 z∈Zn ν(z + ◻ n , p)ν(◻ m , p).

Lemma 6 . 3 . 2 . 2 L 2 ( 2 L 2 (

 6322222 Let X be as in the statement of Theorem 2.1. There exists C(d, Λ) < ∞ such that, for every p ∈ R d , m, n, k ∈ N with m ≥ n ≥ k and X ≤ 3 k , we have the estimate(6.7) ∇v(⋅, ◻ m , p) -∇v(⋅, ◻ n , p) L 2 (◻ k ) ◻ l+1 , p) + 3 -d z∈3 l Z d ∩◻ l+1 ν(z + ◻ l , p)Proof. Fix k ∈ N with X ≤ 3 k . By Theorem 2.1 (here we just use the Lipschitz estimate) and the previous lemma, for every l ∈ N, l ≥ k, we have∇v(⋅, ◻ l+1 , p) -∇v(⋅, ◻ l , p) ◻ k ) ≤ C ∇v(⋅, ◻ l+1 , p) -∇v(⋅, ◻ l , p) ◻ l ) ≤ C ⎛ ⎝ -ν(◻ l+1 , p) + 3 -d z∈3 l Z d ∩◻ l+1 ν(z + ◻ l , p) ⎞ ⎠ .

Lemma 6 . 4 . 1 n=0 3 n 2 .

 64132 Let X be as in the statement of Theorem 2.1. Then, for every m ∈ N with 3 m ≥ X ,Φ(⋅, p) -(Φ(⋅, p)) ◻ m L 2 (◻m) + ∇Φ(⋅, p) H -1 (◻m) (6.8) ≤ C p + C m--ν(◻ m , p) + Z n ◻ l+1 , p) + 3 -d z∈3 l Z d ∩◻ l+1 ν(z + ◻ l , p)Proof. Letting m → ∞ in (6.7) yields that, for every m ∈ N with 3 m ≥ X , (6.9) ∇Φ(⋅, p) -∇v(⋅, ◻ m , p) L 2 (◻m)≤ C ∞ l=m ⎛ ⎝ -ν(◻ l+1 , p) + 3 -d z∈3 l Z d ∩◻ l+1ν(z + ◻ l , p)

2 L 2 3 n

 223 , ◻ m , p) -∇v(⋅, z + ◻ n , p) (z+◻n) ≤ C -ν(◻ m , p) + Z n -1 z∈Zn ν(z + ◻ n , p) .Applying Proposition 6.1, we obtainv(⋅, ◻ m , p) L 2 (◻m) + ∇v(⋅, ◻ m , p) H -1 (◻m) ≤ C ∇v(⋅, ◻ m , p) L 2 (◻m) + C -ν(◻ m , p) + Z n -1 z∈Zn ν(z + ◻ n , p) 1 2

  ◻ l+1 , p) + 3 -d z∈3 l Z d ∩◻ l+1 ν(z + ◻ l , p)

2 ≤ C p + C m- 1 n=0 3 n

 213 [START_REF] Maso | Nonlinear stochastic homogenization and ergodic theory[END_REF]) and obtain the estimate for Φ(x, p) -(Φ(⋅, p)) ◻ m L 2 (◻m) , we argue similarly, but there is an extra step in which we use the Poincaré inequality:Φ(x, p) -(Φ(⋅, p)) ◻ m L 2 (◻m) ≤ 2 v(⋅, ◻ m , p) L 2 (◻m) + Φ(x, p) -v(⋅, ◻ m , p) -(Φ(⋅, p)v(⋅, ◻ m , p)) ◻m L 2 (◻m) ≤ C v(⋅, ◻ m , p) L 2 (◻m) + C3 m ⨏ ◻m ∇Φ(x, p) -∇v(x, ◻ m , p) 2 dx 1 -ν(◻ m , p) + Z n ◻ l+1 , p) + 3 -d z∈3 l Z d ∩◻ l+1 ν(z + ◻ l , p)

For the sum over- 1 n=0 3

 13 n < m, we use Jensen's inequality and Theorem 1.1 to get log E exp λ3 2nβ ′ 1 2 p ⋅ ap -Z n -1 z∈Zn ν(z + ◻ n , p) ≤ C(1 + λ 2 ). Therefore m-n(β ′ -β) E exp λ3 2mβ+2n(β ′ -β) s 2

Acknowledgments. The second author was supported by the Academy of Finland project #258000.