Symmetry vs Wedge vs HP equation
Question :
You should have results for the pipe flow case with the wedge boundary condition.
Now, your task is to make the simulation work with Symmetry Boundary Condition.
- You should have results obtained from solving the Hagen Poisuelle's equation
- You should have results obtained from simulating the pipe flow with the wedge BC
- You should obtain results for the pipe flow with symmetry BC
- 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.
Answer :
Given :-
- Wedge angles 10,25 & 45 degree.
1] Assumptions:-
- D = 0.01m .....Diameter
- R = 0.005m .....Radius
- mu = 1.002e-3 .....In SI units (water at 20 degree celsius)
- rho = 998.21 .....In SI units (water at 20 degree celsius)
- Re = 2100 .....Reynolds number
- Assumed Temperature of water = 20 degree C
- Here we are calculating results for various angles which are 10,25 and 45. degree.
When the fluid enters the pipe its boundary layer increases till the centre the pipe and fills the pipe,this region is entrance region where the flow is developing gradually and the remaining length is fully developed length.
Entrance length is a length which a flow travels before it becomes fully developed flow. When a boundary layer expands to fill pipe, the flow starts to develop. In the the above figure we can see that the flow is developing in the Entrance region and after the entrance region it gets fully developed and thereafter there are no changes in the flow characteristics. Viscous forces are resposible for the velocity profiles.
Entrance Length of Laminar flow pipe is calculated by (Le/D) = (0.06*Re)
2] Analytical values:-
- L_e = 0.06*Re*D ..... Entrance Length
- L_t = L_e + 1.24 ..... Total Length
- v_avg = (Re*mu)/(rho*D) ..... Velocity
- v_max = 2*v_avg ...... Max Velocity
- Pressure Drop = ( 32*mu*v_avg*L_t ) / (D^2)
- Kinematic pressure = Pressure drop/rho.
Below are the matlab codes for the blockmeshDict file and blockMeshDict file for all given angles and also updated p, U, controlDict, Transportproperties files.
3] Code
- Matlab code for BlockMeshDict file for any angle
clear all
close all
clc
% Inputs
Re = 2100;
D = 0.01;
R = D/2;
theta = input('Enter value for Wedge angle =');
mu = 1.002e-3;
rho = 998.21;
nu = mu/rho;
% Analytical Values
L_e = 0.06*Re*D;
L_t = L_e + 1.24;
v_avg = (Re*nu)/D;
v_max = 2*v_avg;
Delta_P = (32*mu*v_avg*L_t)/D^2;
kinematic_P = Delta_P/rho;
% For finding Vertices
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_t 0 0];
v4 = [L_t R*cosd(theta/2) -R*sind(theta/2)];
v5 = [L_t R*cosd(theta/2) R*sind(theta/2)];
% Script for blockmeshDict Directory
f1 = fopen('blockMeshDict.txt','w')
fprintf(f1,'%s','/*--------------------------------*- C++ -*----------------------------------*\');
fprintf(f1,'\n========= |\n');
fprintf(f1,' \\ / F ield | OpenFOAM: The Open Source CFD Toolbox\n');
fprintf(f1,' \\ / O peration | Website: https://openfoam.org\n');
fprintf(f1,' \\ / A nd | Version: 7\n');
fprintf(f1,' \\/ M anipulation |\n');
fprintf(f1,'%s','\*---------------------------------------------------------------------------*/');
fprintf(f1,'\nFoamFile\n');
b1 = blanks(4);
b2 = blanks(5);
b3 = blanks(6);
b4 = blanks(7);
fprintf(f1,'{\n');
fprintf(f1,'%s',' version',b2,'2.0;');
fprintf(f1,'\n');
fprintf(f1,'%s',' format',b3,'ascii;');
fprintf(f1,'\n');
fprintf(f1,'%s',' class',b4,'dictionary;');
fprintf(f1,'\n');
fprintf(f1,'%s',' object',b3,'blockMeshDict;');
fprintf(f1,'\n');
fprintf(f1,'}\n');
fprintf(f1,'%s','// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //');
fprintf(f1,'\n\nconvertToMeters 1;\n\n');
% for Vertices
fprintf(f1,'%s','vertices');
fprintf(f1,'\n');
fprintf(f1,'(\n');
fprintf(f1,b1,'%s','(');
v_0 = fprintf(f1,'(%f %f %f)\n',v0(1),v0(2),v0(3));
fprintf(f1,b1,'%s','(');
v_1 = fprintf(f1,'(%f %f %f)\n',v1(1),v1(2),v1(3));
fprintf(f1,b1,'%s','(');
v_2 = fprintf(f1,'(%f %f %f)\n',v2(1),v2(2),v2(3));
fprintf(f1,b1,'%s','(');
v_3 = fprintf(f1,'(%f %f %f)\n',v3(1),v3(2),v3(3));
fprintf(f1,b1,'%s','(');
v_4 = fprintf(f1,'(%f %f %f)\n',v4(1),v4(2),v4(3));
fprintf(f1,b1,'%s','(');
v_5 = fprintf(f1,'(%f %f %f)\n',v5(1),v5(2),v5(3));
fprintf(f1,'\n);\n');
fprintf(f1,'\n');
fprintf(f1,'%s','blocks');
fprintf(f1,'\n');
fprintf(f1,'(\n');
fprintf(f1,b1,'%s','(');
fprintf(f1,'%s','hex (0 1 2 0 3 4 5 3)');
% Mesh
m1 = input('Enter no.of mesh elements along X-axis =');
m2 = input('Enter no.of mesh elements along Y-axis =');
m3 = input('Enter no.of mesh elements along Z-axis =');
fprintf(f1,blanks(1),'%s','(');
fprintf(f1,'(%2.0f %1.0f %3.0f)', m1, m2, m3);
fprintf(f1,'%s',' simpleGrading (1 1 1)');
fprintf(f1,'\n);');
% edges
fprintf(f1,'\n');
fprintf(f1,'\nedges\n');
fprintf(f1,'(\n');
fprintf(f1,'%s arc 1 2 (0 %f 0)',blanks(1),R);
fprintf(f1,'\n%s arc 4 5 (%f %f 0)',blanks(1),L_t,R);
fprintf(f1,'\n);\n\n');
% boundary conditions
% For Inlet
fprintf(f1,'boundary\n');
fprintf(f1,'(\n');
fprintf(f1,'%s',' inlet');
fprintf(f1,'\n');
fprintf(f1,' {\n');
fprintf(f1,'%s',blanks(10),'type patch;');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'faces');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'(');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(13),'(0 1 2 0)');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),');');
fprintf(f1,'\n }\n');
% for Outlet
fprintf(f1,'%s',' outlet');
fprintf(f1,'\n');
fprintf(f1,' {\n');
fprintf(f1,'%s',blanks(10),'type patch;');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'faces');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'(');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(13),'(3 5 4 3)');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),');');
fprintf(f1,'\n }\n');
% for front
fprintf(f1,'%s',' front');
fprintf(f1,'\n');
fprintf(f1,' {\n');
fprintf(f1,'%s',blanks(10),'type symmetry;');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'faces');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'(');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(13),'(1 0 3 4)');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),');');
fprintf(f1,'\n }\n');
% for back
fprintf(f1,'%s',' back');
fprintf(f1,'\n');
fprintf(f1,' {\n');
fprintf(f1,'%s',blanks(10),'type symmetry;');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'faces');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'(');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(13),'(2 5 3 0)');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),');');
fprintf(f1,'\n }\n');
% for fixedWalls
fprintf(f1,'%s',' fixedWalls');
fprintf(f1,'\n');
fprintf(f1,' {\n');
fprintf(f1,'%s',blanks(10),'type wall;');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'faces');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'(');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(13),'(2 1 4 5)');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),');');
fprintf(f1,'\n }\n');
% for axis
fprintf(f1,'%s',' axis');
fprintf(f1,'\n');
fprintf(f1,' {\n');
fprintf(f1,'%s',blanks(10),'type empty;');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'faces');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),'(');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(13),'(0 3 3 0)');
fprintf(f1,'\n');
fprintf(f1,'%s',blanks(10),');');
fprintf(f1,'\n }\n');
fprintf(f1,');\n');
fprintf(f1,'\nmergePatchPairs');
fprintf(f1,'\n(\n);\n\n');
fprintf(f1,'%s','// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //');
fclose(f1)
- blockMeshDict of 10 degree Wedge angle
/*--------------------------------*- C++ -*----------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Version: 7
\\/ M anipulation |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
convertToMeters 1;
vertices
(
(0.000000 0.000000 0.000000)
(0.000000 0.004981 -0.000436)
(0.000000 0.004981 0.000436)
(2.500000 0.000000 0.000000)
(2.500000 0.004981 -0.000436)
(2.500000 0.004981 0.000436)
);
blocks
(
hex (0 1 2 0 3 4 5 3) (50 1 100) simpleGrading (1 1 1)
);
edges
(
arc 1 2 (0 0.005000 0)
arc 4 5 (2.500000 0.005000 0)
);
boundary
(
inlet
{
type patch;
faces
(
(0 1 2 0)
);
}
outlet
{
type patch;
faces
(
(3 5 4 3)
);
}
front
{
type symmetry;
faces
(
(1 0 3 4)
);
}
back
{
type symmetry;
faces
(
(2 5 3 0)
);
}
fixedWalls
{
type wall;
faces
(
(2 1 4 5)
);
}
axis
{
type empty;
faces
(
(0 3 3 0)
);
}
);
mergePatchPairs
(
);
// ************************************************************************* //- blockMeshDict of 25 degree Wedge angle
/*--------------------------------*- C++ -*----------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Version: 7
\\/ M anipulation |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
convertToMeters 1;
vertices
(
(0.000000 0.000000 0.000000)
(0.000000 0.004881 -0.001082)
(0.000000 0.004881 0.001082)
(2.500000 0.000000 0.000000)
(2.500000 0.004881 -0.001082)
(2.500000 0.004881 0.001082)
);
blocks
(
hex (0 1 2 0 3 4 5 3) (50 1 100) simpleGrading (1 1 1)
);
edges
(
arc 1 2 (0 0.005000 0)
arc 4 5 (2.500000 0.005000 0)
);
boundary
(
inlet
{
type patch;
faces
(
(0 1 2 0)
);
}
outlet
{
type patch;
faces
(
(3 5 4 3)
);
}
front
{
type symmetry;
faces
(
(1 0 3 4)
);
}
back
{
type symmetry;
faces
(
(2 5 3 0)
);
}
fixedWalls
{
type wall;
faces
(
(2 1 4 5)
);
}
axis
{
type empty;
faces
(
(0 3 3 0)
);
}
);
mergePatchPairs
(
);
// ************************************************************************* //- blockMeshDict of 45 degree Wedge angle
/*--------------------------------*- C++ -*----------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Version: 7
\\/ M anipulation |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
convertToMeters 1;
vertices
(
(0.000000 0.000000 0.000000)
(0.000000 0.004619 -0.001913)
(0.000000 0.004619 0.001913)
(2.500000 0.000000 0.000000)
(2.500000 0.004619 -0.001913)
(2.500000 0.004619 0.001913)
);
blocks
(
hex (0 1 2 0 3 4 5 3) (50 1 100) simpleGrading (1 1 1)
);
edges
(
arc 1 2 (0 0.005000 0)
arc 4 5 (2.500000 0.005000 0)
);
boundary
(
inlet
{
type patch;
faces
(
(0 1 2 0)
);
}
outlet
{
type patch;
faces
(
(3 5 4 3)
);
}
front
{
type symmetry;
faces
(
(1 0 3 4)
);
}
back
{
type symmetry;
faces
(
(2 5 3 0)
);
}
fixedWalls
{
type wall;
faces
(
(2 1 4 5)
);
}
axis
{
type empty;
faces
(
(0 3 3 0)
);
}
);
mergePatchPairs
(
);
// ************************************************************************* //Now all remaining files like controlDict, p, U, and transport properties will be same for all angles and are shown below.
- ContolDict
/*--------------------------------*- C++ -*----------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Version: 7
\\/ M anipulation |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "system";
object controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
application icoFoam;
startFrom startTime;
startTime 0;
stopAt endTime;
endTime 0.5;
deltaT 0.005;
writeControl timeStep;
writeInterval 20;
purgeWrite 0;
writeFormat ascii;
writePrecision 6;
writeCompression off;
timeFormat general;
timePrecision 6;
runTimeModifiable true;
// ************************************************************************* //- p
/*--------------------------------*- C++ -*----------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Version: 7
\\/ 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.1693;
}
front
{
type symmetry;
}
back
{
type symmetry;
}
fixedWalls
{
type zeroGradient;
}
axis
{
type empty;
}
}
// ************************************************************************* //- U
/*--------------------------------*- C++ -*----------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Version: 7
\\/ 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.2108 0 0);
}
outlet
{
type zeroGradient;
}
front
{
type symmetry;
}
back
{
type symmetry;
}
fixedWalls
{
type noSlip;
}
axis
{
type empty;
}
}
// ************************************************************************* //- TransportProperties
/*--------------------------------*- C++ -*----------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Version: 7
\\/ M anipulation |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "constant";
object transportProperties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
nu [0 2 -1 0 0 0 0] 1.00379e-6;
// ************************************************************************* //2] Outputs
- For 10 degree
i) Mesh Output
ii) Velocity Distribution at inlet
iii) Velocity Distribution at Outlet
iv) Velocity profile at 0.01m from inlet
v) Velocity profile at 0.05m from inlet
vi) Velocity profile at 1m from inlet
vii) Velocity profile at 2.5m from inlet
viii) Pressure Distribution
ix) Pressure Profile
- for 25 degree
i) Mesh Output
ii) Velocity Distribution at inlet
iii) Velocity Distribution at Outlet
iv) Velocity profile at 0.01m from inlet
v) Velocity profile at 0.05m from inlet
vi) Velocity profile at 1m from inlet
vii) Velocity profile at 2.5m from inlet
viii) Pressure Distribution
ix) Pressure Profile
- For 45 degree
i) Mesh Output
ii) Velocity Distribution at inlet
iii) Velocity Distribution at Outlet
iv) Velocity profile at 0.01m from inlet
v) Velocity profile at 0.05m from inlet
vi) Velocity profile at 1m from inlet
vii) Velocity profile at 2.5m from inlet
viii) Pressure Distribution
ix) Pressure Profile
Result :
|
Sr. No. |
Angle (deg.) |
Max Velocity (m/s) |
|
1. |
4 |
0.292585 |
|
2. |
10 |
0.293146 |
|
3. |
25 |
0.295000 |
|
4. |
45 |
0.300778 |
Conclusion :
- We can see that the Velocity profiles get developed same as Wedge Bc solutions.
- From the above results we can conclude that as the wedge angle increases, max velocity also increases making the velocity to get closer to the analytical velocity value.
- As the angle increases max velocity value gets close to analytical results.
- Percentage error decreases with increase in wedge angle for symmetry boundary conditions.
- Also the simulation time required increases with increase in wedge angle.
- from the above results we ca see that the pressure Distribution is almost same for all the angles with a slight changes in the values.
- Symmetry Bc are better than wedge BC.