BlogFor BusinessTalk to a counsellor
CoursesIIT ProgramsPlacementsLog in
BlogFor Business
Talk to a counsellor
Challenge

Week 12 - Validation studies of Symmetry BC vs Wedge BC in OpenFOAM vs Analytical H.P equation

KKiran More20 Feb 2021
CFDMATLAB

Aim: Validation studies of Symmetry BC vs Wedge BC in OpenFOAM vs Analytical H.P equation

Theory:

For a flow through circular pipe velocity profile is parabolic. Velocity is minimum at surface and maximum at centre of pipe. Similarly shear stress variation is linear. Shear stress is zero at centre of pipe and maximum at surface.

Following parameters were assumed

Reynolds number, Re= 2100

Diameter of pipe = 0.01m

Density of water = 1000 kg/m^3.

Dynamic viscosity of water = 8.90e-4 Pa.s

Entrance is a length within which velocity profile is developing. For a length greater or equal to entrance length, velocity profile is fully developed.

Entrance length Le is given by

Le = 0.05 Re D=1.05 m

Average velocity over a pipe cross section can be calculated as

`V_avg=(R_e μ)/(ρD)=0.1869 m/sec`

Also for fully developed flow, pressure drop over a length L is

`∆p=(32 μ V_avg L)/D^2 =79.84 Pa`

Kinematic pressure drop is`79.84/1000=0.07984`

In this work, flow throgh circular pipe is silumated for a wdge angle of 10,25 and 45 degree respectively with symmetry and wedge boundary condition.

MATLAB code to to generate blockMeshDict file (Symmetry BC)

clear all
close all
clc

% Pipe dimentions
D = 0.01;
r = D/2;
L = 1.5;  

% Input parameters
Re = 2100;
mu = 8.9e-4;
rho = 1000;
v_avg = (Re*mu)/(rho*D);
Le = 0.05*Re*D;
V_max = 2*v_avg;
dp = (32*mu*v_avg*L)/(D^2);

% Wedge angle input
theta = input("Input wedge angle in degrees:");
grading = input("Enter the grading factor");

% To calcuculte vertice of wedge
v0 = [0 0 0];
v1 = [r*sind(theta/2) r*cosd(theta/2) 0];
v2 = [-r*sind(theta/2) r*cosd(theta/2) 0];
v3 = [0 0 L];
v4 = [r*sind(theta/2) r*cosd(theta/2) L];
v5 = [-r*sind(theta/2) r*cosd(theta/2) L];

% To define the block
blk = sprintf("hex (0 0 1 2 3 3 4 5) (1 50 250) simpleGrading (1 %f 1)", grading);
edg1 = "arc 1 2";
edg2 = "arc 4 5";
edg1_mid = [0 r 0];
edg2_mid = [0 r L];

inlet = "(0 0 2 1)";
outlet = "(3 4 5 3)";
wedge_1 = "(0 1 4 3)";
wedge_2 = "(3 5 2 0)";
wall = "(4 1 2 5)";
axis = "(0 3 3 0)";

% Creating blockMeshDict file
h1='/*--------------------------------*- C++ -*----------------------------------*\';
h2='=========                 |';
h3='\\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox';
h4=' \\    /   O peration     | Website:  https://openfoam.org';
h5='  \\  /    A nd           | Version:  7';
h6='   \\/     M anipulation  |';
h7='\*---------------------------------------------------------------------------*/';
h8='// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //';
h9='// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //';
h10='// ************************************************************************* //';

%final = '//****************************************//';

%To write blockmeshdict file

% To write header file
f1 = fopen('blockMeshDict45.txt','w');
fprintf(f1, '%s\n',h1);
fprintf(f1, '%s/\n',h2);
fprintf(f1, '%s\n',h3);
fprintf(f1, '%s\n',h4);
fprintf(f1, '%s\n',h5);
fprintf(f1, '%s\n',h6);
fprintf(f1, '%s\n',h7);
fprintf(f1, '%s\n',h8);

fprintf(f1, 'FoamFile \n {');
fprintf(f1, '%15s %10s \n','version','2.0;');
fprintf(f1, '%16s %10s \n','format',' ascii;');
fprintf(f1, '%15s %10s \n','class','dictionary;');
fprintf(f1, '%17s %10s \n','object ',' blockMeshDict;');
fprintf(f1, '}\n');
fprintf(f1, '%s\n',h9);
fprintf(f1, '%2s %2s \n','convertToMeters','1;');

% To define vertices
fprintf(f1, 'vertices \n (\n');
fprintf(f1, '\t(%d %d %d)\n',v0);
fprintf(f1, '\t(%d %d %d)\n',v1);
fprintf(f1, '\t(%d %d %d)\n',v2);
fprintf(f1, '\t(%d %d %d)\n',v3);
fprintf(f1, '\t(%d %d %d)\n',v4);
fprintf(f1, '\t(%d %d %d)\n',v5);
fprintf(f1, ' \n );\n');

% To define blocks
fprintf(f1, 'blocks \n (\n');
fprintf(f1, '%10s',blk);
fprintf(f1, ' \n ); \n');

% To define edges
fprintf(f1, 'edges \n (\n');
fprintf(f1,'\t %s (%f %f %f)\n','edg1',edg1_mid);
fprintf(f1,'\t %s (%f %f %f )\n','edg2',edg2_mid);
fprintf(f1, ' \n );\n');

% To define boundaries
fprintf(f1, 'boundary \n (\n');
fprintf(f1,'\t %5s \n','inlet','{');
fprintf(f1, '\t %10s %s \n','type','patch;');
fprintf(f1, '\t %11s\n %18s\n','faces','(');
fprintf(f1, '\t %18s %s', inlet);
fprintf(f1, ' \n );\n');
fprintf(f1, '}\n');


fprintf(f1, 'outlet \n {\n');
fprintf(f1, '%s %s \n','type','patch;');
fprintf(f1, 'faces \n (\n');
fprintf(f1, '%s', outlet);
fprintf(f1, ' \n );\n');
fprintf(f1, '}\n');

fprintf(f1, 'front \n {\n');
fprintf(f1, '%s %s \n','type','symmetry;');
fprintf(f1, 'faces \n (\n');
fprintf(f1, '%s', wedge_1);
fprintf(f1, ' \n );\n');
fprintf(f1, '}\n');

fprintf(f1, 'back \n {\n');
fprintf(f1, '%s %s \n','type','symmetry;');
fprintf(f1, 'faces \n (\n');
fprintf(f1, '%s', wedge_2);
fprintf(f1, ' \n );\n');
fprintf(f1, '}\n');

fprintf(f1, 'wall \n {\n');
fprintf(f1, '%s %s \n','type','wall;');
fprintf(f1, 'faces \n (\n');
fprintf(f1, '%s', wall);
fprintf(f1, ' \n );\n');
fprintf(f1, '}\n');

fprintf(f1, 'axis \n {\n');
fprintf(f1, '%s %s \n','type','empty;');
fprintf(f1, 'faces \n (\n');
fprintf(f1, '%s', axis);
fprintf(f1, ' \n );\n');
fprintf(f1, '}\n');
fprintf(f1, ');\n');

fprintf(f1, 'mergePatchPairs \n (\n); \n');
fprintf(f1, '%s\n',h10);

MATLAB code to to generate blockMeshDict file (Wedge BC)

clear all
close all
clc

% Pipe dimentions
D = 0.01;
r = D/2;
L = 1.5;  

% Input parameters
Re = 2100;
mu = 8.9e-4;
rho = 1000;
v_avg = (Re*mu)/(rho*D);
Le = 0.05*Re*D;
V_max = 2*v_avg;
dp = (32*mu*v_avg*L)/(D^2);

% Wedge angle input
theta = input("Input wedge angle in degrees:");
grading = input("Enter the grading factor");

% To calcuculte vertice of wedge
v0 = [0 0 0];
v1 = [r*sind(theta/2) r*cosd(theta/2) 0];
v2 = [-r*sind(theta/2) r*cosd(theta/2) 0];
v3 = [0 0 L];
v4 = [r*sind(theta/2) r*cosd(theta/2) L];
v5 = [-r*sind(theta/2) r*cosd(theta/2) L];

% To define the block
blk = sprintf("hex (0 0 1 2 3 3 4 5) (1 50 250) simpleGrading (1 %f 1)", grading);
edg1 = "arc 1 2";
edg2 = "arc 4 5";
edg1_mid = [0 r 0];
edg2_mid = [0 r L];

inlet = "(0 0 2 1)";
outlet = "(3 4 5 3)";
wedge_1 = "(0 1 4 3)";
wedge_2 = "(3 5 2 0)";
wall = "(4 1 2 5)";
axis = "(0 3 3 0)";

% Creating blockMeshDict file
h1='/*--------------------------------*- C++ -*----------------------------------*\';
h2='=========                 |';
h3='\\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox';
h4=' \\    /   O peration     | Website:  https://openfoam.org';
h5='  \\  /    A nd           | Version:  7';
h6='   \\/     M anipulation  |';
h7='\*---------------------------------------------------------------------------*/';
h8='// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //';
h9='// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //';
h10='// ************************************************************************* //';

%final = '//****************************************//';

%To write blockmeshdict file

% To write header file
f1 = fopen('blockMeshDict45.txt','w');
fprintf(f1, '%s\n',h1);
fprintf(f1, '%s/\n',h2);
fprintf(f1, '%s\n',h3);
fprintf(f1, '%s\n',h4);
fprintf(f1, '%s\n',h5);
fprintf(f1, '%s\n',h6);
fprintf(f1, '%s\n',h7);
fprintf(f1, '%s\n',h8);

fprintf(f1, 'FoamFile \n {');
fprintf(f1, '%15s %10s \n','version','2.0;');
fprintf(f1, '%16s %10s \n','format',' ascii;');
fprintf(f1, '%15s %10s \n','class','dictionary;');
fprintf(f1, '%17s %10s \n','object ',' blockMeshDict;');
fprintf(f1, '}\n');
fprintf(f1, '%s\n',h9);
fprintf(f1, '%2s %2s \n','convertToMeters','1;');

% To define vertices
fprintf(f1, 'vertices \n (\n');
fprintf(f1, '\t(%d %d %d)\n',v0);
fprintf(f1, '\t(%d %d %d)\n',v1);
fprintf(f1, '\t(%d %d %d)\n',v2);
fprintf(f1, '\t(%d %d %d)\n',v3);
fprintf(f1, '\t(%d %d %d)\n',v4);
fprintf(f1, '\t(%d %d %d)\n',v5);
fprintf(f1, ' \n );\n');

% To define blocks
fprintf(f1, 'blocks \n (\n');
fprintf(f1, '%10s',blk);
fprintf(f1, ' \n ); \n');

% To define edges
fprintf(f1, 'edges \n (\n');
fprintf(f1,'\t %s (%f %f %f)\n','edg1',edg1_mid);
fprintf(f1,'\t %s (%f %f %f )\n','edg2',edg2_mid);
fprintf(f1, ' \n );\n');

% To define boundaries
fprintf(f1, 'boundary \n (\n');
fprintf(f1,'\t %5s \n','inlet','{');
fprintf(f1, '\t %10s %s \n','type','patch;');
fprintf(f1, '\t %11s\n %18s\n','faces','(');
fprintf(f1, '\t %18s %s', inlet);
fprintf(f1, ' \n );\n');
fprintf(f1, '}\n');


fprintf(f1, 'outlet \n {\n');
fprintf(f1, '%s %s \n','type','patch;');
fprintf(f1, 'faces \n (\n');
fprintf(f1, '%s', outlet);
fprintf(f1, ' \n );\n');
fprintf(f1, '}\n');

fprintf(f1, 'front \n {\n');
fprintf(f1, '%s %s \n','type','wedge;');
fprintf(f1, 'faces \n (\n');
fprintf(f1, '%s', wedge_1);
fprintf(f1, ' \n );\n');
fprintf(f1, '}\n');

fprintf(f1, 'back \n {\n');
fprintf(f1, '%s %s \n','type','wedge;');
fprintf(f1, 'faces \n (\n');
fprintf(f1, '%s', wedge_2);
fprintf(f1, ' \n );\n');
fprintf(f1, '}\n');

fprintf(f1, 'wall \n {\n');
fprintf(f1, '%s %s \n','type','wall;');
fprintf(f1, 'faces \n (\n');
fprintf(f1, '%s', wall);
fprintf(f1, ' \n );\n');
fprintf(f1, '}\n');

fprintf(f1, 'axis \n {\n');
fprintf(f1, '%s %s \n','type','empty;');
fprintf(f1, 'faces \n (\n');
fprintf(f1, '%s', axis);
fprintf(f1, ' \n );\n');
fprintf(f1, '}\n');
fprintf(f1, ');\n');

fprintf(f1, 'mergePatchPairs \n (\n); \n');
fprintf(f1, '%s\n',h10);

Mesh

1) 10 degree wedge angle

2) 25degree wedge angle

3) 45 degree wedge angle

Pressure Boundary and initial conditions (Symmetry BC)

/*--------------------------------*- 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;  //Neumann 
    }

    outlet
    {
        type            fixedValue;
        value           uniform 0.07984;
    }

    wall
    {
        type            zeroGradient; //Neumann
    }
    
    front
    {
        type           symmetry;
    }
    
    back
    {
        type           symmetry;
    }
    
    axis
    {
        type           empty;
    }
}

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

Pressure Boundary and initial conditions (Wedge BC)

/*--------------------------------*- 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;  //Neumann 
    }

    outlet
    {
        type            fixedValue;
        value           uniform 0.07984;
    }

    wall
    {
        type            zeroGradient; //Neumann
    }
    
    front
    {
        type           wedge;
    }
    
    back
    {
        type           wedge;
    }
    
    axis
    {
        type           empty;
    }
}

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

Velocity Booundary and initial conditions (Symmetry BC)

/*--------------------------------*- 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 (0 0 0.1869); //Dirichlet
    }

    outlet
    {
        type            zeroGradient;  //Neumann
    }

    wall
    {
        type           noSlip;
    }
    
    front
    {
        type           symmetry;
    }
    
    back
    {
        type           symmetry;
    }
    
    axis
    {
        type           empty;
    }
}

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

Velocity Booundary and initial conditions (Wedge BC)

/*--------------------------------*- 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 (0 0 0.1869); //Dirichlet
    }

    outlet
    {
        type            zeroGradient;  //Neumann
    }

    wall
    {
        type           noSlip;
    }
    
    front
    {
        type           wedge;
    }
    
    back
    {
        type           wedge;
    }
    
    axis
    {
        type           empty;
    }
}

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

ControllDict file

/*--------------------------------*- 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;
    location    "system";
    object      controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

application     icoFoam;

startFrom       startTime;

startTime     0;

stopAt          endTime;

endTime         30;

deltaT          0.01;

writeControl    timeStep;

writeInterval   20;

purgeWrite      0;

writeFormat     ascii;

writePrecision  6;

writeCompression off;

timeFormat      general;

timePrecision   6;

runTimeModifiable true;


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

Results (For symmetry Boundary condition)

1) Wedge angle 10 degree

Velocity profile at L=0 m

Velocity profile at L= 0.1 m

Velocity profile at L= 0.5 m

Velocity profile at L= 0.8 m

Velocity profile at L= 1.05 m

Velocity profile at L= 1.2 m

Shear stress profile at L= 1.2 m

2) Wedge angle 25 degree

Velocity profile at L=0 m

Velocity profile at L= 0.1 m

Velocity profile at L= 0.5 m

Velocity profile at L= 0.8m

Velocity profile at L= 1.05 m

Velocity profile at L=1.2 m

Shear stress profile at L= 1.2 m

3) Wedge angle 45 degree

Velocity profile at L=1.2 m

Shear stress profile at L= 1.2 m

Similar velocity and shear stress profile can be obtained for wedge boundary conditions

Discussion

It is found that till L=1.05 m (Entrance length), velocity profile was changing continuously. And after L=1.05 m velocity profile is constant and it is not changing with length of pipe.

Also theoretical value of maximum velocity and value obtained using numerical simulation is matching.

Type of BC

1) Symmetry

Maximum velocity (Analytical) = 0.3738

Wedge angle (Degree) Maximum velocity Simulation time (Sec)
10 0.371545 435
25 0.371866 605
45 0.372546 583

2) Wedge

Maximum velocity (Analytical) = 0.3738

Wedge angle (Degree) Maximum velocity Simulation time (Sec)
10 0.371545 616
25 0.371866 531
45 0.372546 501
Your portfolio, next

Ship work like this
in your own course.

Talk to a career counsellor

Free 20-minute call. We'll map the right course to where you want to work.

By submitting you agree to our Terms and Privacy Policy.

Talk to a career counsellor

Free 20-minute call. We'll map the right course to where you want to work.

By submitting you agree to our Terms and Privacy Policy.