Week 11 - Simulation of Flow through a pipe in OpenFoam
AIM : In this challenge, we will simulate an axi-symmetric flow by applying the wedge boundary condition. The flow is laminar and flows through a constant cross-section pipe at an inlet Reynold's number of 2100. You can do a transient or steady state simulation and you will validate your results with the Hagen- Poiseuille's equations for the fully developed flow.
You will also write code in Matlab to automate the generation of blockMeshDict file. The code should take inputs like pipe length, diameter, wedge angle, etc and output a fully functional blockMeshDict file.
Expected Results
- Validate Hydro-dynamic length with the numerical result
- Validate the fully developed flow velocity profile with its analytical profile
- Validate maximum velocity and pressured drop for fully developed flow
- Post-process Shear stress and validate wall shear stress for fully developed flow
Expected Velocity Profile :
- The distance the fluid has to travel become fully developed is given by the entry length formula.
- Given Reynolds number is 2100.Tis implies flow is laminar.
- Hence entry length formula for laminar is given by L/d = 0.06Re.
- Where L is effective length; d is dia of pipe and Re is reynolds number.
- Assuming the dia of pipe as 10mm; substituting the values.
EQUATION :
`△ p = 8 μ π L Q/ π R ^4 = 8 π μ L Q/ A^ 2`
where,
- triangle p = pressure difference between the two ends
- L = is the length of pipe
- mu = is the dynamic viscosity
- Q = is the volumetric flow rate
- R = is the pipe radius
- A = is the cross section of pipe.
Procedure :
- First the geometry is created in the blockmesh dict file,this is done by creating the file through the matlab code.
- The mesh is generated according to requirement.
- The velocity and pressure are applied.
- The controlDict file is edited for simulation run ttime and timestep.
- After the geometry is input the simulation can be started by calling 'icoFoam' command in terminal.
- 'paraviewWindows' command helps to open paraview.The shear stress diagram is obtained by post-processing the velocity results obtained from openfoam.
Matlab Code For Generating BlockmeshDict file :
clear
close all
clc
Re = 2100;
Ld = 2;
nu = 8.917e-7;
rho = 997;
D = 0.01;
mu = nu*rho;
le = 0.06*Re*D;
l = le +Ld;
V_avg = mu*Re/(rho*D);
V_max = 2*V_avg;
r = D/2 ;
% Hagen Poisulle equation
delta_p = 8*mu*V_avg*(l/r^2);
kinematic_p = delta_p/rho;
% Wedge angle
theta = 4;
% calculating shear stress
R = linspace(0,D/2,1000);
for i = 1:length(R)
tau(i) = 2*mu*V_max*(abs(R(i)))/r^2;
end
% plot
plot(R,tau)
xlabel('Radius in m')
ylabel('shear stress')
title('variation of shear stress along radius cylinder')
grid on
hold on
% Defining the vertices of block
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)];
% header file
h1 ='/*--------------------------------*- C++ -*----------------------------------*\'
h2 ='| ========= | |'
h3 ='| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |'
h4 ='| \\ / O peration | Version: v2112 |'
h5 ='| \\ / A nd | Website: www.openfoam.com |'
h6 ='| \\/ M anipulation | |'
h7 ='\*---------------------------------------------------------------------------*/'
h8 ='FoamFile'
h9 ='{'
h10 =' version 2.0;'
h11 =' format ascii;'
h12 =' class dictionary;'
h13 =' object blockMeshDict;'
h14 ='}'
h15 =' // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //'
f1 = fopen('blockMeshDict','w')
fprintf(f1,'%s \n',h1);
fprintf(f1,'%s \n',h2);
fprintf(f1,'%s \n',h3);
fprintf(f1,'%s \n',h4);
fprintf(f1,'%s \n',h5);
fprintf(f1,'%s \n',h6);
fprintf(f1,'%s \n',h7);
fprintf(f1,'\n');
fprintf(f1,'%s \n',h8);
fprintf(f1,'%s \n',h9);
fprintf(f1,'%s \n',h10);
fprintf(f1,'%s \n',h11);
fprintf(f1,'%s \n',h12);
fprintf(f1,'%s \n',h13);
fprintf(f1,'%s \n',h14);
fprintf(f1,'%s \n',h15);
fprintf(f1 ,'\n');
fprintf(f1,'%2s \n','convertToMeters 1;');
fprintf(f1,'\n');
% Defining vertices
fprintf(f1,'%s \n','Vertices');
fprintf(f1,'%s \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',v5(1),v5(2),v5(3));
fprintf(f1,'%s \n',');');
% Defining blocks
fprintf(f1,'\n');
fprintf(f1,'%s \n','blocks');
fprintf(f1,'%s \n','(');
fprintf(f1,'\t hex (0 3 4 1 0 3 5 2) (%d %d %d) simpleGrading (1 0.1 1)\n',300,30,1);
fprintf(f1,'%s \n',');');
% Defining edges
fprintf(f1,'\n');
fprintf(f1,'%s \n','edges');
fprintf(f1,'%s \n','(');
fprintf(f1,' \t arc 1 2 (%d %d %d)\n',0,r,0);
fprintf(f1,' \t arc 4 5 (%d %d %d)\n',l,r,0);
fprintf(f1,'%s \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, '); \n\n')
fprintf(f1, 'mergePatchPairs\n')
fprintf(f1, '(\n')
fprintf(f1, ');\n')
header9='// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //'
fprintf(f1, '%s\n',header9)
fclose(f1)
Generated BlockmeshDict File :
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2112 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
convertToMeters 1;
vertices
(
(0 0 0)
(0 4.996954e-03 -1.744975e-04)
(0 4.996954e-03 1.744975e-04)
(3.260000e+00 0 0)
(3.260000e+00 4.996954e-03 -1.744975e-04)
(3.260000e+00 4.996954e-03 1.744975e-04)
);
blocks
(
hex (0 3 4 1 0 3 5 2) (300 30 1) simpleGrading (1 0.1 1)
);
edges
(
arc 1 2 (0 5.000000e-03 0)
arc 4 5 (3.260000e+00 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
(
);
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //ContriolDict file :
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2112 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
application icoFoam;
startFrom startTime;
startTime 0;
stopAt endTime;
endTime 3;
deltaT 0.005;
writeControl timeStep;
writeInterval 20;
purgeWrite 0;
writeFormat ascii;
writePrecision 6;
writeCompression off;
timeFormat general;
timePrecision 6;
runTimeModifiable true;
// ************************************************************************* //Pressure File :
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2112 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ 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 0;
}
front
{
type wedge;
}
back
{
type wedge;
}
}
// ************************************************************************* //Velocity File :
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2112 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ 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
{
axis
{
type empty;
}
top
{
type noSlip;
}
inlet
{
type fixedValue;
value uniform (0.187 0 0);
}
outlet
{
type zeroGradient;
}
front
{
type wedge;
}
back
{
type wedge;
}
}
// ************************************************************************* //TransportProperties File :
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2112 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object transportProperties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
nu [0 2 -1 0 0 0 0] 8.9e-7;
// ************************************************************************* //Plot :
The plot shows the variation of shear stress along radius cylinder ,it generated by matlab.
Meshed Geometry :
The above shown image is the sector of pipe which is assumed .It is generated in openFoam.
Velocity Profile :
The above plot is generated in the openfoam by using plotoverline .
The velocity profile changes until the flow is fully developed.Distance from inlet ,the velocity attains its characteristics parabola shape
The plot of velocity profile is parabolic in nature with max velocity attained at the center of the pipe.
Shear -Stress Profile :
Conclusion:
Shear stress increases the friction on the walls of the pipe because of which the velocity is maximum at the center of the pipe.
The shear stress plot is obtained after post-processing and it is found that the shear stress is maximum at the walls and zero at the center of the pipe.