Challenge

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 wedge and symmetry
boundary conditions in OpenFOAM. Consider the flow to be laminar with Reynolds number
2100 for any diameter of the pipe. Compare the results with Hagen-Poiseuille's equations
for flow through circular pipe.

Expected Results:

  • Validation of hydrodynamic length with Numerical results.
  • Validation of velocity profile of fully developed flow with its analytical profile.
  • Validation of max velocity and pressure drop for fully developed flow.
  • Post processing of shear stress and validation of wall shear stress for fully developed flow.

Parameters:

  1. Radius of Pipe (`R`) `= 0.0075 "m"` i.e. `D=0.015 m` .
  2. Reynolds number (`Re`) `= 2100` .
  3. Dynamic Viscosity (`mu`) of water `= 0.001 Pa-s` .
  4. Density of water (`rho`) = `1000 "kg"/m^3` .

Calculations:

The calculations are based on the Hagen Poiseuille's equation for flow through circular pipe, which are as follows:

  1. Hydrodynamic entrance length (`Le`) `= 0.06*"Re"*D =1.89 m` .
  2. Length of pipe (`L`) `= Le +0.5 = 2.39 m` .
  3. Average flow velocity using `Re` , `V_"avg" = (mu*"Re")/(rho*D) = 0.14 m/s` .
  4. Friction factor for Laminar flow (`f`) `= 64/"Re" = 0.0305 ` .
  5. Pressure gradient (`- (delP)/(delx)`) `= (f*rho*(V_"avg")^2)/(2*D) = 19.91 N/m^3`
  6. Maximum velocity (`V_"max"`) `= 2 * V_"avg" = 0.28 m/s` .
  7. Max Shear stress occurs at wall (`tau_"wall"`) `= 8*mu*V_"avg"/D = 0.0747 N/m^2` .
  8. Dynamic Pressure drop at outlet (`DeltaP`) `= (- (delP)/(delx))*"Le"/rho =0.0476 m^2/s^2` .
  9. Radial velocity profile is given by `V_r = (1/(4mu))(- (delP)/(delx))*(R^2 -r^2)`
    where `r `- radial distance (`m`) .

Workflow:

(1) MATLAB code for automatic generation of fully functional the blockMeshDict file with necessary input parameters for wedge and symmetry boundary conditions.

(a) Code to generate blockMeshDict file for wedge BC for 25 degree wedge angle:

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

convertToMeters 1;

vertices
(
    (0.000000 0.000000 0.000000)       //0 
    (0.000000 0.007322 0.001623)       //1 
    (0.000000 0.007322 -0.001623)       //2 
    (2.390000 0.000000 0.000000)       //3 
    (2.390000 0.007322 0.001623)       //4 
    (2.390000 0.007322 -0.001623)       //5 

  );
blocks 
(
     hex (0 2 1 0 3 5 4 3) (50 1 100) simpleGrading (0.100000 1.000000 1.000000) 

  );
edges 
 (
     arc 1 2 (0.000000 0.007500 0.000000) 
     arc 4 5 (2.390000 0.007500 0.000000) 
  );
boundary 
(
    inlet 
     {
        type patch; 
        faces 
        (
             (0 1 2 0) 
        );
    }
    outlet 
    {
        type patch; 
        faces 
        ( 
            (3 4 5 3) 
        ); 
    }
    upperWall 
    {
        type wall; 
        faces 
        (
            (2 1 5 4) 
       ); 
    }
    front 
    {
        type wedge; 
        faces 
        (
            (0 3 4 1) 
        ); 
    }
    back 
    {
        type wedge; 
        faces 
        (
            (0 2 5 3) 
         ); 
    }
);
// ************************************************************************* //



(b) Code to generate blockMeshDict file for symmetry BC for 25 degree angle:

close all
clear all
clc

% input parameters
theta = 12.5;                                 % halfangle of wedge in degrees
radius = 0.0075;                           % radius of pipe in m 
Re = 2100 ;                                % Reynolds number
hydlength = 0.06*Re*2*radius ;             % Hydrodynamic Entrance length
length = hydlength + 0.5;                 % Pipe length
p1 = radius*cosd(theta);                   % calculation for point 1
p2 = radius*sind(theta);                   % calculation for point 2 

% Vertex inputs
xmin= 0;                          
xmax= length;
ymin= 0;
ymax= p1; 
zmax= p2; 
zmin=-1.0*p2;

% Mesh gridpoints and simple grading inputs

nx = 50;
ny = 3;         
nz = 100;
gx = 0.1;
gy = 1;
gz = 1;

%% Writing blockMeshDict file

f1=fopen('blockMeshDict.txt','w');
fprintf(f1,'\n/*--------------------------------*- C++ -*----------------------------------*\\');
fprintf(f1,'\n  =========                 |');
fprintf(f1,'\n  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox ');
fprintf(f1,'\n   \\    /   O peration     | Website:  https://openfoam.org ');
fprintf(f1,'\n    \\  /    A nd           | Version:  9 ');
fprintf(f1,'\n     \\/     M anipulation  |\n');
fprintf(f1,'\n\\*---------------------------------------------------------------------------*/');
fprintf(f1,'\nFoamFile \n {\n    format      ascii;\n    class       dictionary;\n    object      blockMeshDict;\n  }');
fprintf(f1,'\n// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //');
fprintf(f1,'\n\nconvertToMeters 1;\n');
fprintf(f1,'\nvertices\n(\n    (%f %f %f)       //0 \n',xmin,ymin,ymin);
fprintf(f1,'    (%f %f %f)       //1 \n',xmin,ymax,zmax);
fprintf(f1,'    (%f %f %f)       //2 \n',xmin,ymax,zmin);
fprintf(f1,'    (%f %f %f)       //3 \n',xmax,ymin,ymin);
fprintf(f1,'    (%f %f %f)       //4 \n',xmax,ymax,zmax);
fprintf(f1,'    (%f %f %f)       //5 \n\n  );\n',xmax,ymax,zmin);
fprintf(f1,'blocks \n(\n     hex (0 2 1 0 3 5 4 3) (%d %d %d) simpleGrading (%f %f %f) \n\n  );\n',nx,ny,nz,gx,gy,gz);
fprintf(f1,'edges \n (\n     arc 1 2 (%f %f %f) \n     arc 4 5 (%f %f %f) \n  );\n',xmin,radius,0,xmax,radius,0);
fprintf(f1,'boundary \n(\n    inlet \n     {\n        type patch; \n        faces \n        (\n             (0 1 2 0) \n        );\n    }\n');
fprintf(f1,'    outlet \n    {\n        type patch; \n        faces \n        ( \n            (3 4 5 3) \n        ); \n    }\n');
fprintf(f1,'    upperWall \n    {\n        type wall; \n        faces \n        (\n            (2 1 5 4) \n       ); \n    }\n');
fprintf(f1,'    front \n    {\n        type symmetry; \n        faces \n        (\n            (0 3 4 1) \n        ); \n    }\n');
fprintf(f1,'    back \n    {\n        type symmetry; \n        faces \n        (\n            (0 2 5 3) \n         ); \n    }\n);\n');
fprintf(f1,'// ************************************************************************* //');   
fclose(f1);


Output(a): Fully functional blockMeshDict file for Wedge angle 25 with wedge BC.

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

convertToMeters 1;

vertices
(
    (0.000000 0.000000 0.000000)       //0 
    (0.000000 0.007322 0.001623)       //1 
    (0.000000 0.007322 -0.001623)       //2 
    (2.390000 0.000000 0.000000)       //3 
    (2.390000 0.007322 0.001623)       //4 
    (2.390000 0.007322 -0.001623)       //5 

  );
blocks 
(
     hex (0 2 1 0 3 5 4 3) (50 1 100) simpleGrading (0.100000 1.000000 1.000000) 

  );
edges 
 (
     arc 1 2 (0.000000 0.007500 0.000000) 
     arc 4 5 (2.390000 0.007500 0.000000) 
  );
boundary 
(
    inlet 
     {
        type patch; 
        faces 
        (
             (0 1 2 0) 
        );
    }
    outlet 
    {
        type patch; 
        faces 
        ( 
            (3 4 5 3) 
        ); 
    }
    upperWall 
    {
        type wall; 
        faces 
        (
            (2 1 5 4) 
       ); 
    }
    front 
    {
        type wedge; 
        faces 
        (
            (0 3 4 1) 
        ); 
    }
    back 
    {
        type wedge; 
        faces 
        (
            (0 2 5 3) 
         ); 
    }
);
// ************************************************************************* //



Output(b): Fully functional blockMeshDict file for Wedge angle 25 with symmetry BC.

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

convertToMeters 1;

vertices
(
    (0.000000 0.000000 0.000000)       //0 
    (0.000000 0.007322 0.001623)       //1 
    (0.000000 0.007322 -0.001623)       //2 
    (2.390000 0.000000 0.000000)       //3 
    (2.390000 0.007322 0.001623)       //4 
    (2.390000 0.007322 -0.001623)       //5 

  );
blocks 
(
     hex (0 2 1 0 3 5 4 3) (50 3 100) simpleGrading (0.100000 1.000000 1.000000) 

  );
edges 
 (
     arc 1 2 (0.000000 0.007500 0.000000) 
     arc 4 5 (2.390000 0.007500 0.000000) 
  );
boundary 
(
    inlet 
     {
        type patch; 
        faces 
        (
             (0 1 2 0) 
        );
    }
    outlet 
    {
        type patch; 
        faces 
        ( 
            (3 4 5 3) 
        ); 
    }
    upperWall 
    {
        type wall; 
        faces 
        (
            (2 1 5 4) 
       ); 
    }
    front 
    {
        type symmetry; 
        faces 
        (
            (0 3 4 1) 
        ); 
    }
    back 
    {
        type symmetry; 
        faces 
        (
            (0 2 5 3) 
         ); 
    }
);
// ************************************************************************* //

Similarly, for wedge angle 10 and 45 the blockMesh files are generated with wedge and symmetry BC.

(2) OpenFOAM case set up:

BlockMeshDict file is as generated above.

(a)Initial Velocity boundary conditions (U file in 0 folder) with wedge BC:

/*--------------------------------*- 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
{
    inlet
    {
        type            fixedValue;
        value           uniform (0.14 0 0);
    }

    outlet
    {
        type            zeroGradient;
    }

    upperWall
    {
        type            noSlip;
    }

    front
    {
        type wedge;
       
    }
     
    back
    {
        type wedge;
     }

}

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

(b)Initial Velocity boundary conditions (U file in 0 folder) with symmetry BC :

/*--------------------------------*- 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
{
    inlet
    {
        type            fixedValue;
        value           uniform (0.14 0 0);
    }

    outlet
    {
        type            zeroGradient;
    }

    upperWall
    {
        type            noSlip;
    }

    front
    {
        type symmetry;
       
    }
     
    back
    {
        type symmetry;
     }

}

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

(c)Initial pressure boundary conditions (p file in 0 folder) with wedge BC:

/*--------------------------------*- 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;
    }
    
    upperWall
    {
        type zeroGradient;
       
    }
    
    front
    {
        type wedge;
       
    }
     
    back
    {
        type wedge;
     }

    
}

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

(d)Initial pressure boundary conditions (p file in 0 folder) with symmetry BC:

/*--------------------------------*- 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;
    }
    
    upperWall
    {
        type zeroGradient;
       
    }
    
    front
    {
        type symmetry;
       
    }
     
    back
    {
        type symmetry;
     }

    
}

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

(e)Thermophysical properties (Transportproperties file from constant folder):

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

transportModel  Newtonian;

nu              [0 2 -1 0 0 0 0] 1e-06;

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

(f)Solution controls (controlDict file from system folder):

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

deltaT          0.01;

writeControl    timeStep;

writeInterval   100;

purgeWrite      0;

writeFormat     ascii;

writePrecision  6;

writeCompression off;

timeFormat      general;

timePrecision   6;

runTimeModifiable true;


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

(3) Mesh is generated using blockMesh command and it has been checked by checkMesh command. For wedge BC the error regarding nonalignment of edges has occurred for greater than 1 number of cell in Z direction. Also for Wedge angle 10 with symmetry BC , 3 and 2 number of cells in Z direction caused cell aspect ratio high error so its number of cells in Z direction reduce to 1 (default).

Mesh quality is found to be OK for all meshes. Meshes has also been inspected in paraview using paraFoam command.

(a)Mesh For 10 degree wedge angle with wedge and symmetry BC respectively:




(b)Mesh for 25 degree wedge angle with wedge and symmetry BC respectively:



(c)Mesh for 45 degree wedge angle with wedge and symmetry BC respectively:

(4) Simulation is done by using icoFoam command and results are obtained. Same are post processed in the paraview using paraFoam command. As we simulated the pipe flow for the small sector of pipe with different wedge angles, to obtain the results for full cross section of pipe the reflect filter from paraview is used. Then using slice command the diametric plane from the geometry is obtained to do further post processing. The velocity and pressure contours for various time steps are observed found acceptable with the analytical results from Hagen Poiseuille's equations.

Velocity profile of fully developed flow over the entire pipe length:

(i) For Wedge angle 10 degree with wedge BC:


(ii) For Wedge angle 10 degree with Symmetry BC:

(iii) For Wedge angle 25 degree with wedge BC:


(iv) For Wedge angle 25 degree with Symmetry BC:

(v) For Wedge angle 45 degree with wedge BC:


(vi) For Wedge angle 45 degree with Symmetry BC:





Velocity profile over cross section at a distance of 1.89m are inspected for various wedge angles :

(i) For Wedge angle 10 degree with wedge BC:

(ii) For Wedge angle 10 degree with Symmetry BC:

(iii) For Wedge angle 25 degree with wedge BC:

(iv) For Wedge angle 25 degree with Symmetry BC:

(v) For Wedge angle 45 degree with wedge BC:

(vi) For Wedge angle 45 degree with Symmetry BC:

The Pressure drop over the entire length of pipe are also obtained using line over plot filter. The graphs for various wedge angles with wedge and symmetry BC are as follows:

(i) For Wedge angle 10 degree with wedge BC:

(ii) For Wedge angle 10 degree with Symmetry BC:

(iii) For Wedge angle 25 degree with wedge BC:


(iv) For Wedge angle 25 degree with Symmetry BC:



(v) For Wedge angle 45 degree with wedge BC:

(vi) For Wedge angle 45 degree with Symmetry BC:



(5) As shear stress is not a primitive variable that has to be obtained through calculations of primitive variables as follows:

(i) First of all using calculator filter the magnitude of velocity (Mag_U) is obtained for entire domain.

(ii) Then gradient of magnitude of velocity (Grad_Mag_U) is obtained using filter gradient of unstructured dataset.

(iii) As as shear stress (Tau) is gradient of velocity is coefficient of viscosity times gradient of velocity. It is obtained using calculator filter.

(iv) As shear stress is tensor its magnitude (Mag_Tau) is obtained using calculator filter.

(v) The shear stress profiles are obtained at the distance of 1.89m from the inlet. The respective graphs for different wedge angles with wedge and symmetry BC are as follows:

(a)For Wedge angle 10 degree with wedge BC:


(b)For Wedge angle 10 degree with Symmetry BC:



(c)For Wedge angle 25 degree with wedge BC:


(d)For Wedge angle 25 degree with symmetry BC:

(e)For Wedge angle 45 degree with wedge BC:

(f)For Wedge angle 45 degree with Symmetry BC:

(6)Afterwards, all of the required data from the results at the distance of 1.89 m from inlet are obtained through paraview using spreadsheet view (CSV file and text files are attached herewith). The same data is used in MATLAB by file parsing technique and the result data is compared with analytical data by using following MATLAB code.

close all
clear all
clc
 
%% Inputs for Analytical solution of Hagen Poiseuille
theta = [10,25,45];                        % halfangle of wedge in degrees
radius = 0.0075;                           % radius of pipe in m 
Re = 2100 ;                                % Reynolds number
hydlength = 0.06*Re*2*radius ;             % Hydrodynamic Entrance length
length = hydlength + 0.5;                  % Pipe length
mu=0.001;                                  % Dynamic viscosity of water
rho=1000;                                  % Desity of water 
f=64/Re;                                   % Friction factor for laminar flow 
Vavg=Re*mu/(rho*2*radius);                 % Avg velocity computation using Re number 

y =linspace(0,radius,1000);                        % Radial distance
V =(0.0625*f*rho*Vavg^2/(radius*mu))*(radius^2-y.^2); % Velocity profile
Vmax =(0.0625*f*rho*Vavg^2/(radius*mu))*(radius^2);   % Max velocity
tau =(0.125*f*rho*Vavg^2/radius)*y;                % Shear stress
deltap_l =(0.25*f*rho*Vavg^2/radius);              % Pressure gradient per unit length
tau_max =max(tau);                                 % Calculated max shear stress
DeltaP_o =deltap_l*(2.39)/rho;                     % Dynamic pressure diffence at outlet


%% Data Reading from Paraview data file for 10/25/45 degree wedge and data plotting
f1=fopen('data10.txt','r');     
fgetl(f1)
A=fscanf(f1,'%f',[3 inf]);      % Reading data from 10 degree datafile  
A=A';
fclose(f1);

VmaxA=max(A(:,2));              % Maximum velocity from result data
tau_maxA=max(A(:,3));           % Maximum shear stress from result data

f11=fopen('data10w.txt','r');     
fgetl(f11)
A1=fscanf(f11,'%f',[3 inf]);      % Reading data from 10 degree datafile  
A1=A1';
fclose(f11);

VmaxA1=max(A1(:,2));              % Maximum velocity from result data
tau_maxA1=max(A1(:,3));           % Maximum shear stress from result data

f2=fopen('data25.txt','r');
fgetl(f2)
B=fscanf(f2,'%f',[3 inf]);      % Reading data from 25 degree datafile    
B=B';
fclose(f2);

VmaxB=max(B(:,2));              % Maximum velocity from result data
tau_maxB=max(B(:,3));           % Maximum shear stress from result data

f21=fopen('data25w.txt','r');
fgetl(f21)
B1=fscanf(f21,'%f',[3 inf]);      % Reading data from 25 degree datafile    
B1=B1';
fclose(f21);

VmaxB1=max(B1(:,2));              % Maximum velocity from result data
tau_maxB1=max(B1(:,3));           % Maximum shear stress from result data

f3=fopen('data45.txt','r');     % Reading data from 45 degree datafile
fgetl(f3)
C=fscanf(f3,'%f',[3 inf]);
C=C';
fclose(f3);

VmaxC=max(C(:,2));               % Maximum velocity from result data
tau_maxC=max(C(:,3));            % Maximum shear stress from result data

f31=fopen('data45w.txt','r');     % Reading data from 45 degree datafile
fgetl(f31)
C1=fscanf(f31,'%f',[3 inf]);
C1=C1';
fclose(f31);

VmaxC1=max(C1(:,2));               % Maximum velocity from result data
tau_maxC1=max(C1(:,3));            % Maximum shear stress from result data

figure(1)
plot(A(:,2),A(:,1),LineWidth=1)
hold on
plot(A1(:,2),A1(:,1),LineWidth=1)
hold on
plot(B(:,2),B(:,1),LineWidth=1)
hold on
plot(B1(:,2),B1(:,1),LineWidth=1)
hold on
plot(C(:,2),C(:,1),LineWidth=1)
hold on
plot(C1(:,2),C1(:,1),LineWidth=1)
hold on
plot(V,y,'r',LineWidth=1)
hold on
plot(V,-y,'r',LineWidth=1)
grid on
grid minor
title('Velocity Profile Comparison')
xlabel('Velocity (m/s)')
ylabel('Radial distance (m)')
legend ('Comp-Vel-10SBC','Comp-Vel-10WBC','Comp-Vel-25SBC','Comp-Vel-25WBC','Comp-Vel-45SBC','Comp-Vel-45WBC','Analytical Vel',Location='west')
saveas(gcf,'Velocity Profile study-1','png')


figure(2)
plot(A(:,3),A(:,1),LineWidth=1)
hold on
plot(A1(:,3),A1(:,1),LineWidth=1)
hold on
plot(B(:,3),B(:,1),LineWidth=1)
hold on
plot(B1(:,3),B1(:,1),LineWidth=1)
hold on
plot(C(:,3),C(:,1),LineWidth=1)
hold on
plot(C1(:,3),C1(:,1),LineWidth=1)
hold on
plot(tau,y,'r',LineWidth=1)
hold on
plot(tau,-y,'r',LineWidth=1)
grid on
grid minor
title('Shear stress Profile Comparison')
xlabel('Shear stress (N/m^2)')
ylabel('Radial distance (m)')
legend ('Comp-Shear-10SBC','Comp-Shear-10WBC','Comp-Shear-25SBC','Comp-Shear-25WBC','Comp-Shear-45SBC','Comp-Shear-45WBC','Analytical Shear Stress',Location='east') 
saveas(gcf,'Shear stress Profiles study-1','png')

Workspace screenshot:

Velocity Profile Comparison:

Shear stress profile comparison:

Table of Observations:


Conclusions:

  1. By observing the graphs of velocity profile over the length and velocity profiles at 1.89m, it can be concluded that the flow through pipe is fully developed at the 1.89m from the inlet of the pipe. From the observation table the value of Vmax for computational results is more than 99% of the analytical value of Vmax. Thus the hydrodynamic entrance length of 1.89m is validated.
  2. Numerical Velocity profiles are compared with Analytical velocity profile and the results are stated in velocity profile comparison. The velocity profile wedge boundary condition for larger wedge angles show some distortion in the profile towards the wall while having closer values of Vmax with Vmax analytical. Close match of the velocity profile for symmetry BC indicates that the numerical results are very close to the analytical solution. Thus the velocity profiles are validated.
  3. Similarly Shear Stress are compared with analytical shear stress profile and the results are plotted in graph of Shear stress profile comparison. The wall shear stress values for wedge BC have greater error than the corresponding values with symmetry BC. The shear stress profile for larger wedge angles with wedge BC have more deviation from the analytical shear stress profile, it is evident from Shear stress profile comparison. Thus the shear stress profiles are validated.
  4. Numerical Values of `DeltaP` are calculated from the data and compared with the analytical value for `DeltaP` . Same is stated in the observation table. Values are very close to analytical values thus the `DeltaP` are validated. It can be observed that the error in the `DeltaP` values with wedge BC is larger than the error in the `DeltaP` values with symmetry BC conditions. Hence the symmetry BC should be preferred for larger wedge angles.
  5. Thus the Hagen Poiseuille's equations for circular pipe are verified through numerical simulation of laminar flow through pipe with wedge and symmetry BCs for various wedge angles.
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.