Week 11 - Simulation of Flow through a pipe in OpenFoam
Considering a cylindrical pipe with the working fluid water entering the pipe with laminar flow. The aim is to simulate the flow through a section of pipe and compare the flow with the Hagen Poiseuille’s equation.
Hagen Poiseuuile’s equations
In non ideal fluid dynamics, the Hagen-Poiseuille equation, also known as the hagen-poiseuille law, is a physical law that gives the pressure drop in an incompressible and Newtonian fluid in laminar flow flowing through a long cylindrical pipe of constant cross-section.
For a pipe with a laminar flow there are basically two regions.Entrance lenght developing region where the flow is not fully developed 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 fluid layers in the center of the pipe increases to compensate for reduced speed of the fluid layers ner the surface of the pipe.It creates a gradient of velocity across the pipe cross-section.
Entry length for laminar flow is
L=0.06*Ry*D;
let
Diameter,D=0.004;
Radius,r=0.002;
Radial Disatance, R=0 to r;
Dynamic viscosity,mu=0.00089N/m^2;
Density,rho=997kg/m^3;
Reylonds number,Re=2100;
Formulaes
BlockMeshDict
/*--------------------------------*- 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 1.998782e-03 6.979899e-05)
(0 1.998782e-03 -6.979899e-05)
(1.204000e+00 0 0)
(1.204000e+00 1.998782e-03 6.979899e-05)
(1.204000e+00 1.998782e-03 -6.979899e-05)
);
blocks
(
hex (0 3 5 2 0 3 4 1) (300 30 1)SimpleGrading (1 0.1 1)
);
edges
(
arc 1 2 (0 2.000000e-03 0)
arc 4 5 (1.204000e+00 2.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
(
);
// ************************************************************************* //ControlDict
/*--------------------------------*- 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 1;
deltaT 0.001;
writeControl timeStep;
writeInterval 20;
purgeWrite 0;
writeFormat ascii;
writePrecision 6;
writeCompression off;
timeFormat general;
timePrecision 6;
runTimeModifiable true;
// ************************************************************************* //
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.4686 0 0);
boundaryField
{
inlet
{
type fixedValue;
value uniform (0.4686 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 1.0073;
}
axis
{
type empty;
}
top
{
type zeroGradient;
}
front
{
type wedge;
}
back
{
type wedge;
}
}
// ************************************************************************* //
Dynamic viscoscity
/*--------------------------------*- 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] 0.892e-6;
// ************************************************************************* //
Meshing of the sector of pipe
Matlab code
clear all
close all
clc
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,'%s\n',header1);
fprintf(f1,'%s\n',header2);
fprintf(f1,'%s\n',header3);
fprintf(f1,'%s\n',header4);
fprintf(f1,'%s\n',header5);
fprintf(f1,'%s\n',header6);
fprintf(f1,'%s\n',header7);
fprintf(f1,'FoamFile\n{\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');
fprintf(f1,'%s\n\n','// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //')
fprintf(f1,'convertToMeters 1;\n\n')
fprintf(f1,'vertices\n')
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,');\n\n');
fprintf(f1,'blocks\n');
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)\n\n');
fprintf(f1,');\n\n');
fprintf(f1,'edges\n');
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,');\n\n');
fprintf(f1,'boundary\n');
fprintf(f1,'(\n');
fprintf(f1,'\t inlet \n');
fprintf(f1,'\t {\n');
fprintf(f1,'\t\t type patch;\n');
fprintf(f1,'\t\t faces\n');
fprintf(f1,'\t\t (\n');
fprintf(f1,'\t\t\t (0 1 2 0)\n');
fprintf(f1,'\t\t );\n');
fprintf(f1,'\t }\n');
fprintf(f1,'\t outlet\n');
fprintf(f1,'\t {\n');
fprintf(f1,'\t\t type patch;\n');
fprintf(f1,'\t\t faces\n');
fprintf(f1,'\t\t(\n');
fprintf(f1,'\t\t\t (3 5 4 3)\n');
fprintf(f1,'\t\t );\n');
fprintf(f1,'\t }\n');
fprintf(f1,'\t top\n');
fprintf(f1,'\t {\n');
fprintf(f1,'\t\t type wall;\n');
fprintf(f1,'\t\t faces\n');
fprintf(f1,'\t\t(\n');
fprintf(f1,'\t\t\t (1 4 5 2)\n');
fprintf(f1,'\t\t );\n');
fprintf(f1,'\t }\n');
fprintf(f1,'\t front\n');
fprintf(f1,'\t {\n');
fprintf(f1,'\t\t type wedge;\n');
fprintf(f1,'\t\t faces\n');
fprintf(f1,'\t\t(\n');
fprintf(f1,'\t\t\t (0 3 4 1)\n');
fprintf(f1,'\t\t );\n');
fprintf(f1,'\t }\n');
fprintf(f1,'\t back\n');
fprintf(f1,'\t {\n');
fprintf(f1,'\t\t type wedge;\n');
fprintf(f1,'\t\t faces\n');
fprintf(f1,'\t\t(\n');
fprintf(f1,'\t\t\t (0 2 5 3)\n');
fprintf(f1,'\t\t );\n');
fprintf(f1,'\t }\n');
fprintf(f1,'\t axis\n');
fprintf(f1,'\t {\n');
fprintf(f1,'\t\t type empty;\n');
fprintf(f1,'\t\t faces\n');
fprintf(f1,'\t\t(\n');
fprintf(f1,'\t\t\t (0 3 3 0)\n');
fprintf(f1,'\t\t );\n');
fprintf(f1,'\t }\n');
fprintf(f1,'mergePatchPairs\n');
fprintf(f1,'\t\t );\n');
fprintf(f1,'\t }\n');
fprintf(f1,'%s\n','// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //');
fclose(f1);
Velocity profile
2.Velocity profile at 0.35m
3.Velocity profile at 0.5m
At the length of 0.504,the flow gets fully developed because of maximum velocity obtained here.And it cen be seen in the plot for the velocity profile at the length 0.504
For understanding below velocity profile plots for 0.7 and 1.2m are shown below
4.Velocity profile at 0.7
5.Velocity profile at 1.2m
shear stress code:
clear all
close all
clc
D=0.004;
R=D/2;
r=linspace(-0.01,0.01,1000);
Re=2100;
mu=8.9e-4;
rho=997;
V=(Re*mu)/(rho*D);
u_max=2*V;
Tau=zeros(length(r));
for i=1:length(r)
Tau(i)=(2*mu*u_max*abs(r(i)))/(R^2);
end
figure(1)
plot(t,Tau);
xlabel('Radial distance in m');
yalel('Shear strees(Tau)in N/m^2');
grid minor
If we analyse the graph shear stress near the walls is maximum and at the center of the pipe in mimimum(i.e 0).