Week 11 - Simulation of Flow through a pipe in OpenFoam
Objective : To do the simulation of flow through a pipe where using the wedge boundary conditions and to make a study over the 'Hagen poiseuille's equation'.
Hydrodynamic entrance length:
The hydrodynamic entrance region refers to the area of a pipe where fluid entering a pipe develops a velocity profile due to viscous forces spreading from a pipe's interior wall. This region has a non-uniform flow. The fluid enters a pipe at a constant speed, so fluid particles in the layer in contact with the pipe's surface come to a complete stop because of the no-slip rule. The layer in contact with the pipe surface is resistant to the acceleration of neighboring layers due to viscous forces within the fluid and slowly slows down adjacent layers of fluid, creating a velocity profile. To keep the mass true, the speed of the fluid layers in the center of the pipe increases to compensate for the reduced speed of the fluid layers near the surface of the pipe. It creates a gradient of velocity across the pipe cross-section.
Entry length for Laminar flow is `L_e=0.06*Re*D`
Hagen poiseuille's equation :
In non ideal fluid dynamics, Hagen-poiseuille equation is a physical law that gives pressure drop in an incompressible and newtonian fluid in a laminar flow through a long cylinder pipe of constant cross section.

THEORY :
Assumptions :
- The flow is incompressible, laminar ,steady and axisymmetric.
- The working fluid is Newtonian (water).
- The pipe-diameter is constant.
Water flows through a constant-diameter circular pipe of length 2 m and diameter 10 mm. It enters the pipe with Reynold's number, Re = 2100. Due to the no-slip condition at the walls of the pipe, a hydrodynamic boundary layer begins to develop near the walls.The hydrodynamic entrance region refers to the area of a pipe where fluid entering a pipe develops a velocity profile due to viscous forces spreading from a pipe's interior wall. This region has a non-uniform flow. The fluid enters a pipe at a constant speed, so fluid particles in the layer in contact with the pipe's surface come to a complete stop because of the no-slip rule. The layer in contact with the pipe surface is resistant to the acceleration of neighboring layers due to viscous forces within the fluid and slowly slows down adjacent layers of fluid, creating a velocity profile. To keep the mass true, the speed of the fluid layers in the center of the pipe increases to compensate for the reduced speed of the fluid layers near the surface of the pipe. It creates a gradient of velocity across the pipe cross-section.

In non ideal fluid dynamics, Hagen-poiseuille equation(analytical solution) is a physical law that gives pressure drop in an incompressible and newtonian fluid in a laminar flow through a long cylinder pipe of constant cross section, but however, it applies only to the fully developed flow (i.e., x>Le">x > `L_e`) .
THE WEDGE SECTION:

MATLAB CODE:
clear all
close all
clc
% The Reynolds number
Re = 2100
% The Diameter of pipe
D = 0.01
% The Entry length of pipe
L =0.05*Re*D
% The total length of pipe
L_p = L+0.95
% Radius of pipe
r = D/2
% Angle of wedge (less than 5)
theta = 2
% Dynamic viscosity
mu = 0.89e-3
% Density
rho = 997
% Hagen poiseuille flow equations
v_avg = (Re*mu)/(rho*D)
v_max = 2*(v_avg)
del_p = (8*mu*L_p*v_avg)/(r^2)
kinematic_p = del_p/rho
header1='/*----------------------------------*- C++ -*------------------------------------*'
header2=' ========= | '
header3=' \ / F ield | OpenFOAM: The Open Source CFD Toolbox '
header4=' \ / O peration | Website: https://openfoam.org '
header5=' \ / A nd | Version: 8 '
header6=' \/ M anipulation | '
header7='*----------------------------------------------------------------------------------*/'
f1 = fopen('blockMeshDict','w')
fprintf(f1, '%sn', header1)
fprintf(f1, '%sn', header2)
fprintf(f1, '%sn', header3)
fprintf(f1, '%sn', header4)
fprintf(f1, '%sn', header5)
fprintf(f1, '%sn', header6)
fprintf(f1, '%sn', header7)
fprintf(f1, 'FoamFilen{n')
fprintf(f1, '%12s t 2.0;n','version')
fprintf(f1, '%12s t ascii;n','format')
fprintf(f1, '%12s t dictionary;n','class')
fprintf(f1, '%12s t blockMeshDict;n','object')
fprintf(f1, '}n')
header8='// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //'
fprintf(f1, '%snn',header8)
fprintf(f1, 'convertToMeters 1;nn')
fprintf(f1, 'verticesn')
fprintf(f1, '(n')
fprintf(f1, 't (%d %d %d)n',0,0,0)
fprintf(f1, 't (%d %d %d)n',0,r*cosd(theta/2),r*sind(theta/2))
fprintf(f1, 't (%d %d %d)n',0,r*cosd(theta/2),-r*sind(theta/2))
fprintf(f1, 't (%d %d %d)n',L_p,0,0)
fprintf(f1, 't (%d %d %d)n',L_p,r*cosd(theta/2),r*sind(theta/2))
fprintf(f1, 't (%d %d %d)n',L_p,r*cosd(theta/2),-r*sind(theta/2))
fprintf(f1, ');nn')
fprintf(f1, 'blocksn')
fprintf(f1, '(n')
fprintf(f1, 't hex (0 3 5 2 0 3 4 1) (%d %d %d)',300,30,1)
fprintf(f1, ' SimpleGrading (1 0.1 1)nn')
fprintf(f1, ');nn')
fprintf(f1, 'edgesn')
fprintf(f1, '(n')
fprintf(f1, 't arc %d %d (%d %d %d)n',1,2,0,r,0)
fprintf(f1, 't arc %d %d (%d %d %d)n',4,5,L_p,r,0)
fprintf(f1, ');nn')
fprintf(f1, 'boundaryn')
fprintf(f1, '(n')
fprintf(f1, 't inlet n')
fprintf(f1, 't {n')
fprintf(f1, 'tt type patch;n')
fprintf(f1, 'tt facesn')
fprintf(f1, 'tt (n')
fprintf(f1, 'ttt (0 1 2 0)n')
fprintf(f1, 'tt );n')
fprintf(f1, 't }n')
fprintf(f1, 't outletn')
fprintf(f1, 't {n')
fprintf(f1, 'tt type patch;n')
fprintf(f1, 'tt facesn')
fprintf(f1, 'tt (n')
fprintf(f1, 'ttt (3 5 4 3)n')
fprintf(f1, 'tt );n')
fprintf(f1, 't }n')
fprintf(f1, 't topn')
fprintf(f1, 't {n')
fprintf(f1, 'tt type wall;n')
fprintf(f1, 'tt facesn')
fprintf(f1, 'tt (n')
fprintf(f1, 'ttt (1 4 5 2)n')
fprintf(f1, 'tt );n')
fprintf(f1, 't }n')
fprintf(f1, 't frontn')
fprintf(f1, 't {n')
fprintf(f1, 'tt type wedge;n')
fprintf(f1, 'tt facesn')
fprintf(f1, 'tt (n')
fprintf(f1, 'ttt (0 3 4 1)n')
fprintf(f1, 'tt );n')
fprintf(f1, 't }n')
fprintf(f1, 't backn')
fprintf(f1, 't {n')
fprintf(f1, 'tt type wedge;n')
fprintf(f1, 'tt facesn')
fprintf(f1, 'tt (n')
fprintf(f1, 'ttt (0 2 5 3)n')
fprintf(f1, 'tt );n')
fprintf(f1, 't }n')
fprintf(f1, 't axisn')
fprintf(f1, 't {n')
fprintf(f1, 'tt type empty;n')
fprintf(f1, 'tt facesn')
fprintf(f1, 'tt (n')
fprintf(f1, 'ttt (0 3 3 0)n')
fprintf(f1, 'tt );n')
fprintf(f1, 't }n')
fprintf(f1, '); nn')
fprintf(f1, 'mergePatchPairsn')
fprintf(f1, '(n')
fprintf(f1, ');n')
header9='// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //'
fprintf(f1, '%sn',header9)
fclose(f1)
BLOCKMESH DICT:
/*----------------------------------*- C++ -*------------------------------------*
========= |
\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\ / O peration | Website: https://openfoam.org
\ / A nd | Version: 8
\/ M anipulation |
*----------------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
convertToMeters 1;
vertices
(
(0 0 0)
(0 4.999238e-03 8.726203e-05)
(0 4.999238e-03 -8.726203e-05)
(2 0 0)
(2 4.999238e-03 8.726203e-05)
(2 4.999238e-03 -8.726203e-05)
);
blocks
(
hex (0 3 5 2 0 3 4 1) (300 30 1) SimpleGrading (1 0.1 1)
);
edges
(
arc 1 2 (0 5.000000e-03 0)
arc 4 5 (2 5.000000e-03 0)
);
boundary
(
inlet
{
type patch;
faces
(
(0 1 2 0)
);
}
outlet
{
type patch;
faces
(
(3 5 4 3)
);
}
top
{
type wall;
faces
(
(1 4 5 2)
);
}
front
{
type wedge;
faces
(
(0 3 4 1)
);
}
back
{
type wedge;
faces
(
(0 2 5 3)
);
}
axis
{
type empty;
faces
(
(0 3 3 0)
);
}
);
mergePatchPairs
(
);
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
CONTROL DICT:
/*--------------------------------*- C++ -*----------------------------------*
========= |
\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\ / O peration | Website: https://openfoam.org
\ / A nd | Version: 8
\/ M anipulation |
*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "system";
object controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
application icoFoam;
startFrom startTime;
startTime 0;
stopAt endTime;
endTime 8;
deltaT 4.4469e-04;
writeControl timeStep;
writeInterval 20;
purgeWrite 0;
writeFormat ascii;
writePrecision 6;
writeCompression off;
timeFormat general;
timePrecision 6;
runTimeModifiable true;
// ************************************************************************* //
BOUNDARY CONDITIONS:
VELOCITY:
/*--------------------------------*- C++ -*----------------------------------*
========= |
\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\ / O peration | Website: https://openfoam.org
\ / A nd | Version: 8
\/ M anipulation |
*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volVectorField;
object U;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dimensions [0 1 -1 0 0 0 0];
internalField uniform (0.1874 0 0);
boundaryField
{
inlet
{
type fixedValue;
value uniform (0.1874 0 0);
}
outlet
{
type zeroGradient;
}
axis
{
type empty;
}
top
{
type fixedValue;
value uniform (0 0 0);
}
front
{
type wedge;
}
back
{
type wedge;
}
}
// ************************************************************************* //
PRESSURE:
/*--------------------------------*- C++ -*----------------------------------*
========= |
\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\ / O peration | Website: https://openfoam.org
\ / A nd | Version: 8
\/ 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.10709;
}
axis
{
type empty;
}
top
{
type zeroGradient;
}
front
{
type wedge;
}
back
{
type wedge;
}
}
// ************************************************************************* //
TRANSPORT PROPERTIES:
/*--------------------------------*- C++ -*----------------------------------*
========= |
\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\ / O peration | Website: https://openfoam.org
\ / A nd | Version: 8
\/ M anipulation |
*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "constant";
object transportProperties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
nu [0 2 -1 0 0 0 0] 8.9267e-07;
// ************************************************************************* //
RESULTS:
MESHING:

VELOCITY PROFILE:





PRESSURE PROFILES:

SHEAR STRESS PROFILE:


CONCLUSION:
Comparing plot results with Hagen poiseuille equation for different variables, we get
|
Variables |
Inferred from plot |
Hagen-poiseuille equation |
|
Maximum Velocity |
0.3727 |
0.3749 |
|
Kinematic Pressure |
0.10706 |
0.10709 |
|
Max. Shear Stress |
0.13058 |
0.13342 |
We can that this small variation in values of variables(error) is due to the fact that while one is theoretical , the other is obtained by doing simulation(added to the fact that this particular case is transient state, so the value we get from it is not the perfect answer , but very close enough to it for us to use it as its substitute).