Week 12 - Validation studies of Symmetry BC vs Wedge BC in OpenFOAM vs Analytical H.P equation
Aim :
To simulate the flow through pipe using axis symmetry by applying wede and symmetry boundry conditions in OpenFoam.Consider the flow to be laminar with Reynolds number 2100 for any diameter of pipe.Compare th eresults with the Hagen-Poiseuille's equations for flow through circular pipe.
Theory :
Subject of the boundaries is quite complex because their role in modelling is not simply that of geometric entity but
an integral part f solutions and numerics through boudary conditions or inter boundary connections.A discussion of boundaries sits uncomfortablly between discussion on meshes,fields,discretisation,computational processing etc.We first need to consider that,for the purpose of applying boundary conditionsa boundary is generally broken up into a set of paches.One patch may include one or more area of the boundary surface which do not necessarilly need to be physically connected.A type is assigned to every patch as apart of mesh descriptio,as part of the boundary file.It describes the patch in terms of geometry or a data communicaton link.Below shown is the example of an icoFoam case.In this case we are compairing effects of different boundary conditions i.e. wedge and symmetry.
Consider a incompressible laminar and newtonian flow through a pipe having constant circular cross section
throughout the length(Hagen Poiseuille flow).The diameter of the pipe D = 0.03m and length of the pipe L= 5m and water is the working fluid. First we need to generate the blockmesh dict file,which is later used by openfoam to generate the mesh.Then we need to use the above generated file to set up the physical model and run the simulation using OpenFoam.Geomety,loadings,boundary conditions and materials are symmetric with respect to z axis,thus ,the flow has axisymmetric nature.The suitable wedge and symmetry angles arei considered.
So in the results we need to observe the change in velocity profile and pressure drop across the length of the pipe. We also need to calculate the magnitude of shear stress acting on the walls of the domian.At last,validate results with analytical solution obtained from Hagen Poiseullle equation. Assume constant velocity inlet and constant zero gauge pressure outlet as the initial boundary condition.
Wedge :
For 2D axi-symmetric cases,wedge type boundary condition is used e.g. a cylinder,the geometry is specified as a wedge of small angle and 1 cell thick running along the pannel of symmetry,straddling one of the co-ordinate planes.The axi_symmetric wedge planes must be specified as a sepaate patches of wedge type.The angle of wedge should be <5 degrees.
Symmetry :
It is similar to wedge boundary condition in terms of function in addition it is used for any non planer or 3d patch which uses the symmetr condition.No angular limit between the two symmetry planes in axi symmetry problems as in wedge case.
Hagen_Poiseuille Law :
In nonideal fluid dynamics,the Hagen Poiseuille equation also known as Poiseuille law,is a physical law that gives the pressure drop in an incompressibe and newtonian fluid in laminar flow flowing through a long cylindrical pipe of constant cross section. The assumptions of the equations are that the fluid is incompressible and Newtonian and that the flow is laminar through a pipe of constant circular cross section that is substantially longer than it's diameter and there is no
accceleration of flid in the pipe.H-P Law is valid only in hydrodynamically fully developed region i.e only after the
entrance region. In this project we first need to generate the blockmesh dict file usinf matlab,which is later used by openfoam to generate the mesh.As the flow is axis symmetric only the reduced model is considered instead of the whole pipe,thus reducing computational efforts extensivelly.We will compare the results of analytical solution of Hagen Poiseuille's law with the numerical solution obtained using the wedge and symmetry boundary conditions to assure the fidelity of our numerical results.
Matlab Code For blockMeshDict File For Symmetry Boundary Condition :
Main Code -
% Main Code For 10 Degrees Symmetry Angle
close all
clear all
clc
L = 5; %Length of pipe
theta = 10; %Angle of symmetry wrt axis
r = 0.015; %Radius of pipe
rho = 997; %Density of working fluid
Re =2100; %Reynolds number
sf = 1; %Scaling factor
sgx = 0.3; %Simple grading along X-direction
sgy = 1; %Simple grading along Y-direction
sgz = 1; %Simple grading along Z-direction
dx = 20; %Number of cells in X-direction
dy = 1; %Number of cells in Y-direction
dz = 500; %Number of cells in Z-direction
blockMeshdict_script_function_symmetry_10_degrees(sf,L,r,theta,dx,dy,dz,sgx,sgy,sgz)Function -
function blockMeshdict_script_function_symmetry_10_degrees(sf,L,r,theta,dx,dy,dz,sgx,sgy,sgz ...
)
x1 = r*cosd(theta);
y1 = r*sind(theta);
z1 = 0;
x2 = r*cosd(360-theta);
y2 = r*sind(360-theta);
z2 = L;
% Header lines
line1='/*--------------------------------*-C++ -*--------------------------------*';
line2=' ========= | ';
line3=' / F ield | OpenFOAM: The Open Source CFD Toolbox ';
line4=' / O peration | Website: https://openfoam.org ';
line5=' / A nd | Version: 8 ';
line6=' / M anipulation | ';
line7='*-------------------------------------------------------------------------*/';
line8='FoamFile';
line9='// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //';
% BlockMeshDict
f=fopen('blockMeshDict','w');
fprintf(f,'%s\n',line1,line2,line3,line4,line5,line6,line7,line8);
fprintf(f,'{\n');
fprintf(f,'\tversion\t\t 2.0;\n\tformat \t\tascii;\n\tclass\t\tdictionary;\n\tobject\t\tblockMeshDict;\n}\n');
% Starting Line
fprintf(f,'%s',line9);
% Scaling Factor
fprintf(f,'\n\nconvertToMeters %d;\n\n',sf);
% Vertices
fprintf(f,'vertices\n(\n\n\t(0 0 0)\n\t(%f %f %f)\n\t(%f %f %f)\n\t(0 0 %f)\n\t(%f %f %f)\n\t(%f %f %f)\n);\n',x1,y1,z1,x1,y1,z2,z2,x2,y2,z1,x2,y2,z2);
% Blocks
fprintf(f,'\nblocks\n(\n\thex (0 4 1 0 3 5 2 3) (%d %d %d) simpleGrading (%f %f %f)\n);\n',dx,dy,dz,sgx,sgy,sgz);
% Edges
fprintf(f,'\nedges\n(\n\tarc 1 4 (%f %f %f)\n\tarc 2 5 (%f %f %f)\n);\n',x1/2,y1/2,z1,x2/2,y2/2,z2);
% Boundaries
fprintf(f,'\nboundary\n(\n\n\taxis\n\t{\n\t\ttype empty;\n\t\tfaces\n\t\t(\n\t\t\t(0 3 3 0)\n\t\t);\n\t}\n');
fprintf(f,'\tinlet\n\t{\n\t\ttype patch;\n\t\tfaces\n\t\t(\n\t\t\t(0 4 1 0)\n\t\t);\n\t}\n');
fprintf(f,'\toutlet\n\t{\n\t\ttype patch;\n\t\tfaces\n\t\t(\n\t\t\t(3 5 2 3)\n\t\t);\n\t}\n');
fprintf(f,'\tleftwall\n\t{\n\t\ttype symmetry;\n\t\tfaces\n\t\t(\n\t\t\t(4 0 3 5)\n\t\t);\n\t}\n');
fprintf(f,'\trightwall\n\t{\n\t\ttype symmetry;\n\t\tfaces\n\t\t(\n\t\t\t(0 1 2 3)\n\t\t);\n\t}\n');
fprintf(f,'\touterwall\n\t{\n\t\ttype wall;\n\t\tfaces\n\t\t(\n\t\t\t(1 4 5 2)\n\t\t);\n\t}\n\n\n);\n');
% Patch Pairs
fprintf(f,'\nmergePatchPairs\n(\n);\n');
% Ending Line
fprintf(f,'\n%s',line9);
fclose(f);
endGenerated blockMeshDict file :
For Symmetry Angle = `10^o`
/*--------------------------------*-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 0 0)
(0.014772 0.002605 0.000000)
(0.014772 0.002605 5.000000)
(0 0 5.000000)
(0.014772 -0.002605 0.000000)
(0.014772 -0.002605 5.000000)
);
blocks
(
hex (0 4 1 0 3 5 2 3) (20 1 500) simpleGrading (0.300000 1.000000 1.000000)
);
edges
(
arc 1 4 (0.007386 0.001302 0.000000)
arc 2 5 (0.007386 -0.001302 5.000000)
);
boundary
(
axis
{
type empty;
faces
(
(0 3 3 0)
);
}
inlet
{
type patch;
faces
(
(0 4 1 0)
);
}
outlet
{
type patch;
faces
(
(3 5 2 3)
);
}
leftwall
{
type symmetry;
faces
(
(4 0 3 5)
);
}
rightwall
{
type symmetry;
faces
(
(0 1 2 3)
);
}
outerwall
{
type wall;
faces
(
(1 4 5 2)
);
}
);
mergePatchPairs
(
);
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //Simillarly ,blockMeshDict files are generated for `25^o` & `45^o` symmetry angle boundary condition.Also proper changes are done in velocity and pressure file.
U file :
/*--------------------------------*- 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 volVectorField;
object U;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dimensions [0 1 -1 0 0 0 0];
internalField uniform (0 0 0);
boundaryField
{
outlet
{
type zeroGradient;
}
inlet
{
type fixedValue;
value uniform (0 0 0.070000);
}
axis
{
type empty;
}
leftwall
{
type symmetry;
}
rightwall
{
type symmetry;
}
outerwall
{
type noSlip;
}
}
// ************************************************************************* //
P file :
/*--------------------------------*- 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 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;
}
leftwall
{
type symmetry;
}
rightwall
{
type symmetry;
}
outerwall
{
type zeroGradient;
}
}
// ************************************************************************* //Transport Properties :
/*--------------------------------*- 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;
location "constant";
object transportProperties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
nu [0 2 -1 0 0 0 0] 0.000001;
// ************************************************************************* //
The Pressure And Velocity files will stay same for all the wedge and symmetry boundary conditions.
Control Directory :
It take more iteration time in `25^o and 45^o` symmetry angle for velocity profile to propagate than `10^o` symmetry angle.Thus simulations at `25^o and 45^o` were run for longer .
Control Directory for symmetry angle = `10^o`
/*--------------------------------*- 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;
location "system";
object controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
application icoFoam;
startFrom startTime;
startTime 0;
stopAt endTime;
endTime 80;
deltaT 0.011;
writeControl timeStep;
writeInterval 50;
purgeWrite 0;
writeFormat ascii;
writePrecision 6;
writeCompression off;
timeFormat general;
timePrecision 6;
runTimeModifiable true;
// ************************************************************************* //
Control Directory for symmetry angle = `25^o and 45^o`
/*--------------------------------*- 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;
location "system";
object controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
application icoFoam;
startFrom startTime;
startTime 0;
stopAt endTime;
endTime 200;
deltaT 0.011;
writeControl timeStep;
writeInterval 50;
purgeWrite 0;
writeFormat ascii;
writePrecision 6;
writeCompression off;
timeFormat general;
timePrecision 6;
runTimeModifiable true;
// ************************************************************************* //
Results :
1). For Wedge and Symmetry Angle = `10^o`
Velocity profile at 3.15 m for wedge and symmetry angle of `10^O`
. 
Velocity profile at 4 m for wedge and symmetry angle of `10^o`
Velocity profile at 5 m for wedge and symmetry angle of `10^o`
Pressure drop at Wedge and Symmetry Angle = `10^o`
Data Extracted Using Location Probing in Paraview :
2). For Wedge and Symmetry Angle = `25^o`
Velocity profile at 3.15 m for wedge and symmetry angle of `25^O`
Velocity profile at 4 m for wedge and symmetry angle of `25^O`
Velocity profile at 5 m for wedge and symmetry angle of `25^O`
Pressure Drop At Wedge And Symmetry Angle = `25^o`
Data Extracted Using Location Probing in Paraview :
2). For Wedge and Symmetry Angle = `45^o`
Velocity profile at 3.15 m for symmetry angle of `45^O` Velocity profile at 4 m for symmetry angle of `45^O`
Velocity Profile at 5 m for symmetry angle = `45^o`
Wall Shear Stress At Fully Developed Region Near Pipe Walls:
For Wedge Angle = `45^o`
At the angle of `45^o` openFoam has shown the error saying that the angle should be less than `2.5^o`
Compairing Symmetry BC vs Wedge BC in OpenFOAM vs Analytical H.P equation -
| Angle | U | P | Tau At Fully Developed Region | Error(`DeltaU`) | Error(`DeltaP`) | Error(`Deltatau`) | |
| HP | Analytical | 0.14 | 0.0046 | 0.0186 | |||
| `10^o` | 0.137 | 0.0047 | 0.016 | 0.003 | -0.0001 | 0.0026 | |
| SYMMETRY | `25^o` | 0.138 | 0.0056 | 0.0174 | 0.002 | -0.001 | 0.0012 |
| `45^o` | 0.139 | 0.0091 | 0.022 | 0.001 | -0.0045 | -0.0034 | |
| `10^o` | 0.137 | 0.00478 | 0.0160 | 0.003 | -0.00018 | 0.0026 | |
| WEDGE | `25^o` | 0.138 | 0.0056 | 0.0176 | 0.002 | -0.001 | 0.001 |
| `45^o` | --- | --- | --- | --- | --- | --- |
NOTE :- The discription of geometry and formulae to calculate solutions can be found on the following link : https://skill-lync.com/student-projects/week-11-simulation-of-flow-through-a-pipe-in-openfoam-236
Observations And Conclusions :
1). Aalytical hydrodynamic entrance length is validated with velocity profile obtained at 3.15m and numerical results extracted using probe location in paraview.From the velocity plots we can see that velocity profile is fully parabolic in nature with max.velocity close to 0.14m/s.
2). Numerical values extracted for `V_(max),DeltaP and tau` are closely matched with analytical Hagen-Poiseuilles values.
3). At the angle of `45^o` openFoam has shown the error saying that the angle should be less than `2.5^o`
4). The difference in error at everyangle is almost equal,which is satisfying the axi-symmetry property in case of symmetry boundary condition.
5). Thus Hagen-Poiseuilles flow equations for constant circular cross section are verified and are compared using symmetry and wedge boundary condition at different angles.