Simulation of Flow in a Backward Facing Step using OpenFoam
Objective:
To simulate a backward facing step using icoFoam solver in OpenFoam using different mesh grading.
Introduction:
The geometry that is being analysed is as follows:
The mesh is created with the following details:
Number of elements in the x direction:200
Number of elements in the y direction :20(this is selected to get significant difference in the mesh grading cases).
The icoFoam solver from openFoam is used.This is an incompressible solver. The governing equations solved by the solver are as follows.
Continuity Equation:
`nablau=0`
Momentum Equation:
`(delu)/(delt)+nabla(u times u)-nabla(nu nabla u)=-nabla p`
where u is velocity
`nu`is kinematic viscosity
p is kinematic pressure (Pressure/Density).
Method and Case Files:
The blockmesh utility of openFoam is used to create the mesh. The computational domain is divided into 5 blocks as in the figure below.
The blockmesh file for a non graded mesh is as given below:
/*--------------------------------*- C++ -*----------------------------------*\\
========= |
\\\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\\\ / O peration | Website: https://openfoam.org
\\\\ / A nd | Version: 7
\\\\/ M anipulation |
\\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
convertToMeters 0.01;
vertices
(
(0 1 0)
(8 1 0)
(8 1.5 0)
(0 1.5 0)
(8 2 0)
(0 2 0)
(20 2 0)
(20 1.5 0)
(20 1 0)
(20 0 0)
(8 0 0)
(0 1 0.1)
(8 1 0.1)
(8 1.5 0.1)
(0 1.5 0.1)
(8 2 0.1)
(0 2 0.1)
(20 2 0.1)
(20 1.5 0.1)
(20 1 0.1)
(20 0 0.1)
(8 0 0.1)
);
blocks
(
hex (0 1 2 3 11 12 13 14 ) (80 5 1) simpleGrading (1 1 1)
hex (3 2 4 5 14 13 15 16 ) (80 5 1) simpleGrading (1 1 1)
hex (2 7 6 4 13 18 17 15 ) (120 5 1) simpleGrading (1 1 1)
hex (1 8 7 2 12 19 18 13 ) (120 5 1) simpleGrading (1 1 1)
hex (10 9 8 1 21 20 19 12 ) (120 10 1) simpleGrading (1 1 1)
);
edges
(
);
boundary
(
inlet
{
type patch;
faces
(
(0 11 14 3 )
(3 14 16 5 )
);
}
outlet
{
type patch;
faces
(
(9 8 19 20 )
(8 7 18 19 )
(7 6 17 18 )
);
}
fixedWalls
{
type wall;
faces
(
(4 5 16 15)
(6 4 15 17)
(0 1 12 11)
(10 9 20 21)
(1 10 21 12)
);
}
frontAndBack
{
type empty;
faces
(
(0 3 2 1)
(3 5 4 2)
(4 6 7 2)
(1 2 7 8)
(1 8 9 10)
(11 12 13 14)
(14 13 15 16)
(13 18 17 15)
(12 19 18 13)
(21 20 19 12)
);
}
);
mergePatchPairs
(
);
// ************************************************************************* //The mesh obtained can be seen in the figure below:
The grading of the mesh is changed and three different mesh are generated and all the four cases are simulated with the same boundary conditions:
Mesh grading 0.2:
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 (2 7 6 4 13 18 17 15 ) (120 5 1) simpleGrading (5 0.2 1)
hex (1 8 7 2 12 19 18 13 ) (120 5 1) simpleGrading (5 5 1)
hex (10 9 8 1 21 20 19 12 ) (120 10 1) simpleGrading (5 5 1)
);Mesh grading 0.5:
blocks
(
hex (0 1 2 3 11 12 13 14 ) (80 5 1) simpleGrading (1 2 1)
hex (3 2 4 5 14 13 15 16 ) (80 5 1) simpleGrading (1 0.5 1)
hex (2 7 6 4 13 18 17 15 ) (120 5 1) simpleGrading (2 0.5 1)
hex (1 8 7 2 12 19 18 13 ) (120 5 1) simpleGrading (2 2 1)
hex (10 9 8 1 21 20 19 12 ) (120 10 1) simpleGrading (2 2 1)
);
Mesh grading 0.8:
blocks
(
hex (0 1 2 3 11 12 13 14 ) (80 5 1) simpleGrading (1 1.25 1)
hex (3 2 4 5 14 13 15 16 ) (80 5 1) simpleGrading (1 0.8 1)
hex (2 7 6 4 13 18 17 15 ) (120 5 1) simpleGrading (1.25 0.8 1)
hex (1 8 7 2 12 19 18 13 ) (120 5 1) simpleGrading (1.25 1.25 1)
hex (10 9 8 1 21 20 19 12 ) (120 10 1) simpleGrading (1.25 1.25 1)
);
The boundary and initial conditions for p and U are given in the 0 directory.
The pressure used is kinematic pressure as this is an incompressible solver and the density is considered a constant.
/*--------------------------------*- C++ -*----------------------------------*\\
========= |
\\\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\\\ / O peration | Website: https://openfoam.org
\\\\ / A nd | Version: 7
\\\\/ 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;
}
fixedWalls
{
type zeroGradient;
}
frontAndBack
{
type empty;
}
}
// ************************************************************************* //An inlet velocity of 1m/s is given in at the inlet at the instant 0.
/*--------------------------------*- C++ -*----------------------------------*\\
========= |
\\\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\\\ / O peration | Website: https://openfoam.org
\\\\ / A nd | Version: 7
\\\\/ 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;
}
fixedWalls
{
type noSlip;
}
frontAndBack
{
type empty;
}
}
// ************************************************************************* //In the control dict file the total simulation time is set to 1s and the time step is calculated by keeping the courant number less than 1.
`C=nu(Deltat)/(Deltax)^2`
`nu=0.01` as taken from transport properties in the constant directory.
The following table shows the variations in the timestep fot the different cases.
Results and Plots:
The post processing is done in paraview
The velocity profiles and the plots obtained for velocity bear the step i.e. 0.085 from the entrance.A line probe is placed at the point and the velocity profile along the y direction is plotted.
Results withou mesh grading:
Results with 0.2 Mesh Grading
Results with 0.5 Mesh Grading
Results with 0.8 Mesh Grading
From the above plots and velocity contours it can be seen that the velocity profiles are similar for all the grading. with the curve being sharper and smoother for different cases. It can be observeed that the velocity has a peak at 0.014m from the bottim with the peak value of 1.1m/s.