A higher-order large-scale regularity theory for random elliptic operators

We develop a large-scale regularity theory of higher order for divergence-form elliptic equations with heterogeneous coefficient fields $a$ in the context of stochastic homogenization. The large-scale regularity of $a$-harmonic functions is encoded by Liouville principles: The space of $a$-harmonic functions that grow at most like a polynomial of degree $k$ has the same dimension as in the constant-coefficient case. This result can be seen as the qualitative side of a large-scale $C^{k,\alpha}$-regularity theory, which in the present work is developed in the form of a corresponding $C^{k,\alpha}$-"excess decay"estimate: For a given $a$-harmonic function $u$ on a ball $B_R$, its energy distance on some ball $B_r$ to the above space of $a$-harmonic functions that grow at most like a polynomial of degree $k$ has the natural decay in the radius $r$ above some minimal radius $r_0$. Though motivated by stochastic homogenization, the contribution of this paper is of purely deterministic nature: We work under the assumption that for the given realization $a$ of the coefficient field, the couple $(\phi,\sigma)$ of scalar and vector potentials of the harmonic coordinates, where $\phi$ is the usual corrector, grows sublinearly in a mildly quantified way. We then construct"$k$th-order correctors"and thereby the space of $a$-harmonic functions that grow at most like a polynomial of degree $k$, establish the above excess decay and then the corresponding Liouville principle.


Introduction
We are interested in the regularity of harmonic functions u associated with a uniformly elliptic coe cient eld a in d space dimensions (by which we understand a tensor eld satisfying λ|ξ | 2 ≤ ξ · aξ and |aξ | ≤ |ξ | for some λ > 0 and any ξ ∈ R d ) via the divergence-form equation Without continuity assumptions, the local regularity of (weak nite-energy) solutions can be rather low, in particular in case of systems (see e.g., [18,Example 3] for the scalar case and [9, Section 9.1.1] for De Giorgi's celebrated counterexample in the systems case). Because of their homogeneity, the same examples show that even when the coe cients are uniformly locally smooth, the large-scale behavior of a-harmonic functions can be very di erent from the constant coe cient, i.e, Euclidean case; see e.g., Proposition 21 in the appendix below. Largescale regularity is most compactly encoded in a Liouville statement of the following form: The where ε r := sup This sublinear growth (2) was shown to hold under the assumptions of stationarity and qualitative ergodicity. In a second step, large-scale C 1,α -inner regularity estimates for a-harmonic functions were obtained, where the random constant satis es a stretched exponential bound under mild decay assumptions on the spatial covariance of a. In a later version of [11], the optimal stochastic moments for the random constant were obtained.
In the context of nonlinear elliptic systems in divergence form, the result of Armstrong and Smart [3] on the large-scale C 0,1 -estimate was generalized by Armstrong and Mourrat [2] to nonsymmetric coe cients and well beyond nite range, further con rming that the random large-scale regularity theory holds under just a mild quanti cation of ergodicity, like expressed by standard mixing conditions.
In the present work, we go beyond C 1,α and establish a large-scale C k,α -theory in the form of a corresponding excess decay and Liouville result, see Theorem 3 and Corollary 4. This li s the result of Avellaneda and Lin [5] from the periodic to the random case. To streamline presentation, we rst establish the C 2,α -versions of our theorems, see Theorem 7 and Corollary 8.
Let us clearly state that the contribution of this paper is exclusively on the deterministic side. The large-scale regularity is obtained under the assumption that the given realization a of the coe cient eld is such that the corresponding corrector couple (φ, σ ) satis es the following slight quanti cation of (2), namely, lim r→∞ ε 2,r = 0 (4) with ε 2,r := ∞ m=0 min{1, 2 m+1 /r}ε 2 m .
Note that (4) is equivalent to ∞ m=0 ε 2 m < ∞. In a recent preprint by the authors of the present paper [8], it is shown that (4) holds for almost every realization a in case of a stationary ensemble of coe cient elds under mild quanti cation of ergodicity in the form of an assumption on a mild decay of correlations of a: More precisely, given a stationary centered tensor-valued Gaussian random eldã on R d and a bounded Lipschitz map : R d×d → R d×d taking values in the set of λ-uniformly elliptic tensors, the coe cient eld a := (ã) almost surely admits correctors with the property (4) assuming just decay of correlations in the sense | ã(x)ã(y) | ≤ C|x − y| −β for some C > 0 and some β ∈ (0, c(d, λ)) (where · denotes the expectation). Note that under the assumption of a spectral gap for the ensemble, as far as the corrector φ is concerned (but not the "vector potential" σ ), an estimate like (4) could also be deduced to hold almost surely from [12,Proposition 2], modulo the passage from a discrete to a continuum medium.
The key building block for this large-scale C k,α -theory is the space of a-harmonic functions that grow at most like a polynomial of degree k at in nity. Proposition 2 and Corollary 4 imply that under our assumption (4) this space has the same dimension as in the Euclidean casee.g., for k = 2 the space of a-harmonic functions that grow at most quadratically is spanned by 1+d+ d(d+1) 2 −1 maps -, which partially answers the question by Benjamini et al. [6,Chapter 6]. The kth-order excess (11), by the decay of which we encode the C k,α -theory, measures the distance to this space in terms of the averaged squared gradient. As our construction shows, there is a one-to-one correspondence between the asymptotic behavior of functions in this space and a hom -harmonic polynomials of degree k. However, there is no natural one-to-one correspondence between elements of this space and kth-order a hom -harmonic polynomials.
In a recent preprint by Armstrong, Kuusi, and Mourrat published a er our present work, a higher-order regularity result related to our present results is obtained [1], however, under a much stronger assumption on the decorrelation of coe cient elds (namely, nite range of dependence).
Before stating our results, let us recall the de nition of the correctors (φ, σ ). The corrector φ i satis es the equation: − ∇ · a(e i + ∇φ i ) = 0.
The ux correction q ij is de ned as: where a hom is the homogenized tensor, i.e., a hom e i is the expectation of a(e i + ∇φ i ). In our analysis, we will only use that a hom is some constant elliptic coe cient. We introduce the corresponding vector potential σ ijk (antisymmetric in its last two indices) by requiring that For the actual construction of a σ with stationary gradient, we refer to [11]; in this note, we just use the property (8). In the context of periodic homogenization, both the scalar and the vector potentials φ and σ may be chosen to be periodic. In stochastic homogenization, one cannot always expect to have a stationary (φ, σ ) (for instance in d ≤ 2 even in case of nite range of dependence) but, as mentioned above, we expect sublinear growth in the sense of (4) under mild ergodicity assumptions. Finally, let us give a brief historical overview on stochastic homogenization of elliptic PDEs. The qualitative theory of stochastic homogenization was initiated by Kozlov [13] and Papanicolaou and Varadhan [17]; the rst (nonoptimal) quantitative estimate-derived under the assumption of nite range of dependence-is due to Yurinskiȋ [19]. Naddaf and Spencer introduced spectral gap inequalities to quantify ergodicity in stochastic homogenization [16]. Gloria and Otto [12] were the rst to obtain optimal estimates on the size of the homogenization error in the linear elliptic case, though with nonoptimal stochastic integrability. Optimal stochastic integrability-however, with nonoptimal estimates on the size of the error-was obtained by Armstrong and Smart [3]. Finally, recently optimal error estimates with optimal stochastic integrability were established by Gloria and Otto [10] and Armstrong et al. [1]. For a more probabilistic viewpoint of stochastic homogenization of linear elliptic equations, see [15]. In the case of fully nonlinear elliptic equations, a logarithmic rate of convergence has been established by Ca arelli and Souganidis [7] under a very weak assumption on decorrelation; Armstrong and Smart [4] have obtained a power-law rate of convergence in the case of nite range of dependence.
Notation. Throughout the paper, we use the Einstein summation convention, i.e., we implicitly take the sum over an index whenever this index occurs twice. For example, b i ∂ i v is an alternative notation for (b · ∇)v and b i ∇v i is an alternative notation for d i=1 b i ∇v i . By C, we denote a generic constant whose value may be di erent in each appearance of the expression C; similarly, by e.g., C(d, λ), we denote a generic constant depending only on d and λ whose value again may be di erent for every use of the expression C(d, λ).
By E := {E ∈ R d×d : (E ij + E ji )(a hom ) ij = 0}, we denote the space of matrices E ij for which E ij x i x j is an a hom -harmonic second-order polynomial. The notation P (or P(x)) generally refers to a polynomial. By P k , we denote the space of homogeneous polynomials of degree k. By P k a hom , we denote the space of homogeneous polynomials of degree k which are a hom -harmonic. On the space P k , we introduce the norm ||P|| := sup x∈B 1 |P(x)|; note that any other norm on this nite-dimensional space would do as well, since we do not care for C(k)-constants.

Main results
The proof of our large-scale C k,α regularity theory relies in an essential way on the existence of kth-order correctors for the homogenization problem, which enable us to correct a hom -harmonic polynomials of degree k by adding a small (in the L 2 -sense) perturbation.
The ansatz for the deformation of an a hom -harmonic polynomial P, homogeneous of degree k (i.e., P ∈ P k a hom ), into an a-harmonic function u with the same growth behavior is motivated by homogenization: We consider P as the "homogenized solution of the problem solved by u, " so that we think in terms of the two-scale expansion u ≈ P + φ k ∂ k P and have that the error ψ P := u − (P + φ k ∂ k P) satis es −∇ · a∇ψ P = ∇ · ((φ k a − σ k )∇∂ k P). To construct u, we reverse the logic and rst construct a solution ψ P to the above elliptic equation and then set u := P + φ k ∂ k P + ψ P .
Our ψ P indeed enable us-in conjunction with the rst-order correctors φ i -to correct a hom -harmonic kth-order polynomials. Proposition 2. Let d ≥ 2, k ≥ 2, and let P ∈ P k a hom . Suppose that ψ P satis es (9). We then have −∇ · a∇(P + φ i ∂ i P + ψ P ) = 0.
Let us now state our C k,α large-scale regularity result.
Theorem 3 (C k,α large-scale excess-decay estimate). Let d ≥ 2, k ≥ 2, and suppose that (4) holds. Let u be an a-harmonic function. Let ψ P ≡ 0 for linear polynomials P (in order to simplify notation) and let ψ P be the functions constructed in Theorem 1 for higher-order polynomials.
Consider the kth-order excess Let 0 < α < 1 and let r 0 be large enough so that ε 2,r 0 ≤ ε 0 holds (the existence of such r 0 is ensured by (4)), where ε 0 = ε 0 (d, k, λ, α) > 0 is a constant de ned in the proof below. Then for all r, R ≥ r 0 with r < R the C k,α excess-decay estimate is satis ed.
Our large-scale C k+1,α excess-decay estimate entails the following kth-order Liouville principle.
Corollary 4 (kth-order Liouville principle). Let d ≥ 2, k ≥ 2, and suppose that the assumption (4) is satis ed. Then the following property holds: Any a-harmonic function u satisfying the growth condition is of the form with some a ∈ R, b ∈ R d , and P κ ∈ P κ a hom for 2 ≤ κ ≤ k (i.e., P κ is a homogeneous a hom -harmonic polynomial of degree κ). Here, the ψ P denote the higher-order correctors whose existence is guaranteed by Theorem 1.
In particular, the space of all a-harmonic functions satisfying (13) has the same dimension as if a was replaced by a constant coe cient, say a hom .
Note that the de ning Eq. (9) and the growth condition together determine the corrector of order k only up to a-harmonic "polynomials" of order k − 1: The rst-order corrector φ i is determined only up to an additive constant; the secondorder corrector ψ P (for a quadratic polynomial P) is determined only up to corrected a ne functions of the form x → ξ · (x + φ) + c with ξ ∈ R d and c ∈ R, and so on. Let us denote by P k a the space of solutions to the problem −∇ · a∇v = 0 which satisfy the growth condition With this notation, our higher-order correctors yield a canonical isomorphism of the quotient spaces P k a hom / P k−1 a hom ∼ = P k a / P k−1 a de ned by: for any P ∈ P k a hom . Note that this isomorphism is independent of the particular choice of the correctors φ and ψ P .
The basic strategy of the proof of Theorems 1 and 3 is as follows: • First, under the assumption that we already have constructed an appropriate kth-order corrector on a ball B R , we show a C k,α excess-decay estimate on large scales within this ball for a-harmonic functions (Lemma 14). This result directly implies Theorem 3 as soon as we have proven the existence of a corrector on R d (i.e., as soon as we have established Theorem 1). The basic idea for this rst part of the proof is a standard approach from regularity theory: We transfer the regularity properties of the constant-coe cient equation −∇ · a hom ∇u hom = 0 to the equation −∇ · a∇u = 0. To accomplish this, we employ an error estimate for the homogenization error. • Our C k,α estimate implies a C k−1,1 theory for a-harmonic functions on balls B R , provided that we have already constructed an appropriate kth-order corrector on B R . This is done in Lemma 17. • At last, we are able to build our corrector, starting from small balls and iteratively doubling the size of our balls: We decompose the right-hand side of Eq. (9) into contributions from dyadic annulli. In each step, we add the contribution from the next larger scale ξ 2 m r 0 P determined as the Lax-Milgram solution to the problem −∇ · a∇ξ 2 m r 0 to the corrector on the old scale ψ 2 m r 0 P . At this point, we make use of the C k−1,1 theory to show that a er possibly subtracting an appropriate k−1-th order a-harmonic "polynomial, " the new contribution ξ 2 m r 0 P displays kth-order decay in the interior {|x| < 2 m r 0 }, down to the ball {|x| < r 0 }. This ensures that on a ball of a given xed size r with r < 2 m r 0 , the contribution from the next larger scale does not destroy the smallness of the corrector. We are therefore able to construct the corrector on the next larger scale ψ 2 m+1 r 0 P as the sum of the corrector on the old scale ψ 2 m r 0 P and the new contribution ξ 2 m r 0 P minus the aforementioned a-harmonic "polynomial. " This iterative enlargement is carried out in Lemma 18 and nally enables us to prove Theorem 1.
• The kth-order Liouville principle stated in Corollary 4 is an easy consequence of our C k+1,α large-scale excess-decay estimate.

A C 2,α large-scale regularity theory for homogeneous elliptic equations with random coe cients
For the reader's convenience, we shall rst provide a proof for the C 2,α case of our theorems, as in this case, the proofs are less technical while already containing the key ideas. In particular, the overall structure of our proofs is the same as in the C k,α case. Since we shall use a somewhat simpli ed notation in the C 2,α case, let us reformulate the C 2,α case of our theorems using this notation.
Theorem 5 (Existence of second-order correctors). Let d ≥ 2 and suppose that the corrector φ and the ux-correction potential σ satisfy the growth assumption (4). Let r 0 be large enough so that ε 2,r 0 ≤ ε 0 holds [the existence of such r 0 is ensured by (4)], where ε 0 = ε 0 (d, λ) > 0 is a constant de ned in the proof below. Given any E ∈ R d×d , there exists a second-order corrector ψ E satisfying as well as for any r ≥ r 0 . Moreover, the corrector ∇ψ E depends linearly on E.
Due to the linear dependence of ψ E on E, below we shall also write E ij ψ ij in place of ψ E . Note that our second-order correctors indeed enable us-in conjunction with the rstorder correctors φ i -to correct a hom -harmonic second-order polynomials. Proposition 6. Let d ≥ 2 and let E ∈ E (i.e., assume that the polynomial E ij x i x j is a homharmonic). Suppose that ψ E satis es (14). We then have Our C 2,α large-scale regularity theorem reads as follows.
Theorem 7 (C 2,α large-scale excess-decay estimate). Let d ≥ 2 and suppose that (4) holds. Let u be an a-harmonic function. Let ψ E be the second-order corrector constructed in Theorem 5. Consider the second-order excess Let 0 < α < 1 and let r 0 be large enough so that ε 2,r 0 ≤ ε 0 holds [the existence of such r 0 is ensured by (4)], where ε 0 = ε 0 (d, λ, α) > 0 is a constant de ned in the proof below. Then for all r, R ≥ r 0 with r < R the C 2,α excess-decay estimate is satis ed.
Our large-scale excess-decay estimate entails the following C 2,α Liouville principle.
Corollary 8 (C 2,α Liouville principle). Let d ≥ 2 and suppose that the assumption (4) is satis ed. Then the following property holds: Any a-harmonic function u satisfying the growth condition lim inf with some a ∈ R, b ∈ R d , and E ∈ E (i.e., some E ∈ R d×d for which E ij x i x j is an a hom -harmonic polynomial).
Let us start with the proof of Proposition 6, which only requires a simple computation.
Proof of Proposition 6. Making use of the fact that E ij ((a hom ) ij + (a hom ) ji ) = 0 (in the third step below), we compute We therefore obtain which together with (14) implies our proposition.

The C 2,α excess-decay estimate
To establish our C 2,α excess-decay estimate, we make use of the following lemma, which essentially generalizes Theorem 7 to correctors which are only available on balls B R .
For any E ∈ E, denote byψ E a solution to the equation of the second-order corrector (14) on the ball B R (without boundary conditions); assume thatψ E depends linearly on E. Set εψ ,r,R := sup For an a-harmonic function u in B R , consider the second-order excess For any 0 < α < 1, there exists a constant ε min > 0 depending only on d, λ, and α such that the following assertion holds: Suppose that r 0 > 0 satis es ε r 0 +εψ ,r 0 ,R ≤ ε min . Then for all r ∈ [r 0 , R] the C 2,α excess-decay estimate is satis ed. Note that the in mum in (18) is actually attained, as the average integral in the de nition of Exc 2 (ρ) is a quadratic functional of b and E. Denote by b ρ,min and E ρ,min a corresponding optimal choice of b and E in (18). We then have the estimates and Proof of Theorem 7. Theorem 7 obviously follows from Lemma 9 by settingψ E := ψ E , with ψ E being the second-order corrector whose existence is guaranteed by Theorem 5.
The following lemma is essentially a special case of our C 2,α large-scale excess-decay estimate Lemma 9; it entails the general case of Lemma 9 (see below).

Lemma 10.
Let d ≥ 2 and let R, r > 0 satisfy r < R/4 and ε R ≤ 1. For any E ∈ E, denote byψ E a solution to the equation of the second-order corrector (14) on the ball B R (without boundary conditions); assume thatψ E depends linearly on E. For an a-harmonic function u in B R , consider again the second-order excess (18). Then the excess on the smaller ball B r is estimated in terms of the excess on the larger ball B R and our quantities ε R and ∇ψ E : We have Before proving Lemma 10, we would like to show how it implies Lemma 9.
Proof of Lemma 9. First choose 0 < θ ≤ 1/4 so small that the strict inequality C(d, λ)θ 4 < θ 2+2α is satis ed (with C(d, λ) being the constant from Lemma 10). Then, choose the threshold ε min for ε r 0 + εψ ,r 0 ,R so small that the estimate Let M be the largest integer for which θ M R ≥ r holds. Applying Lemma 10 inductively with R m := θ m−1 R, r m := θ m R for 1 ≤ m ≤ M, we infer Since we have trivially and since by de nition of M, we have r > θr M and thus It remains to show the estimates for |b r,min − b R,min | and |E r,min − E R,min | as well as the bounds for |b r,min | and |E r,min |. To do so, let us rst estimate the di erences |b R m ,min − b r m ,min | and From Lemma 11 below, we thus obtain Note that a similar estimate for the last increment |b r M ,min − b r,min | + R|E r M ,min − E r,min | can be derived analogously. Taking the sum with respect to m and recalling that R 1 = R and r m = R m+1 , we nally deduce It only remains to establish the last estimate for |b r,min | and |E r,min |. By the previous estimate, it is su cient to prove the corresponding bound for b R,min and E R,min . This in turn is a consequence of the inequality together with Lemma 11 below.
The following lemma quanti es the linear independence of the corrected polynomials where ε 0 = ε 0 (d) is to be de ned in the proof below. Then for any b ∈ R d and any E ∈ E, we have the estimate Proof. Poincaré's inequality (with zero mean) and the triangle inequality imply On the one hand, by transversality of constant, linear, and quadratic functions, we have On the other hand, we have by the triangle inequality and Poincaré's inequality, Putting these estimates together, by boundedness of the integrals in the previous line by ε 2 0 ρ 2 our assertion is established.
Proof of Lemma 10. In the proof of the lemma, we may assume that To see this, recall that the in mum in the de nition of Exc 2 (R) is actually attained. Denote the corresponding choices of b and E by b min and E min .
, we see that we may indeed assume (23): The new function is also a-harmonic due to (6) and Proposition 6.
We then apply Lemma 20 below to our function u. This yields an a hom -harmonic function u hom close to u which in particular satis es By inner regularity theory for elliptic equations with constant coe cients, the a hom -harmonic function u hom satis es Since −∇ · a hom ∇u hom = 0 holds, we infer E R,Taylor ij (a hom ) ij = 0 and therefore E R,Taylor ∈ E (note that E R,Taylor ij = E R,Taylor ji ). By Taylor's expansion of ∇u hom around x = 0, we deduce for any x ∈ B R/4 the bound Making use of the identity the previous estimate yields in connection with the bound for |∇ 3 u hom | and r < R/4 By the Caccioppoli inequality for the a-harmonic function The approximation property of u hom Combining the last three estimates and the equality This nally yields in connection with the above bounds on ∇ 2 u hom in B R/4 (recall that where in the last step we have used the inequality ε 2 The new bound directly implies the desired estimate.

The C 1,1 excess-decay estimate
We now show how our C 2,α excess-decay estimate for the second-order excess Exc 2 from Lemma 9 entails a C 1,1 excess-decay estimate for the rst-order excess Exc. Lemma 12. Let d ≥ 2 and R > 0. For any E ∈ E, denote byψ E a solution to the equation of the second-order corrector (14) on the ball B R (without boundary conditions); assume thatψ E depends linearly on E. There exists a constant ε min > 0 depending only on d and λ such that the following assertion holds: Suppose r 0 ∈ (0, R] is so large that ε r 0 ≤ ε min and hold. Let u be an a-harmonic function on B R . Then there exists b R ∈ R d for which the estimate holds for any r ∈ [r 0 , R]. Furthermore, b R depends linearly on u and satis es Proof. In Lemma 9, x α := 1/2. We then easily verify that Lemma 9 is applicable in our situation. Set b R := b r 0 ,min and E R := E r 0 ,min ; this implies that b R depends linearly on u. The estimate (21) takes the form Furthermore, applying Lemma 9 with r 0 playing the role of r and r playing the role of R, we deduce from (20) In conjunction with the two previous estimates, we infer Our lemma is therefore established.

Construction of second-order correctors
Using the C 1,1 theory established in the previous subsection, we now proceed to the construction of our second-order corrector. The following lemma provides the inductive step; starting from a function which acts as a corrector on a ball B R , we construct a function acting as a corrector on the ball B 2R .
Lemma 13. Let d ≥ 2 and let r 0 > 0 satisfy the estimate ε 2,r 0 ≤ ε 0 , where ε 0 = ε 0 (d, λ) is to be chosen in the proof below. Then the following implication holds: Let R = 2 M r 0 for some M ∈ N 0 . Suppose that for every E ∈ R d×d we have a solution ψ R E to the equation subject to the growth condition is a su ciently large constant to be chosen in the proof below. Assume furthermore that ψ R E depends linearly on E. Then for every E ∈ R d×d there exists a solution ψ 2R E to the equation min{1, 2 m r 0 /r}ε 2 m r 0 for all r ≥ r 0 . Furthermore, ψ 2R E depends linearly on E and we have Proof. To establish the lemma, we rst note that the assumptions of the lemma ensure that the C 1,1 excess-decay lemma (Lemma 12) is applicable on B R withψ E := ψ R E . To see this, we estimate for any r ∈ [r 0 , R] By choosing ε 0 > 0 small enough depending only on d and λ and C 1 (which is to be chosen at the end of this proof), we can ensure that the assumption of Lemma 12 regarding smallness of εψ ,r 0 ,R is satis ed. Let now ξ R E be the weak solution on R d with square-integrable gradient, which is unique up to additive constants and whose existence follows from the Lax-Milgram theorem, to the problem Obviously, ∇ξ R E depends linearly on E; a er xing the additive constant, for e.g., by requiring As ξ R E is a-harmonic in B R , Lemma 12 now implies the existence of some b R E ∈ R d for which the estimates and hold for all r ∈ [r 0 , R] and which linearly depends on E.
. The combination of both r-ranges yields In total, we see that is the desired function (note in particular that the last term is a-harmonic), provided we choose C 1 to be the constant appearing in (27).
We now establish existence of second-order correctors by means of the previous lemma.
Proof of Theorem 5. We just need to construct an "initial" second-order corrector ψ r 0 E subject to the properties of Lemma 13; then Lemma 13 yields a sequence (ψ 2 m r 0 E ) m which is a Cauchy sequence in H 1 (B R ) for every R > 0 due to the last estimate in the lemma and our assumption (4) which implies summability of ε 2 m r 0 . Thus, the limit ψ E satis es the Eq. (14) in the whole space, depends linearly on E, and satis es the estimate min{1, 2 m r 0 /r}ε 2 m r 0 for any r ≥ r 0 .
To construct ψ r 0 E , just use Lax-Milgram to nd the solution ψ r 0 E on R d with squareintegrable gradient (unique up to an additive constant) to the equation Obviously, a er xing the additive constant appropriately ψ r 0 E depends linearly on E. Furthermore, we have the energy estimate i.e., for any r ≥ r 0 We note that this provides the starting point for Lemma 13, possibly a er enlarging the constant C 1 in the statement thereof.

Proof of the C 2,α Liouville principle
The C 2,α Liouville principle (Corollary 8) is an easy consequence of our large-scale excessdecay estimate (Theorem 7).
Proof of Corollary 8. Let α ∈ (0, 1) be such that holds. By the Caccioppoli estimate, we deduce Fix r ≥ r 0 . The excess-decay estimate from Theorem 7 yields together with the trivial bound Exc 2 (R) ≤ − B R |∇u| 2 dx that Passing to the limit R → ∞, we deduce that Exc 2 (r) = 0 holds for every r ≥ r 0 . Therefore, on every B r with r ≥ r 0 , ∇u can be represented exactly as the derivative of a corrected polynomial of second order (since the in mum in the de nition of Exc 2 is actually attained, as noted at the beginning of the proof of Lemma 10), i.e., we have in B r for some b r ∈ R d and some E r ∈ E. It is not di cult to show that for r large enough, the b r and E r are actually independent of r and de ne some common b ∈ R d and E ∈ E: For example, one may use Lemma 9 to compare the b r , E r for two di erent radii r 1 , r 2 ≥ r 0 ; the estimate for |b r 1 − b r 2 | and |E r 1 − E r 2 | then contains the factor Exc 2 (max(r 1 , r 2 )) and is therefore zero. Moreover, the gradient ∇u determines the function u itself up to a constant, i.e., we have for some a ∈ R, some b ∈ R d , and some E ∈ E ⊂ R d×d .

A C k,α large-scale regularity theory for elliptic equations with random coe cients
We now generalize our proofs from the C 2,α case in order to correct polynomials of order k and obtain our C k,α large-scale regularity theory. We proceed by induction in k.
To establish our C k,α regularity theory, let us rst show Proposition 2, which -like the proof of Proposition 6 in the C 2,α case -only requires a simple computation.

This yields
which together with (9) implies our proposition.

The C k,α excess-decay estimate
To establish our C k,α excess-decay estimate, we make use of the following lemma, which essentially generalizes Theorem 3 to correctors that are only available on balls B R . Lemma 14. Let d ≥ 2 and k ≥ 2. Suppose that Theorem 1 holds for orders 2, . . . , k − 1, and set ψ P ≡ 0 for rst-order polynomials P to simplify notation. For any P ∈ P k a hom , denote byψ P a solution to the Eq. (9) on the ball B R (without boundary conditions); assume that theψ P depend linearly on P. Set εψ ,r,R := sup For an a-harmonic function u in B R , consider the kth-order excess For any 0 < α < 1, there exists a constant ε min > 0 depending only on d, k, λ, and α such that the following assertion holds: Suppose that r 0 > 0 satis es ε 2,r 0 + εψ ,r 0 ,R ≤ ε min . Then for all r ∈ [r 0 , R] the C k,α excessdecay estimate is satis ed. Note that the in mum in (29) is actually attained, as the average integral in the de nition of Exc 2 (ρ) is a quadratic functional of P κ . Denote by P ρ,min κ a corresponding optimal choice of P κ in (29). We then have the estimates: and k κ=1 R 2(κ−1) ||P r,min Proof of Theorem 3. Once we have shown Theorem 1, Theorem 3 obviously follows from Lemma 14 by settingψ P k := ψ P k , with ψ P k being the kth-order corrector whose existence is established in Theorem 1.
The following lemma is essentially a special case of our C k,α large-scale excess-decay estimate Lemma 14; it entails the general case of Lemma 14 (see below).

Lemma 15.
Let d ≥ 2, k ≥ 2, and let R, r > 0 satisfy r < R/4 and ε 2,R ≤ ε 0 (d, k − 1, λ), with ε 0 (d, k − 1, λ) being the constant from Theorem 1 for the orders 2, . . . , k − 1. Assume that Theorem 1 holds for orders 2, . . . , k − 1, and let ψ P ≡ 0 for linear polynomials P in order to simplify notation. For any P ∈ P k a hom , denote byψ P a solution to the Eq. (9) on the ball B R (without boundary conditions); assume thatψ P depends linearly on P. For an a-harmonic function u on B R , consider again the kth-order excess (29). Then the excess on the smaller ball B r is estimated in terms of the excess on the larger ball B R and our quantities ε 2,R and ∇ψ P : We have Before proving Lemma 15, we would like to show how it implies Lemma 14.
Proof of Lemma 14. First choose 0 < θ ≤ 1/4 so small that the strict inequality C(d, k, λ)θ 2k < θ 2(k−1)+2α is satis ed (with C(d, k, λ) being the constant from Lemma 15). Then, choose the threshold ε min for ε 2,r 0 + εψ ,r 0 ,R so small that the estimate Let M be the largest integer for which θ M R ≥ r holds. Applying Lemma 15 inductively with R m := θ m−1 R, r m := θ m R for 1 ≤ m ≤ M, we infer Since we have trivially From Lemma 16 below, we thus obtain A similar estimate for the last increment k κ=1 R κ−1 ||P r M ,min κ − P r,min κ || can be derived analogously. Taking the sum with respect to m and recalling that R 1 = R and r m = R m+1 , we nally deduce It only remains to establish the last estimate for ||P r,min κ ||. By the previous estimate, it is su cient to prove the corresponding bound for ||P R,min κ ||. This in turn is a consequence of the obvious inequality in conjunction with Lemma 16 below.
The following lemma quanti es the linear independence of the corrected polynomials P κ + φ i ∂ i P κ + ψ P κ (with 1 ≤ κ ≤ k); it is needed for the previous proof.

Lemma 16.
Suppose that the functions φ andψ P κ (2 ≤ κ ≤ k) satisfy where ε 0 = ε 0 (d, k) is to be de ned in the proof below. Setψ P ≡ 0 for linear polynomials P in order to simplify notation. Then for any P κ ∈ P κ a hom (1 ≤ κ ≤ k), we have the estimate Proof. Poincaré's inequality (with zero mean) and the triangle inequality imply On the one hand, by transversality of constant, linear, homogeneous second-order, …, and homogeneous kth-order polynomials, we have On the other hand, we have by the triangle inequality and Poincaré's inequality, Putting these estimates together, by boundedness of the integrals in the previous line by ε 2 0 ρ 2(κ−1) , our assertion is established.
Proof of Lemma 15. In the proof of the lemma, we may assume that To see this, recall that the in mum in the de nition of Exc k (R) is actually attained. Denote the corresponding choices of P κ by P min κ . Replacing u by u − k−1 κ=1 (P min κ + φ i ∂ i P min κ + ψ P min κ ) − (P min k + φ i ∂ i P min k +ψ P min k ), we see that we may indeed assume (34): The new function is also a-harmonic due to (6) and Proposition 2.
We then apply Lemma 20 below to our function u. This yields an a hom -harmonic function u hom close to u which in particular satis es By inner regularity theory for elliptic equations with constant coe cients, the a hom -harmonic function u hom satis es Let P R,Taylor κ (for 1 ≤ κ ≤ k) be the term of order κ in the Taylor expansion of u hom at x 0 = 0. We now show (for κ ≥ 2, as for κ = 1 this assertion is trivial) that P R,Taylor κ ∈ P κ a hom . The term-wise Hessian of the Taylor series of u hom yields the Taylor series of ∇ 2 u hom . We now know that a hom : ∇ 2 u hom = 0; thus, the Taylor series of a hom : ∇ 2 u hom is identically zero and by equating the coe cients, we deduce a hom : ∇ 2 P R,Taylor κ = 0 for 2 ≤ κ ≤ k. As the term-wise derivative of the Taylor series of u hom yields the Taylor series of ∇u hom , we obtain by the standard error estimate for the Taylor expansion of ∇u hom at x 0 = 0 for any x ∈ B R/4 the estimate Making use of the identity the previous estimate yields in connection with the bound for |∇ k+1 u hom | and r < R/4 By the Caccioppoli inequality for the a-harmonic function x i + φ i (6), we have The approximation property of u hom + φ i ∂ i u hom in B R/2 from Lemma 20 below implies Combining the last three estimates and the equality This nally yields in connection with the bounds on ∇ κ u hom in B R/4 (35) which in particular where in the last step, we have used the inequality ε 2 and ε R ≤ ε 2,R as well as (10) for 2 ≤ κ ≤ k − 1. Our new estimate now implies the desired bound.

The C k−1,1 excess-decay estimate
Like in the C 2,α case, we now show how the C k,α excess-decay estimate for the kth-order excess Exc k (in Lemma 14) entails a C k−1,1 excess-decay estimate for the (k − 1)th-order excess Exc k−1 .
Lemma 17. Let d ≥ 2, k ≥ 2, and R > 0. Assume that Theorem 1 holds for the orders 2, . . . , k − 1, and let ψ P ≡ 0 for linear polynomials P in order to simplify notation. For any P ∈ P k a hom , denote byψ P a solution to the Eq. (9) on the ball B R (without boundary conditions); assume that theψ P depend linearly on P. Then there exists a constant ε min > 0 depending only on d, k, and λ such that the following assertion holds: Suppose r 0 ∈ (0, R] is so large that ε 2,r 0 ≤ ε min and hold. Let u be an a-harmonic function on B R . Then there exists P R κ ∈ P κ a hom (1 ≤ κ ≤ k − 1) for which the estimate holds for any r ∈ [r 0 , R]. Furthermore, the P R κ depend linearly on u and satisfy Proof. In Lemma 14, x α := 1/2. We then easily verify that Lemma 14 is applicable in our situation. Set P R κ := P r 0 ,min κ ; this implies that the P R κ depend linearly on u. The estimate (32) takes the form: Furthermore, applying Lemma 14 with r 0 playing the role of r and r playing the role of R, we deduce from (31) We now estimate Exc k (R) Exc k (R) + C(d, k, λ)||P r,min k || 2 r 2(k−1) In conjunction with the two previous estimates, we infer Our lemma is therefore established.

Construction of correctors of order k
Using the C k−1,1 theory established in the previous subsection, we now proceed to the construction of our kth-order corrector. The following lemma provides the inductive step; starting from a function which acts as a kth-order corrector on a ball B R , we construct a function acting as a kth-order corrector on the ball B 2R .
Lemma 18. Let d ≥ 2, k ≥ 2, and assume that Theorem 1 holds for the orders 2, . . . , k − 1. Let r 0 > 0 satisfy the estimate ε 2,r 0 ≤ ε 0 , where ε 0 = ε 0 (d, k, λ) is to be chosen in the proof below. Then the following implication holds: Let R = 2 M r 0 for some M ∈ N 0 . Suppose that for every P ∈ P k , we have a solution ψ R P to the equation −∇ · a∇ψ R P = ∇ · (χ B R (φ i a − σ i )∇∂ i P) subject to the growth condition min{1, 2 m r 0 /r}ε 2 m r 0 for all r ≥ r 0 , where C 1 (d, k, λ) is a su ciently large constant to be chosen in the proof below. Assume furthermore that ψ R P depends linearly on P. Then for every P ∈ P k there exists a solution ψ 2R P to the equation min{1, 2 m r 0 /r}ε 2 m r 0 for all r ≥ r 0 . Furthermore, ψ 2R P depends linearly on P and we have Proof. To establish the lemma, we rst note that the assumptions of the lemma ensure that the C k−1,1 excess-decay lemma (Lemma 17) is applicable on B R withψ P := ψ R P . To see this, we estimate for any r ∈ [r 0 , R] By choosing ε 0 > 0 small enough depending only on d, k, λ, and C 1 (which is to be chosen at the end of this proof), we can ensure that the assumption of Lemma 17 regarding smallness of εψ ,r 0 ,R is satis ed. We now turn to the construction of ψ 2R P − ψ R P and to that purpose denote by ξ R P the weak solution on R d with zero mean in B 2R and square-integrable gradient, whose existence and uniqueness follows by the Lax-Milgram theorem, to the problem . Obviously, ξ R P depends linearly on P. Furthermore, by ellipticity, we have the estimate The last estimate in turn implies We now obtain ψ 2R P −ψ R P by modifying ξ R P by an a-harmonic function of degree k−1. As ξ R P is a-harmonic in B R , Lemma 17 now implies the existence of some P R κ,P ∈ P κ for 1 ≤ κ ≤ k − 1 which depend linearly on P and for which the estimates ≤ C(d, k, λ)||P|| 2 r 2(k−2) R 2 ε 2 2R . The combination of both r-ranges yields In total, we see that is the desired function (note in particular that the last term is a-harmonic), provided we choose C 1 to be the constant appearing in (39).
We now establish existence of kth-order correctors by the previous lemma.
Proof of Theorem 1. We just need to construct an "initial" kth-order corrector ψ r 0 P subject to the properties of Lemma 18; then Lemma 18 yields a sequence (ψ 2 m r 0 P ) m which (a er subtracting appropriate constants) is a Cauchy sequence in H 1 (B R ) for every R > 0 due to the last estimate in the lemma and our assumption (4) which implies summability of ε 2 m r 0 . Thus, the limit ψ P satis es the Eq. (9) in the whole space, depends linearly on P, and satis es the estimate min{1, 2 m r 0 /r}ε 2 m r 0 ≤ C 1 (d, k, λ)||P||ε 2,r for any r ≥ r 0 .
To construct ψ r 0 P , we use Lax-Milgram to nd the (unique) solution ψ r 0 P on R d with squareintegrable gradient and zero mean on B r 0 to the equation: Obviously, ψ r 0 P depends linearly on P. Furthermore, we have the energy estimate We therefore get This yields in particular for any r ≥ r 0 ≤ C(d, k, λ)||P|| 2 r 2(k−1) min{1, (r 0 /r) 2 }ε 2 r 0 . We note that this provides the starting point for Lemma 18, possibly a er enlarging the constant C 1 in the statement thereof.

Proof of the kth-order Liouville principle
Like in the C 2,α case, the C k,α Liouville principle (Lemma 19 below) is an easy consequence of our large-scale excess-decay estimate (Theorem 3). The kth-order Liouville principle (Corollary 4) in turn is an easy consequence of the C k+1,α Liouville principle. Lemma 19. Let d ≥ 2, k ≥ 2, and suppose that the assumption (4) is satis ed. Then the following property holds: Any a-harmonic function u satisfying the growth condition with some a ∈ R, b ∈ R d , and P κ ∈ P κ a hom for 2 ≤ κ ≤ k (i.e., P κ is a homogeneous a hom -harmonic polynomial of degree κ). Here, the ψ P denote the higher-order correctors whose existence is guaranteed by Theorem 1.
Proof of Corollary 4. Obviously, (13) entails (40) with k + 1 in place of k and e.g., α := 1 2 . By Lemma 19, any a-harmonic function u subject to condition (13) must be of the form: with some a ∈ R, b ∈ R d , and P κ ∈ P κ a hom for 2 ≤ κ ≤ k + 1. Our stronger growth condition (13) however shows that we have P k+1 ≡ 0: Since the φ i grow sublinearly (2) and since ψ P k+1 grows slower than a polynomial of degree k + 1 (10), we see that for large |x| the term P k+1 would be the dominating term in (41) if it was nonzero, contradicting our growth condition (13).
Proof of Lemma 19. Let α ∈ (0, 1) be such that holds. By the Caccioppoli estimate, we deduce lim inf Fix r ≥ r 0 . The excess-decay estimate from Theorem 3 together with the trivial bound Exc k (R) ≤ − B R |∇u| 2 dx yields Passing to the lim inf R → ∞, we deduce that Exc k (r) = 0 holds for every r ≥ r 0 . Therefore, on every B r with r ≥ r 0 , ∇u can be represented exactly as the derivative of a corrected polynomial of kth order (since the in mum in the de nition of Exc k is actually attained, as noted at the beginning of the proof of Lemma 15), i.e., we have in B r for some b r ∈ R d and some P r κ ∈ P κ a hom (2 ≤ κ ≤ k); recall that we have used the convention ψ P ≡ 0 for linear polynomials P. It is not di cult to show that for r large enough, the b r and P r κ are actually independent of r and de ne some common b ∈ R d and P κ ∈ P κ a hom : For example, one may use Lemma 14 to compare the b r , P r κ for two di erent radii r 1 , r 2 ≥ r 0 ; the estimate for |b r 1 − b r 2 | and ||P r 1 κ − P r 2 κ || then contains the factor Exc k (max(r 1 , r 2 )) and is therefore zero. Moreover, the gradient ∇u determines the function u itself up to a constant, i.e., we have for some a ∈ R, b ∈ R d , and P κ ∈ P κ a hom (2 ≤ κ ≤ k).

Appendix A. Approximation of a-harmonic functions by corrected a homharmonic functions
Our proofs make use of the following lemma, which is implicitly derived in the course of the proof of Lemma 2 in [11]. For the reader's convenience, we recall its proof here. The lemma essentially states that an a-harmonic function u on a ball B R may be approximated on the ball B R/2 up to a small error (of order ε 1/(d+1) 2 R ) by an appropriate a homharmonic function u hom and correcting this function u hom using the rst-order corrector φ i .
The purpose of the lemma is the same as in classical elliptic regularity theory: The function u hom satis es an elliptic equation with constant coe cients, i.e., it is smooth and good estimates for its higher derivatives are available. In our proof above, we show by the present lemma that this high regularity of u hom transfers (in an appropriate sense) to u itself.
Lemma 20. Let R > 0 and let u be a-harmonic on B R . Suppose that ε R ≤ 1 [with ε R as de ned in (3)]. Then there exists an a hom -harmonic function u hom on B R/2 satisfying the following two properties: First, we have the energy estimate Second, the "corrected" function u hom + φ i ∂ i u hom is a good approximation for u in the sense that Proof. Choose some R ′ ∈ [ 3 4 R, R] for which holds. Let u hom be the a hom -harmonic function in B R ′ which coincides with u on ∂B R ′ . Testing the equation −∇ · a hom ∇u hom = 0 with u hom − u (note that this test function is admissible since we have u hom − u = 0 on ∂B R ′ ), we infer by ellipticity of a and (in the second step) Young's inequality which because of R/2 ≤ R ′ ≤ R gives the desired energy estimate. It remains to establish the approximation property of u hom + φ i ∂ i u hom . Denote by η 0 : R → R a smooth function with η 0 (s) = 1 for s ≥ 1 and η 0 (s) = 0 for s ≤ 0. Let 0 < ρ < R/4 and set η(x) := η 0 (2(R ′ − ρ/2 − |x|)/ρ). Note that we have |∇η| ≤ C(d)/ρ as well as η ≡ 0 outside of B R ′ −ρ/2 and η ≡ 1 in B R ′ −ρ . Due to ρ ≤ R/4, we also have R ′ − ρ ≥ R/2. We will optimize in this "boundary layer thickness" ρ at the end of the proof.
Let us abbreviate where the purpose of η is to have v ≡ 0 on ∂B R ′ . The desired approximation property of u hom +φ i ∂ i u hom as stated in the lemma will be a consequence of an appropriate energy estimate for v (recall that we have η ≡ 1 in B R/2 since ρ < R/4 and R ′ > 3R/4).
Testing the weak formulation of this equation with v (recall that v ≡ 0 on ∂B R ′ ) and using the ellipticity of a, we deduce using Young's inequality and the properties of η Since our function u hom is a hom -harmonic, we have the regularity estimates sup where p := 2d/(d − 1): The rst estimate is a standard constant coe cient interior regularity estimate (which is a consequence example of an iterative application of Theorem 4.9 in [9] and the Sobolev embedding). The second estimate follows by combining 1) the existence of an extensionū of u hom subject to the estimate ||∇ū|| L p (B R ′ ) ≤ C(d)||∇ tan u hom || L 2 (∂B R ′ ) and 2) the Calderon-Zygmund estimate on B R ′ , which reads ||∇w|| L p (B R ′ ) ≤ C(d, λ)||∇ū|| L p (B R ′ ) for any solution w ∈ H 1 (B R ′ ) with w −ū ∈ H 1 0 (B R ′ ) to the equation −∇ · a hom ∇w = 0. For the latter estimate, see Theorem 7.1 in [9].

Appendix B. Failure of Liouville principle for smooth uniformly elliptic coecient elds
We now provide the argument that smoothness of a uniformly elliptic coe cient eld does not prevent Liouville's theorem from failing: Even for smooth uniformly elliptic coe cient elds, sublinearly growing harmonic functions are not necessarily constant, implying a failure even of the zeroth-order Liouville theorem.
Proposition 21. For any α ∈ (0, 1) there exists a smooth, bounded, and uniformly elliptic symmetric coe cient eld a on R 2 such that the following holds: There exists a smooth function u which is a-harmonic and satis es Proof. By a classical example in dimension d = 2 [18], for any exponent α ∈ (0, 1), there exists a uniformly elliptic, symmetric coe cient eld a 0 of a scalar equation, and a weakly a 0harmonic function u 0 (in particular, it is locally integrable and of locally integrable gradient) whose modulus on average grows like |x| α , for instance as expressed by: Moreover, in this example a 0 and u 0 are homogeneous and smooth outside the origin.
We now argue that this example may be post-processed to an example of an everywhere smooth uniformly elliptic symmetric coe cient eld a and a smooth a-harmonic function u such that still (45) holds. Indeed, because of (47), we can easily construct a uniformly elliptic coe cient eld a that agrees with a 0 outside of B 1 and is smooth. Next, we observe that (47) also implies (using d = 2 and α > 0) that ∇u 0 is locally square integrable, so that by Riesz' representation theorem, there exists a weak solution of − ∇ · a∇w = ∇ · (a − a 0 )∇u 0 (48) in the sense that w and its gradient are locally integrable and that Equation (48) is made such that u = u 0 + w is a weak solution (i.e., locally integrable with locally integrable gradient) of −∇ · a∇u = 0, and thus smooth since a is smooth by classical uniqueness and regularity results. It remains to give the argument in favor of (45), which in view of (46) follows once we show that (49) implies in particular for large R − B R w 2 dx This is a well-known argument related to "bounded mean oscillation": By Poincaré's estimate with mean value zero, we have on every dyadic ball around the origin Proposition 22. There exists a smooth, bounded, and uniformly elliptic symmetric coe cient eld a on R 3 such that the following holds: There exists a smooth map u : R 3 → R 3 which is a-harmonic and satis es where α = 1 2 (1 − 3
Proof. By a classical example of De Giorgi in dimension d = 3 (Chapter 9.1.1, [9]), there exists a bounded, symmetric, and uniformly elliptic coe cient eld a 0 which is radial and smooth away from the origin, for which the map ) is a 0 -harmonic. Choose a to be a smooth, bounded, and uniformly elliptic coe cient eld which agrees with a 0 outside of the unit ball B 1 .
We now show that the a 0 -harmonic map u 0 may be modi ed to yield an a-harmonic map u with the same decay properties on large scales. To construct the di erence u − u 0 , let w be the Lax-Milgram solution (which is unique up to a constant) to the problem − ∇ · a∇w = ∇ · (a − a 0 )∇u 0 . (54) Since a−a 0 is supported in B 1 , since a and a 0 are bounded, and since ∇u 0 belongs to L 2 loc (R 3 ), we deduce by the standard energy estimate |∇w| 2 dx ≤ C |(a − a 0 )∇u 0 | 2 dx ≤ C.
Poincaré's inequality now implies for any R > 0 which entails We therefore deduce that the sequence − B 2 n w dx is Cauchy: We have for any N > n ≥ 0 Possibly adding a constant to w (to ensure that the limit of the above sequence is zero), we therefore may assume that − B 2 n w dx ≤ C2 −n/2 .