Explicit Algebraic Stress Model for Three-Dimensional Turbulent Buoyant Flows Derived Using Tensor Representation

Ronald M. C. So

Department of Mechanical Engineering, The Hong Kong Polytechnic University, Hong Kong, China.

DOI: 10.4236/jamp.2022.104082                

Abstract

An explicit algebraic stress model (EASM) has been formulated for two-dimensional turbulent buoyant flows using a five-term tensor representation in a prior study. The derivation was based on partitioning the buoyant flux tensor into a two-dimensional and a three-dimensional component. The five-term basis was formed with the two-dimensional component of the buoyant flux tensor. As such, the derived EASM is limited to two-dimensional flows only. In this paper, a more general approach using a seven-term representation without partitioning the buoyant flux tensor is used to derive an EASM valid for two- and three-dimensional turbulent buoyant flows. Consequently, the basis tensors are formed with the fully three-dimensional buoyant flux tensor. The derived EASM has the two-dimensional flow as a special case. The matrices and the representation coefficients are further simplified using a four-term representation. When this four-term representation model is applied to calculate two-dimensional homogeneous buoyant flows, the results are essentially identical with those obtained previously using the two-dimensional component of the buoyant flux tensor. Therefore, the present approach leads to a more general EASM formulation that is equally valid for two- and three-dimensional turbulent buoyant flows.

Keywords

Explicit Algebraic Stress Model, Buoyant Flows, Tensor Representation

1. Introduction

Previously, tensor representation theory has been used by So et al. (2002) Ref. [1], to derive an explicit algebraic stress model (EASM) for homogeneous turbulent buoyant shear flows. Two different approaches have been adopted: one based on the proposal of Gatski and Speziale (1993) Ref. [2], and another employed the more general approach of Junger and Gatski (1998) Ref. [3]. The starting point of their analysis is the normalized anisotropic Reynolds stress equation simplified assuming homogeneous equilibrium turbulence. Thus derived, the equation is quadratic in the anisotropic stress tensor $b=b_{ij}=\left({\tau }_{ij}-\left(2/3\right)k{\delta }_{ij}\right)/2k$ which is symmetric and traceless. Here, k = t_ii/2 is the turbulent kinetic energy, ${\tau }_{ij}=\overline{u_i{u}_j}$ is the kinematic Reynolds stress tensor, u_i is the i^th component of the fluctuating velocity and the overbar is used to denote ensemble average. The resulting equation is implicit in b. However, tensor representation theory can be used to render this equation explicit; thus allowing explicit algebraic stress models to be developed for buoyant shear flows.

Up to this point, the derivation is quite general because no assumptions have been made to limit the model for two-dimensional (2-D) flows. As for the assumptions of homogeneity and equilibrium, they were made to simplify the b equation so that it can be reduced to an algebraic equation. Gatski and Speziale (1993) Ref. [2], and Jongen and Gatski (1998) Ref. [3] also invoked these assumptions in their derivation of the non-buoyant EASM. In principle, the representation theory can be used to formulate the tensor representation for b if the basis tensors chosen are traceless. For incompressible flows, the kinematic shear strain rate tensor $S=S_{ij}=\left(\partial U_i/\partial x_j+\partial U_j/\partial x_i\right)/2$ and the rotation rate tensor $W=W_{ij}=\left(\partial U_i/\partial x_j-\partial U_j/\partial x_i\right)/2$ are traceless. Since the buoyant flux tensor, $\Gamma ={\Gamma }_{ij}=\left(G_{ij}-2G{\delta }_{ij}/3\right)/2G$, that appears naturally in the b equation for buoyant flows is also traceless, it would appear that S, W and $\Gamma$ could be used as basis tensors to form the tensor representation for b. Here, $G=G_{ij}=-\beta g_i\overline{u_j\theta }-\beta g_j\overline{u_i\theta }$ is the buoyant production tensor, G = G_ii/2 is the buoyant production of k, g_i is the ith component of the gravitational vector, U_i is the ith component of the mean flow velocity vector, x_i is the ith component of the space vector, θ is the fluctuating temperature, β is the coefficient of thermal expansion of the fluid and δ_ij = 1 for i = j and 0 for $i\ne j$. The instantaneous velocity vector ${\tilde{u}}_i$ and temperature $\tilde{\theta }$ are decomposed into ensemble mean and fluctuating parts, i.e., ${\tilde{u}}_i=U_i+u_i$ and $\tilde{\theta }=\Theta +\theta$.

If the assumption of 2-D flow were invoked, the properties of the tensors S, W and $\Gamma$ will become quite different. While the tensors S and W become 2-D with elements in one row and one column being identically zero, $\Gamma$ remains a 3-D tensor and its use to form the basis tensors for b with S and W becomes very complicated. Therefore, an assumption was made to split the buoyant flux tensor into two parts, a symmetric traceless tensor f = f_ij that has the same properties as S and W for 2-D homogeneous shear flows, and a three-dimensional tensor N = N_ij. This was accomplished through the introduction of a 2-D tensor ${\delta }_{ij}^{\left(2d\right)}$, first proposed by Pope (1975) Ref. [4], to define $f_{ij}=\left(G_{ij}-{\delta }_{ij}^{\left(2d\right)}G\right)/G$ and $N_{ij}={\delta }_{ij}^{\left(2d\right)}/2-{\delta }_{ij}/\text{3}$. The 2-D tensor ${\delta }_{ij}^{\left(2d\right)}$ is defined as

${\delta }_{ij}^{\left(2d\right)}=\left|\begin{array}{ccc}1& 0& 0\\ 0& 0& 0\\ 0& 0& 1\end{array}\right|$. (1)

Thus partitioned, f has properties similar to that of S and W in the case of 2-D homogeneous shear flows, i.e., elements in one row and one column are identically zero, and the theoretical derivation of the EASM can be carried out to give an explicit expression for b.

The tensor representation thus proposed limits the applicability of the EASM to 2-D flows only. In the EASM derived by So et al. (2002) Ref. [1], a five-term representation was assumed. Among the five basis tensors, three involve S and W, while two involve S, W and f. It is the presence of these last two basis tensors that render the representation not applicable to 3-D flows. Therefore, the derived EASM needs further extension to 3-D flows. It is recognized that the limitation stems from the partitioning of $\Gamma$ into f and N. Therefore, if the derived EASM is to be equally applicable to 3-D flows, a way has to be found so that $\Gamma$ can be chosen as one of the tensors used to form the basis tensors. This suggests that the basis tensors should be formed from S, W and $\Gamma$.

The present paper is an extension of the model of So et al. (2002) Ref. [1] to 3-D flows. However, the analysis will invoke the incompressibility assumption because S and W will remain as traceless tensors even for 3-D flows. Instead of splitting $\Gamma$ into a 2-D symmetric traceless tensor f and a 3-D tensor N, the present approach uses S, W and $\Gamma$ as the basis tensors to formulate a representation for b. A seven-term representation was derived for the general 3-D case. For pure shear flows without buoyancy, the EASM reduces to the form proposed by Gatski and Speziale (1993) Ref [2], i.e. a three-term representation. The seven-term representation of b is simplified to a four-term representation involving S, W and $\Gamma$ as basis tensors, and its performance is found to be essentially identical to the EASM of So et al. (2002) Ref. [1] for 2-D homogeneous buoyant shear flows with a five-term representation involving S, W and f as basis tensors.

2. Mathematical Formulation

The EASM of So et al. (2002) Ref. [1] was derived from the Reynolds stress equation assuming the validity of Boussinesq approximation and equilibrium turbulence. The governing equation of the Reynolds stress anisotropy tensor as given in So et al. (2002) Ref. [1] is

$-\frac{1}{g\lambda }b-a_3\left(bS+Sb-\frac{2}{3}\left\{bS\right\}I\right)+a_2\left(bW-Wb\right)+\frac{a_4}{\lambda }\left(b^2-\frac{1}{3}\left\{b^2\right\}I\right)=a_{1}S+L$,(2)

where $L=d/\lambda +\left[a_{\text{6}}\left(G/\epsilon \right)/\lambda \right]\Gamma$ and $d=d_{ij}$ is the anisotropic dissipation rate tensor. In Equation (2), $1/g=\left(C_1^{*}/2+1\right)\left(\tilde{P}/\epsilon \right)+C_1/2-1+\left(G/\epsilon \right)={\alpha }^{\prime }\left(\tilde{P}/\epsilon \right)+{\beta }^{\prime }+\left(G/\epsilon \right)$, $a_1=\left(4/3-C_2\right)/2$, $a_2=\left(2-C_4\right)/2$, $a_3=\left(2-C_3\right)/2$, $a_4=C_5/2$, $a_6=C_6-1$, $\lambda =k/\epsilon$, ${\alpha }^{\prime }=C_1^{*}/2+1$, and ${\beta }^{\prime }=C_1/2-1$. The model constants are specified as $C_1=3.4$ , $C_2=0.36$, $C_3=1.25$, $C_4=0.40$, $C_5=4.2$, $C_1^{*}=1.8$ and $C_6=0.3$. Also, $\tilde{P}=P_{ii}/2$ is the shear production of k, e is the dissipation rate of k, and $P_{ij}=-\overline{u_i{u}_k}\left(\partial U_j/\partial x_k\right)-\overline{u_j{u}_k}\left(\partial U_i/\partial x_k\right)$ is the shear stress production tensor. Equation (2) is implicit in b and further derivation is needed to render it explicit. In deriving Equation (2), the pressure-strain model of Speziale et al. (1991) Ref. [5] has been assumed. It can be written here as:

$\begin{array}{c}{\Pi }_{ij}=-\left(C_1+C_1^{*}\tilde{P}\right)\epsilon b_{ij}+C_{2}k{S}_{ij}+C_{3}k\left(b_{ik}{S}_{jk}+b_{jk}{S}_{ik}-\frac{2}{3}{b}_{mn}{S}_{mn}{\delta }_{ij}\right)\\ -C_{4}k\left(b_{ik}{W}_{jk}-b_{jk}{W}_{ik}\right)+C_5\epsilon \left(b_{ik}{b}_{kj}-\frac{1}{3}{b}_{mn}{b}_{nm}{\delta }_{ij}\right)-C_6\left(G_{ij}-\frac{2}{3}G{\delta }_{ij}\right)\end{array}$ (3)

Following the work of Rivlin and Ericksen (1955) Ref. [6], let $T^{\left(\text{1}\right)},T^{\left(2\right)},\cdots ,T^{\left(N\right)}$ be any symmetric 3 ´ 3 traceless tensors formed from the tensors S, W and $\Gamma$. A linear relation of the form,

$b={\displaystyle \underset{n=1}{\overset{N}{\sum }}{Q}_n{T}^{\left(n\right)}}$, (4)

can be established between b and the tensors $T^{\left(n\right)}\left(n=1,2,\cdots ,N\right)$, which are not necessarily linearly independent. With b being a function of S, W and $\Gamma$, there will be 41 tensors forming this basis altogether. According to Jongen and Gatski (1998) Ref. [3], for pure shear flows, three of the basis tensors are sufficient to give a relatively good approximation for 3-D flows. In view of this, a seven-term representation for b with four terms involving $\Gamma$ is assumed in the present analysis. These seven basis tensors are given by

$\begin{array}{l}{T}^{\left(1\right)}=S\text{,}{T}^{\left(2\right)}=SW-WS,T^{\left(\text{3}\right)}=S^{\text{2}}-\frac{1}{3}\left\{S^2\right\}I,\\ T^{\left(4\right)}=\Gamma \text{,}{T}^{\left(5\right)}=\Gamma W-W\Gamma ,T^{\left(6\right)}={\Gamma }^{\text{2}}-\frac{1}{3}\left\{{\Gamma }^2\right\}I,\\ T^{\left(7\right)}=\Gamma S+S\Gamma -\frac{2}{3}\left\{\Gamma S\right\}I.\end{array}$ (5)

Assuming d_ij= 0 and a₄= 0 and forming the scalar product of Equation (2) with each of the tensors $T^{\left(m\right)},m=1,2,\cdots ,7$, and using (4), the following equation is obtained,

$\begin{array}{l}-\frac{1}{g\lambda }{\displaystyle \underset{n=1}{\overset{7}{\sum }}{Q}_n\left(T^{\left(n\right)},T^{\left(m\right)}\right)}-2a_3{\displaystyle \underset{n=1}{\overset{7}{\sum }}{Q}_n\left(T^{\left(n\right)}S,T^{\left(m\right)}\right)}+2a_2{\displaystyle \underset{n=1}{\overset{7}{\sum }}{Q}_n\left(T^{\left(n\right)}W,T^{\left(m\right)}\right)}\\ =a_1\left(S\text{,}{T}^{\left(m\right)}\right)+\frac{a_6\left(G/\epsilon \right)}{\lambda }\left(\Gamma \text{,}{T}^{\left(m\right)}\right)\end{array}$ (6)

The above equation can be written in a compact form as,

${\displaystyle \underset{n=1}{\overset{7}{\sum }}{A}_{mn}^{\left(7\right)}{Q}_n}=\left(R,T^{\left(m\right)}\right)$. (7)

where $A_{mn}^{\left(7\right)}$ and (R, T^(m)) are given by,

$A_{mn}^{\left(7\right)}=-\frac{1}{g\lambda }\left(T^{\left(n\right)},T^{\left(m\right)}\right)-2a_3\left(T^{\left(n\right)}S,T^{\left(m\right)}\right)+2a_2\left(T^{\left(n\right)}W,T^{\left(m\right)}\right)$, (8a)

$\left(R,T^{\left(m\right)}\right)=a_1\left(S\text{,}{T}^{\left(m\right)}\right)+\frac{a_6\left(G/\epsilon \right)}{\lambda }\left(\Gamma \text{,}{T}^{\left(m\right)}\right)$. (8b)

The solution of (7) can be obtained by determining the matrices (T⁽ⁿ⁾S, T^(m)), (T⁽ⁿ⁾W, T^(m)), (T⁽ⁿ⁾, T^(m)) and (R, T^(m)). Here, a distinction can be made depending on whether the solution is sought for 2-D or 3-D flow. In the following, two different solutions are presented: one for 3-D flows with a 7-term basis and another for 2-D and 3-D flows with a 4-term basis. The 7-term representation gives rise to matrices that are very complicated, and it is not possible to simplify them to give analytical expressions for the coefficients Q_n. The choice of a 4-term representation, on the other hand, allows simplifications to be made to deduce analytical expressions for Q_n, and yet permits b to reduce to the buoyant case without shear identically. That way, the solution can be compared with the 5-term representation previously derived by So et al. (2002) Ref. [1] with S, W and f as the basis tensors.

3. Solution with a 7-Term Basis

According to Clapham and Nicholson (2009) Ref. [7], in linear algebra, the Cayley-Hamilton theorem states that every square matrix over a commutative ring satisfies its own characteristic equation. In view of this property, it is helpful to invoke the Cayley-Hamilton theorem in the present analysis. Therefore, using the Cayley-Hamilton theorem identity for their entries, the matrices (T⁽ⁿ⁾S, T^(m)), (T⁽ⁿ⁾W, T^(m)), (T⁽ⁿ⁾, T^(m)) and (R, T^(m)) can be reduced to:

(9a)

(9b)

$\left(T^{\left(n\right)},T^{\left(m\right)}\right)=\left(\begin{array}{ccccccc}{\eta }_1& 0& {\eta }_3& {\eta }_8& 2{\eta }_{14}& {\eta }_9& 2{\eta }_{10}\\ 0& {\eta }_1{\eta }_2-6{\eta }_5& 0& -2{\eta }_{14}& {\eta }_2{\eta }_8-6{\eta }_{15}& -2{\eta }_{16}& 2{\eta }_{17}\\ {\eta }_3& 0& \frac{1}{6}{\eta }_1^2& {\eta }_{10}& -2{\eta }_{17}& {\eta }_{11}-\frac{1}{3}{\eta }_1{\eta }_6& \frac{1}{3}{\eta }_1{\eta }_8\\ {\eta }_8& -2{\eta }_{14}& {\eta }_{10}& {\eta }_6& 0& {\eta }_7& 2{\eta }_9\\ 2{\eta }_{14}& {\eta }_2{\eta }_8-6{\eta }_{15}& -2{\eta }_{17}& 0& {\eta }_2{\eta }_6-6{\eta }_{13}& 0& 2{\eta }_{16}\\ {\eta }_9& -2{\eta }_{16}& {\eta }_{11}-\frac{1}{3}{\eta }_1{\eta }_6& {\eta }_7& 0& \frac{1}{6}{\eta }_6^2& \frac{1}{3}{\eta }_6{\eta }_8\\ 2{\eta }_{10}& 2{\eta }_{17}& \frac{1}{3}{\eta }_1{\eta }_8& 2{\eta }_9& 2{\eta }_{16}& \frac{1}{3}{\eta }_6{\eta }_8& \frac{2}{3}{\eta }_8^2+{\eta }_1{\eta }_6-2{\eta }_{11}\end{array}\right)$, (9c)

$\left(R,T^{\left(m\right)}\right)=\left(\begin{array}{c}{a}_1{\eta }_1+\frac{a_6\left(G/\epsilon \right)}{\lambda }{\eta }_8\\ -2\frac{a_6\left(G/\epsilon \right)}{\lambda }{\eta }_{14}\\ a_1{\eta }_3+\frac{a_6\left(G/\epsilon \right)}{\lambda }{\eta }_{10}\\ a_1{\eta }_8+\frac{a_6\left(G/\epsilon \right)}{\lambda }{\eta }_6\\ 2a_1{\eta }_{14}\\ a_1{\eta }_9+\frac{a_6\left(G/\epsilon \right)}{\lambda }{\eta }_7\\ 2a_1{\eta }_{10}+\frac{2a_6\left(G/\epsilon \right)}{\lambda }{\eta }_9\end{array}\right)$, (9d)

where the scalar invariants h_i’s are given by

$\begin{array}{l}{\eta }_1=\left\{S^2\right\},{\eta }_2=\left\{W^2\right\},{\eta }_3=\left\{S^3\right\},{\eta }_4=\left\{W^{2}S\right\},{\eta }_5=\left\{W^2{S}^2\right\},{\eta }_6=\left\{{\Gamma }^2\right\},\\ {\eta }_7=\left\{{\Gamma }^3\right\},{\eta }_8=\left\{\Gamma S\right\},{\eta }_9=\left\{{\Gamma }^{2}S\right\},{\eta }_{10}=\left\{S^2\Gamma \right\},{\eta }_{11}=\left\{S^2{\Gamma }^2\right\},{\eta }_{12}=\left\{W^2\Gamma \right\},\\ {\eta }_{13}=\left\{W^2{\Gamma }^2\right\},{\eta }_{14}=\left\{WS\Gamma \right\},{\eta }_{15}=\left\{W^{2}S\Gamma \right\},{\eta }_{16}=\left\{WS{\Gamma }^2\right\},{\eta }_{17}=\left\{W\Gamma S^2\right\},\\ {\eta }_{18}=\left\{W{S}^2\Gamma S\right\},{\eta }_{19}=\left\{W{\Gamma }^{2}S\Gamma \right\},{\eta }_{20}=\left\{W{S}^2{\Gamma }^2\right\},{\eta }_{21}=\left\{W^2{S}^2\Gamma \right\},\\ {\eta }_{22}=\left\{W^2{\Gamma }^{2}S\right\},{\eta }_{23}=\left\{W^2\Gamma WS\right\}.\end{array}$ (10)

Matrix A⁽⁷⁾ in Equation (8a), its inverse, and the scalar coefficients of Equation (4) are deduced by making use of the MAPLE software. The coefficients thus determined will allow Equation (4) to reduce correctly to the pure shear flow limit when $\Gamma$ is set to zero and the pure buoyant flow limit when S and W vanish in Equation (2).

4. Solution with a 4-Term Basis

In the case of a 4-term representation, solution can be sought depending on whether 2-D or 3-D flows are considered. For 3-D flows, only the first four rows and columns of Equation (9) have to be retained and no further simplifications to the scalar coefficients h_i can be made. On the other hand, for 2-D flows, the 4 * 4 matrices thus obtained are reduced to simpler forms as some of the h_i’s from Equation (10) in this case simplify to

$\begin{array}{l}{\eta }_3=0,{\eta }_4=0,{\eta }_5=\frac{1}{2}{\eta }_1{\eta }_2,{\eta }_{15}=\frac{1}{2}{\eta }_2{\eta }_8,\\ {\eta }_{17}=0,{\eta }_{18}=-\frac{1}{2}{\eta }_1{\eta }_{14},{\eta }_{20}=0,{\eta }_{21}=\frac{1}{2}{\eta }_2{\eta }_{10},\\ {\eta }_{22}=\frac{1}{2}{\eta }_2{\eta }_9,{\eta }_{23}=\frac{1}{2}{\eta }_2{\eta }_{14},\end{array}$ (11)

while others are as given in Equation (10). Thus simplified, the matrices, (T⁽ⁿ⁾, T^(m)), (R, T^(m)), (T⁽ⁿ⁾S, T^(m)) and (T⁽ⁿ⁾W, T^(m)), in the 2-D case become

${\left(T^{\left(n\right)},T^{\left(m\right)}\right)}^{\left(2d\right)}=\left(\begin{array}{cccc}{\eta }_1& 0& 0& {\eta }_8\\ 0& -2{\eta }_1{\eta }_2& 0& -2{\eta }_{14}\\ 0& 0& \frac{1}{6}{\eta }_1^2& {\eta }_{10}\\ {\eta }_8& -2{\eta }_{14}& {\eta }_{10}& {\eta }_6\end{array}\right)$, (12a)

${\left(R,T^{\left(m\right)}\right)}^{\left(2d\right)}=\left(\begin{array}{c}{a}_1{\eta }_1+\frac{a_6\left(G/\epsilon \right)}{\lambda }{\eta }_8\\ -2\frac{a_6\left(G/\epsilon \right)}{\lambda }{\eta }_{14}\\ \frac{a_6\left(G/\epsilon \right)}{\lambda }{\eta }_{10}\\ a_1{\eta }_8+\frac{a_6\left(G/\epsilon \right)}{\lambda }{\eta }_6\end{array}\right)$, (12b)

${\left(T^{\left(n\right)}S,T^{\left(m\right)}\right)}^{\left(2d\right)}=\left(\begin{array}{cccc}0& 0& \frac{1}{6}{\eta }_1^2& {\eta }_{10}\\ 0& 0& 0& 0\\ \frac{1}{6}{\eta }_1^2& 0& 0& \frac{1}{6}{\eta }_1{\eta }_8\\ {\eta }_{10}& 0& \frac{1}{6}{\eta }_1{\eta }_8& {\eta }_9\end{array}\right)$, (12c)

${\left(T^{\left(n\right)}W,T^{\left(m\right)}\right)}^{\left(2d\right)}=\left(\begin{array}{cccc}0& -{\eta }_1{\eta }_2& 0& -{\eta }_{14}\\ {\eta }_1{\eta }_2& 0& 0& {\eta }_2{\eta }_8\\ 0& 0& 0& 0\\ {\eta }_{14}& -{\eta }_2{\eta }_8& 0& 0\end{array}\right)$. (12d)

The matrix A⁽⁴⁾ in the 2-D case is reduced to

${\text{A}}^{\left(4\right)}=\left(\begin{array}{cccc}-\frac{{\eta }_1}{g\lambda }& -2a_2{\eta }_1{\eta }_2& -\frac{1}{3}{a}_3{\eta }_1^2& -\frac{{\eta }_8}{g\lambda }-2a_3{\eta }_{10}-2a_2{\eta }_{14}\\ 2a_2{\eta }_1{\eta }_2& 2\frac{{\eta }_1{\eta }_2}{g\lambda }& 0& 2\frac{{\eta }_{14}}{g\lambda }+2a_2{\eta }_2{\eta }_8\\ -\frac{1}{3}{a}_3{\eta }_1^2& 0& -\frac{1}{6}\frac{{\eta }_1^2}{g\lambda }& -\frac{{\eta }_{10}}{g\lambda }-\frac{1}{3}{a}_3{\eta }_1{\eta }_8\\ -\frac{{\eta }_8}{g\lambda }-2a_3{\eta }_{10}+2a_2{\eta }_{14}& 2\frac{{\eta }_{14}}{g\lambda }-2a_2{\eta }_2{\eta }_8& -\frac{{\eta }_{10}}{g\lambda }-\frac{1}{3}a{}_3\eta {}_1{\eta }_8& -\frac{{\eta }_6}{g\lambda }-2a_3{\eta }_9\end{array}\right)$, (13)

The inverse of A⁽⁴⁾ and the scalar coefficients in this case are again obtained using the MAPLE software and the resultant coefficients Q₁ to Q₄ can be simplified to

$Q_1=\frac{3a_{1}g\lambda }{D_1}+\frac{Q_{12}g\lambda a_6\left(G/k\right)}{D_1{D}_2{\eta }_1}$, (14a)

$Q_2=\frac{3a_1{a}_2{g}^2{\lambda }^2}{D_1}+\frac{Q_{22}{g}^2{\lambda }^2{a}_6\left(G/k\right)}{D_1{D}_2{\eta }_1}$, (14b)

$Q_3=-\frac{6a_1{a}_3{g}^2{\lambda }^2}{D_1}-\frac{6Q_{32}{a}_3{g}^2{\lambda }^2{a}_6\left(G/k\right)}{D_1{D}_2{\eta }_1}$, (14c)

$Q_4=\frac{Q_{42}g\lambda a_6\left(G/k\right)}{D_2}$, (14d)

where

$D_1=6a_2^2{g}^2{\lambda }^2{\eta }_2+2a_3^2{g}^2{\lambda }^2{\eta }_1-3$, (15a)

$D_2=2a_{3}g\lambda {\eta }_1^2{\eta }_2{\eta }_9-4a_{3}g\lambda {\eta }_1{\eta }_2{\eta }_8{\eta }_{10}-{\eta }_1{\eta }_2{\eta }_8^2+{\eta }_2{\eta }_6{\eta }_1^2-6{\eta }_2{\eta }_{10}^2+2{\eta }_1{\eta }_{14}^2$, (15b)

and Q₁₂, Q₂₂, Q₃₂ and Q₄₂ are given in Appendix I.

Once derived, b depends on $\Gamma$, λ, $\tilde{P}/\epsilon$, and G. A model for $\overline{u_i\theta }$ is required if $\Gamma$ is to be evaluated. In So et al. (2002) Ref. [1], $\overline{u_i\theta }$ was evaluated by solving its transport equation with a suitable model for the pressure scrambling term and the modeled transport equations for $\overline{{\theta }^2}$ and ${\epsilon }_{\theta }$, the temperature variance and its dissipation rate. The variable λ can be determined by solving the modeled transport equations for k and ε. Therefore, $\Gamma$, G and λ are all known, and the only unknown left is $\tilde{P}/\epsilon$. An equation for $\tilde{P}/\epsilon$ can be derived following the approach taken in So et al. (2002) Ref. [1]. Again, the equation is found to be cubic, and it can be written as

${{\alpha }^{\prime }}^2{\left(\frac{\tilde{P}}{\epsilon }\right)}^3+A_1{\left(\frac{\tilde{P}}{\epsilon }\right)}^2+A_2\left(\frac{\tilde{P}}{\epsilon }\right)+A_3=0$, (16)

$A_1=2{\alpha }^{\prime }\left(\frac{G}{\epsilon }\right)+2{\alpha }^{\prime }{\beta }^{\prime }$, (17a)

$A_2={\left(\frac{G}{\epsilon }\right)}^2+\left(2{\beta }^{\prime }-2{\alpha }^{\prime }\lambda {\eta }_8{a}_6\right)\frac{G}{\epsilon }+{{\beta }^{\prime }}^2-2a_2^2{\lambda }^2{\eta }_2-\frac{2}{3}{a}_3^2{\lambda }^2{\eta }_1-2a_1{\alpha }^{\prime }{\lambda }^2{\eta }_1$, (17b)

$\begin{array}{c}{A}_3=-2\lambda {\eta }_8{a}_6{\left(\frac{G}{\epsilon }\right)}^2-\left[2\lambda a_6\left({\beta }^{\prime }{\eta }_8-2a_3\lambda {\eta }_{10}+2a_2\lambda {\eta }_{14}\right)+2a_1{\lambda }^2{\eta }_1\right]\left(\frac{G}{\epsilon }\right)\\ -2a_1{\beta }^{\prime }{\lambda }^2{\eta }_1.\end{array}$ (17c)

Details of the derivation of Equation (16) are given in So et al. (2002) Ref. [1]. Since the roots of a cubic equation can be found in any mathematical handbook, such as Ref. [7], they are not given here. If b is to be explicitly determined, the heat flux model of So et al. (2002) Ref. [1] can again be adopted.

5. Discussion

The results given in Equation (4), Equation (14) and Equation (15) reduce to the pure shear flow and the pure buoyant flow limit when $\Gamma$ and S are set to zero individually. When $\Gamma$ is set to zero, Equation (4) yields a 3-term representation with Q₄ = 0. The other scalar coefficients reduce to $Q_1=3a_{1}g\lambda /D_1$, $Q_2=a_{2}g\lambda Q_1$ and $Q_3=-6a_1{a}_3{g}^2{\lambda }^2/D_1$, which are the same as those given by Gatski and Speziale (1993) Ref. [2] and Jongen and Gatski (1998) Ref. [3]. Furthermore, when Equation (4), Equation (14) and Equation (15) are substituted into Equation (2), the equation is satisfied identically. As for the pure buoyant limit, setting S to zero reduces Equation (4) to a 1-term representation and the result is given by $b=-a_6\left(G/k\right)g\lambda \Gamma$, because W is also zero. This is identical to the expression obtained from Equation (2) when S = W = 0 is substituted. In other words, Equation (4), Equation (14) and Equation (15) reduce to the two limiting forms correctly. Therefore, it remains to be shown that when Equation (4), Equation (14) and Equation (15) are used to calculate 2-D homogeneous buoyant shear flows the results are essentially identical to those given previously by So et al. (2002) Ref. [1].

The EASM previously derived by So et al. (2002) Ref. [1] and denoted as EASM/GS or EASM/BJG, has been validated against the direct numerical simulation (DNS) data of Gerz et al. (1989) Ref. [8] and also against calculations using a two-equation turbulence model, i.e., an EASM without accounting for buoyancy effects and a full Reynolds stress model assuming the SSG model of Speziale et al. 1991 Ref. [5] for the pressure strain term. Since EASM/GS and EASM/BJG give identical results, from this point on, the model of So et al. (2002) Ref. [1] will be referred to as EASM/BJG. This comparison reveals that, while accounting for buoyancy effects is important in an EASM, the full Reynolds stress model gives the best prediction of the DNS data because of its ability to model the history of the flow. The present EASM, denoted as EASM/Γ, is essentially a more general version of EASM/BJG. Therefore, its validity can be verified by comparing with the calculations of EASM/BJG and the DNS data. Calculated results from the two-equation and the Reynolds stress model will not be

Figure 1. Comparison of the calculated k with DNS data. (a) is for the Ri = 0.1 case, and (b) is for the Ri = 1.0 case.

shown because they have already been compared to EASM/BJG in So et al. (2002) Ref. [1]. The details of the calculations, including the assumption and the normalization used to reduce the governing equations to ordinary differential equations, the numerical method, and the specification of the initial conditions for the Gerz et al. (1989) Ref. [8] case, have been fully discussed in So et al. (2002) Ref. [1]. These details will not be repeated here, the interested readers could consult the reference. Calculations have been carried out for a range of Richardson numbers ranging from $\text{Ri}=\beta g_c\left(\text{d}\Theta /\text{d}z\right)/{\left(\text{d}U/\text{d}z\right)}^2=0.1$ to 1.0.

Figure 2. Comparison of the calculated temperature variance with DNS data. (a) is for the Ri = 0.1 case, and (b) is for the Ri = 1.0 case.

Here, S = dU/dz is the shear gradient, dΘ/dz is the temperature gradient, g_c is the gravitational constant and z is the vertical coordinate aligned with the gravitational direction g₃.

The plots k/k_o, $\overline{{\theta }^2}/\overline{{\theta }_o^2}$, $\overline{u\theta }/u^{\prime }{\theta }^{\prime }$ and $-\overline{w\theta }/w^{\prime }{\theta }^{\prime }$ versus τ = St are shown in Figures 1-4, respectively. The subscript “o” is used here to denote the initial value and τ is the dimensionless time. Here, the prime is used to denote the root mean square value of the fluctuating quantities. Only the results for Ri = 0.1 and

Figure 3. Comparison of the calculated streamwise heat flux with DNS data. (a) is for the Ri = 0.1 case, and (b) is for the Ri = 1.0 case.

Figure 4. Comparison of the calculated normal heat flux with DNS data. (a) is for the Ri = 0.1 case, and (b) is for the Ri = 1.0 case.

1.0 are shown in Figures 1-4. Essentially the same predicted behavior is seen for other values of Ri. The predictions of EASM/BJG and EASM/Γ are essentially identical, thus showing that a one-term explicit representation of the buoyant flux is sufficient to describe the buoyant behavior as that given by a 2-term representation. Therefore, partitioning of $\Gamma$ into f and N is restrictive. It is demonstrated here that such a partitioning is not necessary. As a result, a fairly general form of EASM for turbulent buoyant shear flows has been derived. The resulting EASM is valid for both 2-D and 3-D flows.

6. Conclusion

This paper shows a more general way to derive an EASM for 3-D buoyant shear flows. The approach reveals that it is not necessary to partition the buoyant flux tensor into a 2-D symmetric traceless tensor f and a 3-D tensor N to deduce analytical expressions for the coefficients of the representation for b. The derivation can be carried out directly using $\Gamma$, even though $\Gamma$ is a symmetric traceless 3-D tensor. Thus, the methodology established here can serve as a general approach to derive EASM for a variety of flows where external body forces give rise to a symmetric traceless 3-D tensor in the b equation.

Acknowledgements

The author wishes to acknowledge support given him by the Mechanical Engineering Department, The Hong Kong Polytechnic University, through administrative help, and library and computing services.

Appendix

$\begin{array}{c}{Q}_{12}:=-6{\eta }_8^2{\eta }_2{\eta }_{10}{a}_3{\eta }_1{a}_4-2{\eta }_8^3{\eta }_2{a}_3^2{\eta }_1^2{a}_4^2-6{\eta }_8^3{\eta }_1{a}_2^2{\eta }_2^2{a}_4^2+6{\eta }_8{\eta }_2{\eta }_1^2{a}_3{\eta }_9{a}_4\\ +2{\eta }_1^3{\eta }_8{\eta }_2{a}_3^2{a}_4^2{\eta }_6+6{\eta }_1^3{\eta }_2{\eta }_{14}{a}_2{a}_4{\eta }_6-6{\eta }_1^2{\eta }_2{\eta }_{10}{a}_3{a}_4{\eta }_6+36{\eta }_2{\eta }_{10}^3{a}_3{a}_4\\ +4{\eta }_1^2{\eta }_8{a}_3^2{a}_4^2{\eta }_{14}^2+6{\eta }_1^2{\eta }_8{\eta }_2^2{a}_2^2{a}_4^2{\eta }_6-36{\eta }_8{\eta }_2^2{a}_2^2{a}_4^2{\eta }_{10}^2-12{\eta }_9{\eta }_1^2{\eta }_2{a}_3^2{a}_4^2{\eta }_{10}\\ +12{\eta }_9{\eta }_1^2{\eta }_2{\eta }_{14}{a}_2{a}_4^2{a}_3+12{\eta }_1{\eta }_8{\eta }_2{a}_2^2{a}_4^2{\eta }_{14}^2-6{\eta }_1{\eta }_8^2{\eta }_2{\eta }_{14}{a}_2{a}_4\\ -24{\eta }_1{\eta }_8{\eta }_2{\eta }_{14}{a}_2{a}_4^2{\eta }_{10}{a}_3+12{\eta }_1{\eta }_8{\eta }_2{\eta }_{10}^2{a}_3^2{a}_4^2+12{\eta }_1{\eta }_{14}^3{a}_2{a}_4\\ -12{\eta }_1{\eta }_{14}^2{\eta }_{10}{a}_3{a}_4-36{\eta }_2{\eta }_{14}{a}_2{a}_4{\eta }_{10}^2\end{array}$

$\begin{array}{c}{Q}_{22}:=-4{\eta }_9{\eta }_1^3{\eta }_{14}{a}_3^3{a}_4^2-6{\eta }_1{\eta }_8^2{\eta }_2{a}_2{a}_4{\eta }_{10}{a}_3+36{\eta }_2{a}_2{a}_4{\eta }_{10}^3{a}_3+6{\eta }_9{\eta }_1^2{\eta }_8{\eta }_2{a}_2{a}_4{a}_3\\ -12{\eta }_9{\eta }_1^2{\eta }_2{a}_2{a}_4^2{a}_3^2{\eta }_{10}-12{\eta }_1{\eta }_8{\eta }_{14}{\eta }_{10}{a}_3-12{\eta }_1{a}_2{a}_4{\eta }_{10}{a}_3{\eta }_{14}^2+6{\eta }_1{\eta }_8{a}_2{\eta }_{14}^3\\ +12{\eta }_1{a}_2^2{a}_4{\eta }_{14}^3+3{\eta }_1^2{\eta }_8{\eta }_2{a}_2{\eta }_6+6{\eta }_1^2{\eta }_2{a}_2^2{a}_4{\eta }_{14}{\eta }_6-6{\eta }_1^2{\eta }_2{a}_2{a}_4{\eta }_6{\eta }_{10}{a}_3\\ -36{\eta }_2{a}_2^2{a}_4{\eta }_{14}{\eta }_{10}^2-6{\eta }_1{\eta }_8^2{\eta }_2{a}_2^2{a}_4{\eta }_{14}+8{\eta }_1^2{\eta }_8{\eta }_{14}{a}_3^3{a}_4^2{\eta }_{10}\\ +24{\eta }_1{\eta }_8{\eta }_2{\eta }_{10}^2{a}_3^2{a}_4^2{a}_2-3{\eta }_1{\eta }_2{a}_2{\eta }_8^3-18{\eta }_8{\eta }_2{a}_2{\eta }_{10}^2+6{\eta }_9{\eta }_1^2{\eta }_{14}{a}_3\end{array}$

$\begin{array}{c}{Q}_{32}:=4{\eta }_9{\eta }_1^2{\eta }_2{\eta }_{14}{a}_2{a}_4^2{a}_3+2{\eta }_8{\eta }_2{\eta }_1^2{a}_3{\eta }_9{a}_4+12{\eta }_9{\eta }_1{\eta }_2^2{a}_2^2{a}_4^2{\eta }_{10}-6{\eta }_9{\eta }_1{\eta }_2{\eta }_{10}\\ -2{\eta }_1^2{\eta }_2{\eta }_{10}{a}_3{a}_4{\eta }_6+2{\eta }_1^2{\eta }_2{\eta }_{14}{a}_2{a}_4{\eta }_6+{\eta }_1^2{\eta }_8{\eta }_2{\eta }_6-8{\eta }_1{\eta }_8{\eta }_2{\eta }_{14}{a}_2{a}_4^2{\eta }_{10}{a}_3\\ -2{\eta }_8^2{\eta }_2{\eta }_{10}{a}_3{\eta }_1{a}_4-2{\eta }_1{\eta }_8^2{\eta }_2{\eta }_{14}{a}_2{a}_4-{\eta }_1{\eta }_2{\eta }_8^3-4{\eta }_1{\eta }_{14}^2{\eta }_{10}{a}_3{a}_4+2{\eta }_1{\eta }_8{\eta }_{14}^2\\ +4{\eta }_1{\eta }_{14}^3{a}_2{a}_4-24{\eta }_8{\eta }_2^2{a}_2^2{a}_4^2{\eta }_{10}^2+12{\eta }_2{\eta }_{10}^3{a}_3{a}_4+6{\eta }_8{\eta }_2{\eta }_{10}^2-12{\eta }_2{\eta }_{14}{a}_2{a}_4{\eta }_{10}^2\end{array}$

$Q_{42}:={\eta }_1{\eta }_2{\eta }_8^2-2{\eta }_1{\eta }_{14}^2+6{\eta }_2{\eta }_{10}^2-{\eta }_2{\eta }_1^2{\eta }_6$

Conflicts of Interest

The author declares no conflicts of interest regarding the publication of this paper.

References

[1]	So, R.M.C., Vimala, P., Jin, L.H., Zhao, C.Y. and Gatski, T.B. (2002) Accounting for Buoyancy Effects in the Explicit Algebraic Stress Model: Homogeneous Turbulent Shear Flows. Theoretical and Computational Fluid Dynamics, 15, 283-302. ttps://doi.org/10.1007/s00162-002-0057-x
[2]	Gatski, T.B. and Speziale, C.G. (1993) On Explicit Algebraic Stress Models for Complex Turbulent Flows. Journal of Fluid Mechanics, 254, 59-78. https://doi.org/10.1017/S0022112093002034
[3]	Jongen, T. and Gatski, T.B. (1998) General Explicit Algebraic Stress Relations and Best Approximation for Three-Dimensional Flows. International Journal of Engineering Science, 36, 739-763. https://doi.org/10.1016/S0020-7225(97)00122-5
[4]	Pope, S.B. (1975) A More General Effective Viscosity Hypothesis. Journal of Fluid Mechanics, 72, 331-340. https://doi.org/10.1017/S0022112075003382
[5]	Speziale, C.G., Sarkar, S. and Gatski, T.B. (1991) Modelling the Pressure-Strain Correlation of Turbulence: An Invariant Dynamical Systems Approach. Journal of Fluid Mechanics, 227, 245-272. https://doi.org/10.1017/S0022112091000101
[6]	Rivlin, R.S. and Ericksen, J.L. (1955) Stress-Deformation Relations for Isotropic Materials. Archive of Rational Mechanics and Analysis, 4, 323-425. https://doi.org/10.1512/iumj.1955.4.54011
[7]	Clapham, C. and Nicholson, J. (2009) The Concise Oxford Dictionary of Mathematics. Oxford University Press, Oxford, England. https://doi.org/10.1093/acref/9780199235940.001.0001
[8]	Gerz, T., Schumann, U. and Elghobashi, S.E. (1989) Direct Numerical Simulation of Stratified Homogeneous Turbulent Shear Flow. Journal of Fluid Mechanics, 200, 563-594. https://doi.org/10.1017/S0022112089000765