OpenFOAM is an open source CFD package written in c++ that I have been interested in trying out for some time. I finally sat down and tried my first model using it - a sort of "hello world" example. My goal was to set up my simple dune geometry (dune of 400mm length and 20mm height)and run a simple flow over it. My objective was just to better understand how to define a model geometry in OpenFOAM and the basics for setting up a model run.
OpenFOAM has excellent documentation and also includes a complete set of tutorial models that are great places to start modelling from. I decided to generate this model using the icoFoam solver which solves the incompressible laminar Navier-Stokes equation (
more here). Again, in reality this case will be turbulent, but icoFoam is much simpler to set up than a turbulent model.
First step was to copy over and set up the required folder structure. I copied the files from the Cavity example which is one of the tutorials under the icoFoam tutorial directory.
Defining the Mesh
OpenFOAM does not have a GUI which can make things difficult and adds to the learning curve. This is particularly true for building the mesh. For this example I decided to use the BlockMesh mesh generation utility which creates parametric meshes (see
here for a further discussion on setting up the mesh).
Before I entered in my mesh I drew it out on paper. This was a key first step since there is no GUI and it can be difficult to visualize what you are doing otherwise. To define my mesh, I first determined all of the vertexes in my problem.
vertices
(
(0 20 0 ) //Vertex 0
(10 20 0 ) //Vertex 1
(10 150 0 ) //Vertex 2
(0 150 0 ) //Vertex 3
(0 20 10 ) //Vertex 4
(10 20 10 ) //Vertex 5
(10 150 10 ) //Vertex 6
(0 150 10 ) //Vertex 7
...
Next, I broke down my model geometry into "blocks" with each block containing 8 vertices. For each block you also define the mesh spacing:
hex (0 1 2 3 4 5 6 7) (5 20 1) simpleGrading (1 1 1)
In total I had 6 blocks that defined my model geometry. Once the blocks are defined, I need to identify patches (boundaries). Patches are defined by a collection of block-faces (e.g. 4 vertices) and have a defined type (e.g. wall or patch). This can be a little confusing - i suggest following through the tutorials in the OpenFoam documentation which gives a better understanding of this - once you get the hang of it, it is quite simple. The patch for the top of my model is given by (note, the name "topWall" is defined by the user):
wall topWall
(
(3 7 6 2)// watch the order of these RHR
(2 6 11 9)
(9 11 15 13)
(13 15 19 17)
(17 19 23 21)
(21 23 27 25)
(25 27 31 29)
)
The definition for the inlet to my model is as follows - note it is a type
patch as opposed to a
wall :
patch inlet
(
(0 4 7 3)
)
The most important thing to remember is that the vertices that define the boundary need to be in a specific order. The order is based on the right-hand-rule - using your right hand curl your fingers in the direction that you ordered the vertex, your thumb should be pointing out of the model. My full
blockMeshDict file looks like this:
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 1.5 |
| \\ / A nd | Web: http://www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
convertToMeters 0.001;
vertices
(
(0 20 0 ) //0
(10 20 0 )
(10 150 0 )
(0 150 0 ) //3
(0 20 10 )
(10 20 10 )
(10 150 10 )
(0 150 10 ) //7
(50 0 0 ) //8
(50 150 0 )
(50 0 10 )
(50 150 10 ) //11
(75 0 0 ) // 12
(75 150 0 )
(75 0 10 )
(75 150 10 ) //15
(137.5 2 0 )//16
(137.5 150 0 )
(137.5 2 10 )
(137.5 150 10 )//19
(325 18.40 0 )// 20
(325 150 0 )
(325 18.40 10 )
(325 150 10 )// 23
(390 20 0 ) //24
(390 150 0 )
(390 20 10 )
(390 150 10 ) //27
(400 20 0 )// 28
(400 150 0 )
(400 20 10 )
(400 150 10 )// 31
);
blocks
(
hex (0 1 2 3 4 5 6 7) (5 20 1) simpleGrading (1 1 1) // dx dy dz
hex (1 8 9 2 5 10 11 6) (30 20 1) simpleGrading (1 1 1)
hex (8 12 13 9 10 14 15 11) (20 20 1) simpleGrading (1 1 1)
hex (12 16 17 13 14 18 19 15) (20 20 1) simpleGrading (1 1 1)
hex (16 20 21 17 18 22 23 19) (40 20 1) simpleGrading (1 1 1)
hex (20 24 25 21 22 26 27 23) (10 20 1) simpleGrading (1 1 1)
hex (24 28 29 25 26 30 31 27) (5 20 1) simpleGrading (1 1 1)
);
edges
(
);
patches
(
wall topWall
(
(3 7 6 2)// watch the order of these RHR
(2 6 11 9)
(9 11 15 13)
(13 15 19 17)
(17 19 23 21)
(21 23 27 25)
(25 27 31 29)
)
wall bottomWall
(
(0 1 5 4)// watch the order of these RHR
(1 8 10 5)
(8 12 14 10)
(12 16 18 14)
(16 20 22 18)
(20 24 26 22)
(24 28 30 26)
)
patch inlet // note these are patches!!
(
(0 4 7 3)
)
patch outlet
(
(28 29 31 30)
)
empty frontAndBack
(
(0 3 2 1)// watch the order of these RHR
(4 5 6 7)
(1 2 9 8)
(5 10 11 6)
(8 9 13 12)
(10 14 15 11)
(12 13 17 16)
(14 18 19 15)
(16 17 21 20)
(18 22 23 19)
(20 21 25 24)
(22 26 27 23)
(24 25 29 28)
(26 30 31 27)
)
);
mergePatchPairs
(
);
// ************************************************************************* //
Once defined, I ran the blockMesh utility which generated the mesh. Just another note. When I was generating the mesh, I defined it using one block at a time and ran the meshing utility which made it much easier to identify problems. I then viewed the mesh in ParaView:
My mesh looks ok, however there is still work to be done to adjust the mesh. Once the mesh was built, I set up the initial/boundary conditions. For icoFaom, the velocity and pressure fields need to be defined. I used a simple inlet velocity and defined the outlet pressure.For all of the walls I defined them as being no-slip including the top wall (I will need to look up how to simulate an open channel type flow). Since this is a 2D case, the front and back walls are defined as being empty.
// * * * * * * * Velocity Field * * * * * * * * * * * * * * * * * //
dimensions [0 1 -1 0 0 0 0];
internalField uniform (0 0 0);
boundaryField
{
inlet
{
type fixedValue;
value uniform (0.5 0 0);
}
outlet
{
type zeroGradient;
}
topWall
{
type fixedValue;
value uniform (0 0 0);
}
bottomWall
{
type fixedValue;
value uniform (0 0 0);
}
frontAndBack
{
type empty;
}
}
// ************************************************************************* //
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 1.5 |
| \\ / A nd | Web: http://www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
object p;
}
// * * * * * * * * * * * Pressure * * * * * * * * * * * * * //
dimensions [0 2 -2 0 0 0 0];
internalField uniform 0;
boundaryField
{
inlet
{
type zeroGradient;
}
outlet
{
type fixedValue;
value uniform 0;
}
topWall
{
type zeroGradient;
}
bottomWall
{
type zeroGradient;
}
frontAndBack
{
type empty;
}
}
// ************************************************************************* //
I used the model run settings (defined in the
controlDict file) the same as for the cavity case. Other than that there was little else to do. I ran the model and viewed the results in ParaView:
|
Velocity Field |
|
Pressure Field |
|
Vertical Velocity (v) |
Again, I an not interested so much in the results of the model, but a good first step towards understanding how to set up a geometry and run the model.