Challenge

Week 11 - Simulation of Flow through a pipe in OpenFoam

In nonideal flueid dynamics,the higen_poiseuille equation, also known as Hagen-poiseuille law, iis a pysical law that gives the pressuree drop in an incompressible and Newtonian fluid in laminar flow flowing through a long cylindrical pipe of constant cross section

EQUATION:

In standard fluid-kinetics notation.

where:

Δp is the pressure dsifference between the two ends

L is the length of pipe

μ is the dynamic viscosity

Q is the volumetric flow rate

R is the pipe radius

A is the cross section of pipe

The equation does not hold close to p[ipe entrance.

Given parameters.

  1. Wedge angle less than 5 degrees.
  2. Flow is axi symetrical.
  3. Fluid is incompressible and the flow is laminar.

Expected Velocity Profile.

Assuming the following:

  1. The distance the fluid has to travel to become fully developed is given by the Entry Length formula.
  2. Given Reynolds number is 2100.This implies flow is laminar.
  3. Hence Entry length formula for Laminar is given by:L/d=0.06.Re
  4. where L is effective length; d is dia of pipe and Re is reynolds number.
  5. Assuming the dia of pipe as 10mm;substituting the values,we get L=1.2m;

Steps to follow

1. First, the geometry is created in the blockmesh Dict file, this is done by creating the file through the matlab code.

2. The mesh is generated according to requirement.

3. The velocity and preassure bc are applied

4. The controlDict file is edited for Simulation run time and timestep;

5. After the geometry is input, the simulation can be started by calling 'icoFoam' command in the terminal.

6. 'paraFoam' command helps to open paraView

7. The shear stress diagram is obtained by post-processing the velocity results obtained from Openfoam.

Matlab programm for automatically generate blockmeshDict file for a wedge angle provided by user.

1. Main Code

clear all
close all
clc
% problem initialization
Re=2100;
rho=1000;
mu=0.89e-3;
nu=mu/rho;

D=0.02;
A_c=(pi/4)*D^2;
p=pi*D;
D_h=4*A_c/p;
r_h=D_h/2;

v_avg=(Re*mu)/(rho*D_h);
v_max=2*v_avg;
tau_max=2*mu*v_max/r_h;

if Re<2300
   L_e=0.06*Re*D_h;
elseif Re==2300
    L_e=115*D_h;
elseif Re>4000
    L_e=1.359*D_h*Re^(1/4);
end

L=L_e+0.18; % length of pipe

%% Pressure loss
f=65/Re;

% Haigen Poiseuille Eqn
P=(32*mu*(L-L_e)*v_avg)/(D^2);
g=9.81;
h_loss=P/(rho*g);
P_kinematicPress=P/rho;
theta=3;
r=r_h;

%% Print out
system(r,L,theta)

2. system()

function system(r,L,theta)

%% creating blockmesh dict file
%%opening part
h1='/=------------------------*-C++ -*------------------------*';
h2=' ========             |';
h3=' \      / F ield      | opennFoam: The open source CFD toolbox';
h4='  \    /  O peration  | website:   http://openfoam.org';
h5='   \  /   A nd        | Version:   9';
h6='    \/    M anipulation |';
h7='--------------------------------------*/';
h8='FoamFile';
h9='{';
h10=' format   ascii;';
h11=' class    dictionary;';
h12=' object   blockMeshDict;';
h13='}';
h14='// *********************************//';
f1=fopen('blockMeshDict.txt','wt');

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',h14);

%%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)];
mid1=[0 r 0];
mid2=[L r 0];

%%write vertices
h15='convertToMeters 1;';
fprintf(f1,'%s\n',h15);

h16='vertices';
fprintf(f1,'%s\n',h16);
h17='(';
fprintf(f1,'%s\n',h17);
fprintf(f1,'(%d %d %d)\n',v0);
fprintf(f1,'(%d %d %d)\n',v2);
fprintf(f1,'(%d %d %d)\n',v3);
fprintf(f1,'(%d %d %d)\n',v4);
fprintf(f1,'(%d %d %d)\n',v5);
h18=');';
fprintf(f1,'%s\n',h18);

%% write blocks
h19='blocks';
h20='(';
h21='  hex (0 3 4 1 0 3 5 2) (200 50 1) simpleGrading (1 0.2 1)';
h22=');';
fprintf(f1,'%s\n',h19);
fprintf(f1,'%s\n',h20);
fprintf(f1,'%s\n',h21);
fprintf(f1,'%s\n',h22);

%% write edges
h23='edges';
h24='(';
h25='  arc 2 1';
h26='  arc 5 4';
h27=');';
fprintf(f1,'%s\n',h23);
fprintf(f1,'%s\n',h24);
fprintf(f1,'%s (%d %d %d)\n',h25,0, r, 0);
fprintf(f1,'%s (%d %d %d)\n',h26,L, r, 0);
fprintf(f1,'%s\n',h27);

%% write boundary
inlet= [0 1 2 0];
outlet=[3 5 4 3];
surface=[2 1 4 5];
axis=[0 3 3 0];
front=[1 0 3 4];
back=[0 2 5 3];

fprintf(f1,'%s\n','boundary');
fprintf(f1,'%s\n','(');
%inlet patch
fprintf(f1,'%s\n','inlet');
fprintf(f1,'%s\n','{');
fprintf(f1,'%s\n','type patch');
fprintf(f1,'%s\n','faces');
fprintf(f1,'%s\n','(');
fprintf(f1,'(%d %d %d)\n',inlet);
fprintf(f1,'%s\n',');');
fprintf(f1,'%s\n','}');

%outlet patch
fprintf(f1,'%s\n','outlet');
fprintf(f1,'%s\n','{');
fprintf(f1,'%s\n','type patch');
fprintf(f1,'%s\n','faces');
fprintf(f1,'%s\n','(');
fprintf(f1,'(%d %d %d)\n',outlet);
fprintf(f1,'%s\n',');');
fprintf(f1,'%s\n','}');

%surface
fprintf(f1,'%s\n','surface');
fprintf(f1,'%s\n','{');
fprintf(f1,'%s\n','type wall');
fprintf(f1,'%s\n','faces');
fprintf(f1,'%s\n','(');
fprintf(f1,'(%d %d %d)\n',surface);
fprintf(f1,'%s\n',');');
fprintf(f1,'%s\n','}');

%axis
fprintf(f1,'%s\n','axis');
fprintf(f1,'%s\n','{');
fprintf(f1,'%s\n','type empty');
fprintf(f1,'%s\n','faces');
fprintf(f1,'%s\n','(');
fprintf(f1,'(%d %d %d)\n',axis);
fprintf(f1,'%s\n',');');
fprintf(f1,'%s\n','}');

%front
fprintf(f1,'%s\n','front');
fprintf(f1,'%s\n','{');
fprintf(f1,'%s\n','type wedge');
fprintf(f1,'%s\n','faces');
fprintf(f1,'%s\n','(');
fprintf(f1,'(%d %d %d)\n',front);
fprintf(f1,'%s\n',');');
fprintf(f1,'%s\n','}');

%back
fprintf(f1,'%s\n','back');
fprintf(f1,'%s\n','{');
fprintf(f1,'%s\n','type wedge');
fprintf(f1,'%s\n','faces');
fprintf(f1,'%s\n','(');
fprintf(f1,'(%d %d %d)\n',back);
fprintf(f1,'%s\n',');');
fprintf(f1,'%s\n','}');
fprintf(f1,'%s\n',')');

fprintf(f1,'%s\n','mergePatchPairs');
fprintf(f1,'%s\n','(');
fprintf(f1,'%s\n',')');
fprintf(f1,'%s\n','// ******************************************* //');

end

OpenFoam Files:

1. System Folder:

a. blockMeshDict:

/*--------------------------------*- 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 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 4 1 0 3 5 2) (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 2 1 0)
		 );
	}
	outlet
	{
		type patch;
		faces
		(	
			(3 4 5 3) 
		);
	}

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

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

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

mergePatchPairs

(
);

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

controlDict:

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

deltaT          0.01;

writeControl    timeStep;

writeInterval   200;

purgeWrite      0;

writeFormat     ascii;

writePrecision  4;

writeCompression off;

timeFormat      general;

timePrecision   6;

runTimeModifiable true;


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

velocity (U):

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

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

    outlet
    {
        type            zeroGradient;
    }
    
    surface
    {
    type               noSlip;
    //type              fixedValue;
    //value            uniform (0 0 0);
    }
    
    axis
    {
    type               empty;
    }
    
    front
    {
        type            wedge;
    }
    
    back
    {
        type            wedge;
    }
}

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

pressure (P):

/*--------------------------------*- 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
{
    inlet
    {
        type            zeroGradient;
    }
   
   outlet
    {
        type            fixedValue;
        value           uniform 0;
    }
    
    surface
    {
         type            zeroGradient;
    }
    
    axis
    {
         type            empty;
    }
    
    	front
    {
        type            wedge;
    }
	
	back
    {
        type            wedge;
    }
}

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

Results and graph:

For wedge angle 3degree

1.Velocity line plot

velocity line plot across the radial distance X=2.5m

velocity line plot across the radial distance X=0.05m

maximum shear stress plot:

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.

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.