Simulation Of Flow Through a Pipe-OpenFoam
PROJECT
INTRODUCTIONS
The purpose of this project is to simulate the laminar incompressible flow through a pipe. The project is divided into two part, in the first part the objective is to make simulations. With the given Reynold’s number and pipe Diameter calculation of inner velocity and the entrance length of the pipe has to be done.
We are considering symmetry conditions for this problem else it would be computationally expensive and minimum wedge angle for proper simulations of flow.
Assume that the loud liquid exhibits laminar flow in a pipe. In the pipe there is bunch of circulation layers (lamina) of liquid each having a velocity determined by the radial distance from the center of the tube .
Fluid velocity in a pipe change from zero at the surface because of the no slip boundary condititons to a maximum at the center.
For laminar flows the Re <= 2100.
Entrance length
The entrance regions refers to that portions of a pipe untill the velocity profile is fully developed.
Hagen Posiseuille Flow ;
It applies to laminar flow of an incompressible fluid of viscosity through a tube of length and radius .
Flow rate is directly proportional to the pressure difference, and inverserly proportional to the length of the tube and viscosity of the fluid.
Geometry and Meshing ;
This project requires the use of MatLab to code and generate the blockMeshDict file for openFoam.
Matlab script :
Matlab code to automatically generate blockMeshDict file for any given wedge angle
wedge angle 4
clear all
close all
clc
% inputs
theta = 4; % angle of wedge
Re = 2100;
D=0.02;
L = 0.06*Re*D; % lenght of the pipe
R = D/2;
mu = 0.89e-3; % viscosity of water at 25 deg C
rho = 997.0; % Density of water at 25 deg C
nu = mu/rho; % kinematic Viscosity of water at 25 deg C
% Calculatic of pressure and velocity :
v_avg = Re*mu/(D*rho); % Average Velocity
V_max = 2*v_avg; % Max velocity
del_P = (32*mu*v_avg*L)/(D^2); % pressure Drop
kin_P = del_P/rho; % kinematic Presure Drop
n1 = 50 % Number of cells along the pipe
n2 = 1 % Number of cells radially
n3 = 300 % number of cells on the curve
gfactor = 0.2;
% vertices of the wedge
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)];
% printing the file
text1 = \'/*-------------------------------*-C++ -*--------------------------------*\';
text2 = \' ========= \';
text3 = \' \\ / F ield [ OpenFOAM: The open Source CFD Toolbox \';
text4 = \' \\ / O perations [ version : 4.1 \';
text5 = \' \\ / A nd [ Web : www.openFOAM.org \';
text6 = \' \\/ M anipulation \';
text7 = \'*-------------------------------------------------------------------------\';
text8 = \'// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *//\';
% creating file
f1 = fopen(\'blockMeshDict.txt\',\'wt\');
fprintf(f1,\'%sn\',text1);
fprintf(f1,\'%sn\',text2);
fprintf(f1,\'%sn\',text3);
fprintf(f1,\'%sn\',text4);
fprintf(f1,\'%sn\',text5);
fprintf(f1,\'%sn\',text6);
fprintf(f1,\'%sn\',text7);
fprintf(f1,\'FoamFilen{n\');
fprintf(f1,\'%16s t 2.0;n\',\'versions\');
fprintf(f1,\'%16s t ascii;n\',\'format\');
fprintf(f1,\'%16s t dictionary;n\',\'class\');
fprintf(f1,\'%16s t blockMeshDict;n}n\',\'object\');
fprintf(f1,\'%snn\', text8);
fprintf(f1,\'convertToMeters 1;nn\');
fprintf(f1,\'verticesn(n\');
fprintf(f1,\'t(%d %d %d)n\',v0(1),v0(2),v0(3));
fprintf(f1,\'t(%d %d %d)n\',v1(1),v1(2),v1(3));
fprintf(f1,\'t(%d %d %d)n\',v2(1),v2(2),v2(3));
fprintf(f1,\'t(%d %d %d)n\',v3(1),v3(2),v3(3));
fprintf(f1,\'t(%d %d %d)n\',v4(1),v4(2),v4(3));
fprintf(f1,\'t(%d %d %d)n);nn\',v5(1),v5(2),v5(3));
fprintf(f1,\'blocksn(n\');
fprintf(f1,\'t hex (0 1 2 0 3 4 5 3)(%d %d %d)\',n1,n2,n3);
fprintf(f1,\' simpleGrading (%0.1f 1 1)n);nn\',gfactor);
fprintf(f1,\'edgesn(ntarc %d %d (%d %0.2f %d)ntarc %d %d (%0.2f %0.2f %d)n);nn\', 1, 2, 0, R, 0, 4, 5, L, R, 0);
fprintf(f1,\'boundaryn(n\');
fprintf(f1,\'t%snt{ntt type %s,ntt Faces ntt (ntt (%d %d %d %d) ntt );nt}nn\',\'axis\',\'empty\',0,3,3,0);
fprintf(f1,\'t%snt{ntt type %s,ntt Faces ntt (ntt (%d %d %d %d) ntt );nt}nn\',\'inlet\',\'patch\',0,1,2,0);
fprintf(f1,\'t%snt{ntt type %s,ntt Faces ntt (ntt (%d %d %d %d) ntt );nt}nn\',\'pipe\',\'wall\',1,4,5,2);
fprintf(f1,\'t%snt{ntt type %s,ntt Faces ntt (ntt (%d %d %d %d) ntt );nt}nn\',\'outlet\',\'patch\',3,4,5,3);
fprintf(f1,\'t%snt{ntt type %s,ntt Faces ntt (ntt (%d %d %d %d) ntt );nt}nn\',\'front\',\'wedge\',0,3,5,2);
fprintf(f1,\'t%snt{ntt type %s,ntt Faces ntt (ntt (%d %d %d %d) ntt );nt}nn\',\'back\',\'wedge\',0,1,4,3);
fprintf(f1,\');nn\');
fprintf(f1,\'mergePatchPairsn(n);nn\');
fprintf(f1,\'%s\',text8);
fclose(f1);main Geomatric type :
Patch: generic type containing no geometric or topological information about the Mesh, eg. Used for an inlet or an outlet.
Wall: for patch that coincides with a solid wall, required for some physical modelling, e.g . wall functions in turbulance modelling.
Empty: dimension for which no solution is required.
Wedge : for 2 dimensions axi-symmetric cases, e.g. a cylinder, the geometry is specfied as a wedge of small angle ( e.g .1) and 1 cell thick, running along the centre line, stadding one of the co-ordinate planes as shown in figure the axi-symmetric wedge planes must be specified as seperated patches of wedge type,
BlockMeshDict file generated for wedge angle 4
/*--------------------------------*- 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 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
(
axis
{
type empty;
faces
(
(0 3 3 0)
);
}
inlet
{
type patch;
faces
(
(0 1 2 0)
);
}
top
{
type wall;
faces
(
(1 4 5 2)
);
}
outlet
{
type patch;
faces
(
(3 5 4 3)
);
}
front
{
type wedge;
faces
(
(0 3 4 1)
);
}
back
{
type wedge;
faces
(
(0 2 5 3)
);
}
);
mergePatchPairs
(
);
// ************************************************************************* //
pressure Boundary condition
/*--------------------------------*- 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;
}
top
{
type zeroGradient;
}
inlet
{
type zeroGradient;
}
outlet
{
type fixedValue;
value uniform 1.0073;
}
front
{
type wedge;
}
back
{
type wedge;
}
}
// ************************************************************************* //
velocity boundary condition
/*--------------------------------*- 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.4686 0 0);
boundaryField
{
axis
{
type empty;
}
inlet
{
type fixedValue;
value uniform (0.04686 0 0);
}
top
{
type fixedValue;
value uniform (0 0 0);
}
outlet
{
type zeroGradient;
}
front
{
type wedge;
}
back
{
type wedge;
}
}
// ************************************************************************* //
ControlDict 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 1;
deltaT 0.001;
writeControl timeStep;
writeInterval 20;
purgeWrite 0;
writeFormat ascii;
writePrecision 6;
writeCompression off;
timeFormat general;
timePrecision 6;
runTimeModifiable true;
// ************************************************************************* //
Transport Properties ;
/*--------------------------------*- 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] 0.892e-6;
// ************************************************************************* //
Mesh Obtained
velocity profile
inlet
middle
outlet
As it can be seen that the velocity profile is a branch of a parabola. the analytical velocity is 0.1875 and we have obtained o.1647 with some relative error.
Pressure Drop
from the graph it can be seen that the pressure is non linear at the beginning and linear onwards. this is because at the beginning the flow is non- developed while later it is developed . the pressure drop analytical is 0.1682 and obtained is 0.021 with some relative error.
Wall shear stress
to obtain this graph we have to selet the patch pipe .
from the graph it can be seen that the wall shear stress is constant at the letter part of the length of the pipe. this is beacause of the fact that the velocity becomes constant