Week 11 - Simulation of Flow through a pipe in OpenFoam
Aim-To simulates the flow through the pipe in OpenFOAM.
- The purpose of this project is to simulate axisymmetric flow by applying the wedge boundary condition. The flow is laminar and flows through a constant cross-section pipe.
- The Reynolds number is given as 2100.
- The laminar flow through a constant cross-sectional pipe is governed by the Hagen poiseuille equation.
- The assumption of this equation is fluid should be incompressible and Newtonian and there is no acceleration in the pipe.
- Consider the pipe with a uniform cross-section of length L and diameter D.The free stream velocity of the fluid is V.
- There is a no-slip condition at the pipe wall, so the velocity of the fluid in contact with the pipe wall is zero.
The flow at the inlet enters with free stream velocity and velocity gradient develops along the y-direction, the flow is fully developed after a certain length of the pipe that length is called the entrance length of the pipe.
The entrance length is calculated by the following formula
Le=0.06*Re*D
Where
Le = entrance length
Re = Reynolds number
D= diameter of the pipe
Properties have been given for the simulation.
Re=2100;
D=0.015 m(Diameter of the pipe)
Le=2.52 m(Entrance length)
`rho=1000 (kg)/m^3`
`mu=8.926*10^-4 pa*s`
Analytical solution according to hagen poiseuille equation
`Vavg=(Re*mu)/(rho*D)`
=0.1249 m/s
Now,according to hagen poiseuille equation
Vmax=2*Vavg
=0.2499 m/s
`Deltap=(32*mu*Vavg*L)/(D^2)`
=31.71 pa
Terms used in openFoam for the simulation of the flow.
Patch=no geometric or topological information
Wall=solid wall
Empty=dimension for which no solution is needed
Wedge=used in the simulation of the axis symmetric flow.
Geometry created in the openfoam with the help of blockMeshdict file.
0 1 2 0 =inlet patch
3 5 4 3=outlet patch
0 3 4 1=Front
0 2 5 3=back
1 4 5 2=Top
matlab code for the blockmeshDict file
clear all
close all
clc
%% creating the geometry file for the openFoam in matlab
d=0.015 ; % diameter of the pipe in meter
theta=4 % wedge angle
L=2.5;
v0=[0 0 0];
v1=[0 r*sind(theta/2) -r*cosd(theta/2)];
v2=[0 r*sind(theta/2) r*cosd(theta/2)];
v3=[L 0 0];
v4=[L r*sind(theta/2) -r*cosd(theta/2)];
v5=[L r*sind(theta/2) r*cosd(theta/2)];
n1=50; %number of cells along the pipe
n2=1 %number of cells radially
n3=300; %number of cells along the curve
gfactor=0.2;
f1=fopen('blockMeshDict.txt','w')
f2=fopen('text.txt','r')
for i=1:17
A=fgetl(f2);
h_line=convertCharsToStrings(A);
fprintf(f1,'%s\n',h_line)
end
fprintf(f1,'vertices\n\n');
fprintf(f1,'(\n');
fprintf(f1,'\t(%f %d %d)\n',v0(1),v0(2),v0(3));
fprintf(f1,'\t(%f %f %f)\n',v1(1),v1(2),v1(3));
fprintf(f1,'\t(%f %f %f)\n',v2(1),v2(2),v2(3));
fprintf(f1,'\t(%f %f %f)\n',v3(1),v3(2),v3(3));
fprintf(f1,'\t(%f %f %f)\n',v4(1),v4(2),v4(3));
fprintf(f1,'\t(%f %f %f)\n\n',v5(1),v5(2),v5(3));
fprintf(f1,');\n');
fprintf(f1,'\n');
fprintf(f1,'%s\n','blocks');
fprintf(f1,'(\n');
%fprintf(f1,'blocksn(n');
fprintf(f1,' \thex (0 1 2 3 4 5 3)(%d %d %d)',n1,n2,n3);
fprintf(f1,'simpleGrading(%0.1f 1 1)\n',gfactor);
fprintf(f1,');\n\n');
fprintf(f1,'edges\n(\n\tarc %d %d (%d %0.2f %d)\n\tarc %d %d(%0.2f %0.2f %d)\n);\n\n',1,2,0,r,0,4,L,r,0,3);
fprintf(f1,'boundary\n(');
fprintf(f1,'\t\n\t%s\n\t{\n\t\t type %s ;\n\t\t faces \n\t\t (\n\t\t\t (%d %d %d %d) \n\t\t );\n\t}\n\n','axis','empty',0,3,3,0);
fprintf(f1,'\t%s\n\t{\n\t\t type %s ;\n\t\t faces \n\t\t (\n\t\t\t (%d %d %d %d) \n\t\t );\n\t}\n\n','inlet','patch',0,2,1,0);
fprintf(f1,'\t%s\n\t{\n\t\t type %s ;\n\t\t faces \n\t\t (\n\t\t\t (%d %d %d %d) \n\t\t );\n\t}\n\n','pipe','wall',2,5,4,1);
fprintf(f1,'\t%s\n\t{\n\t\t type %s ;\n\t\t faces \n\t\t (\n\t\t\t (%d %d %d %d) \n\t\t );\n\t}\n\n','outlet','patch',3,4,5,3);
fprintf(f1,'\t%s\n\t{\n\t\t type %s ;\n\t\t faces \n\t\t (\n\t\t\t (%d %d %d %d) \n\t\t );\n\t}\n\n','front','wedge',0,3,5,2);
fprintf(f1,'\t%s\n\t{\n\t\t type %s ;\n\t\t faces \n\t\t (\n\t\t\t (%d %d %d %d) \n\t\t );\n\t}\n\n','back','wedge',0,1,4,3);
fprintf(f1,');\n\n');
fprintf(f1,'mergePatchPairsn(n);\n\n');
fprintf(f1,'(\n');
fprintf(f1,');');
Blockmesh dict 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;
object blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
convertToMeters 1;
vertices
(
(0.000000 0 0)
(0.000000 7.495e-3 0.2617e-3)
(0.000000 7.495e-3 -0.2617e-3)
(2.5 0 0)
(2.5 7.495e-3 0.2617e-3)
(2.5 7.495e-3 -0.2617e-3)
);
blocks
(
hex (0 3 5 2 0 3 4 1) (200 20 1) SimpleGrading (1 0.1 1)
);
edges
(
arc 1 2 (0 7.5e-3 0)
arc 4 5 (2.5 7.5e-3 0)
);
boundary
(
inlet
{
type patch ;
faces
(
(0 1 2 0)
);
}
outlet
{
type patch;
faces
(
(3 4 5 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 file of openFoam
/*--------------------------------*- 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 50;
deltaT 0.02;
writeControl timeStep;
writeInterval 20;
purgeWrite 0;
writeFormat ascii;
writePrecision 6;
writeCompression off;
timeFormat general;
timePrecision 6;
runTimeModifiable true;
// ************************************************************************* //here simulation time= 50 seconds
time step =0.02
pressure 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;
}
axis
{
type empty;
}
front
{
type wedge;
}
back
{
type wedge;
}
top
{
type zeroGradient;
}
}
// ************************************************************************* //at inlet velocity condition is given so the pressure is given as the zeroGradient and outlet is given as the zero pressure
velocity 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.119 0 0);
}
outlet
{
type zeroGradient;
}
axis
{
type empty;
}
front
{
type wedge;
}
back
{
type wedge;
}
top
{
type fixedValue;
value uniform (0 0 0);
}
}
// ************************************************************************* //Here at inlet 0.119 is selected for the velocity and at the outlet, it is given as zeroGradient because the pressure is given at the outlet.
The top wall condition is given as the fixedValue and the value is given as zero, so at the top no-slip condition is satisfied
post processing
velocity profile with y-direction
at inlet
at 1.25 meter from origin
at 2 meter from origin
From the velocity plot after the 2 meter from the origin we can see that the maximum velocity is at the centre and at the top wall the velocity is zero.
The magnitude of the maximum velocity is 0.2499 m/s from the analytical method which we can conform from the velocity plot.
pressure variation with x direction
shear stress distribution for y directioin
