Open access peer-reviewed chapter

Two Different Formulations for Solving the Navier-Stokes Equations with Moderate and High Reynolds Numbers

Written By

Blanca Bermúdez, Alejandro Rangel-Huerta, Wuiyevaldo Fermín Guerrero-Sánchez and José David Alanís

Submitted: 22 April 2017 Reviewed: 24 October 2017 Published: 20 December 2017

DOI: 10.5772/intechopen.71921

From the Edited Volume

Computational Fluid Dynamics - Basic Instruments and Applications in Science

Edited by Adela Ionescu

Chapter metrics overview

1,002 Chapter Downloads

View Full Metrics

Abstract

In this work, we discuss the numerical solution of the Taylor vortex and the lid-driven cavity problems. Both problems are solved using the Stream function-vorticity formulation of the Navier-Stokes equations in 2D. Results are obtained using a fixed point iterative method and working with matrixes A and B resulting from the discretization of the Laplacian and the advective term, respectively. We solved both problems with Reynolds numbers in the range of 3200 ≤ Re ≤ 7500. Results are also obtained using the velocity-vorticity formulation of the Navier-Stokes equations. In this case, we are using only the fixed point iterative method. We present results for the lid-driven cavity problem and for the Stream function-vorticity formulation with Reynolds numbers in the range of 3200 ≤ Re ≤ 7500. As the Reynolds number increases, the time and the space step size have to be refined. We show results for 3200 ≤ Re ≤ 20,000. The numerical scheme with the velocity-vorticity formulation uses a smaller step size for both time and space. Results are not as good as with the Stream function-vorticity formulation, although the way the scheme behaves gives us another point of view on the behavior of fluids under different numerical schemes and different formulation.

Keywords

  • Navier-Stokes equations
  • velocity-vorticity formulation
  • Stream function-vorticity formulation
  • Reynolds number
  • fixed point iterative process

1. Introduction

The fixed point iterative method has already been used for solving the Navier-Stokes equations and the Boussinesq system under different formulations, see [1, 2, 3, 4].

The idea behind this iterative method was to work with a symmetric and positive definite matrix. The scheme worked very well, as shown in [1, 2, 3, 4], but the processing time was, in general, very large especially for high Reynolds numbers. Working with matrixes A and B, we are dealing with a non-symetric matrix, but it can be proved that it is strictly diagonally dominant for Δt sufficiently small. The processing time with the second method was reduced in approximately 30 or 35%.

Additionally, in order to show that the fixed point iterative method works well for moderate and high Reynolds numbers, we report results for the lid-driven cavity problem and Reynolds numbers in the range of 3200 ≤ Re ≤ 100,000 using the Stream function-vorticity formulation and also the velocity-vorticity formulation, but in the case of the velocity-vorticity formulation, we just arrive to Re = 20,000, because of computing time and memory requirements.

Results, in both formulations, are obtained using the fixed point iterative method reported in [5], which is applied to a non-linear elliptic system resulting after time discretization. The method has shown to be robust enough to handle moderate and high Reynolds numbers, which is not an easy task, see [1, 2].

As the Reynolds number increases, the mesh has to be refined and a smaller time step has to be used, in order to capture the fast dynamics of the flow and, numerically, because of stability reasons, as mentioned in [6], although, with the velocity-vorticity formulation [7, 8], a finer mesh has to be used, both in time and in space.

The computing time is, in general, very large with this numerical scheme and for both formulations, so that is why we are looking forward to reduce computing time working with both matrixes A and B resulting from the discretization of the Laplacian and the advective term, respectively, instead of working just with matrix A, which is symmetric and positive definite.

With the Stream function-vorticity formulation, and for moderate and high Reynolds numbers, the second scheme was faster than the fixed point iterative method (see [9, 10]). With respect to the velocity-vorticity formulation, we are just showing results using the fixed point iterative method, and for lower Reynolds numbers, but we are looking forward to modify the scheme also in order to reduce computing time.

Advertisement

2. Mathematical model

Let Ω R N N = 2 3 be the region of a viscous, incompressible, non-stationary flow and Γ is its boundary

u t 1 Re Δ u + p + u u = f , a u = 0 b E1

These are the Navier-Stokes equations in primitive variables. This system must be provided with appropriate initial conditions u x 0 = u 0 x in Ω and boundary conditions u = g on Γ .

When working in a two-dimensional region Ω , taking the curl in both sides of (Eq. (1a)) and taking into account that

u 1 = ψ y , u 2 = ψ x , E2

followed by (Eq. (1b)), with ψ the Stream function, u 1 and u 2 the two components of the velocity, we arrive to the Stream function-vorticity formulation of the Navier-Stokes equations.

The following system of equations is then obtained:

Δ ψ = ω a ω t v Δ ω + u ω = f ω , b E3

where ω is the vorticity given by ω = u 2 x u 1 y , and v = 1 Re .

This system represents the Navier-Stokes equations in the Stream function-vorticity formulation. The incompressibility condition (Eq. (1b)), because of (Eq. (2)), is automatically satisfied and the pressure does not appear any more.

Advertisement

3. The velocity-vorticity formulation

Taking the curl in

ω = × u , E4

and using the identity × × a = 2 a + a and (Eq. (1b)), we obtain the following Poisson type equation for the velocity:

Δ u = × ω E5

Two Poisson type equations for the velocity are obtained, which together with the equation for the vorticity give us:

Δ u 1 = ω y a Δ u 2 = ω x b ω t v Δ ω + u ω = f ω c E6

These are the Navier-Stokes equations in the velocity-vorticity formulation.

Advertisement

4. Numerical method for the Stream function-vorticity formulation

The following second-order approximation for the time derivative is used:

ω t x n + 1 t 3 ω n + 1 4 ω n + ω n 1 2 Δt , E7

where n 1 , x Ω and t > 0 is the time step.

The resulting discretization system reads:

Δ ψ n + 1 = ω , ψ Γ = ψ bc a α ω n + 1 v Δ ω n + 1 + u n + 1 ω n + 1 = f ω , ω Γ = ω bc , b E8

where α = 3 2 t and f ω = 4 ω n ω n 1 2 t .

At each time level, the following non-linear system has to be solved.

Δ ψ = ω , ψ Γ = ψ bc a αω v Δ ω + u ω = f ω , ω Γ = ω bc , b E9

To obtain ( ψ 1 , ω 1 ) in (Eq. (8)), any second-order strategy using a combination of one step can be applied and steady systems of the form (Eq. (9)) are also obtained.

For solving this system of equations, two strategies were used in this work: First, we used the fixed point iterative method described in [5]:

Denoting R ω ω ψ by

R ω ω ψ = αω v Δ ω + u ω f ω in Ω , E10

our system is equivalent to

Δ ψ = ω , ψ Γ = ψ bc a R ω ψ = 0 , ω Γ = ω bc b E11

Then, at each time level, the following the fixed point iterative process (see [5]) is used:

Given ω n , 0 = ω n and ψ n , 0 = ψ n solve until convergence in ω and ψ

Δ ψ n , m + 1 = ω n , m in Ω , ψ n , m + 1 = ψ bc on Γ ω n , m + 1 = ω n , m ρ ω αI v Δ 1 R ω ω n , m ψ n , m + 1 in   Ω , ω n , m + 1 = ω bc n , m + 1 on Γ , ρ ω > 0 ; E12

and then take ω n + 1 ψ n + 1 = ω n , m + 1 ψ n , m + 1 .

To reduce the computing time, we worked in solving the system by the following method at each time step:

Δ ψ n + 1 = ω n , ψ n + 1 Γ = ψ bc a αI v h 2 A ω n + 1 + 1 2 h B ω n + 1 = f ω , ω n + 1 Γ = ω bc , b E13

where A and B are the matrixes resulting from the discretization of the Laplacian and the advective term respectively, and (Eq. (13b)) is solved using Gauss-Seidel method.

Advertisement

5. Numerical method for the velocity-vorticity formulation

The second-order approximation (Eq. (7)) for the time derivative is used and the following non-linear system is obtained in Ω

Δ u 1 = ω y Δ u 2 = ω x , u Γ = u bc R ω ω u = 0 , ω Γ = ω bc , E14

where

R ω ω u αω v Δ ω + u ω f ω . E15

Using again the fixed point iterative method previously described, we have:

Given ω n , 0 = ω n , u 1 n , 0 = u 1 n , u 2 n , 0 = u 2 n solve until convergence in ω , u 1 and u 2

Δ u 1 n , m + 1 = ω n , m y Δ u 2 n , m + 1 = ω n , m x , u n , m + 1 Γ = u b , c n , m + 1 αI v Δ ω n , m + 1 = αI v ω m ρ ω > 0 , ω n , m + 1 Γ = ω bc n , m E16

and then take ω n + 1 u 1 n + 1 u 2 n + 1 = ω n , m + 1 u 1 n , m + 1 u 2 n , m + 1 .

Advertisement

6. Numerical experiments

With respect to the lid-driven cavity problem and using the Stream function-vorticity formulation Ω = 0 1 × 0 1 , the top wall is moving with a velocity given by (1, 0) and for the other walls, the velocity is given by (0, 0). Ψ is over-determined at the boundary ( ψ n Γ is also known) and there is no boundary condition for ω. In our case, we have followed the alternative proposed by Goyon [11]. Ψ = 0 is chosen over Г. A translation of the boundary condition in terms of the velocity (primitive variable) has to be used. By Taylor series expansion of (Eq. (3a)), we obtained:

ω 0 y t = 1 2 h x 2 8 ψ h x y t ψ 2 h x y t ω a y t = 1 2 h x 2 8 ψ a h x y t ψ a 2 h x y t ω x 0 t = 1 2 h y 2 8 ψ x h y t ψ x 2 h y t ω x b t = 1 2 h y 2 8 ψ x b h y t ψ x b 2 h x t E17

where h x , h y denote the spatial step size in the directions of x and y, respectively.

In Figures 1 and 2 , we show results for the lid-driven cavity problem with Re = 5000 and Re = 7500 , with h x = h y = 1 / 64 .

Figure 1.

Streamlines (left) and isovorticity contours (right) for Re = 5000, hx  = hy  = 1/64.

Figure 2.

Streamlines (left) and isovorticity contours (right) for Re = 7500, hx  = hy  = 1/64.

For the Taylor vortex problem, results are shown in Figures 3 and 4 for Re = 5000 and Re = 7500 , with h x = h y = 2 π / 64 and t = 10 .

Figure 3.

Stream function and vorticity for Re = 5000 hx  = hy  = 2π/64 and t = 10.

Figure 4.

Exact stream function and vorticity for Re = 5000 hx  = hy  = 2π/64 y t = 10.

The exact Stream function and the vorticity are also shown in Figure 5 , for Re = 5000 . For this problem, Ω = 0 2 π × 0 2 π the exact solution is known and is given by:

u 1 x y t = cos x sin y e 2 t Re u 2 x y t = sin x cos y e 2 π Re E18

Figure 5.

Stream function and vorticity for Re = 7500 hx  = hy  = 2π/64 y t = 10.

In the primitive variables formulation, we have as initial conditions:

u 1 x y t = cos x sin y u 2 x y t = sin x cos y E19

For the Stream function, the boundary conditions are:

ψ x 0 t = ψ x π t = cos x e 2 t Re ψ 2 π y t = cos y e 2 t Re E20

For the vorticity, the boundary conditions are:

ω x 0 t = ω x 2 π t = 2 cos x e 2 t Re ω 0 y t = ω 2 π y t = 2 cos y e 2 t Re E21

In Tables 1 and 2 , we show computing times for the above-mentioned problems with both the methods; the Fixed Point Iterative Method and working with matrixes A and B.

Re Fixed point iterative method (s) Working with A and B (s)
5000 153 120
7500 801 610.25

Table 1.

Time in seconds, for both Reynolds numbers and the two methods described for the lid-driven cavity problem.

Re Fixed point iterative method (s) Working with A and B (s)
5000 15.5 12.75
7500 15.5 12.75

Table 2.

Time in seconds, for both Reynolds numbers and the two methods for the Taylor vortex problem.

In Figure 6 , we show the streamlines and isovorticity contours for Re = 25,000 , with h = h y = 1 / 728 . In Figure 7 , we show results for Re = 50,000 , with h = h x = h y = 1 / 1024 . For these values of the Reynolds number, since there is no steady state, we show results for T final = 5 .

Figure 6.

Streamlines (left) and isovorticity contours (right) for Re = 25,000, h = hx  = hy = 1/728 y dt = 0.00025, Tfinal  = 5.

Figure 7.

Streamlines (left) and isovorticity contours (right) for Re = 50,000, h = hx  = hy  = 1/1024 y dt = 0.00025, Tfinal  = 5.

Then, in Figure 8 , we show just the isovorticity contours for Re = 75,000 , with h = h x = h y = 1 / 1024 and for Re = 100,000 , with h = h x = h y = 1 1280 and T final = 5 .

Figure 8.

Isovorticity contours (left) for Re = 75,000, h = hx  = hy  = 1/1024 y dt = 0.00025, Tfinal  = 5, isovorticity contours (right) for Re = 100,000, h = hx  = hy = 1/1280 y dt = 0.00025, Tfinal  = 5.

In the case of the velocity-vorticity formulation and the lid-driven cavity problem, the boundary condition for u is given by u = (1, 0) in the moving boundary y = b and u = (0, 0) anywhere else at the boundary.

Not all the results were obtained using the second-order discretization. In some cases, a fourth-order discretization has to be used, using the fourth-order option of Fishpack [12] (used in this work for solving the elliptic problems appearing).

In Figure 9 , we show the streamlines and the isovorticity contours for Re = 3200 , h = h x = h y = 1 / 512 , T final = 50 .

Figure 9.

Streamlines (left) and isovorticity contours (right) for Re = 3200, h = hx  = hy  = 1/512 y dt = 0.0001, Tfinal  = 50.

In Figure 10 , we show the isovorticity contours for Re = 20,000 with (a) h = h x = h y = 1 / 1512 , T final = 5 , obtained using the velocity-vorticity formulation and (b) with the Stream function-vorticity formulation with h = h x = h y = 1 / 768 , T final = 5 .

Figure 10.

Isocontours for the vorticity for Re = 20,000, (a) velocity-vorticity formulation with h = hx  = hy = 1/1512, dt = 0.0001, Tfinal  = 5, (b) Stream function-vorticity formulation with h = hx  = hy = 1/768, dt = 0.0001, Tfinal  = 5.

As can be noticed with the Stream function-vorticity formulation, we are using a value of h half of the size of the one used with the velocity-vorticity formulation. We assume the results obtained with the first-mentioned formulation are more reliable. Computing time for the velocity-vorticity formulation was much larger. We think there are still some numerical problems with this formulation and for very high Reynolds numbers.

Advertisement

7. Conclusions

For the lid-driven cavity problem results agree very well with those reported in the literature [1, 2, 3, 4, 13, 14], and by working with matrixes A and B it was possible to reduce computing time between a 30 and 35%.

As can be seen in Figures 1 and 2 , numerical oscillations occurred, given the high Reynolds numbers used in such a way that it is necessary to use smaller values of h [6], numerically because of stability of the method and physically in order to capture the fast dynamics of the flow.

For high Reynolds numbers and small values of h the computational work takes a lot of time, so reducing computing time becomes a very important fact. For the Taylor Vortex Problem [8, 15], processing time was also reduced between 30 and 35%.

With the velocity-vorticity formulation, as already mentioned, we only show results using the Fixed Point Iterative Method, and we are looking forward working with both matrixes A and B, in order to reduce computing time also with this formulation. This is the reason why we only show results till Re = 20,000 and not for higher Reynolds numbers.

In conclusion, the numerical scheme applied with the stream function-vorticity formulation is not as good with the velocity-vorticity formulation, although, the way it behaves with some values of the parameters and the order of discretization, gives us another point of view about the behavior of the fluids under different numerical methods and different formulations.

We must also say that our code has not been parallelized since it is difficult to do this. It must be taken into account that the equations, in both formulations, are coupled. We are looking forward to use a solver for the system of linear equations that can be parallelized. This can be viewed as a future work.

Advertisement

Acknowledgments

The authors would like to acknowledge the support received by the National Laboratory of Supercomputing from of the southeast of Mexico BUAP-INAOE-UDLAP (Laboratorio Nacional de Supercómputo del Sureste de México).

References

  1. 1. Bermúdez B, Nicolás A, Sánchez FJ, Buendía E. Operator splitting and upwinding for the Navier–Stokes equations. Computational Mechanics. 1997;20(5):474-477
  2. 2. Nicolás A, Bermúdez B. 2D incompressible viscous flows at moderate and high Reynolds numbers. CMES. 2004;6(5):441-451
  3. 3. Nicolás A, Bermúdez B. 2D thermal/isothermal incompressible viscous flows. International Journal for Numerical Methods in Fluids. 2005;48:349-366
  4. 4. Bermúdez B, Nicolás A. Isothermal/thermal incompressible viscous fluid flows with the velocity–vorticity formulation. Información Tecnológica. 2010;21(3):39-49
  5. 5. Nicolás A. A finite element approach to the Kuramoto Sivashinski equation. In: Advances in Numerical Equations and Optimization. SIAM; 1991
  6. 6. Nicolás Carrizosa A, Bermúdez Juárez B. Onset of two dimensional turbulence with high Reynolds numbers in the Navier‐Stokes equations. Computational Methods for Coupled Problems in Science and Engineering IV M. Papadrakakis, E. Oñate and B. Schrefler (Eds.) First edition, June 2011:996-1006
  7. 7. Nicolás A, Bermúdez B. Viscous incompressible flows by the velocity–vorticity Navier–Stokes equation. CMES. 2007;20(2):73-83
  8. 8. Bermúdez B, Nicolás A. The Taylor vortex and the driven cavity problems by the velocity–vorticity formulation. In: Proceedings of the 7th International Conference on Heat Transfer, Fluid Mechanics and Thermodynamics; 2010. pp. 629-632
  9. 9. Bermúdez B, Juárez L. Numerical solution of an advection-diffusion equation. Información Tecnológica. 2014;25(1):151-160
  10. 10. Bermúdez B, Posadas R. The Taylor vortex and the driven cavity problems in the stream functions-vorticity formulation. In: Proceedings of the 11th World Congress on Computational Mechanics (WCCM XI); 5th European Congress on Computational Mechanics (ECCM V); 6th European Congress on Computational Fluid Dynamics (ECDF VI). A Publication of CIMNE; 2014. pp. 3746-3755
  11. 11. Goyon O. High Reynolds number solutions of Navier–Stokes equations using incremental unknowns. Computer Methods in Applied Mechanics and Engineering. 1996;130:319-335
  12. 12. Adams J, Swarztrauber P, Sweet R. FISHPACK: A Package of FORTRAN Subprograms for the Solutions Elliptic PDE's. Boulder, Colorado, USA: The National Center for Atmospheric Research; 1980
  13. 13. Glowinski R. Finite element methods for the numerical simulations of incompressible viscous flow. In: Introduction to the Control of the Navier‐Stokes Equations, Handbook of Numerical Analysis, Vol. II, North-Holland, Amsterdam, 1991, 17-351
  14. 14. Ghia U, Guia KN, Shin CT. High-re solutions for incompressible flow using the Navier–Stokes equations and a multigrid method. Journal of Computational Physics. 1982;48:387-411
  15. 15. Anson DK, Mullin T, Cliffe KA. A numerical and experimental investigation of a new solution in the Taylor vortex problem. Journal of Fluid Mechanics. 1988:475-487

Written By

Blanca Bermúdez, Alejandro Rangel-Huerta, Wuiyevaldo Fermín Guerrero-Sánchez and José David Alanís

Submitted: 22 April 2017 Reviewed: 24 October 2017 Published: 20 December 2017