Challenge

Week 8 - Simulation of a backward facing step in OpenFOAM

AAmal T15 Oct 2022

We can notice that in our case, the flow is laminar and incompressible and the solution is transient as well.

I will be solving for 2 cases, one without Grading Factor and the other with a Grading Factor.

Grading Factor and its importance in Meshing

When solving a problem numerically, it is essential to do meshing. When less elements are generated in the mesh, we trade off some accuracy of the true solution for faster runtimes and more stability. A Graded mesh essentially reduces the size of the elements in the mesh in one part of the domain while increasing the size in another part. If done correctly, a graded mesh can improve the accuracy of a solution. A graded mesh allows us to increase the accuracy of the solution in the part of the domain we are interested in (near the step walls, in this case), while sacrificing some accuracy in other parts we are not interested in. The advantage is that, since the number of elements remains constant, the stability (courant number) of the solution is not affected.

Case 1: No Grading Factor

First I start off by navigating to the tutorials folder in OpenFOAM, then I find the icoFoam solver under incompressible flow. I decide to use the cavity solver, and I proceed to copy the cavity solver from the tutorials folder to the week8 folder under the name case_1.

We can notice that there are 3 folders inside. The 0 folder is for the initial conditions, it contains the U(velocity) folder and P(Pressure) folder. The system folder contains BlockMeshDict, fvSchemes, fvSolutions and controlDict. The Constant folder contains transportProperties.

These are the files that I changed from their default contents which came from the tutorial.

System Folder:

blockMeshDict - This file contains the geometry of the domain, the faces, edges, etc. It also contains the mesh properties. The mesh used in this case is a uniformly sized mesh containing 200 elements in the longest X direction and 20 in the longest Y direction.

/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2012 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

convertToMeters 1;

vertices
(
// vertices from point 0 to 10
(0 0 0)
(0.08 0 0)
(0.08 0.005 0)
(0 0.005 0)
(0.08 0.01 0)
(0 0.01 0)
(0.2 0.01 0)
(0.2 0.005 0)
(0.2 0 0)
(0.2 -0.01 0)
(0.08 -0.01 0)

// vertices from point 11 to 21
(0 0 0.1)
(0.08 0 0.1)
(0.08 0.005 0.1)
(0 0.005 0.1)
(0.08 0.01 0.1)
(0 0.01 0.1)
(0.2 0.01 0.1)
(0.2 0.005 0.1)
(0.2 0 0.1)
(0.2 -0.01 0.1)
(0.08 -0.01 0.1)
);

blocks
(
// block 1
hex (0 1 2 3 11 12 13 14) (80 5 1) simpleGrading (1 1 1)
//block 2
hex (3 2 4 5 14 13 15 16) (80 5 1) simpleGrading (1 1 1)
//block 3
hex (2 7 6 4 13 18 17 15) (120 5 1) simpleGrading (1 1 1)
//block 4
hex (1 8 7 2 12 19 18 13) (120 5 1) simpleGrading (1 1 1)
//block 5
hex (10 9 8 1 21 20 19 12) (120 10 1) simpleGrading (1 1 1)
);

edges
(
);

boundary
(
inlet
{
type patch;
faces
(
(11 14 3 0)
(14 16 5 3)
);
}
outlet
{
type patch;
faces
(
(6 17 18 7)
(8 7 18 19)
(9 8 19 20)
);
}
otherwall //no slip walls
{
type wall;
faces
(
//top
(15 4 5 16)
(15 17 6 4)
//step
(21 12 1 10)
//bottom left
(11 0 1 12)
//bottom right
(21 10 9 20)
);
}
frontAndBack
{
type empty;
faces
(
//back
(0 3 2 1)
(3 5 4 2)
(2 4 6 7)
(1 2 7 8)
(10 1 8 9)
//front
(13 15 16 14)
(12 13 14 11)
(20 19 12 21)
(19 18 13 12)
(18 17 15 13)
);
}

);
mergePatchPairs
(
);

// ************************************************************************* //

controlDict - This folder contains some essential values which controls the simulation (like timeperiod, timestep, precision, etc.) In this case, the simulation was run for 5 seconds (5 cycles in total) and with a timestep of 0.0005s.

/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2012 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "system";
object controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

application icoFoam;

startFrom startTime;

startTime 0;

stopAt endTime;

endTime 5;

deltaT 0.0005;

writeControl runTime;

writeInterval 0.5;

purgeWrite 0;

writeFormat ascii;

writePrecision 6;

writeCompression off;

timeFormat general;

timePrecision 6;

runTimeModifiable true;

//adaptive timestep
/*
adjustTimeStep true;

maxCo 1;*/

// ************************************************************************* //

fvSchemes:

/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2012 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "system";
object fvSchemes;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

ddtSchemes
{
default Euler;
}

gradSchemes
{
default Gauss linear;
grad(p) Gauss linear;
}

divSchemes
{
default none;
div(phi,U) Gauss linear;
}

laplacianSchemes
{
default Gauss linear orthogonal;
}

interpolationSchemes
{
default linear;
}

snGradSchemes
{
default orthogonal;
}

// ************************************************************************* //

fvSolution:

/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2012 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "system";
object fvSolution;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

solvers
{
p
{
solver PCG;
preconditioner DIC;
tolerance 1e-06;
relTol 0.05;
}

pFinal
{
$p;
relTol 0;
}

U
{
solver smoothSolver;
smoother symGaussSeidel;
tolerance 1e-05;
relTol 0;
}
}

PISO
{
nCorrectors 2;
nNonOrthogonalCorrectors 0;
pRefCell 0;
pRefValue 0;
}

// ************************************************************************* //

folder fvSolutions and fvSchemes were not changed for this simulation.

0 folder

U (Velocity) - This contains the initial and boundary velocity conditions for the domain. The initial conditions are 0.1 m/s across the domain. The boundary condition at the inlet is 0.2 m/s.

/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2012 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volVectorField;
object U;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions [0 1 -1 0 0 0 0];

internalField uniform (0.1 0 0);

boundaryField
{
inlet
{
type fixedValue;
value uniform (0.2 0 0);
}

outlet
{
type zeroGradient;
}

frontAndBack
{
type empty;
}
otherwall
{
type noSlip; // zero at walls
}
}

// ************************************************************************* //

P(pressure):

This contains the initial and boundary pressure conditions for the domain. The initial conditions are 0 Pa across the domain. The boundary condition at the outlet is 0 Pa.

/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2012 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
object p;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions [0 2 -2 0 0 0 0];

internalField uniform 0;

boundaryField
{
inlet
{
type zeroGradient;
}

outlet
{
type fixedValue;
value uniform 0;
}

frontAndBack
{
type empty;
}

otherwall
{
type zeroGradient;
}
}

// ************************************************************************* //

constant Folder:

Upon running the blockMesh command, that script also generates a polyMesh folder in this folder, which contains details of the mesh.

transportProperties - This contains other values required to run the simulation but remain constant through the solution (such as Dynamic Viscosity ν">νν)

Dynamic viscosity was chosen to be 1e-4.

/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2012 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "constant";
object transportProperties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

nu [0 2 -1 0 0 0 0] 1e-4;

// ************************************************************************* //

Now, I begin the simulation process.

1. Running the blockMesh command

blockMesh:
/*---------------------------------------------------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2012 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
Build : _7bdb509494-20201222 OPENFOAM=2012
Arch : "LSB;label=32;scalar=64"
Exec : blockMesh
Date : Oct 14 2022
Time : 22:29:45
Host : MY-PC
PID : 6288
I/O : uncollated
Case : /home/amal/OpenFOAM/amal-v2012/run/week8/case_1/cavity/cavity
nProcs : 1
trapFpe: Floating point exception trapping enabled (FOAM_SIGFPE).
fileModificationChecking : Monitoring run-time modified files using timeStampMaster (fileModificationSkew 5, maxFileModificationPolls 20)
allowSystemOperations : Allowing user-supplied system call operations

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time

Creating block mesh from "system/blockMeshDict"
--> FOAM IOWarning :
Found [v1012] 'convertToMeters' entry instead of 'scale' in dictionary "/home/amal/OpenFOAM/amal-v2012/run/week8/case_1/cavity/cavity/system/blockMeshDict"

This keyword is deemed to be 120 months old.

Creating block edges
No non-planar block faces defined
Creating topology blocks
Creating topology patches

Creating block mesh topology

Check topology

Basic statistics
Number of internal faces : 5
Number of boundary faces : 20
Number of defined boundary faces : 20
Number of undefined boundary faces : 0
Checking patch -> block consistency

Creating block offsets
Creating merge list (topological search)...
Deleting polyMesh directory "constant/polyMesh"

Creating polyMesh from blockMesh
Creating patches
Creating cells
Creating points with scale 1
Block 0 cell size :
i : 0.001 .. 0.001
j : 0.001 .. 0.001
k : 0.1 .. 0.1

Block 1 cell size :
i : 0.001 .. 0.001
j : 0.001 .. 0.001
k : 0.1 .. 0.1

Block 2 cell size :
i : 0.001 .. 0.001
j : 0.001 .. 0.001
k : 0.1 .. 0.1

Block 3 cell size :
i : 0.001 .. 0.001
j : 0.001 .. 0.001
k : 0.1 .. 0.1

Block 4 cell size :
i : 0.001 .. 0.001
j : 0.001 .. 0.001
k : 0.1 .. 0.1

There are no merge patch pairs

Writing polyMesh with 0 cellZones
----------------
Mesh Information
----------------
boundingBox: (0 -0.01 0) (0.2 0.01 0.1)
nPoints: 6842
nCells: 3200
nFaces: 13020
nInternalFaces: 6180
----------------
Patches
----------------
patch 0 (start: 6180 size: 10) name: inlet
patch 1 (start: 6190 size: 20) name: outlet
patch 2 (start: 6210 size: 410) name: otherwall
patch 3 (start: 6620 size: 6400) name: frontAndBack

End

The following commands are executed

icoFoam
paraFoam

This means the solution is now ready, paraFoam opens the paraView window which contains our solution.

Observations:

We can notice that the flow through the inlet section is laminar until the step wall is encountered. The velocity of the flow in the inlet region is high at the centre and progressively decreases going towards the walls. When the flow encounters the step wall for the first time, there are minor vortices formed near the step wall and at the top of the wall which dissipates as the flow approaches steady state. Slowly the flow approaches steady state, at around 125th timestep (2.5s) into the simulation. There is no significant change in flow after it attains steady state.

Case 2 (Grading factor of 0.2) :

These are the files that I changed from their default contents which came from the tutorial.

System Folder:

blockMeshDict - This file contains the geometry of the domain, the faces, edges, etc. It also contains the mesh properties. The mesh used in this case is a graded mesh containing 200 elements in the longest X direction and 20 in the longest Y direction. The grading factor is 0.2 and in such a way that elements near the step wall is smaller than elements on the other side of the domain.

convertToMeters 1;

blocks
(
// block 1
hex (0 1 2 3 11 12 13 14) (80 5 1) simpleGrading (0.2 5 1)
//block 2
hex (3 2 4 5 14 13 15 16) (80 5 1) simpleGrading (0.2 0.2 1)
//block 3
hex (2 7 6 4 13 18 17 15) (120 5 1) simpleGrading (5 0.2 1)
//block 4
hex (1 8 7 2 12 19 18 13) (120 5 1) simpleGrading (5 5 1)
//block 5
hex (10 9 8 1 21 20 19 12) (120 10 1) simpleGrading (5 1 1)
);

edges
(
);

mergePatchPairs
(
);

// ************************************************************************* //

controlDict:

This folder contains some essential values which controls the simulation (like timeperiod, timestep, precision, etc.) In this case, the simulation was run for 5 seconds (5 cycles in total) and with a timestep of 0.0005s.

/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2012 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "system";
object controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

application icoFoam;

startFrom startTime;

startTime 0;

stopAt endTime;

endTime 5;

deltaT 0.0005;

writeControl timeStep;

writeInterval 0.5;

purgeWrite 0;

writeFormat ascii;

writePrecision 6;

writeCompression off;

timeFormat general;

timePrecision 6;

runTimeModifiable true;

// ************************************************************************* //

folder fvSolutions and fvSchemes were not changed for this simulation.

0 folder

U (Velocity):

This contains the initial and boundary velocity conditions for the domain. The initial conditions are 0.1 m/s across the domain. The boundary condition at the inlet is 0.2 m/s.

/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2012 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volVectorField;
object U;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions [0 1 -1 0 0 0 0];

internalField uniform (0.1 0 0);

boundaryField
{
inlet
{
type fixedValue;
value uniform (0.2 0 0);
}

outlet
{
type zeroGradient;
}

frontAndBack
{
type empty;
}
otherwall
{
type noSlip; // zero at walls
}
}

// ************************************************************************* //

p (pressure):

This contains the initial and boundary pressure conditions for the domain. The initial conditions are 0 Pa across the domain. The boundary condition at the outlet is 0 Pa.

/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2012 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
object p;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions [0 2 -2 0 0 0 0];

internalField uniform 0;

boundaryField
{
inlet
{
type zeroGradient;
}

oUtlet
{
type fixedValue;
value uniform 0;
}

frontAndBack
{
type empty;
}

otherwall
{
type zeroGradient;
}
}

// ************************************************************************* //

constant Folder

Upon running the blockMesh command, that script also generates a polyMesh folder in this folder, which contains details of the mesh.

transportProperties:

This contains other values required to run the simulation but remain constant through the solution (such as Dynamic Viscosity `nu`)

Dynamic viscosity was chosen to be 1e-5.

/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2012 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "constant";
object transportProperties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

nu [0 2 -1 0 0 0 0] 1e-5;

// ************************************************************************* //

Now, I begin the simulation process.

blockMesh

/*---------------------------------------------------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2012 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
Build : _7bdb509494-20201222 OPENFOAM=2012
Arch : "LSB;label=32;scalar=64"
Exec : blockMesh
Date : Oct 15 2022
Time : 19:57:26
Host : MY-PC
PID : 1590
I/O : uncollated
Case : /home/amal/OpenFOAM/amal-v2012/run/week8/case_2/cavity/cavity
nProcs : 1
trapFpe: Floating point exception trapping enabled (FOAM_SIGFPE).
fileModificationChecking : Monitoring run-time modified files using timeStampMaster (fileModificationSkew 5, maxFileModificationPolls 20)
allowSystemOperations : Allowing user-supplied system call operations

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time

This keyword is deemed to be 120 months old.

Creating block edges
No non-planar block faces defined
Creating topology blocks
Creating topology patches

Creating block mesh topology

Check topology

Basic statistics
Number of internal faces : 5
Number of boundary faces : 20
Number of defined boundary faces : 20
Number of undefined boundary faces : 0
Checking patch -> block consistency

Creating block offsets
Creating merge list (topological search)...

Creating polyMesh from blockMesh
Creating patches
Creating cells
Creating points with scale 1
Block 0 cell size :
i : 0.00200653 .. 0.000401307
j : 0.000382406 .. 0.00191203
k : 0.1 .. 0.1

Block 1 cell size :
i : 0.00200653 .. 0.000401307
j : 0.00191203 .. 0.000382406
k : 0.1 .. 0.1

Block 2 cell size :
i : 0.00040166 .. 0.0020083
j : 0.00191203 .. 0.000382406
k : 0.1 .. 0.1

Block 3 cell size :
i : 0.00040166 .. 0.0020083
j : 0.000382406 .. 0.00191203
k : 0.1 .. 0.1

Block 4 cell size :
i : 0.00040166 .. 0.0020083
j : 0.001 .. 0.001
k : 0.1 .. 0.1

There are no merge patch pairs

End

The following commands are executed

icoFoam

paraFoam

This means the solution is now ready, paraFoam opens the paraView window which contains our solution.

Observations: We can notice that the flow through the inlet section is laminar until the step wall is encountered. The velocity of the flow in the inlet region is high at the centre and progressively decreases going towards the walls. When the flow encounters the step wall for the first time, there are minor vortices formed near the step wall and at the top of the wall which dissipates as the flow approaches steady state. Slowly the flow approaches steady state, at around 125th timestep (2.5s) into the simulation. There is no significant change in flow after it attains steady state.

Velocity Profile at x= 0.085 for both Cases

Plot over line tool was used to generate a csv file for both cases, the file contains the distribution of velocity for both cases at x=0.085m. A simple plot was generated from those values.

case 1:

case 2:

Observations: We can notice that case 1's curve is quite smoother than case 2's curve at all points. Additionally, the maximum value for Case 1: 0.285561 is slightly higher than that of Case 2: 0.279766

Week 8 - Simulation of a backward facing step in OpenFOAM

Ship work like this
in your own course.

Talk to a career counsellor

Ship work like thisin your own course.

Talk to a career counsellor

Ship work like this
in your own course.