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

Week 12 - Symmetry vs Wedge vs HP equation

RRam Prasanth Anjaneyulu21 Jun 2020
  1. AIM

To Program and generate a computational mesh for any wedge angle and grading schemes. The following are the problem objectives:

  • Perform simulation work with Symmetry Boundary Condition.
  • You should have results obtained from solving the Hagen Poisuelle's equation
  • You should have results obtained from simulating the pipe flow with the wedge BC
  • You should obtain results for the pipe flow with symmetry BC
  • Write a Matlab program that takes an angle as input and generates a blockMesh file for the given angle.
  • Angles to test 10,25,45
  • Compare the above results and discuss your findings.
  1. INTRODUCTION

Hagen Poiseuille’s Law gives the pressure drop in incompressible and non Newtonian fluid in a laminar flow through a long cylindrical pipe of constant cross sectional area. In this challenge we would simulate a wedge section of a pipe to study the flow through pipe.

Assumptions

  • Fluid is incompressible and flow is laminar
  • Cross sectional changes are negligible
  • Reynolds’s Number = 2100
  • Wedge angle = 4 degree
  • Diameter of pipe = 0.02 m
  • Average Velocity remains constant in incompressible flow
  • Working fluid is water 25
  • Density of water = `997.0 kgm^-3`
  • dynamic viscosity of water `mu = 0.89*10^-3 Pa.S`
  • kinematic viscosity of water ` nu = 8.927*10^-7 m^2s^-1`
  • Fluid velocity in a pipe is zero at surface and increases to maximum at the centre of pipe.

FORMULAE:

  • Reynolds Number `Re = (V_(avg)D)/nu = (rho*V_(avg)*D)/mu`
  • Laminar flow entrance length ` ~= 0.06*Re*D`
  • The velocity profile for a fully developed laminar flow in a pipe is parabolic with a maximum at the centreline and minimum at R = 0.
  • `u(r) = R^2/(4mu)*((dp)/dx)(1-r^2/R^2)`
  • `V_(avg) = R^2/(8mu)*((dp)/dx)`
  • Velocity Profile à `u(r) = 2V_(avg)(1-r^2/R^2)`
  1. MATLAB CODE

3.1 blockMeshDict file generation using matlab for wedge angle =10 

%Matlab program to study flow through pipe
clear all
close all
clc

%% Input  %%
Re= 2100;     % Reynold's Number 
D= 0.02;       % Diameter of pipe
R=D/2;        % Radius 
theta= input('Enter  Wedge angle =');      % Wedge angle 
mu=0.89e-3;  %Dynamic viscocity @ 25 degree Celsius 
rho=997;   %Density of water  @ 25 degree Celsius 
nu= mu/rho;   %Kinematic viscocity  @ 25 degree Celsius 

% Analytical Velocity, pressure

L_e= 0.06*Re*D; % Entrance length 
L_t= L_e + 0.3; %  Total length 
v_avg= (Re*mu)/(rho*D); % Average velocity 
v_max= 2*v_avg; % Maximum velocity 
Delta_P= (32*mu*v_avg*L_t)/(D^2); % Pressure Drop
kin_P= Delta_P/rho; % Kinematic Pressure

% Finding vertices
v0= [0 0 0];
v1= [0 R*cosd(theta/2) -R*sind(theta/2)];
v2= [0 R*cosd(theta/2) R*sind(theta/2)];
v3= [L_t 0 0];
v4= [L_t R*cosd(theta/2) -R*sind(theta/2)];
v5= [L_t R*cosd(theta/2) R*sind(theta/2)];


%blockMeshDict file.

f1= fopen('blockMeshDictwege10.txt','w'); % create blockmesh file

fprintf(f1,'%s','/*--------------------------------*- C++ -*----------------------------------*');
fprintf(f1,'\n=========                 |\n');
fprintf(f1,'\n      /  F ield         | OpenFOAM: The Open Source CFD Toolbox');
fprintf(f1,' \n    /   O peration     | Website:  https://openfoam.org');
fprintf(f1,'  \n  /    A nd           | Version:  7');
fprintf(f1,'   \n/     M anipulation  |\n');
fprintf(f1,'%s','*---------------------------------------------------------------------------*/');
fprintf(f1,'\nFoamFile\n');

b1=blanks(5);
b2=blanks(6);
b3=blanks(7);
b4=blanks(4);

fprintf(f1,'{\n');
fprintf(f1,'%s','    version',b1,'2.0;');
fprintf(f1,'\n');
fprintf(f1,'%s','    format',b2,'ascii;');
fprintf(f1,'\n');
fprintf(f1,'%s','    class',b3,'dictionary;');
fprintf(f1,'\n');
fprintf(f1,'%s','    object',b2,'blockMeshDict;');
fprintf(f1,'\n');
fprintf(f1,'}\n');
fprintf(f1,'%s','// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //');

fprintf(f1,'\n\nconvertToMeters 1;\n\n');

% Finding vertices

fprintf(f1,'vertices\n');

fprintf(f1,'(\n');
fprintf(f1,b4,'%s','(');
v_0=fprintf(f1,'(%f %f %f)\n',v0(1),v0(2),v0(3));
fprintf(f1,b4,'%s','(');
v_1=fprintf(f1,'(%f %f %f)\n',v1(1),v1(2),v1(3));
fprintf(f1,b4,'%s','(');
v_2=fprintf(f1,'(%f %f %f)\n',v2(1),v2(2),v2(3));
fprintf(f1,b4,'%s','(');
v_3=fprintf(f1,'(%f %f %f)\n',v3(1),v3(2),v3(3));
fprintf(f1,b4,'%s','(');
v_4=fprintf(f1,'(%f %f %f)\n',v4(1),v4(2),v4(3));
fprintf(f1,b4,'%s','(');
v_5=fprintf(f1,'(%f %f %f)\n',v5(1),v5(2),v5(3));
fprintf(f1,'\n);\n');

%  blocks creation

fprintf(f1,'\nblocks\n');

fprintf(f1,'(\n');
fprintf(f1,b4,'%s','(');
fprintf(f1,'%s','hex (0 1 2 0 3 4 5 3)');

%  mesh generation
m1= input('Enter no. of mesh element along X-axis:');
m2= input('Enter no. of mesh element along Y-axis:');
m3= input('Enter no. of mesh element along Z-axis:');

fprintf(f1,blanks(1),'%s','(');
fprintf(f1,'(%2.0f %1.0f %3.0f)',m1,m2,m3);
fprintf(f1,'%s',' simpleGrading (1 1 1)');
fprintf(f1,'\n);');

%  edges creation

fprintf(f1,'\nedges\n');
fprintf(f1,'(\n');
fprintf(f1,'%s   arc 1 2 (0 %f 0)',blanks(1),R);
fprintf(f1,'\n%s   arc 4 5 (%f %f 0)',blanks(1),L_t,R);
fprintf(f1,'\n);\n\n');

% Apply boundary conditions

%  axis
fprintf(f1,'boundary\n(\n');
fprintf(f1,'%s','   axis');
fprintf(f1,'\n   {\n');
fprintf(f1,'%s',blanks(10),'type empty;');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'faces');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'(');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(13),'(0 3 3 0)');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),');');
fprintf(f1,'\n   }\n');

% fixedwalls
fprintf(f1,'%s','   fixedWalls');
fprintf(f1,'\n   {\n');
fprintf(f1,'%s',blanks(10),'type wall;');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'faces');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'(');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(13),'(2 1 4 5)');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),');');
fprintf(f1,'\n   }\n');

% front
fprintf(f1,'%s','   front');
fprintf(f1,'\n   {\n');
fprintf(f1,'%s',blanks(10),'type wedge;');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'faces');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'(');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(13),'(1 0 3 4)');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),');');
fprintf(f1,'\n   }\n');

% back
fprintf(f1,'%s','   back');
fprintf(f1,'\n   {\n');
fprintf(f1,'%s',blanks(10),'type wedge;');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'faces');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'(');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(13),'(2 5 3 0)');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),');');
fprintf(f1,'\n   }\n');

% inlet
fprintf(f1,'%s','   inlet');
fprintf(f1,'\n   {\n');
fprintf(f1,'%s',blanks(10),'type patch;');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'faces');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'(');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(13),'(0 1 2 0)');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),');');
fprintf(f1,'\n   }\n\n');

% outlet
fprintf(f1,'%s','   outlet');
fprintf(f1,'\n   {\n');
fprintf(f1,'%s',blanks(10),'type patch;');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'faces');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'(');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(13),'(3 5 4 3)');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),');');
fprintf(f1,'\n   }\n\n');

fprintf(f1,'\n);\n');

fprintf(f1,'\nmergePatchPairs');
fprintf(f1,'\n(\n);\n\n');

fprintf(f1,'%s','// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //');
fclose(f1);

3.2 BlockMeshDict File for wedge BC and wedge angle =10

Wedge boundary condition use the same Pressure and velocity initial conditions as done for Challenge 11.

/*--------------------------------*- 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 1;

vertices
(
    (0.000000 0.000000 0.000000)
    (0.000000 0.009962 -0.000872)
    (0.000000 0.009962 0.000872)
    (2.820000 0.000000 0.000000)
    (2.820000 0.009962 -0.000872)
    (2.820000 0.009962 0.000872)

);

blocks
(
    hex (0 1 2 0 3 4 5 3) (30 1 100) simpleGrading (1 1 1)
);
edges
(
    arc 1 2 (0 0.010000 0)
    arc 4 5 (2.820000 0.010000 0)
);

boundary
(
   axis
   {
          type empty;
          faces
          (
             (0 3 3 0)
          );
   }
   fixedWalls
   {
          type wall;
          faces
          (
             (2 1 4 5)
          );
   }
   front
   {
          type wedge;
          faces
          (
             (1 0 3 4)
          );
   }
   back
   {
          type wedge;
          faces
          (
             (2 5 3 0)
          );
   }
   inlet
   {
          type patch;
          faces
          (
             (0 1 2 0)
          );
   }

   outlet
   {
          type patch;
          faces
          (
             (3 5 4 3)
          );
   }


);

mergePatchPairs
(
);

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

3.3 BlockMeshDict File for symmetry BC and wedge angle =10

/*--------------------------------*- 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 1;

vertices
(
    (0.000000 0.000000 0.000000)
    (0.000000 0.009962 -0.000872)
    (0.000000 0.009962 0.000872)
    (2.820000 0.000000 0.000000)
    (2.820000 0.009962 -0.000872)
    (2.820000 0.009962 0.000872)

);

blocks
(
    hex (0 1 2 0 3 4 5 3) (30 1 100) simpleGrading (1 1 1)
);
edges
(
    arc 1 2 (0 0.010000 0)
    arc 4 5 (2.820000 0.010000 0)
);

boundary
(
   axis
   {
          type empty;
          faces
          (
             (0 3 3 0)
          );
   }
   fixedWalls
   {
          type wall;
          faces
          (
             (2 1 4 5)
          );
   }
   front
   {
          type symmetry;
          faces
          (
             (1 0 3 4)
          );
   }
   back
   {
          type symmetry;
          faces
          (
             (2 5 3 0)
          );
   }
   inlet
   {
          type patch;
          faces
          (
             (0 1 2 0)
          );
   }

   outlet
   {
          type patch;
          faces
          (
             (3 5 4 3)
          );
   }


);

mergePatchPairs
(
);

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

3.4 ControlDict symmetry BC and wedge angle =10

/*--------------------------------*- 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         10;

deltaT          0.001;

writeControl    timeStep;

writeInterval   20;

purgeWrite      0;

writeFormat     ascii;

writePrecision  6;

writeCompression off;

timeFormat      general;

timePrecision   6;

runTimeModifiable true;


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

3.5 Pressure 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
{
    axis
    {
        type            empty;
    }
    fixedWalls
    {
        type            zeroGradient;
    }
    front
    {
        type            symmetry;
    }
    back
    {
        type            symmetry;
    }
    inlet
    {
        type            zeroGradient;
    }
    outlet
    {
        type            fixedValue;
        value           uniform 0.0189;
    }
 
    
}

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

3.6 Velocity 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.0937312 0 0);

boundaryField
{
    inlet
    {
        type            fixedValue;
        value           uniform (0.0937312 0 0);
    }
    fixedWalls
    {
        type            fixedValue;
        value           uniform (0 0 0);
    }
    front
    {
        type            symmetry;
    }
    back
    {
        type            symmetry;
    }
    axis
    {
        type            empty;
    }
    outlet
    {
        type            zeroGradient;
    }
}

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

3.6 Transport Properties 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       dictionary;
    location    "constant";
    object      transportProperties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

nu              [0 2 -1 0 0 0 0] 8.9268e-07;


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

  1. RESULTS AND DISCUSSION

4.1 Results for Wedge Angle 10 degree – Wedge BC

Inference:

  1. b) The velocity profile at the entrance has three distinct regions:

Viscous Sub layer: The effect of velocity boundary layer is evident in this region

Transitional Buffer layer: It gradually decreases as the viscous sub layer merges with the turbulent core with increasing distance from the entrance.

Turbulent Core: Velocity profile is almost flat and comprised of the turbulent components of the flow.

  1. c) It is envisaged from the results that velocity profile almost resembles a fully developed laminar flow with only meagre presence of turbulent component at the pipe centre line
  2. d) The velocity profile shows a fully developed laminar flow.
  3. e) The shear stress plot matches with the conventional image in literature having zero at centre line and maximum at the walls due to boundary layer effect.

4.2 Results for Wedge Angle 10 degree – Symmetry BC

Inference:

b) The velocity profile at the entrance has three distinct regions:

Viscous Sub layer: The effect of velocity boundary layer is evident in this region

Transitional Buffer layer: It gradually decreases as the viscous sub layer merges with the turbulent core with increasing distance from the entrance.

Turbulent Core: Velocity profile is almost flat and comprised of the turbulent components of the flow.

c) It is envisaged from the results that velocity profile almost resembles a fully developed laminar flow with only meagre presence of turbulent component at the pipe centre line

d) The velocity profile shows a fully developed laminar flow.

e) The shear stress plot matches with the conventional image in literature having zero at centre line and maximum at the walls due to boundary layer effect.

4.3 Table depicts the Comparison of Results Wedge and Symmetry BC

Wedge Angle (Degrees)

Hagen Poiseuille Maximum Velocity (m/s)

Wedge Boundary Condition

Symmetry Boundary Condition

Computational Time (sec)

Maximum Velocity (m/s)

Courant Number Max

Computational Time (sec)

Maximum Velocity (m/s)

Courant Number Max

10

0.1875

464.49

0.1758

0.00926297

582.29

0.1756

0.00926563

25

0.1875

497.72

0.1767

0.00935719

642.68

0.1768

0.00935755

45

0.1875

507.01

0.1791

0.00964046

693.81

0.1795

0.00964062

Inference:

  • Symmetry boundary condition results are closer to the analytical solution and thus it is slightly better when compared to Wedge BC.
  • CFL number increases with the increase in wedge angle
  • Simulation Execution time is larger for Symmetry BC in comparison to Wedge BC.
  • As Wedge Angle increases the results are closer to the analytical solution ie., Accuracy increases with increase in wedge angle
  1. CONCLUSION

The program to compare Wedge and symmetry BC of the laminar flow through pipe using OpenFoam & Matlab was developed and executed successfully.

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.