BlockMesh Drill down challenge- openFOAM

Objective:  To use the icoFOAM solver to simulate the flow through a backward-facing step. To create multiple meshes and compare the results obtained.

Desired output

• How does the velocity magnitude profile change as a function of the mesh grading factor? Use factors, 0.2, 0.5,0.8
• Measure the velocity profile at 0.085 m from the inlet of the geometry
• A plot must be a line plot
• Compare velocity magnitude contours near the step region for the above-mentioned grading factors. 

GIVEN:

Specifications to follow for the figure:

Phases to Procced:

1. Pre-processing:
     The phase where problem statement is converted into Idealised or discretized computer mode
     Assumptions made based on the type of flow to be modeled. Another process involves mesh generation, application of initial & boundary condition, etc.

2. Solving:
    The actual computation is done by solvers. For this solving phase, actual computational power is required. Multiple solvers are available varying in efficiency and capability of solving certain physical phenomena.

3. Post-processing:
    Obtained results are analyzed and visualized in this phase. At this stage, we can verify results and conclusions can be drawn based on obtained results.

Mesh specification is as follows:

1. Number of cells along the x-direction (longer dimension) = 200
2. Number of cells along the y-direction = 10
3. Use a graded mesh. Grading factor near all walls should be 0.2

icoFOAM solver:
icoFoam solves the incompressible laminar Navier-Stokes equations using the PISO algorithm. The code is inherently transient, requiring an initial condition (such as zero velocity) and boundary conditions. The icoFOAM code can take mesh non-orthogonality into account with successive non-orthogonality iterations. The number of PISO corrections and non-orthogonality corrections is controlled through user input.

Path to find solver:
`tut->`Incompressible `-> icoFOAM` OR `sol ->` Incompressible `-> icoFOAM`

Mesh creation:

In order to create a given figure as per given mesh specifications, we need to edit a blockMeshDict file located inside the system folder of the case we are going to run.

blockMeshDict file:

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

convertToMeters 1;

vertices
(
    (0 0 0) 		// Ver 0
    (0.08 0 0)		// Ver 1
    (0.08 0.005 0)	// Ver 2
    (0 0.005 0)		// Ver 3
    (0.08 0.01 0)	// Ver 4
    (0 0.01 0)		// Ver 5
    (0.2 0.01 0)	// Ver 6
    (0.2 0.005 0)	// Ver 7
    (0.2 0 0)		// Ver 8
    (0.2 -0.01 0)	// Ver 9
    (0.08 -0.01 0)	// Ver 10
    (0 0 0.001)		// Ver 11
    (0.08 0 0.001)	// Ver 12
    (0.08 0.005 0.001)	// Ver 13
    (0 0.005 0.001)	// Ver 14
    (0.08 0.01 0.001)	// Ver 15
    (0 0.01 0.001)	// Ver 16
    (0.2 0.01 0.001)	// Ver 17
    (0.2 0.005 0.001)	// Ver 18
    (0.2 0 0.001)	// Ver 19
    (0.2 -0.01 0.001)	// Ver 20
    (0.08 -0.01 0.001)	// Ver 21
);

blocks
(
 	hex (0 1 2 3 11 12 13 14) (130 10 1) simpleGrading (1 5 1) 	//Block 1

	hex (3 2 4 5 14 13 15 16) (130 10 1) simpleGrading (1 0.2 1) 	//Block 2

	hex (2 7 6 4 13 18 17 15) (200 10 1) simpleGrading (5 0.2 1) 	//Block 3

	hex (1 8 7 2 12 19 18 13) (200 10 1) simpleGrading (5 5 1) 	//Block 4

	hex (10 9 8 1 21 20 19 12) (200 20 1) simpleGrading (5 5 1) 	//Block 5
);

edges
(
);

boundary
(
    inletWall
    {
        type patch;
        faces
        (
            (0 11 14 3)
	    (3 14 16 5)  
        );
    }
    outletWalls
    {
        type patch;
        faces
        (
            (9 8 19 20)
            (8 7 18 19)
            (7 6 17 18)
        );
    }
    frontAndBack
    {
        type empty;
        faces
        (
		(0 3 2 1)
		(3 5 4 2)
		(2 4 6 7)
		(1 2 7 8)
		(10 1 8 9)
		(11 12 13 14)
		(14 13 15 16)
		(13 18 17 15)
		(12 19 18 13)
		(21 20 19 12)
        );
    }
    Noslipwalls
    {
      type wall;
      faces 
      (
	(5 16 15 4)
	(4 15 17 6)
	(0 1 12 11)
	(10 9 20 21)
	(1 10 21 12)


      );
    }
);

mergePatchPairs
(
);

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

Along Y-direction only to create finner mesh towards the wall.
using Grading Factor: 0.2

Closer look at created mesh

Grading Factor: 0.5

Keeping the remaining blockMeshDict file as same just changes the values of simpleGrading.
The script will be:

blocks
(
 	hex (0 1 2 3 11 12 13 14) (130 10 1) simpleGrading (1 2 1) 	//Block 1

	hex (3 2 4 5 14 13 15 16) (130 10 1) simpleGrading (1 0.5 1) 	//Block 2

	hex (2 7 6 4 13 18 17 15) (200 10 1) simpleGrading (2 0.5 1) 	//Block 3

	hex (1 8 7 2 12 19 18 13) (200 10 1) simpleGrading (2 2 1) 	//Block 4

	hex (10 9 8 1 21 20 19 12) (200 20 1) simpleGrading (2 2 1) 	//Block 5
);

A closer look at created mesh

Grading Factor: 0.8

The script will be:

blocks
(
 	hex (0 1 2 3 11 12 13 14) (130 10 1) simpleGrading (1 1.25 1) 	//Block 1

	hex (3 2 4 5 14 13 15 16) (130 10 1) simpleGrading (1 0.8 1) 	//Block 2

	hex (2 7 6 4 13 18 17 15) (200 10 1) simpleGrading (1.25 0.8 1) 	//Block 3

	hex (1 8 7 2 12 19 18 13) (200 10 1) simpleGrading (1.25 1.25 1) 	//Block 4

	hex (10 9 8 1 21 20 19 12) (200 20 1) simpleGrading (1.25 1.25 1) 	//Block 5
);

A closer look at created mesh

As we can see from the above picture lesser the grading factor finner will be mesh as the grading factor is given by

where, `delta_e=`ending cell and `delta_s=`starting cell.

Applying Initial And Boundary Conditions:

Setting appropriate boundary conditions is vital for a successful simulation. Ill-posed boundary conditions will lead to physically incorrect predictions, and in many cases solver failure. Users must specify the boundary conditions for each solved field. The tutorials provided with OpenFOAM show examples of good practice in terms of selection and application for various cases.

Here, Initial and Boundary conditions can be applied by editing pressure and velocity script files contained inside 0 folder of case.

Pressure script file:

This is Kinematic Pressure where we are going to apply the necessary boundary condition.

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

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

internalField   uniform 1;

boundaryField
{
    inletWall
    {
        type            zeroGradient;
    }

    outletWalls
    {
        type           fixedValue;
	value	       uniform 0;

    }

    frontAndBack
    {
        type            empty;
    }
    Noslipwalls
    {
        type            zeroGradient;
    }
}

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

Velocity script file:

This is Velocity where we are considering 3 m/s velocity at inlet.

/*--------------------------------*- C++ -*----------------------------------*\\
  =========                 |
  \\\\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\\\    /   O peration     | Website:  https://openfoam.org
    \\\\  /    A nd           | Version:  6
     \\\\/     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
{
    inletWall
    {
        type            fixedValue;
        value           uniform (3 0 0);
    }

    outletWalls
    {
        type            zeroGradient;
    }

    frontAndBack
    {
        type            empty;
    }
    Noslipwalls
    {
        type            noSlip;
	value		uniform (0 0 0);
    }
}

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

ControlDict file:

This is where we control our simulation by defining timestep, end time, etc. properly.

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

application     icoFoam;

startFrom       startTime;

startTime       0;

stopAt          endTime;

endTime         0.02;

deltaT          1e-5;

writeControl    timeStep;

writeInterval   15;

purgeWrite      0;

writeFormat     ascii;

writePrecision  6;

writeCompression off;

timeFormat      general;

timePrecision   6;

runTimeModifiable true;


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

POST-PROCESSING results:

The velocity contour plot is as below based on different grading factors:

Grading Factor: 0.2

Grading Factor: 0.5

Grading Factor: 0.8

The above velocity contour plot shows the velocity distribution along with the domain in X-axis, here data capturing capacity of grading factor 0.2 is more comparatively than that of other grading factors i.e. 0.5 and 0.8. This will be observed if we take a closer look at the backward-facing step.

Velocity Plot:

Velocity Profile at X=0.085m from the inlet of geometry based on grading factors are as below:

Grading Factor: 0.2

Grading Factor: 0.5

Grading Factor: 0.8

From the above graphs, we can make a suitable conclusion:

1. For all grading factors, we have a maximum velocity at a particular region of the domain, the only difference is that their velocity distribution along the Y-axis of the domain.

2. For grading factor 0.8 we are getting smooth velocity distribution but it is lacking in capturing velocity near walls. In contrast, grading factor 0.2 has maximum potential to capture all velocity data near-wall also that\'s why it has a dome-shaped like plot.

3. The maximum velocity value for grading factor 0.2 is less because of the distribution of information to maximum mesh elements compare to other grading factors.