Week 11 - Simulation of Flow through a pipe in OpenFoam
The laminar flow through a pipe is constrained by the walls and the viscous effects grow and meet and permeate the entire flow.
The internal flow through a long duct shows there is an entrance region where an inviscid upstream flow converges and enters the tube. Viscous boundary layers grow downstream and retard the axial flow u(r,x) at the wall, thereby accelerating the center core flow to maintain continuity.
At some finite distance from entrance, boundary layers merge and the flow is then entirely viscous. The axial velocity adjusts slightly further until x = Le, it no longer changes with x and is said to be fully developed.
For laminar flow the entrance length correlation with Reynolds number is
`L_e/d = 0.06*Re_d`
The value of Reynolds for transition from laminar to turbulent is,
`Re_(d,crit) = 2300`
For laminar fully developed Poiseuille flow in a round pipe of diameter d, the pressure drop is given by the Hagen-Poiseuille equation.
`Deltap = (32*mu*V*L)/d^2 `
where V is the average velocity which is one-half of the maximum velocity. The velocity follows a paraboloid profile.
`u = u_max(1-r^2/R^2)`
`V = u_max/2`
`Re_d = (rho*V*d)/mu` Reynolds number is less than about 2300.
FLOW DOMAIN
The flow geometry is choosen to be a wedge of angle 4 degrees. This is done since the flow is axi-symmetric and simulating a flow through an entire pipe would require more resources and computation time.
The angle for wedge is choosen as such because the surface normal gradient scheme requires a correction for a non-orthogonality greater than 5 degrees(i.e the angle between the cell-to-cell vector and the normal-to-face vector should be close to 90 degrees).
The number of cells in each directions are, x = 50, y = 1, z = 300.
Using the wedge boundary condition, the flow is simulated and the results are compared with the analytical values.
For this case, the values have been defined as follows,
d = 0.01m
`Re_d = 2100`
The values of dynamic and kinematic viscosities, density are taken for water at 1 atm, 20 degrees Celsius.
`mu = 1.00*10^-3 N.s/m^2`
`nu = 1.00*10^-6 m^2/s`
`V = (Re_d*nu)/d = 0.21 m/s`
`u_max = 0.42 m/s`
`Deltap = (32*mu*V*L)/d^2 = 134.4 Pa`
kinematic pressure drop = `(Deltap)/rho = 0.1347 m^2/s^2`
`L_e = 0.06*Re_d*d = 1.26m`
1] The hydrodynamic entry length from numerical results is,
From point A(1.276,0.41681) to exit point B(1.998,0.41943) the change in velocity is only 0.6%. Hence, we can say the flow is now fully developed. Therefore, the entry length from the numerical result can be assumed to be Le = 1.276m, which is fairly close to the analytical value of 1.26m.
2] Comparison of fully developed velocity profile, from numerical results at 1.5m v/s the analytical profile
As seen, there is not much difference between the fully developed numerical and analytical results, slightly different near the central part of the pipe.
3] Maximum velocity and pressure drop
| Numerical | Analytical | Difference(numerical -analytical) | |
| Velocity m/s | 0.4183 | 0.42 | -0.0017 |
| Pressure drop (kinematic) m^2/s^2 | 0.16 | 0.1347 | +0.0253 |
The pressure drops linearly after the flow has fully developed along the length of pipe, which validates with the theory for Poiseuille flow through pipe,
`-(dP)/dx = DeltaP/L = 32*mu*V/d^2 = constant`
3]The shear stress profile after processing it from the results is,
It is zero at the center and maximum at the walls. The variation is linear along the radius which is consistent with the theoretical results.
| Theoretical `tau = 4*mu*V/R` | Numerical | Difference | |
| Shear stress at the walls`tau` Pa | 0.192 | 0.164 | -0.028 |
Table of comparison for different flow properties.
| Theoretical | Numerical | Difference(num - theor) | |
| Hydrodynamic entry length in meters | 1.26 | 1.276 | 0.016 |
| Maximum velocity in m/s | 0.42 | 0.4183 | -0.0017 |
| Kinematic pressure drop in m^2/s^2 | 0.1347 | 0.16 | 0.0253 |
| Shear stress in Pa | 0.192 | 0.164 | -0.028 |
CODE & SIMULATION FILES
The code for creating the blockMeshDict file.
The code consists a function which takes input of the file name in the '.txt','.cpp' format and outputs a text file.
The code also takes values for diameter, wedge angle, pipe length. It also computes the entrance length based on the given values and compares it with the given length and prints an error if the assumed length is shorter than entrance length.
function blockmeshdict(filename)
d = 0.01; % diameter of pipe
l = 2; %length of pipe
theta = 2*(pi/180); %half angle of wedge
nu = 1E-06;
Re = 2100;
Le = 0.06*Re*d;
Vavg = (Re*nu)/d;
Vmax = 2*Vavg;
fid = fopen(filename,'w');
if fid<0
fprintf('error opening file\n');
return;
end
if l<Le
disp('length of pipe is insufficient');
fprintf('%f %f',l ,Le);
return;
end
Title_01 = "FoamFile";
data_01 = {'version','format','class','object'};
data_02 = {'2.0;','ascii;','dictionary;','blockMeshDict;'};
Title_02 = "convertToMeters 0.01;";
Title_03 = "vertices";
C = [0 d*0.5*100*cos(theta) d*0.5*100*sin(theta) -d*0.5*100*sin(theta) l*100];
V = [C(1) C(1) C(1); C(2) C(3) C(1); C(2) C(3) C(5); C(1) C(1) C(5); C(2) C(4) C(1); C(2) C(4) C(5)];
Title_04 = "blocks";
sub04 = "hex (0 4 1 0 3 5 2 3) (20 1 10) simpleGrading (1 1 1)";
Title_05 = "edges";
Title_06 = "boundary";
front_f = "(0 1 2 3)";
back_f = "(0 3 5 4)";
tankWall_f = "(1 4 5 2)";
inlet_f = "(0 4 1 0)";
outlet_f = "(3 5 2 3)";
axis_f = "(0 3 3 0)";
Title_07 = "mergePatchPairs";
fprintf(fid,'%s\n',"/*--------------------------------*- C++ -*----------------------------------*\");
fprintf(fid,'%s\n'," ========= ");
fprintf(fid,'%s\n'," \\ / F ield | OpenFOAM: The Open Source CFD Toolbox ");
fprintf(fid,'%s\n'," \\ / O peration | Website: https://openfoam.org ");
fprintf(fid,'%s\n'," \\ / A nd | Version: 8 ");
fprintf(fid,'%s\n'," \\/ M anipulation | ");
fprintf(fid,'%s\n',"\*---------------------------------------------------------------------------*/");
fprintf(fid,'%s\n',Title_01);
fprintf(fid,'{\n');
for i = 1:4
fprintf(fid,'%10.8s',data_01{i});
fprintf(fid,'%20.14s\n',data_02{i});
end
fprintf(fid,'}\n\n');
fprintf(fid,'%s\n',"// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //");
fprintf(fid,'%s\n\n',Title_02);
fprintf(fid,'%s\n',Title_03);
fprintf(fid,'(\n');
for i = 1:6
fprintf(fid,'(');
fprintf(fid,'%6.2f,%6.5f,%2.1f',V(i,:));
fprintf(fid,')\n');
end
fprintf(fid,');\n');
fprintf(fid,'%s\n',Title_04);
fprintf(fid,'(\n');
fprintf(fid,'%s\n\n',sub04);
fprintf(fid,');\n');
fprintf(fid,'%s\n',Title_05);
fprintf(fid,'(\n');
fprintf(fid,');\n');
fprintf(fid,'%s\n',Title_06);
fprintf(fid,'(\n');
fprintf(fid,'%8.5s\n','front');
fprintf(fid,'{\n');
fprintf(fid,'%15.11s\n','type wedge;');
fprintf(fid,'%8.5s\n','faces');
fprintf(fid,'(\n');
fprintf(fid,'%15.9s\n',front_f);
fprintf(fid,');\n');
fprintf(fid,'}\n\n');
fprintf(fid,'%8.5s\n','back');
fprintf(fid,'{\n');
fprintf(fid,'%15.11s\n','type wedge;');
fprintf(fid,'%8.5s\n','faces');
fprintf(fid,'(\n');
fprintf(fid,'%15.9s\n',back_f);
fprintf(fid,');\n');
fprintf(fid,'}\n\n');
fprintf(fid,'%8.8s\n','tankWall');
fprintf(fid,'{\n');
fprintf(fid,'%15.11s\n','type wall;');
fprintf(fid,'%8.5s\n','faces');
fprintf(fid,'(\n');
fprintf(fid,'%15.9s\n',tankWall_f);
fprintf(fid,');\n');
fprintf(fid,'}\n\n');
fprintf(fid,'%8.5s\n','inlet');
fprintf(fid,'{\n');
fprintf(fid,'%15.11s\n','type patch;');
fprintf(fid,'%8.5s\n','faces');
fprintf(fid,'(\n');
fprintf(fid,'%15.9s\n',inlet_f);
fprintf(fid,');\n');
fprintf(fid,'}\n\n');
fprintf(fid,'%8.6s\n','outlet');
fprintf(fid,'{\n');
fprintf(fid,'%15.11s\n','type patch;');
fprintf(fid,'%8.5s\n','faces');
fprintf(fid,'(\n');
fprintf(fid,'%15.9s\n',outlet_f);
fprintf(fid,');\n');
fprintf(fid,'}\n\n');
fprintf(fid,'%8.5s\n','axis');
fprintf(fid,'{\n');
fprintf(fid,'%15.11s\n','type empty;');
fprintf(fid,'%8.5s\n','faces');
fprintf(fid,'(\n');
fprintf(fid,'%15.9s\n',axis_f);
fprintf(fid,');\n');
fprintf(fid,'}\n\n');
fprintf(fid,');\n\n');
fprintf(fid,'%s\n',Title_07);
fprintf(fid,'(\n');
fprintf(fid,');\n\n');
fprintf(fid,'%s\n',"// ************************************************************************* //");
end
blockMeshDictFile
/*--------------------------------*- 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 0.01;
vertices
(
( 0.00,0.00000,0.0)
( 0.50,0.01745,0.0)
( 0.50,0.01745,200.0)
( 0.00,0.00000,200.0)
( 0.50,-0.01745,0.0)
( 0.50,-0.01745,200.0)
);
blocks
(
hex (0 4 1 0 3 5 2 3) (20 1 10) simpleGrading (1 1 1)
);
edges
(
);
boundary
(
front
{
type wedge;
faces
(
(0 1 2 3)
);
}
back
{
type wedge;
faces
(
(0 3 5 4)
);
}
tankWall
{
type wall;
faces
(
(1 4 5 2)
);
}
inlet
{
File: /home/shritej/Downloads/new.txt Page 2 of 2
type patch;
faces
(
(0 4 1 0)
);
}
outlet
{
type patch;
faces
(
(3 5 2 3)
);
}
axis
{
type empty;
faces
(
(0 3 3 0)
);
}
);
mergePatchPairs
(
);
// ************************************************************************* //-
Initial pressure openFoam 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
{
front
{
type wedge;
}
back
{
type wedge;
}
tankWall
{
type zeroGradient;
}
inlet
{
type zeroGradient;
}
outlet
{
type fixedValue;
value uniform 0;
}
axis
{
type empty;
}
}
// ************************************************************************* //-
Initial 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
{
front
{
type wedge;
}
back
{
type wedge;
}
tankWall
{
type fixedValue;
value uniform (0 0 0);
}
inlet
{
type fixedValue;
value uniform (0 0 0.21);
}
outlet
{
type zeroGradient;
}
axis
{
type empty;
}
}
// ************************************************************************* //controlDict 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 "system";
object controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
application icoFoam;
startFrom startTime;
startTime 0;
stopAt endTime;
endTime 10;
deltaT 0.01;
writeControl timeStep;
writeInterval 200;
purgeWrite 0;
writeFormat ascii;
writePrecision 6;
writeCompression off;
timeFormat general;
timePrecision 6;
runTimeModifiable true;
// ************************************************************************* //