Week 11 - Simulation of Flow through a pipe in OpenFoam
Objective of the project:-
1. Calculation of flow variables.
2. Generation of the blockMeshDict using the wedge boundary condition. create the blockMeshDict file in MATLAB such that the mesh is autoatically created for various wedge angles.
3. Simulating the pipe flow in openFoam using icoFoa solver.
4. Post process the velocity profile and shear strass and comparfe the simulation results with analytical results obtained using the hagen poiseuille's equation.Post processing is carried out in the paraVIEW.
Assumptions:-
Fluid obey's the neuton's law of viscosity`(tau = mu * (dudy))`
No slip of the fluid particles at the boundary.
Flow variables:-
The working fluid is water.
Density =1000 kg/m^3
Dynamic viscosity = `8.9e-4 pa s`
Reynolds number= 2100
diameter of pipe=0.004 m
Calculated results:-
L=0.42 m
V=0.4686 m/s
`V_(max) = 0.9373 ms`
Kinematic pressure = 0.3493161
`tau_(max) = 0.831705Nm^2`
since the minimum length required for fully developed flow is 0.42m, hence we take the pipe length as 1m.
Hydrodynamic Entrance Length :-
It is the distance a flow travels in the pipe before it becomes fully developed. Entrance length refers to the length of entry region, the area following the pipe entrance where effects originating from the interior wall of the pipe propagate into the flow as an expanding boundary layer. When the boundary layer expands to fill the entire pipe, the developing flow becomes a fully developed flow, where flow characteristics no longer changes with increased distance along the pipe.
For laminar flow, the hydrodynamic entrance length is given as-
`L = 0.05 Re*D or L = .05m`
Hence the length of the pipe for simulation is taken as 1.2m.
Analytical solution by Hagen-Poiseuille equation:-
`average velocity = V_(avg) = Re. mu rho D = 0.46725 m/s`
`V_(max) = 2 V_(avg) = 0.9345 m/s`
pressure drop = `delta P = 32 mu L V_(avg) D2 = 349.316`
Kinematic pressure drop = `delta P_k = delta P rho = 0.3493161`
Matlab code for automatic generation of the blockMeshDict file:-
clear all
close all
clc
L=1;
D=0.004;
r=D/2;
theta=3;
nx=500;
ny=50;
nz=1;
line1='/*----------------------------------------*- C++ -*-------------------------------------*\';
line2='| ========= | |';
line3='| \\ / F ield | openFoam: The open source CFD Toolbox |';
line4='| \\ / O peration | Version: 4.1 |';
line5='| \\ / A nd | Web: www.openFoam.org |';
line6='| \\/ M anipulation | |';
line7='\*--------------------------------------------------------------------------------------*/';
bmd=fopen('blockMeshDict','w')
fprintf(bmd,'%s\n',line1);
fprintf(bmd,'%s\n',line2);
fprintf(bmd,'%s\n',line3);
fprintf(bmd,'%s\n',line4);
fprintf(bmd,'%s\n',line5);
fprintf(bmd,'%s\n',line6);
fprintf(bmd,'%s\n',line7);
fprintf(bmd,'Foamfile\n{\n');
fprintf(bmd,'%12s\t 2.0;\n','version');
fprintf(bmd,'%11s\t ascii;\n','format');
fprintf(bmd,'%10s\t dictionary;\n','class');
fprintf(bmd,'%11s\t blockMeshDict;\n','object');
fprintf(bmd,'}\n');
line8='//**************************************************************************************//';
fprintf(bmd,'%s\n\n',line8);
fprintf(bmd,'convertToMeters 1;\n\n');
fprintf(bmd,'vertices\n')
fprintf(bmd,'(\n')
fprintf(bmd,'\t (%d %d %d)\n',0,0,0); %point 0
fprintf(bmd,'\t (%d %d %d)\n',0,r*cosd(theta/2),r*sind(theta/2)); %point 1
fprintf(bmd,'\t (%d %d %d)\n',0,r*cosd(theta/2),-r*sind(theta/2)); %point 2
fprintf(bmd,'\t (%d %d %d)\n',L,0,0); %point 3
fprintf(bmd,'\t (%d %d %d)\n',L,r*cosd(theta/2),r*sind(theta/2)); %point 4
fprintf(bmd,'\t (%d %d %d)\n',L,r*cosd(theta/2),-r*sind(theta/2)); %point 5
fprintf(bmd,');\n\n');
fprintf(bmd,'blocks\n');
fprintf(bmd,'(\n');
fprintf(bmd,' hex (0 3 5 2 0 3 4 1) (%d %d %d) ',nx,ny,nz);
fprintf(bmd,'simpleGrading(1 0.1 1)\n\n');
fprintf(bmd,');\n\n');
fprintf(bmd,'edges\n(\n\tarc %d %d (%d %d %d)\n',1,2,0,r,0);
fprintf(bmd,'\tarc %d %d (%d %d %d)\n',4,5,L,r,0);
fprintf(bmd,');\n\n');
fprintf(bmd,'boundary\n(\n');
fprintf(bmd,'\t inlet\n');
fprintf(bmd,'\t{\n');
fprintf(bmd,'\t\t type patch;\n');
fprintf(bmd,'\t\t faces\n');
fprintf(bmd,'\t\t(\n');
fprintf(bmd,'\t\t\t (0 1 2 0)\n');
fprintf(bmd,'\t\t );\n');
fprintf(bmd,'\t}\n');
fprintf(bmd,'\t outlet\n');
fprintf(bmd,'\t{\n');
fprintf(bmd,'\t\t type patch;\n');
fprintf(bmd,'\t\t faces\n');
fprintf(bmd,'\t\t (\n');
fprintf(bmd,'\t\t\t (3 5 4 3)\n');
fprintf(bmd,'\t\t );\n');
fprintf(bmd,'\t}\n');
fprintf(bmd,'\t top\n');
fprintf(bmd,'\t{\n');
fprintf(bmd,'\t\t type wall;\n');
fprintf(bmd,'\t\t faces\n');
fprintf(bmd,'\t\t (\n');
fprintf(bmd,'\t\t\t (1 4 5 2)\n');
fprintf(bmd,'\t\t );\n');
fprintf(bmd,'\t}\n');
fprintf(bmd,'\t back\n');
fprintf(bmd,'\t{\n');
fprintf(bmd,'\t\t type wedge;\n');
fprintf(bmd,'\t\t faces\n');
fprintf(bmd,'\t\t (\n');
fprintf(bmd,'\t\t\t (0 2 5 3)\n');
fprintf(bmd,'\t\t );\n');
fprintf(bmd,'\t}\n');
fprintf(bmd,'\t front\n');
fprintf(bmd,'\t{\n');
fprintf(bmd,'\t\t type wedge;\n');
fprintf(bmd,'\t\t faces\n');
fprintf(bmd,'\t\t (\n');
fprintf(bmd,'\t\t\t (0 3 4 1)\n');
fprintf(bmd,'\t\t );\n');
fprintf(bmd,'\t}\n');
fprintf(bmd,'\t axis\n');
fprintf(bmd,'\t{\n');
fprintf(bmd,'\t\t type empty;\n');
fprintf(bmd,'\t\t faces\n');
fprintf(bmd,'\t\t (\n');
fprintf(bmd,'\t\t\t (0 3 3 0)\n');
fprintf(bmd,'\t\t );\n');
fprintf(bmd,'\t}\n');
fprintf(bmd,');\n\n')
fprintf(bmd,'mergePatchPairs\n(\n');
fprintf(bmd,');\n')
fprintf(bmd,'%s\n',line8);
fclose(bmd)blockMeshDict file generated for wedge angle =`3 degree`
/*--------------------------------*- C++ -*----------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Version: 9
\\/ M anipulation |
\*---------------------------------------------------------------------------*/
FoamFile
{
format ascii;
class dictionary;
object blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
convertToMeters 0.1;
vertices
(
(0 0 0)
(0.6 0 0)
(0 0.4 0)
(0.6 0.4 0)
(1 0.4 0)
(0 1 0)
(0.6 1 0)
(1 1 0)
(0 0 0.1)
(0.6 0 0.1)
(0 0.4 0.1)
(0.6 0.4 0.1)
(1 0.4 0.1)
(0 1 0.1)
(0.6 1 0.1)
(1 1 0.1)
);
blocks
(
hex (0 1 3 2 8 9 11 10) (12 8 1) simpleGrading (1 1 1)
hex (2 3 6 5 10 11 14 13) (12 12 1) simpleGrading (1 1 1)
hex (3 4 7 6 11 12 15 14) (8 12 1) simpleGrading (1 1 1)
);
edges
(
);
boundary
(
lid
{
type wall;
faces
(
(5 13 14 6)
(6 14 15 7)
);
}
fixedWalls
{
type wall;
faces
(
(0 8 10 2)
(2 10 13 5)
(7 15 12 4)
(4 12 11 3)
(3 11 9 1)
(1 9 8 0)
);
}
frontAndBack
{
type empty;
faces
(
(0 2 3 1)
(2 5 6 3)
(3 6 7 4)
(8 9 11 10)
(10 11 14 13)
(11 12 15 14)
);
}
);
mergePatchPairs
(
);
// ************************************************************************* //
Velocity boundary condition:-
/*--------------------------------*- C++ -*----------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Version: 9
\\/ 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.46725 0 0);
}
outlet
{
type zeroGradient;
}
top
{
type noSlip;
}
back
{
type wedge;
}
front
{
type wedge;
}
axis
{
type empty;
}
}
// ************************************************************************* //Pressure boundary condition:-
/*--------------------------------*- C++ -*----------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Version: 9
\\/ 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.3493161;
}
top
{
type zeroGradient;
}
back
{
type wedge;
}
front
{
type wedge;
}
axis
{
type empty;
}
}
// ************************************************************************* //OpenFoam case setUp:-
A wedge boundary condition is used in this case as it is cylinder and flow is axi-symmetric along the axis.
After updating the blockMeshDict and pressure and velocity boundary conditions file simulation is ran using the icoFoam solver which is used to simulate the turbulent, incompressible flows.
Base mesh for wedge angle=`3 degree`
Velocity contours:-
Pressure contours:-
Velocity profile At 0.01m:-
the velocity profile near the entrance of the pipe is flat. It means the flow is still developing.
Velocity profile At 0.1m:-
velocity profile At 0.43m:-
velocity profile At 0.75m:-
A parabolic velocity profile is obtained after a certain distance in the pipe. This profiles remains unchanged thereafter. it concludes that the length of the pipe as calculated theoritically is suffiecient to get a fully developed flow at the outlet of the pipe.
The shear stress is calculated from the given formula
`tau = mu dudy Pa`