BlogFor BusinessTalk to a counsellor
CoursesIIT ProgramsPlacementsLog in
BlogFor Business
Talk to a counsellor
Challenge

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

SSridhar Kumar20 Jun 2023
CFDMATLAB

Validation studies of Symmetry BC vs Wedge BC in OpenFOAM vs Analytical Hagen

AIM:

To Simulate the flow through pipe with wedge BC and symmetry BC and to validate the results in MATLAB for (10 ,25 & 45 degrees).

OBJECTIVE:

  1. Validate Hydro-dynamic length with the numerical result
  2. Validate the fully developed flow velocity profile with its analytical profile
  3. Validate maximum velocity and pressured drop for fully developed flow
  4. Post-process Shear stress and validate wall shear stress for fully developed flow
  5. Highligh the differences in both boundary conditions in terms of their usage and aplication

THEORY:

Wedge BC:

Wedge boundary conditions are used in mathematical models to solve boundary problems in a wedge-shaped domain. The wedge boundary condition is used in computational fluid dynamics to define an axisymmetric situation, for example, a cylinder. The model and flow must be axisymmetric along a central line such that all physical variables of the flow have the same value and distribution at a given radius for all angles.In OpenFOAM, wedge-type boundary conditions are used for purely axisymmetric calculations. The geometry is specified as a wedge of small angle (<5°) and 1 cell thick running along the plane of symmetry, straddling one of the coordinate planes. Wedge symmetry is also used as a boundary condition in some computational fluid dynamics simulations.

Symmetry BC:

Symmetry boundary condition implies that the normal velocities are set to zero and the scalar flux across the boundary is zero. In other words, there is no flow across the symmetry plane and no heat transfer or diffusion of scalars. This means that the values of properties just adjacent to the solution domain are taken as values at the nearest node just inside the domain.

INPUTS:

Diameter (D) 0.02 m
Radius (R) D/2
Reynolds Number(Re) 2100
Dynamic Viscosity (`mu`) 8.9e-4 (Ns/m2)
Hydrodynamic length (l_e) 0.05*Re*D m = 2.1m
Total length of the pipe (l) l_e+0.3 m = 2.4 m
Density (`rho`) 997 (kg/m3)
Velocity Average (Vavg) `(Re*mu)/(rho*D)` =9.373e-2 (m/s2)
maximum velocity (Umax) vavg*2 = 1.875e-4 (m/s2)
pressure drop (`DeltaP`) `Deltap = (8*mu*V_(avg)*l)/R^2` =16.0168
kiematic viscosity (`nu`) `mu/rho` = 8.9268e-7

solving Procedure:

In order to solve this case in openfoam it is necessary to create the blockMeshDict file. blockMesh Dict file is created same as the week 11 challenge and the same has been modified for different boundary conditions and the wedge angle (theta).

Software OpenFoam
Flow Incompressible laminar flow
Time Transient
Solver icoFoam
Result visualization Paraiew

Below image shows the representation of the vertices of the wedge;

MATLAB Code to create blockMeshDict file: (for Wedge Boundary Condition)

clear all
close all
clc

% Input Parameters

D = 0.02;   %Diameter of the pipe in meters
R = D/2; %Radius of the wedge
Re = 2100;  %Reynolds Number
l_e = 0.05*Re*D;  %Hydrodynamic entry Length
theta = input('Wedge Angle in deg = ');    %Wedge Angle
rho = 997; %kg/m3
mu = 0.00089; %dynamic viscosity
nu = mu/rho; %kinematic iscosity

% Boundary Conditions
v_avg = ((Re*mu)/(rho*D));
v_max = v_avg*2;
l = l_e+0.3; %Total Length
delp = (8*mu*v_avg*l)/(R^2);

scale = 1; %scaling the values to meters

%  Defining the vertices
v0 = [0 0 0];
v1 = [0 R*cosd(theta/2) -R*sind(theta/2)];
v2 = [0 R*cosd(theta/2) R*sind(theta/2)];
v3 = [l 0 0];
v4 = [l R*cosd(theta/2) -R*sind(theta/2)];
v5 = [l R*cosd(theta/2) R*sind(theta/2)];

% Defining blocks and grading factor
cellx = 50; %cell count along x axis
celly = 1;  %cell count along y axis
cellz = 400; %cell count along z axis
gradx = 1; %grading factor along x axis
grady = 1; %grading factor along y axis
gradz = 1; %grading factor along z axis

% Arc Inputs
arc1 = [1,2];
arc2 = [4,5];

% creating the blockMesh structure
f1 = fopen('blockMeshDict.txt','w');
fprintf(f1,'%s\n','/*--------------------------------*- C++ -*----------------------------------*\');
fprintf(f1,'%s\n','========                |');
fprintf(f1,'%s\n','\      /  Field          | OpenFOAM: The Open Source CFD Toolbox');
fprintf(f1,'%s\n',' \    /   O peration     | Website:  https://openfoam.org');
fprintf(f1,'%s\n','  \  /    A nd           | Version:  6');
fprintf(f1,'%s\n','   \/     M anipulation  |');
fprintf(f1,'%s\n',' *---------------------------------------------------------------------------*/');
fprintf(f1,'%s\n','FoamFile');
fprintf(f1,'%s\n','{');
fprintf(f1,'%s\n','    version     2.0;');
fprintf(f1,'%s\n','    format      ascii;');
fprintf(f1,'%s\n','    class       dictionary;');
fprintf(f1,'%s\n','    object      blockMeshDict;');
fprintf(f1,'%s\n','}');
fprintf(f1,'%s\n','// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //');
fprintf(f1,'%s %d %s\n','ConvertToMeters',scale,';');

% Defining Vertices
fprintf(f1,'\n');
fprintf(f1,'%s\n','vertices');
fprintf(f1,'(\n');
fprintf(f1,'     (%d %d %d)\n',v0(1),v0(2),v0(3));
fprintf(f1,'     (%d %d %d)\n',v1(1),v1(2),v1(3));
fprintf(f1,'     (%d %d %d)\n',v2(1),v2(2),v2(3));
fprintf(f1,'     (%d %d %d)\n',v3(1),v3(2),v3(3));
fprintf(f1,'     (%d %d %d)\n',v4(1),v4(2),v4(3));
fprintf(f1,'     (%d %d %d)\n',v5(1),v5(2),v5(3));
fprintf(f1,'%s\n',');');


% Defining Blocks
fprintf(f1,'\n');
fprintf(f1,'%s\n','blocks');
fprintf(f1,'(\n');
fprintf(f1,'hex (0 1 2 0 3 4 5 3)(%d %d %d) %s (%d %d %d)\n',cellx,celly,cellz,'simpleGrading',gradx,grady,gradz);
fprintf(f1,'%s\n',');');

% Defining Edges
fprintf(f1,'\n');
fprintf(f1,'%s\n','edges');
fprintf(f1,'%s\n','(');
fprintf(f1,'   arc %d %d (%d %d %d)\n',arc1(1),arc1(2),0,R,0);
fprintf(f1,'   arc %d %d (%d %d %d)\n',arc2(1),arc2(2),l,R,0);
fprintf(f1,'%s\n',');');

% Defining Boundaries
fprintf(f1,'\n');
fprintf(f1,'%s\n','boundary');
fprintf(f1,'(\n');
fprintf(f1,'   inlet\n');
fprintf(f1,'   {\n');
fprintf(f1,'      type patch;\n');
fprintf(f1,'      faces\n');
fprintf(f1,'      (\n');
fprintf(f1,'          (0 2 1 0)\n');
fprintf(f1,'      );\n');
fprintf(f1,'   }\n');

fprintf(f1,'\n');
fprintf(f1,'   front\n');
fprintf(f1,'   {\n');
fprintf(f1,'      type wedge;\n');
fprintf(f1,'      faces\n');
fprintf(f1,'      (\n');
fprintf(f1,'          (0 3 5 2)\n');
fprintf(f1,'      );\n');
fprintf(f1,'   }\n');

fprintf(f1,'\n');
fprintf(f1,'   back\n');
fprintf(f1,'   {\n');
fprintf(f1,'      type wedge;\n');
fprintf(f1,'      faces\n');
fprintf(f1,'      (\n');
fprintf(f1,'          (0 1 4 3)\n');
fprintf(f1,'      );\n');
fprintf(f1,'   }\n');

fprintf(f1,'\n');
fprintf(f1,'   top\n');
fprintf(f1,'   {\n');
fprintf(f1,'     type wall;\n');
fprintf(f1,'      faces\n');
fprintf(f1,'      (\n');
fprintf(f1,'          (1 2 5 4)\n');
fprintf(f1,'      );\n');
fprintf(f1,'   }\n');

fprintf(f1,'\n');
fprintf(f1,'   axis\n');
fprintf(f1,'   {\n');
fprintf(f1,'     type empty;\n');
fprintf(f1,'      faces\n');
fprintf(f1,'      (\n');
fprintf(f1,'          (0 3 3 0)\n');
fprintf(f1,'      );\n');
fprintf(f1,'   }\n');

fprintf(f1,'\n');
fprintf(f1,'   outlet\n');
fprintf(f1,'   {\n');
fprintf(f1,'      type patch;\n');
fprintf(f1,'      faces\n');
fprintf(f1,'      (\n');
fprintf(f1,'          (3 4 5 3)\n');
fprintf(f1,'      );\n');
fprintf(f1,'   }\n');
fprintf(f1,');\n');
fprintf(f1,'\n');

% Writing the mergepach pairs
fprintf(f1,'%s\n','mergePatchPairs');
fprintf(f1,'(\n');
fprintf(f1,');\n');
fprintf(f1,'%s\n','// ************************************************************************* //');
fclose(f1);

Output of the blockMeshDict file:

For Wedge 10°

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

vertices
(
     (0 0 0)
     (0 9.961947e-03 -8.715574e-04)
     (0 9.961947e-03 8.715574e-04)
     (2.400000e+00 0 0)
     (2.400000e+00 9.961947e-03 -8.715574e-04)
     (2.400000e+00 9.961947e-03 8.715574e-04)
);

blocks
(
hex (0 1 2 0 3 4 5 3)(50 1 400) simpleGrading (1 1 1)
);

edges
(
   arc 1 2 (0 1.000000e-02 0)
   arc 4 5 (2.400000e+00 1.000000e-02 0)
);

boundary
(
   inlet
   {
      type patch;
      faces
      (
          (0 2 1 0)
      );
   }

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

   back
   {
      type wedge;
      faces
      (
          (0 1 4 3)
      );
   }

   top
   {
     type wall;
      faces
      (
          (1 2 5 4)
      );
   }

   axis
   {
     type empty;
      faces
      (
          (0 3 3 0)
      );
   }

   outlet
   {
      type patch;
      faces
      (
          (3 4 5 3)
      );
   }
);

mergePatchPairs
(
);
// ************************************************************************* //

For Wedge 25°;

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

vertices
(
     (0 0 0)
     (0 9.762960e-03 -2.164396e-03)
     (0 9.762960e-03 2.164396e-03)
     (2.400000e+00 0 0)
     (2.400000e+00 9.762960e-03 -2.164396e-03)
     (2.400000e+00 9.762960e-03 2.164396e-03)
);

blocks
(
hex (0 1 2 0 3 4 5 3)(50 1 400) simpleGrading (1 1 1)
);

edges
(
   arc 1 2 (0 1.000000e-02 0)
   arc 4 5 (2.400000e+00 1.000000e-02 0)
);

boundary
(
   inlet
   {
      type patch;
      faces
      (
          (0 2 1 0)
      );
   }

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

   back
   {
      type wedge;
      faces
      (
          (0 1 4 3)
      );
   }

   top
   {
     type wall;
      faces
      (
          (1 2 5 4)
      );
   }

   axis
   {
     type empty;
      faces
      (
          (0 3 3 0)
      );
   }

   outlet
   {
      type patch;
      faces
      (
          (3 4 5 3)
      );
   }
);

mergePatchPairs
(
);
// ************************************************************************* //

For Wedge 45°;

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

vertices
(
     (0 0 0)
     (0 9.238795e-03 -3.826834e-03)
     (0 9.238795e-03 3.826834e-03)
     (2.400000e+00 0 0)
     (2.400000e+00 9.238795e-03 -3.826834e-03)
     (2.400000e+00 9.238795e-03 3.826834e-03)
);

blocks
(
hex (0 1 2 0 3 4 5 3)(50 1 400) simpleGrading (1 1 1)
);

edges
(
   arc 1 2 (0 1.000000e-02 0)
   arc 4 5 (2.400000e+00 1.000000e-02 0)
);

boundary
(
   inlet
   {
      type patch;
      faces
      (
          (0 2 1 0)
      );
   }

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

   back
   {
      type wedge;
      faces
      (
          (0 1 4 3)
      );
   }

   top
   {
     type wall;
      faces
      (
          (1 2 5 4)
      );
   }

   axis
   {
     type empty;
      faces
      (
          (0 3 3 0)
      );
   }

   outlet
   {
      type patch;
      faces
      (
          (3 4 5 3)
      );
   }
);

mergePatchPairs
(
);
// ************************************************************************* //

MATLAB Code to create blockMeshDict file: (for Symmetry Boundary Condition)

clear all
close all
clc

% Input Parameters

D = 0.02;   %Diameter of the pipe in meters
R = D/2; %Radius of the wedge
Re = 2100;  %Reynolds Number
l_e = 0.05*Re*D;  %Hydrodynamic entry Length
theta = input('Wedge Angle in deg = ');    %Wedge Angle
rho = 997; %kg/m3
mu = 0.00089; %dynamic viscosity
nu = mu/rho; %kinematic iscosity

% Boundary Conditions
v_avg = ((Re*mu)/(rho*D));
v_max = v_avg*2;
l = l_e+0.3; %Total Length
delp = (8*mu*v_avg*l)/(R^2);

scale = 1; %scaling the values to meters

%  Defining the vertices
v0 = [0 0 0];
v1 = [0 R*cosd(theta/2) -R*sind(theta/2)];
v2 = [0 R*cosd(theta/2) R*sind(theta/2)];
v3 = [l 0 0];
v4 = [l R*cosd(theta/2) -R*sind(theta/2)];
v5 = [l R*cosd(theta/2) R*sind(theta/2)];

% Defining blocks and grading factor
cellx = 50;      %cell count along x axis
celly = 1;       %cell count along y axis
cellz = 400;     %cell count along z axis
gradx = 1;       %grading factor along x axis
grady = 1;       %grading factor along y axis
gradz = 1;       %grading factor along z axis

% Arc Inputs
arc1 = [1,2];
arc2 = [4,5];

% creating the blockMesh structure
f1 = fopen('blockMeshDict.txt','w');
fprintf(f1,'%s\n','/*--------------------------------*- C++ -*----------------------------------*\');
fprintf(f1,'%s\n','========                |');
fprintf(f1,'%s\n','\      /  Field          | OpenFOAM: The Open Source CFD Toolbox');
fprintf(f1,'%s\n',' \    /   O peration     | Website:  https://openfoam.org');
fprintf(f1,'%s\n','  \  /    A nd           | Version:  6');
fprintf(f1,'%s\n','   \/     M anipulation  |');
fprintf(f1,'%s\n',' *---------------------------------------------------------------------------*/');
fprintf(f1,'%s\n','FoamFile');
fprintf(f1,'%s\n','{');
fprintf(f1,'%s\n','    version     2.0;');
fprintf(f1,'%s\n','    format      ascii;');
fprintf(f1,'%s\n','    class       dictionary;');
fprintf(f1,'%s\n','    object      blockMeshDict;');
fprintf(f1,'%s\n','}');
fprintf(f1,'%s\n','// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //');
fprintf(f1,'%s %d %s\n','ConvertToMeters',scale,';');

% Defining Vertices
fprintf(f1,'\n');
fprintf(f1,'%s\n','vertices');
fprintf(f1,'(\n');
fprintf(f1,'     (%d %d %d)\n',v0(1),v0(2),v0(3));
fprintf(f1,'     (%d %d %d)\n',v1(1),v1(2),v1(3));
fprintf(f1,'     (%d %d %d)\n',v2(1),v2(2),v2(3));
fprintf(f1,'     (%d %d %d)\n',v3(1),v3(2),v3(3));
fprintf(f1,'     (%d %d %d)\n',v4(1),v4(2),v4(3));
fprintf(f1,'     (%d %d %d)\n',v5(1),v5(2),v5(3));
fprintf(f1,'%s\n',');');


% Defining Blocks
fprintf(f1,'\n');
fprintf(f1,'%s\n','blocks');
fprintf(f1,'(\n');
fprintf(f1,'hex (0 1 2 0 3 4 5 3)(%d %d %d) %s (%d %d %d)\n',cellx,celly,cellz,'simpleGrading',gradx,grady,gradz);
fprintf(f1,'%s\n',');');

% Defining Edges
fprintf(f1,'\n');
fprintf(f1,'%s\n','edges');
fprintf(f1,'%s\n','(');
fprintf(f1,'   arc %d %d (%d %d %d)\n',arc1(1),arc1(2),0,R,0);
fprintf(f1,'   arc %d %d (%d %d %d)\n',arc2(1),arc2(2),l,R,0);
fprintf(f1,'%s\n',');');

% Defining Boundaries
fprintf(f1,'\n');
fprintf(f1,'%s\n','boundary');
fprintf(f1,'(\n');
fprintf(f1,'   inlet\n');
fprintf(f1,'   {\n');
fprintf(f1,'      type patch;\n');
fprintf(f1,'      faces\n');
fprintf(f1,'      (\n');
fprintf(f1,'          (0 2 1 0)\n');
fprintf(f1,'      );\n');
fprintf(f1,'   }\n');

fprintf(f1,'\n');
fprintf(f1,'   front\n');
fprintf(f1,'   {\n');
fprintf(f1,'      type symmetry;\n');
fprintf(f1,'      faces\n');
fprintf(f1,'      (\n');
fprintf(f1,'          (0 3 5 2)\n');
fprintf(f1,'      );\n');
fprintf(f1,'   }\n');

fprintf(f1,'\n');
fprintf(f1,'   back\n');
fprintf(f1,'   {\n');
fprintf(f1,'      type symmetry;\n');
fprintf(f1,'      faces\n');
fprintf(f1,'      (\n');
fprintf(f1,'          (0 1 4 3)\n');
fprintf(f1,'      );\n');
fprintf(f1,'   }\n');

fprintf(f1,'\n');
fprintf(f1,'   top\n');
fprintf(f1,'   {\n');
fprintf(f1,'     type wall;\n');
fprintf(f1,'      faces\n');
fprintf(f1,'      (\n');
fprintf(f1,'          (1 2 5 4)\n');
fprintf(f1,'      );\n');
fprintf(f1,'   }\n');

fprintf(f1,'\n');
fprintf(f1,'   axis\n');
fprintf(f1,'   {\n');
fprintf(f1,'     type empty;\n');
fprintf(f1,'      faces\n');
fprintf(f1,'      (\n');
fprintf(f1,'          (0 3 3 0)\n');
fprintf(f1,'      );\n');
fprintf(f1,'   }\n');

fprintf(f1,'\n');
fprintf(f1,'   outlet\n');
fprintf(f1,'   {\n');
fprintf(f1,'      type patch;\n');
fprintf(f1,'      faces\n');
fprintf(f1,'      (\n');
fprintf(f1,'          (3 4 5 3)\n');
fprintf(f1,'      );\n');
fprintf(f1,'   }\n');
fprintf(f1,');\n');
fprintf(f1,'\n');

% Writing the mergepach pairs
fprintf(f1,'%s\n','mergePatchPairs');
fprintf(f1,'(\n');
fprintf(f1,');\n');
fprintf(f1,'%s\n','// ************************************************************************* //');
fclose(f1);

Output of the blockMeshDict file:

For Symmetry BC & theta - 10°

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

vertices
(
     (0 0 0)
     (0 9.961947e-03 -8.715574e-04)
     (0 9.961947e-03 8.715574e-04)
     (2.400000e+00 0 0)
     (2.400000e+00 9.961947e-03 -8.715574e-04)
     (2.400000e+00 9.961947e-03 8.715574e-04)
);

blocks
(
hex (0 1 2 0 3 4 5 3)(50 1 400) simpleGrading (1 1 1)
);

edges
(
   arc 1 2 (0 1.000000e-02 0)
   arc 4 5 (2.400000e+00 1.000000e-02 0)
);

boundary
(
   inlet
   {
      type patch;
      faces
      (
          (0 2 1 0)
      );
   }

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

   back
   {
      type symmetry;
      faces
      (
          (0 1 4 3)
      );
   }

   top
   {
     type wall;
      faces
      (
          (1 2 5 4)
      );
   }

   axis
   {
     type empty;
      faces
      (
          (0 3 3 0)
      );
   }

   outlet
   {
      type patch;
      faces
      (
          (3 4 5 3)
      );
   }
);

mergePatchPairs
(
);
// ************************************************************************* //

For Symmetry BC & theta - 25°;

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

vertices
(
     (0 0 0)
     (0 9.762960e-03 -2.164396e-03)
     (0 9.762960e-03 2.164396e-03)
     (2.400000e+00 0 0)
     (2.400000e+00 9.762960e-03 -2.164396e-03)
     (2.400000e+00 9.762960e-03 2.164396e-03)
);

blocks
(
hex (0 1 2 0 3 4 5 3)(50 1 400) simpleGrading (1 1 1)
);

edges
(
   arc 1 2 (0 1.000000e-02 0)
   arc 4 5 (2.400000e+00 1.000000e-02 0)
);

boundary
(
   inlet
   {
      type patch;
      faces
      (
          (0 2 1 0)
      );
   }

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

   back
   {
      type symmetry;
      faces
      (
          (0 1 4 3)
      );
   }

   top
   {
     type wall;
      faces
      (
          (1 2 5 4)
      );
   }

   axis
   {
     type empty;
      faces
      (
          (0 3 3 0)
      );
   }

   outlet
   {
      type patch;
      faces
      (
          (3 4 5 3)
      );
   }
);

mergePatchPairs
(
);
// ************************************************************************* //

For Symmetry BC & theta - 45°;

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

vertices
(
     (0 0 0)
     (0 9.238795e-03 -3.826834e-03)
     (0 9.238795e-03 3.826834e-03)
     (2.400000e+00 0 0)
     (2.400000e+00 9.238795e-03 -3.826834e-03)
     (2.400000e+00 9.238795e-03 3.826834e-03)
);

blocks
(
hex (0 1 2 0 3 4 5 3)(50 1 400) simpleGrading (1 1 1)
);

edges
(
   arc 1 2 (0 1.000000e-02 0)
   arc 4 5 (2.400000e+00 1.000000e-02 0)
);

boundary
(
   inlet
   {
      type patch;
      faces
      (
          (0 2 1 0)
      );
   }

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

   back
   {
      type symmetry;
      faces
      (
          (0 1 4 3)
      );
   }

   top
   {
     type wall;
      faces
      (
          (1 2 5 4)
      );
   }

   axis
   {
     type empty;
      faces
      (
          (0 3 3 0)
      );
   }

   outlet
   {
      type patch;
      faces
      (
          (3 4 5 3)
      );
   }
);

mergePatchPairs
(
);
// ************************************************************************* //

ControlDict:

Here we define the timestep and no of iterations to solve. ControlDict file can be found in the system subfolder of the respective case files.

This is same for all the cases.

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

deltaT          0.01;

writeControl    timeStep;

writeInterval   10;

purgeWrite      0;

writeFormat     ascii;

writePrecision  6;

writeCompression off;

timeFormat      general;

timePrecision   6;

runTimeModifiable true;


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

Physical (Transport Properties) properties:

Here we define the physical properties of the file which is kinematic viscosity in this case and the physical properties file can be found in the constant folder in the respective case files. this is the same for all the cases.

/*--------------------------------*- 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.926e-7;


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

Defining the Pressure;

For Pressure distribution around the wedge profile,

Here we define the named selection of the geometry and the behaviour of pressure acting on the surfaces.

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

    front
    {
        type            wedge;
    }

    back
    {
        type            wedge;
    }

    top
    {
        type            zeroGradient;
    }

    axis
    {
        type            empty;
    {

    outlet
    {
        type            fixedValue;
        value           uniform 0;
    }
}

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

Defining the Velocity;

For Velocity distribution around the wedge profile,

Here we define the named selection of the geometry and the behaviour of Velocity acting on the surfaces.

/*--------------------------------*- 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.09373 0 0);

boundaryField
{
    inlet
    {
        type            fixedValue;
        value           uniform (0.09373 0 0);
    }

    front
    {
        type            wedge;
    }

    back
    {
        type            wedge;
    }

    top
    {
        type            noSlip;
    }

    axis
    {

        type            empty;
    }

    outlet
    {
        type            zeroGradient;
    }
}

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

After setting the necessary data for the simulation, we need to build the mesh for all the cases;

To build the mesh use - blockMesh

To check the mesh - checkMesh

After checking the mesh next step is to run the solver use - icoFoam

COMPARING AND ANALYSING THE RESULTS IN THE MATLAB ENVIRONMENT

VALIDATING WITH MATLAB CODE:

clear all 
close all
clc

% Required variables
% Reynold's number
D = 0.02;
R = D/2;
R_e = 2100;
L_e = 0.05*R_e*D;
L = 0.3+L_e;
theta = 3;
mu = 0.00089;  % Dynamic viscosity (N.s/m^2)
rho = 997; % (kg/m3)
nu = mu/rho;

% Analytical Method
r = linspace(0,R,1000);   % Radial Distance
r_t = r';
% Hagen Poiseuille's solution
v_avg = (R_e*mu)/(rho*D); %Average velocity
u_max = 2*v_avg %maximum Velocity
v_r = u_max*(1-(r.^2/R^2));
v_rt = v_r';

% Hagen poiseille's Equation for Pressure Drop
P_len = linspace(0,L,1000)
P_len1 = P_len';
for i = 1:length(P_len1)
    delta_p(i) = (8*mu*P_len1(i)*v_avg)/R^2;
    k_pre(i) = delta_p(i)/rho;  %kinematicc pressure
end

k_pre = fliplr(k_pre);

% Analytical shear stress calculation
for i = 1:length(r)
    tau_an(i) = (2*mu*u_max*r(i))/R^2;
    tau_an1=tau_an';
end

% Extracting data from openfoam
% For wedge10
f1 = fopen('v10.csv','r');
fgetl(f1);
for i = 1:1000
    B = fgetl(f1);
    C = strsplit(B,',');
    v10(i) = C(1);
    L10(i) = C(2);
end
v10 = str2double(v10);
L10 = str2double(L10);

% Reading and Extracting shear stress 
f2 = fopen('t10.csv','r');
fgetl(f2);
for i = 1:1000
    B = fgetl(f2);
    C = strsplit(B,',');
    t10(i) = C(1);
    Lw10(i) = C(2);
end
t10 = str2double(t10);
Lw10 = str2double(Lw10);


% For Wedge25
% Extracting data from openfoam
f1 = fopen('v25.csv','r');
fgetl(f1);
for i = 1:1000
    B = fgetl(f1);
    C = strsplit(B,',');
    v25(i) = C(1);
    L25(i) = C(2);
end
v25 = str2double(v25);
L25 = str2double(L25);

% Reading and Extracting shear stress 
f2 = fopen('t25.csv','r');
fgetl(f2);
for i = 1:1000
    B = fgetl(f2);
    C = strsplit(B,',');
    t25(i) = C(1);
    Lw25(i) = C(2);
end
t25 = str2double(t25);
Lw25 = str2double(Lw25);


% For Wedge45
% Extracting data from openfoam
f1 = fopen('v45.csv','r');
fgetl(f1);
for i = 1:1000
    B = fgetl(f1);
    C = strsplit(B,',');
    v45(i) = C(1);
    L45(i) = C(2);
end
v45 = str2double(v45);
L45 = str2double(L45);

% Reading and Extracting shear stress 
f2 = fopen('t45.csv','r');
fgetl(f2);
for i = 1:1000
    B = fgetl(f2);
    C = strsplit(B,',');
    t45(i) = C(1);
    Lw45(i) = C(2);
end
t45 = str2double(t45);
Lw45 = str2double(Lw45);


% For Symmetry10
% Extracting data from openfoam
f1 = fopen('vs10.csv','r');
fgetl(f1);
for i = 1:1000
    B = fgetl(f1);
    C = strsplit(B,',');
    vs10(i) = C(1);
    Ls10(i) = C(2);
end
vs10 = str2double(vs10);
Ls10 = str2double(Ls10);

% Reading and Extracting shear stress 
f2 = fopen('ts10.csv','r');
fgetl(f2);
for i = 1:1000
    B = fgetl(f2);
    C = strsplit(B,',');
    ts10(i) = C(1);
    Lts10(i) = C(2);
end
ts10 = str2double(ts10);
Lts10 = str2double(Lts10);

% For Symmetry25
% Extracting data from openfoam
f1 = fopen('vs25.csv','r');
fgetl(f1);
for i = 1:1000
    B = fgetl(f1);
    C = strsplit(B,',');
    vs25(i) = C(1);
    Ls25(i) = C(2);
end
vs25 = str2double(vs25);
Ls25 = str2double(Ls25);

% Reading and Extracting shear stress 
f2 = fopen('ts25.csv','r');
fgetl(f2);
for i = 1:1000
    B = fgetl(f2);
    C = strsplit(B,',');
    ts25(i) = C(1);
    Lts25(i) = C(2);
end
ts25 = str2double(ts25);
Lts25 = str2double(Lts25);


% For Symmetry45
% Extracting data from openfoam
f1 = fopen('vs45.csv','r');
fgetl(f1);
for i = 1:1000
    B = fgetl(f1);
    C = strsplit(B,',');
    vs45(i) = C(1);
    Ls45(i) = C(2);
end
vs45 = str2double(vs45);
Ls45 = str2double(Ls45);

% Reading and Extracting shear stress 
f2 = fopen('ts45.csv','r');
fgetl(f2);
for i = 1:1000
    B = fgetl(f2);
    C = strsplit(B,',');
    ts45(i) = C(1);
    Lts45(i) = C(2);
end
ts45 = str2double(ts45);
Lts45 = str2double(Lts45);

% Plotting and comparing the results

figure(1)
plot(r_t,v_rt,'r','linewidth',1); % Anaytical Data
hold on;
plot(L10,v10,'--','linewidth',1); % Numerical Data
plot(L25,v25,'--','linewidth',1);
plot(L45,v45,'--','linewidth',1);
plot(Ls10,vs10,':','linewidth',1);
plot(Ls25,vs25,':','linewidth',1);
plot(Ls45,vs45,':','linewidth',1);
xlabel('Location of Radial Dist. from the center to wall of pipe (in m)', ...
    'FontWeight','bold','FontSize',12);
ylabel('Velocity (in m/s)','FontWeight','bold','FontSize',12);
title('Velocity Profile Comparison of Analytical and Numerical values','FontSize',14);
legend('Analytical Solution','For 10 Deg Wedge','For 25 Deg Wedge','For 45 Deg Wedge', ...
   'For 10 Deg Symmetry','For 25 Deg Symmetry','For 45 Deg Symmetry');
grid on


figure(2)
plot(r_t,v_rt,'r','linewidth',1); % Anaytical Data
hold on;
plot(Lw10,t10,'--','linewidth',1); % Numerical Data
plot(Lw25,t25,'--','linewidth',1);
plot(Lw45,t45,'--','linewidth',1);
plot(Lts10,ts10,':','linewidth',1);
plot(Lts25,ts25,':','linewidth',1);
plot(Lts45,ts45,':','linewidth',1);
xlabel('Location of Radial Dist. from the center to wall of pipe (in m)', ...
    'FontWeight','bold','FontSize',12);
ylabel('Shear stress (in Pa)','FontWeight','bold','FontSize',12);
title('Shear Stress Profile Comparison of Analytical and Numerical values','FontSize',14);
legend('Analytical Solution','For 10 Deg Wedge','For 25 Deg Wedge','For 45 Deg Wedge', ...
   'For 10 Deg Symmetry','For 25 Deg Symmetry','For 45 Deg Symmetry','Location','northwest');
grid on

OUTPUT FOR VELOCITY DISTRIBUTION:

OUTPUT FOR SHEAR STRESS DISTRIBUTION:

CONCLUSION:

Hence the Simulation of the flow through pipe with wedge BC and symmetry BC is carried out successfully and the results are validated in MATLAB for (10 ,25 & 45 degrees). BlockMesh creation from MATLAB is done to simulate the flow in Openfoam and the results are imported in txt or csv and the results are plotted (analytical vs Numerical).From the results shown above, we can see that the values with the Wedge Boundary condition is slightly accurate than the Symmetry boundary condition.However, we can notice that the accuracy difference between the Wedge and Symmetry boundary condition gets lower as the angle increases, therefore it is recommended to use Wedge boundary condition with angle less than5°.

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.