OpenFOAM

Pipe flow-Symmetry BC using OpenFOAM

Objective– To use the icoFOAM solver and simulate the flow through a pipe by creating a wedge and providing symmetry boundary condition.

To write a program in Matlab that can generate the computational mesh automatically for any wedge angle and grading schemes.

Compare the results for symmetry,wedge and HP equation.

Dimensions of Pipe:

Diameter(D)=10mm

To calculate length of pipe entry length formula is used

Le=0.06⋅Re⋅D

Re=2100 

   =0.06*2100*10=1260mm=12.6m

So the length is taken as

    L=1500mm=1.5m

Calculation for inlet velocity:

Re=VDν

ν=8.9⋅10-7

V=2100⋅8.9⋅10-70.01=0.1869m/s

blockMeshDict file

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

vertices
(
    (0 0 0)
    (0 4.98097 -0.435779)
    (0 4.98097 0.435779)
    (1500 0 0)
    (1500 4.98097 -0.435779)
    (1500 4.98097 0.435779)
);
blocks
(
    hex (0 1 2 0 3 4 5 3) (50 1 400) simpleGrading (0.1 1 10)
);
edges
(
    arc 1 2 (0 5 0)
    arc 4 5 (1500 5 0)
);
boundary
(
inlet
   {
     type patch;
     faces
(
      (0 2 1 0)
);
   }
outlet
{
     type patch;
     faces
   (
     (3 4 5 3)
);
}
upperWall
{
       type wall;
      faces
   (
      (1 2 5 4)
);
}
lowerWall
{
       type empty;
      faces
(
      (0 3 3 0)
);
}

front
{
       type symmetry;
      faces
(
      (0 3 5 2)
);
}

   Back
{
       type symmetry;
      faces
(
       (0 1 4 3)
);
}
);
mergePatchPairs
(
);
// ************************************************************************* //

Automated blockmesh generation:

close all;
clear all;
clc;
theta=10;
L=1500;
dia=10;
r=dia/2;
y=r*cosd(theta/2);
z=r*sind(theta/2);


f1=fopen('blockmeshdict.txt','wt');

fprintf(f1,'/*--------------------------------*- C++ -*----------------------------------*\n');
 fprintf (f1,'=========                |n')
 fprintf (f1,'\      /  F ield         | OpenFOAM: The Open Source CFD Toolboxn');
  fprintf (f1,'\    /   O peration     | Website:  https://openfoam.orgn');
  fprintf  (f1,' \  /    A nd           | Version:  7n');
   fprintf  (f1,'  \/     M anipulation  |n');
fprintf(f1, '')
fprintf(f1,'*---------------------------------------------------------------------------*/n');
fprintf(f1,'FoamFilen{n')
   fprintf(f1,' version     2.0;n');
   fprintf(f1,' format      ascii;n');
    fprintf(f1,' class       dictionary;n');
    fprintf(f1,' object      blockMeshDict;n}n');

fprintf(f1,'// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //n');
fprintf(f1,'convertToMeters 0.001;nn');
fprintf(f1,'verticesn');
fprintf(f1,'(n');
fprintf(f1,'    (0 0 0)n')
fprintf(f1,'    (0 %g -%g)n',y,z)
fprintf(f1,'    (0 %g %g)n',y,z)
fprintf(f1,'    (%g 0 0)n',L)
fprintf(f1,'    (%g %g -%g)n',L,y,z)
fprintf(f1,'    (%g %g %g)n',L,y,z)
fprintf(f1,');n');
fprintf(f1,'blocksn');
fprintf(f1,'(n');
fprintf(f1,'    hex (0 1 2 0 3 4 5 3) (50 1 400) simpleGrading (0.1 1 10)n')
fprintf(f1,');n');
fprintf(f1,'edgesn');
fprintf(f1,'(n');
fprintf(f1,'    arc 1 2 (0 %g 0)n',r)
fprintf(f1,'    arc 4 5 (%g %g 0)n',L,r)
fprintf(f1,');n');
fprintf(f1,'boundaryn');
fprintf(f1,'(n');
fprintf(f1,'inletn   {n');
        fprintf(f1,'     type patch;n');
        fprintf(f1,'     facesn(n');
            fprintf(f1,'      (0 2 1 0)n);n   }n');

       
   fprintf(f1,'outletn{n');
        fprintf(f1,'     type patch;n');
        fprintf(f1,'     facesn   (n');
             fprintf(f1,'     (3 4 5 3)n);n}n');
   fprintf(f1,'upperWalln{n');
      fprintf(f1,'       type wall;n');
       fprintf(f1,'      facesn   (n');
            fprintf(f1,'      (1 2 5 4)n);n}n');
     fprintf(f1,'lowerWalln{n');
      fprintf(f1,'       type empty;n');
       fprintf(f1,'      facesn(n');
            fprintf(f1,'      (0 3 3 0)n);n}nn');
fprintf(f1,'frontn{n');
      fprintf(f1,'       type symmetry;n');
       fprintf(f1,'      facesn(n');
            fprintf(f1,'      (0 3 5 2)n);n}nn');
            fprintf(f1,'   Backn{n');
      fprintf(f1,'       type symmetry;n');
       fprintf(f1,'      facesn(n');
            fprintf(f1,'       (0 1 4 3)n);n}n);n');

fprintf(f1,'mergePatchPairsn(n);n');



fprintf(f1,'// ************************************************************************* //');

controlDict file

In this file time step and time of simulation is provided.The time-step was taken as 1e-3 which gives a courant number below 1.

The simulation was run for 9 seconds.

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

application     icoFoam;

startFrom       startTime;

startTime       0;

stopAt          endTime;

endTime         9;

deltaT          0.001;

writeControl    timeStep;

writeInterval   20;

purgeWrite      0;

writeFormat     ascii;

writePrecision  6;

writeCompression off;

timeFormat      general;

timePrecision   6;
runTimeModifiable true;


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

Transport Properties:

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

nu              [0 2 -1 0 0 0 0] 8.9e-7;


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

fvShemes

The discretization scheme is provided here.

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


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

fvSolution

The iterative solver is included here.

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


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

Initial and Boundary conditions:

Initial and boundary conditions are given in p and U file located in the 0 folder.

The initial condition is indicated by the name Internal Field and is set to 0,which means pressure is 0 at the start.

Boundary condition is given in the boundary field

Intlet velocity is given as 0.1869 m/s and pressure is extrapolated from inner nodes

Oultet velocity is extrapolated from inner nodes and pressure is set to 0

No-slip boundary condition at top wall and empty at bottom walls.

Front and back wedge is set to symmetry boundary conditions

p file 

 

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

    front
    {
        type            symmetry;
    }
    Back
    {
        type            symmetry;
    }
}

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

U file

/*--------------------------------*- C++ -*----------------------------------*
  =========                 |
  \      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \    /   O peration     | Website:  https://openfoam.org
    \  /    A nd           | Version:  7
     \/     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
{
    inlet
    {
        type            fixedValue;
        value           uniform (0.1869 0 0);
    }
    outlet
    {
        type            zeroGradient;
    }
    upperWall
    {
        type            noSlip;
    }
    lowerWall
    {
        type            empty;
    }

    front
    {
        type            symmetry;
    }
    Back
    {
        type            symmetry;
    }
}

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

Case-1 : Wedge angle=10 Degrees

Mesh specification

  1. Number of cells along the x direction (longer dimension) = 400
  2. Number of cells along the y direction = 50

 

Output:

Velocity contour

Velocity magnitude at 1.5m from inlet (Fully Developed)

Pressure Drop across pipe:

Case-2 : Wedge angle=25 Degrees

Mesh specification

  1. Number of cells along the x direction (longer dimension) = 400
  2. Number of cells along the y direction = 50

Velocity contour

Velocity magnitude at 1.5m from inlet (Fully Developed)

Pressure Drop across pipe:

Case-3 : Wedge angle=45 Degrees

Mesh specification

  1. Number of cells along the x direction (longer dimension) = 400
  2. Number of cells along the y direction = 50

Velocity contour

Velocity magnitude at 1.5m from inlet (Fully Developed)

Pressure Drop across pipe:

 Haigen poissule equation:

Analytical solution,

According to Haigen poissule equation

Δp=32μLVD2

           =32⋅8.9⋅10-4⋅0.1869⋅1.50.012=79.84 pa

           =0.0798 m2s2

And

Vmax=2*Vavg=2*0.1869=0.3738 m/s

Simulation results from Wedge boundary condition:

`Deltap=0.1 m^2/s^2`

And

Vmax=0.375 m/s

Results in Table:

Conclusion:

So it can be seen that the maximum velocity is same in each case but the pressure drop is also close to analytical solution.So it can be concluded that both BCs give similar results.

Leave a Reply

Your email address will not be published. Required fields are marked *