BlogFor BusinessTalk to a counsellor
CoursesIIT ProgramsPlacementsLog in
BlogFor Business
Talk to a counsellor
Project

Symmetry vs Wedge vs HP equation

AAnam Ali10 May 2020
CFD

Objective

The objective of this project is to

  • Perform the simulation with Symmetry Boundary Condition of the wedge geometry.
  • 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

Note:

All the steps are same as we did in week 11 challenge only change occur in the fornt back boundary condition which we take symmetry in this case. And we perform simulation with three different wedge angles.

Procedure

1. In matlab blockMeshDict file we change only the front and back boundary conditions wedge to symmetry.

2. And also change the wedge angle 10,25 and 45 respectively.

3. After getting the blockMeshDict file for symmetry boundary conditions and changed wedge angle we copy the file in the system folder which is located in the cavity folder.

4. After copying the blockMeshDict file we have to edit the velocity and pressure file.

5. Now we edit the controlDict file and put the value of delta t and set the simulation time.

6. Then go to system folder open the transport properties file and put the value of kinematic viscosity.

7.After setting all the parameters run the blockMesh command in the terminal.

8.Then run the \'icoFoam\' solver.

9.After the solver stops, run the paraFoam for postprocessing the results.

Matlab code 

clear all
close all
clc
%%%inputs%%%
Re = 2100;      %Reynold numbers 
D = 0.004;      %Diameter of pipe
L = 0.06*Re*D;  %Entrance length of pipe
L_p = L+0.7;    %total length of pipe
r = D/2;        %radius of pipe
theta =input(\'value of angle\');         %wedge angle, it should be less than 5.
mu = 0.89e-03;  %Dynamic viscosity
rho = 997;      %density
%%%Analytical values%%%
% Hagen poiseuille flow equations
v_avg = (Re*mu)/(rho*D);
v_max = 2*(v_avg);
del_p = (8*mu*L*v_avg)/r^2;
kinematic_p = del_p/rho;
% Calculating the shear stress
R = linspace(0,D/2,1000)
for i= 1:length(R)
    tau(i) = 2*mu*v_max*(abs(R(i)))/r^2;
end
plot(R,tau)
xlabel(\'Radius\')
ylabel(\'shear stress\')
grid on
legend(\'shear stress\')
header1=\'/*--------------------------------*- C++ -*----------------------------------*\\\';
header2=\'=========                 |                                                    \'
header3=\'\\\\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox              \'
header4=\' \\\\    /   O peration     | Website:  https://openfoam.org                     \'
header5=\'  \\\\  /    A nd           | Version:  7                                        \'
header6=\'   \\\\/     M anipulation  |                                                    \'
header7=\'\\*---------------------------------------------------------------------------*/\'
f1 = fopen(\'blockMeshDict.txt\',\'w\')
fprintf(f1, \'%s\\n\', header1)
fprintf(f1, \'%s\\n\', header2)
fprintf(f1, \'%s\\n\', header3)
fprintf(f1, \'%s\\n\', header4)
fprintf(f1, \'%s\\n\', header5)
fprintf(f1, \'%s\\n\', header6)
fprintf(f1, \'%s\\n\', header7)
fprintf(f1, \'FoamFile\\n{\\n\')
fprintf(f1, \'%12s \\t 2.0;\\n\', \'version\')
fprintf(f1, \'%12s \\t ascii;\\n\', \'format\')
fprintf(f1, \'%12s \\t dictionary;\\n\', \'class\')
fprintf(f1, \'%12s \\t blockMeshDict;\\n\', \'object\')
fprintf(f1,\'}\\n\')
header8=\'// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //\'
fprintf(f1, \'%s\\n\\n\', header8)

  fprintf(f1, \'convertToMeters 1;\\n\\n\')
  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_p,0,0)
  fprintf(f1, \'\\t (%d %d %d)\\n\', L_p,r*cosd(theta/2),r*sind(theta/2))
  fprintf(f1, \'\\t (%d %d %d)\\n\', L_p,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\')

fprintf(f1,\'boundary\\n\')
fprintf(f1,\'(\\n\')
fprintf(f1, \'\\t inlet \\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 (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 topWall\\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 }\\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,\'\\t }\\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\\t );\\n\')
fprintf(f1,\'\\t }\\n\')
fprintf(f1, \');\\n\\n\')

fprintf(f1, \'mergePatchPairs\\n\')
fprintf(f1, \'(\\n\')
fprintf(f1, \');\\n\')

header9=\'// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //\'
fprintf(f1, \'%s\\n\', header9)

fclose(f1)

blockMeshDict file for 10 degree wedge

/*--------------------------------*- 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 0 0)
	 (0 1.992389e-03 1.743115e-04)
	 (0 1.992389e-03 -1.743115e-04)
	 (1.204000e+00 0 0)
	 (1.204000e+00 1.992389e-03 1.743115e-04)
	 (1.204000e+00 1.992389e-03 -1.743115e-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)
		 );
	 }
	 topWall
	{
		 type wall;
		 faces
		 (
			 (1 4 5 2)
		 );
	 }
	 front
	{
		 type symmetry;
		 faces
		 (
			 (0 3 4 1)
		 );
	 }
	 back
	{
		 type symmetry;
		 faces
		 (
			 (0 2 5 3)
		 );
	 }
	 axis
	{
		 type empty;
		 faces
		 (
			 (0 3 3 0)
		 );
	 }
);

mergePatchPairs
(
);
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

blockMeshDict file for 25 degree wedge

/*--------------------------------*- 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 0 0)
	 (0 1.952592e-03 4.328792e-04)
	 (0 1.952592e-03 -4.328792e-04)
	 (1.204000e+00 0 0)
	 (1.204000e+00 1.952592e-03 4.328792e-04)
	 (1.204000e+00 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)
		 );
	 }
	 topWall
	{
		 type wall;
		 faces
		 (
			 (1 4 5 2)
		 );
	 }
	 front
	{
		 type symmetry;
		 faces
		 (
			 (0 3 4 1)
		 );
	 }
	 back
	{
		 type symmetry;
		 faces
		 (
			 (0 2 5 3)
		 );
	 }
	 axis
	{
		 type empty;
		 faces
		 (
			 (0 3 3 0)
		 );
	 }
);

mergePatchPairs
(
);
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

BlockMeshDict file for 45 degree wedge

/*--------------------------------*- 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 0 0)
	 (0 1.847759e-03 7.653669e-04)
	 (0 1.847759e-03 -7.653669e-04)
	 (1.204000e+00 0 0)
	 (1.204000e+00 1.847759e-03 7.653669e-04)
	 (1.204000e+00 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)
		 );
	 }
	 topWall
	{
		 type wall;
		 faces
		 (
			 (1 4 5 2)
		 );
	 }
	 front
	{
		 type symmetry;
		 faces
		 (
			 (0 3 4 1)
		 );
	 }
	 back
	{
		 type symmetry;
		 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:  7
     \\\\/     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.001;

writeControl    timeStep;

writeInterval   20;

purgeWrite      0;

writeFormat     ascii;

writePrecision  6;

writeCompression off;

timeFormat      general;

timePrecision   6;

runTimeModifiable true;


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

Pressure boundary condition

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

    topWall
    {
        type            zeroGradient;
    }
    front
    {
        type            symmetry;
    }
    back
    {
        type            symmetry;
    }
    Axis
    {
        type            empty;
    }
}
// ************************************************************************* //

velocity boundary condition

/*--------------------------------*- 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.4686 0 0);

boundaryField
{
    inlet
    {
        type            fixedValue;
        value           uniform (0.4686 0 0);
    }

    outlet
    {
        type            zeroGradient;
    }

    topWall
    {
        type            noSlip;
    }
    front
    {
        type            symmetry;
    }
    back
    {
        type            symmetry;
    }
    Axis
    {
        type            empty;
    }
}

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

Transport properties

/*--------------------------------*- 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] 0.000000893;


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

Generated Mesh for 10 degree wedge

\"\"

\"\"

\"\"

Postprocessing

Velocity profile at various lengths

1.Velocity profile at entrance length

\"\"

2.Velocity profile at 0.01m length

\"\"

3.Velocity profile at 0.1m length

\"\"

3.Velocity profile at 0.2m length

\"\"

4.velocity profile at 0.35 m length

\"\"

5.Velocity profile at 0.504 m length

\"\"

6.velocity profile at 0.7m length

\"\"

7.velocity profile at 1m length

\"\"

8.velocity profile at 1m length

\"\"

9.velocity profile at 1.2m length

\"\"

\"\"

Generated Mesh and geometry at 25 degree wedge

\"\"

\"\"

Velocity profiles for 25 degree wedge

1. Velocity profiles at 0.01m length

\"\"

2.Velocity profiles at 0.3m length

\"\"

3.Velocity profile at 0.504m length

\"\"

4. Velocity profile at 1m length

\"\"

5.Velocity profile at 1.2m length

\"\"

\"\"

Generated Mesh and geometry at 45 degree Mesh

\"\"

\"\"

\"\"

Velocity profiles at 45 degree angle

1.velocity profile at 0.01m length

\"\"

2. Velocity profile at 0.1m length

\"\"

3.Velocity profile at 0.3m length

\"\"

4.Velocity profile at 0.504m length

\"\"

5.Velocity profile at 1m length

\"\"

6.velocity pofile at 1.2m length

\"\"

\"\"

Results

\"\"

Parameter Wedge wedge symmetry_10 symmetry_25 symmetry_45

N.Velocity 0.925 0.9280 0.9291 0.9295

A.Velocity 0.9373 0.9373 0.9373 0.9373

Error 0.0123 0.0093 0.0082 0.0078

Conclusions

1.Parabolic profile is maintained in all the wedge angles.

2.The velocity is almost constant in all wedge angle.

3.Error is decreasing with the wedge angle increases

4.Symmetry boundary conditions are gave the more accurate results as compared to wedge boundary condition.

5. Execution time is greater for symmetry boundary conditions than wedge boundary conditions.

6.Behaviour of velocity profiles is almost same in symmetry and wedge type boundary conditions.

Shear Stress

For all Newtonian fluids in a laminar fow , the shear stress is proportional to the strain rate in the fluid, where the viscosity is the constant of proportionality.

Formula of shear stress is

`tau=mufrac{partialu}{partialy}`

Steps involve calculating the shear stress in a ParaView

1. In step1, we are calculating the magnitude of the velocity field using calculator tool.

2. In step2, we are calculating the gradient of the velocity field using the \'Gradients of unstructure data set tool\'.

3. In step 3, we are calculating the shear stress using the calculator tool.

4. In step 4, we are plotting the shear stress using the plot over line command. And we are plotting the shear stress at radial distance of 1.2m which is the fully developed region.

Numerical results

\"\"

\"\"

\"\"

Shear stress Analytical code of matlab

clear all
close all
clc
D=0.004;
r=D/2;
R = linspace(-0.01,0.01,1000)
Re=2100;
mu=8.9e-04;
rho=997;
V= (Re*mu)/(rho*D);
v_max =2*V;
for i= 1:length(R)
    tau(i) = 2*mu*v_max*(abs(R(i)))/r^2;
end
figure(1)
plot(R,tau)
xlabel(\'Radius\')
ylabel(\'shear stress\')
grid on
legend(\'shear stress\')

Graph of analytical code

\"\"

Results comparison

Parameter Values Analytical Value Error

Wedge angle 10 Shear stress 0.81 0.87 0.06%

Wedge angle 25 Shear stress 0.83 0.87 0.04%

Wedge angle 45 Shear stress 0.86 0.87 0.01%

Conclusions

1. In the above numerical and analytical graphs we observe that behaviour of the shear stress graphs is linear.

2. At the center of the wall the value of the shear stress is zero and slowly slowly when it is near about the walls its values becomes maximum.

3. Approximately the analytical and numerical values of shear stress are matched.

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.