Week 11 - Simulation of Flow through a pipe in OpenFoam

Aim :- Simulation of Flow through a pipe in OpenFoam.

Objective :- Simulate an axi-symmetric flow by applying the wedge boundary condition. Validate results with the Hagen- Poiseuille's equations for the fully developed flow. Also write code in Matlab to automate the generation of blockMeshDict file. 

Theory :- 

Laminar Flow,

Laminar or streamline flow in pipes (or tubes) occurs when a fluid flows in parallel layers, with no disruption between the layers.  At low velocities, the fluid tends to flow without lateral mixing, and adjacent layers slide past one another. There are no cross-currents perpendicular to the direction of flow, nor eddies or swirls of fluids.  In laminar flow, the motion of the particles of the fluid is very orderly with all particles moving in straight lines parallel to the pipe walls.   Any lateral mixing (mixing at right angles to the flow direction) occurs by the action of diffusion between layers of the liquid.  Diffusion mixing can be slow however if the diameter of the pipe of the tube is small then this diffusive mixing can be very significant.

Turbulent Flow,

Turbulent flow is a flow regime characterized by chaotic property changes. This includes a rapid variation of pressure and flows velocity in space and time. Turbulent flow is a type of fluid (gas or liquid) flow in which the fluid undergoes irregular fluctuations, or mixing, in contrast to laminar flow, in which the fluid moves in smooth paths or layers. In turbulent flow, the speed of the fluid at a point is continuously undergoing changes in both magnitude and direction. 

Reynolds Number,

The Reynolds number is the ratio of inertial forces to viscous forces within a fluid that is subjected to relative internal movement due to different fluid velocities. A region where these forces change behaviour is known as a boundary layer, such as the bounding surface in the interior of a pipe.

The Reynolds number (`R_e`) helps predict flow patterns in different fluid flow situations. At low Reynolds numbers, flows tend to be dominated by laminar (sheet-like) flow, while at high Reynolds numbers flows tend to be turbulent.

Reynolds Number for a internal flow through circular pipe is given as, `R_e = (rhoV_(avg)D)/mu`

Reynolds number for flow through pipes 

`R_e <= 2300`  Laminar flow

`2300 <= R_e <= 4000`Transitional flow

`R_e >= 4000`  Turbulent Flow

Entrance Length

The entrance length is the distance a flow travels after entering a pipe before the flow becomes fully developed. Entrance length refers to the length of the entry region, the area following the pipe entrance where effects originating from the interior wall of the pipe propagate into the flow as an expanding boundary layer. When the boundary layer expands to fill the entire pipe, the developing flow becomes a fully developed flow, where flow characteristics no longer change with increased distance along the pipe. Many different entrance lengths exist to describe a variety of flow conditions. 

The region from the pipe inlet to the point at which the velocity profile is fully developed is called the hydrodynamic entrance region, and the length of this region is called the hydrodynamic entry length `L_h`.

 Velocity profile for Laminar flow

Hagen- Poiseuille's Equation,

`DeltaP = (8muLQ)/(piR^4)`where,

 `DeltaP` = pressure difference between the two ends

`L` = length of pipe

`mu` = dynamic viscosity

`Q` = volumetric flow rate                    

`R` = pipe radius                                 

Hagen Poiseuille states that for a specified flow rate, the pressure drops and thus the required pumping power is proportional to the length of the pipe and the viscosity of the fluid, but it is inversely proportional to the fourth power of the radius (or diameter) of the pipe. Therefore, the pumping power requirement for a laminar-flow piping system can be reduced by a factor of 16 by doubling the pipe diameter.

blockMeshDict File in MATLAB,

close all
clear all
clc

%Given Inputs
re=2100;    %Reynolds Number

%Assumptions
d=0.0254/2;   %Diamter of the Pipe
r=d/2;      %Radius of the Pipe

%theta=input("Enter Wedge Angle Less Than 5°=");
theta=2;

%Entry Length for Laminar Flow
l_entry=0.05*re*d;  %x-direction or flow direction
l=l_entry+1; %Assuming 1 extra for fully developed flow profile

%Finding the coordinates for the wedge
z=r*sind(theta/2);  %z-direction
y=r*cosd(theta/2);  %y-direction or radial incase of circular pipe
point0=[0 0 0];
point1=[0 y z];
point2=[0 y -z];
point3=[l 0 0];
point4=[l y z];
point5=[l y -z];

%Creating the blockMeshDict File
f= fopen('blockMeshDict','w');

fprintf(f,'%sn','/*--------------------------------*- C++ -*----------------------------------*');
fprintf(f,'%sn','  =========                 |');
fprintf(f,'%sn','  \      /  F ield         | OpenFOAM: The Open Source CFD Toolbox');
fprintf(f,'%sn','   \    /   O peration     | Website:  https://openfoam.org');
fprintf(f,'%sn','    \  /    A nd           | Version:  8');
fprintf(f,'%sn','     \/     M anipulation  |');
fprintf(f,'%sn','*---------------------------------------------------------------------------*/');
fprintf(f,'%sn','FoamFile');
fprintf(f,'%sn','{');
fprintf(f,'%sn','    version     2.0;');
fprintf(f,'%sn','    format      ascii;');
fprintf(f,'%sn','    class       dictionary;');
fprintf(f,'%sn','    object      blockMeshDict;');
fprintf(f,'%sn','}');
fprintf(f,'%sn','// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //');
fprintf(f,'%sn','');
fprintf(f,'%sn','convertToMeters 1;');
fprintf(f,'%sn','');
fprintf(f,'%sn','vertices');
fprintf(f,'%sn','(');
fprintf(f,'t(%f %f %f)n',point0(1),point0(2),point0(3));
fprintf(f,'t(%f %f %f)n',point1(1),point1(2),point1(3));
fprintf(f,'t(%f %f %f)n',point2(1),point2(2),point2(3));
fprintf(f,'t(%f %f %f)n',point3(1),point3(2),point3(3));
fprintf(f,'t(%f %f %f)n',point4(1),point4(2),point4(3));
fprintf(f,'t(%f %f %f)n',point5(1),point5(2),point5(3));

fprintf(f,'%sn',');');
fprintf(f,'%sn','');
fprintf(f,'%sn','blocks');
fprintf(f,'%sn','(');
fprintf(f,'t%sn','hex (0 3 5 2 0 3 4 1) (1000 20 1) simpleGrading (1 0.5 1)');
fprintf(f,'%sn',');');
fprintf(f,'%sn','');
fprintf(f,'%sn','edges');
fprintf(f,'%sn','(');
fprintf(f,'tarc 1 2 (0 %f 0)n',r);
fprintf(f,'tarc 4 5 (%f %f 0)n',l,r);
fprintf(f,'%sn',');');
fprintf(f,'%sn','');
fprintf(f,'%sn','boundary');
fprintf(f,'%sn','(');
fprintf(f,'t%sn','inlet');
fprintf(f,'t%sn','{');
fprintf(f,'tt%sn','type patch;');
fprintf(f,'tt%sn','faces');
fprintf(f,'tt%sn','(');
fprintf(f,'ttt%sn','(0 1 2 0)');
fprintf(f,'tt%sn',');');
fprintf(f,'t%sn','}');
fprintf(f,'t%sn','outlet');
fprintf(f,'t%sn','{');
fprintf(f,'tt%sn','type patch;');
fprintf(f,'tt%sn','faces');
fprintf(f,'tt%sn','(');
fprintf(f,'ttt%sn','(3 5 4 3)');
fprintf(f,'tt%sn',');');
fprintf(f,'t%sn','}');
fprintf(f,'t%sn','top');
fprintf(f,'t%sn','{');
fprintf(f,'tt%sn','type wall;');
fprintf(f,'tt%sn','faces');
fprintf(f,'tt%sn','(');
fprintf(f,'ttt%sn','(1 4 5 2)');
fprintf(f,'tt%sn',');');
fprintf(f,'t%sn','}');
fprintf(f,'t%sn','front');
fprintf(f,'t%sn','{');
fprintf(f,'tt%sn','type wedge;');
fprintf(f,'tt%sn','faces');
fprintf(f,'tt%sn','(');
fprintf(f,'ttt%sn','(0 2 5 3)');
fprintf(f,'tt%sn',');');
fprintf(f,'t%sn','}');
fprintf(f,'t%sn','back');
fprintf(f,'t%sn','{');
fprintf(f,'tt%sn','type wedge;');
fprintf(f,'tt%sn','faces');
fprintf(f,'tt%sn','(');
fprintf(f,'ttt%sn','(0 3 4 1)');
fprintf(f,'tt%sn',');');
fprintf(f,'t%sn','}');
fprintf(f,'t%sn','bottom');
fprintf(f,'t%sn','{');
fprintf(f,'tt%sn','type empty;');
fprintf(f,'tt%sn','faces');
fprintf(f,'tt%sn','(');
fprintf(f,'ttt%sn','(0 3 3 0)');
fprintf(f,'tt%sn',');');
fprintf(f,'t%sn','}');
fprintf(f,'%sn',');');
fprintf(f,'%sn','');
fprintf(f,'%sn','mergePatchPairs');
fprintf(f,'%sn','(');
fprintf(f,'%sn',');');
fprintf(f,'%sn','');
fprintf(f,'%sn','// ************************************************************************* //');

fclose(f);

 BlockMeshDict Flle,

/*--------------------------------*- C++ -*----------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | Website:  https://openfoam.org
    \\  /    A nd           | Version:  10
     \\/     M anipulation  |
\*---------------------------------------------------------------------------*/
FoamFile
{
    format      ascii;
    class       dictionary;
    object      blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

convertToMeters 1;

vertices
(
    /*0*/(0 0 0)
    /*1*/(0 0.009994 0.000349)
    /*2*/(0 0.009994 -0.000349)
    /*3*/(2.500000 0 0)
    /*4*/(2.500000 0.009994 0.000349)
    /*5*/(2.500000 0.009994 -0.00349)
    
);

blocks
(
    hex (0 2 1 0  3 5 4 3) (30 1 500) simpleGrading (0.2 1 1)
);

edges
(
    arc 2 1 (0 0.010000 0)
    arc 5 4 (2.500000 0.010000 0)
);

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

Simulation Run Step,

ControlDict File,

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

application     icoFoam;

startFrom       startTime;

startTime       0;

stopAt          endTime;

endTime         0.5;

deltaT          0.0002;

writeControl    timeStep;

writeInterval   5;

purgeWrite      0;

writeFormat     ascii;

writePrecision  6;

writeCompression off;

timeFormat      general;

timePrecision   6;

runTimeModifiable true;


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

Initial and Boundary Conditions,

Pressure :

/*--------------------------------*- C++ -*----------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | Website:  https://openfoam.org
    \\  /    A nd           | Version:  10
     \\/     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;
    }

    front
    {
        type            wedge;
    }
    
    back
    {
        type            wedge;
    }
    
    axis
    {
        type            empty;
    }
    
    pipewall
    {
        type            zeroGradient;
    }
}

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

Velocity :

/*--------------------------------*- C++ -*----------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | Website:  https://openfoam.org
    \\  /    A nd           | Version:  10
     \\/     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.08927 0 0);
    }

    outlet
    {
        type            zeroGradient;
    }

    front
    {
        type            wedge;
    }
    
    back
    {
        type            wedge;
    }
    
    axis
    {
        type            empty;
    }
    
    pipewall
    {
        type            noSlip;
    }
}
// ************************************************************************* //

Transport Properties File,

/*--------------------------------*- C++ -*----------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | Website:  https://openfoam.org
    \\  /    A nd           | Version:  10
     \\/     M anipulation  |
\*---------------------------------------------------------------------------*/
FoamFile
{
    format      ascii;
    class       dictionary;
    location    "constant";
    object      physicalProperties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

nu              [0 2 -1 0 0 0 0] 8.90e-4;


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

FVSchemes & FVSolution :

Along with the blockMeshDict and controlDict files there are also fvschemes and fvSolution files in the system folder.

fvSchemes script file sets the numerical schemes for terms, such as derivatives in equations, that are calculated during a simulation. 

FVschemes,

/*--------------------------------*- C++ -*----------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | Website:  https://openfoam.org
    \\  /    A nd           | Version:  10
     \\/     M anipulation  |
\*---------------------------------------------------------------------------*/
FoamFile
{
    format      ascii;
    class       dictionary;
    location    "system";
    object      fvSchemes;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

ddtSchemes
{
    default         Euler;
}

gradSchemes
{
    default         Gauss linear;
    grad(p)         Gauss linear;
}

divSchemes
{
    default         none;
    div(phi,U)      Gauss linear;
}

laplacianSchemes
{
    default         Gauss linear orthogonal;
}

interpolationSchemes
{
    default         linear;
}

snGradSchemes
{
    default         orthogonal;
}


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

FVSolutions,

/*--------------------------------*- C++ -*----------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | Website:  https://openfoam.org
    \\  /    A nd           | Version:  10
     \\/     M anipulation  |
\*---------------------------------------------------------------------------*/
FoamFile
{
    format      ascii;
    class       dictionary;
    location    "system";
    object      fvSolution;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

solvers
{
    p
    {
        solver          PCG;
        preconditioner  DIC;
        tolerance       1e-06;
        relTol          0.05;
    }

    pFinal
    {
        $p;
        relTol          0;
    }

    U
    {
        solver          smoothSolver;
        smoother        symGaussSeidel;
        tolerance       1e-05;
        relTol          0;
    }
}

PISO
{
    nCorrectors     2;
    nNonOrthogonalCorrectors 0;
    pRefCell        0;
    pRefValue       0;
}


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

Velocity profile,

Velocity profile and plot,