Challenge

Week 11 - Simulation of Flow through a pipe in OpenFoam

HHarsha Villuri15 Sept 2020

The length of the entrance is the distance that fluid travels after entering a pipe before the full development of the flow. Entrance length refers to the length of the entry region, effected the pipe's interior wall spread as an expanding boundary layer into the flow. When the boundary layer extends to fill the entire pipe, the flow grows into a fully developed flow where the flow characteristics no longer change with increased distance along the pipe.

Hydrodynamic entrance:

Hydrodynamic entrance length describes the formation of a velocity profile caused by viscous forces propagating from the pipe wall. This region has a non-uniform flow. The fluid enters a pipe at uniform velocity, then the fluid particles make contact with the surface of the pipe come to a complete stop due to no-slip condition. The contact layer with the pipe surface is resistant to the acceleration of neighboring layers due to viscous forces within the fluid and slowly slows down the adjacent layer of the fluid. The velocity of the fluid in the center of the pipe increases to compensate for the reduced velocity of the fluid near the pipe. The boundary layer divides the flow in a pipe into 2 regions. They are the boundary layer region and the irrotational flow region.

The length of the hydrodynamic entry region along the pipe is known as the hydrodynamic entry length. This length for laminar flow is calculated by:

Le = 0.05*Re*D

where Re is Reynold's number and D is the diameter of the pipe.

Reynold's number, Re = 2100

Diameter, D = 0.005m

Length of the pipe, L = 1m

Radius, r = 0.0025m

Dynamic viscosity of water at 30c, mu = 0.798e-3 N s/m2

kinematic viscosity, v = 0.801e-6 m2/s

Hagen-poiseuille's law:

In nonideal fluid dynamics, the poiseuille equation is a physical law that gives the pressure drop in an incompressible and Newtonian fluid in laminar flow through a long cylindrical pipe of a constant cross-section.

Average velocity, `V_(avg) = (R_e*mu)/(rho*D)` = 0.33642 m/s

Maximum velocity, `V_max = 2*V_(avg)` = 0.67284 m/s

Pressure drop, `DeltaP = (32*mu*L*V_(avg))/D^2` = 0.6872 Pa/m

Matlab Code:

clear all
close all
clc

Re = 2100;   % reynolds number
theta = 4,   % angle of wedge
D = 0.005;   % Diameter of the pipe 
L = 1;       % length of the pipe 
Le = 0.05*Re*D;  % length of hydrodynamic enterance length  
r = D/2;      % radius of pipe
mu = 0.798e-3; % dynamic viscosity 
v = 0.801e-6;  % kinematic viscosity
V = (v*Re)/D;  % velocity
R = 0.0025;    % radius of pipe 
theta = 4;     % wedge angle
v_max = 2*0.33642;

% calculating shear stress
R = linspace(0,D/2,500);

for i = 1:length(R)
   tau(i) = 2*mu*v_max*(abs(R(i)))/r^2;
endfor

figure 1
plot(R,tau)
xlabel('radius')
ylabel('shear stress')
grid on
legend('shear stress')
title('shear stress')


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

fprintf(f1,'%s\n','/*--------------------------------*- C++ -*----------------------------------*\');
fprintf(f1,'%s\n','  =========                 |');
fprintf(f1,'%s\n','  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox');
fprintf(f1,'%s\n','   \\    /   O peration     | Website:  https://openfoam.org');
fprintf(f1,'%s\n','    \\  /    A nd           | Version:  8');
fprintf(f1,'%s\n','     \\/     M anipulation  |');
fprintf(f1,'%s\n','\*---------------------------------------------------------------------------*/');
fprintf(f1,'%s\n','FoamFile');
fprintf(f1,'%s\n','{');
fprintf(f1,'%s\n','    version     2.0;');
fprintf(f1,'%s\n','    format      ascii;');
fprintf(f1,'%s\n','    class       dictionary;');
fprintf(f1,'%s\n','    object      blockMeshDict;');
fprintf(f1,'%s\n','}');
fprintf(f1,'%s\n','// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //');
fprintf(f1,'%s\n','convertToMeters 1;');

fprintf(f1,'%s\n','vertices');
fprintf(f1,'%s\n','(');
fprintf(f1,'\t(%d %d %d)  //V1\n', 0,0,0);
fprintf(f1,'\t(%d %d %d)  //V2\n', 0, R*cosd(theta/2), R*sind(theta/2));
fprintf(f1,'\t(%d %d %d)  //V3\n', 0, R*cosd(theta/2), -R*sind(theta/2));
fprintf(f1,'\t(%d %d %d)  //V4\n', 1, 0, 0);
fprintf(f1,'\t(%d %d %d)  //V5\n', 1, R*cosd(theta/2), R*sind(theta/2));
fprintf(f1,'\t(%d %d %d)  //V6\n', 1, R*cosd(theta/2), -R*sind(theta/2));
fprintf(f1,'%s\n',');');
fprintf(f1,'%s\n','blocks');
fprintf(f1,'%s\n','(');
fprintf(f1,'%s\n','    hex (0 0 2 1 3 3 5 4) (1 20 200) simpleGrading (1 0.2 1)');
fprintf(f1,'%s\n',');');
fprintf(f1,'%s\n','edges');
fprintf(f1,'%s\n','(');
fprintf(f1,'%s\n',');');
fprintf(f1,'%s\n','boundary');
fprintf(f1,'%s\n','(');
fprintf(f1,'%s\n','   inlet');
fprintf(f1,'%s\n','    {');
fprintf(f1,'%s\n','        type patch;');
fprintf(f1,'%s\n','        faces');
fprintf(f1,'%s\n','        (');
fprintf(f1,'%s\n','            (0 1 2 0)');
fprintf(f1,'%s\n','        );');
fprintf(f1,'%s\n','    }');
fprintf(f1,'%s\n','    outlet');
fprintf(f1,'%s\n','    {');
fprintf(f1,'%s\n','        type patch;');
fprintf(f1,'%s\n','        faces;');
fprintf(f1,'%s\n','        (');
fprintf(f1,'%s\n','            (3 5 4 3)');
fprintf(f1,'%s\n','        );');
fprintf(f1,'%s\n','    }');
fprintf(f1,'%s\n','    frontwedge');
fprintf(f1,'%s\n','    {');
fprintf(f1,'%s\n','        type wedge;');
fprintf(f1,'%s\n','        faces');
fprintf(f1,'%s\n','        (');
fprintf(f1,'%s\n','            (1 0 3 4)');
fprintf(f1,'%s\n','        );');
fprintf(f1,'%s\n','    }');
fprintf(f1,'%s\n','    backwedge');
fprintf(f1,'%s\n','    {');
fprintf(f1,'%s\n','        type wedge;');
fprintf(f1,'%s\n','        faces');
fprintf(f1,'%s\n','        (');
fprintf(f1,'%s\n','            (0 2 5 3)');
fprintf(f1,'%s\n','        );');
fprintf(f1,'%s\n','    }');
fprintf(f1,'%s\n','    topwall');
fprintf(f1,'%s\n','    {');
fprintf(f1,'%s\n','        type wall;');
fprintf(f1,'%s\n','        faces');
fprintf(f1,'%s\n','        (');
fprintf(f1,'%s\n','            (1 4 5 2)');
fprintf(f1,'%s\n','        );');
fprintf(f1,'%s\n','    }');
fprintf(f1,'%s\n','    bottomwall');
fprintf(f1,'%s\n','    {');
fprintf(f1,'%s\n','        type empty;');
fprintf(f1,'%s\n','        faces');
fprintf(f1,'%s\n','        (');
fprintf(f1,'%s\n','            (0 3 3 0)');
fprintf(f1,'%s\n','        );');
fprintf(f1,'%s\n','    }');
fprintf(f1,'%s\n',');');
fprintf(f1,'%s\n','');
fprintf(f1,'%s\n','mergePatchPairs');
fprintf(f1,'%s\n','(');
fprintf(f1,'%s\n',');');
fprintf(f1,'%s\n','// ************************************************************************* //');

fclose(f1);

Mesh:

blockMeshDict:

/*--------------------------------*- 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) //V1
    (0 2.498477e-3 8.724874e-5) //V2
    (0 2.498477e-3 -8.724874e-5) //V3
    (1 0 0) //V4
    (1 2.498477e-3 8.724874e-5) //V5
    (1 2.498477e-3 -8.724874e-5) //V6
);

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

edges
(
);

boundary
(
   inlet
    {
        type patch;
        faces
        (
            (0 1 2 0)
        );
    }
    outlet
    {
        type patch;
        faces
        (
            (3 5 4 3)          
        );
    }
    frontwedge
    {
        type wedge;
        faces
        (
            (1 0 3 4)
        );
    }
    backwedge
    {
        type wedge;
        faces
        (
            (0 2 5 3)
            
        );
    }
    topwall
    {
        type wall;
        faces
        (
            (1 4 5 2)
        );
    }
    bottomwall
    {
        type empty;
        faces
        (
            (0 3 3 0)          
        );
    }    
);

mergePatchPairs
(
);

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

controlDict:

/*--------------------------------*- 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.5;

stopAt          endTime;

endTime         10;

deltaT          0.001;

writeControl    timeStep;

writeInterval   20;

purgeWrite      0;

writeFormat     ascii;

writePrecision  6;

writeCompression off;

timeFormat      general;

timePrecision   6;

runTimeModifiable true;


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

pressure:

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

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

internalField   uniform 0;

boundaryField
{
    inlet
    {
        type            zeroGradient;
    }
    outlet
    {
        type            fixedValue;
        value           uniform 0;   
    }
    frontwedge
    {
        type            wedge;
    }
    backwedge
    {
        type            wedge;
    }
    topwall
    {
        type            zeroGradient;          
    }
    bottomwall
    {
        type            empty;
    }
}

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

velocity:

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

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

internalField   uniform (0 0 0);

boundaryField
{
    inlet
    {
        type            fixedValue;
        value           uniform (0.33642 0 0);
    } 
    outlet
    {
        type            zeroGradient;   
    }
    frontwedge
    {
        type            wedge;
    }
    backwedge
    {
        type            wedge;
    }
    topwall
    {
        type            noSlip;          
    }
    bottomwall
    {
        type            empty;
    }
}

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

Thermoproperties:

/*--------------------------------*- 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    "constant";
    object      transportProperties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

nu              [0 2 -1 0 0 0 0] 0.801e-6;


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

Results:

Mesh:

Pressure mesh: