Solution of Incompressible Viscous Flow in a Lid-Driven Cavity using Crank-Nicholson Scheme

Solution of Incompressible Viscous Flow in a Lid-Driven Cavity using Crank-Nicholson Scheme

Solution of Incompressible Viscous Flow in a Lid-Driven Cavity using Crank-Nicholson Scheme

Banamali Dalai

Centre for Advanced Post-Graduate Studies, Mechanical Engineering, Biju Patanaik University of Technology,

Odisha, Rourkela, India,

Abstract-The stream function-vorticity form of the unsteady Navier-Stokes equation is solved using second order accurate in space and first order accurate in time Crank-Nicholson scheme in a finite difference mesh. The scheme shows implicit in nature. The geometry of the problem is a lid-driven cavity. The solution is obtained up to maximum Reynolds number 22,500 in the grid size 129×129 and 257×257. The flow pattern is studied in the form of stream function & vorticity contours, central velocity profiles and location of the eddies.

Keywords – Unsteady, Incompressible, Crank-Nicholson

In this paper, the Navier-Stokes equation is solved in a lid-driven cavity using Crank-Nicholson scheme. The Crank- Nicholson scheme has the flexibility of taking the average value at n+1 and n time steps. This scheme is implicit in nature.

1. FORMULATION

   The governing equation in stream function-vorticity form is expressed as:

   scheme, Finite difference scheme.

   1

   2

   2

   1. INTRODUCTION

      • u v

   t x y

   Re x2

   2

   y

   (1)

   The lid-driven cavity is a square cavity in which all the walls except the top lid are fixed. The top lid is allowed to

   The stream function equation is expressed as:

   2 2

   move towards right with non-dimensional velocity unity. Initially, the cavity is filled with fluid which is at rest. When the lid starts moving towards right the fluid flow in the cavity

   x2

   y 2

   (2)

   is set up. The fluid is assumed to be viscous and incompressible.

   In 1966, Burggraf[1] studied the viscous flow in two- dimensional lid-driven cavity. He solved the stream function and vorticity form of the Navier-Stokes equation up to

   Where , and Re represent the stream function,

   vorticity and Reynolds number respectively. u and v represent the components of velocities along x and y-directions respectively.

   The boundary conditions for the lid-driven cavity are:

   Reynolds number 400 in the grid size 40×40 and 60×60. He has proved analytically that the vorticity value at the centre of the primary vortex will not exceed -1.8859. Later Erturk et

   At y 0, 0 x 1, u v 0; At y 1, 0 x 1, u 1, v 0;

   (3)

   al.[2] solved the stream function-vorticity form of unsteady

   Navier-Stokes equation up to maximum Reynolds number

   At x 0, 1; 0

   y 1, u v 0

   21000 in the grid sizes 401×401, 501×501 and 601×601. Their solution for vorticity value at the centre of the primary

   Eq.(1) & (2) are discretised using Crank-Nicholson

   scheme in finite difference uniform mesh.

   vortex crosses the theoretical limit of -1.8859 within

   n1 n

   un n1 n1

   n n

   Reynolds number 21000. Erturk & Gockol[3] solved the

   ij ij ij ij 1 ij 1 ij 1 ij 1

   stream function-vorticity form of the Navier-Stokes equation using fourth order accurate alternating direction implicit

   t 2

   2x

   2x

   vn n1 n1

   n n

   method upto Reynolds number 20,000 in the grid sizes

   ij i 1 j i 1 j

   i 1 j i 1 j

   2

   401×401, 501×501 and 601×601. Here Erturk and

   Gockols[3] solution for vorticity value at the centre of the

   primary vortex does not cross the theoretical limit of -1.8859. Erturk[4] solved the stream function-vorticity form of the Navier-Stokes equation using Gauss-Seidel iteration techniques in a grid size 1025×1025. He observed that the vorticity value at the centre of the primary eddy does not cross the theoretical limit of -1.8859 within Reynolds number 20000.

   2y

   2y

   n1

   n1 n1 n

   n n

   various Reynolds numbers are shown in Table.1. There are

   1

   ij 1

   2.0 ij

   2

   ij 1

   ij 1 2.0

   ij

   2

   ij 1

   small eddies at the corners along the lower wall of the cavity. As the Reynolds number increases, the stream function

   Re

   n

   2x

   2.0 n1 n1

   2x

   n 2.0 n n

   (4)

   contours in the primary eddy become circular in nature and its centre moves towards the geometric centre of the cavity.

   i 1 j ij i 1 j i 1 j ij i 1 j

   As the Reynolds number increases more and more, the stream

   2y 2

   2y 2

   functions in the primary eddy become more circular in nature

   and the centre of the primary eddy moves towards the

   The discretised equation (3) is first order accurate in time and second order accurate in space.

   Similarly, the stream function equation (2) is discretised using Crank-Nicholson scheme as:

   1 n1 2.0 n1 n1

   geometric centre of the cavity as shown in Fig.1. The number of secondary eddies increase in the lower corner of the cavity. The number of secondary eddies become more than the right corner as the Reynolds number increases. The number of secondary eddies also increase in the left vertical wall of the cavity after Reynolds number 2500 as shown in Fig.1. The

   n1

   ij 1 ij ij 1

   2

   ij

   n 2.0 n n

   x 2

   n1 2.0 n1 n1

   vorticity contours also show large variation near to the walls

   of the cavity. The centre of the primary eddy becomes constant vorticity core with increase of Reynolds number as

   ij 1 ij ij 1 i 1 j ij i 1 j

   shown in Fig.2. The radius of the constant vorticity core

   x 2

   n n n

   y 2

   increases with increase of Reynolds number.

   i1 j 2.0 ij

   y 2

   i1 j

   (5)

   TABLE I: Location and strength of the centre of the primary eddy

   Re	Grid	x	y			Ref.
   100	129×129	0.617	0.734	-0.100	-3.046
   	257×257	0.617	0.738	-0.102	-3.121
   400	129×129	0.563	0.609	-0.107	-2.156
   	257×257	0.555	0.606	-0.111	-2.231
   1000	129×129	0.531	0.563	-0.107	-1.855
   	257×257	0.531	0.566	-0.113	-1.966
   	601×601	0.530	0.565	-0.119	-2.067	Ref[2]
   2500	129×129	0.523	0.539	-0.101	-1.632
   	257×257	0.523	0.543	-0.112	-1.814
   	601×601	0.520	0.543	-0.121	-1.974	Ref[2]
   5000	129×129	0.524	0.531	-0.101	-1.632
   	257×257	0.511	0.535	-0.108	-1.706
   	601×601	0.515	0.536	-0.122	-1.936	Ref[2]
   7500	129×129	0.523	0.531	-0.084	-1.312
   	257×257	0.516	0.531	-0.104	-1.633
   	601×601	0.513	0.532	-0.122	-1.919	Ref[2]
   10000	129×129	0.523	0.531	-0.078	-1.211
   	257×257	0.516	0.527	-0.101	-1.573
   	601×601	0.512	0.530	-0.122	-1.908	Ref[2]
   12500	129×129	0.523	0.523	-0.073	-1.129
   	257×257	0.516	0.527	-0.098	-1.522
   	601×601	0.512	0.528	-0.121	-1.899	Ref[2]
   15000	129×129	0.523	0.523	-0.069	-1.060
   	257×257	0.516	0.527	-0.095	-1.477
   	601×601	0.510	0.528	-0.121	-1.893	Ref[2]
   17500	129×129	0.523	0.523	-0.065	-1.001
   	257×257	0.516	0.527	-0.093	-1.436
   	601×601	0.510	0.527	-0.121	-1.887	Ref[2]
   20000	129×129	0.523	0.523	-0.062	-0.095
   	257×257	0.516	0.523	-0.090	-1.399
   	601×601	0.510	0.527	-0.121	-1.881	Ref[2]
   22500	129×129	0.523	0.523	-0.059	-0.906

   The equation (4) and (5) are second order accurate in space. Equation (4) and (5) produces the vorticity and stream function values respectively in new time steps (n+1). Applying Eq.(3) in Thoms formula[1], the first order accurate boundary values for vorticity can be obtained. The discretised boundary conditions are represented as:

   ij

   2

   p

   ij 1

   ij

   O(h)

   for j 0 and all i;

   ij

   2

   h 2 ij

   ij 1

   O(h)

   for j N and all i;

   ij

   2

   p

   i1 j

   ij

   O(h)

   for i 0 and all j &

   ij

   (6)

   2

   h 2

   ij 1

   ij

   O(h)

   for i N and all j

   Where N is the total number of grid points along i and j directions respectively. Here the boundary conditions are first order accurate in space. Since the discretised equations are implicit in nature, the change of time step does not have any effect to the convergence. The solution of the discretised equation (4) and (5) are obtained using Gauss-Seidel iteration techniques with relaxation parameter. The convergence of the solution is assumed to be the residue of discretised equation

   (1) and (2) approaches towards 10-10. The grid meshes used for the solution are 129×129, 257×257.

2. RESULTS & DISCUSSION

The solution for the incompressible viscous flow in the lid-driven cavity is obtained up to Reynolds number 22,500 in both grid sizes 129×129 and 257×257. It is observed that at lower Reynolds number 100, the stream function contour forms closed structures and its centre is towards right near the lid. This is named as primary eddy of the cavity. The location and strength of the primary eddy for

The u and v-velocity profiles along the centre of the cavity at different Reynolds numbers are shown in Fig.3 and 4 respectively. Both the velocity profiles match very well up to Reynolds number 5000. At Reynolds number 10,000 and 20,000 the computed velocity profiles do not match perfectly with Erturk et al[2]. This happens due to large grid size difference with Erturk et al[2]. The computed velocity profiles are at grid size 257×257 whereas Erturk et als.[2] grid size is 401×401. Here the important point is to be

observed that the vorticity contour plot at Reynolds number 20,000 show zig-zag nature near the top right corner of the cavity in the grid size 257×257 which shows the insufficient grid resolution at that Reynolds number because the u & v- velocity profiles do not match with Erturk et al.[2]. So the Crank-Nicholson scheme seems to produce better result than Erturk et als[2] result but accuracy is the main point of concern.

Figure: 1 Stream function contour plots at various Reynolds numbers

Figure: 2 Vorticity contour plots at various Reynolds numbers

Figure: 3 u-velocity profiles along the centerline of the cavity at different Reynolds numbers

Figure: 4 v-velocity profiles along the centerline of the cavity at different Re

From Fig.5 it is seen that the vorticity value in the grid size 257×257 crosses the theoretical limit value of – 1.8859 near to Re=2000 which is much earlier than Erturk et al[2]. So, it seen that though the high Reynolds number

solution is obtained at lower grid size but the accuracy is very poor. From Fig.6 comparing the magnitude of the velocity values along the centre of the cavity, it seen that as the Reynolds number increases the length of the central portion of linear profile increases. The linear nature of the velocity profile shows the constant vorticity region which increases with increase of Reynolds number. The velocity profiles become closer to the wall with increase of Reynolds number which indicates the presence of boundary layers on the walls. This type of flow feature can also be verified from the vorticity contour plot in Fig.2.

Figure: 5 Vorticity value at the centre of the primary eddy with Re

Figure: 6 u and v-velocity profiles along the centreline of the cavity with Re

IV CONCLUSION

The stream function vorticity form of the Navier- Stokes equation is solved using Crank-Nicholson scheme. The solution is obtained up to maximum Reynolds number 22,500 using the grid size 129×129 and 257×257. The Crank- Nicholson scheme produces better result at lower grid size than that of Erturk et al but the accuracy of the former scheme is very poor.

REFERENCES

[1] R Burggraf, Analytical and numerical studies of the structure of steady separated flows, Journal of Fluid Mechanics 24 pp.113-151 (1966).

[2] E. Erturk, T. C. Corke and C Gockol, Numerical solutions of 2-D steady incompressible driven cavity flow at highe Reynolds numbers, Int. Journal of Numerical Methods in Fluids 48 pp. 747- 774 (2005).

[3] E Erturk and C. Gokcol, Fourth-order compact formulation of Navier-Stokes equations and driven cavity flow at high Reynolds numbers, Int. J. for Num. Methods in Fluids. 50 pp. 421-436 (2006).

[4] E Erturk, Discussions on driven cavity flow. Int. J. for Num.

Methods in Fluids. 60 pp. 275-294 (2009).