Simulating the flow over a backward facing step using the icoFoam solver in OpenFOAM and studying the effect of mesh grading factor
Objective
To simulate the flow through a backward facing step using the icoFoam solver in OpenFOAM and to study the effect of mesh grading factor on the velocity magnitude profile.
Mesh specification
- Number of cells along x-direction (longer dimension) = 200
- Number of cells along y-direction = 20
Creating the mesh
The geometry of the backward facing step as given in the problem is as shown below;
This geometry is divided into 5 blocks with the order of vertices as shown below;
In the above geometry, the top wall consisting of faces (16 15 4 5), (15 21 10 4), the bottom wall consisting of faces (0 1 12 11), (6 7 18 17) and the left step face (1 6 12 17) are considered no-slip walls and the grading factor is to be applied near these walls. The front and back faces are not considered and denoted by empty as this is a 2D simulation. All this information is provided in the blockMeshDict file located inside the system folder of the case.
blockMeshDict file
/*--------------------------------*- C++ -*----------------------------------*\\
| ========= | |
| \\\\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\\\ / O peration | Version: 4.1 |
| \\\\ / A nd | Web: www.OpenFOAM.org |
| \\\\/ M anipulation | |
\\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
convertToMeters 0.001;
vertices
(
(0 0 0)
(80 0 0)
(80 5 0)
(0 5 0)
(80 10 0)
(0 10 0)
(80 -10 0)
(200 -10 0)
(200 0 0)
(200 5 0)
(200 10 0)
(0 0 1)
(80 0 1)
(80 5 1)
(0 5 1)
(80 10 1)
(0 10 1)
(80 -10 1)
(200 -10 1)
(200 0 1)
(200 5 1)
(200 10 1)
);
blocks
(
hex (0 1 2 3 11 12 13 14) (80 5 1) simpleGrading (1 5 1)
hex (3 2 4 5 14 13 15 16) (80 5 1) simpleGrading (1 0.2 1)
hex (1 8 9 2 12 19 20 13) (120 5 1) simpleGrading (5 5 1)
hex (2 9 10 4 13 20 21 15) (120 5 1) simpleGrading (5 0.2 1)
hex (6 7 8 1 17 18 19 12) (120 10 1) simpleGrading (5 5 1)
);
edges
(
);
boundary
(
inlet
{
type patch;
faces
(
(0 11 14 3)
(3 14 16 5)
);
}
outlet
{
type patch;
faces
(
(8 9 20 19)
(9 10 21 20)
(7 8 19 18)
);
}
noSlipWalls
{
type wall;
faces
(
(0 1 12 11)
(6 7 18 17)
(16 15 4 5)
(15 21 10 4)
(1 6 17 12)
);
}
frontAndBack
{
type empty;
faces
(
(11 12 13 14)
(14 13 15 16)
(12 19 20 13)
(13 20 21 15)
(17 18 19 12)
(0 3 2 1)
(3 5 4 2)
(1 2 9 8)
(2 4 10 9)
(6 1 8 7)
);
}
);
mergePatchPairs
(
);
// ************************************************************************* //Time-step
The time control information is provided inside the \'controlDict\' file inside the system folder. For this case, the time-step size is chosen such that the Courant number remains well below 1 throughout the geometry. The value for time-step size is given 1e-4 and simulation time is given 0.1s. The \'controlDict\' file is shown below;
controlDict file
/*--------------------------------*- C++ -*----------------------------------*\\
| ========= | |
| \\\\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\\\ / O peration | Version: 4.1 |
| \\\\ / A nd | Web: www.OpenFOAM.org |
| \\\\/ M anipulation | |
\\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location \"system\";
object controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
application icoFoam;
startFrom startTime;
startTime 0;
stopAt endTime;
endTime 0.01;
deltaT 0.0001;
writeControl timeStep;
writeInterval 20;
purgeWrite 0;
writeFormat ascii;
writePrecision 6;
writeCompression off;
timeFormat general;
timePrecision 6;
runTimeModifiable true;
// ************************************************************************* //Discretization scheme
The discretization information is contained inside the \'fvSchemes\' file in the systems folder. The first-order accurate Forward Euler scheme is applied for time marching and Linear Gauss scheme, that is second-order accurate is applied for space marching. The file is shown below;
/*--------------------------------*- C++ -*----------------------------------*\\
| ========= | |
| \\\\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\\\ / O peration | Version: 4.1 |
| \\\\ / A nd | Web: www.OpenFOAM.org |
| \\\\/ 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;
}
// ************************************************************************* //Iterative Solver
This data is contained inside the \'fvSolution\' file inside the system folder. For this case, we\'re using a \'Pre-Conjugate Gradient\' solver with \'Diagonal-base Incomplete Cholesky\' pre-conditioner. The file is as follows;
/*--------------------------------*- C++ -*----------------------------------*\\
| ========= | |
| \\\\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\\\ / O peration | Version: 4.1 |
| \\\\ / A nd | Web: www.OpenFOAM.org |
| \\\\/ 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;
}
// ************************************************************************* //Assigning Initial and boundary conditions
The initial and boundary conditions for pressure and velocity are available inside the \'p\' and \'U\' files inside the \'0\' folder. The initial condition is denoted by \'internalField\' value inside these files. For this case, both the velocity and pressure are considered to be zero, initially.
The boundary conditions are given under \'boundaryField\' inside \'p\' and \'U\' files. At the inlet, the velocity is fixed and pressure floated. At the outlet, the pressure is fixed and the velocity is floated. Floating means applying a Neumann boundary condition where the quantity at the outer nodes is calculated by linear extrapolation of inner nodes. Providing a fixed value means applying a Dirichlet BC. The top and bottom walls are applied no-slip boundary condition for velocity and \'zeroGradient\' (Neumann BC) for pressure.
The value of the velocity at inlet is given 1m/s and the value for pressure is given 0 Pa. The respective files are shown below;
p file
/*--------------------------------*- C++ -*----------------------------------*\\
| ========= | |
| \\\\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\\\ / O peration | Version: 4.1 |
| \\\\ / A nd | Web: www.OpenFOAM.org |
| \\\\/ 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;
}
noSlipWalls
{
type zeroGradient;
}
frontAndBack
{
type empty;
}
}
// ************************************************************************* //U file
/*--------------------------------*- C++ -*----------------------------------*\\
| ========= | |
| \\\\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\\\ / O peration | Version: 4.1 |
| \\\\ / A nd | Web: www.OpenFOAM.org |
| \\\\/ M anipulation | |
\\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volVectorField;
object U;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dimensions [0 1 -1 0 0 0 0];
internalField uniform (0 0 0);
boundaryField
{
inlet
{
type fixedValue;
value uniform (1 0 0);
}
outlet
{
type zeroGradient;
}
noSlipWalls
{
type noSlip;
}
frontAndBack
{
type empty;
}
}
// ************************************************************************* //Result
The geometry mesh for three different grading factors is as follows;
Grading Factor = 0.2
Grading Factor = 0.5
Grading Factor = 0.8
The velocity distribution is as follows;
Grading Factor = 0.2
Grading Factor = 0.5
Grading Factor = 0.8
The velocity magnitude profile is as follows;
Grading Factor = 0.2
Grading Factor = 0.5
Grading Factor = 0.8
Velocity profile comparison
The velocity profile data for each mesh was imported into an excel file and compared using plot function in Matlab.
From the velocity profile comparison, it can be observed that velocity profiles of mesh with higher grading factors appear to be smoother than the one with lower grading factor. This is because a mesh with a lower grading factor depicts more accurate results.
Therefore with grading, we are able to refine the mesh in key regions without increasing the number of grid points or cells inside a block. This implies that, for the same number of cells, we can manipulate the mesh size at the desired region to get a more accurate result.