19 KiB
simulation_workflow
To do a simulation of a reactor, that has already been simulated once, one can use this workflow to aqcuire the same results. For other reactor geometriesd the workflow hasd to be adjusted accordingly. The Following tools are needed for that: freeCAD, GMSH, openFoam and OpenMPI.
You start the workflow with freecad
Freecad
Freecad is only needed to construct the circles with the diameters of the pipereactor, because everything else is extruded in gmsh. For that use the tools to draw the circles with thew sketch tool. If you want to build the reactor tubes with membrane and swap as whole regions you need the three circles depicting the crosscut of the pipes. See picture below.
If you have drawn the circles you must split the the circles into equivalent arc sections to get the arch points. those will be needed in gmsh to refine the mesh. You can split the arch with the *split edge command:
If you split it equaly it should look something like this example depicted below.
Save the drawing as .brep and put it into the gmsh folder.
gmsh
Make a .geo file that adapts the name of your .brep file or choose another descriptive name
You can do that wit the command: touch <name_of_file>.geo
Importing the brep and setting up the geo kernel
open the file and write the import command at the top line of the file
Merge "<name_of_brep_file>.brep";This comnmand imports the curves and points from the brep file specified. It is important to include a ; at the end of every line. This implies the end of the command, not the end of the line, because you can theoretically write the whole geo file in a single line separating every command with a ;. But for the sake of code clarity we give every command its own line. At least for this workflow. It is strictly advised to also do so in every other project.
In the next the next step you choose the Kernel taht should be used for the Computation. Therefore write the command: The best option is the OpenCascade geo Kernel. This one should always be used. The alternative is the gmsh kernel
SetFactory("OpenCASCADE");Setting up Mesh refinement
This section can be inserted at any point before actually running th emesh command. But for simplicity we define the refinement variables at the top of the document.
//+ Mesh refinement
// Variables
i_ms = 3.5; // inner mesh size
i_mr = 0.15; // inner mesh refinement
m_mr = 0.15; // membrane mesh refinement
o_mr = 0.15: // outer mesh refinementWe set the refinement to some arbitrary values, those can be changed according to the mesh at hand. In addition to defining the meshing variables we have to define a point in the middle of the pipes. In our example we used the center of the coordinate system as geometric centre as well. So we define the Point as follows:
Point (25) = {0,0,0, i_ms}; //Add Point to control inner mesh size
The syntax for creating the point is as follows
<Geometric entity (Point)> (<ref No>) = {x,y,z,<meshsize at point>};
The reference number has to be unique for each geometric entity. That means there can be a Point(1) as well as a Line(1) but no two points with number (1). There are 24 Points in the brep file so to define a new point we need to assign it the number 25
Now we have to set up the points from the brep file with their intended Mesh size. Because the Points are already defined inside the brep file, we can not create them new like we did with Point 25. That has the advantage that we do not need to now the exact position of the points. But now we need to assign the Mesh size to the Points.
//membrane
MeshSize {1,2,3,4,5,6,7,8} = m_mr;
// outer_circle
MeshSize {9,10,11,12,13,14,15,16} = o_mr;
// inner_circle
MeshSize {17,18,19,20,21,22,23,24} = i_mr;
This command assigns the mesh Size of i_mr m_mr and o_mr to the points listed between the curly brackets. As we did define the Meshsize of Point 25 while creating the point, we do not need to do that anymore. The order in which the arches are numberded depends on the order in which we defined the arches in the CAD file. But the arches of a whole circle are always numbered in order.
Boolean addition of the arches
If you now make a Line Loop out of those 8th circles, it would create cylindric walls based on the length of those arches instead of a whole cylindric wall at the interface region. To prevent this behaviour, we connect the arches beforehand and create a whole circle again. This sounds contra intuitive, because in the previous step we split up the circles and now we add them back together. This is because of the way freecad handles the export of points into the brep file format. Only points belonging to a curve are exported into the file, even if you fuse them to the curve inside the CAD file they wont be exported. Those Points are needed for the mesh refinement, so thy have to be included.
To fuse the lines together, we use the Boolean Union command and add them part by part. This is done as follows in the .geo file (explained with the middle curve as an example):
BooleanUnion{Line{1 };}{Line{2};}
BooleanUnion{Line{25};}{Line{3};}
BooleanUnion{Line{26};}{Line{4};}
BooleanUnion{Line{27};}{Line{5};}
BooleanUnion{Line{28};}{Line{6};}
BooleanUnion{Line{29};}{Line{7};}
BooleanUnion{Line{30};}{Line{8};} // to generate Line 31The first command gets the first and second eights arches as input and fuses them together to create a new arch-line. There are already 24 Arches in the .brep file so the next one to be created gets the consecutive number (25). To add another arch to that line, we add this line and the next arch in the series to the Boolean union, that creates another longer arch. This is repeated until the circle is full and the last line (full circle) is created (31). The commands inside the curly brackets have to be escaped (ended) with a ; but not the boolean statement, because the way it is intended to be written is as follows:
BooleanUnion{
Line{1 };}{
Line{2};
}The other circles are done accordingly. This then creates the circle lines 38 and 45.
creating Surfaces
To transform those circles to pipes we first have to create surfaces out of thos lines. We cannot create a surface by using the curves (lines) we just created because for gmsh those lines are not closed (unclosed lines are called wires). To close them to form a loop that can wrap a surface the command Line Loop() is used. It takes an ordered list of the lines to loop together and creates with them a line loop that needs a reference number.
Curve Loop(1) = {31}; //middle
Curve Loop(2) = {45}; //inner
Curve Loop(3) = {38}; //outerIf it is done correctly it should look like the picture below.
With those loops loops we can create the three surfaces for the pipes. For that to work, we heve to visualize the cylindrical structure of our parts. The first surface is simple. It only wraps around the inner loop (2). The other surfaces need to be created by subtracting the curves for the hollow cylinders. For example the membrane hollow cylindre is created by subtracting the inner curve loop from the middle loop (1-2 => {1,2}). The surface is created by wrapping those loops with in a plane surface.
Plane Surface(1) = {2}; //inner
Plane Surface(2) = {1,2}; //membrane
Plane Surface(3) = {3,1}; //outerIf the visibility of the surface labels is toggled on, it should look something like this:
In the last surface creation step, the middle point (25) is fused to the surface. The others do not have to be fused, because they already are inside the .brep file.
Point {25} In Surface {1};creating boundary layers at the interfaces of the surfaces
In most cases there is friction on the walls or other wall dependend regimes that could impact the overall flow simulatiuon. therefore in those regions we need so called boundary cells or a boundary layer. Doing that is not easy in a 3 dimensional body. that is the reason we use the 2D cylindrical surface to define the boundary layer end extrude the surface to a volume afterwards. To simplify the mesh creation for every boundary layer to be created, we define some variables to be used. Those must be changed according to your geometry.
//+ Boundarylayer-fluid
// Variables
nb_layer = 1;
layerthickness = 0.2;
hfar_l = 0.2;
hwall_nl = 0.2;
layersize = 0.175;
ratio_l = 1.1;
q = 1;
im = 1;
//+ Boundarylayer-solid
// Variables
nb_layer_s = 1;
layerthickness_s = 0.2;
hfar_ls = 0.2;
hwall_nls = 0.2;
layersize_s = 0.1;
ratio_ls = 1.1;
qs = 1;
ims = 1;Those variables are set for the solid and fluid regions accordingly. The Boundary layers in gmsh are set via the Field function. So to generate those layers we first have to define the field function as a boundary field. In the next step we define the edges of the surface, were the layer should be deployed. It is important to know that you cannot list line loops or curve loops in this option. In our example we set the edges of every surface respectively because we do not want to define the same boundary twice. So when we list the edges we have to exclude the adjacent other surface of the curve the boundary layer is deployed on. The code below shows an example of defining the boundary layer for the inner pipe.
// inner_circle
Field[1] = BoundaryLayer;
Field[1].EdgesList = {45};
Field[1].hfar = hfar_l;
Field[1].hwall_n = hwall_nl;
Field[1].thickness = layerthickness;
Field[1].ratio = ratio_l;
Field[1].IntersectMetrics = im;
Field[1].ExcludedSurfacesList = {2};
Field[1].Quads = q;
BoundaryLayer Field = 1;
Field[1].NbLayers = nb_layer;
Field[1].Size = layersize;
As can be seen here, we list the edge 45 as the curve were the boundary layer should be positioned and exclude the adjacent surface (2) to that edge to prohibit the boundary layer on the other surface. The other defined variables are set do their associated field. the other two boundaries are listed in the example below.
// outer_circle
Field[3] = BoundaryLayer;
Field[3].EdgesList = {38,31};
Field[3].hfar = hfar_l;
Field[3].hwall_n = hwall_nl;
Field[3].thickness = layerthickness;
Field[3].ratio = ratio_l;
Field[3].IntersectMetrics = im;
Field[3].ExcludedSurfacesList = {2};
Field[3].Quads = q;
BoundaryLayer Field = 3;
Field[3].NbLayers = nb_layer;
Field[3].Size = layersize;
// membrane
Field[2] = BoundaryLayer;
Field[2].EdgesList = {31,45};
Field[2].hfar = hfar_ls;
Field[2].hwall_n = hwall_nls;
Field[2].thickness = layerthickness_s;
Field[2].ratio = ratio_ls;
Field[2].IntersectMetrics = ims;
Field[2].ExcludedSurfacesList = {1,3};
Field[2].Quads = qs;
BoundaryLayer Field = 2;
Field[2].NbLayers = nb_layer_s;
Field[2].Size = layersize_s;Extrusion of surfaces
Now that we have defined the boundary layers we can extrude the volume of thy cylinders from the circular surfaces. WE extrude not only the geometrical volumes but also opt to extrude all future mesh layers that are created after the 3D meshing. The commands to extrude are as follows:
First we set up the variable that represents the length of the extrusion into the Z dimension.
// Extrusion
z_dim = 320; //length in mm
This variable is used in the actual command to extrude. As we only want extrusion in the z doirection we set the x and y dimension to 0. Inside the curly brackets we define the surfaces to be extruded. As we want to choose all circular surfaces we list all of those in the command and end every surface with a ;. The Layers specification sets the number of mesh layers that are extruded from the surface. We set this up to not be too many layers because this will be only a flow region and nothing complicated will happen there. The Recombine parameter adds the triangular mesh cells to rectangular cells together. Those are better for structured meshes in the Z direction. The Extrude command on its own ends without a semicolon after the curly bracket but we added a variable declaration at the forefront of that command (extr[]=). This assigns all curves, surfaces and volumes that are created with the extrusion of the chosen surfaces to on cell of the extr array. which surface and volume are written in which cell is set and can't be changed. In Our example no curves are extruded, as that can only happen when points are specifically chosen for an extruseion. The extr[0] goes to the surface that is created from the first surface chosen for extrusion and is on the exact opposite site of the extrusion axis (mirror suface). The second entry of the extr vector is assigned to the extruded volume of the first chosen surface. the third entry (extr[2]) is set to the mantle surface that surrounds the struded surface (cylinder mantle surface). The other surfaces and volumes are set according to this specification with increasing numbers. This setup is specifically useful when extruding multiple pipes or reactor parts is the goal, because the numbering inside the extr[] vector is consecutive and can therefore be easily applied to a for-loop.
//+
// 320 for 32cm in z-Direction (10 lengthunits = 1cm)
// 160 for layers in z-Direction (note that after 3D Meshing the layers will be doubled, 10 layer per 10 lenghtunits)
extr[] = Extrude {0, 0, z_dim} {
Surface{1};
Surface{2};
Surface{3};
Layers {z_dim/48};
Recombine;
};Setup of physical groups
After all the surfaces are extruded, we close the setup of the .geo file with the definition of the physical groups. Only curves, surfaces and volumes defined in those grpups get exported alongside the mesh to later be defined as walls, inlets outlets and other boundaries. We use the varaible vector that we created in the last step to assign the physical groups. We extruded two flow pipes (inner and outer pipe) and a membrane layer so we only need to define two inlets and two outlets. the other top and bottom durfaces of the membrane are later set to walls. The mantle surface and the three volumes are defined as well.
// grouped surfaces and Volumina (note that you have three volumina)
//+
Physical Surface("inlet_inner", 9) = {1};
//+
Physical Surface("wall_bottom_membrane", 10) = {2};
//+
Physical Surface("inlet_outer", 11) = {3};
Physical Volume("inner_zone", 55) = {extr[1]};
Physical Volume("membrane", 56) = {extr[4]};
Physical Volume("outer_zone", 57) = {extr[8]};
//+
Physical Surface("outlet_inner", 26) = {extr[0]};
//+
Physical Surface("wall_top_membrane", 27) = {extr[3]};
//+
Physical Surface("outlet_outer", 28) = {extr[7]};
//+
Physical Surface("mantle_inner_membrane_complete", 1001) = {extr[2]};
//+
Physical Surface("mantle_membrane_outer_complete", 1002) = {extr[5]};
//+
Physical Surface("mantle_outer_complete", 1003) = {extr[9]};Meshing
If you did the forementioned steps successfully, your geometry in gmsh should look like this when you load the .geo file:
Meshing steps can be included here, but if you did not already create a mesh for this geometry, it is better to make the leftover steps with the GUI of gmsh as the mesh is visualized there. The options for the meshing algorithm are as follows:
//Mesh-Options
Mesh.Algorithm = 6; // Frontal-Delaunay for 2D meshes
Mesh.Algorithm3D =10; // HXT Parallel (if neccessary)
Mesh.RecombinationAlgorithm = 1; //blossom
Mesh.SubdivisionAlgorithm = 2; //all hexahedral
Mesh.Smoothing = 10; //smoothing steps
Mesh.MeshSizeFactor = 0.6; //element size factor
Mesh.MeshSizeMin = 0.7; //min mesh size
Mesh.MeshSizeMax = 2.2; //max mesh size
//initiate Meshing
Mesh 3; //3 for creating 3D Mesh
All that is left is to sve the file:
Save "<name_of_mesh_file>.msh2";
The .msh2 extention specifies that the msh extension should be used but with the ASCII 2 support for the file export.
Now the .geo file needs to be computed with the gmsh programme. There are many flags that can be invoked here, see the documentation on GMSH for more details. But we will only set the "parse file flag "and the number of threads flag for this operation to cut the computation time in half.
gmsh -parse_and_exit -nt 40 /path/to/geo/file/<name_of_geo_file>.geoThe important part to remember is, that the .geo file has to be in the same folder as the corresponding .brep file or the path has to be given while importing the .brep file.
This concludes the gmsh part of the workflow
converting the mesh
To Convert a GMSH mesh to an OpenFoam campatible format we have to transform it with the following openfoam command:
gmshToFOAM <name_of_mesh>.msh2This will create a lot of warnings because we do a multibody Mesh and openfoam on its own does not recognize that the cells in between two volumes belong to both volumes. This warning can be ignored. So in the next step: If the Mesh in question is made to be applied to multibody simulation we need to split the mesh after conversion. The mesh is split at the boundary edges automatically.
splitMeshRegions -overwrite -cellZonesOnly
This command does not need a mesh file as input parameter. It takes the foam mesh created in the previos step and splits it accordingly. The -overwrite flag deletes all previous intries in the occuring files. The flag -cellZonesOnly does not do a walk and uses the cellZones only. Use this if you don't mind having disconnected domains in a single region. This option requires all cells to be in one (and one only) cellZone.