The Navier-Stokes equations in primitive variables
F. Lam


Navier-Stokes Equations]The Navier-Stokes equations in primitive variables

F. Lam]F. Lam

 



[
    [
    July 9, 2019
================



Navier-Stokes Equations; Leray-Hopf Local Solution; Diffusion; Linearisation; Perturbation; Couette; Hagen-Poiseuille; Primitive Formulation

The Navier-Stokes equations in the primitive formulation for incompressible flow describe the evolution of velocity and pressure, without recourse to vorticity. We show that, beyond the finite Leray-Hopf regularity interval, every postulated strong solution is accompanied by infinitely many diffusion-dominated percolations of arbitrary size, while the momentum deficit caused by the non-linearity is compensated by the pressure gradient. In the upper half space, we demonstrate how sequences of these collective companions can be re-scaled into an absurd singularity. Owning to the passive nature of the pressure, there exist no essential a priori bounds for establishing the uniqueness of primitive solutions. With the illustration of well-exploited examples of closed-form basic flows, we elucidate the reason why perturbations, infinitesimal or finite, instigate indeterminate states that render the concept of flow instability inadmissible. An effort has also been made to reappraise a number of important issues in fluid dynamics. Unfortunately, the primitive theory cannot serve as a reliable tool for prediction. 

Nevertheless, a dedicated effort has been made to elaborate a priori bounds for vorticity dynamics for ideal as well as real fluids. As a result, we are able to establish long-time regularity of the Cauchy problem for incompressible flows. A main conclusion is that no events of finite-time blow-up can ever occur in the Euler or Navier-Stokes equations for initial data of finite energy.






§ INTRODUCTION



The Navier-Stokes equations of motion are derived from the conservation principle of linear momentum, and the continuity (Navier 1823; Stokes 1845). For viscous incompressible flows, they are 

    ∂_t  + (·∇)   = - ρ^-1∇ p + νΔ,   ∇·=0,

where the velocity vector (,t)=(u,v,w)(,t), the coordinates =(x_1,x_2,x_3), and the scalar p(,t) denotes the pressure, and the kinematic viscosity, ν=μ/ρ, is the ratio of the dynamic viscosity and the density. We assume that body force has a potential which has been absorbed into the pressure. The dynamics is invariant in any homogeneous direction. In ^3, we must observe

    (,t)=(-,t),   ∇ p(,t)=∇ p(-,t),

for arbitrary finite vector, =(a_1,a_2,a_3). In cylindrical polar co-ordinates (r,θ,z), the axial z-direction is homogeneous in a number of specific cases. 

The initial solenoidal velocity is assumed to be smooth with compact support

    (,0)=_0() ∈  C_c^∞(^3).
 
Specifically, initial data, _0=(u_0,v_0,w_0), can be an element of the Schwartz class of functions (denoted by 𝒮). Every component of _0 is a subset of L^1(^3), defined to be smooth functions, where the derivatives of all orders decay faster than any polynomial,

    𝒮() = { v ∈ C^∞():         α,  β∈ℕ,     _0, 
       sup_∈|(-_0)^β  D^α v(-_0)| < ∞      },

where constant _0 is specified with the data. Any decay property of v at large distances is with respect to _0, in view of the invariance (<ref>).
In many applications, the requirement that initial data are elements of the Schwartz class (S) can be unnecessarily restrictive. Instead, smooth functions of

    |_0| ≤B/(A+|-_0|)^n,    (  A > 0,  n ≥ 1),

are admissible alternatives. The nominator B is a constant or specific bounded functions, |B(,_0)| ≤ 1, for instance. Initial vorticity is found from 

    _0()=∇×_0().


For any arbitrary smooth function, =(t), (0)='(0)=0, spontaneous solutions to the Navier-Stokes or the Euler equations are given by (,∇ p)=(,- ρ'). In order to avoid under-determinacy of the Cauchy problem, we need to impose a `boundary condition' at infinity or an equivalent constraint. Physics suggests that the energy must be finite,

    ∫_ ||^2  < ∞,

for t ∈_T, where _T=[0, T< ∞). In the Sobolev space, the finiteness translates into u_i ∈ W_0^0,2() or  having compact support. Instead of (<ref>), we require that the velocity tends to zero at large distances,

    → 0        || →∞.

At this stage, we can only state the decay loosely, as we do not have more knowledge about the precise decay rate. In subsequent discussions, all similar expressions of attenuation are understood in the relative sense, |-_0| →∞. Conversely, this property implicitly expresses a bound that ∈ W_0^0,1 on . In an extreme scenario, for a source singularity in potential flow, it is known that the velocity decays at infinity ∼ O(||^-2). 
By virtue of Gauss's divergence theorem, we have 

    ∫_∇·  = lim_r→∞∫_S ·n⃗  = 0,

where, S=4 π r^2, stands for the surface area of a ball centred at the origin _0, the radius r of which tends to large values. It follows that, as one possibility, the continuity is valid provided ∼ O(||^-2) as ||→∞.[Condition (<ref>) merely states that the sum of the radial velocity must be zero. There remains a possible scenario that some fluid enters one portion of the spherical surface while the exact same amount leaves at another portion. Moreover, we do not have complete knowledge about local tangential components; the whole ball might be under a uniform rotation, for instance. As r tends to be unbounded, these non-zero velocities cannot be sustained for finite energy.] Similarly, the solenoidal vorticity implies that, at least, ∼ O(||^-2) at infinity. Hence `boundary condition' (<ref>) does not compromise analytical scopes of flow-fields consisting of - and, in fact, it serves our purposes well, compared with periodic flows.

In bounded domains with C^2 boundaries, we relax the data 

    (,0)=_0() ∈  C^∞(),

as we are dealing with flows of finite energy. To be definitive, we also impose the no-slip condition, (,t)=0,∀∈.
Taking divergence on (<ref>), and using the continuity, we obtain the Poisson equation for the pressure 

    Δ p()= - ρ (  |∇|^2 - |∇×|^2   )().

In bounded domains, this is a Neumann boundary value problem.

The conservation laws are universal and, indeed, they do not depend on the measuring units in any particular problem. We view the equations of motion (<ref>) (and associated initial and boundary conditions) as inherently dimension-independent with respect to the SI standard references of length, time, and mass. To be consistent, the viscosity is a physical property of fluids, hence it is made non-dimensional in reference to unit kinematic viscosity ν_0=1  ([ ^2  ^-1 ]). In fact, a reference kinematic viscosity rather than unity can be chosen, as long as the ambient parameters are standardised, for instance, at a fixed set of the SI temperature and density.

The theoretical framework, (<ref>) to (<ref>), is known as the initial-boundary value problem in the primitive variables, (,p). In a nutshell, the primitive equations admit local-in-time smooth solutions, depending on the size of  the initial data (Leray 1934a; Hopf 1951); the global regularity is out of the question, due to lack of suitable a prioir bounds. The popularly-adopted classification of two distinct modes of flow, streamlined laminar and fluctuating turbulent, is completely based on the primitive setting (Reynolds 1883). Perturbation approaches on shear flow instabilities, linear or non-linear, are all formalised in terms of the -p variables. Assisted by modern computers and algorithms, numerical solutions are readily available in applications, despite the absence of rigorous mathematical justifications. It has been well-known that the computations themselves often stimulate new analytical problems, such as error-controls in discretisation, convergence and stability of numerics. In the present note, we examine the compatibility and the global regularity of the primitive formulation. 

By examining Helmholtz's decomposition of vector fields on , we establish the vorticity-velocity compatibility,  Δ = - ∇×. The importance lies in the fact that fluid dynamics is tantamount to solving vorticity equation (see equation (<ref>)) in conjunction with this kinematic constraint. Then, we make an effort on derivation of a priori bounds for the vorticity equation. These estimates are crucial in determining the global regularity of the Cauchy problem.




§ NON-UNIQUENESS



Let the Leray-Hopf regular solution be (^*,p^*), and its existence time interval 0 < t ≤ t^*. 
Suppose that there is an extended local solution, or even a global solution of the Navier-Stokes equations, denoted by (,p)(, t > t^*) with the starting data (<ref>). We consider the case that the solution is smooth, or sufficiently differentiable. The extended solution may be obtained by either an analytic construction or a numerical procedure. We are particularly interested in the numerical solutions established from a projection method, or a fractional step technique, where the primitive variables, (,t) and p(,t), are iterated in parallel. It is evident that the solutions (,p) and (^*,p^*) must coincide over the interval 0 ≤ t ≤ t^*. 




 §.§ Weak formulation and buffer functions



All the derivatives in the equations of motion are interpreted in the sense of distributions. Let the test functions,

    (,t)=(φ_1, φ_2, φ_3)(,t),   (∈ C^∞_c(^3)×)

be vectors. Multiplying the first equation of (<ref>) by , integrating the result over space and time, followed by integration by parts, 
we obtain the weak formulation

    ∫_∫_^3 ( ∂/∂ t + νΔ )  t = ∑_i,j=1^3∫_∫_^3 u_i ∂ u_j/∂ x_iφ_j   t - ∫_∫_^3 p ∇.  t.

This is Leray's definition of weak solutions (Leray 1934a), where the velocity  is square-integrable, and the pressure p locally integrable. If the flow field is found to be incompressible and sufficiently regular, the continuity (<ref>) makes sense. Leray also showed that weak solutions (<ref>) satisfy energy inequality,

    1/2∫_^3 |(,t)|^2  + ν∫_0^t∫_^3 |∇(,t)|^2  t  ≤ 1/2∫_^3 |(,0)|^2 .


As proposed by Hopf (1951), we may emphasise the role of the initial condition,  as long as the test functions are re-defined as

    (,t)=(φ_1, φ_2, φ_3)(,t),   (∈ C^∞_c(^3)× [0,∞))

on a set of measure zero (time-wise). The formulation (<ref>) is extended to 

    ∫_0^∞∫_^3 ( ∂/∂ t   + νΔ )  t + ∫_^3 _0  (0)  
        - ∑_i,j=1^3∫_0^∞∫_^3 u_i ∂ u_j/∂ x_iφ_j   t + ∫_0^∞∫_^3 p ∇·  t=0.

Apparently, the revised formula does not contain noticeable differences compared to Leray's. However, given every initial velocity ∇._0=0, the linear kernel, ∂_t - νΔ, is data-preserving in the sense that the solenoidal field is maintained throughout the motion. Thus the continuity of (<ref>) makes sense and is readily satisfied. This observation suggests that the formulation (<ref>) can be interpreted as a composition of two degenerated dynamics,

    ∫_0^∞∫_^3 ( ∂/∂ t + νΔ )  t + ∫_^3 _0  (0)    = ∫_0^∞∫_^3  t,

and

    ∑_i,j=1^3∫_0^∞∫_^3 u_i ∂ u_j/∂ x_iφ_j   t - ∫_0^∞∫_^3 p ∇·  t=∫_0^∞∫_^3  t,

if there exists a function =(,t) mediating the parities. Although the Cauchy problem (<ref>) or (<ref>) shares identical initial data with formulation (<ref>), their solutions do not overlap in general, except over a very short interval [0,ε], as the former setting aims at wholly non-local and non-linear solutions. 

To determine such a function , we first state the precise definition of the solutions of (<ref>). It is known that the regular Leray-Hopf solutions do exist over a finite period of time. As strong solutions are also weak solutions, we shall work with regular functions and tactically put aside the regularity issue which entails a priori bounds.  We call (,p) a Navier-Stokes solution if the pair satisfy all of the following conditions: 

    .
        1.  , ∇·=0;    t > 0, 
        _0; 
        2.   t > 0;
        3.  ,
       Δ p = - ρ∇·(  (.∇) )   
    ( p_0(t));
        4.  ,   ν>0,
       1/2∫_^3 ||^2  + ν∫_0^t∫_^3 |∇|^2  t   =  1/2∫_^3 |_0|^2 .
    }

In addition, the Neumann boundary value problem for the pressure is self-consistent. Conceptually, the last two constraints are de facto the derivatives of the first two.

Let ∈ C^∞_c(^3) or ∈ C^∞() in bounded domains. The second condition implies that the momentum equation in (<ref>) or its weak forms, (<ref>) and (<ref>)-(<ref>), can be interpreted as one of the following equivalent statements:

    ∂_t    + (·∇)   = - ρ^-1∇ p + νΔ       ⟺
    ∂_t    - νΔ= - (·∇)  - ρ^-1∇ p =        ⟺
    ∂_t    - νΔ -  = - (·∇)  - ρ^-1∇ p -        ⟺
    ∂_t    + (·∇)   + ρ^-1∇ p - νΔ=0,

where each form respects the conservation principle for every solenoidal field . We single out the `forced' diffusion dynamics, 

    ∂_t  - νΔ = ,

subject to the solenoidal initial data _0. Its solution can be expressed as

    (,t)     = ∫_^3 H(-,t)_0()  + ∫_0^t ∫_^3 H(-,t-t')(,t')   t' 
          = ∫_^3 H(,t)_0(+)  + ∫_0^t ∫_^3 H(,t-t')(+,t')   t',

where H is the free-space heat kernel,

    H(-,t)=( 4 πν t )^-3/2exp( - | -  |^2/(4 ν t) ).

It is this property of continuous dependence on initial data that motivates us to split the pressure-related momentum equation into two separate parts. From (<ref>), we see that the pressure can be found from the Neumann problem

    Δ p = - ρ∇·(  (·∇) +  ),

for a proper choice of the function , because the non-linear term is known, though the convection now becomes moderated. To meet every criterion of (<ref>), we specify that the buffer function is calorically mollifiable, and  within free multiples,

    (,t)=(t) (),

where time-dependent  is any regular function, and  satisfies one of the following conditions:


  * ∇·=0, and ·=0 (orthogonal for energy conservation);

  * =∇×, and ·=0;

  * =_1×_2, _1 and _2 being irrotational, and ·=0;

  * =(∇ψ) ×(∇ϕ), for scalars ψ and ϕ, and ·=0;

  * =Δ, for solenoidal , and ∇·∇=0;

  * = ∇ϕ(), and ϕ is harmonic.

Choosing a buffer from the list, we see that the Poisson's equation (<ref>) is reduced to the third condition of (<ref>). To balance the momentum, we simply set

    ∇ p = - ρ(  (·∇) +  ),

where the right-hand side is completely given. For bounded domains, the evaluation is slightly more complicated but is tractable; the dependence of the pressure on  can be found from the Neumann boundary conditions, unless we further impose the no-slip on . These relations rectify the auxiliary nature of the pressure, or more precisely its gradient ∇ p. 

We point out that list (<ref>) is by no means comprehensive. Linear superpositions of these classes are also good candidates. Specifically, the Helmholtz's vector decomposition is included. The last category represents weak solutions of (<ref>) for any solenoidal test function (Serrin 1962). 




 §.§ Percolation and pressure compensation



Now, we return to the issues raised at the beginning of the present section. Given a regular solution (,p) of (<ref>), or (<ref>), which has been postulated to exist beyond the Leray-Hopf time t>t^*. Consider the following superimposed pair

    ( +α,  p+ /β )(,t),

where the parameters, α and β=β(t), are arbitrary but finite. The pair (,) are regarded as perturbations to the postulated solution. We substitute the superposition (<ref>) into the Navier-Stokes equations (<ref>). Since solution (,p) satisfies the equations of motion, we are left with

    ∂_t  + α (·∇) + (·∇) + (·∇) = -(αβρ)^-1∇ + νΔ,   ∇·=0.

To solve these exactly, we characterise the following diffusion problem as part of the solution:

    ∂_t  - νΔ = ,

subject to the initial conditions

    _0()=^*(,t^*)=(,t^*).

An extension of the Leray-Hopf solution is feasible, as the diffusion equation has solution 

    (,t)=∫_^3 H(,t)  _0(+)   + ∫_t^*^t ∫_^3 H(,t-t')(+,t')   t'.

Since the initial data ∇·(,t^*)=∇·^*=0, the solution  remains solenoidal for all times, as indicated by the expression in (<ref>), and the properties of the buffers. Moreover,  is globally regular in space and in time, by virtue of the well-known smoothing properties of the heat kernel. In fact, we only need 

    _0(,t^*) ∈ L^1_loc() ∀∈^3.


The analytical structure of equation (<ref>) suggests that we can compensate the momentum deficiency due to (,t) exactly, by a complement of pressure gradient ∇. Analytically, equations (<ref>) are fully satisfied as long as the pressure  is found from

    ∇= -   A_0  (  α (·∇) +(·∇) + (·∇) ) - ρβ = - ,

where A_0=αβρ. In practice, only the gradient rather than the actual pressure matters. Nevertheless, we proceed to describe a general approach. Let us introduce the abbreviations:

    [_x]=(∂ u/∂ x    ∂ v/∂ x   ∂ w/∂ x)',      
    S[u]=(∂ u/∂ x   ∂ u/∂ y   ∂ u/∂ z),

where the symbol  stands for a 3×1 column matrix operating on vector , S for a 1×3 row matrix on a scalar. Taking divergence on equation (<ref>), or equivalently on (<ref>), we readily establish the following equation:

    Δ()= - A_0 (    α(  S[û] [_x]+S[] [_y]+S[] [_z]  )   + 
          2 (   S[u] [_x]+S[v] [_y]+S[w] [_z]  ) )() = -F().

Evidently, the pressure  is determined, up to a constant, for fixed α and β because each of the 18 velocity terms is known, or

    = 1/4 π∫_^31/|-| F()   + _0.

For every postulation (,p), the full Navier-Stokes system (<ref>) admits the companion solutions (<ref>) over the time interval t > t^*. We call the pair (,) a (viscous) percolation, because they diffuse everywhere at all times.  

As a well-known practice, the energy conservation is obtained by integrating complete momentum (<ref>) or the diffusion (<ref>)  

    1/2∫_^3 ||^2  + ν∫_0^t∫_^3 |∇|^2  t= 1/2∫_^3 |_0|^2 ,

which is scale-invariant to parameter α. We draw our attention to the converse: if we rely on the law as a bound, we cannot distinguish a percolation from a complete non-linear solution. In particular, the percolation is arbitrary, while the full solution, or a result of numerical computation, is expected to be unique.

In numerical schemes, such as a fractional step algorithm or a projection method, the appearance of non-zero buffer function  is inevitable. Over the mesh and time, discretisation errors of the operator ∂_t -νΔ may well be in the same order of those on the non-linearity and the pressure (cf. (<ref>)), implying the existence of buffers. If we strictly solve the primitive equations, we have no reliable means to detect or filter out the differences of the errors which evolve with the computations. In the case of zero buffer, the structure of the primitive formulation still allows the flow-field partition, according to the criteria (<ref>).

It has been suggested that, the pressure plays the role of a Lagrangian multiplier which ensures incompressibility. This is a misleading view, given the partition. The gradient is determined by (<ref>), in a passive manner, resulting in the superficial Poisson equation (<ref>). One implication is that use of the Helmholtz-Leray projection, which eliminates the pressure in conjunction with the continuity, lacks generality.




 §.§ Multiplicity inside a cube



Let the domain be the cubic box with all sides of unity, 
_b:   0 ≤ x,y,z ≤ 1.
We impose the no-slip condition on the walls of the cube.[It has been known that the initial-boundary value problems for incompressible flows (ν≥ 0) are ill-posed in periodic domains with periodic boundary conditions (Lam 2018b). We emphasise the fact that the conjecture of anomalous energy dissipation for inviscid-limit continuum turbulence by Onsager (1949) was made on an assumption of the periodicity for the primitive variables. In the last paragraph on p.286, Onsager discussed his hypothesised breakdown of energy conservation in terms of the Fourier formulation of the governing equations on the macroscopic scales (cf. his equations (15) and (17) of p.284).] We choose =0 so that the velocity is defined by pure diffusion

    ∂_t  - νΔ = 0.

By analogy, the companion velocity to any extended Leray-Hopf flow is still evaluated by the first integral of (<ref>). Instead of (<ref>), the heat kernel now reads

    H_f(x,x',t) × H_f(y,y',t) × H_f(z,z',t),

where the finite-line heat kernel has the form

    H_f(x,x',t)     =2  ∑_k=1^∞ sin(kπ x)sin(kπ x')exp(-ν k^2 π^2  t )
          = 1/√(4 πν t)∑_k=-∞^∞[ exp( -(x-x'-2k)/4 ν t)-exp(-(x-x'+2k)/4 ν t) ].

Since =  = 0 on the walls, the zero Neumann boundary condition is found from (<ref>), namely ∂/∂n⃗|__b=0 in the present case. The pressure is governed by a compatible system because, from (<ref>) and (<ref>), we have

    ∫__b F()  = - ∫__bΔ ()   = ∫__b (∇·)()   = ∫__b (·n⃗)   s=0,

by virtue of Gauss' divergence theorem. For the demonstration of the existence, it is preferable to express the pressure solution (<ref>) in terms of the Neumann function, N_b, which is defined by 

    Δ N_b=-δ(,')+1,   ∂ N_b/∂n⃗|__b=0,

where δ=δ(x) is Dirac's delta. We start from a cosine expansion for the delta in one dimension,

    δ(x-x')=1/L+2/L∑_l=1^∞cos(lπ x/L)cos(lπ x'/L) = 1/L+2/L∑_l c_l^x c_l^x',

for 0 < x,x'< L. Then it is easy to establish, with justifications, see, for instance, Chapter 9 of Roach (1982), 

    -δ(-')+1/L^3   =-2/L^3(∑_l c_l^x c_l^x'+∑_m c_m^y c_m^y'+∑_n c_n^z c_n^z'+2∑_l,m c_l^x c_l^x'c_m^y c_m^y'
          + 2 ∑_m,n c_m^y c_m^y'c_n^z c_n^z'
    +2 ∑_n,l c_n^z c_n^z'c_l^x c_l^x' + 4∑_l,m,n c_l^x c_l^x'c_m^y c_m^y'c_n^z c_n^z').

By symmetry, the Neumann function (L=1) must have the neat expansion,

    N_b(,')=8∑_l,m,n=0 A(l,m,n,x',y',z')  c_l^x c_m^y c_n^z,

which ensures the derivative homogeneous data on the boundary. From the differential equation (<ref>), we find the required Fourier expansions in cosines

    N_c(,')   =-8/π^2∑_l,m,n=0^∞( A_lmncos(l π x)cos(l π x') cos(m π y)cos(m π y')
         cos(n π z)cos(n π z') ) / ( l^2+m^2+n^2 ),

where the normalising factors are

    A_000   =0;    A_00n=A_0m0=A_l00=1/4; 
    
    A_0mn   =A_l0n=A_lm0=1/2;      A_lmn=1,      l,m,n≥1.


An alternative way to determine the pressure is as follows. The Neumann boundary value problem,

    Δ = - F,     ∇|__b=0,

can be reduced to a Dirichlet problem, as -∇|_∂_b is found from momentum (<ref>). This approach may have an advantage in the circumstances where we have computed the velocities as well as their derivatives. By virtue of (<ref>), what we need to analyse is a vector Poisson equation for  with homogeneous boundary data

    -Δ(∇)=∇ F,     |__b=-∇|__b=0.

Solving (<ref>) without specifying the unknown gauge, we find the gradient 

    ∇() = ∫__b G_b(,')∇ F (')  ',

where the Green's function G_b equals

    8/π^2∑_l,m,n=1^∞sin(l π x)sin(l π x') sin(m π y)sin(m π y')sin(n π z)sin(n π z') /( l^2+m^2+n^2 ).


Our partition proposition, as afforded by (<ref>), (<ref>) or (<ref>), tactically bypasses the task of handling the non-linearity. It is instructive to recall that Stokes' system,

    νΔ - ∇ p /ρ = σ;   ∇·=0,   (|__b=0)

has been extensively analysed for the spectra of eigen-values σ_k, and the normalised eigen-functions. A weak solution of (<ref>) is represented, in a Faedo-Galerkin approximation, by the complete orthonormal basis, in analogous ways to Fourier series expansions. Technically, this is a linear strategy as well. However, it by no means inspires the `lower-order' view that prioritises the dynamic role of viscous diffusion.




 §.§ Artificial singularity in upper half space



Denote the flow domain by _h: -∞ < x,y < ∞; 0 ≤ z < ∞. Let (^*,p^*) be a Leray-Hopf regular solution in 0≤ t≤ t^*, where ^*(x,y,0,t)=0. Suppose that there is a mean solution, (,), which extends (^*,p^*) beyond t^*. At t=t^*, these two solutions coincide. The no-slip on the plane z=0 holds for the mean motion t>t^*. Consider a bounded time interval [t^*≤ t ≤ T], which T lies within the existing time of the mean motion. We then sub-divide it into n equal parts, Δ t= (T-t^*)/n. Denote t_k=t^*+k Δ or

    t^*< t^* + Δ t <  t^* + 2Δ t < ⋯ < t^* + k Δ t < ⋯ < T.


We introduce the first percolation (^(1),^(1)) on the mean motion over [t^*,t_1], starting from

    _0() = ^*(,t^*)=(,t^*).

The solution of (<ref>) is expressed in the integral

    ^(1)(,t)=∫__h H_s(,',t) _0(')  ',   (t^* < t < t_1)

where the kernel H_s is familiar and vanishes on the plane z=0,

    H_s(,',t)=H_1(x,x',t) × H_1(y,y',t)×(H_1(z,z',t)-H_1(z,-z',t)),

and the function H_1 is the free-line heat kernel

    H_1(x,x',t)=1/√(4 πν t) exp( - (x-x')^2/(4 ν t) ).

Given the no-slip on z=0, the Neumann function can be found by the method of images. The formula for the solution of (<ref>) reads

    ^(1)()= 1/4 π∫__h( 1/r_-+1/r_+ )F^(1) (')  ' + _0,

where F=F^(1) contains the first percolation ^(1), and

    r_∓(,')=( (x-x')^2+(y-y')^2+(z∓z')^2 )^1/2.

The first companion solution is

    (_c, p_c) = (+α^(1)^(1), +^(1)/β^(1))

over [t^*,t_1]. Next we introduce a second percolation, (^(2),^(2)) on (_c,p_c) over [t_1,t_2] with starting data

    _0()=_c(,t_1).

We repeat our analysis to fix the second percolation (^(2),^(2)), and renew the forcing function F^(2)=F^(2)(^(2)) in (<ref>) to settle ^(2). Now the companion is updated to

    (_c, p_c) = (+α^(1)^(1)+α^(2)^(2), +^(1)/β^(1)+^(2)/β^(2))

over [t^*,t_2]. Assume that, after k-time percolations, we obtain the sum

    +∑_i=1^kα^(i)^(i), +∑_i=1^k^(i)/β^(i)

over [t^*,t_k]. To determine the (k+1)th percolation, we substitute (_c+α^(k+1)^(k+1), p_c+^(k+1)/β^(k+1)) into the Navier-Stokes and subtract the companion, we obtain the analogous equations (<ref>) and (<ref>) for ^(k+1) and ^(k+1) respectively. We solve 

    ∂_t^(k+1)-Δ^(k+1)=0,

subject to data _0=_c(,t_k), and its solution is expressed in the integral (<ref>). The pressure, ^(k+1), is given in the formula (<ref>) with F^(k+1)=F^(k+1)(^(k+1)). Evidently, all of the the percolations are exact and smooth. Specifically, the companion solution to (,) over [t_k-1,T] is the superposition of the individual percolations

    _c(,t)=∑_k=0^n α^(k)^(k)(,t),     p_c(,t)=∑_k=0^n^(k)()/β^(k)(t),

where ^(0)=, α^(0)=1; ^(0)=, β^(0)=1.
 
The localisation and smoothing properties of the heat kernel (<ref>) show that every ^(k) decays in time and in space, relative to its initial data _0. Moreover, all the percolations vanish on the plane z=0. Choose =_s ∈_h, t=T so that (_s,T)≠0. We fix our α's by

    α^(k)^(k)(_s,T) = (_s,T)=A_s.

Thus, the velocity sum becomes _c(_s,T)=(n'+1)A_s, where n'≤ n, taking into account a finite number of (possible) zero .
By increasing the number of intervals over [t^*,T] and by Cantor's diagonal process, we assert that

    _c(_s,T)=lim_n→∞(n'+1)A_s  ⟶∞.

Furthermore, we choose a sub-interval, [t_k,t_k+1], and divide it into m equal parts, and repeat the whole analysis. The result is

    '_c(,t)=∑_i=0^m α^(i)'^(i)(,t),   t_k≤ t ≤ t_k+1,

which diverges at a point, (_s,t_s), as m →∞, as long as (_s,t_s) does not vanish.

There is a considerable leeway to construct β(t) in the sum (<ref>), in order to deliver a divergent gradient in time, perhaps as an algebraic or exponential function. The perturbations given in (<ref>) do not have to be singular, nevertheless, they do imply the existence of multiple backdrop flow-fields, starting from the identical initial conditions. 
Therefore, it makes no sense in talking about an a priori bound on the mean (,). The essence is that the Leray energy inequality (<ref>) or similarity scalings are simply off-target, if not irrelevant. Without a background of the unpredictable percolated flows, claims of finite-time blow-up in the vicinity of the plane z=0 represent a misunderstanding of fluid motion. 




§ FLOW BETWEEN PLATES



We consider two steady planar flows between two parallel plates, namely plane Couette and Poiseuille flows. These flows have a long history of study, mainly in  stability theory. Without causing confusion, we will denote _2=(x,y) and =(u,v). The channel is unbounded: 
_0: -∞ < x  < ∞;  0 ≤  y ≤  1.
The governing equation for real fluids μ>0 has a simple form,
μ^2 u/ y^2 = p/ x. 
Plate boundary conditions yield a couple of basic profiles

    ()   :    u(x,0)   =0,     u(x,1)   =1;    ∂_x p   =0,       u(x,y)   =y;
    
    ()   :    u(x,0)   =0,     u(x,1)   =0;    ∂_x p   =-2μ,     u(x,y)   =y(1-y).

We prevent ourselves from challenging the time-independence. To investigate the flow development with respect to perturbations, we consider the initial-boundary value problem of the full equations of motion (<ref>). Specify the initial perturbation, _0, by

    ( [ _0; _0 ])(x,y)
    	=1/4 π∫__0log(  cosh(x_-)-cos(y_+)/cosh(x_-)-cos(y_-) ) ( [  ∂_y'ω_0; -∂_x'ω_0 ])	  x'  y',

where x_-=π(x-x'), y_±=π(y ± y'), see, for instance, Carslaw & Jaeger (1986). The properties of the logarithmic Green's function show _0 satisfies the no-slip on the plates. The smooth function, ω_0(x,y), is arbitrary. We require it to be any localised and bounded function around the origin. Since both the planar flows are time-independent, we denote the starting instant when the perturbation is first introduced by t=0. 




 §.§ Plane Couette and plane Poiseuille profiles





Consider the perturbation on the linear profile, 
(+α, /β), =(,). In the governing equations (<ref>), we replace  by , and set   the solutions of the planar diffusion equation, 

    ∂_t  - ν(∂_xx + ∂_yy) =0,

subject to starting data _0, and homogeneous boundary value, |__0=0. As the heat kernel takes the form, H_2= H_1(x,x',t) × H_f(y,y',t), where H_f and H_1 are the kernels (<ref>) and (<ref>), we find the solution of (<ref>)

    (_2,t)= ∫__0 H_2(_2,_2,t) _0(_2)  _2.

The size of the percolation flow-field is much larger than unity for some values of α, so that the rectilinear flow disappears soon after the perturbation is set-in. In two space dimensions, the right-hand function in (<ref>) takes a simpler form: 

    Δ= - A_0  (    α(  S_2[û] _2[_x]+S_2[] _2[_y] )   + 
          2 (   S_2[u] _2[_x]+S_2[v] _2[_y] ) )   (= -F_2),

where, by analogy, _2[_x]=(∂_x u  ∂_x v)', S_2[u]=(∂_x u ∂_y u). Because of the linearity, the last term in (<ref>) is reduced to 2 ∂_x. Given , the gradient counterpart of (<ref>) actually upholds the momentum principle. In fact, the gradient is homogeneous on _0, as seen by inspecting (<ref>). Pressure  is evaluated from

    ∫__0 N_2(_2,_2) F_2(_2, )  _2 + _0,

where, after an easy construction of the images, Neumann's function is 

    N_2(_2,_2)=-1/4 πlog(( cosh(x_-)-cos(y_-) ) ( cosh(x_-)-cos(y_+) ) ).
 



In parallel to the linear profile, we seek a percolation of the form (+α, + /β) on the parabolic mean flow. Given the class of initial data (<ref>), the expression (<ref>) remains unchanged for the description of the diffusive perturbations (_2,t). Now, the forcing function for the pressure reads

    F_2=- A_0 ( α(  S_2[û] _2[_x]+S_2[] _2[_y] ) +2(1 -2y) ∂_x  ),

while the pressure has the identical form of (<ref>).

By concise computations, our analysis of the primitive Navier-Stokes dynamics shows that there exist no definitive disturbed states in each of these two shear flows, for any choice of ω_0 of (<ref>). Once initiated, the basic profiles are completely outweighed by the imposing arbitrary perturbations which may be measured as infinitesimal or finite, albeit vaguely. The concept of disturbing established steady flows in the primitive setting is poised to be paradoxical.




 §.§ Ill-posed initial-boundary value problem



Instead of calling _0 and _0 of (<ref>) the perturbations, we set them as the initial data. Then one set of velocity solution  has the components

    (x,y,t)=_0(x,y)α̂(t),   (x,y,t)=_0(x,y)α̂(t),

where α̂ is finite, regular but arbitrary, and α̂(0)=1. Evidently, the flow-field remains solenoidal for all times t > 0 for every time-dependent function α̂. The following specifications of the pressure gradients ensure that the planar primitive equations are satisfied,

    ∂_x     = -ρ( _0α̂' + α̂^2 (_0 ∂_x _0 +_0 ∂_y _0) - α̂ν (∂_xx + ∂_yy) _0)  ),
    ∂_y     = -ρ( _0α̂'  + α̂^2 (_0 ∂_x _0  +_0 ∂_y _0)  - α̂ν (∂_xx + ∂_yy) _0)  ),

for given ω_0.

Our result does not contradict the global regularity of the two-dimensional equations established by Leray (1934b) because his proof of the well-posedness is based a maximum principle which bounds the scalar vorticity in ^2. In this respect, the primitive formulation has been avoided.




 §.§ Spurious instability by normal-modes



In the present section, we generalise the two basic uni-directional flows so that they depend on x_3 only, (x_3), and they are bounded between two infinite planes:

    _∥:     -∞ < x_1,x_2 < +∞,     0 ≤ x_3 ≤ 1.

By rewriting the basic profiles as =x_3, and =x_3(1-x_3), we look for connections of the preceding analysis to the linear stability theory which postulates that any viscous flow may be decomposed into a basic profile plus perturbation

    =+',   p=+p'.

The normal-modes analysis in the instability theory is concerned with the excitation and attenuation of the wave-like perturbations

    '(x_1,x_2,x_3,t)∝ ϕ(x_3,t) exp( α̃ x_1 + β̃ x_2   ),

where α̃ and β̃ stand for the respective wave numbers, and =√(-1). Only the real part of the complex representation in the eigen space is relevant to the motions. These wave-trains are known as instability waves which, presumably, exist and develop in space or in time, depending on whether the wave numbers, α̃,β̃, are complex or real.
 
In the domain of the parallel planes _∥, the complete Navier-Stokes equations (<ref>) are invariant under the space translation in either the x_1-direction or the x_2-direction or both, namely

    (x_1,x_2,x_3,t)=(x_1-a_1,x_2-a_2,x_3,t).

It follows that the invariance principle must also apply to the combined flow, +'. Because the basic flow itself is translation-invariant along the x_1,x_2 axes, the disturbances must preserve their analytical equivalence:

    A ϕ(x_3,t)   exp( α̃ x_1 + β̃ x_2  ) ≡  
    A  ϕ(x_3,t) exp(  α̃ (x_1-a_1)+β̃ (x_2-a_2)  ),

where A is the proportional constant, and is independent of x_1 and x_2. For arbitrary choice of the translation scales a_1 or a_2, we must have 

    α̃=0,      β̃=0,

for consistency of (<ref>). At the planes, the boundary conditions such as the no-slip, are applied to the eigen part ϕ. The non-existence deduction still holds regardless the presence of solid surfaces. As the perturbations are solenoidal, ∇·'=0, the wave-free constraints (<ref>) show that ∂ϕ/∂ x_3=0, which imposes strong restrictions on the eigen functions. Furthermore, the principle (<ref>) applies equally well to inviscid flows which are assumed to exist by formally setting ν=0 in (<ref>). Hence, the postulated sinuous instability waves (<ref>) are ill-defined in the framework of the exact incompressible Navier-Stokes-Euler dynamics, as their introduction implies incompatible physics.[The invariance law (<ref>) also rules out the following mean profiles which depend linearly on the co-ordinates, 

    _i(x_j,t)=σ_ij(t)x_j+v̅_i(t),    x_j ∈,    (i,j=1,2,3).

A corresponding pressure can be found so that the equations of motion are satisfied. In ^3, the proposed mean flows contain an infinite amount of energy for bounded σ's. For arbitrary , the principle demands

    _i(x_j,t)-_i(x_j-a_j,t)=σ_ij(t)a_j=0,

which stipulates the discriminant: spatial strain rates σ_ij≡ 0 at every time t. In rigorous fluid dynamics, we never recommend any exploration of spatially periodic perturbations on these profiles, as nothing insightful can be gained. In any event, these pathological flows would motivate false interpretations on fluid flows, such as the `generalised Kelvin modes'. There have been some confusions about Kelvin's paper (1887) in technical literature. Kelvin had worked exclusively on bounded uni-directional plane Couette flow between 0 ≤ y ≤ b, though he mentioned an extension, b →∞, on p.191 of his paper. Apparently, it was Hopf (1952) who first considered the boundary-free turbulence (see his 7) as a limit of spatially periodic flow. However, he expressed strong reservations about the unbounded behaviours at large distances, when the period tends to infinity. Given our current knowledge, we assert that, as a result of careful analyses, there are no infinite-energy viscous motions. Consequently, the invariance principle is an equivalent statement to the non-existence of the Navier-Stokes solutions in periodic domains.]




 §.§ Unpredictability due to disturbances



By the decompositions (<ref>), the superimposed flow of a basic profile and a perturbation satisfies the equations of motion

    ∂_t ' + ('·∇) '+ (·∇)' + ('·∇) = - ∇ p' /ρ+ νΔ',   ∇·'=0,

because the basic flow, (,)(,t), is a solution. As in many practical applications, it has been claimed that equations (<ref>) can be solved by numerical methods, for given initial data. Let us call these (primitive) solutions (,), and they are valid over time interval [t_0,T^*]. 

Choose a sub-interval [t_0,T], T<T^*. Divide the time interval [t_0,T] into n equal parts: Δ t = (T-t_0)/n or t_k= kΔ t,  k=1,2,⋯,n. We consider the first percolation on the mean flow field at t_0,

    (+α^(1)^(1), +^(1)/β^(1)).

Let the operator H()=(∂_t  -  νΔ ).
Then the percolation is determined by solving the heat equation with a buffer source

    H(^(1))= ,

subject to (arbitrary) solenoidal initial data, _p. Because of the sandwiched geometry, there are useful symmetries in the Green's function H_3 which is governed by 

    H(H_3)=-δ(x-x')δ(y-y')δ(z-z')δ(t).

Taking into account the non-slip on ∂_∥, we find that H_3 has a composite form

    H_3(, ',t)= H_1(x,x',t) × H_1(y,y',t) × H_f(z,z',t).

Thus the solenoidal diffusion is guaranteed,

    ^(1)(,t)=∫__∥ H_3 (-',t) _p(') '+∫_t_0^t_1∫__∥ H_3 (-',t-s) (') '  s,

where _p ∈ L^1_loc(_∥), for any choice of the source  within the constraints imposed on the buffer functions (<ref>). One of the important properties of the kernel H_3 is that it is capable of mollifying a large class of the initial data. In experimental work, various methods have been developed to disturb the mean flow, including wave generators, acoustic excitations by loud speakers, flow impulses, controlled directional jets, and others. All these artificial excitations are well-covered in the theory of L^1_loc space. Computationally, the diffusive component ^(1) can be recognised as the numerical errors from discretisation.

The first equality of (<ref>) is fully maintained if

    ∇^(1) = -  q^(1)(^(1),),

where the non-linear terms are grouped into 

    q^(i)(,) = ρα^(i)β^(i)( α^(i) (·∇) + (·∇) + (·∇) ) + ρβ^(i).

The actual  can be found by solving the problem
Δ^(i) = - (∇· q^(i)),
subject to inhomogeneous Neumann data of (<ref>). Some algebra gives

    ^(i)()   =∫__∥ N_3(,') (∇· q^(i))(') '  
           + ρβ^(i)∫_^2( N_3(,_2',0)(_2',0) -N_3(,_2',1)(_2',1) ) _2'+ _0,

where

    N_3=1/4 π∑_k=-∞^∞( 1/R_++1/R_- )

with the image locations at R^2_±= (x-x')^2+(y-y')^2+(z± z' - 2k)^2.
Hence our first companion solution deviates the mean motion

    (_c,p_c)(,t)=( +α^(1)^(1), +^(1)/β^(1) )(,t).

In general, the postulated solution (,) and the companion do not coincide over [t_0,t_1], even though both are perturbed by the same initial data _p.

Next we consider a series of consecutive percolations

    (α^(k)^(k), ^(k)/β^(k)),   k=2,3,⋯,

on the companion flow

    (^(k-1)_c,p^(k-1)_c)=(+∑_i=1^k-1α^(i)^(i), +∑_i=1^k-1^(i)/β^(i)),

over t_k-1≤ t ≤ t_k.
The k-th percolation is found from H(^(k))=, subject to data _c^(k-1), and it has the form (<ref>). By the same token, we must have 

    ∇^(k) = -  q^(k)( ^(k),^(k-1)_c ).

Over the last interval [t_n-1,t_n], the basic profile (,) has an arbitrary companion

    (^(n)_c,∇ p^(n)_c)(,t)=(+∑_k=1^nα^(k)^(k), ∇+∑_k=1^n∇^(k))(,t).

By increasing the number of sub-intervals, n →∞, we envisage an experimental condition of continuously introducing perturbation _p into the mean flow. A typical example of this situation is found in all wind tunnels where the free-stream contains high levels of turbulence or acoustic noise. Another example includes the macroscopic realisation of thermal fluctuations or Brownian motion which is an inherent part in all flow states. Evidently, the claimed primitive solution, (,∇), can never be unique over [t_0,T], in spite of the origin of the perturbations. By deduction, the primitive theory implies that any mean flow will be strongly distorted immediately after the start. This theoretical account does not reconcile with experimental facts. 

Our evaluations so far have been carried out by exact computations. We have kept all the terms in equation (<ref>). In linear stability theory, it has been assumed that the term, ('. ∇) ', is small compared to other terms in (<ref>). Neglecting this non-linear term, we have 

    ∂_t ' + (·∇)' + ('·∇) = - ∇ p'+ νΔ',   ∇·'=0.

In principle, the linearisation is made at the expense of the momentum conservation; there are no justifications for ignoring the non-linearity unless one tolerates absurd physics. Nevertheless, let the solution of (<ref>) be (_L,_L) over its validity interval [t_0,T]. This basic flow may have been obtained, for instance, by a normal-modes analysis of linear instability or a procedure of asymptotic expansion at high Reynolds numbers. 
Starting with the first percolation (<ref>), H(_L^(1))= remains unchanged. Clearly, the linearisation does not modify our diffusion problem for _L^(k). Only minor adjustments are required for the compensating pressure _L. We need

    q'^(i)(,) = ρα^(i)β^(i)(  (·∇) + (·∇) )+ρβ^(i),

instead of q. Without repeating the whole analysis for consecutive percolations, we deduce that the `linearised' companion to (_L,_L) contains the similar arbitrariness

    (^(n),∇ q^(n))(,t)=(_L+∑_k=1^nα^(k)_L^(k), ∇_L+∑_k=1^n∇^(k)_L)(,t),

where
∇^(k)_L = -  q'^(k)(_L^(k),^(k-1)).
As in the case of the upper half space, infinitely many consecutive percolations are able to saturate any basic profile with the aid of selective α's and β's. What we are trying to emphasise is the fact that, the connection of the sinusoids (<ref>) to fluid motions are only found in the linearised theory. The key issue is that these wave-like instabilities comprise just a small subset of the allowable L^1_loc perturbations which take any physically-relevant forms, for instance, spatiotemporal aperiodic fluctuations. To justify specific normal modes by evoking Fourier transforms is a half-done argument. Indeed counter examples to normal-modes analyses are well-known, as there are flows whose spectra are ostensibly incomplete. In such circumstances, the expansion theory of instability modes is at least inapplicable, as the continuum flow scales do not amend to wave descriptions.

It is an inconvenient reality that our description of fluid motions solely in terms the primitive variables precipitates the inevitability of the momentum partition, which gives rise to successive percolations under the auspices of the buffer . The grave consequence is that any disturbed state, (+', +p'), amounts to indeterminacy. Hence, those methods of solving the primitive (<ref>), such as linear or weakly non-linear instability analyses, non-modal optimised energy approaches, Fourier or eigen-function expansions, high Reynolds-number asymptotic approximations, can hardly be promising.




 §.§ Boundary layers



A thin viscous layer develops on the surface of a plate when it is suddenly set into motion. Theoretically, the thickness of the layer can be estimated from the impulsively-started Stokes' flow on an infinite plane, 

    ∂_t u-ν∂_yy u=0,

subject to u_0(y,0)=0, u(0,t≥0)=U_∞, and u(∞,t)=0. The solution of this equation is given by

    u(y,t)=U_∞∫_0^t  T(t-s) s = y U_∞/√(4 πν)∫_0^t exp(-y^2/(4ν(t-s))) /(t-s)^3/2 s =y U_∞/√(4 πν) I(t),

where, from (<ref>),

    T(t-s)=∂( H_1(y,y',t-s) - H_1(y,-y',t-s))/∂ y'|_y'=0.

The integral, I(t), is a measure of time scale. A fluid particle on the plate impulsively moves with the plate to a location L in time L/U_∞=t_L. In fluids ν≪ 1, the exponential function is ∼ O(1) for y^2/(ν t) ≪ 1. Asymptotically, t_L is the characteristic time 

    I(t) ∼∫_0^t (t-s)^-3/2 s ≈ t^-3/2∫_0^t(1 - 3s/(2t))  s ≈ t_L^-1/2 = √(U_∞/L).

Close to the plate, u/U_∞∼ O(1), the wall layer, δ∼ y, is called a boundary layer. For an observer moving with the plate, the thin viscous shear  

    δ∼ O(√(ν))           δ/L∼√(ν/L U_∞) = Re^-1/2_L,

for Reynolds' number Re_L ≫ 1.

The boundary layer theory of Prandtl (1904) is an approximation of the planar Navier-Stokes dynamics. For a steady boundary layer over a semi-infinite flat-plate, the governing equations are 

    ∂ u/∂ t+u∂ u/∂ x + v∂ u/∂ y = - 1/ρ∂ p/∂ x + ν∂^2 u/∂ y^2,      ∂ u/∂ x +∂ v/∂ y=0.

The boundary conditions are u(x,0)=v(x,0)=0; u → U_∞, v → 0,y →∞,
and upstream compatibility u(0,y)=U_∞, v(0,y)=0. In the scaling U^2_∞/L, v-momentum is of small order, and ∂ p/∂ y ∼ O(δ/L). Since there are three unknowns in (<ref>), a supplementary condition must be imposed. It is argued that, at the edge of the layer, velocity u presumes the value of the inviscid free-stream U_∞. Thus the stream-wise pressure satisfies

    p/ x = - ρ U_∞  U_∞/ x.

Effectively, the gradient in (<ref>) is fixed for every given U_∞=U_∞(x). The well-studied case is a Blasius profile on a flat-plate, U_∞=U_0= In steady flows, a number of free-streams have been proposed, for example, U_∞(x) = U_0 (x/L)^m for an integer m. In the flow over a circular cylinder, the initial flow and free-stream are matched, u_0(x,0)=U_∞(x)=U_0sin(x/L), and hence viscous layers develop in the vicinity of the forward stagnation point. The expectation is that, from the solutions of (<ref>), boundary layers in favourable and adverse gradients may be elucidated. However, the approximations (<ref>) are a reductionist's view of fluid motion and are not always reasonable. For instance, elevating passive pressure (<ref>) to regulative control (<ref>) introduces a misconception that an adverse pressure gradient causes a boundary layer to separate.

Incidentally, the supposition of a fixed gradient circumvents the non-uniqueness of the primitive formulation, as it closes the partition loop of the momentum. But, considerations of perturbing a Blasius flow _BL, in the sense of (<ref>), necessarily involve an analogous analysis to (<ref>) or (<ref>). In practice, a typical flat-plate boundary layer has a stream-wise length L ∼ O(1) behind the leading edge at moderate Reynolds numbers (cf. the estimate (<ref>)). Comparisons between measurements and linearised calculations may seem to be impressive. Even so, we must treat any of the apparent successes with caution, as these limited test realisations by no means elude the analytical indeterminacy. Accordingly, the contrast between convective and absolute instabilities is inconsequential.[Our findings of the present section are consistent with experiment and careful computations. The classic experiment on the boundary-layer transition by Schubauer & Skramstad (1947) conclusively demonstrated that perturbations have no correlations with the measured transition location. Instead of the primitive variables, solutions must be sought in the continuum theory of the full vorticity dynamics. In applications, we do not regard unqualified ripples and stitches in vorticity iso-contours as evidence of the flow instabilities. Yet the critical layers are ill-defined and hence dissolved. Highly skewed vortices haul retarded flow elements away from solid surfaces, leading to separations. Convergence in discretisation and successive asymptotes is of importance. We exemplify these issues by a dipole-wall simulation, see figures 7 to 11 of Lam (2018a).]




§ CYLINDRICAL POLAR CO-ORDINATES



Denote the velocity components by =(u_r,u_θ,u_z) in the cylindrical co-ordinates =(r,θ,z). The primitive equations of motion are

    ∂ u_r/∂ t+ (·∇)u_r - u^2_θ/r    =-1/ρ∂ p/∂ r+ ν(Δ_1 u_r -2/r^2∂ u_θ/∂θ),
    ∂ u_θ/∂ t + (·∇)u_θ + u_r  u_θ/r    =-1/ρ  r∂ p/∂θ + ν( Δ_1 u_θ + 2/r^2∂ u_r/∂θ),
    ∂ u_z/∂ t + (·∇)u_z    =-1/ρ∂ p/∂ z +  ν Δ_0 u_z,

and the continuity,

    ∇'· = ∂ u_r/∂ r+u_r/r+1/r∂ u_θ/∂θ + ∂ u_z/∂ z=0,

where the operators are defined by Δ_1=Δ_0-1/r^2,  

    Δ_0 = ∂^2/∂ r^2+ 1/r∂/∂ r + 1/r^2∂^2/∂θ^2+ ∂^2 /∂ z^2.

The gradient operator, ∇ = ∂/∂ r +∂/(r∂θ)+∂/∂ z.
Equations (<ref>) and (<ref>) are made dimensionless by the SI units. To derive the analogous expression of (<ref>) for the pressure, we take into account the fact that the time-dependent terms must vanish in view of the incompressibility hypothesis. In practice, we look for an expression containing Δ_0p or Δ_1 p. The following identities facilitate the subsequent algebraic simplifications:

    ∂/∂ r (Δ_0 u_r)   = Δ_0( ∂ u_r/∂ r )-1/r^2∂ u_r/∂ r - 2/r^3 ∂^2 u_r /∂θ^2 ,
    1/r (Δ_0 u_r)   = Δ_0( u_r/r)+2/r^2 ∂ u_r/∂ r -u_r/r^3 ,
    1/r∂/∂θ(Δ_0 u_θ)   = Δ_0(1/r∂ u_θ/∂θ)+2/r^2∂^2 u_θ/∂ r ∂θ-1/r^3∂ u_θ/∂θ,

where each term of the left-hand side can be found by manipulating the governing equations. By direct computations, the momentum equations in (<ref>) then show that the diffusion-related terms drop out due to the continuity.
Effectively, we need to simplify the expression

    (∂/∂ r+1/r)((·∇)u_r-u^2_θ/r) + 1/r∂/∂θ((·∇)u_θ+u_r u_θ/r)+ ∂/∂ z((·∇)u_z),

which equals to -Δ_0 p/ρ.
It is found that the lengthy identity can be paraphrased by subtracting a variant of the continuity, (.∇)(∇'·),
so that all the terms containing second derivatives cancel out.
After further reductions and some rearrangements, we obtain the Poisson equation for the pressure

    - 1/ρΔ_0 p     = (∂ u_r/∂ r)^2+( 1/r∂ u_θ/∂θ+u_r/r)^2+(∂ u_z/∂ z)^2
          + 2 ( 1/r∂ u_r/∂θ∂ u_θ/∂ r + 1/r∂ u_θ/∂ z∂ u_z/∂θ + ∂ u_r/∂ z∂ u_z/∂ r- u_θ/r∂ u_θ/∂ r )=P().





 §.§ Non-existence of Burgers' vortex



Flows in the domain,

    _z:     0 ≤ r < ∞;     0 ≤θ < 2 π;    - ∞ < z < +∞,

must obey the invariance principle (<ref>) in the z-direction. Burgers (1948) suggested a vortex motion whose velocities are 

    u_r=- γ r;   u_θ=ϕ(r,t);   u_z= 2 γ z,

where γ denotes a constant, representing a strain rate as it has the dimension [^-1]. The invariance law in the z-direction is expressed in the identity,

    u_z(z)-u_z(z-a_3)= 2 γ a_3 =0,

for arbitrary translation a_3. Thus, Burgers' vortex cannot exist for non-zero γ. Few variants of (<ref>), for instance, γ=γ(t) or γ=γ f(r), have been proposed. These revised vortex flows set the scene for a number of strained vorticity fields, which were treated as building elements in simulating turbulence. By inference, the conclusions of these works have marginal validity, and large parts of them, if not all, must be reappraised for important applications. In finite energy flows, the rule is that the Navier-Stokes equations do not admit solutions of power series 

    u_i ∝ (x_i-a_i)^k,   (i=1,2,3)

in homogeneous direction x_i for non-zero integer k.




 §.§ Rotational normal-modes for Hagen-Poiseuille flow



In the infinite circular pipe of unit radius,
_1: 0 ≤ r ≤ 1;  0 ≤θ < 2 π;  z ∈,
we consider a uni-directional flow _z=_z(r) with zero radial and circumstantial pressure gradients. The continuity is satisfied automatically, and the z-momentum of (<ref>) reduces to μ (∂_rr +∂_r/r)_z = ∂q̅/∂ z. Applying the no-slip on wall r=1, we obtain the steady velocity, 

    =(_r, _θ,  _z
    	)=( 0,  0,  1-r^2 ),

which is driven by an axial pressure gradient, ∂_z q̅=-4 μ. In laboratory experiments, the steady-state profile (<ref>) may be realised at locations a few diameters downstream of a streamlined entry, as long as the transient of starting conditions becomes settled. In this way, the comparison with the theory tentatively removes the dependence of initial conditions. Given the subtleties and difficulties in fluid experiments, this unitary case should never be generalised into a conception that initial conditions are of secondary importance. Practically, we may devise many entry conditions, and hence downstream flows are strong functions of the applied external pressure gradient. 

In theory, parabolic distribution (<ref>) is applicable over the entire domain _1 and is thus treated as a mean flow. In parallel to the wave-form (<ref>), stability analysis considers possible amplifying effects of the disturbances 

    (u'_r,u'_θ,u'_z,p')(r,θ,z,t)∝(ũ_r,ũ_θ,ũ_z,p̃)(r,t) exp(  ( γ̃ z + n θ )  ),

where n=0,±1,±2,⋯, so as to ensure rotational periodicity. The principle of translation invariance asserts that the dynamic compatibility must be observed at all times, or the following expression must make sense

    exp(  ( γ̃ z + n θ )  ) 
    = exp(  ( γ̃ (z-a_3) + n θ )  ).

In the exact non-linear Navier-Stokes dynamics, the travelling waves (<ref>) cannot exist for arbitrary a_3 unless γ̃=0. As a result, the expected perturbed structures are less interesting, as they degenerate into fixed periodic functions in θ. 

By integration, the actual pressure driving the flow (<ref>) is found to have an analytic expression q̅=-4 μ z + q̅_0(t). From a strictly analytical point of view, we may articulate that there can be no pressure of this form which satisfies the invariance law. However, we shall ignore this indeterminate quantity, and stick to the fact that the momentum conservation involves only the gradient.




 §.§ Ambiguity in pipe flow



Let us introduce an add-on flow, , superimposed on flow (<ref>), 

    (,q̂)(r,θ,z,t)=(_r,_θ,_z,q̂)(,t),

where the velocity vanishes at the pipe wall, (1,θ,z,t)=0, and decays to zero as z →±∞. We substitute the sum, 
+α, q̅+q̂/β,
into equations (<ref>) and (<ref>). Since (,) is a solution pair, the analogous system to (<ref>) reads

    ∂_t _r + f_r    =-(∂/∂ r)/A_0+ ν(Δ_1 _r -2r^-2∂_θ/∂θ),
    ∂_t _θ + f_θ   =-(∂∂θ)/(r A_0) + ν( Δ_1 _θ + 2r^-2∂_r/∂θ),
    ∂_t _z + f_z    =-(∂/∂ z)/A_0 + ν Δ_0 _z, 
    ∇'.   =0.

where the non-linear terms are grouped into

    f_r     =α (·∇)_r + (1-r^2) ∂_r/∂ z - ^2_θ/r,
     
     f_θ    = α(·∇)_θ + (1-r^2) ∂_θ/∂ z + _r _θ/r, 
    
      f_z     = α(·∇)_z + (1-r^2) ∂_z/∂ z - 2r _r.

The Neumann data can be inferred from (<ref>), the equation gives the pressure 

    Δ_0q̂ = -A_0 (∂ f_r/∂ r+f_r/r + ∂ f_θ/∂θ/r+ ∂ f_z/∂ z)=-Q̂_p(, ,).


Choose any buffer =(b_r,b_θ,b_z)(,t), as postulated. The components of the percolation are diffusions

    h_1(_r)=-2r^-2∂_θ/∂θ+b_r,   h_1(_θ)=2r^-2∂_r/∂θ+b_θ,   h_0(_z)=b_z,

where h_0=∂_t -νΔ_0, and h_1= h_0+1/r^2.
Because the Hagen-Poiseuille flow is steady, we are free to set the switch-on time of percolations, say, t=0. The initial data are

    ( _r,_θ,_z )(,t=0)=( _0, _0, _0 )().

The solution of the third equation in (<ref>) is expressed as

    _z(,t) =∫__1K_0(,',t)   r' _0(')  '+ ∫_0^t ∫__1K_0(,',t-s)   r'b_z(',s)  '  s,

where K_0(,',t)=H_1(z,z',t)G_0(r,θ,r',θ',t), with the cross-section Green function being

    G_0    =1/π ∑_i=1^∞1/(J_1(κ_0i))^2  J_0(κ_0i r)   J_0(κ_0i r')  exp(-κ^2_0iν t) 
          + 2/π∑_n,i=1^∞cos(n(θ-θ'))/(J_n+1(κ_ni))^2 J_n(κ_ni r) J_n(κ_ni r') exp(-κ^2_niν t),

where J_n stands for the Bessel function of order n=0,1,2,⋯. Each of the Bessel functions, J_n(x), has infinitely many zeros in the interval x>0. In every summation over n, κ_ni is the ith positive root of the function J_n(κ)=0. 

Let u_±=_r ±_θ. The first two equations of (<ref>) are reduced to the system

    .
    
    	(∂_t     - νΔ_1 ) u_+= 2ν ∂_θ u_-/r^2+b_r +b_θ,
    
    	(∂_t     - νΔ_1 ) u_-= -2ν ∂_θ u_+/r^2+b_r - b_θ.
    	}

Since the diffusion operator on the left is separable, we readily establish the analogous kernel, K_1(,',t)=H_1(z,z',t) G_1(r,θ,r',θ',t), where the Green function G_1 is

    G_1(r,θ,r',θ',t)=2/π∑_n,i=1^∞cos(n(θ-θ')) /(J_n+1(κ_ni))^2  J_n(κ_ni r)   J_n(κ_ni r')  exp(-κ^2_niν t).

Transforming (<ref>), the solutions are defined in the coupled integral equations:

    .
    
    u_+(,t)    + 2 ν∫_0^t∫__1 K_2(,',t-s)/r'  u_-(',t-s) '  s = F_+(,t), 
    
    u_-(,t)    - 2 ν∫_0^t∫__1 K_2(,',t-s)/r'  u_+(',t-s) '  s = F_-(,t), 
    }

where the integral terms are obtained by integration by parts, taking into account the no-slip condition. Kernel K_2(,',t)=H_1(z,z',t) G_2(r,θ,r',θ',t). The drivers on the right, F_±, are known

    ∫__1 K_1 (,',t)r'   (_0 ±_0) '+ ∫_0^t ∫__1 K_1 (,',t-s)r'   (b_r ± b_θ) '  s.

In view the rotational symmetry of G_1, the modified Green function is easily derived

    G_2(r,θ,r',θ',t)=2/π∑_n,i=1^∞nsin(n(θ-θ')) /(J_n+1(κ_ni))^2  J_n(κ_ni r)   J_n(κ_ni r')  exp(-κ^2_niν t).


By virtue of the asymptotic properties of the Bessel functions (n≥1), we see that the kernel K_2/r is well-behaved at the pipe centre, as

    n J_n(r) /r ∼r^n-1/(2 (n-1)!)        r → 0.

This evaluation ensures that the flow is regular at the axis of symmetry. It follows that both the integral kernels, K_1 and K_2, are bounded in the whole flow domain. By a standard algebraic manipulation to eliminate (say) u_- from the equations (<ref>), the resulting equation becomes a linear Volterra equation of the second kind with a regular kernel. Its solutions can be established by successive approximations of resolvent kernels, see, for instance, Tricomi (1957). Without going into the technical details, we deduce that system (<ref>) is uniquely solvable, and the diffusion (<ref>) admits solutions, (_r,_θ, _z), for every given data (<ref>) and buffer . 

As in the derivation of Poisson equation (<ref>), system (<ref>) always defines a solenoidal velocity field for any incompressible initial data, i.e., the last equation in (<ref>) holds at all times. Any momentum imbalance caused by the diffusions is restored, because the percolation gradients are fixed by the receptive non-linearities

    ∂_r q̂ = - ραβ f_r- ρβ b_r,  ∂_θq̂ /r = - ραβ f_θ-ρβ b_θ,  ∂_z q̂ = - ραβ f_z - ρβ b_z,

where the right-hand sides are known functions of  and .

If necessary, the pressure (<ref>) can be calculated with extra effort. By virtue of the solutions of (<ref>), the Neumann condition is determined and denoted by

    ∂q̂/∂ r|_r=1 = -ρβ b_r(1,θ,z)=Q_0(θ,z).

Both  and  are assumed to decay at |z|→∞.
So the integral for the pressure is

    q̂(r,θ,z) =∫__1N_p(,')  r' Q̂_p(')  ' + ∫_∫_0^2 πN_p(,1,θ',z') Q_0 θ'  z' + q̂_0,

where the Neumann kernel is governed by Δ_0 N_p = - δ(r-r')δ(θ-θ')δ(z-z')/r, subject to homogeneous data. Starting with the expression for Dirac's delta,

    δ(θ-θ')=1/2 π+ 1/π∑_n=1^∞cos(n(θ-θ'))=1/2 π∑_n=-∞^∞cos(n(θ-θ')),

it is straightforward to derive the product for the cross-section Neumann function

    δ(r-r')δ(θ-θ')/r=1/π∑_n=-∞^∞∑_i=1^∞J_n(σ_ni r)J_n(σ_ni r')/(1-n^2/σ^2_ni) (J_n(σ_ni))^2cos(n(θ-θ')),

where σ_ni is the ith positive root of the function J'_n(σ)=0, so as to satisfy the boundary condition at the wall. By inspection, N_p must have the form

    N_p(,')=∑_n=-∞^∞∑_i=1^∞N_niJ_n(σ_nir)cos(n(θ-θ')).

In view of the pressure decay at large |z|, we determine N_ni by solving the governing equation for N_p. The result is

    N_p= 1/π∑_n=-∞^∞∑_i=1^∞   J_n(σ_ni r)J_n(σ_ni r')/(σ_ni-n^2/σ_ni) (J_n(σ_ni))^2 ×
         cos(n(θ-θ'))exp(-σ_ni|z-z'| ).

In (<ref>), the gauge _0 does not affect the evaluation of gradients ∇.

There are no restrictions on the cycles of consecutive percolations. An interval, [0,T], is sub-divided into n equal parts. In each of these sub-intervals, we repeat our calculations for the local percolations. Over the last interval [t_n-1,t_n], the Hagen-Poiseuille profile (,) has companions

    (+∑_k=1^nα^(k)^(k), ∇+∑_k=1^n∇^(k)/β^(k))(x,y,z,t).


The critical issue is that the complete Navier-Stokes equations admit non-linear as well as elementary solutions, though the class of the present solutions is not related to the most interesting phenomenon: turbulence. Our analysis demonstrates the fact that the companions decimate the parabolic Hagen-Poiseuille profile immediately after the commencement of the motion, as the percolations are multivalued at any time. Strictly, there does not exist one single identifiable characteristic velocity in the whole flow-field. Thus, the notion that turbulence in the pipe experiment of Reynolds (1883) is triggered by some form of finite-amplitude perturbation must be fundamentally flawed, as there are no well-defined Reynolds' numbers. Likewise, it is a truism to say that the laminar-turbulent transition in pipe flow spans a complex, non-linear, spatiotemporal process involving an enormous range of space and time scales, if such a rhetoric merely presumes the primitive particulars. In conclusion, flow instability and pleonastic bi-stability have no precise meanings.




§ COUETTE MOTION BETWEEN CYLINDERS



The flow between two rotating cylinders about a common axis has been studied extensively in fluid dynamics. In this section, we consider time-independent flows in the annulus between two rotating co-axial cylinders,

    _i:   0< R_1 ≤ r ≤ R_2,   0≤θ < 2 π,    - ∞ < z < + ∞,

where R_1 and R_2 are the radii of the inner and outer cylinders.
Denote their angular velocities by Ω_1 and Ω_2. The simplicity of the geometry suggests that the flow is characterised by the non-zero component _θ=_θ(r) such that the incompressibility (<ref>) is fulfilled. The governing equations (<ref>) reduce to

    ^2_θ/r=∂_rq̅/ρ,      ν (∂_rr +∂_r/r-1/r^2)_θ = 0,

respectively. The boundary conditions for real fluids (ν >0) are

    _θ(R_1,θ,z)=R_1Ω_1,   _θ(R_2,θ,z)=R_2Ω_2.

Algebra shows that the velocity is the sum of a rigid-body rotation and a viscous shear layer:

    =(  0,    Ar+B/r,   0 ),

or its angular velocity:

    Ω̅=(  0,    A+B/r^2,   0 ),

where the integrating constants are expressed in symmetric forms

    A=(Ω_2 R_2^2- Ω_1 R_1^2)/(R_2^2-R_1^2),   B=-(Ω_2 - Ω_1) R_2^2R_1^2/(R_2^2-R_1^2).

The mean pressure or its gradient is also known. Because of the no-slip condition,  the cylinders' walls drag the fluid into axi-symmetric rotational motion. 




 §.§ Fallacy of Rayleigh's instability criterion



By formally putting ν=0 in (<ref>), we get 

    D_*/D_*t_r-^2_θ/r= -1/ρ∂/∂ r,     D_*/D_*t (r_θ) = 0,     D_*/D_*t_z=-1/ρ∂/∂ z,

where D_*/D_*t is the total derivative with axi-symmetry. The continuity reads,

    ∂_r/∂ r+_r/r+∂_z/∂ z=0.

The governing equations remain unchanged under z-translation. From the middle momentum of (<ref>), Rayleigh (1916) argued that the circulation, r_θ, or the angular momentum, r^2Ω̃, is time-invariant. In the specific case of purely azimuthal flow, _θ = _θ(r,t) where _r = _z = 0 or _r = 0, _z = f(r,t), the θ-momentum reduces to (r_θ)/ t=0. The kinetic energy per unit mass of the motion is E_θ=^2_θ/(2r^2)=(r^2Ω̃)^2/(2r^2). The stability of the motion depends on the sign of the discriminant,

    (r)=/ r(r^2 E_θ) = r^2Ω̃ / r(r^2 Ω̃).

By virtue of the conservation law, the constancy, r^2Ω̃(r,t)=r_0^2Ω̃(r_0,0), does not necessarily imply stability or instability; we must restore to perturbation analyses (see below).

Rayleigh contemplated that the invariance of r_θ still holds in three-dimensional motions with non-trivial radial and axial velocities, _r≠0 and _z≠0, provided the term ^2_θ/r in the first equation is balanced by a radial force. Then the whole motion would proceed as if _θ were absent. His invariance postulation is imprecise without supplementary boundary conditions. Because of z homogeneity, we consider only motions in planes normal to the axis, where _z is independent of z. From axi-symmetry, it is plain to see that continuity (<ref>) implies _r=α̃(t)/r with any finite α̃, and α̃(0)≠0. Then a solution of (<ref>) follows

    _r=α̃(t)/r,    r_θ=exp(-r^2/2+∫_0^t 
    α̃(s)  s),   _z=g(r,t),

as long as there are no boundary conditions. The pressure gradients are determined accordingly. Although the total derivative D_* (r_θ)/D_*t vanishes, it does not imply the time invariance of r_θ, because the circulation depends on the choice of function α̃(t) for fixed initial value _θ(r,0)=exp(-r_0^2/2)/r_0. In an infinitely long hollow circular tube, 0≤ r ≤ a, solution (<ref>) is singular and inadmissible for flows of finite energy. Similarly, it does not exist wherever there is a boundary, on which _r vanishes. At least, the axi-symmetry and the condition _r|_∂=0 are incompatible propositions, unless _r=0.

If the restriction on the axial symmetry is relaxed, Rayleigh's invariance on r_θ is unavailable. Let us return to the inviscid flow in the hollow tube, with the boundary condition, (a,θ,z,t)=0. We find an arbitrary solution to the governing equations (<ref>)-(<ref>) (ν=0),

    _r = (r-a) α̃(t)h(θ),  _θ=-(2r-a)α̃(t) ∫_0^θh(s) s,  _z=Z(r,h(θ))β̃(t).

We assume Z(r,h(θ)) are finite, or Z=0. Function h is periodic, h(θ)=h(2 π m + θ), where m=0,±1,±2,⋯. The time functions, α̃(t) and β̃(t), are arbitrary, smooth, and bounded, and α̃(0) ≠ 0, and β̃(0) ≠ 0, to avoid trivial initial data. In fact, the inviscid momentum conservation is fully respected, once time-dependent  is specified, as each of ∇ in (<ref>) is fixed. 

For the corresponding exterior flow (r≥ a), a comparable solution is given by

    _r = 1/r(1-d_0)α̃(t)h(θ),   _θ=-(k+1)(r-a)^k d_0α̃(t) ∫_0^θh(s) s,

where the decay at large r is controlled by d_0=exp(-(r-a)^k+1), k ≥ 1. If we consider  the particular planar motion in real fluids, i.e., _z=0 and ν>0, the wall no-slip is satisfied by both _r and _θ. In other words, for given viscosity ν, all the velocity terms in the planar Navier-Stokes equations (<ref>) are known so that the momentum conservation is maintained by the pressure gradients. Note that any suitable functions can be assigned to α̃(t) and h(θ) that, potentially, give rise to a family of disparate flow-fields. The novelty of the present example is that we do not have to rely on the use of the diffusive percolations to demonstrate the incompleteness of the primitive dynamics. 




 §.§ Arbitrary perturbations



In the annulus domain _i, the inviscid equations (<ref>) admit a steady mean flow (,)(r),

    _r=_z=0,   _θ=Ṽ(r)=r Ω̃(r),   (r)=∫(Ṽ^2/r)   r,

where the azimuthal velocity is an arbitrary function of r. Following the practice of perturbation theory, we consider the possibility that the mean flow is perturbed as a superposition

    '_r,  Ṽ(r)+'_θ,  '_z,  (r)+'.

Now, the continuity is governed by

    ∂'_r/∂ r+'_r/r+(∂'_θ/∂θ)/r + ∂'_z/∂ z=0.

The perturbation equations are written as

    ∂'/∂ r =-ρ (∂_t '_r + f'_r),  ∂'/∂θ = -ρ r (∂_t '_θ + f'_θ),  ∂'/∂ z= -ρ (∂_t '_z + f'_z),

where the full non-linearities are

    f'_r     =('·∇)'_r + Ω̃ ∂'_r/∂θ - '^2_θ/r - 2Ω̃'_θ,
     
     f'_θ    = ('·∇)'_θ + Ω̃ ∂'_θ/∂θ  + (Ṽ/ r+Ω̃) '_r + '_r '_θ/r, 
    
      f'_z     = ('·∇)'_z + Ω̃ ∂'_z/∂θ.

System (<ref>)-(<ref>) must be solved subject to the wall condition,

    '_r(R_1,θ,z,t)='_r(R_2,θ,z,t)=0.

By observation, one set of the velocities is

    '_r    = (r-R_1)(R_2-r) α(t) h(θ),
    '_θ    =(3r^2+2(R_2-R_1)r+R_1R_2)α(t)∫_0^θh(s) s,

and '_z has the identical form of (<ref>). Another flow-field reads

    '_r = (1/R_1-1/r) (1/r-1/R_2) α(t) h(θ),  '_θ =-1/r(1/r-1/R_1R_2)α(t)∫_0^θh(s) s.

The sum of (<ref>) and (<ref>) constitutes additional solutions. Given an inviscid Couette mean and the perturbations, equations (<ref>) indicate that the gradients are exactly known. Note that these solutions also satisfy the usual linearised equations, where the non-linear terms in the functions, f'_r,f'_θ, and f'_z, are neglected. It becomes clear that the gradient ∇' needs to be updated to reflect the approximations due to linearisation. The point is that our finite-energy velocity fields can be superimposed on any claimed solutions of (<ref>) which establish unstable properties of the mean. It follows that there is no instability as postulated by criterion (<ref>). 




 §.§ Circular shear between infinitely long cylinders



Couette motion (<ref>) or its inviscid counterpart (<ref>) of the unbounded cylindrical gap is translation-invariant in the z-direction. Hence the wave-trains of the specific form (<ref>) are not supported in these flows. 

To establish percolations (<ref>), we solve equations (<ref>) subject to initial data similar to those in (<ref>), and the non-slip condition on ∂_i.
Instead of Green's function G_0 in (<ref>), we need 

    g_0   = π/4∑_i=1^∞ F_0(λ_ni) E_0(λ_ni r) E_0(λ_ni r')exp(-λ^2_niν t)   + 
         π/2∑_n,i=1^∞cos(n(θ-θ'))F_n(λ_ni) E_n(λ_ni r) E_n(λ_ni r') exp(-λ^2_niν t),

where J_n and Y_n stand for the Bessel functions of order n, and F_n denotes the normalisation factor,

    F_n(λ_ni)=λ_ni^2(J_n(λ_ni R_2))^2/(J_n(λ_ni R_1))^2-(J_n(λ_ni R_2))^2,

and E_n the expansion functions

    E_n(λ_ni r)=J_n(λ_ni R_1) Y_n(λ_ni r) - Y_n(λ_ni R_1) J_n(λ_ni r).

In every summation over n=0,1,2,⋯, λ_ni is the ith positive root of the equation

    J_n(λ R_1) Y_n(λ R_2) - Y_n(λ R_1) J_n(λ R_2)=0.

Similarly, we obtain the analogous system to (<ref>) in terms of a regular kernel K_2, where G_2 is replaced by 

    g_2=π/2∑_n,i=1^∞n sin(n(θ-θ'))   F_n(λ_ni)   E_n(λ_ni r)   E_n(λ_ni r')  exp(-λ^2_niν t).

Furthermore, the integral kernel to replace G_1 is

    g_1=π/2∑_n,i=1^∞cos(n(θ-θ'))   F_n(λ_ni)   E_n(λ_ni r)   E_n(λ_ni r')  exp(-λ^2_niν t).

We obtain the gradients directly from the law of momentum conservation

    ∂/∂ r   = - ραβ ( α(·∇)_r-α^2_θ/r+(A+B/r^2 )( ∂_r/∂θ - 2 _θ )) - ρβ b_r, 
    ∂/∂θ   = - ραβ ( α  r (·∇)_θ + ( Ar+B/r ) ( _θ/∂θ +_r ) + ( A-B/r^2 )_r)- ρβ b_θ, 
    ∂/∂ z    = - ραβ (α(·∇)_z +( A+B/r^2 )∂_z/∂θ)- ρβ b_z.

In short, the companions have the form (<ref>).




 §.§ Perturbations to Taylor's linearised solution



If the disturbances are assumed to take a symmetric form about the axial direction, we drop all terms ∂/∂θ in the equations of motion. Taylor (1923) considered the superimposed flow of the Couette profile (<ref>) and linear axi-symmetric disturbances, +^T, +^T. The perturbation flow is governed by 

    ∂_t ^T_r     - 2 (A+B/r^2) ^T_θ  = νΔ_1^* ^T_r-1/ρ∂^T/∂ r, 
    ∂_t ^T_θ + 2A ^T_r     = νΔ_1^* ^T_θ,   ∂_t ^T_z = νΔ^*_0 ^T_z- 1/ρ∂^T/∂ z,

and by the continuity

    ∇^*·^T=∂^T_r/∂ r+^T_r/r+∂^T_z/∂ z=0.

The Laplacians Δ_1^*=Δ_0^*-1/r^2, and 
Δ_0^* = ∂^2/∂ r^2 + (∂/∂ r)/r + ∂^2/∂ z^2.
By making use the expressions in (<ref>) for axi-symmetric flow, we find

    Δ_0^* ^T = 2 ρ ( (   A+B/r^2 ) ∂^T_θ/∂ r + ( A-B/r^2 ) ^T_θ/r ).

In the narrow-gap approximation, Taylor sought the solutions of (<ref>) and (<ref>) by expanding the disturbances in terms of the eigen-functions,

    ^T_r=_r (r) cos(k z),   ^T_θ=_θ(r) cos(k z),   ^T_z=_z(r) sin(k z),

with boundary conditions ^T=0 at r=R_1,R_2. On the end-plates, z=0,2 π/k, the inviscid condition applies to ^T_z=0; ^T_r and ^T_θ are periodic in the axial direction.

Let the annulus be unit height,
_a:R_1 ≤ r ≤ R_2, 0≤θ < 2 π, 0 ≤ z ≤ 1.
Denote the corresponding Taylor's steady solutions by

    +^T==(U,  V,  W)(r,z),   +^T=P̅(r,z).

Now we consider a set of perturbations 

    (+α,  P̅+/β)(r,z,t).

On substituting (<ref>) into governing equations (<ref>), retaining the axi-symmetry and linearising, we get the following dynamic equations:

    ∂_t v_r +(·∇^s)v_r + (·∇^s)U - 2 V v_θ/r    = νΔ_1^* v_r-(∂/∂ r)/(ραβ), 
    ∂_t v_θ +(·∇^s)v_θ + (·∇^s)V + U v_θ/r+V v_r/r    = νΔ_1^* v_θ,
    ∂_t v_z +(·∇^s)v_z + (·∇^s)W    = νΔ^*_0 v_z-(∂/∂ z)/(ραβ),

where the gradient operator ∇^s=∂/∂ r+∂/∂ z, and ∇^*·=0. The no-slip condition on ∂_a reads

    (R_1,z,t)=(R_2,z,t)=0;   (r,0,t)=(r,1,t)=0.

We are interested in the initial-boundary value problem (<ref>) subject to data

    _s=(U_s,  V_s,  W_s)(r,z,t=0).


For any axi-symmetric buffer =(b_r,b_z)(,t), the percolations are defined by

    ∂_t v_r-νΔ_1^* v_r=b_r,       ∂_t v_z - νΔ^*_0 v_z=b_z.

The solution for the z-component is readily found to be

    v_z(,t)=∫__a K_a(,',t) r' W_s(')  ' + ∫_0^t ∫__a K_a(,',t-s) r' b_z(',s)  '  s,

where K_a(,',t)=H_f(z,z',t)M_0(r,r',t), and M_0 is given by 

    π/4∑_i=1^∞  F_0(λ_0i)   E_0(λ_0i r)   E_0(λ_0i r')  exp(-λ^2_0iν t),

(see (<ref>) to (<ref>) at n=0). Similarly,

    v_r(,t)=∫__a H_a(,',t)r' U_s(')  ' +∫_0^t ∫__a H_a(,',t-s) r' b_r(',s)  '  s,

where H_a(,',t)=H_f(z,z',t)M_1(r,r',t), and 

    M_1= π/4∑_i=1^∞  F_1(λ_1i)   E_1(λ_1i r)   E_1(λ_1i r') exp(-λ^2_1iν t).

In comparison with (<ref>), the forcing term on the pressure, -Δ_0^*  /(ραβ), has a more complicated expression

    (∂/∂ r+1/r)( (.∇^s)v_r + (.∇^s)U - 2 V v_θ/r) +  ∂/∂ z( (.∇^s)v_z + (.∇^s)W).

We give this equation for the sake of completeness. The Neumann boundary values are inhomogeneous because the derivatives of  do not vanish on the boundary unless the condition, |__a=0, is imposed. In practice, we need only the gradients. To balance the diffusive v_r and v_z, we equate the z-gradient to 

    - ραβ (  U ∂_r v_z + W ∂_z v_z + v_r ∂_r W + v_z ∂_z W   ) - ρβ b_z.

From (<ref>), the kernel of operator ∂_t-νΔ^*_1 is simply M_1. Then the integral for v_θ is

    v_θ(,t)   = ( ∫__a M_1 V_s  '+ ∫_0^t ∫__a M_1 (  (.∇')V + V v_r/r  ) '  t' )
       + ∫_0^t ∫__a (M_1 U/r)   v_θ '  t' + ∫_0^t ∫__a M_1 (U ∂_r  + W ∂_z)  v_θ '  t' 
        = F_a(,t) + ∫_0^t ∫__a M(,',t-s)   v_θ(',s)  '  s,

where the revised kernel, 
M = M_1 U / r  - ∂ (M_1U)/∂ r - ∂ (M_1W)/∂ z, in view of the bounded U and W on _a, and the no-slip v_θ. It follows that v_θ is known, as equation (<ref>) is a linear Volterra integral equation of the second kind with a bounded regular kernel. Since  has been found, the radial gradient satisfies

    ∂/∂ r = - ραβ (  U ∂_r v_r + W ∂_z v_r + v_r ∂_r U + v_z ∂_z U -2V v_θ/r  )- ρβ b_r,

for the compatibility of the first equation in (<ref>). Instead of 
the linearisation and the axi-symmetry, the analogous full non-linear system exists, and its solution can be sought in the same way as we did for the case of infinitely long cylinders. The multiplicity of percolations (<ref>) implies that the perturbative proposition for the annulus flow between co-axial rotating cylinders, as stipulated by the linearised boundary-value problem (<ref>)-(<ref>), loses its potential for prediction. 




§ REMARK ON FLOW SIMILARITY



If we write equation (<ref>) in dimensionless form

    ∂_t ^0 + (^0 ·∇) ^0  = - ∇ p^0/ρ + Δ^0/Re,   ∇·^0=0,

for Re=UL/ν, we can still exercise the flow partition as follows:

    ∂_t ^0 -Re^-1Δ^0= ^0;   ∇ p^0=-ρ (^0 ·∇) ^0- ρ^0.

Since the partition is driven by arbitrary dimension-free buffer ^0, the Re-specified primitive formulation does not admit unique solutions. As there can be no consensus on well-defined velocity fields, the idea of flow similarity is self-contradictory.

The analysis of the Hagen-Poiseuille flow shows that there is not a single characteristic velocity anywhere inside the pipe. The only sensible choice for a velocity scale lies in its initial value. Yet, different motions can be generated in a pipe apparatus (diameter d). In practice, there are no particular reasons that the following profiles cannot be initiated by dedicated pressure-differentials:


  * A parabola, u_z(r,t=0)=U_0(1-r^2);

  * An oscillation, u_z(r,t=0)=U_0(sin(5π(1/2-r))+1)/2;

  * A precipice, u_z(r,t=0)=U_0tanh(π^2(1-r));

  * A slide, u_z(r,t=0)=U_0(1-exp(-π^2(1-r)^2)),

assuming ∇·_0=0. It is logical to select the maximum velocity U_0 as a characteristic scale. It is plausible that the local flows over downstream distances ∼ O(5d - 10d) must be distinct in these four scenarios. The implication is that similarity descriptions of these motions by one dimensionless number, U_0d/ν, are too simplistic to accommodate the diversity of the initial data. Thus, the nature of flow similarities departs from its original formalism. Possibly, a logical way forward is to investigate classifications of starting flows in order for non-dimensional groups to be consequential.




 §.§ Vorticity setting



In principle, fluid dynamics is governed by the conservation law of angular momentum. Then the vorticity, =∇×=(ξ,η,ζ), is of paramount importance (Helmholtz 1858). Dynamically, the vorticity evolves according to  

    ∂_t  - νΔ = (·∇)   - (·∇ ),

which is to be solved subject to data (<ref>) and (<ref>).
The compatibility, = (Δ)^-1 (∇×), implies that  is driven by . The fact is that only the vorticity is a priori bounded, and it has an invariant property under Galilean transform,
=∇×=∇×(- c). Note that a partition of the equation ceases to work, as the pressure is truly dependent, see (<ref>). Briefly, the continuum description of flow evolution, including turbulence, is tantamount to the interactions of both non-linear terms. They are responsible for instigating multitudinous vortices. 

If we take out the common dimension factor [^-2], equation (<ref>) becomes

    ∂_t ^0 + (^0 ·∇ )^0 -(^0 ·∇) ^0 =  ν^0 Δ^0,

where the viscosity ν^0 is relative to unity. This normalised equation shows that the denominator [^2] is a measure of enstrophy. Hence, a non-dimensional parameter may be found by reference the kinematic viscosity, μ/ρ, to a representative circulation ω L^2.
Plausibly, dynamic similarities may be understood for individual categories of initial data in similar geometries. 

The rate of dissipation of kinetic energy per unit mass (denoted by ε_E) has the dimension [^2/]. As implied in (<ref>), the instantaneous dissipation occurs across the whole range of the vorticity eddies, which measure the local strains and shears. On the contrary, the aggregate effects of all the vortices assemble the local velocity. In the absence of a representative velocity, this rate of dissipation may be understood as the energy induced by the local eddies, ε_E ∼ O(·), where =(u^2,v^2,w^2), and =(√(ξ^2),√(η^2),√(ζ^2)) are the energy and eddy vectors respectively. This estimate is just one of the possible results by dimensional analysis. An analogy to the familiar O(||^3/L_Turb) is an approximate, O(L^2_Turb ||^3), which is interpreted as the circulation of local enstrophy. Given vorticity , the governing equations (<ref>) assert that the energy dissipates according to ε_E= - ν ||^2,
which defines a glocal dynamic process, depending on the viscosity at all times. 




§ HELMHOLTZ'S DECOMPOSITION



The reason why the formulation of the primitive variables is ill-advised is due to the fact that the role of the vorticity has not been accounted for. To analyse the vorticity dynamics, we must first consider the well-known vector decomposition. In his paper on vorticity theory, Helmholtz (1858) shows that any smooth velocity field in the whole Euclidean space  can be decomposed into the gradient of a scalar potential and the curl of a vector potential:

    = ∇ϕ + ∇×,

consult 147-148 of Lamb (1932); 20 of Sommerfeld (1950); 1.5 and 2.3 of Morse & Feshbach (1953). Since  is incompressible, we have Δϕ=0. Taking curl of the decomposition, vector calculus gives Δ = ∇(∇·)-∇×. The usual practice is to set ∇·=0. 
In smooth flows decaying at infinity, Laplace's and Poisson's equations are well-posed, and hence the decomposition is necessarily unique. 

In application, there are two scenarios that call for special attention: (a) Generalised flow fields; (b) Incomplete decomposition. First, the velocity has a potential,

    = ∇ϕ,       (∇×=0)

where the vanishing of the curl is assumed to hold. This is a matter of steady solutions for the Euler's equations. Second, in the vorticity dynamics, it is paramount to recover the velocity field from the vorticity,

    = ∇×,

where, once again, it is assumed that a gradient (∇ϕ) does not contribute. We calculate  by solving either 

    Δ() = - ∇×(),   ∈.

or Δ=-, where =∇×. In the latter approach, the gradient is missing in view of the full decomposition (<ref>). The key is that (<ref>) and (<ref>) are incomplete decompositions as some properties of the vector fields are neglected. There is a compatibility issue in these formulas. 

In flow-fields having distributional derivatives, the continuity is now interpreted as

    ∫_(∇·) ṽ  = - ∫_·∇ṽ  = 0

for every test function ṽ∈ C_0^∞(). The right integral makes sense for ∈ L^1_loc. Likewise, we interpret irrotational vectors as generalised derivatives as follows:

    ∫_(∇×)·  = - ∫_·(∇×)   = 0

for every vector test function . There must exist irrotational test functions. An example is constructed in

    _i()={    ( E_0,   E_0,    E_0 ),        s < 1, 
         0,        s ≥ 1, 
    .

where s=√(x^2+y^2+z^2+2xy+2yz+2zx), and

    E_0()=Cexp( - 1/1 - s^2),

for normalisation constant C.

Let us start from the full Helmholtz decomposition for vorticity,

    =∇× + ∇φ,      (∇·=0).

Multiplying (<ref>) by ∇ṽ for test function ṽ∈ C_0^∞ and integrating by parts, we arrive at

    ∫_ṽ (∇· - Δφ)   = 0.

Since this integral vanishes for every test function ṽ, we must have,

    Δφ = ∇·,

For unique (vector) stream function  in (<ref>), we put ='+∇ϕ'+ c, where Δϕ'=-∇·('+ c), and c is any irrotational constant. Then  can be recovered by solving

    Δ = - ,   ∇·=0,

assuming  decays. Conversely, the vector potential is given by 

    () = 1/4 π∫_1/|-| () =1/4 π∫_1/ ||  (+) .

The Biot-Savart law shows that gauge choice, ∇·=0, holds everywhere, provided vorticity is solenoidal, ∇·=0. By Gauss's divergence theorem (cf. (<ref>)),  decays at infinity. Hence potential φ must be harmonic:

    Δφ =  0,    φ→ 0,     || →∞,

where the decay is necessary to avoid harmonic functions which are unbounded at infinity, such as, φ = C_1 x^2 +C_2 y^2 - (C_1+C_2)z^2. The decomposition components are orthogonal as

    ∫_∇×·∇φ     = [ φ∇×]_[]- ∫_φ (∇·∇×)  
        = [ ·∇φ]_[] - ∫_· (∇×∇φ)   = 0,

where the last two integrands vanish identically by vector identities.

Taking dot product of (<ref>) with vector test function ∇×, and integrating, we get the integral condition

    ∫_·(Δ + ∇×)   = 0,

which shows that velocity  is recovered by solving elliptic equation (<ref>).

Conversely, given a solenoidal vector, ∇·=0, where  decays at large ||, how do we decompose this vector field? 
Let  be a non-vanishing vector on  which differs from . Furthermore, we impose the condition → 0 at large distances. Vector analysis shows that

    Δ = ∇(∇·) - ∇×∇×.

By curl of test function , this differential equation is formulated as 

    ∫_∇×·Δ    = - ∫_·Δ (∇×) 
        = -∫_·∇×∇(∇·)  + ∫_·(∇×∇×∇×) 
        = - ∫_·Δ ,

where we introduce function , and ∇·=0, and ∇×∇×=-Δ. Since it is solenoidal, vector  is better characterised as |_|| →∞→ 0. This integral condition reduces to

    ∫_·Δ (∇× -)  = 0,

that shows that  differs from ∇× within multiples of a harmonic function. We denote it by , where ∇·=0, ∇×=0, giving =∇×+ a.e. As in potential flows,  is understood as a generalised function according to (<ref>) and (<ref>). Without loss of generality, we choose a scalar potential φ such that ∇φ =, and hence φ satisfies (<ref>). 

In identity (<ref>), we replace ∇× by . The solution can be expressed by

    ()= 1/4 π∫_1/|-|∇(∇·)  + 1/4 π∫_∇_(1/|-|)×()  .

For every test function ṽ∈ C_0^∞(), we evaluate

    ∫_·∇ṽ  = - ∫_ṽ ∇·  =  (I_1 + I_2)/(4 π),

where

    I_1     = - ∫_ṽ ∫_∇_·(1/|-|)∇(∇·)  
        = ∫_ṽ ( ∫_( -/|-|^3)·∇(∇·) )  ,

and

    I_2=-∫_ṽ ( ∫_^3 ∇_×∇_(1/|-|)·()  )   =0,

where the last equality is established due to the - symmetry in the Newtonian potential.
A quick recap tells us that the divergence is determined by a homogeneous integro-differential equation with an integrable kernel: 

    g() = -1/4 π∫_ K(-)·∇ g() ,

for shorthand g()=∇·(). Denote the Fourier transform of g by

    ℱ(g)()=ĝ()=1/(√(2 π))^3∫_exp( -  ·)  g()  .   ( = √(-1))

We transform equation (<ref>), as the kernel is a convolution type. Re-arranging the result leads to algebraic formula

    ĝ() ( 1 +  ·ℱ( K)√(π/2) )=0,

which is valid for arbitrary wave-number vector . The only consistency is ĝ≡ 0. It follows that its inverse

    ∇·()= 1/(√(2 π))^3∫_exp( ·)  ĝ()   = 0,

identically. 

In summary, the vector field  characterised by differential constraints, (<ref>) and (<ref>), must be incompressible and tends to zero at infinity. The essential reason for its existence lies in the solenoidal property of . It follows that vector  can be identified with vorticity , and  with velocity . In fact, we have recovered the Helmholtz decomposition: =∇×+∇φ, where the kinematics is understood in the generalised sense. In particular, the presence of harmonic φ is a direct consequence of ∇·=0. Vector algebra suggests that the velocity may be non-unique, as ∇×=∇×(+∇q+ c), where q is any scalar function. This arbitrariness is insignificant, as what is being decomposed is the vorticity vector, i.e., the velocity is computed from (<ref>) for given vorticity. 

In fact, our discussion of full Helmholtz's decomposition links to the properties of Laplace's equation, as defined in (<ref>). Multiplying it by test function ṽ∈ C_0^∞(^3) and integrating, we get

    ∫_Δφ ṽ  = - ∫_∇φ·∇ṽ  = ∫_φ Δṽ  = 0,

where the weak function φ is smooth on  in view of Weyl's lemma on the interior regularity (Weyl, 1940). By (<ref>) and (<ref>), the theorem reiterates the physics that there are no sources generating vorticity in the `interior' of the whole space, hence ruling-out discontinuous derivatives. Nevertheless, Weyl's proof does not include the `boundary': the origin and infinity. 
 
For φ∈ C^2(), Green's identity gives 

    ∫_φΔφ  = - ∫_(∇φ·∇φ)   = - ∫_ |∇φ|^2   =0.

Thus, the only smooth solution is trivial φ≡ 0. This fact has long been established, see p.490 of Helmholtz (1858) for regular simply-connected domains, and 35 of Poincaré (1893) for . The textbook concept, ≡∇×, must have been made on the presumption that flow-fields of , are the classical solutions of their respective governing equations. 

If the kinematic solutions are considered in the functional spaces other than C^λ(), λ≥ 2, well-documented harmonic functions show that either they are unbounded at infinity (an example is φ = Ax + Byz + Cxyz), or they possess a singularity at the origin, as a source or a sink or a doublet. The former case is rejected on the ground of finite energy. It follow that Laplace's equation holds on -_1 only. Specifically, the gradient is well-behaved at locations away from _1,

    ∇φ = - κ/4 π(-_1)/ |-_1|^3,     ≠_1,

for arbitrary κ. Because of the singularity, Helmholtz's decomposition on the vorticity (or the velocity) breaks down at =_1. If  becomes singular at location =_s, this singularity has no significance, as ∇φ can be translated to _s=_1.

In view of (<ref>) and (<ref>), the gauge choice, ∇·=0, implies ∇·=0. It is evident that vector identity, ∇×∇×=-, is inherently solenoidal by formula ∇·∇× A=0. From the sufficiency point of view, an arbitrary gauge of  cannot be verified by the Biot-Savart integral in the absence of  regularity. Let us fix a constant or any continuous function f(t)=C=∇·≠0. Solution (<ref>) is now self-contradictory. 




§ VORTICITY BOUNDS



By vector identity, (·∇) = ∇××  + ∇(||^2)/2, the momentum equation is rewritten as

    ∂_t   - νΔ - × (∇×) = - ∇χ,

where χ = χ(,p)=p/ρ + ||^2/2, is known as the total pressure.

The principle of conservation of angular momentum holds in fluid motions. Applying this law to control volume analysis, vorticity dynamics can be derived. Analytically, the pressure can be eliminated from system (<ref>). Vector algebra shows that compatibility of solenoidal (,t) requires Δ = - ∇×∇× at every given instant. Consider test functions, (,t) ∈ C_0^∞()×(0,∞). We take dot product on momentum (<ref>) by ∇×, integrate the result over space and time to obtain

    ∫_∫_·∂_t (∇×)   t + ν∫_∫_   (∇×) · (∇×∇×)   t 
        + ∫_∫_×(∇×) · (∇×)  t = 0,

where p is eliminated from integration by parts as 

    ∫_∫_ p  ∇·(∇×)   t = 0.

The continuity also holds in the weak sense

    ∫_∫_(∇·)  ∇×  t=-∫_∫_ ∇·(∇×)  t = 0.

We call the quantity, ∇×, vorticity, denoted by , the components of which are (ξ,η,ζ). Vorticity describes the local rotation of fluid elements.
Equations, (<ref>)-(<ref>), make sense for ∈ L^1 ∩ L^2, ∈ L^1 ∩ L^2. Now suppose that  is sufficiently regular, so that we can perform further integrations. We have

    ∫_∫_·( ∂_t (∇×) - νΔ (∇×) + ∇×(∇×)× )  t = 0,

for every test function . Thus the vorticity dynamics is written in (<ref>).
In inviscid flows, the contribution from the viscous term is ignored by setting ν≡0. The advantage of the vorticity dynamics is evident as the pressure gradient now exerts no mechanical moment or couple on fluid elements. Furthermore, the vorticity field inherits the incompressibility constraint ∇·=0.  




 §.§ Preliminaries



Given a multi-index of order n  (≥ 1), α=(α_1, α_2, ⋯, α_n) with positive α's, denote |α|=α_1+α_2+⋯ +α_n. Denote the α derivative by

    D^αf() = ∂^|α|/∂ x_1^α_1∂ x_2^α_2⋯∂ x_n^α_n f()= ∂_x_1^α_1∂_x_2^α_2⋯∂_x_n^α_n f().

For two functions f() and g(), derivatives of their products can be found by Leibniz's rule,

    D^α(fg)=∑_α≥β α!/β! (α-β)!  (D^βf)  (D^α-βg),

where α and β are two multi-indices.

To allow for possible unbounded vorticity or velocity, we will discuss our solutions and their derivatives in the weak or generalised sense. For integers, 1 ≤ p ≤∞ and k ≥ 0, we denote the Sobolev space by

    W^k,p() = { v ∈ L^p():    D^α v ∈ L^p(),      α∈ℕ,  |α| ≤ k },

where D^αv are the distributional derivatives, see, for instance, Chapter 3 of Adams & Fournier (2003) for well-known properties. Specifically,  

    W_0^k,p() = W^k,p().

The Sobolev norm, ·_W^k,p(), is defined by

    v_W^0,p()    = v_L^p()=( ∫_ |  v  |^p )^1/p,     (1≤ p < ∞)
    v_W^k,p()    = ( ∑_|α|≤ kD^αv ^p_W^0,p())^1/p, 
    v_W^0,∞()    = v_L^∞() =  ∈  |v()|,
    v_W^k,∞()    = max_|α| ≤ k D^αv_L^∞().


All bold-faced letters are vector qualities. On occasion, we drop the space dependence and write _L^p() as _L^p. Frequently, we have to deal with integrals with integrand f ∇ g. We introduce notation for the boundary term,

    ∫_ f ∂_x g  = - ∫_ g ∂_x f  +lim_x →±∞∫_^2 f g  y  z =-∫_ g ∂_x f   + [ fg ]_[x].

The meanings of [fg]_[y] and [fg]_[z] are clear by cyclic permutation. All three can be summed and extended in [·]_[] as well as in [×]_[].




 §.§ The continuity



The incompressibility, ∇·_0=∇·=0, is a finite sum of three reals. It follows that every ∂ u_i/∂ x_i must be finite. There is a natural number N_1 such that every one of these three derivatives can be bounded

    |∂ u_i/∂ x_i| < N_1  ∑_n=1^3  sup_∈|∂ u_n(,t=0)/∂ x_n|.

The conservation of mass actually imposes strong restrictions on the magnitude of the velocity or the vorticity. Conceptually, it is difficult to comprehend how the continuity would be maintained somehow when  becomes unbounded, as fluid material would have been annihilated—an absurd physics. 

At t=0, the postulation of bounded data (<ref>) enables us to fix a useful bound

    A^*_0=N_2( sup_∈D^βξ_0_L^1 + sup_∈D^βη_0_L^1 + sup_∈D^βζ_0_L^1 ),

for β=0,1,2,⋯, natural number N_2.




 §.§ Properties of non-linearity



Denote the total velocity and vorticity by 

    U(,t) = u(,t) + v(,t) + w(,t),   ,   Ω(,t) = ξ(,t) +η(,t) +ζ(,t),

respectively. The inviscid vorticity equation (ν=0 in (<ref>)) is

    ∂_t  = (·∇)   - (·∇ ).

Consider the x-component ξ. Integrating over space, we obtain

    / t∫_ξ  = [ u    ξ]_[x] + [ u η]_[y] + [ u ζ]_[z]- ∫_ (∇·)  
        - [ u ξ]_[x] - [ v ξ]_[y] -[ w ξ]_[z] + ∫_ (∇·)  =0,

where all the boundary terms vanish by virtue of the decay (<ref>). We repeat these steps for the other two components, η and ζ. This result states that the total torque on fluid elements equals to zero. The dynamics equation for the total vorticity reads

    ∂_t Ω = (·∇) U  - (·∇ )Ω.

Use of the total vorticity is advantageous in scenarios where a pair of counter-rotating vortices involves highly concentrated strong shears. Integration yields

    / t∫_Ω  = [ U    ξ]_[x] + [ U η]_[y] + [ U ζ]_[z]- ∫_ (∇·)  U  
        - [ u Ω]_[x] - [ v Ω]_[y] -[ w Ω]_[z] + ∫_ (∇·) Ω =0.

Every component, as well as the total, satisfy an invariance principle

    ∫_(,t)  =∫__0() ,     ∫_Ω(,t)  = ∫_Ω_0() .

But 

    ∫_ξ = ∫_max(ξ,0)   + ∫_min(ξ,0)  ,

where the last sum is finite and hence each of the integrals must be bounded. Thus, ξ is integrable, as are η and ζ. We come to the summability of the field

    ξ,  η,  ζ∈ L^1(),   ,   Ω∈ L^1().

By the Archimedean property for two positive reals,

    _L^1() < A^*_0.

As  is a member of the equivalent class W^0,1, it can be well approximated by equivalent class of C_0^∞ functions on , in view of (<ref>). In the sequence of approximations, each dies out like ∼ ||^-ε (ε>0) or ||^-2 as a polynomial in ∇·=0 near infinity. In the absence of precise decay rates, which depend on the initial data, the best we can infer is

    → 0,       || →∞.

We have

    ∈ L^∞(_T)× W^(0,1)_0().





 §.§ Weak derivatives



Denote test function ṽ(,t) ∈ C_0^∞()×(0,∞)
Multiplying the x-component of vorticity equation (<ref>) by D^αṽ, α=1,2,3,⋯ and integrating, we obtain 

    ∫_∫_(  ξ  (∂_t + νΔ ) D^αṽ - ( u -  ξ) ·∇ D^αṽ  )  t =0,

in view of ∇·=0 and ∇·=0. Together with the other components, we see that equations (<ref>) make sense for , ∈ L^1 ∩ L^2. For the inviscid case α=1, we carry out integration by parts, and acquire a conservation law for the variation of the torque

    ∫_∫_ṽ ( ∂_t (D ) - D (  (·∇)   - (·∇ ) ) ) =0,

for every test function. Hence, the differential equation for vorticity derivative reads

    ∂_t (D ) = D (  (·∇)   - (·∇ ) ).

Consider D=∂_x first. Integration over space shows that

    / t∫_∂_x    = [ (·∇)   - (·∇ )]_[x]=0,

where use is made of estimate (<ref>). In particular, no use has been made of the assumed regularity on the right-hand function in (<ref>). Analogous operations on the total vorticity admit the invariance

    / t∫_∂_x Ω   = [ (·∇) U  - (·∇ )Ω]_[x]=0.

Repeating these integrations with D=∂_y and D=∂_z in turn, we assert that

    ξ,  η,  ζ∈ W_0^1,1(),      Ω∈ W_0^1,1(),

respectively. By the rationale leading to (<ref>), the rate of change in vorticity must tend to zero at infinity,

    D→ 0,       || →∞.

The upgrade of (<ref>) is obvious.

Simple algebra gives D=∇×D and ∇·D=D(∇·)=0, implying the velocity-vorticity compatibility, Δ D=-∇×D. The velocity inversion is done by

    D() = -1/4 π∫_(-/|-|^3) × D() ,

where the boundary term vanishes in view of results (<ref>) and (<ref>). We denote

    D⊂B_1(_0,R_1),      R_1=R_1(_0) ≫ 1.

For || ≫ || ∼ O(R_1), the denominator of the kernel can be approximated as

    |-| = √(·) √( 1 - 2  (·)/(·) + (·)/(·))≈ ||.

As ||→∞,

    D() ∼/ ||^3 ×∫_B_1 D() ∼/ ||^3 ×∫_B_0 D_0() → 0,

because the last integral is bounded. The support of initial data _0 is denoted by B_0. This decay implies ≡ 0 at infinity. In addition, |_|| →∞≡ 0 by (<ref>). This kinematic state is a consequence of dynamics (<ref>), the boundary condition (<ref>), and the solenoidal properties. Specifically, no claim has been made on the regularity of D in (<ref>). Applying this procedure of estimation to Δ=- ∇× as a consistency check, we recover the  boundary condition for integrable vorticity. The physics dictates the existence of the stagnation, because the vicinity at infinity must remain undisturbed, and cannot be accessed by any localised incompressible flow. One implication is that all higher derivatives of the velocity have to vanish at large far distances. However, we shall proceed with our inductive reasoning. 

The weak differentiation on (<ref>) gives us an equation for the second-order derivatives or α=2 in (<ref>),

    ∂_t D^2  = D ( D (  (·∇)   - (·∇ ) )).

For ∂_t ∂_x D, explicitly, the boundary term from the right-hand side,

    [ D ((·∇) )  - D ((·∇ )) ]_[x] =0,

by deductions (<ref>) and (<ref>). Note that we have no knowledge of the `interior regularity'. In summary, D^2∈ W_0^0,1, and D^2 decays. 

For α≥ 2, since ∇·D^α=D^α(∇·)=0, it follows that

    Δ D^α = - ∇× D^α.

Inverting the Laplacian yields

    D^α()=-1/4 π∫_(-/|-|^3) × D^α() ,

provided that D^α is integrable and tends to zero,

    D^α→ 0,       || →∞.

Following the analysis leading to (<ref>), we establish

    D^α→ 0,       || →∞,

which are valid regardless of the differentiability of the velocity at locations ≪∞.
We take derivatives of order α on (<ref>) to get

    ∂_t (D^α+1) = D { D^α(  (·∇)   - (·∇ ) )},

as we wish to show that D^α+1 is integrable. In fact, we only need to consider the term in the outer brackets at infinity. Each term in D^α(·) can be expanded according to Leibniz's rule (<ref>). In each of the two expansions, every resulting term consists of products of derivatives of the vorticity and velocity. All these derivatives are equal to zero as x,y,z →±∞, in view of (<ref>), (<ref>), as well as the decays of the lower-orders. Induction shows

    / t∫_ D^α (ξ, η, ζ)  =0,     α=1,2,3,⋯.

An alternative interpretation of this outcome suggests that the integrations and differentiations commute in view of the attenuation at large distances, even though the smoothness of the flow in finite spheres is undetermined. Consequently, the bound analogous to (<ref>) is given by

    D^α_L^1() < A^*_0,

and hence

    ∈ W_0^k,1(),     k=1,2,3,⋯.

By virtue of the Sobolev embedding theorem (see, for instance, Theorem 4.12 of Adams & Fournier, 2003), we arrive at a space-wise bound 

    ∈ C^2_B(),

which, a posteriori, justifies the assumption made in (<ref>), as the partial Helmholtz decomposition = ∇×
enshrines the definition of vorticity at every given instant. The point is that it is best to recover the velocity from the vorticity by solving (<ref>), 
Δ () = - ∇× (). Remarkably, the integrability of vorticity (L^1-norms) and the `boundary' decay (<ref>) are the essential ingredients. The space-time bound (<ref>) is improved to

    ∈ L^∞(_T)× C^2_B().

Additional information of the velocity is known according to the well-developed theory of elliptic partial differential equations (see, for instance, Gilbarg & Trudinger, 1998). For the time being, it is sufficient to notice

    ∈ C^2_B().


For Hölder continuous ∈ C^0,λ(), 0 ≤λ <1, it is well-established that ∇∈ C^1(), which is regular enough to justify the omission of the harmonic component in =∇×+∇φ. The implication is that our integrations of the vorticity field in the Sobolev space W^k,1 are very effective but are not necessarily optimum. As C^0 is embedded in W^3,1, the next regularity upgrade results in ∈ W^4,1 which misses all fractional values of λ.
 
Taking divergence on Euler's momentum,

    ∂_t  = - (·∇) - ∇ p/ρ,

gives us Poisson's equation for the pressure

    Δ p() = - ρ ∇·((·∇)) = - ρ ∇·( (·∇)) = -ρ N_0(),

where every term on the right-hand side has been shown to be continuous and bounded, see (<ref>) and (<ref>). Consequently, the pressure gradient is bounded in space. In particular, the gradient tends to zero at large distances, giving rise to the gauge pressure at infinity. 




 §.§ Vorticity dynamics for real fluids



Let the vorticity equation,

    ∂_t  - νΔ = (·∇)   - (·∇ )≡ V(,t).

This dynamics is to be solved subject to data (<ref>) and decay (<ref>). The total vorticity satisfies the scalar equation

    ∂_t Ω - νΔΩ = (·∇) U  - (·∇ )Ω.

By Duhamel's principle, equation (<ref>) can be transformed into an equivalent integral equation

    (,t)=∫_  H(-,t) _0()  + ∫_0^t ∫_  H(-,t-s)   V(,s)   s,

where H is the diffusion kernel 

    H(,t)= 1/ (4 π t)^3/2 exp( -  ||^2/4 t).

Integration of (<ref>) over space defines a constraint

    ∫_(,t)   = ∫__0()  + ∫_0^t ( ∫_ V(,s)  )  s.

The last integral equals to zero, as analysed in the preceding paragraphs. Thus the invariance,

    / t∫_  = 0,     (ν > 0)

holds, whereby viscous  satisfies bound (<ref>). Evidently, the above derivation applies equally well to the total vorticity governed by (<ref>). 

In view of (<ref>), a sequence of differential equations also works for weak -derivatives,

    (∂_t  - νΔ )D^α = D^α V(,t).

For every α, the analogous equation to (<ref>) reads

    D^α(,t)=∫_ H(-,t) D^α_0()  + ∫_0^t ∫_ H(-,t-s)  D^α V(,s)   s,

and hence

    ∫_ D^α(,t)   = ∫_ D^α_0()  + ∫_0^t ( ∫_ D^α V(,s)  )  s.

Without repeating the analysis on the integral of D^α V, we assert that bounds, such as (<ref>) and (<ref>), are also valid for the viscous equation. 

In summary, we derive a priori bound 

    ||,  |∇| < M_0(_0)< ∞,       (∈).

These space-wise estimates are effective at every given instant. 




 §.§ Enstrophy and energy



The rate of change of enstrophy can be computed as follows. Taking dot product of equation (<ref>) by  and integrating the result, we see that

    1/2/ t∫_||^2   + ν∫_|∇|^2   = ∫_· (·∇)  - ∫_· (·∇) ,

where we have made use of the vorticity decay on the second term on the left. In view of the solenoidal velocity and vorticity, the last term on the right vanishes as a result of integration by parts. The remaining term is estimated as

    ∫_|· (·∇)|  ≤_L^∞∇_L^1≤_L^∞^2_L^2,

by applying Hölder's inequality. Thus, the evolution is governed by

    / t^2_L^2 < 2 M_0 ^2_L^2,

the solution of which bounds the enstrophy,

    (·,t)^2_L^2 < _0^2_L^2exp( 2 M_0 t).

The exponential is likely an over-estimate. Importantly, the law of energy conservation holds at all times t>0,

    1/2∫_ρ|(,t)|^2   + μ∫_0^t ∫_ |(,t)|^2  = 1/2∫_ρ|_0()|^2  .


Finally, we have arrived at the bounds over time (μ≥ 0):

    , ∈ L^∞(_T) ×  L^2() ∩ C^2_0(),   ∇ p () ∈ C^1().

The most important result is written as

    (,t)_L^∞(_T)× L^∞() < ∞.

Admittedly, the present analyses on vorticity dynamics are much more precise and complete when compared with an earlier paper of Lam (2013). Note that the revised bounds are valid from t>0, including the Leray-Hopf regularity interval.




 §.§ Uniqueness



Let (,q) be a second solution of the Navier-Stokes equations in addition to (,p). Denote =-, and ∇·=∇·-∇·=0. The equation governing  reads

    ∂_t  + (·∇) + (·∇) = ∇ (p-q)/ρ + νΔ.

We multiply this equation by  and integrate the result over space. Then we estimate the magnitude of third term

    ∫_|· (·∇)|  ≤∇_L^∞^2_L^2,

by (<ref>). Since the enstrophy ∇×_L^2 is never negative, the energy is determined by

    / t^2_L^2 ≤ 2 ∇(·,t) _L^∞^2_L^2.

The sum of the energy shows that the difference satisfies

    (·,t)-(·,t) ^2_L^2≤  C _0(·)-_0(·) ^2_L^2.

If the initial data of the flows coincide, solution (,t) agrees with (,t) everywhere. In conclusion, given the complete knowledge about the vorticity, the solutions of the Navier-Stokes equations are unique for any localised initial data of finite energy. A corollary proposition is that the Cauchy solutions of Euler's equations are well-posed for t>0.

It follows that the Cauchy problem of the incompressible vorticity equation possesses classical solutions for given smooth initial data having compact support. The solutions of the integral equation,

    (,t)=∫_    H(-,t) _0() 
       + ∫_0^t ∫_ H(-,t-s)  ( (·∇)   - (·∇ )) (,s)   s,

where

    () = -1/4 π∫_-/|-|^3×()  ,

and

    ∇() = 1/4 π∫_∇(1/|-|)∇×()  ,

have been found to be a superposition of multiple convolutions of the initial data. 
The continuum dynamics instigates vortices or eddies of spatio-temporal scales of numerous ranges, see equation (8.11) of Lam (2013). The non-linearity in (<ref>) is solely responsible for the production of the multitudinous structures, giving rise to turbulence. Critically, the field of the velocity-vorticity is self-consistent, and is kinematically regulated by the Helmholtz decomposition. The integral (<ref>) for the velocity gradient does not contain a singular kernel. In practice, it is almost unavoidable to involve heavy numerical computations in fluid simulation. 




§ CONCLUSION



The Cauchy problem of the primitive equations does not have any unique solution beyond the local regularity time. We have shown that the primitive setting for incompressible flow is over-specified in the sense that the dynamics can be split into a solenoidal diffusion and a reactive non-linear convection. The exception is the flows of zero pressure gradient. Our analysis clarifies the non-existence of strong a priori bounds, that are crucial in verifying the theoretical consistency. In essence, the Navier-Stokes equations are incompletely formulated. As asserted, the concept of unstable flow stems from the approximations due to the ad hoc linearisation. Above all, there exist no immediate instability counterparts in the exact vorticity dynamics. One consolation is that we understand why we cannot satisfactorily explain most of the well-established experimental observations about fluid flows, should we squarely cling on to the primitive perspectives.[The author recalled that, in an introductory undergraduate course on fluid mechanics, one of his professors said the question whether the velocity is induced by the pressure or vice versa is a chicken-and-egg argument.]

Fortunately, the global regularity of the vorticity equation for incompressible flows can be proved beyond doubt, thanks to the invariance of the vorticity field. 



tocsectionReferences




  * Adams, R.A. & Fournier, J.J.F. 2003 Sobolev Spaces, 2nd Ed. Amsterdam:Elsevier Academic Press.


  * Burgers, J.M. 1948 A mathematical model illustrating the theory of turbulence. Adv. Appl. Mech. 1, 171-199.


  * Carslaw, H.S. & Jaeger, J.C. 1986 The Conduction of Heat in Solids. Oxford:Clarendon. 


  * Gilbarg, D. & Trudinger, N.S. 1998 Elliptic Partial Differential Equations of Second Order. Berlin:Springer.


  * von Helmholtz, H. 1858 Über Integrale der hydrodynamischen Gleichungen, welche den Wirbelbewegungen entsprechen. J. Reine Angew. Math. 55, 25-55. English translation by P.G. Tait, (1867) On integrals of the hydrodynamical equations, which express vortex motion. Phil. Mag. 33, 485-512.


  * Hopf, E. 1951 Über die Anfangswertaufgabe für die hydrodynamischen Grundgleichungen. Math. Nachr. 4, 213-231.


  * Hopf, E. 1952 Statistical hydromechanics and functional calculus. J. Rational Mech. Anal. 1, 87-123.


  * Kelvin, Lord (Thomson, W.) 1887 Stability of fluid motion: rectilinear motion of viscous fluid between two parallel plates. Phil. Mag. 24(5), 188-196.


  * Lam, F. 2013 .


  * Lam, F. 2018a .


  * Lam, F. 2018b .


  * Lamb, H. 1932 Hydrodynamics. Cambridge Univ. Press.


  * Leray, J. 1934a Sur le Mouvement d'un Liquide Visquex Emplissant l'Espace. Acta Math. 63, 193-248.


  * Leray, J. 1934b Essai sur les mouvements plans d'un liquide visqueux que limitent des parois. J. de Math. XIII (IV),  331-418.


  * Morse, P.M. & Feshbach, H. 1953 Methods of Theoretical Physics I. New York:McGraw-Hill.


  * Navier, C.L.M.H. 1823 Memoire sur les Lois du Mouvement des Fluides. Mem. Acad. Sci. Inst. France 6, 389-440.


  * Onsager, L. 1949 Statistical hydrodynamics. Nuovo Cimento (Supplemento) 6, 279-287.


  * Poincaré, H. 1893 Théorie des Tourbillons. Paris:Georges Carré.


  * Prandtl, L. 1904 Über Flüssigkeitsbewegung bei sehr kleiner Reibung. 
In Proc. Third Intern. Math. Congress, Heidelberg. Reprinted in Vier Abhandlungen zur Hydrodynamik und Aerodynamik, pp. 1-8. Göttingen (1927); English translation: Motion of fluids with very little viscosity, NACA TM-452 (1928).


  * Rayleigh, Lord (Strutt, J.W.) 1916 On the dynamics of revolving fluids. Proc. Roy. Soc. A 93, 148-154. Also Scientific papers (1920), vol. VI, pp.447-453. Cambridge Univ. Press.


  * Reynolds, O. 1883 An experimental investigation of the circumstances which determine whether the motion of water shall be direct or sinuous, and of the law of resistance in parallel channels. Phil. Trans. Roy. Soc. 174, 935-982. Also Scientific Papers (1901), vol. II, pp.51-105. Cambridge Univ. Press.


  * Roach, G.F. 1982 Green's Functions. Cambridge Univ. Press.


  * Schubauer, G.B. & Skramstad, H.K. 1947 Laminar boundary-layer oscillations and transition on a flat plate. J. Res. Nat. Bur. Stand. 38, 251-292.


  * Serrin, J. 1962 On the interior regularity of weak solutions of the Navier-Stokes equations. Arch. Rational Mech. Anal. 9, 187-195.


  * Sommerfeld, A. 1950 Mechanics of Deformable Bodies. New York:Academic Press.


  * Stokes, G.G. 1845 On the theories of the internal friction of fluids in motion. Trans. Cambridge Phil. Soc. 8, p287. Also Mathematical and Physical Papers (1880), vol. I, pp.75-129. Cambridge Univ. Press.


  * Taylor, G.I. 1923 Stability of a viscous liquid containing between two rotating cylinders. Phil. Trans. Roy. Soc. Lond. A 132, 289-343.


  * Tricomi, F.G. 1957 Integral Equations. New York:Interscience Publishers.


  * Weil, H. 1940 The method of orthogonal projection in potential theory. Duke Math. J. 7, 411-444.








28 March 2021





