Steady turbulent flow over a backward-facing step

## 3.1 Steady turbulent flow over a backward-facing step

Tutorial path:

In this example we shall investigate steady turbulent flow over a backward-facing step. The problem description is taken from one used by Pitz and Daily in an experimental investigation [**] against which the computed solution can be compared. This example introduces the following OpenFOAM features for the first time:

• generation of a mesh using blockMesh using full mesh grading capability;

### 3.1.1 Problem specification

The problem is defined as follows:

Solution domain
The domain is 2 dimensional, consisting of a short inlet, a backward-facing step and converging nozzle at outlet as shown in Figure 3.1.

Governing equations

• Mass continuity for incompressible flow  (3.1)

• Steady flow momentum equation  (3.2)

where is kinematic pressure and (in slightly over-simplistic terms) is the viscous stress term with an effective kinematic viscosity , calculated from selected transport and turbulence models.

Initial conditions
, — required in OpenFOAM input files but not necessary for the solution since the problem is steady-state.
Boundary conditions

• Inlet (left) with fixed velocity m/s;
• Outlet (right) with fixed pressure ;
• No-slip walls on other boundaries.
Transport properties

• Kinematic viscosity of air
Turbulence model

• Standard ;
• Coefficients: .
Solver name
simpleFoam: an implementation for steady incompressible flow.
Case name
pitzDaily, located in the \$FOAM_TUTORIALS/incompressible/simpleFoam directory.

The problem is solved using simpleFoam, so-called as it is an implementation for steady flow using the SIMPLE algorithm. The solver has full access to all the turbulence models in the incompressibleTurbulenceModels library and the non-Newtonian models incompressibleTransportModels library of the standard OpenFOAM release.

### 3.1.2 Mesh generation

We expect that the flow in this problem is reasonably complex and an optimum solution will require grading of the mesh. In general, the regions of highest shear are particularly critical, requiring a finer mesh than in the regions of low shear. We can anticipate where high shear will occur by considering what the solution might be in advance of any calculation. At the inlet we have strong uniform flow in the direction and, as it passes over the step, it generates shear on the fluid below, generating a vortex in the bottom half of the domain. The regions of high shear will therefore be close to the centreline of the domain and close to the walls.

The domain is subdivided into 12 blocks as shown in Figure 3.2.

The mesh is 3 dimensional, as always in OpenFOAM, so in Figure 3.2 we are viewing the back plane along . The full set of vertices and blocks are given in the mesh description file below:

1/*--------------------------------*- C++ -*----------------------------------*\
2| =========                 |                                                 |
3| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
4|  \\    /   O peration     | Version:  v1812                                 |
5|   \\  /    A nd           | Web:      www.OpenFOAM.com                      |
6|    \\/     M anipulation  |                                                 |
7\*---------------------------------------------------------------------------*/
8FoamFile
9{
10    version     2.0;
11    format      ascii;
12    class       dictionary;
13    object      blockMeshDict;
14}
15// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
16
17scale   0.001;
18
19vertices
20(
21    (-20.6 0 -0.5)
22    (-20.6 25.4 -0.5)
23    (0 -25.4 -0.5)
24    (0 0 -0.5)
25    (0 25.4 -0.5)
26    (206 -25.4 -0.5)
27    (206 0 -0.5)
28    (206 25.4 -0.5)
29    (290 -16.6 -0.5)
30    (290 0 -0.5)
31    (290 16.6 -0.5)
32
33    (-20.6 0 0.5)
34    (-20.6 25.4 0.5)
35    (0 -25.4 0.5)
36    (0 0 0.5)
37    (0 25.4 0.5)
38    (206 -25.4 0.5)
39    (206 0 0.5)
40    (206 25.4 0.5)
41    (290 -16.6 0.5)
42    (290 0 0.5)
43    (290 16.6 0.5)
44);
45
46negY
47(
48    (2 4 1)
49    (1 3 0.3)
50);
51
52posY
53(
54    (1 4 2)
55    (2 3 4)
56    (2 4 0.25)
57);
58
59posYR
60(
61    (2 1 1)
62    (1 1 0.25)
63);
64
65
66blocks
67(
68    hex (0 3 4 1 11 14 15 12)
69    (18 30 1)
71
72    hex (2 5 6 3 13 16 17 14)
73    (180 27 1)
74    edgeGrading (4 4 4 4 \$negY 1 1 \$negY 1 1 1 1)
75
76    hex (3 6 7 4 14 17 18 15)
77    (180 30 1)
78    edgeGrading (4 4 4 4 \$posY \$posYR \$posYR \$posY 1 1 1 1)
79
80    hex (5 8 9 6 16 19 20 17)
81    (25 27 1)
83
84    hex (6 9 10 7 17 20 21 18)
85    (25 30 1)
87);
88
89edges
90(
91);
92
93boundary
94(
95    inlet
96    {
97        type patch;
98        faces
99        (
100            (0 1 12 11)
101        );
102    }
103    outlet
104    {
105        type patch;
106        faces
107        (
108            (8 9 20 19)
109            (9 10 21 20)
110        );
111    }
112    upperWall
113    {
114        type wall;
115        faces
116        (
117            (1 4 15 12)
118            (4 7 18 15)
119            (7 10 21 18)
120        );
121    }
122    lowerWall
123    {
124        type wall;
125        faces
126        (
127            (0 3 14 11)
128            (3 2 13 14)
129            (2 5 16 13)
130            (5 8 19 16)
131        );
132    }
133    frontAndBack
134    {
135        type empty;
136        faces
137        (
138            (0 3 4 1)
139            (2 5 6 3)
140            (3 6 7 4)
141            (5 8 9 6)
142            (6 9 10 7)
143            (11 14 15 12)
144            (13 16 17 14)
145            (14 17 18 15)
146            (16 19 20 17)
147            (17 20 21 18)
148        );
149    }
150);
151
152// ************************************************************************* //

A major feature of this problem is the use of the full mesh grading capability of blockMesh that is described in section ?? of the User Guide. The user can see that blocks 4,5 and 6 use the full list of 12 expansion ratios. The expansion ratios correspond to each edge of the block, the first 4 to the edges aligned in the local direction, the second 4 to the edges in the local direction and the last 4 to the edges in the local direction. In blocks 4, 5, and 6, the ratios are equal for all edges in the local and directions but not for the edges in the direction that corresponds in all blocks to the global . If we consider the ratios used in relation to the block definition in section ?? of the User Guide, we realize that different gradings have been prescribed along the left and right edges in blocks 4,5 and 6 in Figure 3.2. The purpose of this differential grading is to generate a fine mesh close to the most critical region of flow, the corner of the step, and allow it to expand into the rest of the domain.

The mesh can be generated using blockMesh from the command line and viewed as described in previous examples.

### 3.1.3 Boundary conditions and initial fields

Edit the case files to set the initial and boundary fields for velocity , pressure , turbulent kinetic energy and dissipation rate . The boundary conditions can be specified as: the upper and lower walls are set to Wall, the left patch to Inlet and the right patch to Outlet. These physical boundary conditions require us to specify a fixedValue at the inlet on , and . is given in the problem specification, but the values of and must be chosen by the user in a similar manner to that described in section 2.1.8.1 of the User Guide. We assume that the inlet turbulence is isotropic and estimate the fluctuations to be of at the inlet. We have

 (3.3)

and

 (3.4)

If we estimate the turbulent length scale to be of the width of the inlet then

 (3.5)

At the outlet we need only specify the pressure .

### 3.1.4 Case control

The choices of fvSchemes are as follows: the timeScheme should be steadyState; the gradSchemes and laplacianSchemes should be set as default to Gauss; and, the divSchemes should be set to upwind to ensure boundedness.

Special attention should be paid to the solver settings of the fvSolution dictionary. Although the top level simpleFoam code contains only equations for and , the turbulence model solves equations for , and , and tolerance settings are required for all 5 equations. A tolerance of and relTol of 0.1 are acceptable for all variables with the exception of when and 0.01 are recommended. Under-relaxation of the solution is required since the problem is steady. A relaxationFactor of 0.7 is acceptable for , , and but 0.3 is required for to avoid numerical instability.

Finally, in the controlDict dictionary, the time step deltaT should be set to 1 since in steady state cases such as this is effectively an iteration counter. With benefit of hindsight we know that the solution requires 1000 iterations reach reasonable convergence, hence endTime is set to 1000. Ensure that the writeInterval is sufficiently high, e.g. 50, that you will not fill the hard disk with data during run time.

### 3.1.5 Running the case and post-processing

Run the case and post-process the results. After a few iterations, e.g. 50, a vortex develops beneath the corner of the step that is the height of the step but narrow in the -direction as shown by the vector plot of velocities is shown Figure 3.3(a). Over several iterations the vortex stretches in the -direction from the step to the outlet until at 1000 iterations the system reaches a steady-state in which the vortex is fully developed as shown in Figure 3.3(b-c).