-
Notifications
You must be signed in to change notification settings - Fork 0
Home
The finite elements/volumes are defined as irregular hexahedral cells. A hexahedron is defined by its eight corner nodes. Therefore, the position of such an element is completely known by the coordinates of the corner nodes. For this to be in a logical way, the numbering of the nodes must be done in a fixed sequence. Here we adopt the principle that the numbering is first anti clockwise in the base plane and then, starting above the first node, anti-clockwise in the top plane. However, what is down (the base plane) and up (the upper plane) is free to decide by the user. Since it is rather tedious to compose such information for a large number of elements and nodes, a finite element grid generator is a useful tool.
The grid generator works as follows. First we have to define the coordinate system:
- It is a right-handed 3D cartesian coordinate system;
- The z-axis is vertical and pointing upwards;
- The x and y axes can be chosen arbitrarily in a horizontal plane.
The study area is divided into a structured finite element grid. By structured we mean that the grid can be considered as a 3D matrix, such that nodes can be numbered and identified by three indices, the first one is i, a column number going from 1 to ni, the second one is j, a row number going from 1 to nj, and the third one l is a layer number going from 1 to nl. By convention we consider i to go from left to right, j to go from front to back, and l to go from top to bottom. In a structured mesh ni, nj, and nl are constants independent of the row, column or layer positions.
The fact that structured grids are less adaptable to real physical space dimensions of study areas, is of no concern here, because the user will have the possibility to exclude elements from the mesh in the simulations later.
It is further assumed that projections of all elements into a horizontal plane coincide with each other, and lead to a cubiform draughtboard structure, with rows i and columns j. In the horizontal plane quadrilateral elements can be distinguished such that all rows contain the same number of elements , ni−1, and each column contains nj−1 number of elements. This implies that x and y coordinates are independent of the layer number l.
Each node is defined by an i-index, a j-index and a l-index, indicating respectively, the column, row and layer containing the node. The coordinates of a node are denoted by x(i, j), y(i, j) and z(i, j, l). Notice that the x and y coordinates are independent of the vertical position.
The finite volume/element generator consists of two parts. The first part of the mesh generation constructs a two-dimensional mesh in the horizontal plane, such that all x and y coordinates are defined. There are two types of nodes: nodes with coordinates fixed by the user as desired, and nodes automatically interpolated by the mesh generator.
Fixed nodes can be specified individually one by one, but in order to speed up the procedure, series of nodes lying on a line segment can be fixed simultaneously. In the first case, the user will provide for each fixed node the indices i and j and the coordinates x(i, j) and y(i, j). In the second case, the user gives the indexes and coordinates of the first node and the last node, and the program will fix the coordinates of intermediate nodes, such that they are spaced at equal distances from one another. Notice that this only makes sense when other nodes exist between the two end-nodes.
The rule for identifying intermediate nodes is as follows:
- Let the first node be (i1, j1), and the last node (i2, j2);
- Calculate ia and ja as the absolute values of i2 − i1 and j2 − j1, and md as the max divider of ia and ja (if ia = 0 then md = ja and vice versa);
- Now, there are md − 1 intermediate nodes, which can be denoted by: [i1 + sign(i2 − i1) × n × ia/md, j1 + sign(j2 − j1) × n × mdja ]
For instance between nodes (4,5) and (4,9) there are 3 intermediate nodes (4,6), (4,7) and (4,8); between (4,5) and (5,9) there are no in-between points, between (4,5) and (6,9) there is one node (5,7), between (4,5) and (8,9) there are 3 nodes (5,6), (6,7) and (7,8), etc.
The grid generator takes these rules into consideration, and will fix the positions of all intermediate nodes accordingly. For all nodes for which the coordinates have not been specified, the program will calculate the coordinates by means of a Laplace interpolation. This technique involves iterative updating of the nodal positions with following equations:
The iteration is stopped after differences in positions become smaller than a tolerance specified in the program. The result of the Laplace interpolator is that nodes are as much as possible evenly distributed between the nodes fixed by the user.
For nodes on the boundary, the equations are modified by excluding interior nodes. For instance for node (i, 1) the interpolation equations become:
This makes that the boundary is completely determined by an envelope of straight lines drawn along fixed boundary nodes. For instance, if only the four corner nodes are fixed, the grid boundaries will consist of the quadrilateral, defined by these four points. It is recommended that, as much as possible, boundary node coordinates are fixed in order to prevent unexpected mesh arrangements. However, there is always a trade-off between accuracy and mesh density, and this is free to decide by the user.
It is also recommended to satisfy some requirements on the shape of the elements. In general, the elements should as equiangular as possible, e.g. regular hexahedra. Highly distorted elements, e.g. long, thin quadrilaterals, can lead to numerical stability problems caused by round-off errors. In addition the convergence properties of the finite element method are proved in the limit of small element size only if the maximum angle in an element is bounded away from π. These requirements are modified for boundary layers, where highly stretched elements are desired and the finite element formulation allows to them. Even in this case this property is required.
The second part of the grid generation involves the z coordinates of the nodes. This is done in the following way. For any layer, the user can specify a number of points with x, y and z coordinates. The program uses an inverse distance interpolator to fix the z-values of all the nodes in that specific layer. For the remaining layers, the z-coordinates are linearly interpolated between fixed upper and lower layers. This technique only works when at least the top and base layers have been specified by the user.
GridGen uses the information specified in the data file .grd as explained in the next paragraph. GridGen produces two output files: filename.dat, filename.fem, and filename.plt. The data file filename.dat is a text file with basic information on the considered problem. The data file filename.fem is a text file containing all the information about the finite element mesh. Finally, filename.plt is equivalent to the previous one especially formatted for mesh visualization by Tecplot or Paraview software.
The input for the mesh generator has to be specified by the user in the file project.grd. This is a text file containing all information for fixing and calculating nodal coordinates. The data must be entered separated by blanks or commas in a sequence as explained below.
There are first three lines of titles or comments that will be used for identification of the active study.
On the fourth line, the number of columns, rows and layers is specified, i.e. ni, nj and nl.
On the next lines information concerning fixed x and y coordinates is given. Each line contains an index i, index j, coordinate x(i, j), and coordinate y(i, j), followed by 0 or 1. The zero implies that we are dealing with an individual fixed node, while the one means that the node is a starting point of a line segment, for which the end point is the node given on the next line. In this way, different line segments can be joined together. Obviously the last line should finish with a zero.
The next lines contain information about fixed z coordinates. First the layer number should be given for which z values will be specified. This is done with the statement LAYER followed by the number of the considered layer. On next lines, x, y, and z coordinates are given, one set on each line. This is repeated for every layer for which z coordinates are specified. The sequence should be in increasing layer number, this is from top to bottom.
For example, consider the data file given below:
First line title
Demonstration problem
Some simple situation
11 11 5
1 1 0. 0. 0
11 1 1000. 0. 0
11 11 1000. 1000. 0
1 11 0. 1000. 0
1 1 0. 0. 0
LAYER 1
500. 500. 50.
LAYER 5
500. 500. 10.
It defines an 11 by 11 nodal grid in the horizontal plane with 5 layers vertically. The mesh forms a square of 1000m in the horizontal plane. The top plane is at an elevation of 50m, while the bottom plane is at an elevation of 10m. With this information the grid generator will construct a 11 × 11 × 5 nodal grid (or 10 × 10 × 4 cells), with the nodes nicely distributed at equal distances of 100m in the horizontal plane, and with 10m intervals in the vertical direction.
Data file filename.dat is a text file created automatically by GridGen. This is a virgin version. It contains basic information about the problem case study, and will be used by other simulators such as GwMove. It contains the titles, number of nodes and cells, and number of columns, rows and layers of the finite element grid. Lateron other data will be edited by the user for running other software.
The virgin version of filename.dat that would result from the previous grid generation sample is shown below:
Title 1 : First line title
Title 2 : Demonstration problem
Title 3 : Some simple situation
Number of nodes : 605 11 11 5
Number of elements : 400
Number of soil types : 0
The file lists the titles that define the problem, the total number of nodes (nn = ni × nj × nl) and the dimensions of the grid (ni, nj, nl), the total number of elements (ne = (ni − 1) × (nj − 1) × (nl − 1)), the total number of soil types, which in the virgin version is zero (no soil types have been defined yet).
The generated grid is stored in a text data file called filename.fem. A free text file format is chosen because of reading and writing friendliness and platform portability. The file will be automatically processed by GwMove, whenever necessary.
The data on filename.fem consists of a list of nodes with their coordinates, and a list of the cells with their corner nodes, ordered in the expected direction as discussed previously. In addition for each element an integer code is specified, which for the time being is put to zero. Lateron this code will be used to identify the type of soil for each element. Hence, filename.fem contains the following data:
(x(i), y(i), z(i), i=1, nn)
(node(i,j), j=1, 8), kg(i), i=1, ne)
where:
- x(i), y(i), z(i) are respectively x, y and z coordinates of node i;
- nn is the total number of nodes;
- node(i,j) : 8 nodal numbers, that correspond to the corner nodes of element i, given in the correct order;
- kg(i) : code for each element;
- ne: total number of elements.
Notice that the coordinates are real values, while all other parameters are integers. Also notice, that the nodes and elements are listed in sequential order, one by one, without reference to the structure of the grid. This ordering is achieved as follows. The list starts with the first node with indexes (1,1,1) and proceeds along the columns along increasing i-index, when reaching the end of a column, the next row is taken, i.e increasing the j-index by one; when arriving at the last node of the last row, this is the last node of the first layer, the next layer is considered, and so on. The last node has number nn = ni × nj × nl.
The ordering of the elements is done in a similar fashion, first along the columns, than along the rows, and finally along the layers. Notice that the number of elements along the columns is ni − 1, along the rows nj − 1, and along the layers nl − 1, such that the total number of elements is ne = (ni − 1) × (nj − 1) × (nl − 1).
The code for each element, kg(i), is used lateron as a soil type information. At this stage all these numbers are equal to zero.