Week 12 - Symmetry vs Wedge vs HP equation
Aim:- Symmetry vs Wedge vs HP equation
Objective:-
To Write a Matlab program that takes an angle as input and generates a blockMesh file for the given angle.
Angles to test 10,25,45. Compare the above results and discuss your findings.
Theory:-
Hagen Poiseuille's equation:
In nonideal fluid dynamics, the Hagen–Poiseuille equation, also known as the Hagen–Poiseuille law, Poiseuille law or Poiseuille equation, is a physical law that gives the pressure drop in an incompressible and Newtoian fluid in laminar flow flowing through a long cylindrical pipe of constant cross section.
The pressure drop is given by:
ΔP=32μLVD2ΔP=32μLVD2
For a pipe with laminar flow there are basically two regions. Entrance length developing region where the flow is not fully developed and the laminar flow layers are at equal velocity with respect to each other whereas in fully deveoped region the flow is fully developed and follows the hagen-Poiseuille law.
Assumed parameters for simulation:
1) Geometry consideration - axi-symmeteric
2) Wedge angle(θθ) = 4 degrees
3) Working fluid = water (density - 997kgm3kgm3
4) Temperature of working fluid = 25°C
5) Reynold's number(Re) = 2100
6) Dynamic viscosity of working fluid(μμ) = 8.9 X10−410-4 Pa.sec
7) Kinematic viscosity of working fluid(νν) = 8.926 X10−710-7 m2secm2sec
8) Diameter of pipe(D) = 0.02m
Analytical solution:
1) Average velocity:
Re=Vavg×D×ρμRe=Vavg×D×ρμ
⇒Vavg=8.9×10−4×2100997×0.02⇒Vavg=8.9×10-4×2100997×0.02
⇒Vavg=0.0937msec⇒Vavg=0.0937msec
2) Maximum velocity:
Vmax=2×VavgVmax=2×Vavg
Vmax=2×0.0937=0.187msecVmax=2×0.0937=0.187msec
3) Hydrodynamic Entry Length:
Lh=0.06×Re×DLh=0.06×Re×D
Lh=0.06×2100×0.02=2.52mLh=0.06×2100×0.02=2.52m
This means that minimum length of pipe should be 2.52m. The total length of pipe can be assumed as 2.7m.
4) Pressure Drop:
ΔP=32μLVavgD2ΔP=32μLVavgD2
ΔP=32×8.9×10−4×2.7×0.0930.022=17.878PaΔP=32×8.9×10-4×2.7×0.0930.022=17.878Pa
5) Kinematic pressure drop:
Pk=ΔPρ=17.878997=0.0179m
Matlab program to generate blockMeshDict file for symmetry boundary condition:
clear all
close all
clc
d = 0.02; %diameter
r = d/2;
L = 2.7;
theta = input('Enter the value of wedge angle =');
x = input('Enter the grading factor along x =');
y = input('Enter the grading factor along y =');
z = input('Enter the grading factor along z =');
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)];
A1= '/*--------------------------------*- C++ -*----------------------------------*\';
A2= ' ========= | ';
A3= ' \\ / F ield | OpenFOAM: The Open Source CFD Toolbox ';
A4= ' \\ / O peration | Website: https://openfoam.org ';
A5= ' \\ / A nd | Version: 8 ';
A6= ' \\/ M anipulation | ';
A7= '\*---------------------------------------------------------------------------*/';
A8= 'FoamFile ';
A9= '{ ';
A10=' version 2.0; ';
A11=' format ascii; ';
A12=' class dictionary; ';
A13=' object blockMeshDict; ';
A14='} ';
A15='// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //';
f1 = fopen('blockMeshDict.txt','wt');
fprintf(f1,'%s\n',A1);
fprintf(f1,'%s\n',A2);
fprintf(f1,'%s\n',A3);
fprintf(f1,'%s\n',A4);
fprintf(f1,'%s\n',A5);
fprintf(f1,'%s\n',A6);
fprintf(f1,'%s\n',A7);
fprintf(f1,'%s\n',A8);
fprintf(f1,'%s\n',A9);
fprintf(f1,'%s\n',A10);
fprintf(f1,'%s\n',A11);
fprintf(f1,'%s\n',A12);
fprintf(f1,'%s\n',A13);
fprintf(f1,'%s\n',A14);
fprintf(f1,'%s\n',A15);
fprintf(f1,'\n');
fprintf(f1,'%s\n','convertToMeters 1;');
fprintf(f1,'\n');
fprintf(f1,'%s\n','vertices');
fprintf(f1,'(\n');
fprintf(f1,'\t(%f %f %f)\r',V0(1),V0(2),V0(3));
fprintf(f1,'\t(%f %f %f)\r',V1(1),V1(2),V1(3));
fprintf(f1,'\t(%f %f %f)\r',V2(1),V2(2),V2(3));
fprintf(f1,'\t(%f %f %f)\r',V3(1),V3(2),V3(3));
fprintf(f1,'\t(%f %f %f)\r',V4(1),V4(2),V4(3));
fprintf(f1,'\t(%f %f %f)\r',V5(1),V5(2),V5(3));
fprintf(f1,');\r\n');
fprintf(f1,'blocks\r(\r');
fprintf(f1,'\t hex (0 3 4 1 0 3 5 2) (300 50 1) simpleGrading (%f %f %f)\r',x,y,z);
fprintf(f1,');\r\n');
fprintf(f1,'edges\r(\r');
fprintf(f1,'\t\t arc 1 2 (%f %f %f)\r',0,r,0);
fprintf(f1,'\t\t arc 4 5 (%f %f %f)\r',L,r,0);
fprintf(f1,');\r\n');
fprintf(f1,'boundary\r(\r');
fprintf(f1,'\tinlet\r\t{\r');
fprintf(f1,'\t\ttype patch;\r');
fprintf(f1,'\t\tfaces\r\t\t(\r');
fprintf(f1,'\t\t\t(0 2 1 0)\r');
fprintf(f1,'\t\t);\r');
fprintf(f1,'\t}\r');
fprintf(f1,'\toutlet\r\t{\r');
fprintf(f1,'\t\ttype patch;\r');
fprintf(f1,'\t\tfaces\r\t\t(\r');
fprintf(f1,'\t\t\t(3 4 5 3)\r');
fprintf(f1,'\t\t);\r');
fprintf(f1,'\t}\r');
fprintf(f1,'\ttop\r\t{\r');
fprintf(f1,'\t\ttype wall;\r');
fprintf(f1,'\t\tfaces\r\t\t(\r');
fprintf(f1,'\t\t\t(1 2 5 4)\r');
fprintf(f1,'\t\t);\r');
fprintf(f1,'\t}\r');
fprintf(f1,'\tback\r\t{\r');
fprintf(f1,'\t\ttype symmetry;\r');
fprintf(f1,'\t\tfaces\r\t\t(\r');
fprintf(f1,'\t\t\t(0 1 4 3)\r');
fprintf(f1,'\t\t);\r');
fprintf(f1,'\t}\r');
fprintf(f1,'\tfront\r\t{\r');
fprintf(f1,'\t\ttype symmetry;\r');
fprintf(f1,'\t\tfaces\r\t\t(\r');
fprintf(f1,'\t\t\t(0 3 5 2)\r');
fprintf(f1,'\t\t);\r');
fprintf(f1,'\t}\r');
fprintf(f1,'\taxis\r\t{\r');
fprintf(f1,'\t\ttype empty;\r');
fprintf(f1,'\t\tfaces\r\t\t(\r');
fprintf(f1,'\t\t\t(0 3 3 0)\r');
fprintf(f1,'\t\t);\r');
fprintf(f1,'\t}\r);\r\n');
fprintf(f1,'mergePatchPairs\n(\n);\n');
fprintf(f1,'\n// ************************************************************************* //');
fclose(f1);BlockMeshDict file for 10 degree wedge with symmetry 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 dictionary;
object blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
convertToMeters 1;
vertices
(
(0.000000 0.000000 0.000000)
(0.000000 0.009962 -0.000872)
(0.000000 0.009962 0.000872)
(2.700000 0.000000 0.000000)
(2.700000 0.009962 -0.000872)
(2.700000 0.009962 0.000872)
);
blocks
(
hex (0 3 4 1 0 3 5 2) (300 50 1) simpleGrading (1.000000 1.000000 1.000000)
);
edges
(
arc 1 2 (0.000000 0.010000 0.000000)
arc 4 5 (2.700000 0.010000 0.000000)
);
boundary
(
inlet
{
type patch;
faces
(
(0 2 1 0)
);
}
outlet
{
type patch;
faces
(
(3 4 5 3)
);
}
top
{
type wall;
faces
(
(1 2 5 4)
);
}
back
{
type symmetry;
faces
(
(0 1 4 3)
);
}
front
{
type symmetry;
faces
(
(0 3 5 2)
);
}
axis
{
type empty;
faces
(
(0 3 3 0)
);
}
);
mergePatchPairs
(
);
// ************************************************************************* //BlockMeshDict file for 25 degree wedge with symmetry 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 dictionary;
object blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
convertToMeters 1;
vertices
(
(0.000000 0.000000 0.000000)
(0.000000 0.009763 -0.002164)
(0.000000 0.009763 0.002164)
(2.700000 0.000000 0.000000)
(2.700000 0.009763 -0.002164)
(2.700000 0.009763 0.002164)
);
blocks
(
hex (0 3 4 1 0 3 5 2) (300 50 1) simpleGrading (1.000000 1.000000 1.000000)
);
edges
(
arc 1 2 (0.000000 0.010000 0.000000)
arc 4 5 (2.700000 0.010000 0.000000)
);
boundary
(
inlet
{
type patch;
faces
(
(0 2 1 0)
);
}
outlet
{
type patch;
faces
(
(3 4 5 3)
);
}
top
{
type wall;
faces
(
(1 2 5 4)
);
}
back
{
type symmetry;
faces
(
(0 1 4 3)
);
}
front
{
type symmetry;
faces
(
(0 3 5 2)
);
}
axis
{
type empty;
faces
(
(0 3 3 0)
);
}
);
mergePatchPairs
(
);
// ************************************************************************* //BlockMeshDict file for 45 degree wedge with symmetry 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 dictionary;
object blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
convertToMeters 1;
vertices
(
(0.000000 0.000000 0.000000)
(0.000000 0.009239 -0.003827)
(0.000000 0.009239 0.003827)
(2.700000 0.000000 0.000000)
(2.700000 0.009239 -0.003827)
(2.700000 0.009239 0.003827)
);
blocks
(
hex (0 3 4 1 0 3 5 2) (300 50 1) simpleGrading (1.000000 1.000000 1.000000)
);
edges
(
arc 1 2 (0.000000 0.010000 0.000000)
arc 4 5 (2.700000 0.010000 0.000000)
);
boundary
(
inlet
{
type patch;
faces
(
(0 2 1 0)
);
}
outlet
{
type patch;
faces
(
(3 4 5 3)
);
}
top
{
type wall;
faces
(
(1 2 5 4)
);
}
back
{
type symmetry;
faces
(
(0 1 4 3)
);
}
front
{
type symmetry;
faces
(
(0 3 5 2)
);
}
axis
{
type empty;
faces
(
(0 3 3 0)
);
}
);
mergePatchPairs
(
);
// ************************************************************************* //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 25;
deltaT 0.01;
writeControl timeStep;
writeInterval 20;
purgeWrite 0;
writeFormat ascii;
writePrecision 6;
writeCompression off;
timeFormat general;
timePrecision 6;
runTimeModifiable true;
// ************************************************************************* //Pressure 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
{
inlet
{
type zeroGradient;
}
outlet
{
type fixedValue;
value uniform 0.0179;
}
top
{
type zeroGradient;
}
front
{
type symmetry;
}
back
{
type symmetry;
}
axis
{
type empty;
}
}
// ************************************************************************* //Velocity 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
{
inlet
{
type fixedValue;
value uniform (0.0937 0 0);
}
outlet
{
type zeroGradient;
}
top
{
type noSlip;
}
front
{
type symmetry;
}
back
{
type symmetry;
}
axis
{
type empty;
}
}
// ************************************************************************* //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;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
nu [0 2 -1 0 0 0 0] 0.0000008927;
// ************************************************************************* //Mesh:
1) At 10 degrees
2) At 25 degrees
3) At 45 degrees
velocity distance graph :
1] For a pipe of wedge angle 10 degrees.
i) At distance of 0.001m from inlet.
ii) At distance of 0.05m from inlet.
iii) At distance of 0.15m from inlet
iv) At outlet
v) Shear stress graph
2] for pipe of edge angle 25 dgrees.
i) At distance of 0.001m from inlet
ii) At distance of 0.05m from inlet
iii) At distance of 0.15m from inlet
iv) At oulet
v) Shear stress graph
3] For a pipe of wedge angle 45
i) At distance of 0.001m from inlet
ii) At distance of 0.05m from inlet
iii) At distance of 0.15m from inlet
iv) At outlet
v) Shear stress graph
Comparison:
Wedge angle Max velocity(H.P equation) Max velocity(simulation) Max Shear stress
1) 4 degree(wedge) 0.187 0.172 0.323
2) 10 degree(symmetry) 0.187 0.171 0.325
3) 25 degree(symmetry) 0.187 0.186 0.332
4) 45 degree(symmetry) 0.187 0.187 0.344
Conclusion:-
1) Here we can see as the wedge angle increases the velocity comes closer to the analytical solution.
2) Symmetry BC is better than the wedge BC as it gives accurate results.
3)Symmetry BC takes more computational time as compared to wedge BC.