Flow through Pipe using OpenFoam
Objective : Simulate the flow through the pipe in the openFoam using the tutorials.
Need to show :
1. Shear stress Distribution.
2. Checking for the Hydrodynamic Entrance length Numerically.
3. validate wall shear stress for the Fully developed Flow.
4. Maximum velocity and pressure drop for the fully developed flow.
5. Velocity Profile for the Fully Developed Flow.
6. Comparing all the results with the analyticall solution.
Assumptions :
1. Steady , incompressible , Laminar flow.
2. Axissymmetric flow ( `frac{delU_{theta}}{deltheta}` = 0 ).
3. Constant Fluid Properties (density , viscosity ).
Since the flow is axissymmetric flow , we can consider a part of the pipe (certain degree wedge).
Fluid Considered : ( Water - Newtonian Fluid ) -> At Normal Conditions
1. Density = 996.5 `frac{kg}{m^3}` = 0.9965`frac{g}{cm^3}`
2. Dynamic Viscosity = 0.8509 `*``10^3` pa.s.
3. Kinematic Viscosity = 0.8539 `*10^6 frac{m^3}{s}`.
Hagen-poisuille Flow :
Fluid Flow based on the Hagen-Poiseuille law , gives the pressure drop along a cylindrical pipe of constant cross-section in a Newtonian , incompressible and laminar Flow.
Poiseuille's Equations decribe the pressure drop due to the viscosity of the fuid.
Considered No acceleration along the flow.
characteristics Of pipe :
1. Diameter of pipe(D) = 0.002m.
2. Radius of pipe (R) = 0.001m.
Hagen-Poiseuille relations :
1. Pressure Drop (`gradP`) = `frac{32*mu*V_{avg}*L}{D^2}`
2. Velocity (U) = `-frac{1}{4*mu}*frac{delP}{delX}*(R^2-r^2)`
3. `V_max` = `-frac{1}{4*mu}*frac{delP}{delX}*(R^2)`
4. `V_{avg} = -frac{1}{8*mu}frac{delP}{delX}*R^2`
5. `V_{max} = 2* V_{avg}`
6. Shear Stress (`tau)` = `frac{2*mu*V_{max}*r}{R^2`
Reynolds Number for the Laminar Floe(re) = 2100.
As known that the prandtl explained due to the No slip conditions , fluid adheres the solid surface Due to the cohesive forces between the molecules of the fluids , causes the above Layer to restrict to
Flow the Downstream , and this effect goes on increasing so the Boundary layer thickness Increases . Boundary layer is the region where the velocity gradient present , In other words can say boundary
layer as the locus of the points of the velocity = 0.99*U , U(input velocity). Outside the boundary layer the flow is potential flow( Incompressible + Inviscid + Irrotational ).
A velocity Gradients developes along the axial direction or along the length of the pipe.To make the mass flow rate constant ( Continuity equation satisfied ).
Region in which the inside the Boundary Layer is velocity Boundary Layer , simply the Boundary Layer . No friction effects outside the Boundary Layer. Region of the from the inlet of the pipe , to the point at
which the flow becomes the fully developed Flow called the Hydeodynamic Entrance , and the length is called the Hydrodynamic Entrance Length() . Developing of the Flow is called the Hydrodynamically
developing flow , after the Hydrodynamic entrance length , flow is Fully Developed Flow . Velocity Profile in the Fully developed region is parabolic for the laminar flow , shown below...,
Hyderodynamically Fully Developed Flow ->
Hydrodynamic Entrance Length ->
Since in the Fully developed the velocity profile remains constant , so the wall shear stress remains
constant , as the dynamic viscosity of the fluid considered constant.
solver ::
Solver used for this project -> icoFoam. (Steady + incompressible + Laminar ) , solver.
simpleFoam solver can also be used ( steady + incompressible + Turbulent ).
Matlab File for the BlockMeshDict File ::
Creating a MatlabScriptFile Automatically using the MATLAB , for any model of the cylinder shape.
Here considered ;
re (Reynolds Number) = 2100.
d (Diamter) = 0.02m;
rho (Density) = 997kg/`m^3`.
mu(Dynamic Viscosity) = 0.89`*10^3`.
close all
clear all
clc
re = input("Enter the Reynolds Number for the laminar flow for a flow inside the cylinder");
rho = input("Enter the Density of the fluid");
mu = input("Enter the Dynamic Viscosity of the fuid");
d = input("Enter the diameter of the cylinder");
r = d/2;
tht = input("Enter the required angle");
l_entry = 0.06*re*d;
L = l_entry + 0.5;
% Boundary conditions
v_avg = (re*mu)/(rho*d);
v_max = 2*v_avg;
dp = (32*v_avg*mu*L)/(d)^2;
nu = mu/rho;
dp_k = dp/rho;
%
R = linspace(0 , r , 1000);
for i = 1:length(R)
tau(i) = 2*mu*v_max*(abs(R(i)))/r^2;
end
plot(R,tau,'b')
grid on
xlabel('radial Distance')
ylabel('shear stress(N/m)')
title('shear stress using analytical results')
%%
format short e
f1 = fopen("blockMeshDictFile3.txt",'w');
%Specifying string for header part
h{1}='/*––––––––––––––––––––––––––––––--*- C++ -*–––––––––––––––––––––––––––––––––-* ';
h{2}='| ========= | | ';
h{3}='| / F ield | OpenFOAM: The Open Source CFD Toolbox | ';
h{4}='| / O peration | Version: 4.0 | ';
h{5}='| / A nd | Web: www.OpenFOAM.org | ';
h{6}='| / M anipulation | | ';
h{7}=' *–––––––––––––––––––––––––––––––––––––––––––––––––––––––––––––––––––––––––––*/ ';
h{8}=' FoamFile';
h{9}='{';
h{10}=' version 2.0;';
h{11}=' format ascii;';
h{12}=' class dictionary;';
h{13}=' object blockMeshDict;';
h{14}='}';
h{15}='// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //';
tht = tht/2;
v1 = [0 0 0];
v2 = [0 r*cosd(tht) r*sind(tht)];
v3 = [0 r*cosd(tht) -r*sind(tht)];
v4 = [L 0 0];
v5 = [L r*cosd(tht) r*sind(tht)];
v6 = [L r*cosd(tht) -r*sind(tht)];
%
for i=1:15
fprintf(f1,'%s n',h{i});
end
fprintf(f1,'convertToMeters 1;n');
s1 = blanks(5);
s2 = blanks(10);
s3 = blanks(15);
fprintf(f1,'vertices n(');
fprintf(f1,'n',s1);
fprintf(f1,'%s (0 0 0) // 0 n',s1);
fprintf(f1,'%s (%d %d %d) // 1 n',s1,v2(1),v2(2),v2(3));
fprintf(f1,'%s (%d %d %d) // 2 n',s1,v3(1),v3(2),v3(3));
fprintf(f1,'%s (%d %d %d) // 3 n',s1,v4(1),v4(2),v4(3));
fprintf(f1,'%s (%d %d %d) // 4 n',s1,v5(1),v5(2),v5(3));
fprintf(f1,'%s (%d %d %d) // 5 n); n',s1,v6(1),v6(2),v6(3));
fprintf(f1,'blocks n(n');
fprintf(f1,'%s hex (0 3 5 2 0 3 4 1) (300 50 1) simpleGrading (1 0.2 1) n);n',s1);
fprintf(f1,'edgesn');
fprintf(f1,'%s arc 1 2 (0.000000 %d 0.000000)n',s1,r);
fprintf(f1,'%s arc 4 5 (%d %d 0.00000)n);n',s1,L,r);
fprintf(f1,'boundaryn(n');
fprintf(f1,'%sinflown%s{n%s',s1,s1,s2);
fprintf(f1,'type patch;n');
fprintf(f1,'%s',s2);
fprintf(f1,'facesn%s(n',s2);
fprintf(f1,'%s',s3);
fprintf(f1,'(0 1 2 0)n%s);n%s}n',s2,s1);
fprintf(f1,'%soutflown%s{n%s',s1,s1,s2);
fprintf(f1,'type patch;n');
fprintf(f1,'%s',s2);
fprintf(f1,'facesn%s(n',s2);
fprintf(f1,'%s',s3);
fprintf(f1,'(3 5 4 3)n%s);n%s}n',s2,s1);
fprintf(f1,'%sfrontn%s{n%s',s1,s1,s2);
fprintf(f1,'type wedge;n');
fprintf(f1,'%s',s2);
fprintf(f1,'facesn%s(n',s2);
fprintf(f1,'%s',s3);
fprintf(f1,'(0 3 4 1)n%s);n%s}n',s2,s1);
fprintf(f1,'%sbackn%s{n%s',s1,s1,s2);
fprintf(f1,'type wedge;n');
fprintf(f1,'%s',s2);
fprintf(f1,'facesn%s(n',s2);
fprintf(f1,'%s',s3);
fprintf(f1,'(0 2 5 3)n%s);n%s}n',s2,s1);
fprintf(f1,'%swalln%s{n%s',s1,s1,s2);
fprintf(f1,'type wall;n');
fprintf(f1,'%s',s2);
fprintf(f1,'facesn%s(n',s2);
fprintf(f1,'%s',s3);
fprintf(f1,'(1 4 5 2)n%s);n%s}n',s2,s1);
fprintf(f1,'%saxisn%s{n%s',s1,s1,s2);
fprintf(f1,'type empty;n');
fprintf(f1,'%s',s2);
fprintf(f1,'facesn%s(n',s2);
fprintf(f1,'%s',s3);
fprintf(f1,'(0 3 3 0)n%s);n%s}n',s2,s1);
fprintf(f1,');nn');
fprintf(f1,'mergePatchPairsn');
fprintf(f1,'(n);n');
fprintf(f1,'// *************************************************************************** //');
fclose(f1);
BlockMeshDictFile ::
/*––––––––––––––––––––––––––––––--*- C++ -*–––––––––––––––––––––––––––––––––-*
| ========= | |
| / F ield | OpenFOAM: The Open Source CFD Toolbox |
| / O peration | Version: 4.0 |
| / A nd | Web: www.OpenFOAM.org |
| / M anipulation | |
*–––––––––––––––––––––––––––––––––––––––––––––––––––––––––––––––––––––––––––*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
convertToMeters 1;
vertices
( (0 0 0) // 0
(0 9.996573e-03 2.617695e-04) // 1
(0 9.996573e-03 -2.617695e-04) // 2
(3.020000e+00 0 0) // 3
(3.020000e+00 9.996573e-03 2.617695e-04) // 4
(3.020000e+00 9.996573e-03 -2.617695e-04) // 5
);
blocks
(
hex (0 3 5 2 0 3 4 1) (300 50 1) simpleGrading (1 0.1 1)
);
edges
(
arc 1 2 (0.000000 1.000000e-02 0.000000)
arc 4 5 (3.020000e+00 1.000000e-02 0.000000)
);
boundary
(
inflow
{
type patch;
faces
(
(0 1 2 0)
);
}
outflow
{
type patch;
faces
(
(3 5 4 3)
);
}
front
{
type wedge;
faces
(
(0 3 4 1)
);
}
back
{
type wedge;
faces
(
(0 2 5 3)
);
}
wall
{
type wall;
faces
(
(1 4 5 2)
);
}
axis
{
type empty;
faces
(
(0 3 3 0)
);
}
);
mergePatchPairs
(
);
// *************************************************************************** //ControlDictFile ::
/*--------------------------------*- 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 30;
deltaT 0.01;
writeControl timeStep;
writeInterval 20;
purgeWrite 0;
writeFormat ascii;
writePrecision 6;
writeCompression off;
timeFormat general;
timePrecision 6;
runTimeModifiable true;
// ************************************************************************* //TransportProperties 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 "constant";
object transportProperties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
transportModel Newtonian;
nu [0 2 -1 0 0 0 0] 0.892e-06;
rho [1 3 0 0 0 0 0] 997;
// ************************************************************************* //Initial and Boundary Conditions ::
Pressure Boundary Condition:::
/*--------------------------------*- 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
{
outflow
{
type fixedValue;
value uniform 0;
}
inflow
{
type zeroGradient;
}
wall
{
type zeroGradient;
}
front
{
type wedge;
}
back
{
type wedge;
}
axis
{
type empty;
}
}
// ************************************************************************* //Velocity Boundary Condition::
/*--------------------------------*- 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
{
inflow
{
type fixedValue;
value uniform (0.0937 0 0);
}
outflow
{
type zeroGradient;
}
front
{
type wedge;
}
back
{
type wedge;
}
wall
{
type noSlip;
value uniform (0 0 0);
}
axis
{
type empty;
}
}
// ************************************************************************* //
Post-processing of the Results ::
ParaView is invoked by using the command paraFoam in the command terminal/
This Figure shows the Boundary Profile Formation at the wall of the pipe.
Fully Developed Boundary Layer formed at a distance x = 2.52m.
Velocity Profile over different axial distance along the cylinder ::
V_max from the analytical Result = 0.1875 m/s .
V_max from the CFD Result = 0.186447 m/s .
error = |0.1875 - 0.186447| = 0.001053.
Pressure Distribution along the cylinder ::
Pressure along the length ::
analytical result:
`DeltaP = frac{32*mu*v_{avg}*L}{rho*d^2}`
`Deltap = 0.0202`
CFD result : 0.025798.
Error -> |0.0202 - 0.025798| = 0.00558.
Tau_max = 0.032pa.
Analytical Shear Stress ::
Using the Reflect property along the y-axis ::
Velocity profile for the fully developed flow ::
Profile of the velocity is parabola , velocity is maximum at the center of the cylinder -> As we know form the Hagen-Poisuille theory.
CFD Information ::
v
Velocity profile of analytical and CFD results , along the radial distance for fully Developed Flow ::
Wedge Boundary Condition used in the CFD , for to define the Axisymmetric situation.
Conclusion ::
From the above results , can say that the CFD result are more accurate to the analytical result.
Residual is very less.
As observed the shear stress are maximum at the walls.
Hagen-Poiseuille Flow was simulated by using a reduced model in the form of a wedgw with a specific
angle , instead of a fully cylindrical pipe , since flow is considered as a axis-symmetric flow.
results from the CFD simulation ,in good agreement with the analytical result.