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

Simulation of flow through the pipe using OpenFOAM Wedge type

GGanesh Kumar4 Aug 2019
CFDMATLAB

Simulation of flow through the pipe using OpenFoam

Objective:

The main objective of this project is to simulate flow through the pipe using OpenFOAM software. For that,

  • Write a program in Matlab that can generate the computational mesh automatically for any wedge angle and grading scheme.
  • Calculate the length of the pipe using the entry length formula for laminar flow through a pipe
  • Show that entry length is sufficient to produce a fully developed flow. You will have to compare velocity profiles at different positions along the length of the pipe.
  • Show that the velocity profile is fully developed
  • Post-process velocity and shear stress.

Assumption:

  • Steady flow, Incompressible & Laminar
  • Flow is two dimensional (Flow WRT theta is ignored)
  • Constant fluid properties ( Viscosity, Density, etc.,)

Physics behind this problem

Introduction to Laminar flow through a pipe

In this project, we aim to simulate a section of pipe (Wedge shaped geometry) to understand various aspects of flow through the pipe. The fluid chosen is water & the flow is laminar & incompressible. Also, we choose to neglect any changes in cross-section. The specified Reynold\'s number Re is 2100.

Internal flow occurs when the fluid is filled entirely in a conduit & is driven primarily by a pressure difference.

Fluid velocity in a pipe changes from zero at the surface because of the no-slip boundary condition to a maximum at the pipe center. The average velocity remains constant in an incompressible flow when the cross-section of the pipe is constant.

Reynold\'s Number (Re):

The transition from laminar flow to turbulent flow occurs from laminar (Re`<=`2100) to Transition (2100`<=`Re`<=`4000) to turbulent (Re `>=`4000) for smooth pipe with no viscous dissipation. It depends on the geometry, face roughness, flow velocity, surface temperature, etc. Reynold\'s number is an important non- Dimensionless number which is the ratio of the inertial force to viscous force and is expressed for internal flow in a circular section pipe as,

Hydrodynamic Length (Le) :

The region from the pipe inlet to the point at which the velocity boundary layer merges at the centreline is called the hydrodynamic entrance region in which the velocity profile is fully developed is called the fully developed region.

In laminar flow, the entry length is given approximately as,

Velocity Profile for laminar flow:

The velocity profile for a fully developed laminar flow in a pipe is parabolic with a maximum at the center line and minimum at R=0.

It is obtained for a fluid element at a radial distance from the center. The average velocity can be written as:

Hence the final expression for the velocity profile can be written as:

Hagen poiseuille Flow: It states that \"for a specified flow rate, the pressure drop & thus the required pumping power is inversely proportional to the length of the pipe & the viscosity of the fluid, but it is inversely proportional to the fourth power of the radius (or diameter) of the pipe.

Analytical calculation:

Renold\'s number(Re) - 2100, So the flow is laminar

Diametre : 0.01m

For water, at 25 Degree

Density, rho = 997 kg/m^3

Dynamic viscosity, mu=8.9e-3 Pa.s

Kinematic viscosity,nu=rho/mu

nu=8.927e-7 m^2/s

Average velocity,Vavg=0.1869 m/s

Maximum velocity, Vmax =0.3738 m/s

Pressure drop, dp=67.0687 pa

Kinematic pressure =dp/rho

= 0.0671 MPa

MATLAB PROGRAM:

(blockMeshDict) Computational mesh generated using the MATLAB program.

clear all
close all
clc

th = input(\'Enter the value of th=\');
%Inputs for entry length
Re=2100;
D=0.01;
L=0.06*Re*D; %length of the pipe by entry length formulae
r=D/2;
rho=1000; %Density
mu=0.00089;%Dynamic viscosity
nu=mu/rho; %%Kinematic viscosity

%% Analytical solutions
Vavg=Re*nu/D; %Average velocity
Vmax=2*Vavg; %Maximum velocity
dp=(32*mu*L*Vavg)/(D^2); %Difference pressure
dp_k=dp/rho; %%Kinematic pressure

s1=blanks(4);
s2=blanks(8);
s3=blanks(12);

%Vertice
V0=[ 0 0 0 ]; 
V1=[0 r*cosd(th/2) -r*sind(th/2)];
V2=[0 r*cosd(th/2) r*sind(th/2)];
V3=[L 0 0];
V4=[L r*cosd(th/2) -r*sind(th/2)];
V5=[L r*cosd(th/2) r*sind(th/2)];

f1=fopen(\'blockMeshDict.txt\',\'w\');

%%edge1=arc 1 2(0 r 0);
%%edge2=arc 4 5(L r 0);
%Specifying blockMeshDict file
h1=\"/*--------------------------------*- C++ -*----------------------------------*\\\";
h2=\"| ========= |                                                                 |\";
h3=\"| \\\\      / F ield | OpenFOAM: The Open Source CFD Toolbox                    |\";
h4=\"|  \\\\    /  O peration | Version: 4.0                                         |\";
h5=\"|   \\\\  /   A nd | Web: www.OpenFOAM.org                                      |\";
h6=\"|    \\\\/    M anipulation |                                                   |\";
h7=\"\\*---------------------------------------------------------------------------*/\";

h8=\"FoamFile\";
h9=\"{\";
h10=\"    version     2.0;\";
h11=\"    format      ascii;\";
h12=\"    class       dictionary;\";
h13=\"    object      blockMeshDict;\";
h14=\"}\";
h15=\"// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //\";

fprintf(f1,\'%s \\n\',h1);
fprintf(f1,\'%s \\n\',h2);
fprintf(f1,\'%s \\n\',h3);
fprintf(f1,\'%s \\n\',h4);
fprintf(f1,\'%s \\n\',h5);
fprintf(f1,\'%s \\n\',h6);
fprintf(f1,\'%s \\n\',h7);
fprintf(f1,\'%s \\n\',h8);
fprintf(f1,\'%s \\n\',h9);
fprintf(f1,\'%s \\n\',h10);
fprintf(f1,\'%s \\n\',h11);
fprintf(f1,\'%s \\n\',h12);
fprintf(f1,\'%s \\n\',h13);
fprintf(f1,\'%s \\n\',h14);
fprintf(f1,\'%s \\n\',h15);

conv=\"convertToMeters 1;\";
ver=\"vertices\";
block=\"hex (0 1 2 0 3 4 5 3) (50 1 200) simpleGrading (0.2 1 1)\";
axis = \"(0 3 3 0)\";
wall= \"(2 5 4 1)\";
inlet= \"(0 0 2 1)\";
outlet= \"(3 3 4 5)\";
wedgeFront= \"(0 3 5 2)\";
wedgeBack= \"(0 1 4 3)\";

fprintf(f1,\'%s \\n\',conv);
fprintf(f1,\' %s\\n\',ver);
fprintf(f1,\'%s \\n\',\'(\');

fprintf(f1,\'%s (%d %d %d) \\n\',s1,V0);
fprintf(f1,\'%s (%d %d %d) \\n\',s1,V1);
fprintf(f1,\'%s (%d %d %d) \\n\',s1,V2);
fprintf(f1,\'%s (%d %d %d) \\n\',s1,V3);
fprintf(f1,\'%s (%d %d %d) \\n\',s1,V4);
fprintf(f1,\'%s (%d %d %d) \\n\',s1,V5);
fprintf(f1,\'%s \\n\',\');\');

fprintf(f1,\'\\n%s\\n\',\'blocks\');
fprintf(f1,\'%s \\n\',\'(\');
fprintf(f1,\'%s %s\\n\',s1,block);
fprintf(f1,\'%s \\n\',\');\');

fprintf(f1,\'\\n%s\\n\',\'edges\');
fprintf(f1,\'%s\\n\',\'(\');
fprintf(f1,\'\\narc %d %d(%d %d %d)\',1,2,0,r,0);
fprintf(f1,\'\\narc %d %d(%d %d %d)\',4,5,L,r,0);
fprintf(f1,\'\\n%s\\n\',\');\');

fprintf(f1,\'\\n%s\\n\',\'boundary\');
fprintf(f1,\'%s\\n\',\'(\');
fprintf(f1,\'%s %s\\n\',s1,\'axis\');
fprintf(f1,\'%s %s\\n\',s1,\'{\');
fprintf(f1,\'%s %s\\n\',s2,\'type empty;\');
fprintf(f1,\'%s %s\\n\',s2,\'faces\');
fprintf(f1,\'%s %s\\n\',s2,\'(\');
fprintf(f1,\'%s %s\\n\',s3,axis);
fprintf(f1,\'%s %s\\n\',s2,\');\');
fprintf(f1,\'%s %s\\n\',s1,\'}\');
fprintf(f1,\'%s %s\\n\',s1,\'wall\');
fprintf(f1,\'%s %s\\n\',s1,\'{\');
fprintf(f1,\'%s %s\\n\',s2,\'type wall;\');
fprintf(f1,\'%s %s\\n\',s2,\'faces\');
fprintf(f1,\'%s %s\\n\',s2,\'(\');
fprintf(f1,\'%s %s\\n\',s3,wall);
fprintf(f1,\'%s %s\\n\',s2,\');\');
fprintf(f1,\'%s %s\\n\',s1,\'}\');
fprintf(f1,\'%s %s\\n\',s1,\'inlet\');
fprintf(f1,\'%s %s\\n\',s1,\'{\');
fprintf(f1,\'%s %s\\n\',s2,\'type patch;\');
fprintf(f1,\'%s %s\\n\',s2,\'faces\');
fprintf(f1,\'%s %s\\n\',s2,\'(\');
fprintf(f1,\'%s %s\\n\',s3,inlet);
fprintf(f1,\'%s %s\\n\',s2,\');\');
fprintf(f1,\'%s %s\\n\',s1,\'}\');
fprintf(f1,\'%s %s\\n\',s1,\'outlet\');
fprintf(f1,\'%s %s\\n\',s1,\'{\');
fprintf(f1,\'%s %s\\n\',s2,\'type patch;\');
fprintf(f1,\'%s %s\\n\',s2,\'faces\');
fprintf(f1,\'%s %s\\n\',s2,\'(\');
fprintf(f1,\'%s %s\\n\',s3,outlet);
fprintf(f1,\'%s %s\\n\',s2,\');\');
fprintf(f1,\'%s %s\\n\',s1,\'}\');
fprintf(f1,\'%s %s\\n\',s1,\'wedgeFront\');
fprintf(f1,\'%s %s\\n\',s1,\'{\');
fprintf(f1,\'%s %s\\n\',s2,\'type wedge;\');
fprintf(f1,\'%s %s\\n\',s2,\'faces\');
fprintf(f1,\'%s %s\\n\',s2,\'(\');
fprintf(f1,\'%s %s\\n\',s3,wedgeFront);
fprintf(f1,\'%s %s\\n\',s2,\');\');
fprintf(f1,\'%s %s\\n\',s1,\'}\');
fprintf(f1,\'%s %s\\n\',s1,\'wedgeBack\');
fprintf(f1,\'%s %s\\n\',s1,\'{\');
fprintf(f1,\'%s %s\\n\',s2,\'type wedge;\');
fprintf(f1,\'%s %s\\n\',s2,\'faces\');
fprintf(f1,\'%s %s\\n\',s2,\'(\');
fprintf(f1,\'%s %s\\n\',s3,wedgeBack);
fprintf(f1,\'%s %s\\n\',s2,\');\');
fprintf(f1,\'%s %s\\n\',s1,\'}\');
fprintf(f1,\'%s %s\\n\',\');\');
fprintf(f1,\'\\n \\n%s\',\'mergePatchPairs\');
fprintf(f1,\'\\n%s\',\'(\');
fprintf(f1,\'\\n%s\',\');\');
fprintf(f1,\'\\n%s\\n\\n\',\'// ************************************************************************* //\');

fclose(f1);

blockMeshDict file:

/*--------------------------------*- C++ -*----------------------------------*\\ 
| ========= |                                                                 | 
| \\\\      / F ield | OpenFOAM: The Open Source CFD Toolbox                    | 
|  \\\\    /  O peration | Version: 4.0                                         | 
|   \\\\  /   A nd | Web: www.OpenFOAM.org                                      | 
|    \\\\/    M anipulation |                                                   | 
\\*---------------------------------------------------------------------------*/ 
FoamFile 
{ 
    version     2.0; 
    format      ascii; 
    class       dictionary; 
    object      blockMeshDict; 
} 
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // 
convertToMeters 1; 
 vertices
( 
     (0 0 0) 
     (0 4.998287e-03 -1.308847e-04) 
     (0 4.998287e-03 1.308847e-04) 
     (1.260000e+00 0 0) 
     (1.260000e+00 4.998287e-03 -1.308847e-04) 
     (1.260000e+00 4.998287e-03 1.308847e-04) 
); 

blocks
( 
     hex (0 1 2 0 3 4 5 3) (50 1 200) simpleGrading (0.2 1 1)
); 

edges
(

arc 1 2(0 5.000000e-03 0)
arc 4 5(1.260000e+00 5.000000e-03 0)
);

boundary
(
     axis
     {
         type empty;
         faces
         (
             (0 3 3 0)
         );
     }
     wall
     {
         type wall;
         faces
         (
             (2 5 4 1)
         );
     }
     inlet
     {
         type patch;
         faces
         (
             (0 0 2 1)
         );
     }
     outlet
     {
         type patch;
         faces
         (
             (3 3 4 5)
         );
     }
     wedgeFront
     {
         type wedge;
         faces
         (
             (0 3 5 2)
         );
     }
     wedgeBack
     {
         type wedge;
         faces
         (
             (0 1 4 3)
         );
     }
); 
 
mergePatchPairs
(
);
// ************************************************************************* //

Boundary conditions:

Pressure (p):

/*--------------------------------*- C++ -*----------------------------------*\\
| =========                 |                                                 |
| \\\\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\\\    /   O peration     | Version:  5                                     |
|   \\\\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\\\/     M anipulation  |                                                 |
\\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      p;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 2 -2 0 0 0 0];

internalField   uniform 0;

boundaryField
{
    axis
    {
        type            empty;
    }

    wall
    {
        type            zeroGradient;
    }

    inlet
    {
        type            zeroGradient;
    }

    outlet
    {
        type            fixedValue;
        value           uniform 0.0671;
    }

    wedgeFront
    {
        type            wedge;
    }

    wedgeBack
    {
        type            wedge;
    }
}

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

Velocity(U):

/*--------------------------------*- C++ -*----------------------------------*\\
| =========                 |                                                 |
| \\\\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\\\    /   O peration     | Version:  5                                     |
|   \\\\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\\\/     M anipulation  |                                                 |
\\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       volVectorField;
    object      U;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 1 -1 0 0 0 0];

internalField   uniform (0.1869 0 0);

boundaryField
{
    
    axis
    {
        type            empty;
    }

    wall
    {
        type            fixedValue;
        value           uniform (0 0 0);
    }

    inlet
    {
        type            fixedValue;
        value           uniform (0.1869 0 0);
    }

    outlet
    {
        type            zeroGradient;
    }

    wedgeFront
    {
        type            wedge;
    }

    wedgeBack
    {
        type            wedge;
    }
}

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

controlDict file

/*--------------------------------*- C++ -*----------------------------------*\\
| =========                 |                                                 |
| \\\\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\\\    /   O peration     | Version:  5                                     |
|   \\\\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\\\/     M anipulation  |                                                 |
\\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    \"system\";
    object      controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

application     icoFoam;

startFrom       startTime;

startTime       0;

stopAt          endTime;

endTime         5;

deltaT          0.01;

writeControl    timeStep;

writeInterval   100;

purgeWrite      0;

writeFormat     ascii;

writePrecision  6;

writeCompression off;

timeFormat      general;

timePrecision   6;

runTimeModifiable true;


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

transportProperties:

/*--------------------------------*- C++ -*----------------------------------*\\
| =========                 |                                                 |
| \\\\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\\\    /   O peration     | Version:  5                                     |
|   \\\\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\\\/     M anipulation  |                                                 |
\\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    \"constant\";
    object      transportProperties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

nu              [0 2 -1 0 0 0 0] 0.00000089267;


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

Grading factor of 0.2, Wedge angle 3° :

Wedge inlet section:

Wedge outlet section:

Velocity distribution:

At x=0.005m

At x=0.01m

At x=0.02m

At x=0.05m

At x=0.25m

At x=0.5m

At x=0.75m

At x=1m

At x=1.25m

Wedge angle 4°

Wedge angle 5°

Fully developed flow, velocity profile:

After getting the single line graph plot, we need to export all the velocity(U_magnitude) to \'.csv\' spreadsheet format.

By using the file reading parse method to plot the fully developed velocity profile & shear stress. 

clear all
close all
clc


f2=fopen(\'shear_stress.txt\',\'r\');
fgetl(f2);

A=fscanf(f2,\'%f\',[2 inf]);
B=A\';

mu=0.00089;
U=B(1:end-1,2);
r=B(1:end-1,1);

tau = zeros(length(U));

for i=2:length(r)
    tau(i)=mu*(U(i-1)-U(i))/(r(i)-r(i-1));
end

tau=tau(:,1);
figure(1)
plot(r,tau,\'b\',\'linewidth\',1.5)
hold on
plot(-r,tau,\'b\',\'linewidth\',1.5)
grid on
xlabel(\'Radius in m\');
ylabel(\'Shear stress in Pa\')
title(\'Shear stress profile at L=1.26m\')

figure(2)
plot(U,r,\'b\',\'linewidth\',1.5)
hold on
plot(U,-r,\'b\',\'linewidth\',1.5)
grid on
xlabel(\'Velocity in m/s\')
ylabel(\'Radius in m\')
title(\'Velocity profile at L=1.26m\')

Fully developed velocity profile:

From the above plot, we can see that the fully developed flow velocity profile is parabolic as per the theoretically discussed.

Shear stress inside pipe flow:

Form the above plot we can observe that the shear stress is increasing at the walls compared with the center of the pipe.

Conclusion:

From the above result, for a wedge angle of 5 degrees show that the less error value, And comparatively low execution time & courant number.

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.