Week 8: Literature review RANS derivation and analysis
OBJECTIVE:
To apply Reynolds decomposition to the NS equations and come up with the expression for Reynolds stress.
- Explain your understanding of the terms Reynolds stress
- What is turbulent viscosity? How is it different from molecular viscosity?
THEORY:
In fluid dynamics, a turbulent regime refers to irregular flows in which eddies, swirls, and flow instabilities occur. It is governed by high momentum convection and low momentum diffusion. It is in contrast to the laminar regime, which occurs when a fluid flows in parallel layers with no disruption between the layers. The turbulence regime is extremely frequent in natural phenomena and human applications; some examples are the rise of cigarette smoke, waterfalls, blood flow in arteries, and most of the terrestrial atmospheric recirculation. In terms of human applications, the turbulent regime occurs in the aerodynamics of vehicles, but also in many industrial applications such as heat exchangers, quenching processes, or continuous casting of steel.
Reynolds showed the transition between a laminar and a turbulent regime through a set of experimental investigations. He also suggested that this transition was directly linked to the ratio between inertial and viscous forces. This ratio was computed by George Gabriel Stokes in 1851 and has been named “Reynolds number” in honor of Osborne Reynolds who popularized it. This dimensionless number is defined as:
where:
- ρ is the density of the fluid
- u is the macroscopic velocity of the flow
- d is the characteristic length (or hydraulic diameter)
- μ is the dynamic viscosity of the fluid
- ν is the kinematic viscosity of the fluid.
Turbulent flows occur when Re exceeds a certain threshold (dependent on the application’s topology and flow physics) called “Critical Reynolds number”.
This principle was motivated by energetical considerations; big eddies are highly inertial and tend to be unstable. Their motion feeds smaller eddies thanks to a local transfer of kinetic energy. These smaller eddies undergo the same process, giving rise to even smaller eddies that inherit the energy of their parent eddy, and so on. This transfer of energy is usually called “energy cascade” and it is mainly inertial, thus almost no energy dissipation occurs until reaching a sufficiently small length scale such that the viscosity of the fluid can effectively dissipate the kinetic energy. This latter scale of the turbulence exhibits a local laminar regime and is characterized by a low value of Re. This process has been depicted in Figure 1. Richardson’s studies highlight an important feature of turbulent flows: they are energy demanding. A turbulent flow will dissipate energy and decay to a laminar flow unless it is fed by an external source of energy.
RANS:
It is normally believed (but not proven) that Navier-Stokes (NS) equations can model any kind of flow, turbulent flows included. The problem is that for very high values of Re, the resolution of the NS equation is very challenging and not stable. Thus a small perturbation in the parameter like the initial or boundary conditions may lead to a completely different solution. This problem is partially overcome by the use of the Reynolds-Averaged Navier-Stokes equations (RANS)3. To understand RANS we need to look into time scales employed in fluid mechanics.
The main idea of RANS is that the fluctuating quantity is not of importance, the transport equations are modeled on the basis of mean quantities of flow parameters, with the time scale adjusted correspondingly.
Reynolds Decomposition:
Following the idea of Reynolds (1895), we assume that the fluid is in a random unsteady turbulent state and work with the time-averaged or mean equations of motion. Any variable 'Q' is resolved into a mean value of `bar plus a fluctuating value Q', where by definition,
`barQ = 1/T int_(t_0)^(t_0+T)Q dt`
Instantaneous quantity at any point in space and time = Time-averaged Quantity + Fluctuating component.
where T is large compared to the relevant period of fluctuations. The mean value `bar itself may vary slowly with time as shown in the image below.
Let us consider only incompressible turbulent flow with constant transport properties but with possible significant fluctuations in velocity, pressure, and temperature.
`u(x,y,t) = baru(x,y)+u'(x,y,t)`
`v(x,y,t) = barv(x,y)+v'(x,y,t)`
`P(x,y,t) = barP(x,y)+P'(x,y,t)`
`T(x,y,t) = barT(x,y)+T'(x,y,t)`
`ϕ(x,y,t) = bar ϕ(x,y)+ϕ'(x,y,t)`
Time averaged quantities:
`barQ(x,y) = 1/T int_(t_0)^(t_0+T)Q(x,y,t) dt`
mean values are not the function of 't' .
Moreover, the averaging of fluctuating terms lead to zero as shown below.
`1/Tint_0^Tu(x,y,t) dt= 1/Tint_0^Tbaru(x,y)dt+1/Tint_0^Tu'(x,y,t)dt`
and
`1/Tint_0^Tu(x,y,t) dt= baru(x,y)`
`1/Tint_0^Tbaru(x,y)dt = baru(x,y)`
implies
`int_0^Tu'(x,y,t)dt = 0`
Before substituting these functions in the NS equations, we can verify from basic integral relation that the following rules of averaging apply for any two turbulent quantities 'f' and 'g'.
`bar(f') = 0`
`bar(barf)=barf`
`bar(fbarg)=barfbarg`
`bar(f'barg) = 0`
`bar(f+g) = barf+barg`
`bar(∂f)/(∂s) = (∂barf)/(∂s)`
`bar(fg)=barfbarg+bar(f'g')`
`bar(intfds) = intbarf ds`
Substituting Averged quantites in NS-equations:
Contnuity Equation:
`(∂u)/(∂x)+(∂v)/(∂y)=0`
substituting the averaged quantities
`(∂(baru+u'))/(∂x)+(∂(barv+v'))/(∂y)=0`
expanding the terms
`(∂baru)/(∂x)+(∂u')/(∂x)+(∂barv)/(∂y)+(∂v')/(∂y)=0`
averaging over time step 'T'
`1/Tint_0^T(∂baru)/(∂x)dt+1/Tint_0^T(∂u')/(∂x)dt+1/Tint_0^T(∂barv)/(∂y)dt+1/Tint_0^T(∂v')/(∂y)dt=0`
since the time integral of fluctuation component is zero and integral limit not a function of space.
`1/Tint_0^T(∂baru)/(∂x)dt+color(red)cancel(1/Tint_0^T(∂u')/(∂x)dt)+1/Tint_0^T(∂barv)/(∂y)dt+color(red)cancel(1/Tint_0^T(∂v')/(∂y)dt)=0`
The equation reduces to
`1/Tint_0^T(∂baru)/(∂x)dt+1/Tint_0^T(∂barv)/(∂y)dt=0`
since the `baru` and `barv` is the function of space only, the above equation can be re-written as:
`(∂baru)/(∂x)1/Tint_0^Tdt+(∂barv)/(∂y)1/Tint_0^Tdt=0`
and
`1/Tint_0^Tdt=(T-0)/T = 1`
implies, the continuity equation after the substitution of averaged quantities depends on only mean values of decomposed quantities as shown below.
`(∂baru)/(∂x)+(∂barv)/(∂y)=0`
x-Momentum Equation:
`ρ[(∂u)/(∂t)+u(∂u)/(∂x)+v(∂u)/(∂y)]=-(∂P)/(∂x)+ μ ((∂^2u)/(∂x^2)+(∂^2u)/(∂y^2))`
multiplying both sides by `1/ρ` and replacing `μ/ρ` with kinematic viscosity `ν`
`(∂u)/(∂t)+u(∂u)/(∂x)+v(∂u)/(∂y)=-1/rho(∂P)/(∂x)+ ν (∂^2u)/(∂x^2)+ ν (∂^2u)/(∂y^2)`
to convert the above equation in suitable form, we multiply x-component of velocity with continuity equation and add it to left hand side and its resultant is zero.
`(∂u)/(∂t)+u(∂u)/(∂x)+v(∂u)/(∂y)+u((∂u)/(∂x)+(∂v)/(∂y))=-1/rho(∂P)/(∂x)+ ν (∂^2u)/(∂x^2)+ ν (∂^2u)/(∂y^2)`
rearranging the terms in the above equation
`(∂u)/(∂t)+2u(∂u)/(∂x)+v(∂u)/(∂y)+u(∂v)/(∂y)=-1/rho(∂P)/(∂x)+ ν (∂^2u)/(∂x^2)+ ν (∂^2u)/(∂y^2)`
From Product rule,
`(∂(uv))/(∂y) =u(∂v)/(∂y)+v(∂u)/(∂y)`
and
`2u(∂u)/(∂x) = (∂u^2)/(∂x)`
substituting the terms of product rule in the reduced equation.
`(∂u)/(∂t)+(∂u^2)/(∂x)+(∂(uv))/(∂y)=-1/rho(∂P)/(∂x)+ ν (∂^2u)/(∂x^2)+ ν (∂^2u)/(∂y^2)`
substituting the averaged quantities
`(∂(baru+u'))/(∂t)+(∂(baru+u')^2)/(∂x)+(∂((baru+u')(barv+v')))/(∂y)=-1/rho(∂(barP+P'))/(∂x)+ ν (∂^2(baru+u'))/(∂x^2)+ ν (∂^2(baru+u'))/(∂y^2)`
expanding the terms
`(∂baru)/(∂t)+(∂u')/(∂t)+(∂baru^2)/(∂x)+(∂(u')^2)/(∂x)+(∂(2baru u'))/(∂x)+(∂(barubarv))/(∂y)+(∂(baru v'))/(∂y)+(∂(u'barv))/(∂y)+(∂(u'v'))/(∂y)=-1/rho(∂barP)/(∂x)-1/rho(∂P')/(∂x)+ ν (∂^2baru)/(∂x^2)+ ν (∂^2u')/(∂x^2)+ ν (∂^2baru)/(∂y^2)+ ν (∂^2u')/(∂y^2)`
averaging over time step 'T'
`1/Tint_0^T(∂baru)/(∂t)dt+1/Tint_0^T(∂u')/(∂t)dt+1/Tint_0^T(∂baru^2)/(∂x)dt+1/Tint_0^T(∂(u')^2)/(∂x)dt+1/Tint_0^T(∂(2baru u'))/(∂x)dt+1/Tint_0^T(∂(barubarv))/(∂y)dt+1/Tint_0^T(∂(baru v'))/(∂y)dt+1/Tint_0^T(∂(u'barv))/(∂y)dt+1/Tint_0^T(∂(u'v'))/(∂y)dt=-1/Tint_0^T1/rho(∂barP)/(∂x)dt-1/Tint_0^T1/rho(∂P')/(∂x)dt+ 1/Tint_0^Tν (∂^2baru)/(∂x^2)dt+ 1/Tint_0^Tν (∂^2u')/(∂x^2)dt+ 1/Tint_0^Tν (∂^2baru)/(∂y^2)dt+ 1/Tint_0^Tν (∂^2u')/(∂y^2)dt`
since the time integral of fluctuation component is zero and the rules mentioned above.
`1/Tint_0^T(∂baru)/(∂t)dt+color(red)cancel(1/Tint_0^T(∂u')/(∂t)dt)+1/Tint_0^T(∂baru^2)/(∂x)dt+1/Tint_0^T(∂(u')^2)/(∂x)dt+color(red)cancel(1/Tint_0^T(∂(2baru u'))/(∂x)dt)+1/Tint_0^T(∂(barubarv))/(∂y)dt+color(red)cancel(1/Tint_0^T(∂(baru v'))/(∂y)dt)+color(red)cancel(1/Tint_0^T(∂(u'barv))/(∂y)dt)+1/Tint_0^T(∂(u'v'))/(∂y)dt=-1/Tint_0^T1/rho(∂barP)/(∂x)dt-color(red)cancel(1/Tint_0^T1/rho(∂P')/(∂x)dt)+ 1/Tint_0^Tν (∂^2baru)/(∂x^2)dt+ color(red)cancel(1/Tint_0^Tν (∂^2u')/(∂x^2)dt)+ 1/Tint_0^Tν (∂^2baru)/(∂y^2)dt+ color(red)cancel(1/Tint_0^Tν (∂^2u')/(∂y^2)dt)`
The above equation reduces to
`1/Tint_0^T(∂baru)/(∂t)dt+1/Tint_0^T(∂baru^2)/(∂x)dt+1/Tint_0^T(∂(u')^2)/(∂x)dt+1/Tint_0^T(∂(barubarv))/(∂y)dt+1/Tint_0^T(∂(u'v'))/(∂y)dt=-1/Tint_0^T1/rho(∂barP)/(∂x)dt+ 1/Tint_0^Tν (∂^2baru)/(∂x^2)dt+ 1/Tint_0^Tν (∂^2baru)/(∂y^2)dt`
Upon integrating, the equation reduces to
`(∂baru)/(∂t)+(∂baru^2)/(∂x)+1/Tint_0^T(∂(u')^2)/(∂x)dt+(∂(barubarv))/(∂y)+1/Tint_0^T(∂(u'v'))/(∂y)dt=-1/rho(∂barP)/(∂x)+ ν (∂^2baru)/(∂x^2)+ ν (∂^2baru)/(∂y^2)`
Rearranging the terms
`(∂baru)/(∂t)+(∂baru^2)/(∂x)+(∂(barubarv))/(∂y)=-1/rho(∂barP)/(∂x)+ ν (∂^2baru)/(∂x^2)+ ν (∂^2baru)/(∂y^2)-color(red)(1/Tint_0^T(∂(u'v'))/(∂y)dt)-color(red)(1/Tint_0^T(∂(u')^2)/(∂x)dt)`
From the above equation it is clear that the momentum equation obtained with mean values of components is similar to momentum equation of actual transport properties except with the red terms.
If we consider that the `u'` value is fluctuating between -1 and 1 then by averaging it the value of `u'` comes equal to 0 but if we consider the same case for `(u')^2` then the square of -1 and 1 will be between 0 and 1 which when integrated will not give zero. Hence,
`1/Tint_0^T(∂(u')^2)/(∂x)dt !=0`
Moreover, the above term can be approximated to zero in the boundary layer where the gradient along y is higher in comparison to the x-direction. Though we are not using this assumption in this case.
Similarly in the term, `1/Tint_0^T(∂(u'v'))/(∂y)dt` the fluctuations along the x-axis `u'` and the fluctuations along the y-axis `v'` are correlated. That is when a low momentum y-component is flowing in the positive y-direction then the high momentum (i.e., high velocity) x-component gets slowed down which is caused by `v'`. So, the positive increase in the `v'` is going to cause a negative result in the `u'`, and the product of two will come out to be negative value only. Hence, this term will also not give zero.
`ubrace(color(blue)((∂baru)/(∂t)))_("change w.r.t time")+ubrace(color(green)((∂baru^2)/(∂x)+(∂(barubarv))/(∂y)))_("convective terms")=ubrace(color(violet)(-1/rho(∂barP)/(∂x)))_("Pressure Gradient")+ ubrace(ν ((∂^2baru)/(∂x^2)+(∂^2baru)/(∂y^2)))_("Shear Stress")-ubrace(color(red)(1/Tint_0^T(∂(u'v'))/(∂y)dt)-color(red)(1/Tint_0^T(∂(u')^2)/(∂x)dt))_("Reynolds Stress")`
In the above equation, integral values are replaced by averages.
`(∂baru)/(∂t)+(∂baru^2)/(∂x)+(∂(barubarv))/(∂y)=-1/rho(∂barP)/(∂x)+ ν (∂^2baru)/(∂x^2)+ ν (∂^2baru)/(∂y^2)-(∂(bar(u'v')))/(∂y)-(∂(bar((u')^2)))/(∂x)`
The equation can be further rearranged as follows:
`(∂baru)/(∂t)+(∂baru^2)/(∂x)+(∂(barubarv))/(∂y)=-1/rho(∂barP)/(∂x)+ ν (∂^2baru)/(∂x^2)+ ν (∂^2baru)/(∂y^2)+1/ρ[(∂(-ρbar(u'v')))/(∂y)+(∂(-ρbar((u')^2)))/(∂x)]`
In the compact notation:
`"Stress term" = color(blue)(ubrace(μ∇^2u_i)_("shear stress gradient"))-color(red)(ubrace(ρ(∂bar(u_i^')bar(u_j^'))/(∂x_j))_("Reynolds Stress gradient")`
`= μ∂/(∂x_j)((∂u_i)/(∂x_j))-ρ∂/(∂ x_j)(bar(u_i^')bar(u_j^'))`
`= ∂/(∂x_j)[μ((∂u_i)/(∂x_j))-ρ(bar(u_i^')bar(u_j^'))]`
`"Total Shear Stress" = ubrace(μ((∂u_i)/(∂x_j)))_("Laminar Stress")-ubrace(ρ(bar(u_i^')bar(u_j^')))_("Turbulent Stress")`
In the 3D cartesian cordinate system the above equation will be as follows:
`(∂baru)/(∂t)+(∂baru^2)/(∂x)+(∂(barubarv))/(∂y)+(∂(barubarw))/(∂z)=-1/rho(∂barP)/(∂x)+ ν [(∂^2baru)/(∂x^2)+(∂^2baru)/(∂y^2)+(∂^2baru)/(∂z^2)]+1/ρ[(∂(-ρbar((u')^2)))/(∂x)+(∂(-ρbar(u'v')))/(∂y)+(∂(-ρbar(u'w')))/(∂z)]`
Hence the resulting Reynolds stress in the matrix form, along the x, y, and z axes is as follows:
`τ_(ij)^(Reynolds) = [[τ_(x x), τ_(xy), τ_(xz) ],[τ_(yx), τ_(yy), τ_(yz) ],[τ_(zx), τ_(zy), τ_(zz) ]]=[[-ρbar((u')^2),-ρbar((u'v')),-ρbar((u'w'))], [-ρbar((v'u')),-ρbar((v')^2),-ρbar((v'w'))], [-ρbar((w'u')),-ρbar((w'v')),-ρbar((w')^2)]]`
Reynold's stress: In fluid dynamics, the Reynolds stress is the component of the total stress tensor in a fluid obtained from the averaging operation over the Navier-Stokes equations to account for turbulent fluctuations in fluid momentum This stress offer significantly effects for the complex interactions in turbulent flows. This stress can be used to model the Reynolds Stress Transport Equation which will allow using the comparatively larger grid to flow analysis.
MOLECULAR & TURBULENT VISCOCITY