Example

2D cavity flow

Our first example of simulating Navier-Stokes flow is flow within a 2D cavity. Flow is driven by a kinematic boundary condition at one side, while all other sides are walls with zero velocity. The model setup is shown in Fig. 3 .

../../_images/2d_cavity_setup.png

Fig. 3 Two-dimensional cavity model geometry and boundary conditions.

Preparing the case

Two-dimensional cavity flow is included in the official tutorials of OpenFoam. We get started by copying the tutorial case into our work directory. You need to do this from your docker shell.

cd $HOME/HydrothermalFoam_runs
cp $FOAM_TUTORIALS/incompressible/icoFoam/cavity/cavity ./cavity2D

Check the directory structure:

.
├── 0
│   ├── U
│   └── p
├── a.foam
├── clean.sh
├── constant
│   └── transportProperties
└── system
    ├── blockMeshDict
    ├── controlDict
    ├── fvSchemes
    └── fvSolution

The 0 directory now contains the boundary and initial conditions for pressure and velocity (our primary variables), the constant directory has the transport properties (here viscosity), and the system directory contains information on the mesh, time steps, numerical schemes etc.

The two files run.sh and clean.sh are actually not included and we need to create them. These files are usually part of every OpenFoam case; run.sh includes all commands that need to executed to run the case and clean.sh cleans the case (removes mesh and output directories).

Listing 6 The run.sh file.
 1#!/bin/sh
 2cd ${0%/*} || exit 1    # Run from this directory
 3
 4# Source tutorial run functions
 5. $WM_PROJECT_DIR/bin/tools/RunFunctions
 6
 7application=`getApplication`
 8
 9./clean.sh
10runApplication blockMesh
11runApplication $application
Listing 7 The clean.sh file.
1#!/bin/sh
2cd ${0%/*} || exit 1 # run from this directory
3
4# Source tutorial run functions
5. $WM_PROJECT_DIR/bin/tools/CleanFunctions

Make the scripts executable.

chmod u+x clean.sh run.sh

Making the mesh

We will use OpenFoam’s blockMesh utility to make a simple 2D mesh. The corresponding blockMeshDict file that has all the meshing information is located in the system folder.

../../_images/cavity2d_bm.png

Fig. 4 Structure of the blockMeshDict

First we need to define the vertices of the mesh, the nodes.

../../_images/cavity2d_vertices.png

Fig. 5 Numbering of the vertices.

The next step is define the connectivity between the vertices in order to describe the modeling domain.

../../_images/cavity2d_vorder.png

Fig. 6 The order by which the vertices are passed to the hex command matters!

Order of vertices

The OpenFoam documentation provides a nice description of the vertices ordering.

  • the axis origin is the first entry in the block definition, vertex 0 in our example

  • the x direction is described by moving from vertex 0 to vertex 1

  • the y direction is described by moving from vertex 1 to vertex 2

  • vertices 0, 1, 2, 3 define the plane z = 0

  • vertex 4 is found by moving from vertex 0 in the z direction

  • vertices 5,6 and 7 are similarly found by moving in the z direction from vertices 1,2 and 3 respectively.

Next boundary patches are defined and labeled in the blockMeshDict. Also here care must be take to provide the vertices in a consistent order (right-hand coordinate system). Two easy ways to remember this is to:

  • apply the right-hand rule, which means if the thumb of your right hand points to the outside of a face, the numbering has to follow your fingers.

  • or, looking onto a face and starting from any vertex, the numbering has to be counter-clockwise.

../../_images/cavity2d_bounds.png

Fig. 7 Assigning boundary labels and types.

Now we are ready to run the blockMesh utility and create the mesh

blockMesh

You can visualize the mesh using paraview

touch a.foam
paraview a.foam

Boundary conditions

We now have velocity and pressure as primary variables and need to set initial and boundary conditions for them. First we look at the velocity boundary conditions:

code 0/u
../../_images/cavity2d_u.png

Fig. 8 Velocity boundary conditions. The front and back sides are set to empty because we are doing a 2D calculation.

Next we look into the pressure boundary conditions.

code 0/p
Listing 8 Pressure boundary conditions. Front and back are of type “emtpy” for 2-D runs.
 1/*--------------------------------*- C++ -*----------------------------------*\
 2=========                 |
 3\\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
 4 \\    /   O peration     | Website:  https://openfoam.org
 5  \\  /    A nd           | Version:  7
 6   \\/     M anipulation  |
 7\*---------------------------------------------------------------------------*/
 8FoamFile
 9{
10    version     2.0;
11    format      ascii;
12    class       volScalarField;
13    object      p;
14}
15// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
16
17dimensions      [0 2 -2 0 0 0 0];
18
19internalField   uniform 0;
20
21boundaryField
22{
23    movingWall
24    {
25        type            zeroGradient;
26    }
27
28    fixedWalls
29    {
30        type            zeroGradient;
31    }
32
33    frontAndBack
34    {
35        type            empty;
36    }
37}
38
39// ************************************************************************* //

Tip

One has to be careful about the dimensions of pressure in OpenFoam. In incompressible runs, like we are doing here, the pressure is usually the relative pressure \(\frac{p}{\rho}\) and has units \(\frac{m^2}{s^2}\)

Run controls

The time stepping, run time, and output frequency are again set in system/controlDict. Open it and check that you understand the entires.

In case you wondered how OpenFoam is solving the equations. We will cover the details later in the course, but you can have a preview by opening the system/fvSchemes file. In this dictionary, the various discretization schemes can be set. Fig. 9 gives some further explanations.

../../_images/cavity2_num.png

Fig. 9 The exact discretization schemes can be set in system/fvSchemes.

Time to run the case! Just start the solver

icoFoam

Technical detail

The icoFoam solver is actually using the PISO algorithm and the not the SIMPLE algorithm described above. SIMPLE is a steady-state solver, while PISO can also resolve transient flows. PISO also make some smart corrections for better convergence. We will learn about these technical issues later in this course.

Visualization

Open paraview and look at the results.