Week 12 - Validation studies of Symmetry BC vs Wedge BC in OpenFOAM vs Analytical H.P equation

1. AIM : To simulate the Laminar flow in the Pipe using Symmetry BC and Wedge BC and validate the result with analytical solution of Hagen Poiseuille's equation.

2. SETUP AND SIMULATION

The orginal geometry was reduced and used for simulating as shown in Fig 2.1. In this way the computation time can be reduced and finer mesh can be also used for better result.

Fig 2.1 : Computational Geomtery

Reynolds number was kept below 2100 for a laminar flow as given below.

1. Velocity (V) = 0.1m/s
2. Diameter (D) = 0.015m
3. Kinematic Viscosity (ν">νν) = 8.9x10-7 m2s">m2sm2s

∴Re=VDν">∴Re=VD/ν∴Re=VDν

=(0.1)(0.015)(8.9)⋅10-7">=(0.1)(0.015)(8.9)⋅10−7=(0.1)(0.015)(8.9)⋅10-7

=1685.39">=1685.39=1685.39 .......................... [Laminar Flow]

• Entry length for fully developed Laminar flow = 0.05⋅Re⋅D">0.05⋅Re⋅D0.05⋅Re⋅D =0.05⋅(1685.39)⋅(0.015)">=0.05⋅(1685.39)⋅(0.015)=0.05⋅(1685.39)⋅(0.015) =1.26m">=1.26m=1.26m

Hence, the piple length was kept more than 1.26m to see the effect of fully developed laminar flow.

Pipe Length (L) = 1.5m

2.1 Automated blockmesh generation file using Matlab

Fig 2.2 : Boundary types

The program to generate a blockmesh file is given below. Any value of input varibale as length, diameter and angle of the wedge would automatically generate a blockmesh file which can be used in openfoam. Same code would be used for three different angles as given below.

Angles to test:

• 10o">10o10o
• 25o">25o25o
• 45o">45o45o

For these angles, two simulation were run. One with wedge boundary condition and another with symmetry boundary condtion as shown in Fig 2.2.

Matlab code;

clear all
close all
clc

l = 1500;                %Length of pipe in mm
d = 15;                 %Diameter of pipe in mm
r = d/2;                %radius
theta = 10;             %angle
y = r*cosd(theta/2);    %Y-cordinate
z = r*sind(theta/2);    %X-cordinate

f1 =  fopen('blockMeshDict.txt','w');
fprintf(f1,'FoamFilen{')
fprintf(f1,'n    version     2.0;')
fprintf(f1,'n    format     ascii;')
fprintf(f1,'n    class     dictionary;')
fprintf(f1,'n    object     blockMeshDict;n}')

fprintf(f1,'nnconvertToMeters 0.001;')

fprintf(f1,'nnverticesn(')
fprintf(f1,'n    (0 0 0)')
fprintf(f1,'n    (0 %g %g)',y,z)
fprintf(f1,'n    (0 %g -%g)',y,z)
fprintf(f1,'n    (%g 0 0)',l)
fprintf(f1,'n    (%g %g %g)',l,y,z)
fprintf(f1,'n    (%g %g -%g)',l,y,z)
fprintf(f1,'n);')

fprintf(f1,'nnblocksn(')
fprintf(f1,'n    hex (0 3 5 2 0 3 4 1) (70 20 5) simpleGrading (1 0.1 1)n);')

fprintf(f1,'nnedgesn(')
fprintf(f1,'n    arc 1 2(0 %g 0)',r)
fprintf(f1,'n    arc 4 5(%g %g 0)',l,r)
fprintf(f1,'n);')

fprintf(f1,'nnboundaryn(')
fprintf(f1,'n    inletn    {')
fprintf(f1,'n       type patch;')
fprintf(f1,'n       facesn       (')
fprintf(f1,'n         (0 1 2 0)n       );n    }')
fprintf(f1,'n    outletn    {')
fprintf(f1,'n       type patch;')
fprintf(f1,'n       facesn       (')
fprintf(f1,'n         (3 5 4 3)n       );n    }')
fprintf(f1,'n    upperwalln    {')
fprintf(f1,'n       type wall;')
fprintf(f1,'n       facesn       (')
fprintf(f1,'n         (1 4 5 2)n       );n    }')
fprintf(f1,'n    wedge_1n    {')
fprintf(f1,'n       type wedge;')
fprintf(f1,'n       facesn       (')
fprintf(f1,'n         (0 2 5 3)n       );n    }')
fprintf(f1,'n    wedge_2n    {')
fprintf(f1,'n       type wedge;')
fprintf(f1,'n       facesn       (')
fprintf(f1,'n         (0 3 4 1)n       );n    }')
fprintf(f1,'n);')


fprintf(f1,'nnmergePatchPairs')
fprintf(f1,'n(n);')
fclose(f1)

Output blockMeshDict file

FoamFile
{
    version     2.0;
    format     ascii;
    class     dictionary;
    object     blockMeshDict;
}

convertToMeters 0.001;

vertices
(
    (0 0 0)
    (0 7.49886 0.130893)
    (0 7.49886 -0.130893)
    (1500 0 0)
    (1500 7.49886 0.130893)
    (1500 7.49886 -0.130893)
);

blocks
(
    hex (0 3 5 2 0 3 4 1) (300 40 1) simpleGrading(1 0.1 1)
);

edges
(
    arc 1 2(0 7.49886 0)
    arc 4 5(1500 7.49886 0)
);

boundary
(
    inlet
    {
       type patch;
       faces
       (
         (0 1 2 0)
       );
    }
    outlet
    {
       type patch;
       faces
       (
         (3 5 4 3)
       );
    }
    upperwall
    {
       type wall;
       faces
       (
         (1 4 5 2)
       );
    }
    wedge_1
    {
       type wedge;
       faces
       (
         (0 2 5 3)
       );
    }
    wedge_2
    {
       type wedge;
       faces
       (
         (0 3 4 1)
       );
    }
);

mergePatchPairs
(
);

a) 10o">10o10o b)25o">25o25o c)45o">45o45o

Fig 2.3: Computational Mesh with Grading Factor 0.1 near the wall

2.2 Solver

The steady state result was needed hence time step size was kept larger targeting to obtian the end steady state result faster while maintaining the stability of the solution. The solver was choosen for simulation an incompressible laminar flow.

1. Solver : icoFoam
2. timestep : 0.01
3. Simulation start time : 0
4. Simulation start time : 45

3. VERIFICATION AND VALIDATION

The simualtion results were validated using Hagen Poiseuille equation as shown below.

3.1 Velocity

Fig 3.1 shows a velocity profile of a fully developed laminar flow. As the pipe flow is axis symmetry, the parabolic velocity profile can be seen in half.

Analytical Result,

Vmax=2Vavg">Vmax=2VavgVmax=2Vavg

= 2x0.1

= 0.2m/s

Fig 3.1 : Velocity profile along radial direction of fully developed Laminar flow

3.2 Pressure

Fig 3.2 shows the kinematic pressure along the length of the pipe. The results are as expected as high pressure at inlet and lower pressure at the outlet. Also the pressure drop is linear at the fully developed flow region after x = 1m.

Analytical Result,

Δp=32μLVD2">Δp=32μLVD2Δp=32μLVD2

=(32)⋅(8.9)⋅(10-4)⋅(1.5)⋅(0.1)(0.015)2">=(32)⋅(8.9)⋅(10−4)⋅(1.5)⋅(0.1)(0.015)2=(32)⋅(8.9)⋅(10-4)⋅(1.5)⋅(0.1)(0.015)2

=18.98pa">=18.98pa=18.98pa

= 0.0189 m2/s2

Fig 3.3 : Kinematic Pressure along the Pipe

3.3 Wall Shear Stress

As shown in Fig 3.4, the shear stress is varying linearly along the radial direction which is expected. For laminar flow, the Hagen Poiseuille and Darchy-Weisbach gives same result from where we get the friction factor,

f=64Re">f=64Ref=64Re

= 64/ 1685.39

= 0.038

Analytical Solution

Wall shear stress (τw">τwτw) = 12(f4ρVavg2)">12(f4ρV2avg)12(f4ρVavg2)

= 12⋅0.0384⋅1000⋅(0.1)2">12⋅0.0384⋅1000⋅(0.1)212⋅0.0384⋅1000⋅(0.1)2

= 0.0475pa

Fig 3.4 : Shear Stress along radial direction of fully developed Laminar flow

It is clear that inspite of the wedge angle, the value for all the variables are close to the analytical solution when symmetry boundary conditions were used. The comparision of numerical value and analytical value is shown in Table 1.

Table 1 : Value for Symmetry Boundary Condition

	10o">10o10o	25o">25o25o	45o">45o45o	AnalyticalResult
Maxmimum Velocity(m/s)	0.1961	0.1966	0.1966	0.2
Pressure drop(m2/s2)	0.0252	0.0253	0.0254	0.0189
Wall Shear Stress(pa)	0.0459	0.0459	0.0459	0.0475

When wedge boundary condition were used, the numerical result were not close to the analytical result as shown in Table 2. The values for all the variables were approaching toward the analytical value when the angle was reduced.

Table 2 : Value for Wedge Boundary Condition

	10o">10o10o	25o">25o25o	45o">45o45o	AnalyticalResult
Maxmimum Velocity(m/s)	0.24	0.311	0.319	0.2
Pressure drop(m2/s2)	0.0346	0.0532	0.056	0.0189
Wall Shear Stress(pa)	0.0581	0.0814	0.0836	0.0475

4. CONCLUSION

Wedge BC is more suitable for geometry with wedge angle less than 5 degree. Hence it is used for axis symmetric cases where the flow field variables doesn't change with the angular axis. For the larger wedge angle symmetry boundary condition can be used.

An object that follows a cylindrical cooridnate system and is axis symmetric, both wedge BC and symmetry BC can be used. But in case of object like square, which is not axis symmetric but has symmetric plane, only symmetry BC can be used.