Aim: Validation studies of Symmetry BC vs Wedge BC in OpenFOAM vs Analytical H.P equation.

Objective

• Validate Hydro-dynamic length with the numerical result
• Validate the fully developed flow velocity profile with its analytical profile
• Validate maximum velocity and pressured drop for fully developed flow
• Post-process Shear stress and validate wall shear stress for fully developed flow
• Highlight the differences in both boundary conditions in terms of their usage and application

Preprocessing

BlockMeshDict preparation and analysis

For the block mesh preparation we enter the following values into the matlab script

• 10, 25 and 45 degrees

5 copies of the original analysis folder placed in the $FOAM_RUN to which the following adjustments are made

• Replace “wedge” with “symmetry” for the symmetric boundary condition analyses in the blockmeshDict file.
• Replace the blockmeshDict file for each angle/BC in their respective system analysis folder
• Replace “wedge” with “symmetry” in the U/P initial condition files for the symmetry boundary condition analyses.

The same schemes and solution controls are applied for all analyses.

Only the relevant files for 45deg/Symmetry are shown here.

/*--------------------------------*- 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.000000 0.000000 0.000000)
	 (-0.003827 0.009239 0.000000)
	 (0.003827 0.009239 0.000000)
	 (0.000000 0.000000 2.310000)
	 (-0.003827 0.009239 2.310000)
	 (0.003827 0.009239 2.310000)
);

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

edges
(
);

boundary
(
	 inlet 
	 { 
	 	 type patch; 
	 	 faces 
	 	 ( 
	 	 (0 1 2 0) 
	 	 ); 
	 } 
	 outlet 
	 { 
	 	 type patch; 
	 	 faces 
	 	 ( 
	 	 (3 4 5 3) 
	 	 ); 
	 } 
	 upperwall 
	 { 
	 	 type wall; 
	 	 faces 
	 	 ( 
	 	 (1 4 5 2) 
	 	 ); 
	 } 
	 frontwall 
	 { 
	 	 type symmetry; 
	 	 faces 
	 	 ( 
	 	 (0 3 4 1) 
	 	 ); 
	 } 
	 backwall 
	 { 
	 	 type symmetry; 
	 	 faces 
	 	 ( 
	 	 (0 2 5 3) 
	 	 ); 
	 } 
);

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

stopAt          endTime;

endTime         20;

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:  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.0;//0.014
    }

    upperwall
    {
        type            zeroGradient;
    }
    
    backwall
    {
        type            symmetry;
    }
    
    frontwall
    {
        type            symmetry;
    }
}

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

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

    outlet
    {
        type            zeroGradient;
    }

    upperwall
    {
        type            noSlip;
    }

    backwall
    {
        type            symmetry;
    }

    frontwall
    {
        type            symmetry;
    }
}

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

Analysis

The following parameters are given/calculated for the purpose of comparing the analytical against openfoam values.

Re=2100 %Reynolds number

rho=1000; %density water

mu=0.89e-3; %dynamic viscosity water at 25degC (mPa)

nu=mu/rho; %kinematic viscosity =mu/rho

R=0.01; %radius of pipe (m)

D=2*R;

L=0.05*Re*D=2.1;

L_e=0.05*Re*D +0.1*0.05*Re*D=2.31; %entrance length for laminar flow +10%

L_dp=L_e-L=0.21; %distance at which we evaluate shear stress and pressure drop

theta=22.5; %enter half of total required angle axis symmetric max 5deg

Vavg=Re*nu/D=0.0934 % approximated from Reynolds equation for RE=2100

vmax=2*Vavg=0.1869%max velocity

dp=8*mu*L_dp*Vavg/(R^2)=1.3973; %pressure drop from Hagen's

dp_k=dp/rho=0.0014 %kinematic pressure difference for comparison against OpenFoam

Results/Conclusion

The following plots shows the 3 variants 10, 25 and 45 deg wedge angle

The following plots below show the max velocity, pressure drop and shear stress calculated for the 6 different variants, 3 wedge angles and 2 boundary condition types.

Velocity

• It shows that the wedge and symmetry BC yield the same result for each angle variant
• Max velocity for the 10deg angle and 25deg matches well with the analytical value (max vel=.1869) while the 45deg angle shows a larger error.
• The 10deg/25deg appears to look more parabolic with respect to the “centreline” at radius of the pipe (0,01m)

Pressure

• It shows that the wedge and symmetry BC yield the same result for each angle variant
• The 10deg wedge appears to compare the best with the analytical value of the kinematic pressure drop. This is evaluated along a distance “after” which the flow is fully developed.
• The 45deg variant deviates the largest from the analytical value (Pressure=0,0014).
• Again the 10deg wedge angle appears to predict the value of shear stress closest to the analytical value

Shear Stress

• It shows that the wedge and symmetry BC yield the same result for each angle variant
• The 10deg wedge predicts the most accurate value of shear stress st a radius ~0.009m relative to the analytical value of 0.033Pa
• The 45deg shows the greatest deviation followed by the 25deg variant.
• It is not clear at the present time what causes the sudden upshoot of shear stress closer to the wall. It is believed that this behaviour is spurious for the following reasons
  • A lack of a finer mesh in the region of the wall to correctly approximate the boundary layer
  • The arc is not modelled.

• All in all, it can be said that there is consistency of results in that the 10deg wedge approximates the analytical values better than the 25/45 deg variants.
• Care must be taken in the choice of what boundary conditions to use for axisymmetric problems.
• This requires an understanding of the underlying assumptions of these modelling shortcuts.
• For example, the openfoam literature states a max of 5 degrees for the patch type wedge.
• In this problem, the finer details of the arc which exist at the ends of the model where not modelled.
• This can have an influence on the predicted flow properties.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%defining the verties for blockMesh input file generation
v1=[0 0 0];
v2=[-R*sind(theta) R*cosd(theta) 0];
v3=[+R*sind(theta) R*cosd(theta) 0];
v4=[0 0 L_e];
v5=[-R*sind(theta) R*cosd(theta) L_e];
v6=[+R*sind(theta) R*cosd(theta) L_e];
% mesh discretisation
nx=20
ny=1; % 1 cell thick according to userguide pg U-39
nz=500;
%blockMeshdict file generation
f=fopen('blockMeshDict','w')

fprintf(f,'%s\n','/*--------------------------------*- C++ -*----------------------------------*\');
fprintf(f,'%s\n','=========                 |');
fprintf(f,'%s\n','\\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox');
fprintf(f,'%s\n','\\    /   O peration     | Website:  https://openfoam.org');
fprintf(f,'%s\n','\\  /    A nd           | Version:  8');
fprintf(f,'%s\n','\\/     M anipulation  |');
fprintf(f,'%s\n','\*---------------------------------------------------------------------------*/');
fprintf(f,'%s\n','FoamFile');
fprintf(f,'%s\n','{');
fprintf(f,'%s\n','version     2.0;');
fprintf(f,'%s\n','format      ascii;');
fprintf(f,'%s\n','class       dictionary;');
fprintf(f,'%s\n','object      blockMeshDict;');
fprintf(f,'%s\n','}');
fprintf(f,'%s\n','// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //');
fprintf(f,'\n');
fprintf(f,'%s\n','convertToMeters 1;');
fprintf(f,'\n');
fprintf(f,'%s\n','vertices');
fprintf(f,'%s\n','(');
fprintf(f,'\t (%f %f %f)\n',v1(1), v1(2), v1(3));
fprintf(f,'\t (%f %f %f)\n',v2(1), v2(2), v2(3));
fprintf(f,'\t (%f %f %f)\n',v3(1), v3(2), v3(3));
fprintf(f,'\t (%f %f %f)\n',v4(1), v4(2), v4(3));
fprintf(f,'\t (%f %f %f)\n',v5(1), v5(2), v5(3));
fprintf(f,'\t (%f %f %f)\n',v6(1), v6(2), v6(3));
fprintf(f,'%s\n',');');
fprintf(f,'\n');
fprintf(f,'%s\n','blocks');
fprintf(f,'%s\n','(');
fprintf(f,'\t hex (0 2 1 0 3 5 4 3)\n');
fprintf(f,'\t (%d %d %d)\n',nx,ny,nz);
fprintf(f,'\t %s\n','simpleGrading (1 1 1)');
fprintf(f,'%s\n',');');
fprintf(f,'\n');
fprintf(f,'%s\n','edges');
fprintf(f,'%s\n','(');
fprintf(f,'%s\n',');');
fprintf(f,'\n');
fprintf(f,'%s\n','boundary');
fprintf(f,'%s\n','(');
fprintf(f,'\t inlet \n');
fprintf(f,'\t { \n');
fprintf(f,'\t \t type patch; \n');
fprintf(f,'\t \t faces \n');
fprintf(f,'\t \t ( \n');
fprintf(f,'\t \t (0 1 2 0) \n');
fprintf(f,'\t \t ); \n');
fprintf(f,'\t } \n');
fprintf(f,'\t outlet \n');
fprintf(f,'\t { \n');
fprintf(f,'\t \t type patch; \n');
fprintf(f,'\t \t faces \n');
fprintf(f,'\t \t ( \n');
fprintf(f,'\t \t (3 4 5 3) \n');
fprintf(f,'\t \t ); \n');
fprintf(f,'\t } \n');
fprintf(f,'\t upperwall \n');
fprintf(f,'\t { \n');
fprintf(f,'\t \t type wall; \n');
fprintf(f,'\t \t faces \n');
fprintf(f,'\t \t ( \n');
fprintf(f,'\t \t (1 4 5 2) \n');
fprintf(f,'\t \t ); \n');
fprintf(f,'\t } \n');
fprintf(f,'\t frontwall \n');
fprintf(f,'\t { \n');
fprintf(f,'\t \t type wedge; \n');
fprintf(f,'\t \t faces \n');
fprintf(f,'\t \t ( \n');
fprintf(f,'\t \t (0 3 4 1) \n');
fprintf(f,'\t \t ); \n');
fprintf(f,'\t } \n');
fprintf(f,'\t backwall \n');
fprintf(f,'\t { \n');
fprintf(f,'\t \t type wedge; \n');
fprintf(f,'\t \t faces \n');
fprintf(f,'\t \t ( \n');
fprintf(f,'\t \t (0 2 5 3) \n');
fprintf(f,'\t \t ); \n');
fprintf(f,'\t } \n');
fprintf(f,'%s\n',');');