Week 11 - Simulation of Flow through a pipe in OpenFoam
Abstract
In this project, the laminar flow of an incompressible fluid through a pipe is simulated using OpenFOAM. The pressure-drop across the pipe length, the hydrodynamic entry length and the velocity profile are then compared to the analytical values obtained using the Hagen-Poiseuille’s equation.
Introduction
The velocity of a fluid entering a circular pipe at a uniform rate starts decreasing along the boundary layer region due to ‘no-slip’ condition which in turn increases the velocity at the mid-section of the pipe. Thus, the boundary layer thickness increases gradually in the flow direction. When the boundary layer thickness reaches the centre of the pipe, the flow is said to be fully developed and the velocity profile doesn’t change anymore. The region from the inlet of the pipe to this point is known as the hydrodynamic entry length.
Fig 1: The development of velocity boundary layer in a pipe [1]
In laminar flow, the non-dimensional hydrodynamic entry length is given approximately as [1]:
The maximum velocity which occurs at the centreline of the pipe can be obtained by substituting r=0 in the above equation.
The Hagen–Poiseuille equation gives the pressure-drop of an incompressible and Newtonian fluid in laminar flow which is flowing through a long cylindrical pipe of constant cross section. This pressure-drop is attributed to the viscous effects only.
Methodology:
Simulation using OpenFOAM:
The laminar incompressible flow of water is simulated in OpenFOAM by modelling a wedge of angle 4 degrees. In this way the axi-symmetric flow is modelled using wedge boundary condition.
The blockMeshDict file is created in MATLAB by file parsing as follows:
Matlab Code for generation of blockmeshDict file:
close all
clear all
clc
Re = 2100;
D = 0.0254; % considering 1" pipe dia
R = D/2;
Le = 0.05*Re*D; % entrance length
L = round(Le); % Total length
theta = 4; % Assumed wedge angle
nu = 8e-07; % Kinematic viscosity of water at 25degC
% 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 0 0];
v4 = [L R*cosd(theta/2) -R*sind(theta/2)];
v5 = [L R*cosd(theta/2) R*sind(theta/2)];
% generating blockmesh file
f1 = fopen('blockMeshDict.txt','w');
fprintf(f1,'%s \n','/*--------------------------------*- C++ -*-------------------------------*\');
fprintf(f1,'%s \n',' ========= |');
fprintf(f1,'%s \n',' \\ / F ield | OpenFOAM: The Open Source CFD Toolbox');
fprintf(f1,'%s \n',' \\ / O peration | Website: https://openfoam.org');
fprintf(f1,'%s \n',' \\ / A nd | Version: 8');
fprintf(f1,'%s \n',' \\/ M anipulation |');
fprintf(f1,'%s \r\n', '\*------------------------------------------------------------------------*/');
fprintf(f1, ' \r\n');
fprintf(f1, 'FoamFile \r\n{ \r\n');
fprintf(f1, '%13s \t 2.0; \r\n', 'version');
fprintf(f1, '%12s \t ascii; \r\n', 'format');
fprintf(f1, '%11s \t directory; \r\n', 'class');
fprintf(f1, '%12s \t blockMeshDict;\r\n}\r\n', 'object');
fprintf(f1, '%s \r\n', '// *******************************************//');
fprintf(f1,'convertToMeters 1; \r\n\r\n');
fprintf(f1,'vertices\r\n(\r\n');
fprintf(f1, '\t(%d %d %d)\r\n', v0(1),v0(2),v0(3));
fprintf(f1, '\t(%d %d %d)\r\n', v1(1),v1(2),v1(3));
fprintf(f1, '\t(%d %d %d)\r\n', v2(1),v2(2),v2(3));
fprintf(f1, '\t(%d %d %d)\r\n', v3(1),v3(2),v3(3));
fprintf(f1, '\t(%d %d %d)\r\n', v4(1),v4(2),v4(3));
fprintf(f1, '\t(%d %d %d)\r\n', v5(1),v5(2),v5(3));
fprintf(f1, ');\r\n\r\n');
fprintf(f1,'blocks\r\n(\r\n');
fprintf(f1,'\thex (0 3 4 1 0 3 5 2) (1000 50 1)');
fprintf(f1,' simpleGrading (%d %d %d)\r\n', 1,1,1);
fprintf(f1,');\r\n');
fprintf(f1,'edges\r\n\t(\r\n');
fprintf(f1,'\t\t arc %d %d (%d %d %d)\r\n', 1,2,0,R,0);
fprintf(f1,'\t\t arc %d %d (%d %d %d)\r\n', 4,5,L,R,0);
fprintf(f1,'\t);\r\n');
fprintf(f1,'boundary\r\n(\r\n');
%assign boundary conditions
fprintf(f1,'\tinlet\r\n\t{\r\n');
fprintf(f1,'\t\ttype patch;\r\n');
fprintf(f1,'\t\tfaces\r\n\t\t(\r\n');
fprintf(f1,'\t\t\t(0 2 1 0)\r\n');
fprintf(f1,'\t\t);\r\n\t}\r\n');
fprintf(f1,'\toutlet\r\n\t{\r\n');
fprintf(f1,'\t\ttype patch;\r\n');
fprintf(f1,'\t\tfaces\r\n\t\t(\r\n');
fprintf(f1,'\t\t\t(3 4 5 3)\r\n');
fprintf(f1,'\t\t);\r\n\t}\r\n');
fprintf(f1,'\tfront\r\n\t{\r\n');
fprintf(f1,'\t\ttype wedge;\r\n');
fprintf(f1,'\t\tfaces\r\n\t\t(\r\n');
fprintf(f1,'\t\t\t(0 3 5 2)\r\n');
fprintf(f1,'\t\t);\r\n\t}\r\n');
fprintf(f1,'\tback\r\n\t{\r\n');
fprintf(f1,'\t\ttype wedge;\r\n');
fprintf(f1,'\t\tfaces\r\n\t\t(\r\n');
fprintf(f1,'\t\t\t(0 1 4 3)\r\n');
fprintf(f1,'\t\t);\r\n\t}\r\n');
fprintf(f1,'\twall\r\n\t{\r\n');
fprintf(f1,'\t\ttype wall;\r\n');
fprintf(f1,'\t\tfaces\r\n\t\t(\r\n');
fprintf(f1,'\t\t\t(1 2 5 4)\r\n');
fprintf(f1,'\t\t);\r\n\t}\r\n');
fprintf(f1,'\taxis\r\n\t{\r\n');
fprintf(f1,'\t\ttype empty;\r\n');
fprintf(f1,'\t\tfaces\r\n\t\t(\r\n');
fprintf(f1,'\t\t\t(0 3 3 0)\r\n');
fprintf(f1,'\t\t);\r\n\t}\r\n');
fprintf(f1, ');\r\n\r\n');
fprintf(f1, 'mergePatchPairs\r\n\r\n(\r\n);');blockmeshDict file:
/*--------------------------------*- 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 directory;
object blockMeshDict;
}
// *******************************************//
convertToMeters 1;
vertices
(
(0 0 0)
(0 1.269226e-02 -4.432236e-04)
(0 1.269226e-02 4.432236e-04)
(3 0 0)
(3 1.269226e-02 -4.432236e-04)
(3 1.269226e-02 4.432236e-04)
);
blocks
(
hex (0 3 4 1 0 3 5 2) (1000 50 1) simpleGrading (1 1 1)
);
edges
(
arc 1 2 (0 1.270000e-02 0)
arc 4 5 (3 1.270000e-02 0)
);
boundary
(
inlet
{
type patch;
faces
(
(0 2 1 0)
);
}
outlet
{
type patch;
faces
(
(3 4 5 3)
);
}
front
{
type wedge;
faces
(
(0 3 5 2)
);
}
back
{
type wedge;
faces
(
(0 1 4 3)
);
}
wall
{
type wall;
faces
(
(1 2 5 4)
);
}
axis
{
type empty;
faces
(
(0 3 3 0)
);
}
);
mergePatchPairs
(
);ControlDict file:
/*--------------------------------*- 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 20;
deltaT 0.01;
writeControl timeStep;
writeInterval 200;
purgeWrite 0;
writeFormat ascii;
writePrecision 4;
writeCompression off;
timeFormat general;
timePrecision 6;
runTimeModifiable true;
// ************************************************************************* //
pressure input file:
/*--------------------------------*- 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;
}
front
{
type wedge;
}
back
{
type wedge;
}
wall
{
type zeroGradient;
}
axis
{
type empty;
}
}
// ************************************************************************* //
Velocity input file:
/*--------------------------------*- 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 0 0);
boundaryField
{
inlet
{
type fixedValue;
value uniform (0.073804 0 0);
}
outlet
{
type zeroGradient;
}
front
{
type wedge;
}
back
{
type wedge;
}
wall
{
type noSlip;
}
axis
{
type empty;
}
}
// ************************************************************************* //
Mesh generation:
Velocity distribution:
Pressure distribution:
Results:
1. The below plot shows the velocity profile at different distances from the inlet. It is observed that the velocity profile is fully developed at a distance of 2.667m from the inlet in agreement with the theoretical hydrodynamic entrance lengths.
The plot shown below for the theoretical velocity profile for the fully developed region is matching with the numerical results in the above plot. The maximum velocity obtained theoretically is 0.14761m/s which is matching in both the calculations as well.
2. The plot for the numerically computed pressure loss for the fully developed region is shown below. The theoretical value of the pressure loss is 0.0011 m^2/s^2 which is matching to the numerically obtained value as can be observed from below plot.
3. The below plot shows the shear stress obtained by post processing the velocity data in paraview. The maximum shear stress from the plot is 0.0208 Pa whereas the max shear stress calculated theoretically is 0.02069 Pa.
Conclusion:
The numerical results obtained using icoFoam solver in OpenFOAM are in good agreement with the theoretical results obtained using Hagen-Poiseuilli's equations.
References:
[1] Cengel, Yunus A. and Cimbala. Fluid Mechanics: Fundamentals and Applications, 3rd Edition.