2014 OFoam Tut Complete

Tutorials & Lecture notes for “166.049 Fluiddynamik (CFD) Thermischer Trennverfahren”

2nd edition, Oct. 2014

This offering is not approved or endorsed by ESI® Group, ESI-‐OpenCFD® or the OpenFOAM® Foundation, the producer of the OpenFOAM® software and owner of the OpenFOAM® trademark.

Editors: • Jozsef Nagy • Christian Jordan • Michael Harasek • Bahram Haddadi Cover picture from: • J. Nagy und M. Harasek. Investigation of the aerobreakup of a liquid droplet at high Weber number with different turbulence models. In Proceedings of the Fifth International Conference from Scientific Computing to Computational Engineering, Athens, Greece, 2012. •

J. Nagy, A. Horvath, C. Jordan, und M. Harasek. Turbulent phenomena in the aerobreakup of liquid droplets. CFD Letters, 4(3) 2012.

Attribution-NonCommercial-ShareAlike 3.0 Unported (CC BY-NC-SA 3.0) This is a human-readable summary of the Legal Code (the full license). Disclaimer You are free: to Share — to copy, distribute and transmit the work to Remix — to adapt the work Under the following conditions: Attribution — You must attribute the work in the manner specified by the author or licensor (but not in any way that suggests that they endorse you or your use of the work). Noncommercial — You may not use this work for commercial purposes. Share Alike — If you alter, transform, or build upon this work, you may distribute the resulting work only under the same or similar license to this one. With the understanding that: Waiver — Any of the above conditions can be waived if you get permission from the copyright holder. Public Domain — Where the work or any of its elements is in the public domain under applicable law, that status is in no way affected by the license. Other Rights — In no way are any of the following rights affected by the license: Your fair dealing or fair use rights, or other applicable copyright exceptions and limitations; The author's moral rights; Rights other persons may have either in the work itself or in how the work is used, such as publicity or privacy rights. Notice — For any reuse or distribution, you must make clear to others the license terms of this work. The best way to do this is with a link to this web page.


Tutorial One: elbow Solver: icoFoam Geometry: 2 dimensional Purpose: Different meshes

Tutorial Two: forwardStep Solver: sonicFoam Geometry: 2 dimensional Purpose: Built in mesh

Tutorial Three: shockTube Solver: sonicFoam Geometry: 1 dimensional Purpose: Patching fields

Tutorial Four: shockTube Solver: scalarTransportFoam Geometry: 1 dimensional Purpose: Discretization

Tutorial Five: circle Solver: scalarTransportFoam Geometry: 2 dimensional Purpose: Discretization

Tutorial Six: pitzDaily Solver: scalarTransportFoam Geometry: 2 dimensional Purpose: Physics of the problem

Tutorial Seven: pitzDaily Solver: simpleFoam Geometry: 2 dimensional Purpose: Steady state, Turbulence, Parameter

i


Tutorial Eight: pitzDaily Solver: pisoFoam Geometry: 2 dimensional Purpose: Turbulence, Parameter

Tutorial Nine: damBreak Solver: interFoam Geometry: 2 dimensional Purpose: Multiphase

Tutorial Ten: depthCharge3D Solver: compressibleInterFoam Geometry: 3 dimensional Purpose: Compressible, Incompressible, Parallel

Tutorial Eleven: depthCharge3D Solver: compressibleInterFoam Geometry: 3 dimensional Purpose: Parallel processing, Manual method in parallel processing

Tutorial Twelve: TJunction Solver: simpleFoam, scalarTransportFoam Geometry: 3 dimensional Purpose: Residence Time Distribution

Appendix A: Important commands on Linux

Appendix B: Frequently Asked Questions (FAQ)

Appendix C: paraView

Appendix D: OpenFOAM® solvers

ii


icoFoam – elbow (mesh)1

Simulation Using icoFoam solver simulate 75 s of flow in an elbow for following GAMBIT meshes. • Tri-mesh (comes with OpenFOAM® tutorial) • hex-mesh coarse (created in GAMBIT tutorial) • 2xfiner hex-mesh (refine GAMBIT mesh) Objectives • • • •

Looking at the initial values for p and U Make sure all boundaries in GAMBIT are called accordingly Make sure after importing the mesh no additional surfaces are introduced Make sure the surfaces in constant/polyMesh/boundary are defined correctly.

Post processing Follow step-by-step simulation, import your simulation into paraview, then copy your data from server to your computer and make along this line two diagrams (e.g. by Excel) of pressure and velocity magnitude respectively including the results of all three simulations.

1 Provided by Bahram Haddadi

1


Step by step simulation Log into SSH server: Windows: Run PuTTY. Set the following: Category>Session

Host name: openhost.university.edu Connection type: SSH Category>Connection>SSH>Tunnels

Source port: 5901 Destination: localhost:59**1 Press Add (don’t forget to press Add!) Please don’t use any other Display number, so that you don’t disturb other students. Category>Connection>Data

Auto-login username: openFoamUser2 Category>Session

Saved Sessions: openFoamUser Press Save Now from saved sessions choose your session (openFoamUser) and press Open. In the opened Command Prompt window, it prompts for your password. The password is not echoed to the screen. To log out use whatever command is used to logout from the server you are logged into (typically ctrl + d). Mac OS X and Linux: Open your Terminal application. You will see a window with a $ symbol and a blinking cursor. From here, you may issue the following command to establish the SSH connection to your server (be careful about caps lock ‘L’ in the -gL). >ssh –gL 5901:localhost:59** [email protected]

Immediately after issuing this command, your computer will attempt to establish a connection to your server. If it is your first time connecting to that server, you will see a message asking you to confirm the identity of the machine. Make sure you entered the address properly, and type yes, followed by the return key, to proceed. You will then be prompted to enter your password, Type or copy/paste your SSH user's password into Terminal. You will not see the cursor move while entering your password. This is normal. Once you are finished inputting your password, press return on your keyboard. At that point, you will be connected to your server remotely through SSH. For Linux it’s also the same.

1 Display number (which is already given to you) 2 Session ID (which is already given to you) 2


Run VNC Install the appropriate VNC Viewer and run it: VNC Server: localhost: 01 Press Connect Press Continue Enter your password Press Ok Open terminal Left click on VNC desktop from opened menu choose Xterm, in this toturials openfoam version 2.3.0 is used, for using it, in the opened terminal execute following command: >. ./bashrc_OF230

for older versions according to the version of openfoam you want to use, in the opened terminal execute one of following commands: >. ./bashrc_OF171

for openfoam version 1.7.1, or >. ./bashrc_OF221

for openfoam version 2.2.1. In different versions some solver may not be available or they may need different settings for executing. Open tutorial In the opened window type following command to change your directory to elbow directory: >cd ~/OpenFOAM/OpenFOAM-2.3.0/tutorials/incompressible/icoFoam/ elbow

Copying mesh Copy elbow.msh to this directory: >cp {mesh file address}/elbow.msh .

Note: the mesh which is imported to OpenFOAM® should be a three dimensional mesh. For simulating 2D (also the same 1D) simulations a three-dimensional mesh should be created with just one cell in the third direction (For 1D, one cell in second and third directions). Converting mesh The mesh, which is produced by Gambit, is not appropriate for OpenFOAM® solver, so first the mesh should be converted to OpenFOAM® mesh, for this:

3


>fluentMeshToFoam elbow.msh –scale 1.0

the -scale flag is used for converting the created mesh dimensions from other units to SI1 units, e.g. if the mesh is created in mm by using –scale 0.001 it will be converted to meter and if the flag is omitted, uses 1. Note: if there are internal boundaries in the mesh, there is another command, fluent3DMeshToFoam. Using this command the internal boundaries will be kept during conversion. Checking files If you check inside the elbow tutorial folder, >ls

there are three directories (0, constant, system) and one file (elbow.msh). Initial value The initial values are set in directory 0, checking inside it, >ls 0

there are two files: p and U. p is initial value for pressure, checking it: >nano2 p

it will be like this: /*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 2.3.0 | | \\ / A nd | Web: www.OpenFOAM.com | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class volScalarField; object p; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions

[0 2 -2 0 0 0 0];

internalField

uniform 0;

boundaryField {

1 International System of Units 2 nano is the text editor used in Linux OS (for closing it ctrl+x) 4

wall-4 { type } velocity-inlet-5 { type } velocity-inlet-6 { type } pressure-outlet-7 { type value } wall-8 { type } frontAndBackPlanes { type }


zeroGradient;

zeroGradient;

zeroGradient;

fixedValue; uniform 0;

zeroGradient;

empty;

} // ************************************************************************* //

In the dimensions the dimension of the quantity is defined, for example here it shows that the p dimension is (m/s)2. Note: as you can see here the p unit is not pressure unit (Pa), it is due to the fact that in the incompressible solvers in OpenFOAM®, p is defined as pressure divided by density. Note: in OpenFOAM® all units are SI units. Note: in the dimension matrix the first number presents mass unit power, the second one the length, the third one time and the forth one the temperature. The internalField sets the initial field of that quantity in the solution domain. The type of each of our boundaries and value of that quantity on the boundaries are defined in the boundaryField (see introduction). Note: There is new type of boundary here: empty. This type of boundary is for sides, which are vertical to the direction that are not going to be considered (e.g. in 2D simulations this boundaries are vertical to the third direction). In this boundary type both of the sides vertical to one direction should be selected together and named as one boundary. U is the initial value for velocity field, checking it: >nano U

it’s also the same as p. As velocity is a vector, it has to be defined via three components, where the first one stands for x component, second one for y component, and the third

5


one z component. z component is always zero because it’s a 2D simulation and no calculations will be done in the z direction (The boundaries vertical to z direction have been already set to empty). constant In the constant folder there is a directory and a file. In the directory (polyMesh) the mesh data, which were converted, are stored. >cd polyMesh >nano boundary

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // 6 (

…

wall-4 { type nFaces startFace } velocity-inlet-5 { type nFaces startFace } frontAndBackPlanes { type nFaces startFace }

wall; 100; 1300;

patch; 8; 1400;

empty; 1836; 1454;

) // ************************************************************************* //

Check the boundary names here with what you set in Gambit, they should be the same, and also check the boundary types (Walls should be wall, inlet and outlets should be patch, empties should be empty). Starting cell number and also number of each face cells can also be checked here. By opening the transportProperties file, properties dimensions and also after that property amount can be found and edited, e.g.: nu

nu [ 0 2 -1 0 0 0 0 ] 0.01;

nu is the fluid kinematic viscosity, which is 0.01 for this example. system System settings can be found and changed in this directory. In the fvSchemes file the discretization scheme that is used for each term of the equations can be found and edited. The fvSolution file contains the settings to the coupling method of pressure and velocity

6


and also the solver, which is used for each quantity, and also the final tolerance for convergence of that quantity. Numbers of orthogonal and nonorthogonal corrections for pressure are also set here (see introduction). In the controlDict file, the time, when simulation starts (startFrom), the time when the simulation finishes (stopAt), time step (deltaT), data saving interval (writeInterval), saved data file format (writeFormat), saved file data precision (writePrecision), and also if changing the files during the run can affect the run or not (runTimeModifiable). // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // application

icoFoam;

startFrom

latestTime;

startTime

0;

stopAt

endTime;

endTime

75;

deltaT

0.05;

writeControl

timeStep;

writeInterval

20;

purgeWrite

0;

writeFormat

ascii;

writePrecision

6;

writeCompression uncompressed; timeFormat

general;

timePrecision

6;

runTimeModifiable yes; // ************************************************************************* //

Note: in here the simulation continues from the last time step data which is saved (latestTime), and if there be no saved data it will start from start time (startTime), which is zero here. Running simulation Simulation can be run by type the solver name and execute it: >icoFoam

Note: for running the simulation the solver command (e.g. icoFoam) should be executed from the tutorial main folder, for example in here the command should be executed at icoFoam folder, if it be at some subfolders or somewhere else the simulation will fail.

7


Exporting simulation The data files created by OpenFOAM® should be exported by the appropriate command to the post processing tools data format. for paraview: >foamToVTK

where VTK is paraview data format. Post processing Here, paraview is used as the post processing tool, for running it, >paraview &

Note: By putting & at the end of command, the command line will be also active as that program is running. The comparisons and charts are shown below.

8

Mesh


Pressure

Velocity

Tri

Hex

Hex Fine

Figure 1.1 Comparison of different mesh types results at t=75 s

9


4

Velocity

3.5 3 2.5 U

U(Tri)

2

U(Hex)

1.5

U(FineHex)

1 0.5 0 0 2 5 7 10 12 14 17 19 22 24 27 29 31 34 36 39 41 43 46 48 Arc Length(m) 2

Pressure/Density

1.5 1 0.5 P

P(Tri)

0

P(Hex)

-‐0.5

P(FineHex)

-‐1 -‐1.5 -‐2 0 2 5 7 10 12 14 17 19 22 24 27 29 31 34 36 39 41 43 46 48 Arc Length(m)

Figure 1.2 Pressure and velocity for different meshes at t=75 s, along the arc shown The comparison plots were along a line between the middle of small tube entrance and middle of large tube exit part.

10


sonicFoam – forwardStep1

Simulation Using sonicFoam solver simulate 10 s of flow over a forward step. Objectives • •

Understand blockMesh Find out how to define vertices via coordinates as well as surfaces and volumes via vertices.

Post processing Follow step-by-step simulation, import your simulation into ParaView, and examine the mesh and the results in detail.


11


Step by step simulation Log into SSH server Run VNC Open tutorial >cd ~/OpenFOAM/OpenFOAM-2.3.0/tutorials/compressible/sonicFoam/ laminar/forwardStep

Checking files >ls 0

constant

system

Initial value >ls 0 p

T

U

T is temperature initial value. Internal pressure and temperature fields are set to 1, and the initial velocity in the domain is set to zero except at the inlet boundary, where it is 3. Note: as it can be seen here the p unit is the same as pressure unit, it’s because of that sonicFoam is a compressible solver. Note: Don’t forget, this example is a purely numeric example (you might have noticed that from pressure values). Constants >ls constants polyMesh thermophysicalProperties

turbulenceProperties

In the thermophysicalProperties file gas type is set to ideal gas are set. By opening the turbulenceProperties can be set appropriate turbulent mode (Here it’s laminar): simulationType

laminar;

Opening polyMesh directory, there are files where the mesh properties are stored. In this example the mesh is not imported from other programs (e.g. GAMBIT), and it will be created inside OpenFOAM®, for that the blockMesh command is used. blockMeshDict file is like this: >cd constants >cd polyMesh >ls blockMeshDict

boundary

faces

neighbour

owner

points

>nano blockMeshDict // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // convertToMeters 1;

12


vertices ( (0 0 -0.05) (0.6 0 -0.05) (0 0.2 -0.05) (0.6 0.2 -0.05) (3 0.2 -0.05) (0 1 -0.05) (0.6 1 -0.05) (3 1 -0.05) (0 0 0.05) (0.6 0 0.05) (0 0.2 0.05) (0.6 0.2 0.05) (3 0.2 0.05) (0 1 0.05) (0.6 1 0.05) (3 1 0.05) ); blocks ( hex (0 1 3 2 8 9 11 10) (25 10 1) simpleGrading (1 1 1) hex (2 3 6 5 10 11 14 13) (25 40 1) simpleGrading (1 1 1) hex (3 4 7 6 11 12 15 14) (100 40 1) simpleGrading (1 1 1) ); edges ( ); boundary ( inlet { type patch; faces ( (0 8 10 2) (2 10 13 5) ); } outlet { type patch; faces ( (4 7 15 12) ); } bottom { type symmetryPlane; faces ( (0 1 9 8) ); } top { type symmetryPlane; faces ( (5 13 14 6)

13


(6 14 15 7) );

);

} obstacle { type patch; faces ( (1 3 11 9) (3 4 12 11) ); }

mergePatchPairs ( ); // ************************************************************************* //

As told before units in OpenFOAM® are SI units, if the vertex coordinates are different from SI, they can be converted with the convertToMeters command to SI units, the number in the front of convertToMeters shows the constant, which should be multiplied in the dimensions to change them to meter (SI unit of length). For example convertToMeters

0.001

could show that the dimensions are in millimeter, and by multiplying them into 0.001 they are converted to meters. In the vertices part, the geometry vertices coordinates are defined, the vertices are stored and numbered from zero, e.g. vertex (0 0 -0.05) is numbered zero, and vertex (0.6 1 0.05) points to number 6. In the block part, blocks are defined that from which vertices they are made, e.g. the first block is made of vertices (0 1 3 2 8 9 11 10). After each block the mesh in every direction is defined. e.g. (25 10 1) shows that this block is divided to 25 parts in x direction, 10 parts in y direction and 1 part in z direction . As it was explained before even for 2D simulations the mesh and geometry should be 3D, but with one cell in the direction, which is not going to be simulated, e.g. here number of cells in z direction is one and it’s because of that it’s a 2D simulation in x-y plane. The last part (simpleGrading (1 1 1)) shows the size function. At the patches part each boundary is defined by vertices it is made up, and also its type and name are defined. Note: for creating a face the vertices should be chosen clockwise when look at them from inside of geometry. System >cd ../.. >cd system >ls controlDict

fvScheme

fvSolution

14


Running simulation Before running simulation the mesh has be created. In the previous step the mesh and geometry were set, for creating it the following command should be executed: >blockMesh

and after that the mesh is created in the polyMesh folder. For running simulation type the solver name and execute it: >sonicFoam

Exporting simulation >foamToVTK

The mesh is seen like this in ParaView, and you can easily see the three blocks, which were created.

Figure 2.1 Mesh generated by blockMesh Note: the mesh in here is a triangular mesh, in fact ParaView changes the mesh to triangular mesh for visualization, however it’s a hex mesh (see introduction) where every square is represented by two triangles.

15


The simulation results are as follows: Time

Pressure

Velocity

Temperature

0.5s

1s

10s

Figure 2.2 Pressure, velocity and temperature contours at different time steps

16


sonicFoam – shockTube1

Simulation Using sonicFoam solver simulate 0.007s of flow inside a shock tube, with a mesh with 100, 1000 and 10000 cells in one dimension, for initial values 10bar/1bar and 100bar/1bar. Objectives Understanding setFields Post processing Follow step-by-step simulation, import your simulation into ParaView, and compare results.


17


Step by step simulation Log into SSH server Run VNC Open tutorial >cd ~/OpenFOAM/OpenFOAM-2.3.0/tutorials/compressible/sonicFoam/ laminar/shockTube


Allclean

Allrun

constant

system

Initial value >ls 0 magU p T U

magU is magnitude of velocity vector, which is not used here. Constants >ls constants polyMesh thermophysicalProperties

turbulenceProperties

The model is a laminar model, and thermo physical properties of an ideal gas are applied. >cd constants >cd polyMesh >ls blockMeshDict

boundary

>nano blockMeshDict

By checking geometry and mesh, it is obvious that it’s a 1D mesh, because number of mesh cells in y and z directions is one, and also in the patches plates vertical to these directions are defined as empty boundary condition. The mesh density can be set in the blocks part by changing x direction mesh size (e.g. change it from 1000 to 100 or 10000). System >cd ../.. >cd system >ls controlDict

fvScheme

fvSolution

18

sampleDict

setFieldDict


setFieldDict is used for patching (assign an amount to a region) in the simulation. For example here pressure of 10bar should be patched to half of the region (the geometry is from -5 to 5, so from 0 to 5 will be patched) and 100bar to the other half. >nano setFieldDict

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // defaultFieldValues ( volVectorFieldValue U ( 0 0 0 ) volScalarFieldValue T 348.432 volScalarFieldValue p 1000000 ); regions ( boxToCell { box ( 0 -1 -1 ) ( 5 1 1 ) ; fieldValues ( volScalarFieldValue T 278.746 volScalarFieldValue p 100000 ) ; } ); // ************************************************************************* //

In the defaultFieldValues, a value is assigned to the whole domain, for example here velocity has been set everywhere to zero, temperature 348.432, and pressure 1000000. In the regions at first in the boxToCell part the region to which we want to assign a special amount (with boxToCell the region is chosen by a cube, and the cube is defined by giving the coordinates of one of its diagonals) is defined. After choosing the region, the new amounts are assigned to the parameters (e.g. here temperature 278.746 and pressure 100000). Running simulation >blockMesh

In order to assign the amounts which were set in the setFieldDict: >setFields >sonicFoam


The simulation results are as follows:

19


3.50E+02

Velocity

3.00E+02 2.50E+02 2.00E+02 U(m/s) 1.50E+02 1.00E+02

100 cells 1000 cells 10000 cells(10/1)

5.00E+01 0.00E+00 -‐5.00E+01

0

1

2

3

4

5

6

7

8

9

10

arc length(m)

Figure 3.1 Velocity along tube for 10bar/1bar pressure at t=0.007 s 6.00E+02

Velocity

5.00E+02 4.00E+02 U(m/s) 3.00E+02 2.00E+02 10000 cells(100/1) 1.00E+02 0.00E+00 0

1

2

3

4 5 6 arc length(m)

7

8

9

Figure 3.2 Velocity along tube for 100bar/1bar pressure at t=0.007 s

20

10


1.20E+06

Pressure

1.00E+06 8.00E+05 100 cells 1000 cells 10000 cells(10/1)

P(pa) 6.00E+05 4.00E+05 2.00E+05 0.00E+00 0

1

2

3

4 5 6 arc length(m)

7

8

9

10

Figure 3.3 Pressure along tube for 10bar/1bar pressure at t=0.007 s 1.20E+07

Pressure

1.00E+07 8.00E+06 P(pa) 6.00E+06 10000 cells(100/1)

4.00E+06 2.00E+06 0.00E+00 0

1

2

3

4

5

6

7

8

9

arc length(m)

Figure 3.4 Pressure along tube for 100bar/1bar pressure at t=0.007 s

21

10


450

Temperature

400 350 300 250 T(K) 200 150

100 cells

100

1000 cells 10000 cells(10/1)

50 0 0

1

2

3

4

5

6

7

8

9

10

arc length(m)

Figure 3.5 Temperature along tube for 10bar/1bar pressure at t=0.007 s 600

Temperature

500 400 T(K) 300 200 100 10000 cells(100/1) 0 0

1

2

3

4

5

6

7

8

9

10

arc length(m)

Figure 3.6 Temperature along tube for 100bar/1bar pressure at t=0.007 s

22


scalarTransportFoam – shockTube(discretization)1

Simulation Using scalarTransportFoam solver simulate 5s of flow inside a shock tube, with 1D mesh of 1000 cells (10m long geometry from -5 to 5). Patch with a scalar of 1 from -0.5 to 0.5. Simulate following cases: • Set U to uniform (0 0 0). Set diffusion coefficient to a low a medium and a high value. • Set U to (1 0 0) and try 5 in the case of the pure advection including upwind, linear, linearUpwind Gauss linear, QUICK, SuperBee, vanLeer, cubic. Objectives Understanding different schemes of discretization Post processing Follow step-by-step simulation, import your simulation into ParaView, and make a table of execution times.


23


Step by step simulation Log into SSH server Run VNC Compile tutorial Create directory >cd ~/OpenFOAM/OpenFOAM-2.3.0/tutorials/basic/scalarTransportFoam >mkdir shockTube

copy 0, 0.org, constant and system from ~/OpenFOAM/OpenFOAM-2.3.0/tutorials/compressible/sonicFoam/ laminar/shockTube

to created directory. Delete magU, p, thermophysicalProperties, turbulenceProperties. and Replace transportProperties, controlDict, fvSchemes, fvSolution with those from ~/OpenFOAM/OpenFOAM-2.3.0/tutorials/basic/scalarTransportFoam/ pitzDaily


0.org

constant

system

Initial value >ls 0 T U >ls 0.org T

U

After executing the setFields the patched initial files (e.g. T) are overwritten and a list of the values on the mesh cells are written in them for the available mesh, but for running the simulation with a new mesh (e.g. the same mesh with different number of cells) the original file is needed to be patched again. So the original file in the 0.org folder can be used, just copy and replace it by the one in the 0 folder. >cd 0 >cp ../0.org/{T,U} .

Constants >cd .. >ls constants

24


polyMesh transportProperties

The model is a laminar model, and for thermo physical properties an ideal gas is assumed. In the transportProperties file fluid properties such as diffusivity are set. >cd constants >cd polyMesh >ls blockMeshDict

boundary

faces

neighbour

owner

points

>nano blockMeshDict

System >cd ../.. >cd system >ls controlDict

fvScheme

fvSolution

sampleDict

setFieldDict

As it was mentioned before in fvSchemes, the discretization scheme for each operator of the governing equations can be set. >nano fvScheme

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // ddtSchemes { default Euler; } gradSchemes { default } divSchemes { default div(phi,T) }

Gauss linear;

none; Gauss linearUpwind grad(T);

laplacianSchemes { default none; laplacian(DT,T) Gauss linear corrected; } interpolationSchemes { default linear; } snGradSchemes { default }

corrected;

25


fluxRequired { default no; T ; } // ************************************************************************* //

For each type of operation a default scheme can be set (e.g. for divSchemes is set to no default), and also a special type of discretization for each element can be assigned (e.g. div(phi, T) it is set to upwind). For each element, which a discretization method has not been set, the default method will be applied. Note: different schemes should be applied like this: Gauss + scheme In the setFieldDict the amount of 1 is patched to scalar from -0.5 to 0.5. Running simulation >blockMesh >setFields >scalarTransportFoam



Figure 4.1 Just diffusion with low diffusivity (0.00001)

Figure 4.2 Just diffusion with medium diffusivity (0.01)

26


Figure 4.3 Just diffusion with high diffusivity (1)

1.2

Temperature Transport(@t=4s)

1 upwind

0.8

linear linearUpwind Gauss linear

0.6

QUICK

T(K)

0.4

SuperBee vanLeer

0.2

cubic

0 -‐0.2 0

1

2

3

4

5

6

7

8

9

arc length(m)

Figure 4.4 Temperature along tube at t=4 s Scheme cubic linear linearUpwind QUICK superBee upwind vanLeer

Time (s) 17.2 15.1 20.3 21.1 17.5 18.2 17.5

Figure 4.5 Average run time for different schemes

27


scalarTransportFoam – circle(discretization)1

Simulation Using scalarTransportFoam solver simulate a circle (radius=0.5m) at the middle of a 100x100 cell mesh (10m×10m), then move it to the right, to the top and diagonally.

Figure 5.1 Schematic sketch of the problem Objectives Choosing the best scheme of discretization Post processing Follow step-by-step simulation; examine your simulation in ParaView.


28


Step by step simulation Log into SSH server Run VNC Compile tutorial Set up the case such as previous tutorial. Checking files >ls 0

0.org

constant

system

Initial value >ls 0 T

U

>cd 0 >cp ../0.org/{T,U} .

Modify U, appropriately after the correct time, so there will be a velocity field which will move the circle to right, for moving it up and also diagonally do the same. Constants >cd .. >ls constants polyMesh thermophysicalProperties turbulenceProperties

transportProperties

>cd constants >cd polyMesh >ls blockMeshDict

boundary

faces

neighbour

owner

points

>nano blockMeshDict

Modify blockMeshDict for creating a 2D geometry with 100x100 mesh. System >cd ../.. >cd system >ls controlDict

fvSchemes

fvSolution

sampleDict

setFieldDict

Choose the best discretization scheme from previous example and set the fvSchemes.

29


In the setFieldDict patch a circle to the middle of the geometry (Using following lines). // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // defaultFieldValues (volScalarFieldValue T 0 ); regions

( cylinderToCell { p1 ( 0 p2 ( 0 radius fieldValues (

0 -1 ); 0 1 ); 0.5; volScalarFieldValue T 1 ) ; } );

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

cylinderToCell command is used to patch a cylinder to the region, p1 shows the begging of the cylinder centerline and p2 shows its end, in the radius the radius is set. Running simulation >blockMesh >setFields >scalarTransportFoam

For moving the circle to top, and then diagonally just change the velocity field. Note: after moving the circle to the right and changing the velocity field, and continuing the simulation it will be seen that the circle don’t goes up and it continues its way to right. Solution is very simple, just check inside last time step folder, for example: >cd 4 phi

T

U

uniform

There is a phi (flux) file in here, for starting calculations for next time step OpenFOAM® reads the flux from this file. If it isn’t available, phi will be calculated from U file, so easily delete phi and enjoy! Exporting simulation >foamToVTK


30


0 s

1 s

2 s

3 s

4 s

5 s

6 s

7 s

8 s

9 s

10 s

11 s

Figure 5.2 Position of circle at different time steps

31


scalarTransportFoam – pitzDaily1

Simulation Using scalarTransportFoam solver simulate the following simulation for 0.1s. Objectives Understanding the geometry Understanding how a scalar is transported with a field Post processing Follow step-by-step simulation; examine your simulation in ParaView.


32


Step by step simulation Log into SSH server Run VNC Open tutorial >cd OpenFOAM/OpenFOAM-2.3.0/tutorials/basic/scalarTransportFoam/ pitzDaily


0.org

constant

system

Initial value >ls 0 T

U

Note: In this example the velocity field is a preset field, and no need to modify or change it. The solver uses this field to calculate the scalar transportation. Constants >cd .. >ls constants polyMesh transportProperties >cd constants >cd polyMesh >ls blockMeshDict

boundary

faces

neighbour

System >cd ../.. >cd system >ls controlDict

fvScheme

fvSolution

Running simulation >blockMesh >scalarTransportFoam



33

owner

points


0.005 s

0.025 s

0.075 s

0.1 s

Figure 6.1 Velocity field and temperature fields

34


simpleFoam – pitzDaily (turbulence, stationary)1

Simulation Using simpleFoam solver simulate the following stationary simulation with following turbulence models: • laminar • kEpsilon (RAS) • kOmega (RAS) • LRR (RAS) Objectives Understanding turbulence modeling Understanding steady state simulation Post processing Follow step-‐by-‐step simulation; Show the results of U and the turbulent viscosity in two separate contour plots.


35


Step by step simulation Log into SSH server Run VNC Open tutorial >cd ~/OpenFOAM/OpenFOAM-2.3.0/tutorials/incompressible/ simpleFoam/pitzDaily

Checking files >ls 0 constant

system

Initial value >ls 0 epsilon

k

nut

nuTilda

p

R

U

When a turbulent model is chosen, the constants of that model should be set in the appropriate files, for example in kEpsilon model the k and epsilon files should be edited (see introduction). // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions

[0 2 -3 0 0 0 0];

internalField

uniform 14.855;

boundaryField { inlet { type value } outlet { type } upperWall { type value } lowerWall { type value } frontAndBack { type } }

fixedValue; uniform 14.855;

zeroGradient;

epsilonWallFunction; uniform 14.855;

epsilonWallFunction; uniform 14.855;

empty;

36


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

Note: Here is list files which should be available at 0 directory and modified for each turbulent model: • laminar: no constant • kEpsilon (RAS): k and epsilon • kOmega (RAS): k and omega • LRR (RAS): k, epsilon and R • smagorinsky (LES): nuSgs • oneEqEddy (LES): k and nuSgs • spalartAllmaras (LES): nuSgs and nuTilda Constants >cd .. >ls constants polyMesh RASProperties ransportProperties turbulenceProperties

For choosing turbulent model RASProperties file should be checked (e.g. here kEpsilon). // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // RASModel

kEpsilon;

turbulence

on;

printCoeffs

on;

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

Note: for laminar model both turbulenceProperties and RASProperties should be set to laminar, and also in the RASProperties set turbulence and also printCoeffs to off. System >cd ../.. >cd system >ls controlDict

fvScheme

fvSolution

Note: Here because it’s a steady state simulation in controlDict endTime shows number of iterations and deltaT should be 1, because it’s the amount of increase in the iteration number. Running simulation >blockMesh >simpleFoam

37


Note: When the solution converges, “Solution converged!” message will be displayed in the Shell window. If nothing happens and you see no message after awhile, then you should check the residuals which are displayed in the Shell window (you should check initial residual, it shows the difference between this iteration and the last one), if those are close to amounts you have set in the fvSolution then you can stop simulation (ctrl+c). Time = 2634 DILUPBiCG: Solving for Ux, Initial residual = 1.16764e-07, Final residual = 1.16764e-07, No Iterations 0 DILUPBiCG: Solving for Uy, Initial residual = 8.54427e-07, Final residual = 8.54427e-07, No Iterations 0 DICPCG: Solving for p, Initial residual = 2.42693e-06, Final residual = 9.59852e-07, No Iterations 4 time step continuity errors : sum local = 4.69629e-06, global = 3.77562e-07, cumulative = -9.52898e-06 DILUPBiCG: Solving for epsilon, Initial residual = 9.04529e-06, Final residual = 9.04529e-06, No Iterations 0 DILUPBiCG: Solving for k, Initial residual = 9.07688e-06, Final residual = 9.07688e-06, No Iterations 0 ExecutionTime = 3.12 s ClockTime = 4 s



38


Velocity magnitude

Turbulent viscosity

kEpsilon

kOmega

LRR

Figure 7.1 Comparison of different turbulent models at steady state

39


pisoFoam – pitzDaily (turbulence, transient)1

Simulation Using pisoFoam solver simulate the following simulation for 0.1s with following turbulence models: • Laminar • Smagorinsky (LES) • SpalartAllmaras (LES) • oneEqEddy (LES) • kEpsilon (RAS) Objectives Understanding turbulence modeling Understanding the difference between transient and steady state simulation Find appropriate turbulent model Post processing Follow step-‐by-‐step simulation; Show the results of U and the turbulent viscosity in two separate contour plots at three different time steps. Compare with stationary.


40


Step by step simulation Log into SSH server Run VNC Open tutorial >cd ~/OpenFOAM/OpenFOAM-2.3.0/tutorials/incompressible/pisoFoam/ les/pitzDaily


system

Initial value >ls 0 B

k

nSgs

nuTilda

p

U

Set the turbulence model (set the model constants in here). Note: For different turbulent models, different constant files should be modified (check previous tutorial). Constants >cd .. >ls constants polyMesh LESProperties transportProperties turbulenceProperties

As mentioned before in turbulenceProperties we can set the turbulent model type. // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // simulationType

RASModel;

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

For choosing turbulent model LESProperties file should be checked. Note: If RAS models are being used, in the constant directory there is RASProperties file and we should modify it, but if LES models are used LESProperties file should be found and modified. System >cd ../.. >cd system >ls controlDict

fvScheme

fvSolution

41


Running simulation >blockMesh >pisoFoam


The simulation results are as follows: For kEpsilon model after 0.1 s the results are like stationary simulation, so it seems it has reached steady state. Other models don’t have a steady situation and are changing all the time, so they are not appropriate for stationary simulation.

42


Velocity magnitude Smagorinsky

Turbulent viscosity

0.01s

0.05s

0.1s SpalarAllmaras

0.01s

0.05s

0.1s Velocity magnitude

Turbulent viscosity

43


oneEqEddy

0.01s

0.05s

0.1s

kEpsilon

0.01s

0.05s

0.1s

Figure 8.1 Comparison of different turbulent models for transient simulation.

44


interFoam – damBreak (multiphase)1

Simulation Using interFoam solver simulate breaking of a dam for 2s. Objectives Understanding how to set viscosity, surface tension and density for two phases Post processing Follow step-‐by-‐step simulation; see the results in paraview.


45


Step by step simulation Log into SSH server Run VNC Open tutorial >cd ~/OpenFOAM/OpenFOAM-2.3.0/tutorials/multiphase/interFoam/ laminar/damBreak


system

Initial value >ls 0 alpha1

alpha1.org

p_rgh

U

In the alpha1 and p_rgh files the initial value and also boundary conditions for phase alpha1 and also pressure are set. // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions

[0 0 0 0 0 0 0];

internalField

uniform 0;

boundaryField { leftWall { type } rightWall { type } lowerWall { type } atmosphere { type inletValue value } defaultFaces { type

zeroGradient;

zeroGradient;

zeroGradient;

inletOutlet; uniform 0; uniform 0;

empty;

46


}

} // ************************************************************************* //

Note: the inletOutlet and the outletInlet boundary conditions are used when the flow direction is not known. Infact these are derived types and are a combination of two different boundary types. inletOutlet: when the flux direction is toward outside the domain it works like a zeroGradient boundary condition and when when the flux is toward inside it is like a fixedValue boundary condition. outletInlet: it is the other way round, when the flux direction is toward outside the domain it works like a fixedValue boundary condition and when when the flux is toward inside it is like a zeroGradient boundary condition. e.g. if velocity field outlet is set as inletOutlet and inletValue for is set to (0 0 0), it avoids backflow at the outlet! The “inletValue” or “outletValue” are values for fixedValue type of these boundary conditions and “value” is dummy entery for openfoam for finding the variable type (e.g. using (0 0 0), openfoam understands the variable is a vector). Constants >cd .. >ls constants dynamicMeshDict polyMesh transportProperties turbulenceProperties

In the transport properties file properties of two files can be set:

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // phase1 { transportModel Newtonian; nu nu [ 0 2 -1 0 0 0 0 ] 1e-06; rho rho [ 1 -3 0 0 0 0 0 ] 1000; CrossPowerLawCoeffs { nu0 nu0 [ 0 2 -1 0 0 0 0 ] 1e-06; nuInf nuInf [ 0 2 -1 0 0 0 0 ] 1e-06; m m [ 0 0 1 0 0 0 0 ] 1; n n [ 0 0 0 0 0 0 0 ] 0; } BirdCarreauCoeffs { nu0 nuInf k n }

nu0 [ nuInf k [ 0 n [ 0

0 [ 0 0

2 0 1 0

-1 0 0 0 0 ] 0.0142515; 2 -1 0 0 0 0 ] 1e-06; 0 0 0 0 ] 99.6; 0 0 0 0 ] 0.1003;

} phase2 { transportModel Newtonian; nu nu [ 0 2 -1 0 0 0 0 ] 1.48e-05; rho rho [ 1 -3 0 0 0 0 0 ] 1; CrossPowerLawCoeffs

47


{ nu0 nuInf m n

nu0 [ nuInf m [ 0 n [ 0

0 [ 0 0

2 0 1 0

-1 0 0 0 0 ] 1e-06; 2 -1 0 0 0 0 ] 1e-06; 0 0 0 0 ] 1; 0 0 0 0 ] 0;

nu0 [ nuInf k [ 0 n [ 0

0 [ 0 0

2 0 1 0

-1 0 0 0 0 ] 0.0142515; 2 -1 0 0 0 0 ] 1e-06; 0 0 0 0 ] 99.6; 0 0 0 0 ] 0.1003;

} BirdCarreauCoeffs { nu0 nuInf k n } } sigma

sigma [ 1 0 -2 0 0 0 0 ] 0.07;

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

The model which is used here is Newtonian, so just constants for nu and rho are important; if some other model is used then the constants in the CrossPowerLawCoeffs and BirdCarreauCoeffs also should be modified. In the last line sigma should be set, sigma is the surface tension between two phases, for this example it is surface tension between air and water. System >cd ../.. >cd system >ls controlDict decomposeParDict fvScheme fvSolution setFieldsDict

Running simulation >blockMesh >setFields >interFoam


The simulation results are as follows (these are not results for original mesh, it is a refined mesh):

48


Figure 9.1 Two phase flow simulation.

49


compressibleInterFoam – depthCharge3D1

Simulation Using compressibleInterFoam solver simulate the following simulation for 0.5s. Objectives Understanding the difference between incompressible and compressible solvers Understanding parallel processing Post processing Follow step-‐by-‐step simulation; see the results in paraview.


50


Step by step simulation Log into SSH server Run VNC Open tutorial >cd ~/OpenFOAM/OpenFOAM-2.3.0/tutorials/multiphase/ compressibleInterFoam/laminar/depthCharge3D

Checking files >ls 0 Allclean Allrun constant system

Initial value >ls 0 alpha.water.org p_rgh.org p.org T.org U

Constants >cd .. >ls constants g LESProperties polyMesh thermophysicalProperties thermophysicalProperties.air thermophysicalProperties.water turbulenceProperties

System >cd ../.. >cd system >ls controlDict decomposeParDict fvScheme fvSolution setFieldsDict

In the decomposeParDict number of domains and also how the domain is going to be devided to this subdomains for parallel processing are set. // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // numberOfSubdomains 4; method

hierarchical;

simpleCoeffs { n delta } hierarchicalCoeffs { n

( 1 4 1 ); 0.001;

( 1 4 1 );

51

delta order


0.001; xyz;

} metisCoeffs { } manualCoeffs { dataFile }

"";

distributed

no;

roots

( );

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

For example in here the numberOfSubdomains is set to 4, so there will be 4 sub domains. In the method, the method for dividing domain is chose. n shows number of subdomains in every direction, for example in here there will be 4 sub domains in the y direction. In the setFieldsDict a sphere and also a cube are patched to the domain. Running simulation >blockMesh >setFields

For running the simulation in parallel mode first the computing domain should be divided into subdomains and a processor should be assigned to each subdomain, these processes are done by following command: >decomposePar

After executing this command four new directories will be made in the simulation directory (processor0, processor1, processor2 processor3), and each subdomain calculation will be saved in the respective processor directory. Note: when the domain is divided to subdomains in parallel processing new boundaries are defined, which shows at that boundary the data should be exchanged by the neighbor boundary, which it’s connected to, in the main domain. >mpirun –np 4 compressibleInterFoam –parallel >log.CompressibleInterFoam

The command for running the parallel simulation is long, instead of executing this command; a text file can be easily created and executed, which also run the same command. Use nano editor create a text file (e.g. go): > nano go

add the command to this file: mpirun –np 4 compressibleInterFoam –parallel >log.CompressibleInt

52


-erFoam

exit editor and save the file (ctrl+x , y, enter). For changing this file to an executable file, file permissions should be edited: >ls –la go

by using this command file permissions are displayed: -rw-r--r-- 1 openFoamUser E020D166 73 Aug 23 9:15 go

the first r shows that this text file can be read by user, the w shows that user has the permission to write this file, but the – sign shows that this file is not executable by the user, for changing this permission execute following command: >chmod u+x go

now this file is executable: >ls –la go -rwxr--r-- 1 openFoamUser E020D166 73 Aug 23 9:15 go

Now you can run the simulation by this executable text file: >./go

Note: after running the simulation following errors will occur, just ignore them!

[[email protected]:28107] mca: base: component_find: unable to open /home/openFoamUser/OpenFOAM/ThirdParty-2.3.0/platforms/linux64Gcc/openmpi1.4.1/ lib/openmpi/mca_btl_openib: perhaps a missing symbol, or compiled for a different version of Open MPI? (ignored)

The information about each step will be written to log.compressibleInterFoam. For checking the last information which is written to this file following command can be used during the simulation run: >tail –f

log.compressibleInterFoam

Note: for running a simulation in parallel mode first number of free CPU cores on the server should be examined. And for running simulation on n cores there should be at least 2n free cores (in order not to interrupt the other users). For checking number of free cores the following command can be used: >top

It will show the processes on the server on which you are calculating (You know number of cores on your calculating server, e.g. current server has 32 cores, and by finding number of simulations, which are running on the server and subtracting them, number of free cores can be easily found!). Note: top command execution can be interrupted by: q (or ctrl+c)

53


Exporting simulation For exporting data for post processing, at first all the processors data should be put together, and a unit directory for each time step being created. By executing following command all the processors data will be combined and new directories for each time step will be created in the simulation main directory: >reconstructPar >foamToVTK

Note: to continue the reconstructing or foamToVTK conversion from a special time the following flag can be used: >reconstructPar –time [time name, e.g. 016]:


0s

0.01 s

0.025 s

0.055 s

0.07 s

0.09 s

0.1 s

0.11 s

0.15 s

0.16 s

0.185 s

0.205 s

0.235 s

0.26 s

0.295 s

0.325 s

0.36 s

0.44 s

0.47 s

0.5 s

Figure 10.1 3D depth charge, parallel simulation

54


compressibleInterFoam – depthCharge3D (manual discritization)1 Simulation Using compressibleInterFoam solver simulate the following simulation for 0.5s. Objectives More detailed investigation of parallel processing Using manual method in parallel processing Post processing Follow step-by-step simulation; examine your simulation in ParaView.

1 Provided by Vikram Natarajan

55


Step by step simulation Log into SSH server Run VNC Open tutorial >cd ~/OpenFOAM/OpenFOAM-1.7.1/tutorials/multiphase/ compressibleInterFoam/les/depthCharge3D

Edit the decomposeParDict file, depending on the required specifications. Note: This applies if using simple, hierarchical or scotch methods. The manual method is dealt with later in this tutorial. The decomposeParDict file is shown below: // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // numberOfSubdomains 4; method scotch; simpleCoeffs { n ( 1 4 1 ); delta 0.001; } hierarchicalCoeffs { n ( 1 4 1 ); delta 0.001; order xyz; } manualCoeffs { dataFile ""; } distributed no; roots ( ); // ************************************************************************* //

numberOfSubdomains should show the number of processors used, and method should show the method to be used. In the above example, the case is simulated with the Scotch method and 4 processors. If the simple method is being used, the parameter n must be changed accordingly. The three numbers 1 4 1 indicate the numbe of pieces the mesh is split into in the x, y and z directions respectively. If the hierarchical method is being used, these parameters and also the order in which the mesh should be split up in each direction should be provided. If the scotch method is being used, then no user-supplied parameters are necessary. There is also a parameter delta, known as the cell skew factor. This factor is set to a default value of 0.001, and meaures to what extent skewed cells should be accounted for. In order to check the quality of the mesh, run the following command: >checkMesh -allGeometry

If the message “Mesh OK” is seen – the mesh is fine, and no corrections need to be made. If the mesh fails one or more tests, then a possible course of action is to increase the delta

56


parameter (for example: to 0.01) and then rerun the blockMesh and checkMesh allGeometry commands. Running simulation >cp alpha.water.og alpha.water >cp p.org p >cp T.org T

These commands clone the original files to be used in the simulation, so that they are not altered. Note: Before running any simulation, it is important to run the top command (type the top command in the terminal), to check the number of cores being used on the server. Check the load average, this is on the first line and shows the average number of cores being used. There are three numbers displayed, showing the load averages across three different time scales (one, five and 15 minute respectively). Add the number of cores you plan to use to this number – and you will get the expected load average durirng your simulation. This number should be less than the total number of cores in the server – or the simulation will be slowed causing an inaccurate result. It is recommended to leave at least 6-8 cores free, to allow for any fluctuations in the serve load. >blockMesh

This commands forms a mesh using the instructions contained in the blockMeshDict file. If an error: blockMesh: command not found is thrown, the bashrc file has probably not been sourced. Run the command: . ~/OpenFOAM/OpenFOAM-2.3.0/etc/bashrc and this should fix the problem. >setFields

If this command raises an error, the most likely explanation is that the number of cells in the mesh does not match up with the values in the setFieldDict, this is especially true if the values have been changed from the default 80 160 80. >decomposePar

This decomposes the mesh according to the supplied instructions. If an error is raised, it is possible that the product of the parameters in n do not match up to the number of processors, this goes for the simple and hierarchical methods >mpirun -np compressibleInterFoam parallel > log

is the number of processors being used. For example, if 4 processors are desired, we would run: >mpirun -np 4 compressibleInterFoam -parallel > log

The simulation can take several hours, depending on the size of the mesh, the VNC desktop can be closed during this time, but it is important to NOT close the terminal

57


window, otherwise all progress will be lost.Once the simulation is done, run the following commands: >reconstructPar

This reconstructs the complete solution from the individual processor directories >foamToVTK

This converts the data into VTK files for visualization >paraview &

The solution can now be viewed in paraview. Manual method: The manual method is slightly different from the other three. In order to use it: Set the decomposeParDict file as any other simulation. For decomposition method, choose either simple, hierarchical or scotch. Set the number of processors to the same number which is going to be used for manual. >decomposePar – cellDist

Once the decomposition is done, check the cellDecomposition file in the constant directory. It should have a format similar to: // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // 1024000 ( 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 1 ...) // ************************************************************************* //

58


Note: If the above output is not displayed, but a stream of NUL characters, your text editor is probably printing to binary. To fix this, open system/controlDict, and change the writeFormat field from: writeFormat

binary;

to

writeFormat

ascii;

The first number n after the header, but before the opening brackets, 1024000 in this example, refers to the number of points in the mesh. Within the brackets, n lines follow. Each line contains one number between 0 and n-1, where n is the number of processors to be used for the computation. This number refers to the processor that is to be used to compute the corresponding cell in the points file in the constant directory. For example, if the second line in the points file brackets reads 0.125 0 0 and the second line in the cellDecomposition directoy reads 0, this means that the cell 0.125 0 0 will be processed by processor 0. This cellDecomposition file can now be edited. Although this can be done manually, it is probably not feasible for any sufficiently large mesh. The process must thus be automated, by writing a script to populate the cellDecomposition file according to the desired processor breakdown. When the new file is ready, save it under a different name: >cp cellDecomposition manFile

Now, edit the decomposeParDict file. Select decomposition method manual, and for the dataFile field in the manual coeffs range, specify the path to the file which contains the manual decomposition. Note that OpenFOAM searches in the constant directory by default, in case relative paths are being used: // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // numberOfSubdomains 16; method

manual;

simpleCoeffs { n delta } hierarchicalCoeffs { n delta order } manualCoeffs { dataFile

( 1 4 1 ); 0.001;

( 1 16 1 ); 0.001; xyz;

"manFile";

59


} distributed

no;

roots

( );

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

Run the simulation as normal. Visualiziing the processor breakdown: It may be interesting to visualize how exactly OpenFOAM breaks down the mesh. This can be easily visualized using paraview. After running the simulation, but before running the reconstructPar command, repeat the following for each of the processor directories: >cd processor

where n is the processor number >foamToVTK

convert the individual processor files to VTK, next, open paraview: >paraview &

For each of the processor directories, perform the following steps: - Open the VTK files in the relevant processor directoy - Double click them to open them and click on “Apply” - The part of the mesh decomposed by that processor will appear, in grey. - Change the colour in the drop-down menus in the toolbar. This is to ensure that each individual part can be easily seen Once this is done for all processors, the entire mesh will appear, however, the processor regions can now easily be seen in a different colour. In order to save this, there are two options. The first option is to take a screenshot: File > Save a screenshot

The second option is to save the state as a paraview state file. File > Save State

The current position can then be easily recovered by: File > Load State

Saving state allows changes to be made afterwards, while saving a screenshot just keeps a picture, while losing the ability to make changes after exiting paraview. Doing both is recommended.

60


simpleFoam & scalarTransportFoam – TJunction (residence time distribution)1 Simulation Using simpleFoam and scalarTransportFoam solvers simulate flow through a square cross section T pipe with following values, and measure RTD (residence time distribution) for both inlets using a step function injection: • Pipes cross section: 1×1 m2 • Pipes length: 3 m • Gas in the system: Air at ambient • Oprating pressure: 10e5 pa • Inlet 1: 0.1 m/s • Inlet 2: 0.2 m/s Objectives Understanding RTD calculation using OpenFOAM Post processing Follow step-‐by-‐step simulation; Plot step response function, and also RTD curve.


61


Step by step simulation Log into SSH server Run VNC Open tutorial Copy following tutorial to your working directory: > cp -r ~/OpenFOAM/OpenFOAM-2.3.0/tutorials/incompressible/ simpleFoam/pitzDaily/{0,constant,system} .


system

Creating mesh Edit the bockMeshDict as following for creating appropriate geometry. >nano constant/polyMesh/blockMeshDict

// * * * * * * * * * * convertToMeters 1.0; vertices ( (0 4 0) // 0 (0 3 0) // 1 (3 3 0) // 2 (3 0 0) // 3 (4 0 0) // 4 (4 3 0) // 5 (7 3 0) // 6 (7 4 0) // 7 (4 4 0) // 8 (3 4 0) // 9 (0 4 1) // 10 (0 3 1) // 11 (3 3 1) // 12 (3 0 1) // 13 (4 0 1) // 14 (4 3 1) // 15 (7 3 1) // 16 (7 4 1) // 17 (4 4 1) // 18 (3 4 1) // 19 ); blocks ( hex (0 1 2 9 10 11 hex (9 2 5 8 19 12 hex (8 5 6 7 18 15 hex (2 3 4 5 12 13 ); edges (

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

12 15 16 14

19) 18) 17) 15)

(10 (10 (10 (30

30 10 30 10

10) 10) 10) 10)

simpleGrading simpleGrading simpleGrading simpleGrading

62

(1 (1 (1 (1

1 1 1 1

1) 1) 1) 1)


); patches ( patch inlet_one ( (0 10 11 1) ) patch inlet_two ( (7 6 16 17) ) patch outlet ( (4 3 13 14) ) wall walls ( (0 1 2 9) (2 5 8 9) (5 6 7 8) (2 3 4 5) (10 19 12 11) (19 18 15 12) (18 17 16 15) (15 14 13 12) (0 9 19 10) (9 8 18 19) (8 7 17 18) (2 1 11 12) (3 2 12 13) (5 4 14 15) (6 5 15 16) ) ); mergePatchPairs ( ); // ************************************************************************* //

Initial value >ls 0 epsilon

k

nut

nuTilda

p

R U

Update p, U, nut, nuTilda, k and epsilon files with new boundary conditions. >nano 0/U // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions

[0 1 -1 0 0 0 0];

internalField

uniform (0 0 0);

boundaryField { inlet_one { type

fixedValue;

63

value } inlet_two { type value } outlet { type } walls { type value }


uniform (0.1 0 0)

fixedValue; uniform (-0.2 0 0)

zeroGradient;

fixedValue; uniform (0 0 0)

} // ************************************************************************* //

Constants >ls constants polyMesh RASProperties transportProperties turbulenceProperties

Check RASProperties file for turbulence model (kEpsilon). // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // RASModel

kEpsilon;

turbulence

on;

printCoeffs

on;

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

System >ls system controlDict

fvScheme

fvSolution

Running simulation >blockMesh >simpleFoam

Wait for simulation to converge, after convergency check the results to get sure the solution is converged (?). >foamToVTK


64


Figure 12.1 Simulation results after convergence (114 iterations)

RTD calculation Use the velocity field from last part of simulation to calculate RTD for this geometry. Copy tutorial >cp -r ~/OpenFOAM/OpenFOAM-2.3.0/tutorials/basic/scalarTransportFoam/pitzDaily/{0,constant,system} .


constant

system

Initial value >ls 0 T U

Update T file boundary conditions to match new simulation boundaries, to calculate RTD of the inlet_one, set the internalField value to 0, T value for inlet_one 1.0 and T value for inlet_two 0. no need to modify U! Note: Replace the velocity field with calculated velocity field from first part of tutorial (use the field from last time step), and no need to modify or change it. The solver uses this field to calculate the scalar transportation. Constants >ls constants polyMesh transportProperties >ls constants/polyMesh blockMeshDict

Replace the blockMeshDict file with the one from first part of tutorial.

65


System >ls system controlDict

fvScheme

fvSolution

In the controlDict file change the endTime from 0.1 to 120 (approximately two times ideal resistance time) and also delatT from 0.0001 to 0.1 (Courant number approximately 0.4). Running simulation >blockMesh >scalarTransportFoam >foamToVTK

Simulation results The simulation results are as follows:

Figure 12.2 Simulation result at 120s

Calculating RTD For calculating RTD at first average T concentration at outlet should be calculated. For this purpose integrate variables function of paraview can be used. >foamToVTK

Load the outlet VTK file into paraview using following path: File > Open > VTK > outlet > outlet_..vtk > OK > Apply

Select T from variables menu, then integrate the variables on the outlet: Filters > Data Analysis > Integrate Variables > Apply

66


The values which are given in the opened window are inegrated values in that specific time step, by changing time step values for different time steps are displayed. As mentioned before average value of the property is needed so this values should be devided by outlet area to get average values (1m × 1m). The same procedure should be followed for calculating RTD of inlet_two, except T value for inlet_one should be 0 and for inlet_two it should be 1.0. Calculating RTD 0.8 0.7 0.6 0.5 0.4 0.3 0.2 inlet_one

0.1

inlet_two

0 0

20

40

60

80

100

120

Figure 12.3 Average value of T in the outlet for two inlets versus time

0.06 0.05

RTD1

0.04

RTD2

0.03 0.02 0.01 0 0

20

40

60

80

Figure 12.4 RTD of two inlets

67

100

120


Appendix A Important commands in Linux (Mostly on Unix [IRIX, Alpha Unix… usable]) cat, more, less, vi, vim

cd, cd ..

File viewer with pure read function - in order of ease of operation. In less with pagedown/pageup you can navigate within the file, with / and ? can look for strings, q can be used for closing less. cat is back for universally available on Unix. Changing the directory, cd .. goes one directory up and cd ~ moves to home directory. Important to note is the space between cd and .. as opposed to DOS!

cp, cp -r

Copying files or entire directory trees with -r option. Caution: There is no prompt when overwriting existing files! The important thing is that always a goal has to be given, at least one "." which means, copy the current directory.

Ctrl+r

Reverse search, for searching an already typed command in a terminal window.

du, du -s,

Specifies the amount of space in a directory. For safety reasons you should use the -k option (output in kilobytes), since some systems provide the space in blocks that include only 512 bytes ...

du -k gedit

Text editor with graphical user interface, when working with gedit some temporary files (originalFileName~) are created, they can be deleted after saving.

grep

Search command for plain-text data sets for lines matching a regular expression.

gzip, gunzip

Compression-/decompression program for individual files (as opposed to zip/unzip, this can also edit directory or file lists). The great advantage of gzip: Fluent® and OpenFOAM® are able to read and write gz files directly, which saves about 30-90% space.

kill, kill -9

Stopping processes. For this the process ID is required, which can be

68


found out with top or ps. The "Exit" is irrevocable course - but you cannot shoot processes, if you are not the "owner". ls, ls –la

Lists the contents of a directory, with option -la also "hidden" files are displayed, also the file size and characteristics.

mc

(Navigation in the text window), esc-keys, may be necessary: mc -c, for navigating through mc use function keys or esc+[number] combination, e.g. F9 or esc+9 for moving to the menus at the top

mkdir mv

Creates a new directory Moving files and directories. Caution: There is no prompt when overwriting existing files!

nano

The command to run the nano text editor, a terminal based text editor.

passwd

The command to change the password.

ps, ps –A

Lists all the processes that were started in the respective command window with the options are all running processes on the system display.

ps waux

pwd

Specifies the current working directory.

rm, CAUTION: rm -fr

Deletes files. The option -r also be emptied directories recursively and delete directories, f ("force") prevents any further inquiry. - Incorrectly applied, this command can lead to irreversible loss of all (private) data.

rmdir

Deletes an empty directory

scp

The copy command over the network - as secure FTP replacement. Also dominates the -r (recursive) option. Usage: scp source file destination file with source and the destination format can be USERNAME@ COMPUTER.DOMAIN:PATH/TO/FILE. Source or target can of course also be created locally, then your user name and computer are not required.

69


ssh

Telnet replacement with encryption. On Windows, for example, implemented with putty.

tail, tail -f

File viewer, the default outputs the last 10 lines of a file. With option -n XX can spend the last XX lines, with the -f option, the command is running from those lines, which are attached to a file. The command is therefore perfect for watching log files.

top

Displays a constantly updated list of all running processes, with process ID, memory and CPU usage (for processes of one user “top u openFoamUser” should be used, and for quitting “q” or “ctrl+c” should be applied).

exit

Closing connection (termianl window).

70


Appendix B Frequently Asked Questions (FAQ) Q-

What should I do, in case of GAMBIT failure?

A-

e.g. Program stops responding: • In the command window type "ps", search for of Gambit process number. • "kill -9 PROCESS NUMBER" Enter Gambit creates lock files, which must also be deleted, otherwise no possible opening of the affected files: • "rm *. Lok" Enter Furthermore, "junk" (temporary files from GAMBIT) should be disposed of: • "rm -fr GAMBIT.xxx" erases the complete directory, xxx again is the process number. • If you have forgotten, however, before the crash to save, you should copy the file "jou" (it contains all the commands that have been executed and can be processed automatically in GAMBIT) from the directory, to resume its status before the crash.

Q-

How can I prevent typing long commands in the terminal for couple of times?

A-

By using reverse search, use ctrl+r to search for previous commands typed in the terminal, e.g. typing a part of command show the suggestions and you can navigate through them.

Q-

My VNC is not responding from server side?

A-

First you should kill your VNC server: • vncserver –kill :[YOUR DISPLAY NUMBER] Restart your VNC server: • vncserver:[YOUR DISPLAY NUMBER] -geometry 1600x800 -depth 24

Q-

I’ve deleted some of my files accidently, What should I do?

A-

Sorry, no recycling or undelete in Linux

Q-

I cannot connect to the server?

A-

Check to see if you have an IP address.

Q-

How can I start vnc from my linux computer terminal?

A-

Use command: • vncviewer :[NUMBER OF LOCAL PORT, e.g. 1 or 2]

71

QA-


foamToVTK command don’t work for chtMultiRegionFoam? Use command: • foamToVTK –region[REGION NAME]

QA-

Is it possible to export animations from paraview? Yes, by choosing .ogv file format from “file/save animation” menu. Output will be a video file with .ogv format.

QA-

Is there any tool in linux to convert series of paraview pictures to video? Yes, command line tool ffmpeg: • ffmpeg –r [FRAME PER SECOND RATE] –f image2 –i [images names, e.g. rho.%4d.jpg] [OUTPUT FILE NAME].[OUTPUT FILE FORMAT, e.g avi]

72


Appendix C paraView The visualization application, which is usually used with OpenFOAM® is paraView (Figure C-1) which is a free, open source program. The OpenFOAM® command, foamToVTK, converts OpenFOAM® files to readable formats for paraView.

Figure C-1 The paraView window The tree structure of paraView helps user to easily choose and display suitable submodels for creating the desired image or animation. Adding mesh or velocity vectors to a contour plot of pressure is an example of this functionality. For general operations a selection should be made and then the green Apply bottom should be pressed, the reset bottom is used for resetting the window and delete, deletes the selected operation.

Properties panel Setting for time step, regions and fields can be done in the Properties panel.

73


Figure C-2 The Properties panel for contour plots

Display panel For a given case settings for visualizing the data are in the Display panel. Some important notes: • The max/min data range might not be updated automatically, so check and if needed, rescale the data range after appropriate intervals (e.g. after loading the case). •

• •

Two panels can be accessed by clicking Edit Color Map button: 1. Color Scale panel: scale colors can be choose, for resetting the color to standard blue to red, click choose preset, and from opened window select, Blue to Red HSV. 2. Color Legend panel: legend layout (e.g. font) can be set in this panel. For displaying the mesh select Wireframe from Representation menu of the Style panel. Single color can be used for visualizing the geometry, e.g. a mesh (if Wireframe is selected), by selecting Solid Color from the Color By menu and specifying the color in the Set Ambient Color window. 74


•

The opacity of the image can be set (1 = solid, 0 = invisible) in the Opacity in the Style panel.

Figure C-3 The Display panel

Button toolbars Pull-down menus at the top of the main window and the major panels, in the toolbars below the main pull-down menus increase the functionality of paraView. The function of each button can be easily understood by its icon (Figure C-4), also any button description can be found in the Help menu.

Figure C-4 Toolbars in ParaView

75


Manipulating the view View settings The View menu from Edit menu, contain three items: General, Lights and Annotation. The General panel includes the following items (which are often set at startup): • The background color, from arrow down Choose Color button. • Parallel projection is the usual choice for 2D, CFD simulations. The lighting controls are in Lights panel in the Light Kit panel. For producing images with strong bright colors (e.g. isosurface) Headlight of strength 1 is appropriate. For including annotations in the image Annotation panel should be used. The Orientation Axes feature controls an axes icon in the image window (e.g. to set the color of the axes labels x, y and z). General settings Some default behavior of ParaView can be controlled in the General panel. The Auto Accept button enables, accepting the changes without pressing the Apply button (not a very good option for big cases because re-rendering the image after each change takes lots of time) The Render View panel contains 3 sub-items: General, Camera and Server. The level of detail (LOD) is included in the General panel which controls the rendering of the image while it is being manipulated (e.g. rotated or resized); lower levels, allows cases with large numbers of cells to be re-rendered quickly during manipulation. The Camera panel includes control settings for 3D and 2D movements. User can edit the rotation, translate and zoom map to suit for him, these can be used by a combination of mouse, Shift and Control keys.

Contour plots Selecting Contour from the Filter menu at the top menu bar creates a contour plot. The filter acts on a given module so that, if the module is the If the case is a 3D case module itself, the contours will be a set of 2D surfaces that represent a constant value. There is Isosurface list in the Properties panel that the user can edit by New Range window in the most convenient way. The chosen scalar field is selected from a pull down menu.

Introducing a cutting plane Creating contour plots across a plane is more convenient than isosurfaces; Cutting plates are the tools, which can be used for this purpose, to create surfaces. This can be done by Slice filters, using cutting Plane, Box or Sphere options, a cutting plate can be manipulated like other by mouse. In a similar way contour lines can also be derived out of planes.

76


Vector plots The Glyph filter is used for creating vector plots. Scale Mode menu in the properties panel is used for: setting the length of a vector, weather to be proportional to vector magnitude or not, all with the same length (Vector), or controlling the base length of the glyphs (Set Scale Factor).

Plotting at cell centers Vectors at cell centers are more common, for this purpose the Cell Centers filter should be applied before Glyph filter.

Streamlines Creating tracer lines, using the Stream Tracer filter create streamlines. Tracer points can be along a line or point which are shown in white, and can be choose from Tracer Seed panel. Usually some try and error needed for achieving the desired streamlines, the length of steps tracer takes can be changed in the main Stream Tracer panel, and a smaller length increase calculation time but increase smoothness. For having high quality images Tubes filter can be used after tracer lines have been created. There are also different types for Tubes not just only cylindrical.

Image output For creating a screenshot of the graphs you have made the easiest way is Save Screenshot from File menu. After selecting it in the opened window the picture resolution can be set, and by locking aspect ratio, changing image resolution in one direction cause change in its resolution in the other direction respectively. For high quality images resolution over 1000 pixels is a good choice.

Animation output Some animations can be also created using paraView, in the File menu by choosing animation and setting the resolution, and also frames per time step. You can save your animation by assigning a name and choosing file format, then the captured pictures are saved with this format: “_.”

77


Appendix D OpenFOAM® solvers Here is a list of important solvers, with the Navier-Stokes equation in these solvers (see introduction): •

scalarTransportFoam

Solves a transport equation for a passive scalar. solve (

);

•

fvm::ddt(T) + fvm::div(phi, T) - fvm::laplacian(DT, T)

icoFoam

Transient solver for incompressible and laminar flow fvVectorMatrix UEqn ( fvm::ddt(U) + fvm::div(phi, U) - fvm::laplacian(nu, U) ); solve(UEqn == -fvc::grad(p));

•

simpleFoam

Steady-state solver for incompressible and turbulent flow tmp UEqn ( fvm::div(phi, U) + turbulence->divDevReff(U) ); … ( UEqn() == -fvc::grad(p) )

•

pisoFoam

Transient solver for incompressible and turbulent flow

78


fvVectorMatrix UEqn ( fvm::ddt(U) + fvm::div(phi, U) + turbulence->divDevReff(U) ); UEqn.relax(); if (momentumPredictor) { solve(UEqn == -fvc::grad(p)); }

•

sonicFoam

Transient solver for sonic flow of a compressible gas. fvVectorMatrix UEqn ( fvm::ddt(rho, U) + fvm::div(phi, U) + turbulence->divDevRhoReff(U) ); solve(UEqn == -fvc::grad(p));

•

interFoam

Two phase incompressible flow, using the VOF method. fvVectorMatrix UEqn ( fvm::ddt(rho, U) + fvm::div(rhoPhi, U) - fvm::laplacian(muEff, U) - (fvc::grad(U) & fvc::grad(muEff)) ); UEqn.relax(); if (momentumPredictor) { solve ( UEqn == fvc::reconstruct ( (

79


}

fvc::interpolate(interface.sigmaK())*fvc::snGrad(alpha1) - ghf*fvc::snGrad(rho) - fvc::snGrad(p_rgh) ) * mesh.magSf() ) );

Note: ddt shows time derivatives, div shows divergence terms, laplacian shows laplacian terms and grad shows gradient terms

80

2014 OFoam Tut Complete

Recommend Documents