• Run-time postprocessing
  • Prerequisites
  • 1: Setup of initial case
  • 2: Mesh generation and boundary conditions
  • 3: Defining output data
    • Total filling level
    • Cut plane and surface profile
    • Rate of flow out of domain
  • 4: Running the simulation
  • 5: Results and postprocessing
  • 6: Further work
  • 7: Download

Run-time postprocessing

In the earlier spillway and swak4Foam tutorials Paraview were used a lot for postprocessing. This works rather well on small and moderate cases, but when it comes to really large cases with tenths of millions of cells, this might not be practical. In such cases you should consider what you need, and configure OpenFOAM to give you those data during the simulation (and not try to extract it with Paraview later). Then the data will be extracted from the analysis while it is run in parallel. Most of the common prostprocessing tasks can be done this way. The surface elevation probe in the swak4Foam tutorial is one example of such run-time postprocessing, the forces and lift coefficient calculation in the airfoil is another.

The purpose of this tutorial is to:

• Show how you can use programs such as GNU Octave to create meshes with nonregular shapes by using blockMesh only
• Show a few of the possible run-time postprocessing features of OpenFoam and swak4Foam

The case we will use as an example is the emptying of a bottle mostly filled with a liquid held with the bottom up. The interesting questions are:

• What is the filling level at a certain time?
• At what rate does the liquid flow out?
• How does the liquid-air interface look in the bottle's neck?
• And how does the liquid surface look like?

Prerequisites

Before you start on this tutorial, please make sure that you have

• A decent laptop or workstation
• A working OpenFOAM 2.x.x installation
• GNU Octave (Matlab will probably also work with minor modifications, although not tested)
• Done the basic cavity and damBreak tutorials that is distributed together with OpenFOAM and described in the documentation
• Done the spillway and swak4Foam tutorials provided here (or in any other way installed and familiarized yourself with swak4Foam)

1: Setup of initial case

The infamous damBreak tutorial will be used as a base case. You might choose yourself weather you want to use turbulence modelling or not. In the rest of this guide it is assumed that the simple laminar base case is used. Copy the tutorial into a folder of your choise, and possibly rename the folder to something else, such as beerBottle:

cp -r $FOAM_TUTORIALS/multiphase/interFoam/laminar/damBreak ./beerBottle

You should also switch the direction of gravity from the negative y-direction to the negative z-direction as we have done earlier, i.e. the acceleration field should be (0 0 -9.81) m^2/s.

2: Mesh generation and boundary conditions

The domain will for computational reasons be 2D, and the dimensions are as shown in the figure above. As this might be a bit difficult to implement by hand, we have prepared a small script meshgen.m that writes a constant/polyMesh/blockMeshDict file. This is an extremely fast and efficient way do generate simple meshes compared to the method previous used with snappyHexMesh and extrudeMesh. Please look at the script and the blockMeshDict it produces, and make sure you understand. If in doubt about the syntax of the blockMeshDict-file, please consult the OpenFOAM documentation. The mesh is very easily generated by the commands:

chmod +x meshgen.m
./meshgen.m
blockMesh

You should set the boundary conditions yourself (hint: the inletOutlet patch should be defined with a inletOutlet boundary condition).

3: Defining output data

As previously stated we are interested in some general output variables. Some of this can be extracted by using native OpenFOAM tools, and some output is easiest generated by swak4Foam. Brefore proceeding you should add some libraries to system/controlDict and include some files:

libs (
                "libOpenFOAM.so"
                "libsimpleSwakFunctionObjects.so"
                "libswakFunctionObjects.so"
            );
        
       functions
       {
           #include "sampledSets"
           #include "sampledSurf"
       }

Now create the files system/sampledSets and system/sampledSurf:

touch system/sampledSets
                touch system/sampledSurf

You can of course also just put everything in one file, or even in ciontrolDict itself, but some prefer to put postprocessing configuration in separate files (like done here).

Total filling level

The total liquid filling is easiest calculated as the volume integral of the phase fraction. This is done by a native OpenFOAM volumeIntegrate function. Add the following to the system/sampledSets-file:

totalLiquid
                {
                    type volumeIntegrate;
                    fields ( alpha1 );
                    verbose true;
                }

Cut plane and surface profile

Since we are interested in the air-liquid interface in the bottle's neck, we must define a cut plane through the section we are interested in. For the liquid surface we must calculate an isosurface with alpha1=0.5. Both operations is some by the native OpenFOAM functions cuttingPlane and isoSurface. Add the following to system/sampledSurf:

cuttingPlane
                {
                    type                      surfaces;
                    functionObjectLibs        ("libsampling.so");
                    outputControl             outputTime;
                 
                    surfaceFormat             vtk;
                    fields                    ( p U alpha1 );
                 
                    interpolationScheme        cellPoint;
                 
                    surfaces
                    (
                        neckCut
                        {
                            type                cuttingPlane;
                            planeType           pointAndNormal;
                            pointAndNormalDict
                            {
                                basePoint       (0 0 0.015166667);  // Notet that the plane does not match up with the mesh
                                normalVector    (0 0 1);
                            }
                            interpolate         true;
                        }
                         
                        freeSurface
                        {
                            type                isoSurface;
                            isoField            alpha1;
                            isoValue            0.5;
                            interpolate         true;
                        }
                    );
                }

Unless any modifications are done to the mesh, the cutting plane is not intersecting any internal faces in the mesh. This is done intentionally to show that you can define completely arbitrary planes without concern about the location of the internal faces. It would be quite a lot of work if all of the sample planes were to be relocated each time you changed the mesh resolution (for example in a convergence study).

Interpolation and integration

You should be aware that both integration and interpolation at run-time is costly operations (especially integrating over interpolated planes), and you should carefully consider your needs. If you are careless in your definitions you could end up using more CPU time on postprocessing than the analysis itself (however this might of course actually be what you want).

Rate of flow out of domain

The rate of flow out of the domain is easiest calculated using a swakExpression. Append the following lines to your system/sampledSets-file:

volFlow
                {
                    type                swakExpression;
                    valueType           patch;
                    patchName           inletOutlet;
                     
                    verbose             true;
                    expression          "alpha1*U&Sf()";
                    accumulations       ( sum );
                }

4: Running the simulation

Now you have to configure your decomposition the way you want it (done in system/decomposeParDict), decompose the domain and execute the simulation:

decomposePar
                mpirun -np 4 interFoam -parallel

5: Results and postprocessing

Your output planes, surfaces and variables should show up in folders created in the root directory of the case. If you for example want to plot the filling level together with the flow rate of liquid out of the domain, you can do this in gnuplot with the following commands:

set xlabel 'Time (s)'
set ylabel 'Volume of liquid (m^3)'
set y2label 'Flow rate (m^3/s)'
 
set y2tics
set yrange [0:*]
set y2range [0:*]
set grid
 
set term png size 800,600
set output 'fillLevel.png'
plot 'volumeIntegrate_totalLiquid/0/alpha1' axes x1y1 with lines title 'Volume of liquid', \
     'swakExpression_volFlow/0/volFlow' axes x1y2 with lines title 'Flow rate'

The resulting plot will look like:

Try to look at the other results you get, the volume flow and the isosurfaces.

You can of course still open the case in Paraview. Here is two stills form the discharge process:

6: Further work

If you want to test yourself's blockMeshDict-skills, you can try to make a full 3D grid. This should not be too difficult...

7: Download

You can download all necessary case files in beerBottle.tar.gz. Scripts to run and clean the case is also attached. Note that you need Gnuplot make the fill level plot.