Week 11 - Simulation of Flow through a pipe in OpenFoam

AIM

To simulate an axi-symmetric flow by applying the wedge boundary condition through a constant cross-section pipe.

 

OBJECTIVES

 

1. The main objective of this challenge is to simulate an axi-symmetric flow by applying the wedge boundary condition for a laminar flow through a constant cross-section pipe at an inlet Reynold's number of 2100.
2. Carry out transient or steady state simulation to validate the results with the Hagen- Poiseuille's equations for the fully developed flow.
3. To write a code in Matlab to automate the generation of blockMeshDict file, which takes inputs like pipe length, diameter, wedge angle, etc and outputs a fully functional blockMeshDict file. 
4. To obtain the following results from the simulation,
   
   i.) Validate Hydro-dynamic length with the numerical result.
   
   ii.) Validate the fully developed flow velocity profile with its analytical profile.
   
   iii.) Validate maximum velocity and pressured drop for fully developed flow.
   
   iv.) Post-process Shear stress and validate wall shear stress for fully developed flow.

 

Problem Statement

1. create a wedge from a hex block by collapsing its bottom face (by merging both sets of adjacent vertices)
2. Wedge angle should be less than 5^o
3. You can assume the pipe diameter

 

INTRODUCTION

Flow-through pipes or fluid flow is a type of flow within a closed conduit with a certain pressure. We consider the flow of real fluids i.e. the fluids possessing viscosity. Thus, flow of all real fluids is termed as viscous flows. For real fluids, the property viscosity is characterized by the shear stresses or the frictional forces between the fluid layers and fluid to solid surface.

 

THEORY

 

Axisymmetric Flow

An axisymmetric flow is defined as one for which the flow variables, i.e. velocity and pressure, do not vary with the angular coordinate θ.

            Simulating flow through a cylinder can result to be costly and computationally exhaustive if our computers are not powerful enough to run these simulations successfully. So, for this we are going to split our domain and we are going to simulate just a section of the domain (wedge) by making use of the axi-symmetric nature of the flow by applying a symmetry boundary condition.

 

The cross – section of the wedge is as shown below,

 

Laminar Flow

Laminar flow is a flow where the streamlines are smooth and parallel to each other and flow is a highly ordered motion. At low velocities, the fluid tends to flow without lateral mixing, and adjacent layers slide past one another like playing cards. 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. We see this flow from the below

 

Turbulent Flow

 Turbulent flow is a flow where the streamlines are not smooth, the velocity fluctuations and highly disordered motion. It is a flow regime characterized by chaotic property changes. This includes a rapid variation of pressure and flows velocity in space and time.  In contrast to laminar flow, the fluid no longer travels in layers, and mixing across the tube is highly efficient. We see this flow from the below

      If we have been around smokers, we probably noticed that the cigarette smoke rises in a smooth plume for the first few centimetres and then starts fluctuating randomly in all directions as it continues its rise. So, that is one of the example of laminar and turbulent flow.

 

Reynolds number

The Reynolds number represents the ratio of inertial forces to viscous forces and is a convenient parameter for predicting if a flow condition will be laminar or turbulent. It is defined as a characteristic length multiplied by a characteristic velocity and divided by the kinematic viscosity.

`Re_D="Inertial Forces" /"Viscous Forces"`

`=(ρVD)/μ=(VD)/ν`

 

Where,V = Velocity of the Fluid (m/s)

ρ = Density of the Fluid (kg/m³)

D= Characteristic Length (m)

μ = Dynamic Viscosity (Pa·s  or  N·s/m²   or   kg/(m·s))

ν = Kinematic Viscosity (m²/s).

 

The Reynolds number at which the flow becomes turbulent is called the critical Reynolds number, Re_cr. The value of the critical Reynolds number is different for different geometries and flow conditions. For internal flow in a circular pipe, the generally accepted value of the critical Reynolds number is Re_cr=2300. Under most practical conditions, the flow in a circular pipe is laminar for Re ≲ 2300, turbulent for Re ≳ 4000, and transitional in between. That is,

                                                             (Re ≲ 2300) - laminar flow

                                                             (2300 ≲ Re ≲ 4000) - transitional flow

                                                              (Re ≳ 4000) - turbulent flow

 

Fully Developed Flow

 Fully developed flow occurs when the viscous effects due to the shear stress between the fluid particles and pipe wall create a fully developed velocity profile.  The region from the pipe inlet to the point at which the velocity profile is fully developed is called the hydrodynamic entrance region. Flow in the entrance region is called hydrodynamically fully developed flow since this is the region where the velocity profile develops.

 

Entrance Length

 The entry 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.

The hydrodynamic entry length is usually taken to be the distance from the pipe entrance to where the wall shear stress reaches within about 2 percent of the fully developed value. In laminar flow, the nondimensional hydrodynamic entry length is given approximately as,

`L_(h,"laminar" )/D≅0.05Re`

In turbulent flow, the intense mixing during random fluctuations usually overshadows the effects of molecular diffusion. The nondimensional hydrodynamic entry length for turbulent flow is approximated as,

`L_(h,"turbulent" )/D=1.359(Re)^(1/4)`

Since our project is mainly based on the laminar flow, so we are now going to more discussion on the properties of laminar flow through a circular pipe. 

 

Velocity profile for the laminar flow through a circular pipe

The flow is fully developed if the pipe is sufficiently long (relative to the entry length) so that the entrance effects are negligible. In this section, we consider the steady, laminar, incompressible flow of fluid with constant properties in the fully developed region of a straight circular pipe. We obtain the momentum equation by applying a momentum balance to a differential volume element and we obtain the velocity profile by solving it.

This is a free-body diagram of a ring-shaped differential fluid element of radius r, thickness dr, and length dx oriented coaxially with a horizontal pipe in fully developed laminar flow.

The velocity profile is fully developed laminar flow in a pipe is parabolic with a maximum at the centreline and a minimum (zero) at the pipe wall. After simplifying all the mathematical expressions we get the velocity profile expression which is,

`U(r)=-R^2/(4μ) ((dp)/dx)(1-r^2/R^2 )`

`U(r)=2V_(avg) (1-r^2/R^2 )`

This is a convenient form for the velocity profile since V_avg can be determined easily from the flow rate information. The maximum velocity occurs at the centreline and is determined from Eq.(i) by substituting r = 0,

`U_max = 2V_(avg)`

Therefore, the average velocity in fully developed laminar pipe flow is one-half of the maximum velocity.

 

Pressure Drop

A quantity of interest in the analysis of pipe flow is the pressure drop P since it is directly related to the power requirements of the fan or pump to maintain flow. We note that `(dP)/dx="constant"`

and integrating from   where the pressure is   to  where the pressure is  gives

`∆P=P_2-P_1`

`(dP)/dx=(P_2-P_1)/L`

 

Hagen- Poiseuille's Theory

In nonideal fluid dynamics, the Hagen–Poiseuille equation, also known as the Hagen–Poiseuille law, Poiseuille law or Poiseuille equation, is a physical law that gives the pressure drop in an incompressible and Newtonian fluid in laminar flow flowing through a long cylindrical pipe of constant cross section. 

The assumptions of the equation are that the fluid is incompressible and Newtonian; the flow is laminar through a pipe of constant circular cross-section that is substantially longer than its diameter; and there is no acceleration of fluid in the pipe. 

            The pressure drop in an incompressible and Newtonian fluid in laminar flow through a long cylindrical pipe of constant cross – section can be derived as,

`∆P=P_2-P_1=(8μLV_(avg))/R^2`

`=(32μLV_(avg))/D^2`

Where μ is  the fluid viscosity, L is the length of the pipe(>  entry length), V_avg is the average velocity, R is the radius of the pipe and D is the diameter of the pipe. The Hagen–Poiseuille equation describes the relationship between pressure, fluidic resistance, and flow rate, analogous to voltage, resistance, and current, respectively, in Ohm’s law for electrical circuits (V=RI).

The ratio of length to radius of a pipe should be greater than one forty-eighth of the Reynolds number for the Hagen–Poiseuille law to be valid. If the pipe is too short, the Hagen–Poiseuille equation may result in unphysically high flow rates.

 A few examples in which the Hagen–Poiseuille law is practically applied is the air flow in lung alveoli, or the flow through a drinking straw or through a hypodermic needle. It is also applied for the process of blood flow through a narrow constriction.

 

METHODOLOGY FOLLOWED

For the numerical simulation in the OpenFOAM software, an axisymmetric wedge has been considered as the domain for simulation and respective suitable boundary conditions have been applied in order to simulate the flow as accurately as possible with minimum assumptions.

 

Calculations for Analytical Solution

 

Assumption

• Flow is steady, incompressible (density is constant) and laminar.
• Constant fluid properties like Viscosity and density.
• Working temperature of water 25^oC
• Reynolds Number = 2100 (Given)
• Diameter             = 2cm = 0.02m (Assumed)

 

Analytical Solution

1. Diameter of the pipe `D=0.02m → R=D/2=0.01m`
2. Density of working fluid (water) at 25^oC  `ρ=997 (Kg)/m^3`
3. Wedge angle `θ=4^o`
4. Dynamic Viscosity of water at 25^oC `ν=μ/ρ=(8.90×10^(-4))/997=8.92×10^(-7) m^2/sec`
5. Hydrodynamic Entrance Length  `L_h=0.06⋅R_e⋅D= =0.06×2100×0.02=2.52m`
6. Total Length of the Pipe `L=2.8m`
7. Average Velocity `V_avg=(μ⋅R_e)/(ρ⋅D)==(8.90×10^(-4) × 2100)/(997×0.02)=0.0937 m/s`
8. Maximum Velocity `Δp = (32⋅μ⋅V_(avg)⋅L)/D^2 ==(32×8.90×10^(-4)×0.0937×2.8)/(0.02)^2 =18.680 Pa`
9. Velocity Profile `u(r)=2⋅V_(avg)⋅(1–r^2/R^2 )`
10. Pressure Drop `Δp = (32⋅μ⋅V_(avg)⋅L)/D^2 =(32×8.90×10^(-4)×0.0937×2.8)/(0.02)^2 =18.680 Pa`
11. Kinematic Pressure Drop  `Δp=(Δp)/ρ=18.680/997=0.0187 Pa`

 

SIMULATION PROCEDURE

 

Problem Setup

            We first initialize the problem by navigating to the tutorials directory of OpenFOAM by using the “tut” command in the terminal window. From the tutorials, we have to figure out the suitable case file for our given problem, which is an axi-symmetric laminar flow through a constant cross section pipe. We use the “ls” command to list the contents (files or directories) present inside the directory.

 

As it is an incompressible flow, we navigate into the incompressible directory. We use the “cd” command to navigate and change the working directory using the terminal window.

 

The next step involves navigating into the icoFoam directory, as the type of solver for this problem is icoFoam. The icoFoam command solves the incompressible laminar Navier-Stokes equations. The code is inherently transient, requiring an initial condition (such as zero velocity) and boundary conditions.

 

Here, we use the cavity utility to create our required domain. We copy this directory into the $FOAM_RUN file, which is a default file OpenFOAM creates for us. We use the “cp -r” command to recursively create a copy of the contents of the file or directory into the $FOAM_RUN folder under the file name “Flow Through Pipe”, for personal use. We create a copy of the files as OpenFOAM does not let us modify or make changes to its inbuilt tutorial files.

 

We then navigate to the copied tutorial file in the $FOAM_RUN directory to view the contents inside the folder.

 

PROGRAM CODES AND FILES :

First, we are going to write a MatLab code to automate the generation of the blockMeshDict file, which takes our respective inputs and outputs a fully functional blockMeshDict file.

 

Program Code :

% SIMULATION OF FLOW THROUGH A PIPE
% Code to automate generation of blockMeshDict file
clear all
close all
clc


% INPUTS
% Assumptions
% Diameter of Pipe (m)
D = 0.02; 

% Radius of Pipe (m)
R = D/2;

% Reynold's Number
Re = 2100;

% Wedge Angle (degrees)
theta = 4;

% Total Length of the Pipe (m)
L = 2.8;



% Assuming working fluid is water at 25° C
% Density of water (kg/m³)
rho = 997;

% Dynamic Viscosity of Water (Pa.s)
mu = 8.90 * 10^-4;

% Kinematic Viscosity of Water (m^2/s)
nu = mu/rho;


% Hydrodynamic Entrance Length (m)
L_h = 0.06 * Re * D;



% HAGEN - POISEUILLE'S FLOW EQUATIONS
% Average Velocity (m/s)
V_avg = (mu * Re)/(rho * D);

% Maximum Velocity (m/s)
V_max = 2 * V_avg;

% Pressure Drop (Pa)
delta_p = (32*mu*V_avg*L)/(D^2);

% Kinematic Pressure Drop (Pa)
kinematic_delta_p = delta_p/rho;



% DEFINING THE GEOMETRY
% 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)];


% BLOCKMESHDICT FILE
% Creating the file
% Header
h1 = '/*--------------------------------*- C++ -*----------------------------------*\';
h2 = '  =========                 |';
h3 = '  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox';
h4 = '   \\    /   O peration     | Website:  https://openfoam.org';
h5 = '    \\  /    A nd           | Version:  10';
h6 = '     \\/     M anipulation  |';
h7 = '\*---------------------------------------------------------------------------*/';
h8 = 'FoamFile';
h9 = '{';
h10 = '    format      ascii;';
h11 = '    class       dictionary;';
h12 = '    object      blockMeshDict;';
h13 = '}';
h14 = '// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //';


% Creating BlockMesh File
f1 = fopen('blockMeshDict.txt','w');

% Header
fprintf(f1,'%s\n',h1);
fprintf(f1,'%s\n',h2);
fprintf(f1,'%s\n',h3);
fprintf(f1,'%s\n',h4);
fprintf(f1,'%s\n',h5);
fprintf(f1,'%s\n',h6);
fprintf(f1,'%s\n',h7);
fprintf(f1,'%s\n',h8);
fprintf(f1,'%s\n',h9);
fprintf(f1,'%s\n',h10);
fprintf(f1,'%s\n',h11);
fprintf(f1,'%s\n',h12);
fprintf(f1,'%s\n',h13);
fprintf(f1,'%s\n\n',h14);

% Conversion Factor
fprintf(f1,'%s\n\n','convertToMeters 1;');

% Vertices
fprintf(f1,'%s\n','vertices');
fprintf(f1,'%s\n','(');
fprintf(f1,'\t(%d %d %d)\n',V0);
fprintf(f1,'\t(%d %f %f)\n',V1);
fprintf(f1,'\t(%d %f %f)\n',V2);
fprintf(f1,'\t(%.2f %d %d)\n',V3);
fprintf(f1,'\t(%.2f %f %f)\n',V4);
fprintf(f1,'\t(%.2f %f %f)\n',V5);
fprintf(f1,'%s\n\n',');');

% Blocks
fprintf(f1,'%s\n','blocks');
fprintf(f1,'%s\n','(');
fprintf(f1,'\thex (0 3 4 1 0 3 5 2) (400 30 1) simpleGrading (1 0.1 1)\n');
fprintf(f1,'%s\n\n',');');

% Edges
fprintf(f1,'%s\n','edges');
fprintf(f1,'%s\n','(');
fprintf(f1,'\tarc 2 1 (%d %.3f %d)\n',0,R,0);
fprintf(f1,'\tarc 5 4 (%.3f %.3f %d)\n',L,R,0);
fprintf(f1,'%s\n\n',');');

% Boundary
fprintf(f1,'%s\n','boundary');
fprintf(f1,'%s\n','(');

% Axis
fprintf(f1,'\taxis\n');
fprintf(f1,'\t{\n');
fprintf(f1,'\t\ttype empty;\n');
fprintf(f1,'\t\tfaces\n');
fprintf(f1,'\t\t(\n');
fprintf(f1,'\t\t\t(0 3 3 0)\n');
fprintf(f1,'\t\t);\n');
fprintf(f1,'\t}\n\n');

% Inlet
fprintf(f1,'\tinlet\n');
fprintf(f1,'\t{\n');
fprintf(f1,'\t\ttype patch;\n');
fprintf(f1,'\t\tfaces\n');
fprintf(f1,'\t\t(\n');
fprintf(f1,'\t\t\t(0 1 2 0)\n');
fprintf(f1,'\t\t);\n');
fprintf(f1,'\t}\n\n');

% Outlet
fprintf(f1,'\toutlet\n');
fprintf(f1,'\t{\n');
fprintf(f1,'\t\ttype patch;\n');
fprintf(f1,'\t\tfaces\n');
fprintf(f1,'\t\t(\n');
fprintf(f1,'\t\t\t(3 5 4 3)\n');
fprintf(f1,'\t\t);\n');
fprintf(f1,'\t}\n\n');

% Front Face
fprintf(f1,'\tfront\n');
fprintf(f1,'\t{\n');
fprintf(f1,'\t\ttype wedge;\n');
fprintf(f1,'\t\tfaces\n');
fprintf(f1,'\t\t(\n');
fprintf(f1,'\t\t\t(0 3 4 1)\n');
fprintf(f1,'\t\t);\n');
fprintf(f1,'\t}\n\n');

% Back Face
fprintf(f1,'\tback\n');
fprintf(f1,'\t{\n');
fprintf(f1,'\t\ttype wedge;\n');
fprintf(f1,'\t\tfaces\n');
fprintf(f1,'\t\t(\n');
fprintf(f1,'\t\t\t(0 2 5 3)\n');
fprintf(f1,'\t\t);\n');
fprintf(f1,'\t}\n\n');

% Top Wall
fprintf(f1,'\tpipeWall\n');
fprintf(f1,'\t{\n');
fprintf(f1,'\t\ttype wall;\n');
fprintf(f1,'\t\tfaces\n');
fprintf(f1,'\t\t(\n');
fprintf(f1,'\t\t\t(1 4 5 2)\n');
fprintf(f1,'\t\t);\n');
fprintf(f1,'\t}\n');

fprintf(f1,'%s\n\n',');');
fprintf(f1,'%s\n','mergePatchPairs');
fprintf(f1,'%s\n','(');
fprintf(f1,'%s\n\n',');');
fprintf(f1,'%s\n','// ************************************************************************* //');
From the above code, we obtain our respective blockMeshDict file.
	The vertices of the wedge are formed by calculating the values of  R⋅Cos⁡(θ/2) and  R⋅Sin⁡(θ/2) with respect to the axis of the wedge.

	We use the “fprintf” command in the code to write the contents of the file and thereby obtaining the blockMeshDict text file as shown below.

	We use the hexahedral mesh to create the wedge geometry with the help of the calculated vertices. We use a grading factor of 0.1 for this simulation such that we get finer meshing along the walls.

	The number of cells on x , y and z - axes are defined as 400 , 30 and 1 respectively. 

	We create an arc representing the radius of the wedge and also define the input, output, axis, wall, front and back faces of the wedge.

 

blockMeshDict File :

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

convertToMeters 1;

vertices
(
	(0 0 0)
	(0 0.009994 -0.000349)
	(0 0.009994 0.000349)
	(2.80 0 0)
	(2.80 0.009994 -0.000349)
	(2.80 0.009994 0.000349)
);

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

edges
(
	arc 2 1 (0 0.010 0)
	arc 5 4 (2.800 0.010 0)
);

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

	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)
		);
	}
);

mergePatchPairs
(
);

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

The blockMeshDict file has been written automatically by the MatLab program as shown above, which we use to create the geometry of the wedge and run the simulation for the same.

 

controlDict File :

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

application     icoFoam;

startFrom       startTime;

startTime       0;

stopAt          endTime;

endTime         40;

deltaT          0.01;

writeControl    timeStep;

writeInterval   200;

purgeWrite      0;

writeFormat     ascii;

writePrecision  6;

writeCompression off;

timeFormat      general;

timePrecision   6;

runTimeModifiable true;


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

We run the simulation for a period of 40 seconds in timesteps of 0.01 seconds.

 

transportProperties File :

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

nu              [0 2 -1 0 0 0 0] 8.9268e-07;


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

We input the calculated value of Kinematic Viscosity (ν) of water at 25^oC  as `8.92×10^(-7) m^2/sec`

 

Initial Velocity (U) File :

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

dimensions      [0 1 -1 0 0 0 0];

internalField   uniform (0.0937 0 0);

boundaryField
{
    axis
    {
        type            empty;
    }
    
    inlet
    {
        type            fixedValue;
        value           uniform (0.0937 0 0);
    }
    
    outlet
    {
        type            zeroGradient;
    }

    front
    {
        type            wedge;
    }

    back
    {
        type            wedge;
    }
    
    pipeWall
    {
        type            noSlip;
    }
}

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

We define the velocity condition at the inlet of the pipe V_avg as  0.0937m/s

 

Initial Pressure (p) File :

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

dimensions      [0 2 -2 0 0 0 0];

internalField   uniform 0;

boundaryField
{
    axis
    {
        type            empty;
    }
    
    inlet
    {
        type            zeroGradient;
    }
    
    outlet
    {
        type            zeroGradient;
        value           uniform 0.0187;
    }

    front
    {
        type            wedge;
    }

    back
    {
        type            wedge;
    }
    
    pipeWall
    {
        type            zeroGradient;
    }
}

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

We define the kinematic pressure drop (kinematic  Δp) condition at the outlet of the pipe as 0.01870 Pa

 

 

RESULTS :

By running the simulation with respect to the above mentioned data, we get the outputs as follows,

 

WEDGE
We create a cross section of the pipe, i.e; a wedge geometry, with a wedge angle of 4^o as shown below.
 

Geometry :

 

Meshing (axi - symmetric) :

 

2. VELOCITY PROFILES
The change in velocity of the fluid flowing inside the wedge section is as shown below.

At inlet :

 

At outlet :

 

3. PRESSURE PROFILES
The pressure distribution across the cross sectional pipe is represented as follows.

 

At inlet :

 

At outlet :

 

 

 

 

PLOTS

From the above simulation, we obtain the respective plots as follows.

 

1. VELOCITY PROFILES
Here, we are going to study the velocity profiles at certain points or regions along the pipe.

 

Velocity Profile at x = 0.1m from Inlet :

 

Velocity Profile at x = 0.5m from Inlet :

 

Velocity Profile at x = 1.4m from Inlet :

 

Velocity Profile at x = 2.8m from Inlet (Outlet) :

 

2. PRESSURE DROP
            The pressure drop across the cross section of the wedge is as shown below.

 

3. Shear stress plot

 

Now we check the validation of our results. For that we do numerical analysis in MATLAB

 

clear all;
close all;
clc;

Re = 2100;          % Reynolds number
Rho = 1000;         % Density of water
mu = 0.00089;       % Dynamic viscosity of water
nu = mu/Rho;        % Kinematic viscosity of water
D = 0.02;         % Diameter of the pipe
R = D/2;            % Radius of the pipe
ent_L = 0.05*Re*D;  % Entry length
theta = 4;          % Wedge angle in degree(It should be less than 5)

%Boundery condition
V_avg = (nu*Re)/D;             %Average velocity
V_max = 2*V_avg;               % Maximum velocity
L = ent_L + 0.25;              %Length of the pipe
Dell_P = (32*mu*L*V_avg)/D^2;  % Pressure drop
kin_P = Dell_P/Rho;            % Kinematic pressure drop

r = linspace(0,R,100);
for i = 1:length(r)
    U(i) = 2*V_avg*(1-(r(i)^2)/(R^2));
    tau(i) = 2*mu*V_max*(abs(r(i)))/R^2;

end
%shear stress
figure(1)
plot(r,tau,"r*","LineWidth",1)
xlabel("Radius","FontWeight","bold")
ylabel("Shear stress","FontWeight","bold")
title("Shear stress Analysis")

figure(2)
plot(r,U,"g*","LineWidth",1);
xlabel("Radius","FontWeight","bold")
ylabel("Velocity","FontWeight","bold")
title("Velocity Profile Validation Analysis")

 

Numerical results are

This is Velocity plot for fully developed flow:-

 

This is the Shear stress plot for fully developed flow:-

 

Kinematic Pressure value for our project

 

CONCLUSION

1. From the above study, we have successfully simulated an axi – symmetric flow by applying the wedge boundary condition across a constant cross section pipe, for an inlet Reynold’s Number of 2100 and a wedge angle of 4^o
2. The MatLab code to automate the generation of the blockMeshDict file has been successfully written, which takes the provided inputs and outputs a fully functional blockMeshDict file as shown above.
3. From the simulation run, we obtain our respective wedge geometry and velocity profiles as shown above. To study the velocity profiles, we carry out a test to study the flow characteristics at different regions of the pipe. We study the velocity profile plots at 4 regions of the pipe, at a distance of 1m , 0.5m , 1.4m and 2.8m (outlet) from the inlet of the pipe.
4. From the plots obtained, we can observe that the fluid velocity at the walls of the pipe remains zero, due to the no-slip boundary condition we have provided. We can also observe that the velocity profiles slowly develop and become more and more sharper (parabolic) as the distance from the inlet increases.
5. Due to the no – slip conditions at the walls of the pipe, the adjacent layers of the fluid tend to slow down as the flow progresses, due to friction created between these layers. This leads to an increase in velocity of the fluid at the centre of the pipe. As a result, the velocity profiles keep developing as the distance from the inlet increases and the fluid reaches maximum velocity at the outlet of the pipe, which is referred to as the fully developed flow region.
6. From the pressure plots obtained from the simulation, we can observe that the fluid pressure is the highest at the inlet of the pipe and it gradually decreases throughout the length of the pipe. We can see that the pressure is the least at the outlet and thereby, we can say that the pressure drop obtained from this simulation obeys the Hagen – Poiseuille’s Theory.
7. Shear stresses are developed in the fluid as a result of relative motion of the fluid layers adjacent to the walls of the pipe, which are provided with a no – slip boundary condition, leading to the velocity at the walls of the pipe being zero. From the shear stress plot obtained, we can observe that the shear stresses are maximum at the walls of the pipe and it slowly decreases as it reaches the centre of the pipe, as the relative friction between the layers slowly decrease.