Challenge

Week 12 - Validation studies of Symmetry BC vs Wedge BC in OpenFOAM vs Analytical H.P equation

Introduction:

This report is in continuation of the last project, flow through an pipe. In this case study, symmetry boundary condition is used for front and back faces.The wedge angles used in the case study are 10 degree, 25 degree and 45 degree.

Boundary conditions:

The boundary conditions used for all the three cases are :

Inlet Outlet Top walls front back axis
Patch patch wall symmetry symmetry empty

For wedge angle 45 degrees:

Matlab code for the wedge angle 45 degrees:

clear all
close all
clc
%flow parameters

Re =2100;%reynolds number
D=0.004;%diameter of the pipe
r=D/2;%radius of the pipe
l=0.06*Re*D;%hydrodynamic entrance length
l_p=l+0.7;%length of the pipe
theta=45;%wedge angle
mu=0.89e-03;%dynamic viscosity
rho=997;%density of fluid
v_avg=(Re*mu)/(rho*D);%mean velocity
v_max=2*v_avg;%maximum velocity
del_p=(8*mu*l_p*v_avg)/(r^2);%pressure drop
kin_p=(del_p)/(rho);%kinematic pressure drop
%Writing the blockMeshDict file
f1=fopen('blockMeshdict','w');
%defining the vertioces
fprintf(f1, 'vertices\n');
fprintf(f1, '(\n');
fprintf(f1, '\t (%d %d %d)\n',0,0,0);
fprintf(f1, '\t (%d %d %d)\n',0,r*cosd(theta/2),r*sind(theta/2));
fprintf(f1, '\t (%d %d %d)\n',0,r*cosd(theta/2),-r*sind(theta/2));
fprintf(f1, '\t (%d %d %d)\n',l,0,0);
fprintf(f1, '\t (%d %d %d)\n',l,r*cosd(theta/2),r*sind(theta/2));
fprintf(f1, '\t (%d %d %d)\n',l,r*cosd(theta/2),-r*sind(theta/2));
fprintf(f1, ');\n\n');
fprintf(f1, 'blocks\n');
fprintf(f1, '(\n');
fprintf(f1, '\t hex(0 3 5 2 0 3 4 1)(%d %d %d)',300,30,1);
fprintf(f1, 'SimpleGrading (1 0.1 1)\n\n');
fprintf(f1, ');\n\n');

fprintf(f1, 'edges\n');
fprintf(f1, '(\n');
fprintf(f1, '\t arc %d %d (%d%d%d)\n',1,2,0,r,0);
fprintf(f1, '\t arc %d %d(%d %d %d)\n',4,5,l_p,r,0);
fprintf(f1, ');\n\n');

%boundary conditions
fprintf(f1, 'boundary\n');
fprintf(f1, '(\n');
fprintf(f1, '\t inlet \n');
fprintf(f1, '\t {\n');
fprintf(f1, '\t {\n');
fprintf(f1, '\t\t type patch;\n');
fprintf(f1, '\t\t faces \n');
fprintf(f1, '\t\t\t (0 1 2 0)\n');
fprintf(f1, '\t\t); \n');
fprintf(f1, '\t }\n');

fprintf(f1, '\t outlet\n');
fprintf(f1, '\t {\n');
fprintf(f1, '\t\t type patch ;\n');
fprintf(f1, '\t\t faces\n');
fprintf(f1, '\t\t (\n');
fprintf(f1, '\t\t\t (3 5 4 3)\n');
fprintf(f1, '\t\t); \n');
fprintf(f1, '\t }\n');

fprintf(f1, '\t top\n');
fprintf(f1, '\t {\n');
fprintf(f1, '\t\t type wall;\n');
fprintf(f1, '\t\t faces\n');
fprintf(f1, '\t\t (\n');
fprintf(f1, '\t\t\t (1 4 5 2)\n');
fprintf(f1, '\t\t);\n');
fprintf(f1, '\t\t); \n');
fprintf(f1, '\t} \n');

fprintf(f1, '\t front\n');
fprintf(f1, '\t{\n');
fprintf(f1, '\t\t type symmetry;\n');
fprintf(f1, '\t\t faces\n');
fprintf(f1, '\t\t (\n');
fprintf(f1, '\t\t\t (0 3 4 1) \n');
fprintf(f1, '\t\t );\n');
fprintf(f1, '\ }\n');


fprintf(f1, '\t back\n');
fprintf(f1, '\t{\n');
fprintf(f1, '\t\t type symmetry ;\n');
fprintf(f1, '\t\t faces\n');
fprintf(f1, '\t\t (\n');
fprintf(f1, '\t\t\t (0 2 5 3)\n');
fprintf(f1, '\t\t);\n');
fprintf(f1, '\t }\n');


fprintf(f1, '\t axis\n');
fprintf(f1, '\t{\n');
fprintf(f1, '\t\t type empty;\n');
fprintf(f1, '\t\t faces\n');
fprintf(f1, '\t\t (\n');
fprintf(f1, '\t\t\t (0 3 3 0)\n');
fprintf(f1, '\t}\n');
fprintf(f1, '); \n\n');



blockMeshDict 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;
    object      blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

convertToMeters 1;
vertices
(
	 (0 0 0)
	 (0 1.847759e-03 7.653669e-04)
	 (0 1.847759e-03 -7.653669e-04)
	 (5.040000e-01 0 0)
	 (5.040000e-01 1.847759e-03 7.653669e-04)
	 (5.040000e-01 1.847759e-03 -7.653669e-04)
);

blocks
(
	 hex (0 3 5 2 0 3 4 1) (300 30 1) SimpleGrading (1 0.1 1)

);

edges
(
	 arc 1 2 (0 2.000000e-03 0)
	 arc 4 5 (1.204000e+00 2.000000e-03 0)
);

boundary
(
	 inlet 
	 {
		 type patch ;
		 faces 
		 (
			 (0 1 2 0)
		 ); 
	 }
	 outlet
	 {
		 type patch ;
		 faces
		 (
			 (3 5 4 3)
		 ); 
	 }
	 top
	 {
		 type wall ;
		 faces
		 (
			 (1 4 5 2)
		 );
		 
	} 
	front
	{
		 type wedge ;
		 faces
		 (
			 (0 3 4 1) 
		 );
        }
	back
	{
		 type wedge ;
		 faces
		 (
			 (0 2 5 3)
		 );
	}
	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         1;

deltaT          0.0001;

writeControl    timeStep;

writeInterval   20;

purgeWrite      0;

writeFormat     ascii;

writePrecision  6;

writeCompression off;

timeFormat      general;

timePrecision   6;

runTimeModifiable true;


// ************************************************************************* //

Results and post processing:

  • Surface with edges:

  • Velocity contour:

  • Pressure contour:

  • Velocity profiles for undeveloped and fully developed flows:

  • Shear stress contour:

  • For wedge angle 25 degree:

The matlab code to generate the blockMeshDict file is:

clear all
close all
clc
%flow parameters

Re =2100;%reynolds number
D=0.004;%diameter of the pipe
r=D/2;%radius of the pipe
l=0.06*Re*D;%hydrodynamic entrance length
l_p=l+0.7;%length of the pipe
theta=25;%wedge angle
mu=0.89e-03;%dynamic viscosity
rho=997;%density of fluid
v_avg=(Re*mu)/(rho*D);%mean velocity
v_max=2*v_avg;%maximum velocity
del_p=(8*mu*l_p*v_avg)/(r^2);%pressure drop
kin_p=(del_p)/(rho);%kinematic pressure drop
%Writing the blockMeshDict file
f1=fopen('blockMeshdict','w');
%defining the vertioces
fprintf(f1, 'vertices\n');
fprintf(f1, '(\n');
fprintf(f1, '\t (%d %d %d)\n',0,0,0);
fprintf(f1, '\t (%d %d %d)\n',0,r*cosd(theta/2),r*sind(theta/2));
fprintf(f1, '\t (%d %d %d)\n',0,r*cosd(theta/2),-r*sind(theta/2));
fprintf(f1, '\t (%d %d %d)\n',l,0,0);
fprintf(f1, '\t (%d %d %d)\n',l,r*cosd(theta/2),r*sind(theta/2));
fprintf(f1, '\t (%d %d %d)\n',l,r*cosd(theta/2),-r*sind(theta/2));
fprintf(f1, ');\n\n');
fprintf(f1, 'blocks\n');
fprintf(f1, '(\n');
fprintf(f1, '\t hex(0 3 5 2 0 3 4 1)(%d %d %d)',300,30,1);
fprintf(f1, 'SimpleGrading (1 0.1 1)\n\n');
fprintf(f1, ');\n\n');

fprintf(f1, 'edges\n');
fprintf(f1, '(\n');
fprintf(f1, '\t arc %d %d (%d%d%d)\n',1,2,0,r,0);
fprintf(f1, '\t arc %d %d(%d %d %d)\n',4,5,l_p,r,0);
fprintf(f1, ');\n\n');

%boundary conditions
fprintf(f1, 'boundary\n');
fprintf(f1, '(\n');
fprintf(f1, '\t inlet \n');
fprintf(f1, '\t {\n');
fprintf(f1, '\t {\n');
fprintf(f1, '\t\t type patch;\n');
fprintf(f1, '\t\t faces \n');
fprintf(f1, '\t\t\t (0 1 2 0)\n');
fprintf(f1, '\t\t); \n');
fprintf(f1, '\t }\n');

fprintf(f1, '\t outlet\n');
fprintf(f1, '\t {\n');
fprintf(f1, '\t\t type patch ;\n');
fprintf(f1, '\t\t faces\n');
fprintf(f1, '\t\t (\n');
fprintf(f1, '\t\t\t (3 5 4 3)\n');
fprintf(f1, '\t\t); \n');
fprintf(f1, '\t }\n');

fprintf(f1, '\t top\n');
fprintf(f1, '\t {\n');
fprintf(f1, '\t\t type wall;\n');
fprintf(f1, '\t\t faces\n');
fprintf(f1, '\t\t (\n');
fprintf(f1, '\t\t\t (1 4 5 2)\n');
fprintf(f1, '\t\t);\n');
fprintf(f1, '\t\t); \n');
fprintf(f1, '\t} \n');

fprintf(f1, '\t front\n');
fprintf(f1, '\t{\n');
fprintf(f1, '\t\t type symmetry;\n');
fprintf(f1, '\t\t faces\n');
fprintf(f1, '\t\t (\n');
fprintf(f1, '\t\t\t (0 3 4 1) \n');
fprintf(f1, '\t\t );\n');
fprintf(f1, '\ }\n');


fprintf(f1, '\t back\n');
fprintf(f1, '\t{\n');
fprintf(f1, '\t\t type symmetry ;\n');
fprintf(f1, '\t\t faces\n');
fprintf(f1, '\t\t (\n');
fprintf(f1, '\t\t\t (0 2 5 3)\n');
fprintf(f1, '\t\t);\n');
fprintf(f1, '\t }\n');


fprintf(f1, '\t axis\n');
fprintf(f1, '\t{\n');
fprintf(f1, '\t\t type empty;\n');
fprintf(f1, '\t\t faces\n');
fprintf(f1, '\t\t (\n');
fprintf(f1, '\t\t\t (0 3 3 0)\n');
fprintf(f1, '\t}\n');
fprintf(f1, '); \n\n');



blockMeshDict 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;
    object      blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

convertToMeters 1;
vertices
(
	 (0 0 0)
	 (0 1.952592e-03 4.328792e-04)
	 (0 1.952592e-03 -4.328792e-04)
	 (5.040000e-01 0 0)
	 (5.040000e-01 1.952592e-03 4.328792e-04)
	 (5.040000e-01 1.952592e-03 -4.328792e-04)
);

blocks

(
	 hex (0 3 5 2 0 3 4 1) (300 30 1) SimpleGrading (1 0.1 1)

);

edges
(
	 arc 1 2 (0 2.000000e-03 0)
	 arc 4 5 (1.204000e+00 2.000000e-03 0)
);

boundary
(
	 inlet 
	 {
		 type patch ;
		 faces 
		 (
			 (0 1 2 0)
		 ); 
	 }
	 outlet
	 {
		 type patch ;
		 faces
		 (
			 (3 5 4 3)
		 ); 
	 }
	 top
	 {
		 type wall ;
		 faces
		 (
			 (1 4 5 2)
		 );
		 
	} 
	front
	{
		 type wedge ;
		 faces
		 (
			 (0 3 4 1) 
		 );
        }
	back
	{
		 type wedge ;
		 faces
		 (
			 (0 2 5 3)
		 );
	}
	axis
	{
		 type empty ;
		 faces
		 (
			 (0 3 3 0)
	         );
        }
); 

mergePatchPairs
(
);

// ************************************************************************* //

Mesh statistics:

using the command checkMesh inside the cavity directory, the mesh stats are obtained

Points 18361
faces 36030
Internal faces 17670
cells 9000
boundary patches 6

Cells :

Hexahedra 8700
Prisms 300

Geometry

maximum aspect ratio 102.742
minimum face area 5.97e-09
maximum face area 1.45447e-06
maximum skewness 0.554166
  • Surface with edge:

The wedge with angle 25 degrees

  • Velocity contour:

  • Pressure contour:

  • Velocity plots for undeveloped regime and fully developed hydrodynamic length:

Using plot over line, the velocity plots are created along the x direction while the velocity changes along the y direction

The first two plots are for the undevelo[ed regime. The final plot is the velocity profile in the fully developed regime.

  • Shear stress plot:

For wedge angle 10 degree

clear all
close all
clc
%flow parameters

Re =2100;%reynolds number
D=0.004;%diameter of the pipe
r=D/2;%radius of the pipe
l=0.06*Re*D;%hydrodynamic entrance length
l_p=l+0.7;%length of the pipe
theta=15;%wedge angle
mu=0.89e-03;%dynamic viscosity
rho=997;%density of fluid
v_avg=(Re*mu)/(rho*D);%mean velocity
v_max=2*v_avg;%maximum velocity
del_p=(8*mu*l_p*v_avg)/(r^2);%pressure drop
kin_p=(del_p)/(rho);%kinematic pressure drop
%Writing the blockMeshDict file
f1=fopen('blockMeshdict','w');
%defining the vertioces
fprintf(f1, 'vertices\n');
fprintf(f1, '(\n');
fprintf(f1, '\t (%d %d %d)\n',0,0,0);
fprintf(f1, '\t (%d %d %d)\n',0,r*cosd(theta/2),r*sind(theta/2));
fprintf(f1, '\t (%d %d %d)\n',0,r*cosd(theta/2),-r*sind(theta/2));
fprintf(f1, '\t (%d %d %d)\n',l,0,0);
fprintf(f1, '\t (%d %d %d)\n',l,r*cosd(theta/2),r*sind(theta/2));
fprintf(f1, '\t (%d %d %d)\n',l,r*cosd(theta/2),-r*sind(theta/2));
fprintf(f1, ');\n\n');
fprintf(f1, 'blocks\n');
fprintf(f1, '(\n');
fprintf(f1, '\t hex(0 3 5 2 0 3 4 1)(%d %d %d)',300,30,1);
fprintf(f1, 'SimpleGrading (1 0.1 1)\n\n');
fprintf(f1, ');\n\n');

fprintf(f1, 'edges\n');
fprintf(f1, '(\n');
fprintf(f1, '\t arc %d %d (%d%d%d)\n',1,2,0,r,0);
fprintf(f1, '\t arc %d %d(%d %d %d)\n',4,5,l_p,r,0);
fprintf(f1, ');\n\n');

%boundary conditions
fprintf(f1, 'boundary\n');
fprintf(f1, '(\n');
fprintf(f1, '\t inlet \n');
fprintf(f1, '\t {\n');
fprintf(f1, '\t {\n');
fprintf(f1, '\t\t type patch;\n');
fprintf(f1, '\t\t faces \n');
fprintf(f1, '\t\t\t (0 1 2 0)\n');
fprintf(f1, '\t\t); \n');
fprintf(f1, '\t }\n');

fprintf(f1, '\t outlet\n');
fprintf(f1, '\t {\n');
fprintf(f1, '\t\t type patch ;\n');
fprintf(f1, '\t\t faces\n');
fprintf(f1, '\t\t (\n');
fprintf(f1, '\t\t\t (3 5 4 3)\n');
fprintf(f1, '\t\t); \n');
fprintf(f1, '\t }\n');

fprintf(f1, '\t top\n');
fprintf(f1, '\t {\n');
fprintf(f1, '\t\t type wall;\n');
fprintf(f1, '\t\t faces\n');
fprintf(f1, '\t\t (\n');
fprintf(f1, '\t\t\t (1 4 5 2)\n');
fprintf(f1, '\t\t);\n');
fprintf(f1, '\t\t); \n');
fprintf(f1, '\t} \n');

fprintf(f1, '\t front\n');
fprintf(f1, '\t{\n');
fprintf(f1, '\t\t type symmetry;\n');
fprintf(f1, '\t\t faces\n');
fprintf(f1, '\t\t (\n');
fprintf(f1, '\t\t\t (0 3 4 1) \n');
fprintf(f1, '\t\t );\n');
fprintf(f1, '\ }\n');


fprintf(f1, '\t back\n');
fprintf(f1, '\t{\n');
fprintf(f1, '\t\t type symmetry ;\n');
fprintf(f1, '\t\t faces\n');
fprintf(f1, '\t\t (\n');
fprintf(f1, '\t\t\t (0 2 5 3)\n');
fprintf(f1, '\t\t);\n');
fprintf(f1, '\t }\n');


fprintf(f1, '\t axis\n');
fprintf(f1, '\t{\n');
fprintf(f1, '\t\t type empty;\n');
fprintf(f1, '\t\t faces\n');
fprintf(f1, '\t\t (\n');
fprintf(f1, '\t\t\t (0 3 3 0)\n');
fprintf(f1, '\t}\n');
fprintf(f1, '); \n\n');



Mesh properties:

1) Mesh stats :

points 18361
faces 36030
intternal faces 17670
cells 9000
boundary patches 6

cells:

hex 8700
prism 300

geometry:

maximum aspect ratio 100.402
minimum face area 1.232e-09
maximum face area 2.931e-07
maximum skewness 0.33076
  • Results:

Surface with edge:

Velocity contour:

Pressure contour:

Velocity plots at different x locations:

Shear stress plot:

Computing walls shear stress:

  1. Go to calculator , creat mag(U) and apply
  2. Go to filters-searc -gradient of unstructured data
  3. create the name - and apply in the calculator
  4. `tau=mu*(du)/(dy)` so, in the calculator 3, tau = 0.00089*mag(U) ius calculated.
  5. The shear stress plot is calculated by , using the plot over line option

Discussions:

  • Three wedge angles, 10 degree , 25 degree and 45 degree is studied.
  • For fully developed flows, the maximum velocity obtained in all the three cases are 0.93 m/s.
  • The shear stress plot shows neglizible variation with the change in wedge angle.
  • In all these cases, symmetry boundary condition used. Tha major difference in the use of symmetry and wedge boundary conditio n is that, wedge boundary condition used for theta less than 5 degrees.
Your portfolio, next

Ship work like this
in your own course.

Talk to a career counsellor

Free 20-minute call. We'll map the right course to where you want to work.

By submitting you agree to our Terms and Privacy Policy.

Talk to a career counsellor

Free 20-minute call. We'll map the right course to where you want to work.

By submitting you agree to our Terms and Privacy Policy.