Application of the B-Determining Equations Method to One Problem of Free Turbulence

A three-dimensional model of the far turbulent wake behind a self-propelled body in a passively stratified medium is considered. The model is reduced to a system of ordinary differential equations by a similarity reduction and the B-determining equations method. The system of ordinary differential equations satisfying natural boundary conditions is solved numerically. The solutions obtained here are in close agreement with experimental data.

The far turbulent wake behind an axisymmetric body in a stratified medium is an example of a free shear flow. Sufficiently complete experimental data on the dynamics of turbulent wakes generated by moving bodies in stratified fluids were obtained by Lin and Pao [17].
The far turbulent wake behind an axisymmetric towed body in a linearly stratified medium was numerically simulated in [9]. Chernykh et al. [6] carried out the numerical simulation of the dynamics of turbulent wakes in a stable stratified medium based on hierarchy of second order closure models. Calculation results obtained in [6,9] are in close agreement with experimental data [17].
Similarity solutions for several turbulence models were constructed in [7,13,14,15]. In the current paper we consider three-dimensional semiempirical model of the far turbulent wake behind an axisymmetrical self-propelled body in a passively stratified medium (see [6,4,25] and the references therein). This paper is organized as follows. In Section 3 we determine the most general continuous classical symmetry group of the model and obtain the similarity reduction of the model. In Section 4 we use the B-determining equations (BDEs) method [1,12] to transform the reduced system into a system of ordinary differential equations (ODEs).
In the last section, we consider a boundary value problem for the system of ODEs. We use the modified shooting method and the asymptotic expansion of the solution in the vicinity of the singular point to solve this problem. Finally, computational results are given.

Model
The following three-dimensional semiempirical model of turbulence was constructed in [4,6,25] to calculate characteristics of the far turbulent wake behind an axisymmetric self-propelled body in a passively stratified medium: where e is the turbulent kinetic energy, ε is the kinetic energy dissipation rate, ρ 1 is the averaged density defect, and ρ 2 is the density fluctuation variance. All the unknown functions depend on x, y, and z. The quantities C e = 0.136, C ε = C e /δ, δ = 1.3, C ε2 = 1.92, C ρ = 0.208, C 1ρ = 0.087, C T = 1.25 are generally accepted empirical constant [8,21]. U 0 is the free stream velocity. The marching variable x in the equations (1)-(4) acts as time. This model is based on the three-dimensional parabolized system of averaged Navier-Stokes equations in the Oberbeck-Boussinesq approximation (see [5,26]) where U d = U 0 − U is the defect of the averaged longitudinal velocity component; U , V and W are the mean flow velocity component along x-, y-and z-axes, respectively; p 1 is the deviation from the hydrostatic pressure due to stratification ρ s (z); g is the gravity acceleration; ρ 1 is the averaged density defect: is the undisturbed fluid density assumed to be linear: ρ s (z) = ρ 0 (1 − az), a > 0 is a constant; the prime indicates the pulsating components; indicates averaging. In [9,26] the Reynolds stress tensor components u i u j , the turbulent flows u i ρ , and the density fluctuation variance ρ 2 are defined by the algebraic relations [21]. Since the flow in the far turbulent wake is considered, these relations are simplified as follows Differential transport equations [21] are used in [4,6,25] to determine the kinetic turbulent energy e, the kinetic energy dissipation rate ε, and the shear Reynolds stress v w : The turbulent viscosity coefficients in these equations are K ey = K y , K ez = K z , K εy = K ey /δ, (4) is an analogue of the equations (5)-(22) (for the diffusion approximation in a homogeneous fluid V = 0, W = 0, and g = 0). In what follows, we assume that the free stream velocity U 0 equals unity.

Similarity solution
The infinitesimal symmetry group [18,19] of the model (1)-(4) is spanned by the eight vector fields Available experimental data [10,17] and numerical calculations [9,20,27] show that the flow in the far turbulent wake can be considered to be close to a self-similar flow. We therefore consider the linear combination of scaling vector fields X 6 and X 7 The similarity variables associated with the infinitesimal generator Z are where E, G, H, and R are arbitrary functions. According to physical considerations (the influence of gravity is neglected in this problem), the functions E and G must be of the form We obtain the reduced system by introducing similarity variables. Using (23) and changing to polar coordinates ξ = r cos φ, η = r sin φ the reduced system becomes where E = E(r), G = G(r), H = H(r, φ), and R = R(r, φ). Here, and throughout, subscripts denote derivatives, so H r = ∂H/∂r, etc. Lie's classical method do not provide solution of the reduced system agreed with experimental data. We therefore use the BDEs method.

BDEs method
The concept of BDEs of a system of partial differential equations (PDEs) was introduced in [1,12]. Consider a scalar PDE Ω x, u, u 1 , u 2 , . . . = 0, where x = (x 1 , x 2 , . . . , x n ) denotes n independent variables, u denotes the dependent variable, and u k = ∂u k ∂x i 1 ∂x i 2 · · · ∂x i k denotes the set of coordinates corresponding to all k-th order partial derivatives of u with respect to x. In the BDEs method, an extension of the classical symmetry determining relations is made by incorporating an additional factor B(x, u, u 1 , u 2 , . . . ). For a scalar PDE (28), BDE is where h is a function of x, u, u 1 , u 2 , . . . ; D α = D α 1 x 1 · · · D αn xn , D x i is the total x i derivative; α is a multi-index. Equality (29) must hold for all solutions of (28). Now we use the BDEs method to reduce (26) and (27) to some ODEs. Consider more general equation than (26) where A(r), B(r), and C(r) are arbitrary functions. BDE corresponding to (30) is Here, and throughout D φ , D r are the operators of total differentiation with respect to φ and r.
The function h may depend on r, φ, H and derivatives of H. The functions b 1 (r, φ) and b 2 (r, φ) are to be determined together with the function h. Note that if we let in (31) we obtain the classical determining equation [18,19]. For simplicity, we assume that a solution of (31) is independent of r and partial derivatives of H with respect to r Substituting (32) into (31) gives a polynomial equation for derivatives of the fourth order. This polynomial must identically vanish. We can express the derivatives H rrφφ , H φφφφ , H rφφ , H φφφ , and H φφ using (30). The coefficient of H rrr implies b 1 (r, φ) = A(r). As a result, the left side of (31) is a polynomial in H rr and H rφ . Collecting similar terms we obtain Hence Substituting the functions b 1 , b 2 , and h 1 into the left side of (31) we obtain the polynomial with respect to H r and H φ . This polynomial must identically vanish. Collecting similar terms we have Thus, where h 4 is an arbitrary constant.
Clearly, that the Riccati equation (33) has the partial solution h 2 = tan φ for h 4 = 1/2. Thus we find the solution of (31) The corresponding differential constraint h = 0 has the general solution where H 1 and H 2 are arbitrary functions. Next we use (34) and consider the equation (27) in more general form where K(r), L(r), M (r), N (r), and P (r) are arbitrary functions. The BDEs method applied to the equation (35) gives rise to the following results: Integrating differential constraint h = 0 corresponding to the BDE solution (34), we find where R 1 (r) and R 2 (r) are arbitrary functions. Using (34) and (36) we obtain the following corollary in terms of variables x, y, and z.
Corollary 1. The following expressions for the unknown functions allow us to reduce the model (1)-(4) to a system of ODEs.

Calculation results
The system of ODEs (39)-(43) has to satisfy the conditions Conditions (44) take into account flow symmetry with respect to the OX axis. The boundary conditions (45) imply that all functions take zero values outside the turbulent wake.
The system of ODEs (39)-(43) satisfying boundary conditions (44), (45) was solved numerically. Additional difficulties are caused by the fact that the coefficients of ODEs have singularities. The problem was solved by the modified shooting method and asymptotic expansion of the solution in the vicinity of the singular point [3,15] The value of α is taken to be 0.23 in accordance with experimental data [9,17]. The results for the problem solution are illustrated in Figs. 1 and 2. Fig. 1 shows the profiles of the functions E/E 0 , G/G 0 , H, and R/R 0 as ξ = 0, where subscript 0 denotes the axial value. The functions E/E 0 , G/G 0 , H, and R/R 0 are plotted in Fig. 2. The functions E/E 0 and G/G 0 are bellshaped and determine shapes of the normalized turbulent kinetic energy and the normalized kinetic energy dissipation rate respectively. Similarly, H and R/R 0 determine shapes of the average density defect and the normalized density fluctuation varience respectively.
The function H(0, η) characterizing the degree of fluid mixing in the turbulent wake is given in Fig. 1c. The maximum value of this function equals 0.258. This is in consistent with the present notions of incomplete fluid mixing in the wakes [23,24].
In Fig. 3 adapted from [6] the normalized values of the turbulent energy along the wake axis e 1/2 0 /U 0 = e(x, 0, 0)/U 0 are compared with experimental data [17], computational results [9] and results of numerical simulation based on two semi-empirical turbulence models (Model1 and Model2 in [4,6]). The coordinate x is normalized by the body diameter D. The results obtained here are in close agreement with Lin and Pao's experimental data.