From Jean Leray to the millennium problem: the Navier–Stokes equations

One of the Millennium problems of Clay Mathematics Institute from 2000 concerns the regularity of solutions to the instationary Navier–Stokes equations in three dimensions. For an official problem description, we refer to an article by Charles Fefferman (Existence and smoothness of the Navier–Stokes equation. http://www.claymath.org/sites/default/files/navierstokes.pdf) on the whole space problem as well as the periodic setting in R3\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathbb {R}^3$$\end{document}. In this article, we will introduce the Navier–Stokes equations, describe their main mathematical problems, discuss several of the most important results, starting from 1934 with the seminal work by Jean Leray, and proceeding to very recent results on non-uniqueness and examples involving singularities. On this tour de force we will explain the open problem of regularity.


The model of Navier-Stokes
The Navier-Stokes equations are the standard model to describe the flow of a viscous, incompressible fluid such as water; however, polymeric fluids, blood and oil need more sophisticated models involving further nonlinear terms. Dating back to Newton, Daniel Bernoulli and Euler the laws of conservation of mass and of linear momentum yield for the velocity field u = (u 1 , u 2 , u 3 ) and the hydrostatic pressure p of a viscous and incompressible fluid the set of partial differential equations The only nonlinear term u · ∇u stands for 3 j=1 u j J. Evol. Equ. x(t), of a particle, i.e., the acceleration Based on Newton's law of viscosity and the so-called Stokes' hypothesis (1845) which was more or less already formulated by A.J.C. Barré de Saint-Venant in 1843 the easiest way to model viscosity is given by the term −νΔu + ∇ p where ν > 0 is the constant coefficient of viscosity. For simplicity, we put the constant density ρ to 1. In more complicated models of viscous fluids, the term −νΔu + ∇ p is replaced by nonlinear terms leading to non-Newtonian fluids or by either higher order or nonlocal terms with respect to u. From this point of view, the Navier-Stokes equations are the simplest model for viscous, incompressible fluids. The system (1), (2) is to be completed by an initial value u 0 for the velocity u at time t = 0, say, and boundary conditions for u on the boundary ∂Ω of a domain Ω ⊂ R 3 unless Ω = R 3 or Ω equals the torus T 3 = (R/(2π Z)) 3 with periodic boundary conditions. The most common boundary condition is given by the adhesion or no-slip condition u = 0 on ∂Ω, i.e., particles on ∂Ω at time t = 0 stay on the boundary for all times.

Basic observations
The Navier-Stokes system (1)-(2) together with initial-boundary value conditions define a coupled system of nonlinear partial differential equations so that-as for a nonlinear ordinary differential equation-the question on global-in-time existence and uniqueness of solutions might be nontrivial. Unfortunately, compared with the Poisson problem −Δu = f or the heat equation u t − Δu = f , there exists no maximum principle for the velocity field u. Given boundary data u = 0 on ∂Ω, we cannot expect to get the a priori estimate max x∈Ω |u(x)| ≤ max y∈∂Ω |u(y)| in the stationary case. Actually, considering the flow in a finite channel with a bottleneck between the inlet and outlet the velocity vector will attain its maximum somewhere in the bottleneck due to the condition div u = 0, the incompressibility of the fluid.
Moreover, the pressure plays a hindering but nevertheless crucial role. For the stationary Stokes problem, i.e., for the system − νΔu + ∇ p = 0, div u = 0 in Ω, u = 0 on ∂Ω, the solenoidality condition div u = 0 can be considered as a linear constraint requiring a Lagrange multiplier so that (3) can be solved. This multiplier is given by the pressure. Furthermore, in (1)- (2), there is no time derivative for p so that with a solution p(x, t) (and u(x, t)) also p(x, t) + q(t) for any function q(t) which may be even non-measurable defines a solution.
Usually the pressure will be eliminated by the so-called Helmholtz projection P which is given by the solution to a weak Poisson problem with Neumann boundary condition. The solution operator, a pseudodifferential operator of order 0, acts nonlocally in the domain Ω so that in the Navier-Stokes system the pressure depends nonlocally on the term u · ∇u. The Laplacian −Δ will be replaced by the Stokes operator A = −PΔ which partly has properties as the Laplacian −Δ, in particular in the whole space Ω = R 3 or the periodic setting T 3 ; but it is much more difficult to handle in domains with boundary since the operator P maps vector fields u to solenoidal vector fields with vanishing normal component (Pu) · n = 0 on ∂Ω so that in general −PΔ = −ΔP. In other words, the lack of the maximum principle for the Stokes problem (3) compared to the Laplace equation is related to the nonlocality of the Stokes operator A. 1. 3. The state of the art before 1934 There are only very few explicit solutions to the Navier-Stokes system. These examples are solutions in special geometries like an infinite tube (Hagen-Poiseuille flow) or in an infinite rotating cylinder with inner and outer cylinder (Taylor-Couette flow). The common ingredient is the fact that for these solutions the nonlinear term u · ∇u vanishes. The so-called Jeffery-Hamel flow in an infinite two-dimensional wedge with in-or outflow through the vertex of the wedge admits a radially symmetric solution in the sense that ru(x), r = |x|, solves a nonlinear second-order ordinary differential equation. Solutions to the stationary Stokes system (3) in more general domains, either bounded or exterior domains, were found by the method of singleand double-layer hydrodynamic potentials leading to boundary integral equations and exploiting Fredholm theory. The instationary Navier-Stokes system in the whole space R 3 was solved by Oseen [47, §7.4] using an ansatz by power series. However, the power series will converge only on a small time interval [0, T 0 ], and Oseen remarks [47, §7.5] that a further iterative process on a sequence of consecutive intervals [T 0 , T 1 ], [T 1 , T 2 ], …will probably stop in finite time. The main problem in his and other approaches is the fact that the authors were working in spaces of classical smooth functions for which they could not (and we cannot yet even today) find uniform a priori estimates globally in time. To be more precise, these spaces were not adapted to the underlying physics of the Navier-Stokes system.

Construction of weak solutions
Adequate a priori estimates are given by the energy (in-)equality obtained by testing (1) with the solution u itself, i.e., taking the L 2 (Ω)-scalar product of (1) with u and applying integration by parts. Since u| ∂Ω = 0 and div u = 0 in Ω, (1) yields the equation d dt where Ω ∇ p · u dx = 0 and Ω (u · ∇u) · u dx = 1 2 Ω u · ∇|u| 2 dx = 0. Hence, an integration in time on [s, t] yields the energy equality here, · 2 denotes the L 2 (Ω)-norm of functions, vector fields, etc. With s = 0 and i.e., the flow field u has a kinetic energy 1 2 u(t) 2 2 which is uniformly bounded for all times; the same holds for the dissipation energy, ν t 0 ∇u 2 2 dτ , which is caused by friction of particles close to each other, but with different velocities (strain rate tensor (∇u + (∇u) )(x, t) = 0). The idea to use energy estimates dates back to Leray [41]. He was the first to construct global-in-time weak solutions called solutions turbulentes.
This procedure motivates the definition of weak solutions which-for simplicitywill be stated for bounded domains Ω ⊂ R 3 only. In this definition, L 2 σ (Ω) = P L 2 (Ω) is the set of all weakly solenoidal vector fields v in L 2 (Ω) such that the normal component, v · n, vanishes on ∂Ω in a generalized sense, and C ∞ 0,σ (Ω) denotes the set of all solenoidal test vector fields in C ∞ 0 (Ω). Finally, H 1 0 (Ω) is the Sobolev space of L 2 -vector fields v such that v| ∂Ω = 0 and the weak gradient, ∇v, belongs to L 2 (Ω).

Definition 1. A solenoidal vector field u in the Leray-Hopf class
is called a weak solution in the sense of Leray and Hopf of the Navier-Stokes system (1), (2) with Dirichlet boundary data u| ∂Ω = 0 and initial value u 0 ∈ L 2 σ (Ω) if for all smooth compactly supported and solenoidal (with respect to x) vector fields ϕ and u(t) satisfies for all t ∈ [0, T ) the energy inequality Weak solutions of the Navier-Stokes system in the sense of Leray and Hopf can be constructed for bounded domains and all kinds of unbounded domains, even with rough boundaries. The classical procedure is by the Galerkin approximation in which the nonlinear PDE (1)-(2) is replaced by a sequence of nonlinear (N × N )-ODE systems of first order in time yielding a sequence of approximate solutions (u N ) N ∈N . Using the energy equality (4) for each u N , the sequence (u N ) is shown to be bounded in the Leray-Hopf class LH T . Then, by functional analytic principles, (u N ) possesses a weakly/weakly-( * )-convergent subsequence, called (u N ) again, such that Unfortunately, weak convergence does not allow to pass to the limit in the energy equality (4), but, due its lower semi-continuity with respect to norms, guarantees the energy inequality (7) for the weak limit u. Actually, it is still an open problem whether any weak solution u ∈ LH T of the variational problem (6) satisfies the energy equality (4) or the energy inequality (7). The crucial step in the passage to the limit u N → u is the nonlinear term u N · ∇u N ; for that step a compactness argument is needed to conclude from further properties of the sequence of approximate solutions (u N ) that u N → u in L 2 (0, T ; L 2 (Ω)). The Galerkin approximation method was exploited by Hopf [35] to construct global weak solutions of the Navier-Stokes equations in an arbitrary domain Ω ⊂ R 3 the shape of which may even depend on time. Leray [41] accomplished this analysis for the whole space case Ω = R 3 in which approximate solutions and global a priori estimates are easily obtained by the use of the heat kernel; however, his compactness argument on R 3 was quite involved.

Fundamental properties of weak solutions in R 2 and R 3
Mathematically, it makes sense to discuss the Navier-Stokes system (1)-(2) also in two dimensions. Indeed, a flow field u = (u 1 , u 2 , u 3 ) in a three-dimensional domain Ω = Ω × R ⊂ R 3 such that u 3 ≡ 0 and u 1 , u 2 , p are independent of the third variable x 3 , is a solution of (1)-(2) in the two-dimensional domain Ω ⊂ R 2 .
In the following, let us compare properties of weak solutions in R d , d = 2 or 3. To this aim, we introduce the Serrin number S(s, q) for the Bochner space L s (L q ) := L s (0, T ; L q (Ω)): Then, for u ∈ L ∞ (L 2 ) we have S u = d 2 . Moreover, by Sobolev's embedding theorems, u ∈ L 2 (H 1 0 ) ⊂ L 2 (L 2 * ) where 2 * = 2d d−2 so that S u = d 2 as above (actually, when d = 2, there holds only H 1 ⊂ L q for all 2 ≤ q < ∞). Thus, the Serrin number for weak solutions equals S u = d 2 . Returning to the energy (in-)equality let us test the Navier-Stokes solution with the weak solution u itself. The crucial term is the nonlinear one, u · ∇u, leading to the integral Although it is easily seen by an integration by parts that for almost all t ∈ (0, T ) the integral Ω (u · ∇u) · u dx is well-defined and vanishes, the double integral (8) requires in view of ∇u ∈ L 2 (L 2 ) more or less that u ∈ L 4 (L 4 ). For d = 2, S(4, 4) = 1 = d 2 so that the "necessary" condition u ∈ L 4 (L 4 ) is satisfied and consequently the integral in (8) vanishes. The final conclusion is that weak solutions in the 2d-setting satisfy the energy equality. However, in R 3 , S(4, 4) = 5 4 < 3 2 = d 2 , i.e., the condition u ∈ L 4 (L 4 ) cannot be guaranteed for weak solutions and the validity of the energy equality remains open. On the other hand, Galdi [24] proved that the condition u ∈ L 4 (L 4 ) yields the validity of the energy equality even for weak solutions in the sense of distributions, i.e., u satisfies (6) in Definition 1, but a priori is not necessarily an element of the Leray-Hopf class LH T ; in this case, the condition u ∈ LH T follows automatically.
There is an intrinsic argument why Serrin's condition on u is important. The Navier-Stokes equations have the remarkable property that with a solution u and an associated pressure p also the pair (u λ , p λ ), for any λ > 0 is a solution of (1)-(2) for the same coefficient of viscosity ν. This property guarantees that a "solution" controlled in a wind tunnel experiment for a scaled model is related-together with all its measurands-to a real world solution via the scaling (9). Mathematically, this property would open up the opportunity to work with u λ for sufficiently large or small λ rather than u itself and exploit the smallness of some adequate norm of u λ . However, it is easy to check that We conclude that the analysis in scale invariant function spaces yields in some sense optimal results which cannot be improved by replacing u with u λ .
Looking at uniqueness and regularity of weak solutions, we put the nonlinear term to the right-hand side as a given external force. Then, a rigorous analysis based on Sobolev embeddings and Gagliardo-Nirenberg inequalities shows that u · ∇u can be absorbed by the term u t − νΔu on the left-hand side in the two-dimensional case, but not in the three-dimensional one. A formal check is based on Serrin numbers: Since ∇u ∈ L 2 (L 2 ), we have S ∇u = d 2 + 1 so that with S u = d 2 and Hölder's inequality On the other hand, formally S u t = S Δu = d 2 + 2. In the two-dimensional case where d + 1 = d 2 + 2, the nonlinear term can be controlled by linear ones. Actually, a more careful analysis proves that weak solutions in R 2 are unique and smooth in the classical sense for all times t > 0. However, in the 3d case, where d + 1 > d 2 + 2, the nonlinear term cannot be absorbed by u t − νΔu and uniqueness as well as regularity remain open questions-unless the weak solution u satisfies further conditions. This topic will be discussed in the next section.

Basic results on regularity of weak solutions
Assume in the 3d case that a weak solution u lies in the so-called Serrin class, i.e., Then, u is unique among all weak solutions with the same data satisfying the energy inequality (weak-strong uniqueness, see [51]). The same condition, (10), implies that u is a regular solution. For instance, if Ω ⊂ R 3 is a "uniform C 2 -domain" and u 0 is regular, see (20), (21) in Sect. 2.5 for details, then for each 0 in particular, all derivatives in the Navier-Stokes system exist in the weak sense and the equations are satisfied in the pointwise sense almost everywhere in time and space. Even C ∞ -regularity up to the boundary ∂Ω can be obtained assuming that ∂Ω is of class C ∞ . Precise results and proofs are found in the monograph by H. Sohr [53, Theorems V. 1. 8.1 and V. 1.8.2], whereas the result on uniqueness dates back to Leray [41] and Serrin [50,51].
On the other hand, if u is only a solution in the sense of distributions lying in Serrin's class (10) for each interval (δ, T ), 0 < δ < T , and converges to u 0 weakly in L 2 as t → 0+, then u is actually a Leray-Hopf solution satisfying (5) in Definition 1 and even u(t) → u 0 in L 2 as t → 0, see Galdi [25]. A similar result holds in the limit case when L s (L q ) is replaced by C 0 (L 3 ).
The limit case s = ∞, q = 3 in (10) posed further problems concerning uniqueness as well as regularity since L ∞ is neither reflexive nor separable. The corresponding regularity problem was solved not before 2003 by Escauriaza et al. [18] for the whole space Ω = R 3 ; the main idea is a backward uniqueness property for the heat equation (inequality) for the vorticity ω = rot u: roughly spoken, if |ω t − νΔω| ≤ M(|ω| + |∇ω|) on R 3 ×(−T, 0), ω(t = 0) = 0 on R 3 and ω is restricted by some rough growth condition, then ω vanishes on R 3 × (−T, 0). Subsequent articles deal with the case of domains.
In addition to Serrin's condition S u = 1, there are various other conditions on u, ∇u, on specific components of either u or ∇u, on the vorticity ω = curl u, the pressure p and further physical or mathematical quantities which yield regularity and classical smoothness of weak solutions. Numerous results of this type can be found in the literature under the name conditional regularity; see Sect. 3.1 below for a discussion of selected results.
In particular, the vorticity, has interesting mathematical and physical properties which are different in R 2 and R 3 . For a planar flow, d = 2, the derivatives ∂ 3 and the component u 3 vanish so that only the third component ω 3 of the vorticity ω may be nonzero. In this latter case, the scalar vorticity ω(:= ω 3 ) solves the scalar vorticity transport equation a parabolic equation for ω which admits control of ω by the maximum and minimum principle since in each term ω appears under a differential operator. However, note that in a domain Ω ⊂ R 2 the Dirichlet boundary condition u = 0 on ∂Ω does not yield enough information for ω on the boundary to allow for a rigorous application of the maximum principle. In the three-dimensional case, the vorticity transport equation reads a vectorial equation with the additional term ω · ∇u, called the vortex stretching term. For instance, if locally, (13) ω 3 is likely to increase with respect to t and x 3 and could tend to infinity. For further details on regularity related to vorticity, we refer to Sect. 3.1.

Construction of regular solutions
Although there are several well-established methods to construct global-in-time weak solutions, the problem to prove that they are globally regular is still unsolvedit is one of the Millennium Prize problem.
Another idea is to construct regular solution from the very beginning. Concerning Galerkin's method, we can use smooth ansatz functions u j given by an orthonormal basis of eigenfunctions of the Stokes operator A = −PΔ in a bounded domain. Actually, similarly to the Laplacian on a bounded smooth domain, A has a compact resolvent, and, consequently, there exists a complete set of eigenfunctions u j in the Sobolev space The crucial point to construct a regular solution u on a given interval (0, T ) with the property (11) is to get adequate higher order a priori estimates of u. To this aim, the Navier-Stokes system (1) is tested with Au, rather than with u, and yields, amongst others, the integral This integral is estimated by Hölder's inequality, the interpolation estimate Vol. 21 (2021) From Jean Leray to the millennium problem 3251 and the Sobolev embedding H 1 0 (Ω) ⊂ L 6 (Ω) as follows: The term ν 2 Au 2 2 can be absorbed by the term ν Au 2 2 which appears on the left-hand side of the Navier-Stokes system. However, the function can only be shown to satisfy a differential inequality of the form y ≤ Cy 3 and stays bounded for These formal arguments can be made rigorous and prove that the Galerkin approximation method yields a strong solution at least on the interval [0, T 0 ), the length of which depends on norms of the initial value. By the way, in the planar case, the interpolation estimate (14) can be replaced by the much stronger estimate ∇u 4 ≤ c ∇u On the other hand, Leray [41] has shown that a weak solution (in the whole space case) getting singular as t → T 0 − in the sense that ∇u(t) 2 → ∞ blows up at least with the rate Moreover, since friction (dissipation) leads to a loss of kinetic energy, for any global weak solution u, there exists T ∞ > T 0 such that u is regular on (T ∞ , ∞). These two facts lead to the famous Structure Theorem of Leray. Similar results hold for weak solutions on bounded and exterior domains satisfying a so-called strong energy inequality. To be more precise, for a given global weak solution on [0, ∞) there exist finitely or at most countably many open intervals of regularity, I j , in (0, T ∞ ), such that u is regular on each interval I j as well as on (T ∞ , ∞), and the complement satisfies-in terms of the Hausdorff measure H 1/2 -the condition see, e.g. [23]. In particular, the set of temporal singularities is a Lebesgue set of measure 0, but could contain uncountably many epochs.

Construction of regular solutions by semigroups techniques
Another important method for the construction of regular solutions is based on analytic semigroup theory. It is well-known that the Stokes operator A = −PΔ generates a bounded analytic semigroup e −t A , t ≥ 0, for bounded and exterior domains as well as for layers, tubes and half spaces (and variants of them). Then, Duhamel's formula yields the nonlinear integral equation For the analysis of this equation, we need the maximal regularity property of the Stokes operator, i.e., the a priori estimate for solutions of the instationary Stokes system Under suitable assumptions on Ω, this estimate holds for all 1 < s, q < ∞. The first result of this type is due to Solonnikov [54] for s = q; for a more functional analytic proof see Geissert et al. [28]. Moreover, Noll et al. [46] proved that the Stokes operator has a bounded H ∞ -calculus which, among others, allows to characterize the domains of fractional powers of the Stokes operator; for the latter result, we refer also to Giga [29]. Based on (18) and fractional powers of A, a solution u of (17) can be found by Banach's fixed point theorem. However, in order to get a contraction and self-map on a closed ball in some Banach space smallness assumptions either on the initial value u 0 or on the length of the interval of existence [0, T ) must be posed. For a modern review of the semigroup approach, we refer to a handbook article by Hieber and Saal [34]. Since this method is strongly based on L q -theory, q = 2, a careful analysis of the Helmholtz projection and the Stokes operator on L q -spaces is indispensable. However, due to a famous counter-example by Bogovskiǐ [60] there exist smooth unbounded domains on which the Helmholtz projection fails to satisfy the usual mapping properties of existence and uniqueness. In such a case, the Stokes operator A = −PΔ is not welldefined. On the other hand, Hilbert space methods, e.g., the Lemma of Lax-Milgram, allow for any bounded and unbounded domain, even with rough boundary, to define the operators P and A on L 2 -spaces, see [53]. A combination of local L q -theory and global-in-space L 2 -theory was the key to the analysis of P and A in general smooth unbounded domains in spaces of typẽ see Kozono, Sohr and the author [19,20]. Another interesting observation called "Weak Neumann implies Stokes" is described by Geissert et al. in [27]: For a domain with possibly non-compact smooth boundary, assume that the Helmholtz projection P exists on some L q space, 1 < q < ∞. Then, the Stokes operator generates an analytic semigroup on L q satisfying maximal regularity; hence, a unique local mild solution to the Navier-Stokes equations can be constructed provided q > d.
Looking for strong solutions u in Serrin's class L s (0, T ; L q ) where 2 s + 3 q = 1, 3 < q < ∞, also the initial value u 0 must have more regularity than L 2 σ (Ω). For the Stokes system (19) with f = 0, but u(0) = u 0 = 0, it is obvious that Serrin's condition for its solution u(t) = e −t A u 0 is equivalent to the property which in terms of real interpolation theory can be rewritten in the form where L

Conditional regularity
Besides the Serrin condition (10) for the velocity field u, there are numerous related conditions on components of u or ∇u, etc. A classical generalization by Beirão da Veiga [1] concerns ∇u such that S ∇u = 2. More technical proofs are necessary when only one or two components of u or ∇u are taken into account, and it is not immediate to decide whether the same Serrin condition applies to u 3 as for u, say. Most of the following results are proven for the case Ω = R 3 , but we are not touching on problems and details for domains Ω = R 3 .
In [44], Neustupa et al. proved that if one component of a (suitable) weak solution locally satisfies S u 3 = 1 2 , say, then the solution is locally regular. Recently, an optimal result for one column, ∂ 3 u, of ∇u was found by Z. Skalák [52]: If S ∂ 3 u = 2, then the solution is regular. Looking for the row ∇u 3 of ∇u, the L 2 -case ∇u 3 ∈ L 4 (L 2 ) where S ∇u 3 = 2 was solved by Han et al. [33], together with further results on fractional derivatives. Considering only one component of the tensor ∇u results are less complete in the sense that S ∂ i u j is assumed to be less than 2. Here, for diagonal element ∂ 3 u 3 results are better than for off-diagonal elements. We refer to Cao and Titi [9] for both types of results which are too technical for this short review. If mixed norms are used, i.e., with different integration exponents in x 1 -, x 2 -and x 3 -direction, then the critical Serrin number S = 2 can be approached for diagonal elements ∂ i u i , but not for off-diagonal elements so far; see Guo et al. [32].
Concerning the vorticity ω, it is not the size of ω itself which is responsible for the possible existence of singularities as proposed by the vortex stretching mechanism. Indeed, Constantin and Fefferman [15] found that the crucial term is the angle θ(x, y, t) = (ω(x, t), ω(y, t)) which should not vary too fast for large |ω|; if | sin θ(x, y, t)| ≤ K |x − y| 1/2 for large vorticities |ω(x, t)|, |ω(y, t)| ≥ M, then the solution will not blow up. A substantial improvement was found by Beirão da Veiga and Berselli [2] in the sense that it suffices to get for large ω( Returning to the 2d case in which ω =(0, 0, ω) T , the angle θ(x, y, t) is constant to either 0 or ±π so that the above condition is satisfied everywhere. Actually, if two components of the vorticity in the 3d case can be controlled, then the 3d-situation is close to the 2d one yielding regularity; however, it is an open problem whether the control of only one component of ω is sufficient to get regularity.
For the pressure p, there are less criteria available in the literature. Recall that p is unique only up to functions depending on t so that assumptions like p ≥ 0 or p ≥ −C are meaningless. However, Seregin and Šverák [49] proved that a weak solution is regular provided there exists a function g ≥ 0 such that p(x, t) ≥ −g(x, t) and g satisfies the condition |x−x 0 | dx is assumed to be continuous at t 0 from the left. The same result holds when the total head pressure 1 2 |v| 2 + p is bounded above by g. The condition above generalizes the pointwise condition p ≥ −C and prevents the fluid from developing the phenomenon of cavitation.
A completely different criterion is based on the geometry of the rate of deformation tensor, the symmetrized velocity gradient D(u) = 1 2 ∇u + (∇u) . Since D is symmetric with vanishing trace due to the solenoidality of u, D has three real eigenvalues λ 1 ≤ λ 2 ≤ λ 3 summing up to 0. Obviously, λ 1 ≤ 0 and λ 3 ≥ 0 without any a priori knowledge about the sign of λ 2 . Given an orthonormal basis of R 3 of eigenvectors of D, {e 1 , e 2 , e 3 }, it can be shown from (2) that an infinitely small test volume of the fluid is compressed or stretched in the direction e i if λ i < 0 or λ i > 0, respectively. Hence, in the case when λ 2 ≤ 0 a test cube will be compressed to a thin but long tube, whereas in the second case (λ 2 > 0) the test cube will be pressed flat to a thin sheet. Then, if the positive part of λ 2 as a function of (x, t) can be controlled locally by some integrability condition, the solution is locally regular. This result is due to Neustupa and Penel [45].

Analysis in scale invariant borderline spaces
The classical borderline case C 0 [0, T ]; L 3 (Ω)) ⊂ L ∞ (0, T ; L 3 (Ω)) for strong solutions can neither be reached by the Galerkin method nor the semigroup approach in Sect. 2.5 directly; one reason is the elementary fact that the norm in L ∞ (0, T ) does not necessarily converge to 0 for a bounded function as T → 0. This problem was solved for the first time by Kato [37] by a more involved iteration scheme than the classical fixed point iteration of Banach. The idea is to start a first iteration for (u j ) j∈N based on (17) such that Here, q > 3 is any exponent. Proving the convergence of this scheme, in a second step the convergence of (u j ) in C([0, T ; L 3 (Ω)) and the property are proved. The proof by Kato was designed for the whole space R 3 , but can be adapted to many other domains. Moreover, if u 0 ∈ L 3 (Ω) is sufficiently small, then T = ∞ is possible. A variant of Kato's method was described by Giga [30] where the gradient term in the iteration is no longer needed.
The largest space of initial values allowing for local and/or global in time solutions in scale invariant function spaces was found so far by Koch and Tataru [39]. Assume ∞,∞ is larger thanḂ −1 ∞,∞ , initial values in B −1 ∞,∞ will not be considered, but rather regularity question. E.g., Cheskidov and Shvydkoy [13] proved that a weak Leray-Hopf solution u ∈ C 0 ((0, T ]; B −1 ∞,∞ ) or with sufficiently small jumps in the B −1 ∞,∞ -norm is regular on (0, T ]. On the other hand, weak solutions of Leray-Hopf type may be discontinuous at t = 0 in the topology of B −1 ∞,∞ . Actually, for the torus T 3 , Cheskidov and Shvydkoy [12] construct an initial value u 0 ∈ 1<q≤∞ B 3/q−1 q,∞ ∩ L 2 σ such that each Leray-Hopf type weak solution u with u(0) = u 0 satisfies the discontinuity estimate lim sup here, δ > 0 is independent of u.
The Navier-Stokes system is ill-posed in the spaceḂ −1 ∞,∞ , see Bourgain and Pavlović [3]. Actually, the authors proved that for any δ > 0 there exists an initial value u 0 in the Schwartz class S(R 3 ), a T ∈ (0, δ) and a corresponding solution (u, p) to the Navier-Stokes system such that This result was extended by Yoneda [59] to the spacesḂ −1 ∞,r for r > 2, and by Wang [57] to the homogeneous Besov spacesḂ  x 2 ) and x 3 ∈ R. Note that due to the factor 1 ε the norm u 0,ε Ḃ −1 ∞,∞ increases like 1 ε as ε → 0, see [11, p. 989]. Further examples of highly oscillating initial values For a recent survey on the analysis in critical function spaces, we refer to Gallager [26].

Local regularity and possible singularities
In contrast to global-in-space regularity results for weak solutions u looking for singular epochs as in Leray's Structure Theorem, see (16), local regularity theory is interested in space-time points z 0 = (x 0 , t 0 ) ∈ R 4 such that u / ∈ L ∞ (U δ (z 0 )) for any neighborhood U δ ⊂ R 4 of z 0 . Due to the local character, the weak solution has to satisfy a localized version of the usual energy inequality and is called a suitable weak solution. The proof of existence of suitable weak solutions needs a more careful construction, in particular, the use of cut-off functions leads to pressure terms which must be estimated. Suitable weak solutions are known to exist for bounded, exterior and more general unbounded smooth domains, cf. [7,19]. The celebrated result of Caffarelli et al. [7] states that for a suitable weak solution u the set of singular spacetime points, S ⊂ R 4 , satisfies Due to the different scaling in time and space, the Hausdorff measure H 1 can be replaced by an even stronger parabolic variant, P 1 , of H 1 . This result implies that the set S cannot contain a line or any smooth path; thus, for axially symmetric flows, singularities can occur only on the axis of rotation. The result (22) is an integrated version of a local regularity analysis on space-time cylinders where the smallness parameter ε 0 is an absolute constant. We note that the above integrals are scale-invariant quantities with respect to the scaling (9). Moreover, the term ∇u can be replaced by ω = rot u and the L 3 (L 3 )−norm of u can be replaced by any scale-invariant local, weighted L s (L q )-norm (3 ≤ s, q ≤ ∞) of u. For proofs and further results, we refer to [7,40,42,58].
A classical idea of Leray [41] to construct solutions with singularities uses an ansatz of self-similar solutions of the form with V ∈ L 3 (R 3 ), getting singular as t → 0−. Since u is a weak solution of the Navier-Stokes equations, the vector field V satisfies a t-independent modified Navier-Stokes equations with an additional y-dependent term of the form y · ∇ y V . However, the problem to find a nonzero field V leading to a singular weak solution was left open and answered not before the end of the last century by Nečas et al. [43] in the negative (V ≡ 0). The same negative answer was given by Tsai [56] in case that either V ∈ L q (R 3 ) for any q > 3 or the self-similar solution u has bounded energy on the space-time cylinder (−1, 0) × B 1 (0) in the sense that ess sup In the pointwise sense, the singularities discussed above are of the form |u(x, t)| ∼ c |x|+ √ −t and called singularities of Type I. They can be excluded in the ansatz of self-similar solutions as above and for axisymmetric flows, see [14,38]. All other singularities are called of Type II.

The counter-example by Terence Tao
A solution with smooth initial value and singularity in finite time was constructed by Tao [55] for a modified averaged Navier-Stokes system on R 3 . The counter-example is based on a blow-up solution for the Euler system so that the dissipation term −Δu does act only as a pertubation and finally allows for a singularity of Type II. The first idea is a discrete (dyadic) dynamical system of scalar functions (X n (t)) n∈Z of the form ∂ t X n (t) = −λ 2nα X n + λ n−1 X 2 n−1 − λ n X n X n+1 (24) with parameter λ = (1 + ε 0 ) 5/2 , 0 < ε 0 1. Actually, in the final analysis (X n ) is a family of vector fields in R 4 with pairwise disjoint supports in Fourier space such that energy is either pumped or drained or rotated between different components, either slowly or abruptly. Then, the solution u blowing up in finite time has the form u(x, t) ∼ n X n (t)λ 3n/5 ψ(λ 2n/5 x) and is based on a Schwartz function ψ with Fourier transform supported in an annulus of radii λ −5/2 and λ 5/2 . To control the nonlinear transport of energy in (24) it would be helpful if at a given time only one pair (X n , X n+1 ) interacts. This can be achieved by a introducing a piecewise constant function n = n(t) in (24), but would lead to a t-dependent nonlinear term in the Navier-Stokes system. Here and for the transformation to a continuous model a sophisticated averaging process plays a crucial role. To this aim, we return to the bilinear operator B(u, v) := − 1 2 (u · ∇v + v · ∇u) and the trilinear form which vanishes for solenoidal vector fields u = v = w, see (8). Then, Tao introduces the operators Ru(x) := R ξ u(R −1 ξ x) defining the rotation around the axis ξ ∈ R 3 which map solenoidal vector fields to solenoidal ones. Moreover, classical multiplier operators Mu = F −1 m(·)Fu on L p (R 3 ) are needed. With these definitions, the trilinear form B(u, v), w R 3 is replaced by the averaged form where ω is a random variable defined on a subset Ω of the space The crucial properties ofB are the cancellation property B (u, u), u R 3 = 0 and energy estimates as satisfied by B. The analysis is performed in Fourier space in which the energy cascade is mainly controlled by balls tending to infinity as time t approaches the finite blow-up time T * more or less by jumps on a sequence of epochs (t n ) tending to T * exponentially fast. Tao's example implies that an affirmative answer to the regularity problem must either be based on more sophisticated estimates and hidden properties of B which do not hold forB or on a further balance among the terms u t + u · ∇u, −Δu and the pressure ∇ p.

Recent approaches to non-uniqueness
In [36], Jia and Šverák construct a situation with initial values admitting multiple solutions. Starting point is a self-similar initial value u 0 on R 3 such that the Navier-Stokes system possesses for each σ u 0 , 0 ≤ σ 1, a unique global self-similar solution u σ (x, t) = t −1/2 U σ (xt −1/2 ), cf. (23). Linearizing the modified Navier-Stokes equation for U σ around U σ the linear operator L σ , is the crucial operator for a further bifurcation analysis. It leads to certain assumptions on eigenvalues of L σ such that for σ > σ 0 > 0 a pitchfork bifurcation will occur. The theoretical assumptions could not yet be verified rigorously in concrete cases. However, numerical tests by Guillod and Šverák [31] indicate that the above scenario will develop. The method of convex integration developed about 70 and 50 years ago by Nash and Gromov is the basic tool of De Lellis und Székelyhidi Jr. to prove the existence of weak solutions to the 3D Euler equations in the class C β , 0 < β < 1 3 , in space-time with prescribed kinetic energy E kin (t) ≥ 0, s [6,16,17]. In the last years, these results were generalized to the Navier-Stokes equations where the smoothing property by dissipation (friction) was the crucial problem: there exist "weak" solutions u in the class C 0 ([0, T ]; L 2 (Ω)) such that u will neither necessarily satisfy the energy inequality nor have finite dissipation energy t 0 ∇u(τ ) 2 2 dτ . Buckmaster and Vicol [4] prescribe a smooth function E kin (t) ≥ 0 and construct a weak solution u ∈ C 0 ([0, T ]; L 2 (T 3 )) ∩ W 1,1+β (T 3 )) with β > 2 −16 and 1 2 u(t) 2 2 = E kin (t). Since E kin is not assumed to be monotonically decreasing, u will not be a physically reasonable solution. This result was extended by those authors and Colombo [5] as follows: given two smooth solutions u (1) , u (2) of the Navier-Stokes system on T 3 × [0, T ] with f ≡ 0 there exists a weak solution in the above sense with u(0) = u (1) (0) such that u ≡ u (1) in [0, T /3] und u ≡ u (2) in [2T /3, T ].
In particular, choosing u (1) ≡ 0 and u (2) = 0 there exists a solution which is at rest on some initial time interval, but bifurcates to a nonzero smooth solution. On the other hand, u may start with u ≡ u (1) = 0 on [0, T /3] but come to rest at t = 2T /3. Nevertheless, the set S of singular instants of u is small in the sense H 1−β (S) = 0.
Since these solutions cannot be shown to satisfy the energy inequality on [T /3, 2T /3], the question arises whether there exist weak solutions in the sense of Leray-Hopf with given physically reasonable kinetic energy function. A partial answer is given by Ożański [48] who prescribes an arbitrary continuous, non-increasing energy function E kin (t) and constructs a suitable weak solution u satisfying the localized energy inequality such that 1 2 u(t) 2 2 − E kin (t) is bounded by any given positive constant. However, u satisfies a Navier-Stokes inequality, i.e., u is a solution in the Leray-Hopf class with the property that f := u t − νΔu + u · ∇u + ∇ p, f · u ≤ 0.
In other words, the force f is not given, but is defined a posteriori and satisfies f · u ≤ 0; i.e., formally f acts against the direction of the flow field.

Conclusion
Even after more than 85 years when Leray constructed global weak solutions, the regularity problem for global-in-time weak solutions of the Navier-Stokes equations or the problem of the direct construction of global-in-time regular solutions is still open. After an analysis in classical L 2 , L p −spaces, Lorentz and Sobolev spaces, with and without weights, modern analysis also works in Besov spaces, but there are also results on regularity in spaces of type BMO, Hardy H 1 , Triebel-Lizorkin, Morrey, and others. However, the gap between weak solutions with Serrin number S = 3 2 and strong solutions with S = 1 could not be closed. Further work on the construction of solutions with finite-in-time blow up either failed (e.g., self-similar solutions by Leray) or lead to solutions of a modified Navier-Stokes system (Tao) where it is not clear whether and how the construction can be transferred to the original PDE system. A recent trend is to prove non-uniqueness in a new class of weak solutions. However, these solutions which have striking and unexpected behavior do not satisfy the physical properties (energy inequality, finite dissipation energy) as required from the physical point of view.
Funding Open access funding enabled and organized by Projekt DEAL.
Open Access. This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/ by/4.0/.
Publisher's Note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.