{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# CFD Simulation Case Setup\n", "\n", "**The created CFD domain is now read into the CFD package of interest to setup the CFD simulation. It should be noted that the current tutorial has a significant difference compared to other available CFD tutorials online! This tutorial is structured and developed based on a generic and methodological approach to set up a CFD simulation. The first principals and reasonings for each setting is discussed at each step. Potential alterations and modifications to perform similar analysis for different flow conditions are also addressed and discussed. Hence, in the end user will have the capability of applying potential modifications, improvements or extending the application of the current CFD simulation to a more complex problem of interest, rather than having a one time successful run of a specific simulation with specific and strictly pre-defined boundary conditions.**\n", "\n", "> **_In simple words: Current tutorial teaches users to fish, rather than giving them a fish._**\n", "\n", "## Setting up a CFD simulation has following four steps:\n", "\n", "**0. Import the CFD Domain into Solver:**\n", "\n", "The previously generated CFD domain, fully explained in the Domain section, will be imported properly into the CFD solver.\n", "\n", "**1. Setup Model/s:**\n", "\n", "According to the physics of the flow field user will select required model/s to simulate the flow.\n", "\n", "**2. Setup Working Fluid/s & Solid/s:**\n", "\n", "User will define the physical and thermodynamical properties of the working fluid/s and solid/s in the problem. \n", "\n", "**3. Setup Boundary & Zone Conditions:**\n", "\n", "Solving the governing equations of the flow (i.e. system of partial differential equations) requires well-defined boundary conditions within the CFD domain. These conditions are selected and defined in this step.\n", "\n", "**4. Setup Solution Methods:**\n", "\n", "In CFD simulations the governing equations of the flow are solve numerically. Based on the physics of the problem appropriate numerical schemes and solution methods are selected at this step.\n", "\n", "In the following section the details for the above-mentioned four steps of the CFD simulation setup process for the case study of **2D Laminar Flow Over a Cylinder** are explained in great details. As a general guideline, to setup a CFD simulation in OpenFoam it is recommended to start with files from one of the provided default tutorials in OpenFoam located in `OpenFoam/run/tutorials`. Users should find appropriate tutorial such that the general physics of the chosen tutorial is close to the problem of interest of users. This way, the majority of the required physical settings are pre-set and the others can be changed/modified accordingly. For example for the problem of **2D Laminar Flow Over a Cylinder** the tutorial of **2D over an airfoil** located at `OpenFoam/run/tutorials/incompressible/simpleFoam/airFoil2D` is used. Using the previously developed file, the above-mentioned four steps for CFD simulation development is done via editing the corresponding files, where the settings are defined and are located in the main working directory.\n", "\n", "**0. Importing the Created CFD Domain into Solver:**\n", "\n", "The official mesh generator for OpenFoam is called blockMesh. This is a command base mesher software. Users write/edit the desired settings to generate and discretize the geometry of CFD domain in a dictionary file (i.e. a text file) named blockMeshDict located at \\constant\\\\polyMesh inside the main working directory. In simple words the file is the dictionary, where the geometry and discretization process of it are fully defined. Once edits of `blockMeshDict` file are completed then it is compiled via the command blockMesh. This will generate multiple required files, that include the details of the discretized CFD domain for the use of OpenFoam solver in the directory of the `blockMeshDict` file. Readers can read more about the `blockMesh` utility [here](http://cfd.direct/openfoam/user-guide/blockmesh/).\n", "\n", "Alternatively, OpenFoam has various built-in commands that allow it's users to convert their previously generated CFD domain in different mesher packages into the readable format for the OpenFoam. A list of these commands can be found [here](http://www.openfoam.org/features/mesh-conversion.php). In this tutorial this alternative approach is used to convert the previously generated CFD domain with * .msh format to the readable format in OpenFoam. In order to perform this process users put the * .msh file of the previously generated CFD domain in the main working directory and use the command fluentMeshToFoam with the name of the mesh file after the command to perform the mesh conversion process (i.e. fluentMeshtoFoam [file name]. msh). This will generate multiple required files, that include the details of the discretized CFD domain for the use of OpenFoam in the working directory. Users should make sure that the generated report contains no error or problem. The details of the previously generated CFD domain in * .msh format are fully explained in the Domain (see Common_Files -> Domain) part of this tutorial.\n", "\n", "**1. Setup Model/s:**\n", "\n", "The main assumption in the classical problem of a two dimensional flow over a cylinder is that the free stream flow is uniform and steady. Therefore, the steady model is chosen to develop this CFD simulation. This along with other numerical schemes settings are set in a dictionary file called \"fvSchemes\" located at `\\system` folder. To set the solver to be steady users define `Steadystate` for the `ddtSchemes` and the rest of the settings will remain unchanged:\n", "\n", "```C++\n", "ddtSchemes\n", "{\n", " default Steadystate;\n", "}\n", "```\n", "In cases where the flow properties rapidly changes with respect to time, such as cases where the Reynolds number, based on cylinder diameter (as the length scale), is greater than the critical value of 100 in which unsteady vortex shedding will begin to happen around the cylinder the appropriate Transient model should be chosen in this step. Various available time schemes along with other implementable option within \"fvSchemes\" dictionary file in OpenFoam can be found [here](http://cfd.direct/openfoam/user-guide/fvschemes/).\n", "\n", "In the current problem the focus is to investigate laminar wake and boundary layer behavior when a uniform flow passes over a cylinder. The appropriate model for this problem is set in a dictionary file called `RASPropertise` (RAS stands for Reynolds Averages Stress) located at `\\constant` folder:\n", "\n", "```C++\n", "RASModel laminar;\n", "\n", "turbulence off;\n", "\n", "printCoeffs off;\n", "```\n", "\n", "It is important to note that the critical Reynolds number for turbulent flows over a cylinder is about $10^5$. In cases that the Reynolds number, based on cylinder diameter (as the length scale), is greater than this critical value the appropriate turbulent model should be selected at this step. Various available turbulence models along with other options within `RASPropertise` dictionary file in OpenFoam can be found [here](http://cfd.direct/openfoam/user-guide/turbulence/).\n", "\n", "**2. Setup Working Fluid/s & Solid/s:** \n", "\n", "The choice of working fluid is problem dependent. In the absence of obligations to define a specific working fluid such as air or water, it is suggested that the users define the working fluid such that the important non-dimensional groups of interest, such as Reynolds number in this problem, is matched the desired flow conditions. This strategy removes the uncertainty in the choice of the working fluid and will solidify the bases for expected physical observations. The material properties in OpenFoam can be defined in a dictionary file `transportPropertise` located in `constant` folder. In this file the material is set to be `Newtonian` fluid with a defined kinematic viscosity **\\nu**. It should be noted that the numbers in the [] set the dimension of the flow field variable of interest and the other number sets it's value:\n", "\n", "```C++\n", "transportModel Newtonian;\n", "\n", "nu nu [ 0 2 -1 0 0 0 0 ] 0.04;\n", "```\n", "Considering that the cylinder diameter in 1 [m] and user will define a unit constant uniform velocity (1 [m/s]), the value of **\\nu** will be defined as the inverse of Reynolds number. Therefore, changing this value would set the nominal Reynolds number of the flow based on the cylinder diameter. For other cases users should define the working fluid via this process. Various available transport models along with other options within `transportPropertise` dictionary file in OpenFoam can be found [here] (http://cfd.direct/openfoam/user-guide/transport-rheology/#x40-2210007.3.1).\n", "\n", "**3. Setup Boundary and Zone Conditions:** \n", "In OpenFoam the type of boundaries within the CFD domain are defined in the `boundary` file located at `\\constant\\polyMesh`. Whether the users use blockMesh or the convert command to create the CFD domain and discretize the `boundary` file will be created in this process. Furthermore values for the field variables, such as velocity **U** and pressure **p**, at these boundaries are set in `U` and `p` files located at `0` folder located in the working directory respectively. The `0` folder includes the values of the field variables at boundaries for the time `t=0`.\n", "\n", "As discussed earlier in the CFD domain's creation and discretization section, the CFD domain has four faces as boundaries. The boundary and zone conditions settings to develop CFD simulation for **2D Laminar Flow Over a Cylinder** case in OpenFoam are as follows:\n", "\n", "**Inlet:** The flow enters from the inlet face of the CFD domain with constant velocity and uniform atmospheric pressure (i.e. zero gauge pressure). To set the type of this boundary edit `boundary` file in `\\constant\\polyMesh`. Set the type of face to `patch`, which is a condition that contains no geometric or\n", "topological information about the mesh:\n", "\n", "```C++\n", "inlet\n", "{\n", " type patch;\n", " nFaces 300;\n", " startFace 159200;\n", "}\n", "```\n", "Set the field variables at this boundary edit `U` and `p` files in `\\0` folder to a constant uniform velocity in the streamwise direction (x-direction) with zero pressure gradient:.\n", "\n", "```C++\n", "// In U file at 0 folder:\n", "inlet\n", "{\n", " type fixedValue;\n", " value uniform (1 0 0);\n", "}\n", "\n", "// In p file at 0 folder:\n", "inlet\n", "{\n", " type zeroGradient;\n", "}\n", "```\n", "\n", "**Outlet:** The flow exits the CFD domain from the outlet and top face of the domain and reach atmospheric pressure (i.e. zero gauge pressure). To set the type of this boundary edit `boundary` file in `\\constant\\polyMesh` and set the type of face to `patch`, which is a condition that contains no geometric or topological information about the mesh:\n", "\n", "```C++\n", "outlet\n", "{\n", " type patch;\n", " nFaces 300;\n", " startFace 159600;\n", "}\n", "```\n", "Set the field variables at this boundary edit `U` and `p` files in `\\0` to a constant uniform atmospheric pressure with zero velocity gradient:\n", "\n", "```C++\n", "// In U file at 0 folder:\n", "outlet\n", "{\n", " type zeroGradient;\n", "}\n", "\n", "// In p file at 0 folder:\n", "outlet\n", "{\n", " type fixedValue;\n", " value uniform 0;\n", "}\n", "```\n", "If in a problem of interest, there exist a specific pressure difference between the inlet and outlet or other surfaces, that magnitude can be defined in corresponding faces of the domain.\n", "\n", "**top & bottom:** The CFD domain is bounded via two surfaces on top and bottom named after their location within the domain as **top** and **bottom** respectively. These two surfaces have significant distance from the cylinder, which is located in the center of the domain. Therefore, any potential blockage effect or forced accelerated flow around the cylinder is removed. In order to make sure that these boundaries would replicate a semi-infinite domain in two directions a Symmetry boundary condition is set to these surfaces. This condition would reflect identical flow field to the other side of the surface. To set the type of this boundary edit `boundary` file in `\\constant\\polyMesh` and set the type of face to `symmetryPlane`as follows:\n", "\n", "```C++\n", "top \n", "{ \n", " type symmetryPlane; \n", " inGroups 1(symmetryPlane); \n", " nFaces 300; \n", " startFace 160200; \n", "} \n", "\n", "bottom \n", "{ \n", " type symmetryPlane; \n", " inGroups 1(symmetryPlane); \n", " nFaces 300; \n", " startFace 159900; \n", "} \n", "\n", "```\n", "\n", "Set the field variables at this boundary edit `U` and `p` files in `\\0` to `symmetryPlane` as follows too:\n", "\n", "```C++\n", "// In U file at 0 folder:\n", "top \n", "{ \n", " type symmetryPlane; \n", "} \n", "\n", "bottom \n", "{ \n", " type symmetryPlane; \n", "} \n", "\n", "// In p file at 0 folder:\n", "top \n", "{ \n", " type symmetryPlane; \n", "} \n", "\n", "bottom \n", "{ \n", " type symmetryPlane; \n", "} \n", "```\n", "\n", "**Cylinder:** The flow passes over the cylinder located in the middle of the CFD domain and interacts with it based on the no-slip boundary condition. To set the type of this boundary edit `boundary` field in `\\constant\\polyMesh` and set the type of face to `wall`:\n", "\n", "```C++\n", "cylinder\n", "{\n", " type wall;\n", " inGroups 1(wall);\n", " nFaces 400;\n", " startFace 159200;\n", "}\n", "```\n", "Set the field variables at this boundary edit \"U\" and \"p\" files in `\\0` to a fixed value of zero velocity with zero pressure gradient, which represents the fundamental assumptions of the uniform laminar flow over a flat plate:\n", "\n", "```C++\n", "// In U file at 0 folder:\n", "cylinder\n", "{\n", " type fixedValue;\n", " value uniform (0 0 0);\n", "}\n", "\n", "// In p file at 0 folder:\n", "cylinder\n", "{\n", " type zeroGradient;\n", "}\n", "\n", "```\n", "\n", "If the shear forces and formed boundary layer becomes turbulent in this region user should modify this boundary conditions and provide required mesh resolution to capture the phenomena or set this boundary to free slip condition such that fluid elements would not interact with wall region. For low Reynolds number flow, similar to the current problem of interest, a reasonable mesh resolution is sufficient to capture the laminar boundary layer region.\n", "\n", "Since this problem is two dimensional the boundary type and flow field variables for the `frontToBack` boundary, which defines the third dimension of the flow field are set to `empty`:\n", "\n", "```C++\n", "frontAndBackPlanes\n", "{\n", " type empty;\n", " inGroups 1(empty);\n", " nFaces 20000;\n", " startFace 20200;\n", "}\n", "\n", "// In U and p files at 0 folder:\n", "frontAndBack\n", "{\n", " type empty;\n", "}\n", "```\n", "\n", "For more detail on definitions and settings for defining CFD domain boundary conditions in OpenFoam readers are referred to [here](http://cfd.direct/openfoam/user-guide/boundaries/).\n", "\n", "**4. Setup Solution methods:** \n", "In this step, it is highly recommended to use the default/suggested solution methods set in `fvSolution` file located in `\\systems`, unless based on physics of the problem the user is aware of any specific choices. Upon non-smooth convergence and potential divergence of the CFD simulation user can modify and examine various solution methods.\n", "\n", "The `fvSolution` file includes the solution methods for velocity and pressure fields and the required under relaxation factors for flow field variables.\n", "\n", "Now all boundary conditions and settings for the CFD simulation are fully defined. User can **initialize** the solution through an educated guess to start the iteration process. In OpenFoam the initialization of the velocity and pressure fields is set on top of the `U` and `p` files located in the `0` folder:\n", "\n", "```C++\n", "// In U file at 0 folder:\n", "internalField uniform (1 0 0);\n", "\n", "// In p file at 0 folder:\n", "internalField uniform 0;\n", "```\n", "The solution initialization would incept the flow field variables, such as velocity and pressure, based on the defined values by user. For the current problem the CFD domain is recommended to be initialize by values of velocity and pressure at the inlet.\n", "\n", "Iteration process for solving the flow field governing equation now shall start till converged solution is obtained. In OpenFoam the iteration can initiated by running a command line using the name of the application of choice. The name of the application and set marching time steps to solve the flow field governing equations are set in `controlDict` dictionary file located in `\\system folder`:\n", "\n", "```C++\n", "application simpleFoam;\n", "\n", "startFrom startTime;\n", "\n", "startTime 0;\n", "\n", "stopAt endTime;\n", "\n", "endTime 4000;\n", "\n", "deltaT 0.1;\n", "```\n", "For this problem 4000 seconds time interval is set with time steps of 0.1 seconds. The outputs will be written every 100 second till either the end of time interval or set tolerances for the residual of field variables are reached. It should be noted that these are the corresponding settings for case study with Reynolds number of 25 where the flow field is relatively slow. Therefore, a large time interval is required to get a reasonable converged solution. Users can change the length of the time interval and time steps according to the physics of the flow field and numerical simulation requirements. \n", "\n", "Before start running simulation one additional step for simulation setup should be taken. This is the step of setting up drag force coefficient calculation on the cylinder and reporting it at every single iteration. As discussed earlier in the theory section according to dimensional analysis of flow over a cylinder, Reynolds number and drag force coefficients are the two non-dimensional variables that describe the physics of this problem. Therefore, monitoring the convergence trend of this these two variables is vital. To calculate and report the drag force coefficient for flows with different Reynolds number in every single iteration, user can take advantage of pre-defined functions and libraries developed in Open FOAM. One can simply add the following force calculation function into the end of the `controlDict` dictionary file located in `\\system folder` as follows:\n", "\n", "```C++\n", "functions\n", "{\n", "forceCoeffs1\n", "{\n", " type forceCoeffs;\n", "\n", " functionObjectLibs ( \"libforces.so\" );\n", "\n", " outputControl timeStep;\n", "\n", " log yes;\n", "\n", " patches ( cylinder );\n", " pName p;\n", " UName U;\n", " rhoName rhoInf; // Indicates incompressible\n", " log true;\n", " rhoInf 1; // Redundant for incompressible\n", " liftDir (0 1 0);\n", " dragDir (1 0 0);\n", " CofR (0 0 0); // Cylinder origin\n", " pitchAxis (0 0 1);\n", " magUInf 1;\n", " lRef 1; // Cylinder diameter\n", " Aref 2.8; // Projected area of cylinder face (i.e. diameter x depth)\n", "}\n", "```\n", "\n", "This function is called \"forceCoeffs\". User defines \"cylinder\" walls as the patch to calculate the forces on, in the defined directions. The three available directions are `liftDir`, `dragDir` and `pitchAxis`. Based on the orientation of coordinates in the current CFD-domain the drag force direction is along x-axis (1 0 0), the lift force direction is along y-axis (0 1 0) and pitchAxis is along z-direction (0 0 1). The center of origin called `CodR` is located at (0 0 0). The other inputs for this function are reference values for velocity, density, diameter and area used to normalize the drag force to calculate the drag coefficient. Usually the values for these variables are defined as the free stream flow properties. In this case study reference values of velocity, density and diameter are all equal to 1. However, the value of reference area \"Aref\" should be defined differently! Fundamentally, the reference area to calculate the drag force coefficient is defined as the total \"effective\" area that the force of interest is acting against. In case of the flow over a cylinder the fluid forces on the body of cylinder can be visualized as shown in the Fig. 1:\n", "\n", "\n", "Fig.1 - Schematic of total fluid force and it's decomposition on the body of cylinder.\n", "\n", "Each force component can be decomposed into two components in the horizontal and vertical directions shown with red and green colors respectively. Due to the symmetry of the geometry the vertical components of the forces cancel each other. Therefore, the total drag force is the summation of all horizontal forces that are acting against the projected area along the cylinder's wall as visualized in the zoom-in view of Fig. 1. This area is defined as the reference area `Aref` and is equal to diameter times the depth of the cylinder. For the current CFD domain the diameter is 1 [m] and the created depth of the mesh via \"blockMesh\" is 2.8 [m]. Therefore, the value of the reference area is 2.8 $[m^2]$. Note that for a purely 2D problem (with no depth in mesh) the value of depth is considered to be unity. Therefore, in that case the value of the reference area is 1 $[m^2]$\n", "\n", "At this stage the iterations for this simulation can be started by running the command `simpleFoam` to get outputs of residuals or `foamJob simpleFoam` to just save the residual values and plot them later. There are two reasons behind this recommendation: one is that this command will not print the residual values to the screen and therefore will be faster in long term. Second this command to store history of all residuals and in the end of simulation running the command `foamLog logs` will print them in data sets in a directory called `logs`. If users need to see the iteration progress and residual values the command `tail -f log` can be used." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.5.1" } }, "nbformat": 4, "nbformat_minor": 0 }