3.2 Installation Instructions
4.1.1 Pre-Processing
Point List
Throughout the calculation process it is important to organize all the variables. Misplaced deformations or boundary conditions in relation to the stiffness matrix will result in erro-neous results. Therefore, the first step will be to create a list of all the points (i.e. nodes) from the input geometry of lines. The important thing to note about the point list is that no point occurs twice. This is done deliberately so that the point list will be the ”model order” for assembling the global stiffness matrix from the element stiffness matrices in a later procedure.
Because of lower accuracy in Grasshopper than C#, the input points are only accurate to a certain degree, and tend to have strange numbering for decimals placed after 10−6. Since the Cholesky solve method outlined in Ch. 4.1.2 requires a symmetric matrix, they need to be rounded to stave off errors caused by this phenomenon. The process of creating the point list is presented in Lst. 4.1.
1 List<Point3d> points = new List<Point3d>();
2 foreach (Line line in geometry)
3 Point3d tempFrom = new Point3d(Math.Round(line.From.X, 5), Math.Round(line.From.Y, 5), Math.Round(line.From.Z, 5));
4 Point3d tempTo = new Point3d(Math.Round(line.To.X, 5), Math.Round(line.To.Y, 5), Math.Round(line.To.Z, 5));
5 //adds point unless it already exists in pointlist
6 if (!points.Contains(tempFrom))
7 points.Add(tempFrom);
8 if (!points.Contains(tempTo))
9 points.Add(tempTo);
Listing 4.1: Method CreatePointList for 3D Truss
Having identical points occur more than once could disturb the stiffness relations in the global stiffness matrix and add more equations to the linear system, this would be unnecessary and may cause error in the computation process. The method for creating the point list will therefore skip any point if it already exists in the point list, and add it otherwise. The index for each unique point in the point list will thereafter act as the identifier for each point. It is of no consequence in which order the point occur as long as all points are unique and stays in the same order throughout the computation.
Boundary Conditions
With the arbitrary order of points established, the list of boundary conditions can be constructed. Using theBDC Trusscomponent described in Ch. 4.2, the boundary conditions are given as a list of strings with the format ”x,y,z:fx,fy,fz” as shown on Fig. 4.3. The x,y and z represents the coordinates of a point in three dimensional Cartesian coordinates, given in millimetres. The field ”fx” can be interpreted as the question ”free x?”, and takes the form of an integer, 1 or 0, representing respectivelytrueorfalse, the logic is similar for fy and fz.
An example of how the input is formatted for two nodes is shown in Fig. 4.3, where the nodes coordinates are as given in Eq. 4.1.1 below.
(x, y, z)1= (0,0,0) (x, y, z)2= (2000,0,0) (Eq. 4.1.1)
Figure 4.3: BDC string format.
Here Node 1 is clamped in x-,y-,z-direction (notice the fx,fy,fz = 0,0,0), while Node 2 is clamped in y- and z-direction, but free to move in x-direction.
The information about the points in the inputted string list is used to arrange the boundary conditions according to the order from point list. Thereafter, the boundary conditions is stored as a list with the true/false values for fx, fy and fz separated, resulting in a list with three entries for each point in the point list. For all the other points besides the boundaries, the condition is set to 1, which means it is free to move.
In the 2D Truss software, the fy value is disregarded since the calculations are only performed for two dimensions. It can be specified for testing reasons but it will be dis-regarded in the calculation component. Note that the y axis has been disdis-regarded, which means 2D Truss works in the x and z axes.
Loads
Similarly to the boundary condition input, the supplied load list is a list of strings. The strings is formatted as ”x,y,z:vx,vy,vz”, where x, y and z is the coordinates of the loaded point, and vx, vy and vz is the vector components in hence x-, y- and z-direction. Each vector component has the value of the force in the respective direction, hence the complete vector contains information about direction and the force magnitude. The strings are decoded and transformed into a list of doubles, in much the same manner as for the boundary conditions. The respective force in each direction is set separately and adhering to the ordering given in the point list previously created. For 3D Truss, this results in a list thrice the length of the number of points, where all points without loads are set to zero. For 2D Truss, the length is twice the number of all points, as only two directions are considered.
Lst. 4.2 shows how the method CreateLoadList in 3D Truss parses the text input from the load component and stores them as a List of doubles.
1 for (int i = 0; i < loadtxt.Count; i++)
2 string coordstr = (loadtxt[i].Split(’:’)[0]);
3 string loadstr = (loadtxt[i].Split(’:’)[1]);
4
5 string[] coordstr1 = (coordstr.Split(’,’));
6 string[] loadstr1 = (loadstr.Split(’,’));
7
8 inputLoads.Add(Math.Round(double.Parse(loadstr1[0]), 5));
9 inputLoads.Add(Math.Round(double.Parse(loadstr1[1]), 5));
10 inputLoads.Add(Math.Round(double.Parse(loadstr1[2]), 5));
11
12 coordlist.Add(new Point3d(Math.Round(double.Parse(coordstr1[0]), 5), Math.Round(double.Parse(coordstr1[1]), 5),
Math.Round(double.Parse(coordstr1[2]), 5)));
13
14 //place at load at correct entry in global load list
15 foreach (Point3d point in coordlist)
16 int i = points.IndexOf(point);
17 int j = coordlist.IndexOf(point);
18 loads[i * 3 + 0] = inputLoads[j * 3 + 0];
19 loads[i * 3 + 1] = inputLoads[j * 3 + 1];
20 loads[i * 3 + 2] = inputLoads[j * 3 + 2];
Listing 4.2: Excerpt of method CreateLoadList for 3D Truss
Stiffness Matrices
The element stiffness matrices are established based on the geometry (List of lines), point list (List of points), Young’s Modulus E and cross-sectional area A. The E-modulus and area A has been assumed to apply to all elements. The global stiffness matrix is later assembled by inserting the values from each element stiffness matrix (i.e. the element stiffness matrix of each bar) according to their node numbering, and thus connecting all the elements.
Element Stiffness Matrix
The element stiffness matrix in global coordinatesKeis different for 2D and 3D Truss.
However, the local element stiffness matriceskeare (almost) identical for both, and is similar to Eq. 2.5.59. Since the stiffness matrix from Eq. 2.5.59 is defined for one dimension, and trusses are for two and three dimensions, there is a gradual increase in the matrix size.
The local element stiffness matrix for a 2D truss then becomes
ke2DT russ=
The global element stiffness matricesKefor 2D and 3D becomes rather different as they are multiplied by different transformation matrices. In the 2D case, the transformation matrix will be similar to the one in Eq. 2.7.13, without the rotational dof, and reads
T2DT russ=
By applying Eq. 2.7.6, this results in
Ke2DT russ= EA
Solving forKedirectly as in Eq. 4.1.4 is faster and simpler than first establishingkein local coordinates and then transforming to global coordinates. By this reason the 2D Truss component skips this transformation step and implementsKe. For larger matrices like in 3D Truss, 3D Beam and Shell, the transformation procedure of Eq. 2.7.6 is followed instead, choosing to prioritize readability rather than optimize for time, the time usage of this process will be investigated further in later chapters. In three dimension the only difference between from the 2D local element stiffness matrix is two added rows and column of zeroes, this becomes
Ke3DT russ= EA
The transformation matrixT3DT russfor the 3D Truss elements is found by assembling the directional cosines from Eq. 2.7.23-2.7.24 for each node and is assembled as
T3DT russ=
The resultingKe3DT russwill then be calculated according to Eq. 2.7.6 for each element as the coordinates are obtained.
As an example of how the transformation matrix is used, the 2D truss element matrix from Eq. 4.1.4 is filled for element 1 on Fig. 4.4. The coordinates for the node in the figure is presented in Tab. 4.1.
Figure 4.4: 2D Truss with red indices for elements and white indices for nodes.
Table 4.1: Example 2D Truss. Nodal coordinates for Fig. 4.4.
Node index X-coord Y-coord Z-coord
0 0.0 0.0 0.0
1 1000.0 0.0 1000.0
2 2000.0 0.0 0.0
The angleθis found by taking the arctangent of∆z∆x. θ= arctan 0−0
2000−0 (Eq. 4.1.7)
= arctan 0
2000 (Eq. 4.1.8)
= 0° (Eq. 4.1.9)
The angle is then inserted into the abbreviated sine and cosine from Eq. 2.7.12.
c= cos 0°= 1 s= sin 0°= 0 (Eq. 4.1.10) Material properties are set as
Inserting values from Eq. 4.1.10 to 4.1.11 intoKe2DT russfrom Eq. 4.1.4 results in the complete element stiffness matrix for element 1, in global coordinates.
Ke12DT russ= 210GP a·10000mm2
A remark on the example is that the global element stiffness matrix looks exactly like in Eq. 4.1.4 because it is oriented horizontally. The numbers would be ”messier” for the other diagonal elements.
Global Stiffness Matrix
The element stiffness matrix for bar element 1 is now the 4x4 matrix from Eq. 4.1.13 and next step is to add it to the global stiffness matrix. The element stiffness matrix can be divided into four 2x2 matrices: upper left corner, upper right corner, lower left corner and lower right corner. The placement of each 2x2 matrix in the element stiffness matrix is dependent on which index the start node and end node has in the point list. The four sub-matrices is illustrated in Eq. 4.1.14 with the respective node relation.
K12DT russ= 10.5·105 N
This can be automated by use of for-looping like in Lst. 4.3. Note that this double for-loop is for 3D Truss, while the 2D Truss utilized a more direct placement method as there were just a few term to place. As the element stiffness matrices grow, this process become more comprehensive and a double for-loop seemed like the organized way to accomplish this task.
1 for (int row = 0; row < K_elem.RowCount / ldof; row++)
2 for (int col = 0; col < K_elem.ColumnCount / ldof; col++)
3 //top left 3x3 of K-element matrix
4 K_G[nIndex1 * 3 + row, nIndex1 * 3 + col] += K_elem[row, col];
5 //top right 3x3 of K-element matrix
6 K_G[nIndex1 * 3 + row, node2 * 3 + col] += K_elem[row, col + 3];
7 //bottom left 3x3 of K-element matrix
8 K_G[node2 * 3 + row, nIndex1 * 3 + col] += K_elem[row + 3, col];
9 //bottom right 3x3 of K-element matrix
10 K_G[node2 * 3 + row, node2 * 3 + col] += K_elem[row + 3, col + 3];
Listing 4.3: An automated process forKeplacement intoKGfor the 3D Truss software
For the 2D Truss from Fig. 4.4, element 1 begins at node 0 and ends at node 2, this results in a placement in the global stiffness matrix as shown in Eq. 4.1.15-4.1.16.
KG1 =
An identical procedure is performed for every element, normally starting with element 0 and ending with the last element. A complete global stiffness matrixKGfor Fig. 4.4 will look like in Eq. 4.1.17.
The global stiffness matrix can be preallocated with a dimension of NxN by using the number of nodes multiplied with the number of local degrees of freedom (ldofs). For the 2D Truss there will always be 2 ldofs per node, in respectively x- and z-direction. While for 3D Truss there are 3 ldofs, in respectively x-, y- and z-direction. For the 2D Truss example shown in Fig. 4.4 there are three unique nodes, which gives a global stiffness matrix of 6x6
KG=
When assembled the constant global stiffness matrix can describe the linear static behaviour of the structure, by providing a relation between forces and deformations.
Reduced Global Stiffness Matrix and Reduced Load List
After the global stiffness matrix has been established, the reduced global stiffness matrix (KGr-matrix) must be constructed. In order to create theKGr-matrix, the clamped boundary conditions are removed. It is also necessary to reduce the load list equivalently so that it can be used as the right-hand-side (RHS) of the equation when solving for displacements.
The process of solving is more thoroughly explained in Ch. 2.3 and Ch. 2.9.
For every entry in the boundary list containing zeros, the corresponding index for rows and columns inKGand load list is removed, as illustrated in Eq. 4.1.18. Numbers in black will remain in the new list, while the rest (grey) numbers are removed.
The RHS of Eq. 4.1.18 will be the reaction and loading forces, where the reaction forces at this time are still unknowns and are therefore set as zeroes in the software. The equations with reaction forces as the RHS can not be solved without more information and is therefore not taken into the reduced stiffness matrix. The equations involving reaction forces however, does not relate to any deformations as they relates to clamped dofs, which are not movable. They can therefore be removed and the system can be solved without any complications.
The resultingKGr-matrix becomes as shown in Eq. 4.1.19, along with the reduced load list.
105
While the 3D Truss makes use of inbuilt functions for lists and matrices to remove rows and columns by index, the 2D Truss instead builds the reduced stiffness matrix (and reduced load list) from scratch by adding all the entries which relates to free dof in theKGmatrix and load list. As will be explained in Ch. 5.1.1, the method for 2D Truss is actually superior to the ”improved” 3D Truss method in terms of runtime.
The difference in the algorithms for reducing the global stiffness matrix in 2D Truss and 3D Truss can be seen from respectively Lst. 4.4 and Lst. 4.5. Notice how the 3D Truss reducing method seem simpler because of the methodsRemoveRowandRemoveColumn.
But in fact it is noticeable slower than the 2D Truss method when the matrices grow quite large.
1 int dofs_red = points.Count * 2 - (bdc_value.Count - bdc_value.Sum());
2 double[,] K_redu = new double[dofs_red, dofs_red];
3 List<double> load_redu = new List<double>();
4 List<int> bdc_red = new List<int>();
5 int m = 0;
6 for (int i = 0; i < K_tot.GetLength(0); i++)
7 if (bdc_value[i] == 1)
8 int n = 0;
9 for (int j = 0; j < K_tot.GetLength(1); j++)
10 if (bdc_value[j] == 1)
11 K_redu[m, n] = K_tot[i, j];
12 n++;
13 load_redu.Add(load[i]);
14 m++;
Listing 4.4: CreateReducedGlobalStiffnessMatrix method for 2D Truss
1 K_red = Matrix<double>.Build.SparseOfMatrix(K);
2 List<double> load_redu = new List<double>(load);
3 for (int i = 0, j = 0; i < load.Count; i++)
4 if (bdc_value[i] == 0)
5 K_red = K_red.RemoveRow(i - j);
6 K_red = K_red.RemoveColumn(i - j);
7 load_redu.RemoveAt(i - j);
8 j++;
Listing 4.5: CreateReducedGlobalStiffnessMatrix method for 3D Truss