We will look at the development of the matrix structural analysis method for the simple case of a structure made only out of truss elements that can only deform in one direction. This will allow us to get a taste of how matrix structural analysis works without having to learn about all of the details and complexities that are present in beam and frame systems. As we will see, since we have only one-dimensional truss elements, each node in the system only has one possible degree of freedom (translation along the axis of the structural members).

In previous stiffness methods, each degree of freedom was dealt with separately. For example, in the slope-deflection method, we ended up with one equation for each degree of freedom. Then we had to solve for the unknown deflections and rotations by solving the system of equations.

In matrix structural analysis, we will end up with the same equations. But we will construct those equations using matrices that represent each element stiffness. Putting these stiffness matrices together, we will be able to construct a large matrix for the entire structure, that represents all of the equations that we had previously in the slope-deflection method. Since all of our equations will be in matrix form, we can take advantage of matrix methods to solve the system of equations and determine all of the unknown deflections and forces. Computers are well-adapted to solve such matrix problems.

The One-Dimensional Truss Element

Recall that the deformation of a truss element may be found using the following equation:

\begin{equation} \delta = \frac{FL}{EA} \tag{1} \end{equation}

where $\delta$ is the axial deformation, $F$ is the axial force in the truss element, $L$ is the length of the element, $E$ is the Young's modulus, and $A$ is the cross-sectional area of the element. This equation may be rearranged to find the following relationship between axial force and axial deformation:

\begin{equation} \boxed{ F = \left( \frac{EA}{L} \right) \delta } \label{eq:1D-Truss-Force} \tag{2} \end{equation}

The term that is multiplied by the deformation to get the force is the axial stiffness:

\begin{equation} \boxed{ k = \frac{EA}{L} } \label{eq:1D-Truss-Stiffness} \tag{3} \end{equation}

With this background, we can look at the behaviour of a one-dimensional truss element as shown in Figure 11.1. This truss element has a constant Young's modulus $E$ and cross-sectional area $A$. It has two ends, which we can consider to be connected to two separate nodes in our structure, one labelled '1' and one labelled '2' as shown in the figure. We can represent the complete behaviour of this entire element through the force and displacement of the two nodes. The force at node 1 is labelled $F_{x1}$ and the force at node two is labelled $F_{x2}$. Likewise, the displacement of node 1 (relative to its initial position) is labelled $\Delta_{x1}$ and the displacement of node two is labelled $\Delta_{x2}$.

We can relate the forces to the displacements at the two ends of the member using the stiffness term from equation \eqref{eq:1D-Truss-Stiffness}; however, we need to relate the end displacement to the bar deformation (these are not the same thing). For, example, if both the left and right sides move by 1.0 unit positive (to the right), then the entire bar moves to the right as a rigid body, neither expanding or contracting, so the deformation would be zero. The deformation can be related to the end node displacements as follows:

\begin{align} \delta = \Delta_{x2} - \Delta_{x1} \tag{4} \end{align}

So, the total internal axial force in the bar is equal to:

\begin{align} \boxed{ F = \left( \frac{EA}{L} \right) (\Delta_{x2} - \Delta_{x1}) } \label{eq:truss1D-int-force} \tag{5} \end{align}

Due to horizontal equilibrium, $F_{x1} = -F_{x2}$. The magnitude of these external forces is equal to the internal force in the truss element. If the internal force from equation \eqref{eq:truss1D-int-force} is positive, the bar is in tension, so the force on the left ($F_{x1}$) must point to the left (negative), and the force on the right ($F_{x2}$) must point to the right (positive). So:

\begin{align} F_{x1} &= -\left( \frac{EA}{L} \right) (\Delta_{x2} - \Delta_{x1}) \label{eq:Truss1D-Mat-Line1} \tag{6} \\ F_{x2} &= \left( \frac{EA}{L} \right) (\Delta_{x2} - \Delta_{x1}) \label{eq:Truss1D-Mat-Line2} \tag{7} \end{align}

These two equations define the force/deflection behaviour of the truss at both nodes simultaneously. They are only a function of displacements of the nodes (the nodal displacements) and the forces applied to the nodes (the nodal forces).

We can easily express these two equations in a matrix form as follows:

\begin{align} \begin{Bmatrix} F_{x1} \\ F_{x2} \end{Bmatrix} = \begin{bmatrix} \dfrac{EA}{L} & -\dfrac{EA}{L} \\[10pt] -\dfrac{EA}{L} & \dfrac{EA}{L} \end{bmatrix} \begin{Bmatrix} \Delta_{x1} \\ \Delta_{x2} \end{Bmatrix} \tag{8} \end{align}

where the matrix on the left of the equal sign is called the force vector, the large central matrix is called the stiffness matrix and the smaller matrix on the right with the displacements is called the displacement vector. The term vector just means a matrix with only one column. If we multiply the large central matrix by the vector of displacements on the right, we get:

\begin{align} F_{x1} = \left( \frac{EA}{L} \right) \Delta_{x1} + \left( -\frac{EA}{L} \right) \Delta_{x2} \tag{9} \\ F_{x2} = \left( -\frac{EA}{L} \right) \Delta_{x1} + \left( \frac{EA}{L} \right) \Delta_{x2} \tag{10} \end{align}

which are the same as equations \eqref{eq:Truss1D-Mat-Line1} and \eqref{eq:Truss1D-Mat-Line2}. This matrix equation constitutes a complete model for the behaviour of a one-dimensional truss element. The large matrix in the middle is called the stiffness matrix of the element because it contains all of the stiffness terms. It contains the most important information for the model, and it is useful to think about it as a separate element:

\begin{align} k = \begin{bmatrix} \dfrac{EA}{L} & -\dfrac{EA}{L} \\[10pt] -\dfrac{EA}{L} & \dfrac{EA}{L} \end{bmatrix} \tag{11} \\ k = \frac{EA}{L} \begin{bmatrix} 1 & -1 \\ -1 & 1 \end{bmatrix} \label{eq:1DTruss-Stiffness-Matrix} \tag{12} \end{align}

This is the stiffness matrix of a one-dimensional truss element. Other types of elements have different types of stiffness matrices. For a truss element in 2D space, we would need to take into account two extra degrees of freedom per node as well as the rotation of the element in space. Beam elements that include axial force and bending deformations are more complex still. For real physical systems, stiffness matrices are always square and symmetric about the diagonal axis of the matrix.

Another way to think about the construction of a stiffness matrix is to find the forces at either end of the element if the element experiences a unit deformation at each end (separately). This works because the stiffness is defined as the force per unit deformation. If we don't know the stiffness matrix, we can figure it out by first starting with the general form of the stiffness matrix for our element:

\begin{align} \begin{Bmatrix} F_{x1} \\ F_{x2} \end{Bmatrix} = \begin{bmatrix} k_{11} & k_{12} \\ k_{21} & k_{22} \end{bmatrix} \begin{Bmatrix} \Delta_{x1} \\ \Delta_{x2} \end{Bmatrix} \tag{13} \end{align}

where $k_11$, $k_12$, $k_21$ and $k_22$ are the individual terms within the stiffness matrix that we want to find. The first number in the subscript is the row in the matrix where the stiffness term is located, and the second number is the column in the matrix where it is located.

So let's individually set each displacement to 1.0 while setting the other to zero to calculate the stiffness terms. This process is shown in Figure 11.1. If we set $\Delta_{x1} = 1$ and $\Delta_{x2} = 0$, we get:

\begin{align} \begin{Bmatrix} F_{x1} \\ F_{x2} \end{Bmatrix} = \begin{bmatrix} k_{11} & k_{12} \\ k_{21} & k_{22} \end{bmatrix} \begin{Bmatrix} 1.0 \\ 0 \end{Bmatrix} \tag{14} \end{align}

Multiplying this out, we get:

\begin{align} F_{x1} &= k_{11}(1) + k_{12}(0) \tag{15} \\ F_{x2} &= k_{21}(1) + k_{22}(0) \tag{16} \end{align}

So, when $\Delta_{x1} = 1$ and $\Delta_{x2} = 0$, $F_{x1} = k{11}$ and $F_{x2} = k_{21}$. So, to find the stiffness terms $k_{11}$ and $k_{21}$, we just need to find out what the force is in the truss element at each end when $\Delta_{x1} = 1$ and $\Delta_{x2} = 0$. This situation is shown in the middle of Figure 11.1. When the left side of the truss moves to the right by 1.0 and the right side remains in the same place, the truss element is in compression with a total deformation $\delta = -1.0$. This means that the force at the left end of the bar is:

\begin{align} F_{x1} = \left( \frac{EA}{L} \right) (1) \tag{17} \end{align}

which is positive because it points to the right for compression, as shown in the figure. The force at the right end of the bar is:

\begin{align} F_{x2} = -\left( \frac{EA}{L} \right) (1) \tag{18} \end{align}

which is negative because it points to the left for compression, as shown in the figure. This means that:

\begin{align} k_{11} = F_{x1} = \frac{EA}{L} \tag{19} \\ k_{21} = F_{x2} = -\frac{EA}{L} \tag{20} \end{align}

If we switch the displacements and set $\Delta_{x1} = 0$ and $\Delta_{x2} = 1$, we get:

\begin{align} \begin{Bmatrix} F_{x1} \\ F_{x2} \end{Bmatrix} = \begin{bmatrix} k_{11} & k_{12} \\ k_{21} & k_{22} \end{bmatrix} \begin{Bmatrix} 0 \\ 1.0 \end{Bmatrix} \tag{21} \end{align}

Multiplying this out, we get:

\begin{align} F_{x1} &= k_{11}(0) + k_{12}(1) \tag{22}  \\ F_{x2} &= k_{21}(0) + k_{22}(1) \tag{23} \end{align}

So, when $\Delta_{x1} = 0$ and $\Delta_{x2} = 1$, $F_{x1} = k{12}$ and $F_{x2} = k_{22}$. This situation is shown in the lower diagram in Figure 11.1. When the right side of the truss moves to the right by 1.0 and the left side remains in the same place, the truss element is in tension with a total deformation $\delta = 1.0$. This means that the force at the left end of the bar is:

\begin{align} F_{x1} = -\left( \frac{EA}{L} \right) (1) \tag{24} \end{align}

which is negative because it points to the left for tension, as shown in the figure. The force at the right end of the bar is:

\begin{align} F_{x2} = \left( \frac{EA}{L} \right) (1) \tag{25} \end{align}

which is positive because it points to the right for tension, as shown in the figure. This means that:

\begin{align} k_{12} = F_{x1} = -\frac{EA}{L} \tag{26} \\ k_{22} = F_{x2} = \frac{EA}{L} \tag{27} \end{align}

If we now take all of these solved stiffness terms and construct the stiffness matrix of the element, we get:

\begin{align} k = \begin{bmatrix} \dfrac{EA}{L} & -\dfrac{EA}{L} \\[10pt] -\dfrac{EA}{L} & \dfrac{EA}{L} \end{bmatrix} \tag{28} \end{align}

which is the same stiffness matrix that we derived previously in equation \eqref{eq:1DTruss-Stiffness-Matrix}.

Assembling the Full Stiffness Matrix

If we have a structural analysis problem with multiple one-dimensional truss elements, we must first define the stiffness matrices for each individual element as described in the previous section. After we define the stiffness matrix for each element, we must combine all of the elements together to form on global stiffness matrix for the entire problem. This matrix defines all of the interconnections between the elements and includes all of the information related to the stiffness of each element for each degree-of-freedom.

The resulting global stiffness matrix is put into an equation with the global nodal force vector (which contains all of the forces for each node in each DOF) and the global nodal displacement vector (which contains all of the displacements of each node in each DOF) to get a global system of equations for the entire problem with the following form:

\begin{align} \begin{Bmatrix} F_1 \\ F_2 \\ F_3 \\ \vdots \\ F_n \end{Bmatrix} = \begin{bmatrix} k_{11} & k_{12} & k_{13} & \cdots & k_{1n} \\ k_{21} & k_{22} & k_{23} & \cdots & k_{2n} \\ k_{31} & k_{32} & k_{33} & \cdots & k_{3n} \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ k_{n1} & k_{n2} & k_{n3} & \cdots & k_{nn} \end{bmatrix} \begin{Bmatrix} \Delta_{1} \\ \Delta_{2} \\ \Delta_{3} \\ \vdots \\ \Delta_{n} \end{Bmatrix} \label{eq:truss1D-Full-System} \tag{29} \end{align}

where $F_i$ is the external force on node $i$, $k_{ij}$ is the global stiffness matrix term for the force on node $i$ needed to cause a unit displacement at node $j$, and $\Delta_j$ is the displacement at node $j$. This is for a structure that has $n$ nodes, where you need one row for each nodal DOF. This global stiffness matrix is made by assembling the individual stiffness matrices for each element connected at each node. For example, an element that is connected to nodes 3 and 6 will contribute its own local $k_{11}$ term to the global stiffness matrix's $k_{33}$ term. Likewise, it will contribute its own $k_{12}$ term to the global stiffness matrix's $k_{36}$ term, $k_{21}$ to $k_{63}$ and $k_{22}$ to $k_{66}$. Each element stiffness matrix is added to the global stiffness matrix in this way. Multiple contributions in a single node are added together. This process will be demonstrated using an example.

Solving for Nodal Forces and Displacements

Once the stiffness matrix is formed, the full system of equations in the form shown in equation \eqref{eq:truss1D-Full-System} may be solved. At each nodal DOF (each row), we must either know the external force or the nodal deflection. Then, we can solve only those rows where we don't know the deflection. Once we have all of the nodal deflections, we can solve for the nodal forces. Then, using the individual element stiffness matrices, we can solve for the internal force in each element.

There are numerous different computer algorithms that may be used to solve the matrix of equations, but these are outside the scope of this book.

Example

The full process for a matrix structural analysis for a one dimensional truss will be demonstrated using the simple example shown in Figure 11.2. This is a one dimensional structure, meaning that all of the nodes are only permitted to move in one direction. This structure consists of four different truss elements which are numbered one through four as shown in the figure. These elements are connected at four different nodes, also numbered one through four as shown. The left end of the structure, at node 1 is restrained and cannot move. Nodes 2 and 4 have external loads, and Node 3 has an imposed displacement of $13\mathrm{\,mm}$ to the right (positive).

The truss elements in Figure 11.2 are made of one of two different materials, with Young's modulus of either $E =9000\mathrm{\,MPa}$ or $E = 900\mathrm{\,MPa}$. These are labelled in the figure and are shaded differently as shown. Each element also has its own different cross-sectional area $A$ as shown.

The first step in this analysis is to determine the stiffness matrix for each individual element in the structure. Since these are all one-dimensional truss members, we can use equation \eqref{eq:1DTruss-Stiffness-Matrix}. For element 1 (connected to nodes 1 and 2):

\begin{align*} k_1 &= \frac{E_1 A_1}{L_1} \begin{bmatrix} 1 & -1 \\ -1 & 1 \end{bmatrix} \\ k_1 &= \frac{9000 (50)}{4000} \begin{bmatrix} 1 & -1 \\ -1 & 1 \end{bmatrix} \\ k_1 &= 112.5\mathrm{\,N/mm} \begin{bmatrix} 1 & -1 \\ -1 & 1 \end{bmatrix} \end{align*}

For element 2 (connected to nodes 2 and 3):

\begin{align*} k_2 = \frac{9000 (50)}{5000} \begin{bmatrix} 1 & -1 \\ -1 & 1 \end{bmatrix} = 90.0\mathrm{\,N/mm} \begin{bmatrix} 1 & -1 \\ -1 & 1 \end{bmatrix} \end{align*}

For element 3 (connected to nodes 2 and 4):

\begin{align*} k_3 = \frac{9000 (90)}{8000} \begin{bmatrix} 1 & -1 \\ -1 & 1 \end{bmatrix} = 101.2\mathrm{\,N/mm} \begin{bmatrix} 1 & -1 \\ -1 & 1 \end{bmatrix} \end{align*}

For element 4 (connected to nodes 3 and 4):

\begin{align*} k_4 = \frac{900 (120)}{3000} \begin{bmatrix} 1 & -1 \\ -1 & 1 \end{bmatrix} = 36.0\mathrm{\,N/mm} \begin{bmatrix} 1 & -1 \\ -1 & 1 \end{bmatrix} \end{align*}

Now we can form the global stiffness matrix based on these individual stiffness matrices for each element and the connected node locations for each. For example, element 3, which is connected to nodes 2 and four will contribute stiffness terms to elements 22, 24, 42, and 44 of the global stiffness matrix. The contributions to the global stiffness matrix $[k]$ from each element stiffness look like this:

\begin{align*} [k] = \begin{bmatrix} k1 & -k1 & & \\ -k1 & k1 + k2 + k3 & -k2 & -k3 \\ & -k2 & k2 + k4 & -k4 \\ & -k3 & -k4 & k3 + k4 \end{bmatrix} \end{align*}

Putting in the actual stiffness values from each element, we get:

\begin{align*} [k] = \begin{bmatrix} 112.5 & -112.5 & 0 & 0 \\ -112.5 & 112.5 + 90.0 + 101.2 & -90.0 & -101.2 \\ 0 & -90.0 & 90.0+ 36.0 & -36.0 \\ 0 & -101.2 & -36.0 & 36.0 + 101.2 \end{bmatrix} \end{align*}

and simplifying:

\begin{align*} [k] = \begin{bmatrix} 112.5 & -112.5 & 0 & 0 \\ -112.5 & 303.7 & -90.0 & -101.2 \\ 0 & -90.0 & 126.0& -36.0 \\ 0 & -101.2 & -36.0 & 137.2 \end{bmatrix} \end{align*}

Using this global stiffness matrix, we can now look at the entire system of equations for the entire structure:

\begin{align*} \lbrace F \rbrace &= [k] \lbrace \Delta \rbrace \\ \begin{Bmatrix} F_{1} \\ F_{2} \\ F_{3} \\ F_{4} \end{Bmatrix} &= \begin{bmatrix} 112.5 & -112.5 & 0 & 0 \\ -112.5 & 303.7 & -90.0 & -101.2 \\ 0 & -90.0 & 126.0& -36.0 \\ 0 & -101.2 & -36.0 & 137.2 \end{bmatrix} \begin{Bmatrix} \Delta_{1} \\ \Delta_{2} \\ \Delta_{3} \\ \Delta_{4} \end{Bmatrix} \end{align*}

For an analysis problem to be fully-defined, for each node we always know either the external force on that node, or the deflection of the node. Knowing the force may just mean that we know that the external force is zero on a node, but we don't know the displacement. For this problem, as shown in Figure 11.2, we know that the external force at node 2 is $-350\mathrm{\,N}$ ($F_{2} = -350$) and that the external force at node 4 is $+1100\mathrm{\,N}$ ($F_{4} = 1100$). We know that the displacement at node 1 is zero because it is fixed ($\Delta_{1} = 0$). We also know that there is an imposed displacement at node 3 of $13\mathrm{\,mm}$ ($\Delta_{3} = 13$). This node is forced to move exactly $13\mathrm{\,mm}$. If we weren't given an imposed displacement at node 3 and that node was free to move, then we would know instead that the external force at node 3 is zero. Putting this information into our system of equations, we get:

\begin{align*} \begin{Bmatrix} F_{1} \\ -350 \\ F_{3} \\ 1100 \end{Bmatrix} &= \begin{bmatrix} 112.5 & -112.5 & 0 & 0 \\ -112.5 & 303.7 & -90.0 & -101.2 \\ 0 & -90.0 & 126.0& -36.0 \\ 0 & -101.2 & -36.0 & 137.2 \end{bmatrix} \begin{Bmatrix} 0 \\ \Delta_{2} \\ 13 \\ \Delta_{4} \end{Bmatrix} \end{align*}

This is a system of four equations and four unknowns. This is a system that is easily solved using a computer. For this example, since there are only two free displacement degrees of freedom, we can expand the second and fourth rows (the second and fourth equations) to get:

\begin{align*} -350 &= -112.5 (0) + 303.7(\Delta_{2}) -90.0(13) -101.2(\Delta_{4}) \\ 1100 &= 0 (0) - 101.2(\Delta_{2}) -36.0(13) +137.2(\Delta_{4}) \end{align*}

Solving this system of two equations and two unknowns, we get:

\begin{align*} \Delta_{2} &= 17.79\mathrm{\,mm} \\ \Delta_{4} &= 8.62\mathrm{\,mm} \end{align*}

This solution suggests that both nodes 2 and 4 move towards the right, which makes sense based on the system shown in Figure 11.2.

After solving the displacements at nodes 2 and 4, we now know the displacements at all of the nodes. We can now easily multiply through the first and third rows of the system of equations to get:

\begin{align*} F_{1} &= -970\mathrm{\,N} \\ F_{3} &= +222\mathrm{\,N} \end{align*}

These two forces are equivalent to the reaction forces at the fixed end and the imposed displacement location.

The complete solution for the external forces and displacements of this one-dimensional truss is shown in Figure 11.3.

To find the internal forces in individual elements, we can take the global nodal displacements and use them with the original element stiffness matrices. For element 1:

\begin{align*} \begin{Bmatrix} F_{x1} \\ F_{x2} \end{Bmatrix} &= 112.5\mathrm{\,N/mm} \begin{bmatrix} 1 & -1 \\ -1 & 1 \end{bmatrix} \begin{Bmatrix} \Delta_{x1} \\ \Delta_{x2} \end{Bmatrix} \end{align*}

where $F_{x1}$ and $F_{x1}$ are the local forces at nodes 1 and 2 on element 1, and $\Delta_{x1}$ and $\Delta_{x2}$ are the local displacements of nodes 1 and 2 for element 1. The local displacements are the same as the global displacements, so:

\begin{align*} \Delta_{x1} &= \Delta_1 = 0 \\ \Delta_{x2} &= \Delta_2 = 8.62 \end{align*}

Giving us:

\begin{align*} \begin{Bmatrix} F_{x1} \\ F_{x2} \end{Bmatrix} &= \begin{bmatrix} 112.5 & -112.5 \\ -112.5 & 112.5 \end{bmatrix} \begin{Bmatrix} 0 \\ 8.62 \end{Bmatrix} \end{align*}

Multiplying through we can find the forces at each end of element 1:

\begin{align*} F_{x1} &= 112.5(0) - 112.5 (8.62) \\ F_{x1} &= -970\mathrm{\,N} \\ F_{x2} &= -112.5(0) + 112.5 (8.62) \\ F_{x2} &= 970\mathrm{\,N} \end{align*}

The forces at either end of truss element 1 are equal and opposite, as we would expect.

This process may be repeated for the other elements to get the internal axial force in every one-dimensional truss element.