* Lab 4 - Driven
Cavity Problem*

**Objectives**

This lab explores solutions to the 2-D, incompressible Navier-Stokes CFD equations. Specifically, this lab seeks to solve well-posed, isothermal, laminar, incompressible N-S equations in the streamfunction-vorticity form. The accomplishment of this objective occurs via the following proceedure:

- Execute an accuracy/monotonicity study for the driven cavity benchmmark on16x16 and 32x32 uniform meshes at Re = 100.
- Examine energy norms and integer print fields for solution quality, hence design an improved non-uniform 32x32 mesh based on the data.
- Proceed to Re = 1000, 2000 and repeat the experiments.
- Compare solution energy norms as a function of solution adapted meshings, hence validate the conjecture about existence of a "best" 32x32 meshing.
- Reexecute one of the poorer Re = 1000 or 2000 cases adding TWSh dissipation (0 <= beta <= 1.0). Compare "best" monotone GWSh and TWSh solution energy norms and their distributions.

**Problem Statement **

For the case of 2-D, isothermal, incompressible, laminar, unsteady flow, we can enforce the continuity equation via introduction of the vorticity and the streamfunction. In so doing, the velocities and pressure gradient "disappear" and the system reduces to two equations on the vorticity and the streamfunction.

Following the usual procedure, the unsteady time term is reexpressed in terms of a Taylor series expansion with a Euler-implicitness term.

The vorticity equation may now be expressed in template pseudo-code.

{FOMG} = ()()()(0;1)(B200)(-OMGA) +(REI)()()(1122;-1)(B211)(OMGA) +(REI)()()(3344;-1)(B222)(OMGA) +(REI)()()(1324;-1)(B221)(OMGA) +(REI)()()(1324;-1)(B212)(OMGA) +( )()(PSI)(23;-1)(B3102)(OMGA) +(-)()(PSI)(14;-1)(B3102)(OMGA) +( )()(PSI)(14;-1)(B3201)(OMGA) +(-)()(PSI)(23;-1)(B3201)(OMGA)

The streamfunction equation is similarly derived in the template pseudo-cod where the wall velocity condition is expressed everywhere as zero.

{FPSI} = ()()()(1122;-1)(B211)(PSI) +()()()(3344;-1)(B222)(PSI) +()()()(1324;-1)(B221)(PSI) +()()()(1324;-1)(B212)(PSI) +(-)()()(0;1)(B200)(OMGA)

Again as expected, the non-linear nature of the equations requires solution via an iterative method. We choose the time honored Newton iterative method which then requires the calculation of the Jacobian. The procedure is as follows:

{J OMG/OMG} = ()()()(0;1)(B200)() +(REI)()()(1122;-1)(B211)() +(REI)()()(3344;-1)(B222)() +(REI)()()(1324;-1)(B221)() +(REI)()()(1324;-1)(B212)() +( )()(PSI)(23;-1)(B3102)() +(-)()(PSI)(14;-1)(B3102)() +( )()(PSI)(14;-1)(B3201)() +(-)()(PSI)(23;-1)(B3201)() {J OMG/PSI} = ( )()(OMGA)(23;-1)(B3201)() +(-)()(OMGA)(14;-1)(B3201)() +( )()(OMGA)(14;-1)(B3102)() +(-)()(OMGA)(23;-1)(B3102)() {J PSI/OMG} = (-)()()(0;1)(B200)() {J PSI/PSI} = ()()()(1122;-1)(B211)() +()()()(3344;-1)(B222)() +()()()(1324;-1)(B221)() +()()()(1324;-1)(B212)()

Once a solution is obtained, we must perform a post-processing operation to revert the solution back into the two velocity components and the pressure gradient. The vorticity equation gives a solution for the two velocity components as follows:

{U1} = (-)()()(2;0)(B201)(OMGA) +(-)()()(4;0)(B202)(OMGA) {U2} = ()()()(1;0)(B201)(OMGA) +()()()(3;0)(B202)(OMGA)

The streamfunction equation also gives a solution for the two velocity components:

{UPSI} = (-)()()(2;0)(B201)(PSI) +(-)()()(4;0)(B202)(PSI) {VPSI} = ()()()(1;0)(B201)(PSI) +()()()(3;0)(B202)(PSI)

Finally, we may solve for the pressure gradient in the solution as follows:

{PRES} = ()(U1)()(11;-1)(B211)(U1) +()(U2)()(12;-1)(B211)(U1) +()(U1)()(21;-1)(B211)(U2) +()(U2)()(22;-1)(B211)(U2) +()(U1)()(13;-1)(B212)(U1) +()(U2)()(14;-1)(B212)(U1) +()(U1)()(23;-1)(B212)(U2) +()(U2)()(24;-1)(B212)(U2) +()(U1)()(31;-1)(B221)(U1) +()(U2)()(32;-1)(B221)(U1) +()(U1)()(41;-1)(B221)(U2) +()(U2)()(42;-1)(B221)(U2) +()(U1)()(33;-1)(B222)(U1) +()(U2)()(34;-1)(B222)(U1) +()(U1)()(43;-1)(B222)(U2) +()(U2)()(44;-1)(B222)(U2)

As one final note, it is expected that the solution will experience instabilities due to dispersion error. A modified Taylor Weak Statement (TWS^{h}) can provide added diffusion to help stabalize the system. For a steady-state solution, the following modifications are made to the vorticity equation:

There are three different options to handle this extra diffusion term. The first option simply keeps a time step in the system. This is undesireable since we seek a steady-state solution.

By defining a local time scale via the element length and velocity magnitude, the system can be reformed without a time dependence:

The final option involves performing the Taylor series expansion on x rather than t for the steady-state solution. This leads to:

The second method is chosen for the template pseudo-code as detailed below:

{FOMG} = {FOMG} +(BETA)(UM11)()(11;-0.5)(B211)(OMGA) +(BETA)(UM12)()(12;-0.5)(B211)(OMGA) +(BETA)(UM12)()(21;-0.5)(B211)(OMGA) +(BETA)(UM22)()(22;-0.5)(B211)(OMGA) +(BETA)(UM11)()(13;-0.5)(B212)(OMGA) +(BETA)(UM12)()(14;-0.5)(B212)(OMGA) +(BETA)(UM12)()(32;-0.5)(B212)(OMGA) +(BETA)(UM22)()(24;-0.5)(B212)(OMGA) +(BETA)(UM11)()(31;-0.5)(B221)(OMGA) +(BETA)(UM12)()(32;-0.5)(B221)(OMGA) +(BETA)(UM12)()(41;-0.5)(B221)(OMGA) +(BETA)(UM22)()(42;-0.5)(B221)(OMGA) +(BETA)(UM11)()(33;-0.5)(B222)(OMGA) +(BETA)(UM12)()(34;-0.5)(B222)(OMGA) +(BETA)(UM12)()(43;-0.5)(B222)(OMGA) +(BETA)(UM22)()(44;-0.5)(B222)(OMGA)

**Discussion of Results **

Figures 1-3 show the non-uniform vorticity solution for three Reynolds numbers. It it seen, as expected, that dispersion error becomes a problem as the Reynolds number increases. This is manifested as oscillatory waves around the singularities where the moving top meets the stationary walls.

Figure 1

Figure 2

Figure 3

One way to eliminate the dispersion error is via mesh refinement. Since there are two singularities which exhibit oscillations, it was decided to impliment a symmetrical progression ratio non-uniform mesh as shown in Figure 4. Table 1 details the results of a small optimization study for the symmetric progression ratio. The resulting data does not make much since. The highlighted values select the progression ratio combination that results in the lowest energy norm errors. These values are then used to resolve the system and mitigate for dispersion error as shown in Figures 5-7.

Figure 4

Table 1: Progression ratio optimization.

Figure 5

Figure 6

Figure 7

Here is a representative speed and pressure solution for the case of Re = 1000. The velocity and pressure gradients are simply solved via a post-processing operation.

Figure 8

Figure 9

Finally, the TWS is implimented in an attempt to mitigate dispersion error via the introduction of artificial diffusion. A trial and error approach led to an optimal value for beta of 0.2 as shown in Figure 10 below. When implimented alongside the symmetric progression ratio non-uniform meshing, a superior solution is obtained as seen in Figure 11. Upon inspection, it was seen that the energy norms in both the vorticity and streamfunction solutions were lower in all locations for the TWS versus the like GWS solution.

Figure 10

Figure 11

**Conclusion **

In conclusion, this lab has explored aspects of the streamline-vorticity Navier-Stokes GWS algorithm. The dispersion error that inevitably creeps in at higher Reynolds numbers has been deal with via a combination of non-uniform meshing as well as the addition of artificial diffusion. These combination studies have then led to an optimal TWS solution that, aside from the corner-singularities, accurately expresses the given situation of a steady-state driven cavity.