Quantifying wall turbulence via a symmetry approach. Part I. A Lie group theory

First principle based prediction of mean flow quantities of wall-bounded turbulent flows (channel, pipe, and turbulent boundary layer - TBL) is of great importance from both physics and engineering standpoints. Here (Part I), we present a symmetry-based approach which derives analytic expressions governing the mean velocity profile (MVP) from an innovative Lie-group analysis. The new approach begins by identifying a set of order functions (e.g. stress length, shear-induced eddy length), in analogy with the order parameter in Landau's mean-field theory, which aims at capturing symmetry aspects of the fluctuations (e.g. Reynolds stress, dissipation). The order functions are assumed to satisfy a dilation group invariance - representing the effects of the wall on fluctuations - which allows us to postulate three new kinds of invariant solutions of the mean momentum equation (MME), focusing on group invariants of the order functions (rather than those of the mean velocity as done in previous studies). The first - a power law solution - gives functional forms for the viscous sublayer, the buffer layer, the log-layer, and a newly identified central `core' (for channel and pipe, but non-existent for TBL). The second - a defect power law of form $1-r^{m}$ ($r$ being the distance from the center line) - describes the `bulk zone' (the region of balance between production and dissipation). The third - a relation between the group invariants of the stress length function and its first derivative - describes scaling transition between adjacent layers. A combination of these three expressions yields a multi-layer formula covering the entire flow domain, identifying three important parameters: scaling exponent, layer thickness, and transition sharpness. All three kinds of invariant solutions are validated, individually and in combination, by data from direct numerical simulations (DNS).


I. INTRODUCTION
Wall-bounded turbulent flow has received considerable recent attention because of its relatively simple geometry and its relevance to a large variety of engineering applications [1]. It enables direct insight into fundamental turbulence physics, as well as direct verification by local (such as profiles of mean velocity, Reynolds stress, fluctuation intensities, etc.) and integral measurements (such as skin friction and heat transfer). The fundamental challenge is to predict the mean flow properties, including the mean velocity profile (MVP), the Reynolds stress, the kinetic energy, etc; however, a deductive theory of this kind is still missing.
A recent model by L'vov, Proccacia and Rodenko (cited hereinafter as LPR) [2] is noteworthy, as they deal with the MME and the balance equations for the kinetic energy and the Reynolds stress. Several physical arguments regarding energy cascade, dissipation and balance in the average Reynolds stress equation, enable LPR to prescribe the mixing length profile which involves two artificial functions inspired from the limited DNS data -one connecting the log-layer to the wall (e.g. a wall function) and the other to the center (e.g. a bulk function). This two-layer model achieves a good accuracy in describing the data, but the asymptotic behaviors for the wall and bulk functions are not correct [3], and their extensions to describe other flows and other profiles (energy dissipation, kinetic energy, etc) are not obvious.
Recently, a SED theory [4] has been proposed to reconsider this problem of determining the mean velocity profile of canonical wall-bounded flows from first principles. It is known that the mean momentum equation is unclosed, and the determination of the mean strain requires a modeling of the Reynolds stress. Prandtl's [5] mixing length model provides a way to model the Reynolds stress. In the past, simple models for the mixing length were adopted; the most famous is the Karman's linear law of the distance to the wall, leading to a logarithmic profile for the mean velocity. The log-law has had a tremendous success in predicting the Reynolds number dependence of the friction coefficient and in calibrating the computational fluid dynamics (CFD) models of turbulent flows such as K-ǫ or K − ω [1]. Nevertheless, the validity of the logarithmic profile, even for simple wall-bounded flows, has been questioned [7][8][9], owing to the ad hoc nature of the arguments behind its derivation. In [4], we have extended the classical mixing length theory of Prandtl and von Karman by proposing a multi-layer formula for the entire profile of the mixing length ℓ + M -hence the entire mean velocity profile -written as where y + is the coordinate normal to the wall, r is the radial distance of the pipe, y + sub , y + buf and r core are parameters specifying the thicknesses of the layers (Z c = (1 + r 2 core ) 1/4 is a normalization factor), p 1 and p 2 are parameters describing the sharpness of the transitions between the layers, m is a parameter specifying a geometry-dependent scaling for the bulk flow (empirical study shows that m = 4 for channel and m = 5 for pipe [4]), and the superscript '+' denotes wall unites normalization by the wall friction velocity u τ and viscosity. This formula (1), seemingly complex, is in fact composed of a series of simple functions, namely ℓ + M ≈ ρ/y +3/2 sub y + 3 2 , y + ≪ y + sub , r → 1 (2) ≈ ρ/y +2 sub y +2 , y + sub ≪ y + ≪ y + buf , r → 1 ≈ ρy + buf /y +2 sub y + ≈ κy + , y + buf ≪ y + , r → 1 (4) ≈ κ(1 − r m )Re τ m , y + buf ≪ y + , r core ≪ r ≤ 1 ≈ κr 1/2 core Re τ mZ c r 1/2 , r → 0 (6) Clearly, formula (1) is consistent with the log-law (see (4)), and at the same time describes the finite Re modification (see (5)) taking into account the effect of outer flow geometry (with the parameter m), so it naturally reconciles the log-law and power-law debate. All functional forms from (2) to (6) are analytically given, which is a definite step beyond the LPR model [2]. The purpose of this paper is to present analytic derivation of (1), employing a novel Lie-group formalism.
The main hypothesis in this derivation is that turbulent fluctuations recover dilation symmetry (discussed below) either from the wall or from the center (for channel or pipe), which is believed to be a general feature of wall-bounded turbulent shear flows. Such dilation symmetry manifests in a power law scaling of a relevant length (order) function defined using both fluctuation and mean quantities. The introduction of the order function [4] is inspired by the observation that a real turbulent flow system (such as a turbulent pipe) differs in a fundamental way from idealized system of homogeneous isotropic turbulence or homogeneous shear turbulence by involving a series of statistically distinct states of fluctuations. For example, fluctuations near the wall in the so-called viscous sublayer in a pipe is very different from that in the bulk flow. In each layer/region, the energy dynamics display different dominant mechanisms. For instance, in viscous sublayer, pressure transport balances the viscous dissipation and diffusion, while in buffer layer, production (by the mean shear) balances viscous dissipation. As the number of terms in the energy equation is finite, the number of balance mechanisms is finite, hence the finite number of layers.
This multi-state/multi-layer feature is believed to be generic for all systems at a truly far-from-equilibrium state; wall turbulence is one of them. Any measure that quantifies the symmetry property of each state and also describes transitions between one state to another is thus named order function, in close analogy to 'order parameters' in statistical physics [6], which effectively captures the unknown effects of fluctuation structures on the mean and reveals phase transitions (which means from one macroscopic state to another) at macroscopic scales. Our order function plays exactly the same role -capturing effects of turbulent fluctuations on the mean and displaying varying scaling symmetry in the mean profile. A specific fluctuation field responsible for the vertical momentum transfer is the product of streamwise and vertical fluctuations -the Reynolds stress -which is to be described by the mixing length (order) function. Our central task here is to characterize the statistical state of each layer using the scaling of the mixing length, and then the solution of the mean momentum equation. In the future, a complete deductive theory would explain why each layer is necessary and how its quantitative property is derived.
The paper is organized as follows. In section 2, we establish a formal framework for the Lie-group analysis applied to the MME, where an extended phase space including the mixing length and its gradient is introduced. The analysis is applied in Section 3 to derive a specific solution of the mixing length for turbulent channel and pipe flows. In section 4, we validate the theory using DNS data. Section 5 contains a summary and discussion of future studies.

II. LIE-GROUP ANALYSIS OF THE INNER AND OUTER MME
Lie-group analysis is a well-founded method of finding invariant solutions of differential equations [10]. Recently, Oberlack [11] and Lindgren et al [12] use it to derive local scaling laws in the mean velocity, including log-law and power-law, for unclosed Reynolds-averaged Navier-Stokes (RANS) equations. When studying multi-point statistical symmetry, Oberlack and Andreas [13] assert that there exist a large number of possibilities for defining invariant groups, due to the possible new scaling arising from the fluctuations. On the other hand, empirical data from DNS show well-behaved mean profiles, while finding new governing principles to derive it is a challenge. For a symmetry analysis, this amounts to selecting the right invariant groups, which likely correspond to statistically restored symmetries by fully developed turbulence, postulated by Frisch [14]. Finding these right symmetries embedding true invariant solutions of the RANS equations (i.e. the MME) is an important step forward in the symmetry analysis of turbulence.
A key realization here is that although the RANS equation and its boundary conditions for a channel or pipe possess a dilation group (but not a Galilean group), the empirical mean velocity solution seem to break the dilation invariance. Even locally, the mean velocity does not have a recognized scaling invariance, except near the wall with a rather trivial linear scaling in y + . Specifically, when the RANS equation undergoes a dilation transformation, the group invariant of the mean velocity solved by using the Lie-group method is no longer invariant, as can be easily checked with any empirical mean velocity profile. Classical Lie-group analysis seemed to meet a bottleneck, yielding very few results. This obstacle is overcome here by the introduction of additional variables, namely mixing length and its gradient (both order functions), as group invariants.
Note that in the statistical physics of phase transitions, order parameters are considered to the primary variables describing the macroscopic state; we carry out the same tradition: regarding the order functions (i.e. mixing length and its spatial gradient) as fundamental variables describing the dilation-symmetry of the underlying fluctuation field (i.e. Reynolds stress). In other words, when undergoing dilation transformation, the mixing length and its gradient display specific scaling in each layer, and a generalized symmetry relating the mixing length with its gradient describes the transitions between the layers. This Lie-group formalism thus enables us to show that there exists a dilation-group which leaves the governing equation (i.e. the RANS equation) invariant under the postulated scaling symmetry for the mixing length. While the mean velocity is itself no longer invariant under the dilation transformation, the mixing length maintains the scaling symmetry.
This unconventional use of the Lie-group analysis embodies a physical insight about the problem, namely that the mixing length and its gradient are the right variables to represent the symmetry property of the fluctuations (e.g. the Reynolds stress). As noted by Cantwell [10], group theory begins when all the science needed to deduce the solution has been established. In wall turbulence, the physics that has been overlooked so far is the existence of the dilationgroup symmetry displayed by the mixing length in each of layers like sublayer, buffer, bulk and core. A general feature uncovered here is that distinct layers do possess symmetries, and there exist appropriate order functions to display them. Including order function as additional variables in the symmetry analysis of the RANS equation marks the essential difference between this work and all previous symmetry analysis of the RANS equation.
For a channel/pipe flow, the MME is [15]: where U is the streamwise mean velocity, u and v are streamwise and normal fluctuating velocities, Re τ is wall friction Reynolds number, and other notions are described above. Integrating along the normal direction yields where S + = dU + /dy + , W + = − uv + and r = 1 − y + /Re τ ). Prandtl's [5] mixing length, see Fig.1, defines a length function linking the Reynolds stress and the mean shear stress: Then, the MME becomes: Note that ℓ + M fully determines the mean shear stress, and then the MVP. Prandtl [5] assumed a linear, which leads to  [5]. Inset shows the central power-law scaling with exponent −1/2, where ℓM = ℓ + M /Reτ , and r = 1 − y is the distance to the center. DNS data from Iwamoto et al. [16] at Reτ = 650, Hoyas and Jimenez [17] at Reτ = 940, Wu and Moin [18] at Reτ = 1142.
since W + ≈ 1, and therefore the logarithmic law for MVP with Karman constant κ ≈ 0.41 and B ≈ 5.2 [15]. The logarithmic law is widely used in engineering models, albeit invalid in the near wall region and center flow because ℓ M departs from the linear scaling (13) when the wall distance y is out of the overlap region; see in Fig.1. In fact, in the viscous sublayer, Taylor expansion [15] yields W + ∝ y +3 , S + ≈ 1, thus ℓ + M ∝ y +3/2 . In the core region, according to the central symmetry, the mean velocity profile is parabolic, thus S + ∝ r; which together with W + ≈ r, gives ℓ M ∝ 1/ √ r. Note that these two power-law scalings are sharply different from (13). Although various modifications have been made ever since, such as the damping function by van Driest [19] and the wall function in the LPR model [2], a derivation of mixing length function with correct asymptotic behaviors is not achieved yet.
The Lie-group analysis presented below is novel in two aspects. First, by substituting in W + = (ℓ + M S + ) 2 , we rewrite the governing equation (8) as: with wall conditions: (8) is an ordinary differential equation, but since we consider below a family of solutions for varying Re τ and regard Re τ as an independent variable, the partial differentiation is introduced in (16). This is the so-called inner MME (in wall units), and the primary variables are expanded to include the mixing length ℓ + M and its spatial gradientl + M = ∂ℓ + M /∂y + . Similarly, in central unit r, we have the outer MME (l M = ∂ℓ M /∂r): with center condition U + (0) = U + c , ℓ M (0) = ∞ where U + c is the centerline velocity. Note that in (16) and (17), in addition to independent variable y + (or r) and Re τ , and dependent variable U , we now have an expanded function space including ℓ M andl M . This expansion is essential, since it allows one to formulate specific assumptions about the scaling symmetry of the fluctuation fields (or of some specific measure of it, e.g. the Reynolds stress), which then allows one to specify the solution manifold. In other words, using ℓ M andl M , we quantify the symmetry behavior of the Reynolds stress, and use it to derive a solution for the mean velocity. Our Lie-group analysis not only identifies simple scaling symmetries of ℓ M , locally valid in each of layered structures, but also defines a generalized symmetry in the prolonged space of ℓ M andl M , describing a universal transition between one layer to another, and then fulfills the goal of defining a complete solution for U . The generalized symmetry is a truly hidden symmetry which describes the multi-layer structure of wall-sheared turbulence. Since the multi-layer structures are generic to far-from-equilibrium systems, our treatment is of general interest for many flows.

A. Symmetry transformation
The analysis consists of the following three standard steps.
Step one is to define the symmetry transformations. Here, the Lie group transformation S ε can be expressed as where the superscript * denotes transformed variables, t = (y + , Re τ ) and x = (U + , ℓ + M ) are independent and dependent variables, respectively, and Φ, Ψ are analytic functions of t, x and a continuous parameter ε. In this setting, the symmetry transformation satisfies Note that S ε can also be expressed in an infinitesimal form: where , yields the following equation for the infinitesimals [11]: where X is the so called Lie-group operator, and X p = ηU + ∂U + +ηÜ + ∂Ü + +ηl+ M ∂l+ M is referred to as the 'prolongation'. (21) is an over-determined set of linear homogeneous differential equations for infinitesimals, whose general solutions can be conveniently given by the software Maple13, yielding: where coefficients a i and b i (i = 1, 2) denote dilations and translations, respectively. Now, the no-slip wall condition, i.e. U + * = ℓ + * m = 0 at y + * = 0, leads to b 1 = b 2 = 0; thus, near the wall, only the dilations are permitted. Normalized by the parameter a 1 , we have where α = (−1 + a 2 /a 1 )/2 is the only remaining parameter specifying the local scaling of the mixing length. In addition, the infinitesimals for the derivatives are ηU + = −2αU + , ηÜ + = (−1 − 2α)Ü + and ηl+ M = (α − 1)l + M using the contact condition.
Here, we briefly discuss the meaning of the above symmetry transformation. Consider the transformation from the state characterized by (y + ,Re τ ) to that by (y + * , Re * τ ). The transformed variables following (23) are: We call a symmetry in U + to denote the fact that the actual (empirical) solution of the mean velocity U + at the latter state equals the transformed mean velocity, namely U + * . In other words, if U + * does not equal the actual ) (denoted by a dashed line). Here, the superscript ES denotes empirical solution provided by experimental or DNS data, which is the true solution at (y + * , Re * τ ). The symmetry analysis consists in asking: is there any similarity between A * and B? Previous works ( [11,12]) focus on the symmetry in the mean velocity itself, namely demanding U + * = U +ES ; we argue that this symmetry is broken, at least in most parts of the flow domain, that is, U + * = U +ES . Instead, we postulate the symmetry of the mixing length and its gradient, because they represent the restored symmetry [14] by turbulent fluctuations (e.g. effect of the Reynolds stress on the mean shear).
solution, we say that the symmetry in U + breaks down. Previous works [11,12] applying the Lie-group analysis to the RANS equation were restricted to the case with no symmetry-breaking of U + . But, if one wants to derive solutions of the RANS equation from symmetry analysis, the situation of symmetry-breaking needs to be considered, since the empirical MVP from either DNS or experimental data do not seem to have a simple scaling form.
The new idea introduced here consists in recovering the statistical restored scaling symmetry by turbulent fluctuations in wall-bounded flow with extra variables, namely the mixing length and its gradient, see Fig. 2. By assuming scaling symmetry in ℓ + M , we assert that the transformed mixing length, ℓ + * M , is the actual solution at the latter state corresponding to the independent variables (y + * , Re * τ ). This postulated scaling symmetry leads immediately to a power-law solution of ℓ M , which corresponds to the similarity solution of the first kind described below.
Another similarity solution (of the second kind) corresponds to a case of symmetry-breaking in ℓ M , but a maintained scaling symmetry inl M , which is a solution of defect power-law. A third kind of similarity solution is derived under an assumption of a generalized symmetry relation between ℓ M andl M , when the scaling symmetries of both quantities break down. By including the prolongation of ℓ M , only three possibilities appear; hence, our symmetry consideration yields a fairly complete description of the dilation-invariant solutions to the RANS equation. These details are presented in section 2.2.
Step two is to define group-invariants, which are obtained by solving the characteristic equation: dt i /dξ i = dx i /dη i , written as: The solutions to (25) define following group-invariants: These group-invariants, functions of y+ and Re τ , are the so-called similarity variables [10]. Note that F 1 and F 2 are equivalent to the invariant q = y/ xν/U ∞ defining the similarity variable in a Blasius boundary layer. Finally, in step three, we rewrite C = 0 in terms of the group-invariants, Geometrically, C = 0 defines an invariant-surface in the variable space of t and x, which is invariant under a dilation group of transformation, and the infinitesimals are tangent vectors to the invariant-surface. The so-called invariant solution, for example, Θ(t, x) = 0, satisfying the invariant surface condition, i.e. XΘ| Θ=0 = 0, must be a function of similarity variables. The present work seeks invariant solution in an extended coordinate space including the mixing length and its gradient, namely Θ(F 1 , F 2 ) = 0, (e.g. ξ y + , η ℓ + M and ηl+ M ). Similarly, the transformed variables for the outer MME are r * = e ǫ r, Re * τ = e (1/2+α)ǫ Re τ , and we derive a set of group-invariants for the outer flow, by performing the above analysis to the outer MME (in center units): where · denotes r-gradient; the outer MME in terms of group-invariants is It is interesting to contrast (27) and (30) with similar forms of the Blasius equation for a laminar boundary layer. Following [10], the Blasius equation written in streamfunction ψ is ∂ψ ∂y It can be demonstrated that this equation is invariant under the following dilation symmetry transformations: x * = e ǫ x; y * = e λǫ y; The group parameter λ is determined by the free stream condition, i.e.
which leads to with the following group invariants: Substituting them into (31) yields This equation can be re-written as where X i represents the ith order derivative of the invariant f . The Blasius equation (37) are the equivalent of our equations (27) and (30), but with an important difference: (37) is a closed equation, while (27) and (30) contain two dependent variables (U and ℓ M ). The well-known closure problem of turbulence requires the introduction of closure hypothesis. In a symmetry analysis, closure hypothesis comes from symmetry assumption. In the previous attempts [11,12], symmetry is imposed simply on the mean velocity, that is, the group invariant of U + is assumed to be constant, i.e. G 1 (y + , Re τ ) = G 1 (y + * , Re * τ ) = const.
As shown in ( [11]), with this assumption, candidate solutions of power-law and defect power-law can be obtained for the mean velocity, The power-law solution finds its correspondence with the linear scaling in the viscous sublayer with β = 1; the defect power law solution, U + = U + c − G 1 y β in coordinate y, corresponds to flow solution near the central part of a channel, according to Oberlack [11], with fitting parameters G 1 = 5.83 and β = 1.69. However, it can be shown that this latter solution contradicts the central symmetry condition, U + c − U + ∝ r 2 , hence should be regarded as non-physical. Therefore, we consider that imposing the invariance on U to derive closure solutions of the RANS equation is not successful.
With the introduction of the new dependent variables, ℓ M and its gradient, the situation is different. First, we are given new freedom in constructing the solutions, because the symmetry associated with the fluctuations (included in ℓ M ) can be treated separately. Secondly and most importantly, we can independently verify if the assumed symmetry is true. If ℓ M is not a good variable, we can study a variant of it, until we find the good order function. Thirdly, as shown below, we are able to derive rigorously two local solutions near the wall and around the centerline of a channel/pipe, similarly as deriving the 1/2 exponent in (34) for laminar boundary layer. Finally, this program yields a general procedure, with the concept of the order function, of performing symmetry analysis for unclosed problems involving fluctuations. Below, we show that only three kinds of symmetries are sufficient for determining the full profile of mean shear and hence mean velocity according to (12), which is highly encouraging.

B. Three kinds of similarities for mixing length
When considering the possible symmetries of ℓ + M andl + M , three cases naturally arise.

Power-law
The simplest situation is that the scaling symmetries of ℓ + M andl + M are both maintained, as the group parameter ǫ continuously changes. This can be sketched by the mapping:  where the superscript ES denotes the empirical solution. With this symmetry assumption, (26) yields constants F 1 and F 2 : Solving (26), we thus have a power-law solution for ℓ + M : Note that F 2 /F 1 = α is a key quantity to characterize the solution. Since the scaling symmetry is postulated, the actual solution may realize only a specific value of α. Thus, the scaling symmetry assumption provides a candidate solution to be validated by empirical data. Specifically, the theory proposes to plot the empirical data in terms of a so-called diagnostic function: If the empirical γ displays a plateau in a range of y + (or r), then, a local scaling symmetry of ℓ + M is validated, and the value of the plateau is thus α which guarantees the dilation invariance of the RANS equation. Such tests carried out for DNS channel and pipe flow data indeed show that α = 3/2 in the viscous sublayer, and α = −1/2 in the central core. Both exponents can be derived from physical invariants at the boundaries as (33) for the Blasius equation, see Section 3.

Defect Power-law
The second situation is that the scaling symmetry of ℓ + M breaks down, but that ofl + M is maintained. This takes place in the bulk flow (for r = O(1)), where the mixing length tends to a finite value of the order of the radius of the pipe (or the half width of the channel). This finite length breaks the scaling symmetry and any dilation invariance. Using the outer variable r for the outer MME, we then have the following mapping: It can be expressed as This solution is postulated for describing the bulk flow. Integrating (29) once using (45), we obtain ℓ M = c ′ +F 2 r m /m, where c ′ is an integration constant and m is an exponent characterizing the power-law ofl m . Furthermore, using the boundary condition: ℓ M → 0 as r → 1, then c ′ = −F 2 /m. Thus, Note an interesting asymptotic behavior of (46): as r → 1, ℓ M → −F 2 (1 − r) = κy; thus we obtain κ = −F 2 . In other words, the classical Karman constant is a Lie-group invariant of the gradient of the mixing length which is valid not only for the overlap region, but for the whole bulk region also. Here, we call (45) the similarity of the second kind, and (46) the defect power-law. Note that (46) is obtained by employing a prolongation in the Lie-group analysis, which deserves a further comment. Bluman and Kumei [20] proposed to describe the similarity of the first kind as displaying invariant curve, while the similarity of the second kind (with the symmetry-breaking of the original scaling variable) displaying an invariant family of curves parameterized by ǫ. In other words, the latter case is characterized by an invariant manifold composed of a family of invariant curves. Here, we see a general rule, when an invariance property breaks down (i.e. original invariant curve is no longer invariant), another weaker form of the invariance (i.e. a family of invariant curves or invariant manifold) would be ensured. This feature can be formulated as a "weak-invariance-principle" which may guide the symmetry analysis of complex problems (with substantial symmetry-breaking).

Composite solution
Finally, consider the situation when both symmetries in ℓ + M andl + M break down, as sketched below:  In fact, this happens during a transition from one scaling state to another, when there exist several local scaling states in a system composed of a multi-layer structure, with exponents γ I (I = 1, 2, ...) each characterizing a layer. What would happen now? The weak-invariance-principle formulated above applies here, and suggests finding a new form which is yet left invariant under the group transformation; this would be then a generalized symmetry which describes a continuous transition from one power-law solution, ℓ +I M = F I 1 y +γI to another, ℓ +II M = F II 1 y +γII . Let us sketch this generalized invariance law using an abstract function Θ as: A specific ansatz is proposed as follows: where γ(y + ) = F 2 /F 1 and F 1 (y + ) = ℓ + M (y + )/y +γII are variables, while all other symbols are constants. In particular, n is a parameter specifying the transition sharpness. (48) implies that Θ = 0 is an invariant property under the dilation group transformation; it can be shown that it satisfies the so-called invariant surface condition (see appendix). In other words, Θ = 0 is a generalized symmetry -an invariant surface in the function space of ℓ + M ,l + M that is tangent to the Lie-group operator, in analogy to the simple power-law and defect power-law solutions. The postulated ansatz (49) imposes a specific invariance form, which would be constructive if (49) yields the actual empirical solution. Our subsequent studies, reported here and elsewhere, show that it is indeed the case.
One may wonder how (49) is motivated? In fact, Θ = 0 can be rewritten as For adjacent layers I and II, γ = F 2 /F 1 as a function of y + , must vary smoothly and monotonically from γ I to γ II . As a result, it is easy to verify that (50) has the desired asymptotic behavior: γ goes to γ I at one end and to γ II at the other end. Let y + c = F I 1 /F II 1 1/(γII −γI ) denote such a location in the middle of the transition that . Then, for y + ≫ y + c , F 1 → F II 1 , the r.h.s. of (50) goes to 1, as expected. For y + ≪ y + c , F II 1 /F 1 → (y + /y + c ) γII −γI , and the r.h.s. of (50) goes to zero for p = n(γ II − γ I ) ≫ 1. The latter can always be guaranteed by choosing an appropriate n (thus n specifies the sharpness of the transition between two adjacent layers). Therefore, (50) connects the two asymptotic states of power-law invariance with desired scalings.
From (50), we find an analytical expression for ℓ + M valid across the two layers as: where the transition function B is defined as It works as a blending procedure to switch between layers, which has been used as an interpolation function in many models, e.g. Batchelor [21] used it to describe the transition of the velocity structure function from the dissipative to inertial range. We derive it only by a continuity and smoothness argument on the Lie-group invariants. Here, we show that a dilation-invariant structure is behind B; it would be intriguing to investigate why this form is so generic and what determines the transition sharpness.

III. THE MULTI-LAYER DESCRIPTION OF TURBULENT CHANNEL AND PIPE FLOWS
Now, we apply the above analysis to turbulent channel and pipe flows. The notion of the multi-layer comes from empirical observations well known in the literature [15]. For the inner solution, one first meets with the so-called viscous sublayer very close to the wall; then, a buffer layer where coherent near-wall structures such as streamwise vortices dominate the dynamics [22]. Then is the overlap region between the inner solution and the outer solution. For the outer solution, there is a bulk flow, with a core around the center. With the Lie-group formalism developed above, these layers can be described analytically and quantified with parameters related to local scaling, transition location and transition sharpness. Below, we derive analytical expressions of the mixing length for all of these layers.

A. Inner flow
The viscous sublayer is assumed to hold a similarity of the first kind, and is analytically described by a power-law scaling of the mixing length as: This result can be justified by a near-wall expansion in the small fluctuation amplitude: u ∼ y + and v ∼ y +2 such that W + ∼ y +3 . Since S + ≈ 1 near the wall, ℓ + M ∼ y +3/2 , consistent with (53). Note that the symmetry transformation itself does not select the exponent 3/2; only further physical consideration related to the boundary constraint is able to offer a concrete proposal. Following the argument leading to (34), a similar constraint as (33) is S + = 1 at vanishing y + , which yields W + /ℓ +2 M = 1 = W + * /(ℓ + * M ) 2 ∼ (y + * ) 3 /(ℓ + * M ) 2 . Thus, α = 3/2 in (24). For large y + when (53) ceases to be valid, we postulate new layers with the power-law scaling of a distinct exponent, and identify them with the buffer-layer and the log-layer. The reason behind this assumption lies in the belief that the dilation-group invariance should be a robust property throughout the whole near-wall region. Indeed, an empirical study with the diagnostic function γ of the DNS channel flow data confirms our assumption with the following results: The exponent 2 for the buffer layer is believed to be related to coherent vortical structures, such as streamwise vortices whose core is mostly located at y + ≈ 25 [22], which will be studied in details elsewhere. Then, the generalized symmetry -ansatz (50) -yields a composite formula valid across both the sublayer and buffer layer as where valid across the buffer and log layer, where B b−l = (1 + (y + /y + buf ) p2 ) −1/p2 . The test of these formula is presented below.

B. Outer flow
As mentioned above, the bulk flow is characterized by the finite mixing length which leads to a dilation-symmetrybreaking. By assuming that the dilation-symmetry is still maintained in the gradient of the mixing length, we arrive at a defect power-law solution, Empirical study shows that m = 4 for channel and m = 5 for pipe. We have then a similarity of the second kind. Note that (58) smoothly match the log-layer ℓ log M = κy independent of m. The bulk solution (58) is also a discovery of the present theory. A verification by the DNS channel flow data is presented in Fig.2, where the compensated plot of ℓ + M by (58) clearly shows a plateau. Near the centreline (r → 0), we find an extra layer called core, where the mixing-length diverges as ℓ M ∝ 1/ √ r. This is consistent with a similarity solution of the first kind: where the exponent γ core = −1/2 is derived from an invariant central condition 1 r which, using (28), leads to This determination is again similar to the argument used to derive (34) in the Blasius equation. It is consistent with a central symmetry argument: as r → 0, W + ≈ r and S + ∝ r, hence, ℓ M ∝ r −1/2 . Furthermore, using (50), the Lie-group theory yields a matching from (59) to a constant mixing length in the bulk flow with scaling exponent zero, ℓ M → ℓ 0 = κ/m; and the composite solution from the bulk to the central core is thus where B b−c = (1 + (r/r core ) p3 ) −1/2p3 and Z c = B b−c (r = 1) is a normalization factor.
We have outlined above a novel Lie-group analysis for determining the mixing length profile. This analysis can be subjected to a test by DNS data. In other words, the invariance property of the mixing length function postulated in the analysis can go through scrutiny against DNS data which are the true Navier-Stokes solutions. Currently, there exist several DNS datasets for channel and pipe flows for Re τ up to thousands. We have selected two channel flows from Iwamoto et al. [16] at Re τ = 650 and Hoyas and Jimenez [17] at Re τ = 940, and one pipe flow of Wu and Moin [18] at Re τ = 1142, which are believed to be well-resolved. The data confirm the multi-layer similarity law for the mixing length, as below. Fig. 3 shows a clear evidence for the newly discovered bulk flow structure, 1 − r m (m = 4 for channel and 5 for pipe), using a compensated plot of the observed mixing length (e.g. divided by 1 − r m ), in which each of the inner layers and the core (divergent) layer are also very visible. The existence of the central core is due to a switch of the kinetic energy generation mechanism, from mean shear (production) to turbulent transport. Close to the centreline, the mean shear approaches zero, and finite energy dissipation is balanced by turbulent transport. To this end, a clear four-layer picture for turbulent channel and pipe flows is fully confirmed! The thickness of each layer is an important parameter, sketched in Fig. 3 (the demarcation for the core layer is sketched by y + core = (1 − r core )Re τ , hence increasing with Re), whose accurate determination from DNS data will be discussed elsewhere.

V. CONCLUSION AND DISCUSSIONS
In this paper, we demonstrate that a recently proposed multi-layer formula (1) for the mixing length [4] can be derived by a Lie-group symmetry analysis for the inner and outer MMEs, expressed not only in terms of the mean velocity, but also the mixing length and its gradient -the latter representing effects of fluctuations (Reynolds stress). Under the assumption about the symmetry of the mixing length and its gradient, three kinds of similarity solutions are derived. DNS data soundly validates our description, hence postulates.
A feature of the new theory is that it expresses the closure solutions in terms of three sets of parameters: scaling exponent, transition location (thickness in layered structure) and transition sharpness. These parameters are unique to the dilation-group invariance theory, which are, in general, functions of Re, but at large Re, they are believed to be Re-independent. Regarding the Re-similarity, note that in previous works [11,12], the Re (viscosity) effect is typically neglected by assuming a Re-independent solution of the mean velocity, which is obviously invalid because of the complex matching between the inner and outer solutions. This problem is solved here by writing down separately the inner and outer MME, and the Re-similarity is guaranteed by involving Re explicitly in the transformation. Aa large enough Re, the Re-similarity may be translated to Re-independence; in this case, the Lie-group invariant F 1 becomes truly invariant: where the first equation stems from the symmetry transformation, and the second from Re-independence, as mentioned before. Empirical study in a turbulent pipe shows that the Re independence only holds for Re τ > 5000, while a slow increase of the parameters y + buf and r core is found below this critical Re. The details will be discussed elsewhere. With a model of Re-dependence of the parameters, the theory (1) predicts the mean velocity profile for all Re.
What does the symmetry analysis using the Lie-group formalism add to our understanding more than making a direct postulate of power-law or other ansatz for ℓ M ? The answer is that a Lie-group formalism guarantees that the RANS equation remains invariant under the (dilation) group of transformation, hence any ansatz used here yields a similarity solution. In other words, the family of solutions based on power-law, defect power-law and the generalized invariance ansatz of ℓ M form a set of invariant solutions of the RANS equation, when ℓ M is treated as one dependent variables. Future study would clarify why ℓ M is so special to possess this invariance; once it is done, a new avenue of solving the closure problem of turbulence would be opened.
Thus, we speculate to have formulated a general framework for turbulence closure problem: searching for groupinvariant properties of relevant length order functions representative of the fluctuation structures (e.g. the mixing length and some other length order functions), and then solving for the mean velocity (which has a symmetry-breaking) from the MME. This procedure can be extended to other mean quantities by introducing relevant order functions which possess the required symmetry and guarantee the invariance of the relevant balance equations. Ultimately, whether such closure solutions are true solutions need to be verified against empirical data. If verified, the length order function (and its gradient) is indeed the right quantity, since the order function is, by definition, the quantity displaying distinct symmetry. Hence, we achieve a logically consistent framework.
Note that the Lie-group formalism is a necessary and powerful tool for the present success. Only with the present formalism, an analytical expression for the mixing length is obtained, for the first time, wholly motivated from symmetry consideration, i.e. dilation-group invariance. While power-law form for the mixing length can be motivated from simple scaling arguments, the present formalism accomplishes two additional aspects which would be difficult otherwise: a defect-power law form for the bulk flow and a series of matched forms derived from the generalized invariance ansatz. Note that the ansatz requires only simple assumptions about the continuity and smoothness for describing the variation of the local group invariants. Finally, it is remarkable that there exist only three distinct classes of similarity solutions, associated with the dilation group of the mixing length with its prolongation.
Finally, we note that the mixing length is only one example of the order functions, specifically representing the effect of the Reynolds stress in the mean momentum transfer. The complexity of turbulent flows is that they generally encompass a hierarchy set of order functions, each displaying a specific type of turbulent transfers. One will find other order functions responsible for turbulent transfer of kinetic energy, temperature, etc. Our preliminary study has already identified several other order functions in compressible turbulent boundary layers [24], rough pipe [25], and Reyleigh-Benard convection. Thus, the concept of the order function is general, applicable to many other flow systems, and the Lie-group formalism developed in the present work is also of general interest. The present analysis opens a promising direction for modeling the mean profiles in a wide class of turbulent wall-bounded flows.

VI. ACKNOWLEDGEMENTS
We thank Y. Wu, J. Chen, Y.S. Zhang, Y.Z. Wang and H.Y. Zou for helpful discussions. This work is supported by National Nature Science Fund 90716008 and by MOST 973 project 2009CB724100.