Week 11 - Simulation of Flow through a pipe in OpenFoam
Aim: Simulate an axi-symmetric flow by applying the wedge boundary condition.
Given data: The flow is laminar and flows through a constant cross-section pipe at an inlet Reynold's number of 2100.
Theory:
Pipe flow is a type of liquid flow within a closed conduit, such as a pipe or tube. The other type of flow within a conduit is open channel flow.
Pipe flow does not have a free surface which is found in open-channel flow. Pipe flow, being confined within closed conduit, does not exert direct atmospheric pressure, but does exert hydraulic pressure on the conduit.
The entrance 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 variation of wall shear stress in the flow direction for flow in a pipe from the entrance region into the fully developed region.

Calculations:


Coordinates of the wedge:
P0=(0 0 0)
P1=(0 rcos`alpha/2` rsin`alpha/2`)
P2=(L rcos`alpha/2` -rsin`alpha/2`)
P3=(L 0 0)
P4=( L rcos`alpha/2` rsin`alpha/2`)
p5=(L rcos`alpha/2` -rsin`alpha/2`)
Matlab Code to generate the BlockMesh file:
First of all define the geometry vertices and the define the mesh and then the boundaries.
For help in defining the vertices and boundaries refer to openfoam.rrg
clear all
close all
clc
% Flow variables
rho = 1000; % Density
mu = 0.00089; % Dynamic Viscosity
nu = mu/rho;
re = 2100; % Reynolds Number
d = 0.02; % Dia of Pipe
alpha = 4; % wedge angle <5
a = pi/4*d^2;
p = pi*d; % wetted perimeter
d_h = 4*a/p; % Hydraulic Dia of pipe
r_h = d_h/2; %Hydralic radius of pipe
v_avg = (re*mu)/(rho*d_h); % Avg velocity
v_max = 2*v_avg; %Fully developed V max
tau_max = 2*mu*v_max/r_h;
l_e = 0.06*re*d_h;
l =l_e+0.18; % Meters length of Pipe
%Daracy Friction factory --> for Laminar Circular pipe flow
f = 65/re;
%Haigen Poiseuille Eqn
P = (32*mu*(l-l_e)*v_avg)/(d^2);
g = 9.81;
h_loss = P/(rho*g); %head loss
P_kinematicPress = P/rho;
theta = 3; %3 degree wedge angle
r=r_h;
% Coordinates of the wedge
p0 = [0 0 0];
p1 = [0 r*cosd(alpha/2) r*sind(alpha/2)];
p2 = [0 r*cosd(alpha/2) -r*sind(alpha/2)];
p3 = [l 0 0];
p4 = [l r*cosd(alpha/2) r*sind(alpha/2)];
p5 = [l r*cosd(alpha/2) -r*sind(alpha/2)];
line1 ='/*--------------------------------*- C++ -*----------------------------------*\';
line2 = ' ========= |';
line3 = ' \\ / F ield | OpenFOAM: The Open Source CFD Toolbox';
line4 = ' \\ / O peration | Website: https://openfoam.org';
line5 = ' \\ / A nd | Version: 10';
line6 = ' \\/ M anipulation |';
line7 = '\*---------------------------------------------------------------------------*/';
line8 = 'FoamFile';
line9 = '{';
line10 = 'format ascii;';
line11 = 'class dictionary;';
line12 = 'object blockMeshDict;';
line13 = '}';
line14 = '// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //';
line15 = '';
line16 = 'convertToMeters 1;';
line17 = 'vertices';
line18 = '(';
line19 = ');';
f1=fopen('blockMeshDict','w');
fprintf(f1,'%s\n',line1);
fprintf(f1,'%s\n',line2);
fprintf(f1,'%s\n',line3);
fprintf(f1,'%s\n',line4);
fprintf(f1,'%s\n',line5);
fprintf(f1,'%s\n',line6);
fprintf(f1,'%s\n',line7);
fprintf(f1,'%s\n',line8);
fprintf(f1,'%s\n',line9);
fprintf(f1,'\t %s \n',line10);
fprintf(f1,'\t %s \n',line11);
fprintf(f1,'\t %s \n',line12);
fprintf(f1,'%s\n',line13);
fprintf(f1,'%s\n',line14);
fprintf(f1,'%s\n',line15);
fprintf(f1,'%s\n\n',line16);
fprintf(f1,'%s\n',line17);
fprintf(f1,'%s\n',line18);
% Writing the vertices
fprintf(f1,'\t(%d %d %d)\n',p0(1), p0(2), p0(3));
fprintf(f1,'\t(%d %d %d)\n',p1(1), p1(2), p1(3));
fprintf(f1,'\t(%d %d %d)\n',p2(1), p2(2), p2(3));
fprintf(f1,'\t(%d %d %d)\n',p3(1), p3(2), p3(3));
fprintf(f1,'\t(%d %d %d)\n',p4(1), p4(2), p4(3));
fprintf(f1,'\t(%d %d %d)\n',p5(1), p5(2), p5(3));
fprintf(f1,'%s\n\n',line19);
% blocks
line20 = 'blocks';
line21 = '(';
line22 = 'hex (0 3 5 2 0 3 4 1) (200 50 1) simpleGrading (1 0.2 1)';
line23 = ');';
fprintf(f1,'%s\n',line20);
fprintf(f1,'%s\n',line21);
fprintf(f1,'\t%s\n',line22);
fprintf(f1,'%s\n\n',line23);
% Edges
arc1 = ' arc 2 1' ;
arc2 = 'arc 4 5';
fprintf(f1,'%s\n','edges');
fprintf(f1,'%s\n',line21);
fprintf(f1,'\t%s (%d %d %d)\n',arc1,0,r,0);
fprintf(f1,'\t%s (%d %d %d)\n',arc2,l,r,0);
fprintf(f1,'%s\n\n',line23);
% Boundary
line24 = 'Boundary';
line25 = '(';
line26 = 'Inlet';
line27 = '{';
line28 = '}';
line29 = ');';
inlet = '(0 1 2 0)';
outlet = '(3 5 4 3)';
fprintf(f1,'%s\n',line24);
fprintf(f1,'%s\n',line25);
fprintf(f1,'\t%s\n',line26);
fprintf(f1,'\t%s\n',line27);
fprintf(f1,'\t\t %s \n', 'type patch;');
fprintf(f1,'\t\t %s \n', 'faces');
fprintf(f1,'\t\t %s \n', line25);
fprintf(f1,'\t\t\t %s \n', inlet);
fprintf(f1,'\t\t %s \n', line29);
fprintf(f1,'\t %s \n\n', line28);
% Outlet
fprintf(f1,'\t%s\n','Outlet');
fprintf(f1,'\t%s\n',line27);
fprintf(f1,'\t\t %s \n', 'type patch;');
fprintf(f1,'\t\t %s \n', 'faces');
fprintf(f1,'\t\t %s \n', line25);
fprintf(f1,'\t\t\t %s \n', outlet);
fprintf(f1,'\t\t %s \n', line29);
fprintf(f1,'\t %s \n\n', line28);
% Walls
wall = '(2 1 4 5)';
fprintf(f1,'\t%s\n','fixedWalls');
fprintf(f1,'\t%s\n',line27);
fprintf(f1,'\t\t %s \n', 'type wall;');
fprintf(f1,'\t\t %s \n', 'faces');
fprintf(f1,'\t\t %s \n', line25);
fprintf(f1,'\t\t\t %s \n', wall);
fprintf(f1,'\t\t %s \n', line29);
fprintf(f1,'\t %s \n\n', line28);
% Axis
Axis = '(0 3 3 0)';
fprintf(f1,'\t%s\n','axis');
fprintf(f1,'\t%s\n',line27);
fprintf(f1,'\t\t %s \n', 'type wall;');
fprintf(f1,'\t\t %s \n', 'faces');
fprintf(f1,'\t\t %s \n', line25);
fprintf(f1,'\t\t\t %s \n', Axis);
fprintf(f1,'\t\t %s \n', line29);
fprintf(f1,'\t %s \n\n', line28);
fprintf(f1,'%s \n',');');
fprintf(f1,'\n');
%front
front=[1 0 3 4];
back=[0 2 5 3];
fprintf(f1,'\t%s\n','front');
fprintf(f1,'\t%s\n','{');
fprintf(f1,'\t\t%s\n','type wedge');
fprintf(f1,'\t\t%s\n','faces');
fprintf(f1,'\t\t%s\n','(');
fprintf(f1,'\t\t\t(%d %d %d %d)\n',front);
fprintf(f1,'\t\t%s\n',');');
fprintf(f1,'\t%s\n','}');
%back
fprintf(f1,'\t%s\n','back');
fprintf(f1,'\t%s\n','{');
fprintf(f1,'\t\t%s\n','type wedge');
fprintf(f1,'\t\t%s\n','faces');
fprintf(f1,'\t\t%s\n','(');
fprintf(f1,'\t\t\t(%d %d %d %d)\n',back);
fprintf(f1,'\t\t%s\n',');');
fprintf(f1,'\t%s\n','}');
fprintf(f1,'%s\n',');');
fprintf(f1,'%s \n','mergePatchPairs');
fprintf(f1,'%s \n \n','(');
fprintf(f1,'%s \n \n',');');
fprintf(f1,'\n');
last_line = '// ************************************************************************* //';
fprintf(f1,'%s \n',last_line);
fclose(f1);
BlockMesh Dict 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 1;
vertices
(
(0 0 0)
(0 9.996573e-03 2.617695e-04)
(0 9.996573e-03 -2.617695e-04)
(2.7 0 0)
(2.7 9.996573e-03 2.617695e-04)
(2.7 9.996573e-03 -2.617695e-04)
);
blocks
(
hex (0 3 5 2 0 3 4 1) (200 50 1) simpleGrading (1 0.2 1)
);
edges
(
arc 1 2 (0 1e-02 0)
arc 4 5 (2.7 1e-02 0)
);
boundary
(
inlet
{
type patch;
faces
(
(0 1 2 0)
);
}
outlet
{
type patch;
faces
(
(3 5 4 3)
);
}
surface
{
type wall;
faces
(
(2 1 4 5)
);
}
axis
{
type empty;
faces
(
(0 3 3 0)
);
}
front
{
type wedge; //symmetry;
faces
(
(0 3 4 1)
);
}
back
{
type wedge; //symmetry;
faces
(
(0 2 5 3)
);
}
);
mergePatchPairs
(
);
// ************************************************************************* //
Control Dict
/*--------------------------------*- 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 50;
deltaT 0.01;
writeControl timeStep;
writeInterval 200;
purgeWrite 0;
writeFormat ascii;
writePrecision 6;
writeCompression off;
timeFormat general;
timePrecision 6;
runTimeModifiable true;
// ************************************************************************* //
Constant physical properties
/*--------------------------------*- 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.90000000000000e-07;
// ************************************************************************* //
Expected Results
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

Magnitude of U:

Go to calculator>click on magnitude> select vector as U

Gradient of U

Variation of Tau
