A Simplified Mesh Generation Scheme for 3D Geometries Composed by Planar Faces to be used with BEM

— This paper presents the development of a simplified mesh generation scheme for three-dimensional boundary element method (BEM). The developed program presented in this work takes into account three-dimensional geometries composed by planar faces. This is a simplification regarding the form of the acceptable geometries, but a great variety of geometries and problems can be modelled with sufficient accuracy to perform numerical analysis. The main program’s idea consists of moving each face belonging to the three-dimensional geometry to bi-dimensional space using geometrical transformation matrices. In bidimensional space, then, it is applied a bi-dimensional mesh generation algorithm for element generation. After element generation is done, new geometrical transformations are applied in order to send the meshed face back to its original position in three-dimensional space. The continuity of the final mesh is assured by promoting the discretization of the edges of the geometry only once. The proposed scheme has several positive aspects regarding simplicity, efficiency and robustness. A practical and immediate use of the proposed algorithm can be found for those who already have a bi-dimensional mesh generator code implemented and intent to expand its functionalities to treat simple three-dimensional geometries composed by planar faces. Several examples of meshes in arbitrary three-dimensional geometries composed by planar faces are presented in order to illustrate the capabilities of the developed computational program and some computational simulations have been performed to show the quality of the meshes in problems with specific boundary conditions.


I. INTRODUCTION
The mesh generation process is crucial for computational simulations since it strongly affects the results obtained by the numerical solution of partial differential equations that govern the analyzed problem [1][2][3]. The generated mesh should be able to provide sufficiently accurate results according to the complexity of the geometry and physical phenomena involved. Several complex applications can be analyzed with sufficient accuracy level [4,5]. In computational numerical simulations, the pre-processing phase is one of the "bottle-neck" with respect to the computer processing time. This is due to the fact that in the pre-processing phase the mesh is generated and the problem is modelled with specific boundary conditions. It is very difficult to generate suitable meshes for problems involving complex threedimensional geometries and/or boundary conditions. The algorithms for three-dimensional mesh generation are usually very difficult to be implemented for arbitrary cases and the computational costs are generally expensive [6]. That is why, still nowadays, great efforts and investments are done to develop new techniques and algorithms for improving the mesh generation process on threedimensional geometries. In this paper it is presented a simplified scheme to implement a three-dimensional mesh generator to be used in numerical analyses by the boundary element method. The boundary element method has some particular features that make it more favorable in certain types of applications, when compared to other numerical methods. Among those features, there is the fact that, in boundary element method, only a discretization of the geometry boundaries is necessary [7,8]. In other words, in three-dimensional cases only surface elements must be generated. This is a huge benefit from the mesh generation point of view. The algorithm proposed in this work only takes into account three-dimensional geometries composed by planar faces. This is a simplification regarding the form of the acceptable geometries, but a great variety of geometries and problems can be modelled with sufficient accuracy to perform numerical analysis. The main algorithm's idea consists of moving each face belonging to the three-dimensional geometry to bi-dimensional space using geometrical transformation matrices. In bidimensional space, then, it is applied the Delaunay triangulation method over a grid of points to cover the face domain with triangular-linear elements, but other algorithms and elements could also be used taking into account the general idea. The Delaunay triangulation method is one of the best methods to generate meshes with triangular elements devoted to numerical analysis. It has the property of maximizing the triangles' minimum internal angles which is favorable and suitable for numerical methods in general [9,10]. After element generation is done, new geometrical transformations are applied in order to send the meshed face back to its original position in three-dimensional space. The continuity of the final mesh is assured by promoting the discretization of the edges of the geometry only once. This is essential because each edge belongs to two planar faces of the three-dimensional geometry. Therefore, no duplicate nodes should be generated over an edge and the final mesh is supposed to contain only conformal elements. The proposed scheme has several positive aspects regarding simplicity, efficiency and robustness. Firstly, the code is extremely simplified since the problem is transformed from three-dimensional space to bi-dimensional space by using geometrical transformations matrices. The algorithms for element generation are bi-dimensional, that is to say, they are much simpler, reliable and easier to be implemented [11][12][13]. Regarding efficiency and robustness, as each face of the three-dimensional geometry is moved to bi-dimensional space and the mesh is generated over one face per time, the occurrence of a bad element generation is minimized. As the whole mesh is performed face by face, the problem is subdivided in several small problems and the main algorithm assumes the divide-to-conquer paradigm, which assures high efficiency and computer processing velocity [14,15]. A practical and immediate use of the proposed algorithm can be found for those who already have a bi-dimensional mesh generator code implemented and intent to expand its functionalities to treat simple three-dimensional geometries composed by planar faces. The generated surface meshes over these geometries can be used in threedimensional boundary element analysis or, even though, in numerical analysis using the shell finite element formulation. All the code has been written using the object-oriented paradigm in C++ combined with UML notation. Using this approach the program can be easier modified, the maintenance costs are reduced and new implementations can be carried out as user's and/or programmer's needs [16].

II. THREE-DIMENSIONAL BOUNDARY ELEMENT METHOD
The boundary element method is based upon boundary singular integral equations. The analytical formulation involves the transformation of the governing differential equation applicable to the whole domain into an integral over the boundary [7]. To illustrate the technique it is presented the elasticity problem in the form of the partial differential equation known as the Navier equation of elasticity: , (1) where are body force components, is the shear modulus, is the Young's modulus, is the Poisson's ratio and are displacements components.
Equation 1 can be transformed into an integral equation over the boundary. The displacement boundary integral formulation of elasticity can be derived using the Somigliana's Identity and Betti's reciprocity theorem [8].
The fundamental solutions for displacements and tractions are given respectively by: where is a function of the contour shape at the boundary point . From Equation 4 one can note that there are only boundary dependent terms. Thus, it is not necessary to generate domain elements in the discretization procedure. That is why the threedimensional mesh generator implemented in this work only generates surface elements. For simplicity, the mesh generator was implemented by using a bi-dimensional Delaunay mesh generation algorithm, which gives origin to triangular-linear elements, as shown by Figure 1. Therefore, a boundary element solver with threedimensional elasticity isotropic-linear homogeneous formulation, for example, can be used to perform numerical analyses considering solid mechanics engineering problems. Using the shape functions , , , Equation 4 can be written in a discretized form as [7]: where is the number of nodes per element, is the total number of elements on the mesh, represents the number of the node on the mesh that is being evaluated as source point and is the Jacobian of the transformation.
The Jacobian can be calculated by , where the subscripts  and  denote the derivatives with respect to  and , respectively. Equation 5 can be presented in the matrix form as: , (6) which can be rearranged in and readily solved.

III. MESH GENERATION
This work uses bi-dimensional Delaunay triangulation over a grid of points to generate elements in an arbitrary planar faces of the three-dimensional geometry. There are several ways to perform bi-dimensional mesh generation over an arbitrary planar straight line graph (PSLG). As bidimensional algorithms for mesh generation are extensively know in literature, the main focus of this paper is devoted to present the structure of the program that allows three-dimensional surface mesh generation. To perform this task, geometrical transformation matrices are used to move each planar face from the three-dimensional geometry to bi-dimensional space and vice-versa.

Data Structure
The process to obtain the input file for the mesh generator presented in this work is illustrated by Figure 2. First, the three-dimensional geometry is drawn using CAD software. The three-dimensional geometry must only be constituted by planar faces. Then, the geometry is saved on IGES format and an IGES translator developed during this research project is used to generate the input file accepted by the mesh generator program. This input file contains the geometry information in the format presented by Figure  3(a). A similar data structure for the input file can be found in [11]. The MeshPar parameter is a numerical value provided by the user. It is given in length unit. This parameter affects the level of the mesh refinement. The key-words Vertex, Edge and Facet contain information about the vertexes, edges and planar faces of the 3D geometry, respectively. The format of the output file containing the mesh data is presented by Figure 3

Computational Aspects
In this subsection the main classes of the developed program are presented. The classes were modeled by using the UML notation. A special focus is given to the Node, Segment, Facet, Element, Triag, Mesh and GeoTrafo classes, which are the main classes of the program. In Figure 4 are presented the Node, Element, Triag and Centroid classes. It is also possible to see the Dictionary class, with its hidden attributes and methods. The Dictionary class is very important for the developed program because it makes possible to read the input data file easily and efficiently. It can be seen from the diagram presented in Figure 4 that both the Element and Node classes depend on the Dictionary class. In fact, the Segment and Facet classes also depend on it, but they were not shown in this diagram for simplification purpose.

Fig.4: Class diagram exhibiting the relationship between
Node and Element Classes.
The Element class has a number of attributes and methods, as can be seen. The signature of the methods, i.e., the type of method and the arguments that are passed to it were not shown at the diagram for visual simplicity. The give_CurrentWordInterpretation() method is a static method and, therefore, it can be called without the need to instantiate an object. It is only necessary to use the class name and the "::" operator followed by the name of the method. The Triag class is the only class derived from the Element class, because generated meshes consist only of triangular elements. The give_CurrentWordInterpretation() method is an specialist that has the ability to interpret input file information related to a word previously added to the dictionary. The Node, Segment, and Facet classes also have the give_CurrentWordInterpretation() method to read the data structure shown by Figure 3. The other methods of the Element class, give_Centroid(), give_Jacobian(), and give_Area() are virtual methods, which allow the assignment of Element-type pointers to Triag-type objects.
Thus, the code for calculating the centroid, the Jacobian and the area of a triangular element are implemented in the methods of the Triag class. It is possible to note two associations existing in Figure 4: an association between the Element class and the Node class; and, an association between the Element class and the Centroid class. In the first, a Node-type object is aggregated by reference to the Element class, which implies that every Element-type object will have a Node-type attribute, called pCircumCircle, aggregated by reference. This pCircumCircle attribute will have the coordinates of the center of the circle that circumscribes a Triag element. In the second, an object of the Centroid-type, theCentroid, is added by value to the Element class. Figure 5 shows the Segment, Facet, and Edge classes. Also, the List and Container classes are presented. They are Template classes, that is, they are parameterized. Note that the List class depends on the Node, Segment, Facet, Edge and Element classes, because it is necessary to create lists of Node-type objects, Segment-type objects and so on. The Edge class was abstracted aiming at the improvement of the element generation algorithm. With the Edge class it is possible to know the geometric entities in a certain region of interest more quickly and efficiently. For example, as the Edge class relates to the Element class and the Node class, it is possible to see that each Edgetype object has two Element-type objects (pTriag1 and pTriag2) and two Node-type objects (pNode1 and pNode2). Thus, for a particular Edge-type object, it is immediately known which are the two triangles of the mesh that are part of it, and which are the two nodes that define it. This is extremely beneficial from the processing time saving point of view, because, it is not necessary to search for a whole set of geometric entities, but only in a particular region of interest. Figure 6 shows the GeoTrafo and Mesh classes that are the most important classes of the developed program. Again, to facilitate viewing of the classes, the signatures of the methods were not displayed. Fig.6: GeoTrafo and Mesh classes modeled using UML notation.

Mesh Generator Flowchart
The execution flowchart of the implemented mesh generation program is presented in Figure 7. The mesh

Geraldo Creci
International Journal of Advanced Engineering Research and Science, 8(7)-2021 www.ijaers.com generation process begins with the input file reading. As soon as this file is read, the computer has already stored in memory all lists needed to generate the mesh on the geometry. Three lists are created: a VertexList which contains Vertex-type objects with information about the geometry vertexes; an EdgeList which contains Edge-type objects with information about the geometry edges; and a FacetList which contains Face-type objects with information about the planar faces of the geometry. The MeshPar value is stored in a double type variable. The mesh generation process starts from a loop in the FacetList. For each face belonging to the FacetList, it is applied a geometrical transformation matrix which will move it to bi-dimensional space. After the face is in bidimensional space, a loop in the EdgeList is performed, which defines the current face. Each face owns as attribute an EdgeList which contains all the segments of the face. Each segment of the face owns as attribute a VertexList which contains the Vertex-type objects that define the segment. Initially, it is checked if the VertexList of the current segment has only two Vertex. This means that the segment has not been discretized. If the VertexList of the segment contains more than two Vertex, it means that the segment has already been discretized and should not be discretized again. This verification is crucial because each segment belongs to two faces. So, in order to respect the continuity of the final mesh, the nodes on the edges of the geometry are defined a priori and only once. Afterwards the loop in the EdgeList of the Face object is finished, it is called the method that generates the mesh over the face using a grid of points. Good references to implement a bidimensional mesh generator using the Delaunay triangulation method, can be found in [12,[17][18][19]. However, there are several approaches to generate meshes with different types of elements [13]. An advancing front approach, instead of generating a grid of points over a region, is also very efficient for generating triangular meshes [14]. This technique is particularly suitable for the boundary element method because it starts the element generation from the boundary data. An execution flowchart using the advance front technique for mesh generation, whose structure is similar in some aspects to the proposed flowchart, can be found in [15]. As soon as the mesh on the face is generated in a conformal procedure, the inverse of the transformation matrix, used to move the face to bidimensional space, is applied to send the meshed face back to its original position in 3D space. Then, another face of the geometry is selected and the process goes on until all the faces have been analyzed. At the end of the process, the program will have stored a NodeList containing the nodes of the final mesh and an ElementList containing the nodal connectivity. An output file containing this information is written and used as input file by the boundary element solver.

Geometrical Transformations
The process of moving each planar face of the geometry from three-dimensional space to bi-dimensional plane is obtained by using geometrical transformation matrices [20]. In Figure 8 this process is illustrated step by step. In Figure 8(a) it is presented an arbitrary planar face in three-dimensional space constituted by linear segments. This face can be a non-convex polygon and have holes in its interior. First, it is detected the point of the face which has the minimum -coordinate, in this case . Then it is applied a translation matrix in order to move all the face entities to the origin of the coordinate system, as presented in Figure 8(b). The translation matrix is given by: where is the point used for the face translation. The subscript is the number of the point, represents the cartesian components of the point, and, is the Kronecker delta function. The vector is the face normal vector calculated by the cross product between and . As soon as the face is moved to the origin of the coordinate reference system, it is necessary to apply a geometrical transformation matrix in order to align with -axis and with -axis. This is shown by Figure   8(c). The matrix is given by: where are the cartesian components of the vectors.
The vectors can be calculated by , and . So, the transformation matrix to be implemented is given by the product of matrices with , as: Applying on a set of points , the face is moved to plane as shown in Fig. 8(d) by the transformation: The inverse of the transformation matrix is used to send the planar face back to its original position in threedimensional space.

Generated Meshes
In order to show the capabilities of the developed mesh generator, it is presented four meshes in arbitrary threedimensional geometries composed by planar faces, see  the Delaunay triangulation method, the triangles' minimum internal angles found in every mesh are significantly high, being limited only by the minimum angle of the drawing geometry. Therefore, it can be noticed that all presented meshes shown good quality and are suitable for numerical analysis using the threedimensional boundary element method or shell finite element formulation. The meshes presented in this paper are for illustration purposes and, thus, they have a very small number of nodes and elements. However, the MeshPar parameter can be easily adjusted as user's needs in order to generate more refined meshes, with much more nodes and elements, according to computational hardware availability.

Computational Simulations
Some computational simulations were performed in order to demonstrate practical applications of the generated meshes. The first problem analyzed is a tub with internal pressure. The Poisson ratio used is equal to 0.3 and the Young´s Modulus is 2600 units. The boundary conditions applied are: • -direction displacement restriction on the nodes of the face parallel to plane; • -direction displacement restriction on the nodes of the face parallel to plane; • -direction displacement restriction on the nodes of the faces parallel to plane; • 50 units of internal pressure applied on the elements of the internal wall tub. Figure 13 shows the displacement map for the studied problem considering that the modeled mesh is generated by the mesh generator presented by this work having 200 nodes and 396 triangular-linear elements; and, the solution is calculated by a boundary element solver, ECon-3D, which was developed prior to this research project. In a similar way, the same problem was analyzed again considering a modeled mesh generated by the Ansys software having 212 nodes and 420 Shell63 triangularlinear elements; the solution shown by Figure 14 was also calculated by ECon-3D solver. It can be seen that the mesh modeling strongly influences the results. The analytical solution for this problem can be found in [7] by means of comparison. A good agreement can be observed among all the achieved results.  A second case study was performed considering the structural support whose geometry is presented by Figure  15. To analyze this problem, a Poisson ratio of 0.21 and a Young´s Modulus of 2900 units were adopted. The applied boundary conditions applied are: • -and -direction displacement restriction on the nodes of the two holes; • -direction displacement restriction on the nodes of the two outer side faces containing the holes; • 500 units of pressure applied to the elements at the rod end. Figure 15 shows the displacement map of the solution considering the problem modeled with a mesh generated by the Ansys software having 934 nodes and 1832 Shell63 elements. The respective solution of this problem was calculated out by ECon-3D solver. Figure 16 shows another displacement map of the solution considering the problem modeled with a mesh generated by the mesh generator developed by this work having 936 nodes and 1876 elements. The respective solution of this problem was also performed by ECon-3D solver. A final analysis was performed using the finite element method in the Ansys software by using Solid45 elements. Figure 17 shows the mesh generated by the Ansys software with 1618 nodes and 6013 elements. The Solid45 element was chosen because it has linear interpolation functions, thus allowing a fairer comparison among the results. All three analyzes showed very close results. This demonstrates that the meshes generated by the mesh generator developed by this work are of good quality and suitable for using in numerical analysis.

V. CONCLUSIONS
It was developed a triangular-linear surface mesh generator for arbitrary three-dimensional geometries composed by planar faces. As programming was done using the object-oriented paradigm, it was possible to construct a very clear and efficient program structure. In addition, this structure is extremely flexible and allows for faster and easier code expansion and enhancement. Although the mesh generator developed by this work only accepts geometries composed by planar faces, it can be seen from the analysis of the presented results that a wide variety of three-dimensional problems can be modeled with considerable complexity. All the geometries analyzed by the program lead to the generation of satisfactory meshes which can be progressively refined according to the MeshPar parameter and the processing capacity of the computer hardware. The proposed algorithm makes use of geometrical transformation matrices which allow the development of a fast, robust and efficient program. The program flowchart can be particularly useful for those who already have a bi-dimensional mesh generator code implemented and intend to extend its functionalities to treat simple three-dimensional geometries composed by planar faces. Some computational simulations have been performed to show the quality of the meshes in problems with specific boundary conditions. It could be observed that the developed mesh generator provides triangularlinear elements with good technical features for numerical analysis considering the boundary element method or the shell finite element formulation. Besides, the developed program is an open source engineering program that can be used for research purposes. Future improvements can be done according to user's needs.