Challenge

Week 11 - Simulation of Flow through a pipe in OpenFoam

Hagen- Poiseuille's equations

* Hagen–Poiseuille law, Poiseuille law or Poiseuille equation, is a physic law that gives the pressure drop in an incompressable and newtonian fluid in laminar flow flowing through a long cylindrical pipe of constant cross section

* The assumptions of the equation are that the fluid is incompressable and newtonian the flow is laminar through a pipe of constant circular cross-section that is substantially longer than its diameter; and there is no acceleration of fluid in the pipe. For velocities and pipe diameters above a threshold, actual fluid flow is not laminar but turbulent, leading to larger pressure drops than calculated by the Hagen–Poiseuille equation.

  • Δp is the pressure difference between the two ends,
  • L is the length of pipe,
  • μ is the Dynamic viscosity,
  • Q is the volumetric flow rate,
  • R is the pipe radius,
  • A is the cross section of pipe.

* The layer in which the share viscous flow are significant is called boundary layer this boundary layer is important concept its divide the flow pipe into two region boundary layer region and irrotational flow region

* When the fluid just enter the pipe the thickness of a boundary layer gradually increases from zero moving in direction of flow and eventually reaches the pipe centre and fills the entire pipe the region from the entrance of the pipe to the point where the boundary layer covers the pipe is termed as a hydrodynamic region and length of a pipe in the region is termed as a hydro dynamic and length in this region the velocity profile develops and the floor is called hudrodynamically developing flow

Given

  • Reynolds number :2100
  • Working fluid :water
  • Fluid temperature :25celcus
  • Diameter of pipe :0.02
  • dynamic viscosity : 0.89*10^-3
  • kinematics viscosity :0.89*10^-6
  • difference between pressure :32*L*mu*avg velocity/D^2
  • Theta :4
  • length 2.25m

Analytical solution

  • average velocity :re*mu/rho*D = 0.09373
  • maximum velocity :2*avg velocity = 0.1874ms
  • pressure drop :16.187pa
  • kinematics pressure drop : delP/P = 0.0168m

Matlab Code for creating blockmesh dic:

close all
clear all
clc

s = 1;                      %scaling factor
L = 5;                      %Length of pipe
r = 0.015;                  %Radius of pipe
theta = 2;                  %Angle of wedge wrt axis
rho = 997;                  %density of working fluid
Re =2100;                   %Renolds number
dx = 20;                    %Number of cells in X-direction
dy = 1;                     %Number of cells in y-direction
dz = 500;                   %Number of cells in z-direction
Xsg = 0.3;                   %Simple grading along X-direction
Ysg = 1;                     %Simple grading along Y-direction
Zsg = 1;                     %Simple grading along Z-direction
t = 210;                    %End time
dt = 0.001;                 %Time step
i =50;                      %writing interval
mu = 0.997*10^(-3);         %dynamic viscosity
FOAMscripts(s,L,r,theta,rho,Re,dx,dy,dz,Xsg,Ysg,Zsg,t,dt,i,mu)

Function code to create all required files:

function FOAMscripts(s,L,r,theta,rho,Re,dx,dy,dz,Xsg,Ysg,Zsg,t,dt,i,mu)
nue = mu/rho;               %kinematic viscosity
v = (Re*nue)/(2*r);         %Velocity
x1 = r*cosd(theta);         
y1 = r*sind(theta);
z1 = 0;
x2 = r*cosd(360-theta);
y2 = r*sind(360-theta);
z2 = L;
%common header lines

 line1='/*--------------------------------*-C++ -*--------------------------------*';
 line2='  =========                 |                                               ';
 line3='        /  F ield         |  OpenFOAM: The Open Source CFD Toolbox        ';
 line4='       /   O peration     |  Website: https://openfoam.org                ';
 line5='      /    A nd           |  Version: 8                                   ';
 line6='     /     M anipulation  |                                               ';
 line7='*-------------------------------------------------------------------------*/';
 line8='FoamFile';
 line9='// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //';

 mkdir 'Flow_through_a_pipe';
 mkdir 'Flow_through_a_pipe' 'system';
 mkdir 'Flow_through_a_pipe' 'constant';
 mkdir 'Flow_through_a_pipe' '0';
 %% blockmeshdict
 f=fopen('Flow_through_a_pipesystemblockMeshDict','w');
 fprintf(f,'%sn%sn%sn%sn%sn%sn%sn%sn{n%sn%sn%sn',line1,line2,line3,line4,line5,line6,line7,line8);
 fprintf(f,'tversion t2.0;ntformat tascii;ntclass ttdictionary;ntobject tblockMeshDict;n}n');
 %starting line
 fprintf(f,'%s',line9);
 %defining scaling factor
 fprintf(f,'nnconvertToMeters %d;nn',s);
 %defining vertices
 fprintf(f,'verticesn(nnt(0 0 0)nt(%f %f %f)nt(%f %f %f)nt(0 0 %f)nt(%f %f %f)nt(%f %f %f)n);n',x1,y1,z1,x1,y1,z2,z2,x2,y2,z1,x2,y2,z2);
 %defining blocks
 fprintf(f,'nblocksn(nthex (0 4 1 0 3 5 2 3) (%d %d %d) simpleGrading (%f %f %f)n);n',dx,dy,dz,Xsg,Ysg,Zsg);
 %defining edges
 fprintf(f,'nedgesn(ntarc 1 4 (%f %f %f)ntarc 2 5 (%f %f %f)n);n',x1/2,y1/2,z1,x2/2,y2/2,z2);
 %defining boundaries
 fprintf(f,'nboundaryn(nntaxisnt{ntttype empty;nttfacesntt(nttt(0 3 3 0)ntt);nt}n');
 fprintf(f,'tinletnt{ntttype patch;nttfacesntt(nttt(0 4 1 0)ntt);nt}n');
 fprintf(f,'toutletnt{ntttype patch;nttfacesntt(nttt(3 5 2 3)ntt);nt}n');
 fprintf(f,'tleftwallnt{ntttype wedge;nttfacesntt(nttt(4 0 3 5)ntt);nt}n');
 fprintf(f,'trightwallnt{ntttype wedge;nttfacesntt(nttt(0 1 2 3)ntt);nt}n');
 fprintf(f,'touterwallnt{ntttype wall;nttfacesntt(nttt(1 4 5 2)ntt);nt}nnn);n');
 %defining patch pairs
 fprintf(f,'nmergePatchPairsn(n);n');
 %ending line
 fprintf(f,'n%s',line9);
 fclose(f);
 %% control dict 
 f=fopen('Flow_through_a_pipesystemcontrolDict','w');
 fprintf(f,'%sn%sn%sn%sn%sn%sn%sn%sn{n%sn%sn%sn',line1,line2,line3,line4,line5,line6,line7,line8);
 fprintf(f,'tversion t2.0;ntformat ttascii;ntclass ttdictionary;ntlocationt"system";ntobject ttcontrolDict;n}n');
 %starting line
 fprintf(f,'%sn',line9);
 fprintf(f,'napplicationticoFoam;nnstartFromtstartTime;nnstartTimet0;nnstopAtttendTime;nnendTimett%f;nndeltaTtt%fs;',t,dt);
 fprintf(f,'nnwriteControlttimeStep;nnwriteIntervalt%d;nnpurgeWritet0;nnwriteFormattascii;nnwritePrecisiont6;nnwriteCompressiontoff;nntimeFormattgeneral;nntimePrecisiont6;nnrunTimeModifiablettrue;',i);
 %ending line
 fprintf(f,'nn%s',line9);
 fclose(f);
 %% fvSchemes
 f=fopen('Flow_through_a_pipesystemfvSchemes','w');
 fprintf(f,'%sn%sn%sn%sn%sn%sn%sn%sn{n%sn%sn%sn',line1,line2,line3,line4,line5,line6,line7,line8);
 fprintf(f,'tversion t2.0;ntformat ttascii;ntclass ttdictionary;ntlocationt"system";ntobject ttfvSchemes;n}n');
 %starting line
 fprintf(f,'%sn',line9);
 fprintf(f,'nddtSchemesn{ntdefaultttEuler;n}nngradSchemesn{ntdefaultttGauss linear;ntgrad(p)ttGauss linear;n}nndivSchemesn{ntdefaultttnone;ntdiv(phi,U)tGauss linear;n}nn');
 fprintf(f,'laplacianSchemesn{ntdefaultttGauss linear orthogonal;n}nninterpolationSchemesn{ntdefaultttlinear;n}nnsnGradSchemesn{ntdefaultttorthogonal;n}');
 %ending line
 fprintf(f,'nn%s',line9);
 fclose(f);
 %% fvSolution
 f=fopen('Flow_through_a_pipesystemfvSolution','w');
 fprintf(f,'%sn%sn%sn%sn%sn%sn%sn%sn{n%sn%sn%sn',line1,line2,line3,line4,line5,line6,line7,line8);
 fprintf(f,'tversion t2.0;ntformat ttascii;ntclass ttdictionary;ntlocationt"system";ntobject ttfvSolution;n}n');
 %starting line
 fprintf(f,'%sn',line9);
 fprintf(f,'nsolversn{ntpnt{nttsolverttPCG;nttpreconditionertDIC;ntttolerancet1e-06;nttrelToltt0.05;nt}n');
 fprintf(f,'ntpFinalnt{ntt$p;nttrelToltt0;nt}nntUnt{nttsolverttsmoothSolver;nttsmoothertsymGaussSeidel;ntttolerancet1e-05;nttrelToltt0;nt}n}n');
 fprintf(f,'nPISOn{ntnCorrectorst2;ntnNonOrthogonalCorrectorst0;ntpRefCellt0;ntpRefValuet0;n}');
 %ending line
 fprintf(f,'nn%s',line9);
 fclose(f);
 %% transport properties
 f=fopen('Flow_through_a_pipeconstanttransportProperties','w');
 fprintf(f,'%sn%sn%sn%sn%sn%sn%sn%sn{n%sn%sn%sn',line1,line2,line3,line4,line5,line6,line7,line8);
 fprintf(f,'tversion t2.0;ntformat ttascii;ntclass ttdictionary;ntlocationt"constant";ntobject tttransportProperties;n}n');
 %starting line
 fprintf(f,'%sn',line9);
 %defining kinematic viscosity of water
 fprintf(f,'nnutt[0 2 -1 0 0 0 0]  %f;n',nue);
 %ending line
 fprintf(f,'n%s',line9);
 fclose(f);
 %% pressure
 f=fopen('Flow_through_a_pipep','w');
 fprintf(f,'%sn%sn%sn%sn%sn%sn%sn%sn{n%sn%sn%sn',line1,line2,line3,line4,line5,line6,line7,line8);
 fprintf(f,'tversion t2.0;ntformat ttascii;ntclass ttvolScalarField;ntobject ttP;n}n');
 %starting line
 fprintf(f,'%sn',line9);
 %defining pressure conditions
 fprintf(f,'ndimensionst[0 2 -2 0 0 0 0];nninternalFieldtuniform 0;nnboundaryFieldn{n');
 fprintf(f,'tinletnt{ntttypetzeroGradient;nt}ntoutletnt{ntttypetfixedValue;nttvaluetuniform 0;nt}');
 fprintf(f,'ntaxisnt{ntttypetempty;nt}ntleftwallnt{ntttypetwedge;nt}ntrightwallnt{ntttypetwedge;nt}ntouterwallnt{ntttypetzeroGradient;nt}nn}');
 %ending line
 fprintf(f,'n%s',line9);
 fclose(f);
 %% velocity
 f=fopen('Flow_through_a_pipeU','w');
 fprintf(f,'%sn%sn%sn%sn%sn%sn%sn%sn{n%sn%sn%sn',line1,line2,line3,line4,line5,line6,line7,line8);
 fprintf(f,'tversion t2.0;ntformat ttascii;ntclass ttvolVectorField;ntobject ttU;n}n');
 %starting line
 fprintf(f,'%sn',line9);
 %defining velocity conditions
 fprintf(f,'ndimensionst[0 1 -1 0 0 0 0];nninternalFieldtuniform (0 0 0);nnboundaryFieldn{n');
 fprintf(f,'toutletnt{ntttypetzeroGradient;nt}ntinletnt{ntttypetfixedValue;nttvaluetuniform (0 0 %f);nt}',v);
 fprintf(f,'ntaxisnt{ntttypetempty;nt}ntleftwallnt{ntttypetwedge;nt}ntrightwallnt{ntttypetwedge;nt}ntouterwallnt{ntttypetnoSlip;nt}nn}');
 %ending line
 fprintf(f,'n%s',line9);
 fclose(f);

Excuiting the code blockmesh Dic:

/*--------------------------------*- 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.000000 0.009994 -0.000349)
    (0.000000 0.009994 0.000349)
    (2.520000 0.000000 0.000000)
    (2.520000 0.009994 -0.000349)
    (2.520000 0.009994 0.000349)
  
);

blocks

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

edges
(
    
                  arc 1 2 (0.000000 0.010000 0.000000)
                  arc 4 5 (2.520000 0.010000 0.000000)
                 
                 

);

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

mergePatchPairs
(
);
/*--------------------------------*- 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         10;

deltaT          0.01;

writeControl    timeStep;

writeInterval   20;

purgeWrite      0;

writeFormat     ascii;

writePrecision  6;

writeCompression off;

timeFormat      general;

timePrecision   6;

runTimeModifiable true;


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

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


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

 

/*--------------------------------*- C++ -*----------------------------------*\
/*--------------------------------*- 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.0937312 0 0);
    }

   outlet
    {
        type            zeroGradient;
    }

    pipe
    {
        type           noSlip;
       
    }
   
    axis
    {
    type                empty;
   
   }
  
   front
   {
     type                    wedge;
     }
    
   back
   {
    type                        wedge;
    } 
}

// ************************************************************************* //
 
/*--------------------------------*- 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
{
  axis
    {
        type            empty;
    }

   pipe
    {
        type            zeroGradient;
    }
       
    inlet
    {
      type            zeroGradient;
      }
     
    outlet
    {
       type             fixedValue;
       value            uniform 0.0168682;
       }
    
     front
     {
     type                 wedge;
     }
    
     back
      {
      type                  wedge;
      } 
}

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

* Hear you can notice that the velocity turns out zero at the wall it's mainly because of their is no slip condition as slow progress .

* gradually the boundary layer develops and the adjacent there tends to slow down because the velocity of the centre increases centre because of the pressure difference obtained from the equation velocity profile matches with the help of Hagen- Poiseuille's equations.

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.