Wild solutions of the Navier-Stokes equations whose singular sets in time have Hausdorff dimension strictly less than 1

We prove non-uniqueness for a class of weak solutions to the Navier-Stokes equations which have bounded kinetic energy, integrable vorticity, and are smooth outside a fractal set of singular times with Hausdorff dimension strictly less than 1.


INTRODUCTION
Throughout this paper we consider the incompressible three dimensional Navier-Stokes equations: posed on the torus T 3 = [−π, π] 3 . We consider solutions of zero mean, i.e.´T 3 v(x, t)dx = 0 for all t ∈ [0, T ]. The notion of weak solution of (1.1) that we work with in this paper is that of a distributional solution, which has bounded kinetic energy, and is strongly continuous in time: . Given any zero mean initial datum v 0 ∈ L 2 , we say v ∈ C 0 ([0, T ); L 2 (T 3 )) is a weak solution of the Cauchy problem for the Navier-Stokes equations (1.3) if the vector field v(·, t) is weakly divergence-free for all t ∈ [0, T ), has zero mean, and T 3 v 0 · ϕ(·, 0)dx +ˆT 0ˆT 3 v · (∂ t ϕ + (v · ∇)ϕ + ∆ϕ)dxdt = 0 holds for any test function ϕ ∈ C ∞ 0 (T 3 × [0, T )) such that ϕ(·, t) is divergence-free for all t.

1.2)
Here P H is the Helmholtz projector and e t∆ f is the heat extension of f . Our main result of this paper is as follows. Moreover, for every such v there exists a zero Lebesgue measure set of times Σ T ⊂ (0, T ] with Hausdorff (in fact box-counting) dimension less than 1 − β, such that for every t ∈ [0, T ] \ Σ T there exists ε t > 0 such that v ∈ C ∞ ((t − ε t , t + ε t ) × T 3 ).
That is, the weak solution v is almost everywhere smooth.
The outline of the proof of Theorem 1.1 is given in Section 2, while the detailed estimates are done in Sections 3-5. 1 Remark 1.2 (Non-uniqueness of weak solutions for strong initial datum). Theorem 1.1 immediately implies that weak solutions of the Cauchy problem for the Navier-Stokes equation (1.1), cf. Definition 1.1, are not unique. The cheap way to see this is to take any T > 0, u (1) ≡ 0, and u (2) to be any nontrivial mean-free solution of the Navier-Stokes equation on [0, T ] (e.g. a shear flow). Then the weak solution v given by Theorem 1.1 is nontrivial on [0, T ], and thus 0 is not the only weak solution with 0 initial datum. Conversely, we note that taking u (1) to be any nontrivial solution to the Navier-Stokes equation, and u (2) ≡ 0, Theorem 1.1 gives a counterexample to backwards (in time) uniqueness for weak solutions of (1.1), in the sense of Definition 1.1.
More generally, we emphasize that Theorem 1.1 proves the non-uniqueness of weak solutions to the Cauchy problem for the Navier-Stokes equation (1.1) for any strong initial datum. To see this, consider any v 0 ∈Ḣ 3 and take T = c v 0 −1 H 3 , where c > 0 is a sufficiently small universal constant (cf. Proposition 3.1). Then there exists a unique solution u (1) ∈ C 0 ([0, T ]; H 3 ) to the Cauchy problem (1.3) with datum v 0 . Moreover u (1) (T ) L 2 ≤ v 0 L 2 . However, using Theorem 1.1 one can glue to this solution the shear flow u (2) (x 1 , x 2 , x 3 , t) = (Ae −t sin(x 2 ), 0, 0). Then if A is chosen such that Ae −T > 2 v 0 L 2 , we have v(T ) L 2 = u (2) Therefore v is a weak solution to (1.3) with datum v 0 , but v is not equal to the smooth solution u (1) at time T .
While for the above argument we have considered v 0 ∈Ḣ 3 , it is clear that Theorem 1.1 also implies the nonuniqueness of weak solutions to the Cauchy problem for (1.1) for any initial datum for which one has unique local in time solvability of (1.1) (examples include v 0 ∈Ḣ 1/2 cf. [15], v 0 ∈ L 3 with zero mean cf. [24]; v 0 ∈ BM O −1 which is small and has zero mean cf. [26]; see [29] for further details). Indeed, for any such initial datum the unique local in time solution u (1) is smooth in positive time, and hence for any ε > 0 we have u (1) (·, ε) ∈Ḣ 3 . We then apply Theorem 1.1 on the time interval [ε, T ], rather than [0, T ], in order to glue the strong solution to a shear flow with kinetic energy which is either strictly larger, or strictly less at time T .
1.1. Background. We make a few comments concerning different notions of solutions to the Navier-Stokes equation, other than from those in Definition 1.1 (see [30] for a more detailed discussion). The weakest notion of solution to the Cauchy problem for (1.1) is that of a very weak solution: these are distributional solutions of (1.1) which only lie in C 0 weak (0, T ; L 2 ), and are weakly divergence free. However, one typically proves the existence of solutions which are stronger than this.
Indeed, for any L 2 initial datum v 0 , Leray [31] constructed a distributional solution v ∈ C 0 weak (0, ∞; L 2 ) ∩ L 2 (0, ∞;Ḣ 1 ), and obeys the energy inequality v(t) 2 L 2 + 2´t s ∇v(τ ) 2 L 2 dτ ≤ v(s) 2 L 2 for a.e. s ≥ 0, and all t > s. See also the work of Hopf [18] on bounded domains. These are the Leray-Hopf weak solutions. One nice feature of Leray-Hopf weak solutions is that they possess epochs of regularity, i.e. many time intervals on which they are smooth. In fact, already Leray [31] made the observation that these weak solutions are almost everywhere in time smooth, since the putative singular set of times Σ T has Hausdorff dimension ≤ 1 /2. This fact follows directly from two ingredients: the fact that for v 0 ∈ H 1 the maximal time of existence of a unique smooth solution is bounded from below by c v 0 −4 H 1 , and a Vitali-type covering lemma which may be combined with the L 2 t H 1 x information provided by the energy inequality. Scheffer [42] went further to prove that the 1 /2-dimensional Hausdorff measure of Σ T is 0. These results were strengthened to bounds on the box-counting dimension for Σ T , cf. [41,27]. See [30,40] for further references. Remark 1.3 (Weak solutions with partial regularity in time). We note that while the weak solutions constructed in Theorem 1.1 are not Leray-Hopf, they give the first example of a mild/weak solution to the Navier-Stokes equation whose singular set of times Σ T ⊂ (0, T ] is both nonempty, and has Hausdorff (in fact, box-counting) dimension strictly less than 1. This is in contrast with the prior work [5], where Σ T has dimension 1. It is in an interesting open problem to construct weak solutions to (1.1), in the sense of Definition 1.1, where the 1 /2-dimensional Hausdorff measure of the nonempty set of singular times is 0.
A fundamental step towards understanding the uniqueness and smoothness of weak solutions was to introduce the concept of a suitable weak solution, by Scheffer [42] and Caffarelli-Kohn-Nirenberg [6]. Suitable weak solutions obey a localized in space-time version of the energy inequality, and they have partial regularity in space and time: the putative singular set of points in space-time has 1-dimensional parabolic Hausdorff measure is 0. See the reviews [40,30] for more recent extensions and further references.
The uniqueness of suitable weak solutions or of Leray-Hopf weak solutions is an outstanding open problem. The weak-strong uniqueness result of Prodi-Serrin [39,43] states that if there exists a weak/mild solution v ∈ L ∞ t L 2 x ∩ 2 L 2 tḢ 1 x ∩ L p t L q x of the Cauchy problem for (1.1), with 2 /p + 3 /q ≤ 1 1 and p < ∞, and if u is a Leray-Hopf weak solution with the same initial datum, then u ≡ v. This is a conditional uniqueness result within the class of Leray-Hopf weak solutions. Moreover, the solutions are smooth in positive time [28]. The L ∞ t L 3 x endpoint was established in [21]. Similar weak-strong uniqueness results hold within the class of mild solutions, except the q = 3 endpoint which requires continuity in time [13,16,34]. See [30,Chapter 12] for further references. A very interesting conjecture of Jia-Šverák [22,23] essentially states that the Prodi-Serrin uniqueness criteria are sharp, and that the non-uniqueness of Leray-Hopf weak solutions may already be expected in the regularity class L ∞ t L 3,∞ x . Compelling numerical evidence in support of this conjecture was recently provided by Guillod-Šverák [17]. A related interesting open problem is to establish the non-uniqueness of mild/weak solutions to (1.1) in the regularity class C 0 t L q x ∩ L 2 t H 1 x , for any q ∈ [2, 3). We conclude this subsection by revisiting the non-uniqueness result of Remark 1.2, for rough initial datum: Remark 1.4 (Non-uniqueness of very weak solutions for any L 2 initial datum). If instead of the weak solutions of Definition 1.1 we consider very weak solutions of (1.1), so they only lie in C 0 weak (0, T ; L 2 ), then Theorem 1.1 implies that the non-uniqueness for the Cauchy problem holds for any L 2 initial datum of zero mean, within the class of very weak solutions. Indeed, for any such datum, by the work of Leray there exists at least one very weak solution u to the Cauchy problem for (1.1), which in fact is smooth most of the time. Pick any regular time t 0 > 0 of u, and let v 0 = u(t 0 ) ∈Ḣ 3 . We then apply the argument of Remark 1.2 on the time interval [t 0 , t 0 + T ], with u (1) being the unique local in time smooth solution of (1.1) with initial datum v 0 at time t 0 . Note that by weak-strong uniqueness we in fact have that the Leray solution u is equal to u (1) on [t 0 , t 0 + T ]. In view of Theorem 1.1 we can construct a very weak solution v which is equal to u on [0, t 0 + T /3], and equal to a shear flow of our choice on [t 0 + 2T /3, T ]. This solution v is smooth except for a set of times of Hausdorff dimension < 1, and is different from the Leray solution u.
The system (1.3) was first considered by Lions in [32,33] for α in the critical and subcritical regime α ≥ 5 /4. Lions proved the existence and uniqueness of Leray-weak solutions, for any L 2 initial datum. These solutions are regular in positive time. In [45] it was proven that slightly below the critical threshold α = 5 /4 the existence of a globally regular solution still holds when the right-hand side of the first equation in (1.3) is replaced by a logarithmically supercritical operator. For α ∈ [ 3 /4, 1) and (1, 5 /4) partial regularity resultsà la Caffarelli-Kohn-Nirenberg were established in [25,44] and [8]. These works show the existence of a weak solution whose putative singular set (in space-time) has (5 − 4α)-dimensional Hausdorff measure 0. In the opposite direction, the recent works [7,12] prove the non-uniqueness of Leray-weak solutions to (1.3) in the parameter ranges α < 1 /5, respectively α < 1 /3. The non-uniqueness of weak solutions in the sense of Definition 1.1 is also shown to hold for α < 1 /2.
We note that very recently, by adapting the arguments in [5], Luo and Titi [35] demonstrated the non-uniqueness of very weak solutions for (1.3) in the parameter range α ∈ (1, 5 /4). When compared to [35] the weak solutions constructed in this paper have the additional property that their set of singular times has Hausdorff dimension strictly less than 1. Together, the uniqueness result [33], the non-uniqueness results of [35], and of this work, confirm the well-posedness criticality of the exponent α = 5 /4, within the class of weak solutions defined in Definition 1.1.
We give the proof of Theorem 1.5 for general values of α < 5 /4. Theorem 1.1 follows by restricting to α = 1.

OUTLINE OF THE PROOF
The proof of Theorem 1.5 proceeds via a convex integration scheme based on the scheme introduced in [5], which is itself built on a long line of work initiated by De Lellis and Székelyhidi Jr. [10], culminating in the eventual resolution of Onsager's conjecture by Isett [19] (cf. [11,2,1,3,9,4,20]). Such a scheme is used to inductively define a sequence of approximate solutions, converging to a weak solution of (1.3). The principal new idea of this paper is to create good regions in time where the approximate solutions are strong solutions to (1.3) and are untouched in later inductive steps. This is achieved by employing the method of gluing introduced by Isett [19] (cf. [4]). Taking the countable union of the good regions over each inductive step, one forms a fractal set, whose complement has Hausdorff dimension strictly less than 1. This is explained in detail in Section 2.1 and 2.2 below. The concept of good regions is partially inspired by similar concepts introduced in [1] (cf. [3]). An additional novelty of the present work is the introduction of intermittent jets which replace the intermittent Beltrami flows of [5] as the fundamental building blocks on which the convex integration scheme in based (see Section 2.3 and 4.1).
2.1. Inductive estimates and main proposition. For every index q ∈ N we will construct a solution (v q ,R q ) to the Navier-Stokes-Reynolds system whereR q is a trace-free symmetric matrix. The pressure p q is normalized to have zero mean on T 3 and is explicitly given by the formula Here we use the convention that for a 2 tensor S = (S ij ) 3 i,j=1 the divergence contracts on the second component, i.e. (div S) i = ∂ j S ij . The summation convention on repeated indices is used throughout.
The size of the Reynolds stressR q will be measured in terms of a size parameter where λ q is a frequency parameter defined by where a ≫ 1 is an large real number to be chosen later. Note that δ 1 = λ β 1 = a βb is large if a is sufficiently large. For every q ≥ 0 we assume thatR q obeys the estimates for some ε R > 0 to be chosen later, which depends only on the values of α, β, and b. For the approximate velocity field v q , we assume that it obeys the estimates These inductive estimates will ensure that the approximate solutions v q converge strongly in C 0 (0, T ; L 2 ) to a weak solution v of the Navier-Stokes equation (1.3). Consider T > 0 and fix the parameter sequences {τ q } q≥0 and {ϑ q } q≥1 defined in (2.7) and (2.8) below, which obey the heuristic bounds ϑ q+1 ≪ τ q ≪ ϑ q ≪ 1 . (2.6) In particular, for q ≥ 1 we make the choices and For the special case q = 0 we set For ϑ 0 we do not need to assign a value.
In order to ensure that the singular set of times has Hausdorff dimension strictly less than 1, at every q ≥ 0 we split the interval [0, T ] into a closed good set G (q) and an open bad set which obey the following properties: is a finite union of disjoint open intervals of length 5τ q . 4 (iv) For q ≥ 1, the bad sets have measures which obey (2.10) (vi) The residual Reynolds stress obeys Due to (2.11) and the parabolic regularization of the Navier-Stokes equation (cf. (3.4) below) we have that v q is a C ∞ smooth exact solution of the Navier-Stokes equation on G (q) . In addition, (2.10) implies that v = v q on G (q) \ {0}, and thus the limiting solution v is C ∞ smooth on (G (q) \ {0}) × T 3 . This justifies that the singular set of times Σ T obeys (2.12) It thus follows from (2.9) and the definitions of τ q and ϑ q in (2.7) and (2.8) that (2.13) Here we have also used the definition of λ q , and the fact that b > 2. To estimate the box-counting (Minkowski) dimension of Σ T , we note that for every q ≥ 0, the set Σ T is covered by B (q) , which itself consists of disjoint intervals of length 5τ q . Due to (2.13), the number of such intervals is at most By (2.12), and the super-exponential growth of λ q , we conclude that (2.14) This implies that Σ T also has box-counting dimension strictly less than 1. 4 Observe that this condition is consistent with property (i) and the definition τ 0 = T /15.

Proposition 2.1 (Main Iteration Proposition).
There exists a sufficiently small parameter ε R = ε R (α, b, β) ∈ (0, 1) and a sufficiently large parameter a 0 = a 0 (α, b, β, ε R ) ≥ 1 such that for any a ≥ a 0 satisfying the technical condition (2.24), the following holds: Let (v q ,R q ) be a pair solving the Navier-Stokes-Reynolds system (2.1) in T 3 × [0, T ] satisfying the inductive estimates (2.4a)-(2.5b), and a corresponding set G (q) with the properties (i)-(vi) listed above. Then there exists a second pair (v q+1 ,R q+1 ) solving (2.1) and a set G (q+1) which satisfy (2.4a)-(2.5b) and (i)-(vi) with q replaced by q + 1. In addition we have that Gluing stage. The first stage of proving Proposition 2.1 is to start from the approximate solution (v q ,R q ) which obeys (2.4a)-(2.5b) and (2.11), and construct a new glued pair (v q ,R q ), which solves (2.1), obeys bounds which are the same as (2.4a)-(2.5b) up to a factor of 2, but which has the advantage thatR q ≡ 0 on T 3 × B (q+1) . Specifically, the new velocity field v q is defined as where the η i are certain cutoff functions with support in the intervals [t i , t i+1 + τ q+1 ] (with t i = ϑ q+1 i) that form a partition of unity, and the v i are exact solutions of the Navier-Stokes equation with initial datum given by v i (t i−1 ) = v q (t i−1 ). Due to parabolic regularization, these exact solutions v i are C ∞ smooth in space and time on the support of η i , so that v q inherits this C ∞ regularity. This is in contrast to (v q ,R q ), which is only assumed to be H 3 smooth. Trivially, in the regions where a cut-off η i is identically 1, v q is an exact solution to (1.3).
Observe that property (2.11) ensures that v q is already an exact solution of (1.3) on a large subset of [0, T ], namely the τ q neighborhood of G (q) . In particular if t i−1 and t i both lie within this neighborhood, then by uniqueness of the Navier-Stokes equation in Hence v q = v q is an exact solution here. In order to single out overlapping regions where v q is not necessarily an exact solution of (1.3) we introduce the index set We then define By the discussion above, it will follow that v q is an exact solution on the complement of B (q+1) , namely G (q+1) . We prove in Section 3 below that the above defined good and bad sets at level q + 1 obey the postulated properties (i)-(iv). In Section 3 we prove the following proposition:

18)
and moreover the velocity field v q satisfies: and the stress tensorR q satisfies: for all M, N ≥ 0.

Convex integration stage.
In this step we start from the pair (v q ,R q ), and construct a new pair (v q+1 ,R q+1 ) withR q+1 obeying (2.11) at level q + 1, and which obeys the bounds (2.4a)-(2.5b) at level q + 1. The perturbation w q+1 := v q+1 −v q will be constructed to correct forR q . Moreover, w q+1 will be designed to have support outside a τ q+1 neighborhood of G (q+1) -this ensures properties (v) and (vi) in Section 2.1 will be satisfied. As in [5], the perturbation w q+1 will consist of three parts: the principal part w (p) q+1 , the divergence corrector w (c) q+1 , and the temporal corrector w (c) q+1 . The principal part w (p) q+1 will be constructed as a sum of intermittent jets W (ξ) (defined in (4.4), Section 4.1). The use of intermittent jets replaces the use of intermittent Beltrami waves in [5]. The principal difference of intermittent jets from intermittent Beltrami waves is that their definition is in physical space rather than frequency space. Consequently, intermittent jets are comparatively simpler to define and they can be designed to have disjoint support, mimicking the advantageous support properties of Mikado flows, as introduced in [9]. We note that the intermittent variants of the d − 1 dimensional Mikado flows found in [37,36], lying in d-dimensional space, are insufficiently intermittent to be used as building block for a 3-D Navier-Stokes convex integration scheme. 5 Intermittent jets are inherently 3dimensional (in space), with the trade-off that they are time dependent. We note in passing that utilizing intermittent jets, it is likely that the convex integration results [37,38] on the transport equation may be improved.
In the definition of w (p) q+1 , the intermittent jets W (ξ) will be weighted by functions a (ξ) : w (2.21) for some large parameter µ. As is typical in convex integration schemes, the high frequency error can be ignored since its contribution toR q+1 can be bounded using the gain associated with solving the divergence equation. The temporal corrector w (t) q+1 is then defined to be where P H is the Helmholtz projection, and P =0 is the projection onto functions with mean zero. That is, is divergence free.
The intermittent jets will be defined to have support confined to ∼ (ℓ ⊥ λ q+1 ) 3 many cylinders of diameter ∼ 1 λq+1 and length ∼ ℓ ℓ ⊥ λq+1 . In particular, the support of w should be roughly the size R q 1 /2 L 1 , by Hölder one would expect an L p estimate on w Indeed we will prove estimate (2.22) for p = 2 and prove a slightly weaker estimate for 1 < p < 2 (see Proposition 4.4). Utilizing (2.22), one may heuristically estimate the contribution of (−∆) α w (p) q+1 to the new Reynolds stressR q+1 with p > 1 arbitrarily close to 1. Here we see the necessity of the 3-dimensionality of the intermittent jets. In order to ensure that an identity of the form (2.21) holds, the cylinder supports of the intermittent jets will be shifting at a speed ℓ ⊥ λ q+1 µ. Heuristically, one would then expect that in order to ensure that the contribution of ∂ t w (p) q+1 toR q+1 is small, one would need to impose an upper bound on the choice of µ. One then needs to choose µ carefully in order to balance different contributions to the Reynolds stress error. Explicitly, we will define the parameters µ, ℓ ⊥ and ℓ by . (2.23) With these choices, we have . For technical reasons, we will require that λ q+1 ℓ ⊥ ∈ N. This may be achieved by assuming that where we recall that we have previously assumed that b ∈ N.

2.4.
Proof of Theorem 1.5. Let u (1) and u (2) be as in the statement of the theorem. Let η : where a⊗b denotes the traceless part of the tensor a ⊗ b, and R is a standard inverse divergence operator acting on vector fields v which have zero mean on T 3 as for k, ℓ ∈ {1, 2, 3}. The above inverse divergence operator has the property that Rv(x) is a symmetric trace-free matrix for each x ∈ T 3 , and R is an right inverse of the div operator, i.e. div (Rv) = v. When v does not obeý T 3 vdx = 0, we overload notation and denote Rv := R(v −´T 3 vdx). Note that that ∇R is a Calderón-Zygmund operator, and that R obeys the same elliptic regularity estimates as |∇| −1 .
For q ≥ 1 we inductively apply Proposition 2.1. The bound (2.5b) and (2.15) and interpolation yields Since R q L 1 → 0 as q → ∞, it is straightforward to show that v is a weak solution of the Navier-Stokes equation.
Moreover, as a consequence of properties (i) and (v) from Section 2.1 and the definition of v 0 we have The argument leading to (2.14) implies that the singular set of times of v has box-counting dimension (and hence Hausdorff dimension) less than εR /64. Finally, the claimed C 0 t W 1,1+β ′′ x regularity on v, for some β ′′ > 0, follows from the maximal regularity of the heat equation (fractional heat equation if α > 1), once we note that is chosen suitably small. The theorem then holds withβ = min{β ′′ , β ′ , εR /64} > 0.

GLUING STEP
3.1. Local in time estimates. It is well-known Navier-Stokes equations are locally (in time) well-posed in H 3 , which is a scaling subcritical space. Moreover, away from the initial time, parabolic regularization takes place. We summarize these facts, in version that is suitable for the applications in this paper.
have zero mean on T 3 , and consider the Cauchy problem for (1.3) with this initial condition. There exists a universal constant c ∈ (0, 1] such that if t 1 > t 0 is such that then there exists a unique strong solution to (1.3) on [t 0 , t 1 ), and it obeys the estimates Moreover, assuming that From the Gagliardo-Nirenberg-Sobolev and the Poincaré inequalities, and using that ∇ · v = 0 we obtain x a priori estimate is sufficient to establish the uniqueness of the solution. The higher regularity claimed in (3.4) follows from the mild form of the solution and properties of the fractional heat equation which may be derived from Plancherel. Let us first focus on the case M = 0. For α = 1, estimate (3.4) is well-known, and follows from the instantaneous gain of analyticity of the solution [14], or a small modification of the below argument. For α > 1 we briefly sketch the argument. Using Gallilean invariance, let us only consider the case t 0 = 0. From the inequality the formulation (3.5) and the boundedness of the Leray projector P H on L 2 , we obtain from which (3.4) with N = 1 and M = 0 follows in view of (3.3). In order to treat the case N ≥ 2 and M = 0, we first note that for 1 ≤ n ≤ N − 1 by induction on N we have Using the above estimate with n = N − 1 we obtain that To obtain the desired bounds for M = 1, let us consider the case N = 0 first. Using the equation, the already established bounds for M = 0 and N ≥ 0, the Gagliardo-Nirenberg-Sobolev inequalities, and the fact that the Leray projector is bounded on L 2 , we have that and the desired bound follows from the assumption (3.3). The remaining cases N ≥ 1 are treated in a similar manner, using the Leibniz rule. We omit these details.

Stability estimates.
In this section we estimate the difference between an approximate solution v q and an exact solution of the Navier-Stokes equation. Let R be the inverse divergence operator defined in (2.26). The main result is: ) and an integrability index p 0 ∈ (1, 5 /4). Assuming the parameter δ 0 is sufficiently large, depending on p 0 , the following holds. For Then, in view of (2.5b) and Proposition 3.1 there exists a unique In particular, letting from the bounds (3.8a)-(3.8b) we obtain the following stability estimate: Proof of Corollary 3.3. We show that estimates (3.8a)-(3.8b) imply the bounds (3.10a)-(3.10b). Recall the stressR q has zero mean. For p ∈ (1, 2] and δ ∈ [0, 1] by interpolation we have the inequalities The implicit constant depends only on p and δ. In order to prove (3.10a), we use (3.8a) and apply estimate (3.11) with δ = 1 and p = 2. We obtain from (2.4a)- from which estimate (3.10a) follows, since δ q+1 ≤ λ β 1 , and β is sufficiently small. The leftover power of λ q may be used to absorb any constants.
Similarly, in order to prove (3.10b), we use (3.8b), the bound (3.11) with δ = 0 and p = p 0 , and the embedding In the last inequality above we have used the definitions of δ q+1 and λ q , and the fact that p 0 ≥ 1. Estimate (3.10b) follows from the assumption (3.9) on p 0 , upon using the leftover power of λ q to absorb the implicit constants.
Proof of Proposition 3.2. For simplicity, by temporal translation invariance it is sufficient to consider the case t 0 = 0.
In order to prove (3.8a) we let u = v q − v and q = p q − p. Then div u = 0, u| t=0 = 0, and u obeys the equation where P is the Leray projector. Then, since u(0) = 0 the solution of (3.12) may be written in integral form as Next, we use that that for p ∈ [1, 2], t > 0, and any periodic function φ of zero mean we have that where the implicit constant only depends on α. These estimates follow from L 1 bounds for the Green's function of the fractional heat equation. We will also frequently use the Gagliardo-Nirenberg estimates which hold for zero-mean periodic functions φ. We return to (3.13) and obtain that for a suitable constant C 1 > 0 which only depends on p 0 , since p ∈ [p 0 , 2] and α ∈ [1, 5 /4]. Next, we claim that if t 1 > 0 is chosen sufficiently small, depending on v L ∞ and v q L ∞ , then we have Thus if we ensure that for some universal constant C 1 > 0, and further, using (2.5a), (2.5b), (3.2a) and (3.2b), we obtain that To conclude, we use (3.7), which shows that the left side of (3.19) is bounded from above by by letting a, and hence δ 0 , be sufficiently large. Here we have used that α ≥ 1, and that δ 0 , λ q ≥ 1. Thus, we have shown that (3.17) holds.
In order to prove (3.8b) we denote Note that since div u = 0 we have curl z = −u, and using the Calderón-Zygmund inequality we have Ru(t) L p z(t) L p . Thus our goal is to obtain L p estimates for z(t). We apply ∆ −1 curl to the equation obeyed by u (it is convenient to rewrite (3.12) without Leray projectors, and add a pressure gradient term, which is then annihilated by the curl operator) and obtain For the last term on the right side of in (3.20) we have used the identity which written for the i th component is Here we used that the transposition of two indices in ǫ jkl results in a (−1) factor. Moreover, we have also spelled out the commutator term on the right side of (3.20) as which written for the i th component is Here we have also used that ǫ ijk = 0 if two of the indices i, j, or k repeat, and that ǫ ijk ǫ nlk = δ in δ jl − δ il δ jn , where the δ's refer to the Kronecker symbol.
Using (3.20), upon placing the v · ∇z = div (v ⊗ z) term on the right side, and using that z(t 0 ) = 0, the solution to (3.20) may be written in integral form as (3.21) 13 From (3.14a)-(3.14b) and the boundedness of Calderón-Zygmund operators on L p , similarly to (3.16) we conclude that where C 1 depends only on p 0 and α, since p ∈ [p 0 , α]. Next we claim that if t 1 is chosen sufficiently small, then The argument is similar to the one for the bound for u(t), so we only sketch the details. Let us assume that (3.23) holds. Then from (3.22) we obtain Therefore, if we ensure that t 1 is small enough so that which shows that the bootstrap assumption was justified, and thus (3.23) holds on [0, T ]. Denote by C 1 the universal constant in the Gagliardo-Nirenberg inequalities (3.15a)-(3.15b). By also appealing to (2.5a)-(2.5b), (3.2a)-(3.2b), and our assumption (3.7) for t 1 , we obtain that the left side of (3.25) is bounded from above by once we ensure that a, and hence δ 0 is sufficiently large. This concludes the proof of (3.23).

3.3.
Proof of Proposition 2.2. We first define a C ∞ smooth partition of unity Denoting t i = ϑ q+1 i , this may be achieved by letting η i also have the following properties: where the implicit constant is independent of τ q+1 ϑ q+1 , and i.
As a consequence of the above properties, we have that η i η j = 0 whenever |i − j| > 1, and Having constructed the partition of unity {η i } nq+1 i=0 , we next construct exact solutions v i of the Navier-Stokes equation for suitably defined datum.
For every 1 ≤ i ≤ n q+1 we define v i (x, t) to be the unique smooth solution of the Cauchy problem for Navier-Stokes equation (1.3) with initial condition equal to v q at t i−1 : In view of (2.5a)-(2.5b), and Proposition 3.1 this solution v i is uniquely defined and obeys the estimates where c ∈ (0, 1) is the universal constant from (3.3), and α ≥ 1. Note that the definitions (2.3), (2.7), and the fact that β ≤ 1, imply that Therefore, assuming that δ 0 = λ 3β 1 λ −2β 0 ≥ λ β 0 is sufficiently large, depending on the universal constant c, by (3.31) we have that 3ϑ q+1 ≤ c which is consistent with (3.30). Therefore for all 1 ≤ i ≤ n q+1 the exact solutions v i (x, t) are smooth and well-defined for all t ∈ (t i−1 , t i+2 ] ⊃ sup(η i ). Moreover, since where the implicit constant depends only on N ≥ 0 and M ∈ {0, 1}. At this stage we glue the solutions v i together in order to construct (v q ,R q ). We define the divergence-free (note that the cutoffs η i are only functions of time) velocity and the interpolated pressure as where p i is the pressure associated to the exact solution v i . Also we let Having definedv q , we next prove that (2.18) holds. For t ∈ G (0) , this holds by construction. In view of (3.26), it suffices to show that if for some i ∈ {1, . . . , n q+1 } we have t ∈ supp (η i ) ∩ G (q) , then v i (t) = v q (t). For 15 this purpose recall by (2.5b) and (2.11) that v q is a strong solution of the Navier-Stokes equation for all t such that dist(t, G (q) ) ≤ τ q . Moreover, v i solves the Cauchy problem (3.28), so by the uniqueness of solutions in C 0 t H 3 x of the Navier-Stokes equation, we only need to ensure that dist(t i−1 , G (q) ) ≤ τ q . This follows from the fact that t ∈ G (q) , and 0 < t − t i−1 ≤ 3ϑ q+1 ≤ τ q . The last inequality trivially holds by (2.7) and (2.8) for q ≥ 1, and by taking a sufficiently large for q = 0. Thus, we have proven that (2.18) holds.
We now derive the formula for suppR q . Note that on [0, we have thatv q = v q is a smooth solution of the Navier-Stokes equation, and hence automaticallẙ We observe that v i − v i−1 has zero mean because the exact solutions of the Navier-Stokes equations v i , v i−1 preserve their average in time, and v q has zero mean by assumption. Hence we can apply the inverse divergence operator R to v i − v i−1 and for i ∈ {2, . . . , n q+1 − 1} define the symmetric traceless 2-tensor where we denote by a⊗b the traceless part of the tensor a ⊗ b. We also define the scalar pressure It follows from (3.35) that the pair (v q ,R q ) defined by (3.33a) and (3.36) solves the Navier-Stokes-Reynolds system (2.1) on [0, T ] with associated pressurep q . Next, we prove that (2.20a) holds. Note that by construction we have η i ≡ 1 on [t i + τ q+1 , t i+1 ] for all i ∈ {0, . . . , n q+1 }, and thus on these sets we have ∂ t η i = η i (1 − η i ) = 0. Therefore, by (3.36) we have thatR q (t) = 0 whenever t ∈ [t i + τ q+1 , t i+1 ] for some i. Thus it suffices to consider sets of times of the form (t i , t i + τ q+1 ). If i ∈ C or i − 1 ∈ C, then there is nothing to prove since by definition (2.17), dist((t i , t i + τ q+1 ), G (q+1) ) > 2τ q+1 . 16 Hence, consider the case i, i − 1 / ∈ C. Thus by the definition of C,R q (t) = 0 for all t ∈ [t i−2 , t t+1 + τ q+1 ]. Since , then sinceR q vanishes on [t i−2 , t t+1 + τ q+1 ], it follows by the bounds (3.29b) and (2.5b) and the uniqueness of strong solutions to the Navier-Stokes equations that v i−1 = v i = v q on (t i , t i + τ q+1 ). Thus by (3.36) we have thatR q+1 (t) = 0 for (t i , t i + τ q+1 ).
Here we note that condition (3.7) is satisfied on supp (η i ) due to (3.30). By (3.10a), we obtain that Since at most two terms appear in ( Here we have used that by τ q+1 ≤ ϑ q+1 , the definition (2.7), and the fact that ε R ≤ 1, to conclude and using the leftover term λ −1 q to absorb any implicit constants in (3.39). Combined, (3.38) and (3.39) prove (2.20b).
By an abuse of notation, let us periodize Φ ℓ ⊥ and ψ ℓ so that the functions are treated as functions defined on T 2 and T respectively. For a large real number λ such that λℓ ⊥ ∈ N, we define V ξ,ℓ ⊥ ,ℓ ,λ,µ : where here α ξ ∈ R 3 are shifts that ensure that the set of functions {V ξ,ℓ ⊥ ,ℓ ,λ,µ } ξ have mutually disjoint support. In order for such shifts α ξ to exist, we require that ℓ ⊥ to be sufficiently small, depending on the finite sets Λ α .
Our intermittent jet is then defined to be From the definition, using (4.1) and ℓ ⊥ λ ∈ N, we have that W (ξ) has zero mean, and Moreover, by our choice of α ξ , we have that the W (ξ) have mutually disjoint support, i.e.
Note that the intermittent jets W (ξ) are not divergence free, however assuming ℓ ⊥ ≪ ℓ then they be corrected by a small term, such that the sum with the corrector is divergence free. To see this, let us adopt the shorthand notation and compute . (4.6) Moreover, so long as ℓ ⊥ ≪ ℓ then W (c) (ξ) is comparatively small compared to W (ξ) . Observe that as a consequence of the normalizations (4.2) and (4.3) we have We also note that by definition W (ξ) is mean zero. As a consequence, using Lemma 4.1 we have for every symmetric matrix R satisfying |R − Id| ≤ 1 /2. By scaling and Fubini, we have the estimates where again here we have assumed ℓ −1 ≪ ℓ −1 ⊥ ≪ λ . Finally, we note the essential identity which follows from the fact that by construction we have that W (ξ) is a scalar multiple of ξ, and φ (ξ) is time-independent.
4.2. The perturbation. In this section we will construct the perturbation w q+1 .

Stress cutoffs.
Because the Reynolds stressR q is not spatially homogenous, we introduce stress cutoff functions. We let 0 ≤ χ 0 , χ ≤ 1 be bump functions adapted to the intervals [0, 4] and [1/4, 4], such that together they form a partition of unity: for any y > 0. We then define for all i ≥ 0. Here and throughout the paper we use the notation A = (1 + |A| 2 ) 1/2 where |A| denotes the Euclidean norm of the matrix A. By definition the cutoffs χ (i) form a partition of unity i≥0 χ 2 (i) ≡ 1 (4.14) and we will show in Lemma 4.2 below that there exists an index i max = i max (q), such that χ (i) ≡ 0 for all i > i max , and moreover that 4 imax τ −1 q+1 .

The definition of the velocity increment.
Recall that from Lemma 4.1, the functions γ (ξ) are well-defined and smooth in the 1 /2 neighborhood of the identity matrix. In view of (4.13), this motivates introducing the parameters ρ i by which have the property that on the support of χ (i) , for all i ≥ 0 .

4.2.3.
Estimates of the perturbation. This section closely mirrors [5,Section 4.4], and thus we omit most details where the estimates/proofs are mutatis-mutandis those from [5]. There is an analogy between the mollification parameter ℓ in [5] and the time-scale τ q+1 in this paper, in view of parabolic smoothing. First, similarly to [5, Lemma 4.1 and 4.2] we state a useful lemma concerning the cutoffs χ (i) defined in (4.13), summarizing their size and regularity: 21 For q ≥ 0, there exists i max (q) ≥ 0 such that Moreover, for all 0 ≤ i ≤ i max we have that 29) and the bound holds. Additionally, for 0 ≤ i ≤ i max we have The bound (4.33) follows from (2.20b)-(2.20c) and the Gagliardo-Nirenberg inequality H 3 , which holds for any zero-mean periodic function f ∈ H 3 , and the definition of τ q+1 in (2.8), and the fact that ε R ≤ 1.
then χ (i) (x, t) = 0. Therefore, by the bound (4.33) and the fact that βb ≤ 1 /4, if i ≥ 1 is such that then χ (i) ≡ 0. Therefore, in view of the parameter inequality which holds in view of (2.8) and the fact that βb ≤ 1 /4, upon taking a to be sufficiently large, we may thus define With this choice of i max it from the above argument it follows that (4.28) holds. The bound on i max claimed in the second inequality in (4.29) then follows from the above definition. The bound (4.30) follows from the second estimate in (4.29) which gives an upper bound on i max , the definition (4.15), and using that λ − ε R/16 q log 4 (τ −2 q+1 ) ≤ 8λ − ε R/16 q log 4 (λ q ) can be made arbitrarily small if a is chosen sufficiently large, depending on ε R .

Lemma 4.3. The bounds
hold for all 0 ≤ i ≤ i max and N ≥ 1.   (4.39) for N ∈ N, and p > 1.
From the second estimate in (4.39) it is clear that the incompressibility and temporal correctors obey better estimates than the principal corrector.
In order to establish the bound (4.37), it is essential to use the fact that a (ξ) oscillates at a frequency which is much smaller than that of W (ξ) , which allows us to appeal to the L p de-correlation lemma [5,Lemma 3.6], which we recall here for convenience: Let p ∈ {1, 2}, and let f be a T 3 -periodic function such that there exists a constants C f such that for all 1 ≤ j ≤ M + 4. In addition, let g be a ( T /κ) 3 -periodic function. Then we have that holds, where the implicit constant is universal.
Proof of Proposition 4.4. In order to prove (4.37), we use (4.34) when N = 0, and (4.36) with λ where the implicit constant depends only on N . Since W (ξ) is ( T /λq+1ℓ ⊥ ) 3 periodic, in order to apply Lemma 4.5 with λ = τ −5 q+1 and κ = λ q+1 ℓ ⊥ , we first note that by (2.8) and (2.23), we have , (4.42) by using that β is sufficiently small and b is sufficiently large, depending on α. For instance, we may take In (4.42) we also used that λ − 15β 2 1 2π √ 3 ≤ 1, once a is chosen sufficiently large. Therefore, after a short computation we see that the assumptions of Lemma 4.5 hold with the aforementioned κ and λ, with M = 4 in (4.40). Therefore, we only care about N ≤ 4 in (4.41), which also fixes the implicit constant in this inequality, and we may thus take C f to be proportional to ρ by using the small negative power of λ q to absorb the implicit constants in the first inequality. Consider the estimate (4.38). Observe that by definition (4.21), estimate (4.10) with M = 0, and the bound (4.36), we have Here we have again used (4.29) in order to sum over i, and have used the bound τ −3 q+1 ≤ λ q+1 which holds since β is small and b is large.
For the analogous bound on w (c) q+1 , by (4.8)-(4.10), estimate (4.36), the parameter estimates τ −3 q+1 ≤ λ q+1 and ℓ ⊥ ≤ ℓ , and Fubini (recall that ψ (ξ) and Φ (ξ) are functions of one and respectively two variables which are orthogonal to each other), we have Summing over 0 ≤ i ≤ i max loses an additional factor of τ −1 q+1 , which yields the desired bound for the first term on the left of (4.39). Similarly, to estimate the summands in the definition (4.24) of w (t) q+1 we use (4.8), (4.9), (4.36), the aforementioned parameter inequalities, and Fubini to obtain Summing in i loses a factor of τ −1 q+1 cf. (4.29), and we obtain the bound for the second term on the left of (4.39). For the proof of (4.39), we additionally note that (4.20), and (4.43) imply the parameter inequalities which concludes the proof of the proposition.
The following bound shows that (2.15) holds, and collects a number of useful bounds for the cumulative velocity increment w q+1 , which in turn imply that (2.5a) and (2.5b) hold at level q + 1.
Before turning to the proof of the proposition, we note that estimate (4.47) and the inductive assumption (2.5a) at level q imply that v q+1 L 2 ≤ 2δ which is a consequence of 2λ β q ≤ λ β q+1 . Thus, (2.5a) holds at level q + 1. Similarly, from (2.19b) and (4.48) with s = 3 and p = 2, and the parameter inequality (4.45), we conclude where we have used that b is large and that α ∈ [1, 5 /4). The remaining power of λ − 1 /2 q+1 may be used to absorb the implicit constant, and thus (2.5b) holds also at level q + 1.
Similarly to (4.51), we establish two bounds which will be useful in Section 5 for the proof of Corollary 5.2. First, from (4.48) with s = 9 /2 and p = 2, and (2.19d) with M = 0 and N = 2, it follows that Here we have also used the parameter inequality (4.45). Similarly, by (4.49) and the bound (2.19d) with M = 1 and N = 0 we obtain where in the last inequality we used that (4.45) provides an upper bound for τ −1 q+1 , and that α < 5 /4. In order to estimate ∂ t w (t) q+1 we use estimates (4.8) and (4.9), Fubini, and the bound (4.36) to obtain Here we have used explicitly the parameter choice (2.23), the parameter inequality (4.20), the first bound in (4.45), the bound ℓ −1 ≤ ℓ −1 ⊥ ≤ λ q+1 , and the inequality i max τ −1 q+1 .

CONVEX INTEGRATION STEP: THE REYNOLDS STRESS
The main result of this section may be summarized as: Proposition 5.1. There exists an ε R > 0 sufficiently small, and a parameter p > 1 sufficiently close to 1, depending only on α, b, and β, such that the following holds: There exists a traceless symmetric 2 tensor R and a scalar pressure field p, defined implicitly in (5.5) below, satisfying where the constant depends on the choice of p and ε R , but is independent of q, and R has the support property An immediate consequence of Proposition 5.1 is that the desired inductive estimates (2.4a) and (2.4b) and the support property (2.11) hold for the Reynolds stressR q+1 , which is defined as follows.
Corollary 5.2. There exists a traceless symmetric 2 tensorR q+1 and a scalar pressure field p q+1 such that the triple (v q+1 , p q+1 ,R q+1 ) solves the Navier-Stokes-Reynolds system (2.1) at level q + 1. Moreover, the following bounds hold R q+1 Proof of Corollary 5.2. With R and p defined in Proposition 5.1, we let It follows from (5.1) and the definitions of the inverse-divergence operator R and of the Helmholtz projection P H that the (v q+1 , p q+1 ,R q+1 ) solve the Navier-Stokes-Reynolds system (2.1) at level q + 1. Since the operator RP H div is time-independent, the claimed support property forR q+1 , namely (2.11) at level q + 1, follows directly from (5.3).
With the parameter p > 1 from Proposition 5.1, using that RP H div L p →L p 1, we directly bound The estimate (5.4a) then follows since the residual factor λ −εR q+1 can absorb any constant if we assume a is sufficiently large. In order to prove (5.4b), we use equation (5.1), the support property ofR q+1 which implies that suppR q+1 ⊂ T 3 × [ T /3, 2T /3], and the bounds (4.50)-(4.53). Combining these, we obtain For the dissipative term we have used that α < 5 /4, so that 2α + 2 < 9 /2. Using the residual power of λ − 1 /2 q+1 we may absorb any constants and thus (5.4b) follows. 5.1. Proof of Proposition 5.1. Recall that v q+1 = w q+1 + v q , where v q is defined in Section 3.3 and (v q ,R q ) solve (2.1). Using (4.25) we obtain =: div R linear + R corrector + R oscillation + ∇q . (5.5) Here, the linear error and corrector errors are defined by applying R to the first and respectively second line of (5.5), while the oscillation error is defined in Section 5.1.3 below. The zero mean pressure q is defined implicitly in a unique way.
Besides the already used inequalities between the parameters, ℓ ⊥ , ℓ and λ q+1 , we shall use that if p is sufficiently close to 1 the following bounds hold: To see this, we appeal to the bound (4.45) for τ −1 q+1 and the parameter choices (2.23) to conclude that the left side of (5.6) is bounded from above as where in the last inequality we have chosen p sufficiently close to 1, depending only on α. To conclude the proof of (5.6), note that and therefore if we ensure that ε R and β are sufficiently small, depending on α and b only, such that then the three estimates above imply (5.6).

ACKNOWLEDGMENTS
The work of T.B. has been partially supported by the NSF grant DMS-1600868. V.V. was partially supported by the NSF grant DMS-1652134 and by an Alfred P. Sloan Research Fellowship.