Topics in Computational Mechanics: Part 1

Over the next few posts, I will share the writings and programs that I had created for a class in computational structural analysis—how we efficiently analyze a structure using numerical methods. The field is a subset of computational mechanics, which combines the disciplines of mathematics, computer science, and engineering.

The class consisted of juniors and seniors in mechanical or aerospace engineering. In other words, for brevity, I will assume that you have some familiarity with statics, linear algebra, and calculus. If time permits in future, I will upload my class notes to fill gaps and share more drawings.

1. Definitions

We will be looking at bars, beams, trusses, and frames. These are structures that we classify according to what they can and cannot do.

A bar represents a 1D structural element that can withstand only axial forces—tensile or compressive forces directed along the element’s local (reference) axis. We assume that a bar cannot carry other types of internal forces, such as shear force, bending moment, and twisting moment. A bar is the simplest 1D element, but it can also show complex behavior such as buckling.

A beam represents a 1D element that can withstand not only axial forces (tensile or compressive), but also an internal force that involves transverse shear force, bending moment, and twisting moment. As you can imagine, a beam is more complex than a bar, both physically and mathematically.

We can join 1D elements such as bars and beams to form a complex 3D structure. A truss is a structure where the predominant behavior of each element is axial force transmission (bar-like). If some of the elements in a structure carry loads such that internal forces and/or internal moments exist (beam-like), then we refer to the structure as a frame.

2. Truss Examples

What better way to learn than solving problems by hand?

I had created 3 problems for my recitation. Each is challenging as we allow arbitrary parameters that meet small-strain constraint. The last two problems also mark an important lesson in computational mechanics: whenever possible, choose an easy coordinate system and use symmetry to simplify a problem.

a. Assembly of truss

Consider a truss with two bars, one of length 5L and the other 4L. The bars are joined by a roller at node A, while their other ends are fixed, as shown below:

A two-bar truss with horizontal displacement and vertical load applied
A two-bar truss with horizontal displacement and vertical load applied to node A.

A force of P is applied to node A in the negative y direction, while node A moves by an amount of 0.05L in the positive x direction. Both bars are assumed to be linearly elastic and isotropic, with Young’s modulus E and cross-sectional area A.

Determine the y displacement of node A and the axial force in each bar. For a model problem, we can use the following parameters:

\begin{array}{l} P = 10^{6}\,\mbox{N} \\[16pt] L = 1\,\mbox{m} \\[16pt] E = 210\,\mbox{GPa}\,\mbox{(steel)} \\[16pt] A = 5 \times 10^{-4}\,\mbox{m}^{2}\,\mbox{(}5\,\mbox{cm}^{2}\mbox{)} \end{array}

i. Notations

First, let’s agree on our notations. We will write the displacements of the nodes (also called degrees of freedom) as,

u_{A,x},\,\,\,u_{A,y},\,\,\,u_{B,x},\,\,\,u_{B,y},\,\,\,u_{C,x},\,\,\,u_{C,y}

the bars’ internal forces as,

f_{1},\,\,\,f_{2}

and the reaction forces as,

f_{A,x},\,\,\,f_{A,y},\,\,\,f_{B,x},\,\,\,f_{B,y},\,\,\,f_{C,x},\,\,\,f_{C,y}

ii. Equilibrium equations

The net forces along x and y directions must be 0. To help us derive the equilibrium equations, we will draw free-body diagrams for each node.

Free-body diagram for node A
Free-body diagram for node A
Free-body diagram for node B
Free-body diagram for node B
Free-body diagram for node C
Free-body diagram for node C

For node A, we get the following system of equations:

\begin{array}{rcl} \displaystyle -\frac{3}{5}f_{1} \,+\, f_{A,x} & = & 0 \\[16pt] \displaystyle \frac{4}{5}f_{1} \,+\, f_{2} \,+\, f_{A,y} \,-\, P & = & 0 \end{array}

For node B,

\begin{array}{rcl} \displaystyle \frac{3}{5}f_{1} \,+\, f_{B,x} & = & 0 \\[16pt] \displaystyle -\frac{4}{5}f_{1} \,+\, f_{B,y} & = & 0 \end{array}

Finally, for node C,

\begin{array}{rcl} \displaystyle f_{C,x} & = & 0 \\[16pt] \displaystyle -f_{2} \,+\, f_{C,y} & = & 0 \end{array}

We can combine these equations into a single matrix equation:

\left[\begin{array}{cc} \displaystyle \frac{3}{5} & 0 \\[16pt] \displaystyle -\frac{4}{5} & -1 \\[16pt] \displaystyle -\frac{3}{5} & 0 \\[16pt] \displaystyle \frac{4}{5} & 0 \\[16pt] 0 & 0 \\[16pt] 0 & 1 \end{array}\right] \left[\begin{array}{c} f_{1} \\[16pt] f_{2} \end{array}\right] = \left[\begin{array}{c} f_{A,x} \\[16pt] f_{A,y} - P \\[16pt] f_{B,x} \\[16pt] f_{B,y} \\[16pt] f_{C,x} \\[16pt] f_{C,y} \end{array}\right]

ii. Stress-strain equations (constitutive equations)

Since we are working with small-strain problems, the relation between stress and strain is assumed to be linear:

\left[\begin{array}{c} f_{1} \\[16pt] f_{2} \end{array}\right] = EA \left[\begin{array}{c} \varepsilon_{1} \\[16pt] \varepsilon_{2} \end{array}\right]

Note, \varepsilon_{1} and \varepsilon_{2} denote the axial strains.

iii. Strain-displacement equations

For bar 1,

\begin{aligned} \varepsilon_{1} &\,=\, \displaystyle\frac{u_{A,x} - u_{B,x}}{5L} \cdot \frac{3}{5} + \frac{u_{A,y} - u_{B,y}}{5L} \cdot -\frac{4}{5} \\[16pt] &\,=\, \displaystyle\frac{3}{25L}\bigl(u_{A,x} - u_{B,x}\bigr) - \frac{4}{25L}\bigl(u_{A,y} - u_{B,y}\bigr) \end{aligned}

For bar 2,

\begin{aligned} \varepsilon_{2} &\,=\, \displaystyle\frac{u_{A,x} - u_{C,x}}{4L} \cdot 0 + \frac{u_{A,y} - u_{C,y}}{4L} \cdot -1 \\[16pt] &\,=\, \displaystyle-\frac{1}{4L}\bigl(u_{A,y} - u_{C,y}\bigr) \end{aligned}

Therefore,

\displaystyle\left[\begin{array}{c} \varepsilon_{1} \\[16pt] \varepsilon_{2} \end{array}\right] = \frac{1}{L}\left[\begin{array}{cccccc} \displaystyle\frac{3}{25} & \displaystyle-\frac{4}{25} & \displaystyle-\frac{3}{25} & \displaystyle\frac{4}{25} & 0 & 0 \\[16pt] 0 & \displaystyle-\frac{1}{4} & 0 & 0 & 0 & \displaystyle\frac{1}{4} \end{array}\right] \left[\begin{array}{c} u_{A,x} \\[16pt] u_{A,y} \\[16pt] u_{B,x} \\[16pt] u_{B,y} \\[16pt] u_{C,x} \\[16pt] u_{C,y} \end{array}\right]

iv. Stiffness equation

We combine the equilibrium, stress-strain, and strain-displacement equations to arrive at a single matrix equation:

\underbrace{\displaystyle\frac{EA}{L}\left[\begin{array}{cccccc} \displaystyle\frac{9}{125} & \displaystyle-\frac{12}{125} & \displaystyle-\frac{9}{125} & \displaystyle\frac{12}{125} & 0 & 0 \\[16pt] \displaystyle-\frac{12}{125} & \displaystyle\frac{16}{125} + \frac{1}{4} & \displaystyle\frac{12}{125} & \displaystyle-\frac{16}{125} & 0 & \displaystyle-\frac{1}{4} \\[16pt] \displaystyle-\frac{9}{125} & \displaystyle\frac{12}{125} & \displaystyle\frac{9}{125} & \displaystyle-\frac{12}{125} & 0 & 0 \\[16pt] \displaystyle\frac{12}{125} & \displaystyle-\frac{16}{125} & \displaystyle-\frac{12}{125} & \displaystyle\frac{16}{125} & 0 & 0 \\[16pt] 0 & 0 & 0 & 0 & 0 & 0 \\[16pt] 0 & \displaystyle-\frac{1}{4} & 0 & 0 & 0 & \displaystyle\frac{1}{4} \end{array}\right]}_{\equiv\,K} \underbrace{\left[\begin{array}{c} u_{A,x} \\[16pt] u_{A,y} \\[16pt] u_{B,x} \\[16pt] u_{B,y} \\[16pt] u_{C,x} \\[16pt] u_{C,y} \end{array}\right]}_{\equiv\,\vec{u}} \,=\, \underbrace{\left[\begin{array}{c} f_{A,x} \\[16pt] f_{A,y} - P \\[16pt] f_{B,x} \\[16pt] f_{B,y} \\[16pt] f_{C,x} \\[16pt] f_{C,y} \end{array}\right]}_{\equiv\,\vec{f}}

The matrix K \in \mathbb{R}^{6 \times 6} is called the stiffness matrix. It represents the relation between the nodal displacements and internal forces.

From above, we see that the first two columns of K are linearly independent. On the other hand, each of the remaining four columns can be written as a linear combination of these two. Hence, \mbox{rank}(K) = 2.

This means, we need to specify (at least) 4 displacement constraints in order to get a unique solution.

v. Boundary conditions

From the problem description, we know that,

\begin{array}{l} u_{A,x} = 0.05L \\[16pt] u_{B,x},\,u_{B,y},\,u_{C,x},\,u_{C,y} = 0 \\[16pt] f_{A,y} = 0 \end{array}

We can find the y displacement of node A from the second equation:

\begin{aligned} &\,\, \displaystyle-\frac{12}{125}\bigl(0.05L\bigr) \,+\, \Bigl(\frac{16}{125} \,+\, \frac{1}{4}\Bigr)u_{A,y} \,=\, -\frac{PL}{EA} \\[16pt] \Rightarrow &\,\, u_{A,y} \,=\, \left(\frac{\displaystyle-\frac{P}{EA} \,+\, \frac{12}{125}\bigl(0.05\bigr)}{\displaystyle\frac{16}{125} \,+\, \frac{1}{4}}\right)L \end{aligned}

For the model parameters, this simplifies to,

u_{A,y} \,\approx\, -0.0125\,\mbox{m}

As we might expect, node A moves downwards by a small amount.

The strains are given by,

\begin{array}{l} \begin{aligned} \varepsilon_{1} &\,=\, \displaystyle\frac{3}{25L}\bigl(0.05L - 0\bigr) - \frac{4}{25L}\bigl(u_{A,y} - 0\bigr) \\[16pt] &\,\approx\, 0.0080\,\mbox{(in tension)} \end{aligned} \\[40pt] \begin{aligned} \varepsilon_{2} &\,=\, \displaystyle-\frac{1}{4L}\bigl(u_{A,y} - 0\bigr) \\[16pt] &\,\approx\, 0.0031\,\mbox{(in tension)} \end{aligned} \end{array}

The axial (internal) forces are then,

\begin{array}{l} f_{1} \,=\, EA\,\varepsilon_{1} \,\approx\, 8.3995 \times 10^{5}\,\mbox{N} \\[16pt] f_{2} \,=\, EA\,\varepsilon_{2} \,\approx\, 3.2804 \times 10^{5}\,\mbox{N} \end{array}

b. Effect of coordinate system

Consider a truss with 5 bars. All bars are linearly elastic and isotropic, with Young’s modulus E and cross-sectional area A.

A five-bar truss with horizontal load applied to node A.
A five-bar truss with horizontal load applied to node A.

A force of P is applied to node A in the negative x direction. Notice that node D is on a roller that is inclined at an angle of 45^{\circ}.

Find the unknown displacements. Then, determine which bars increase in length, decrease, or remain the same.

i. Element stiffness matrices

Rather than deriving the equilibrium, stress-strain, and strain-displacement equations again, we can compute the element stiffness matrix K_{e} for each bar and assemble the element matrices (stitch them together, so to speak) to create the stiffness matrix K.

For a bar with nodes i and j, its element stiffness matrix is given by,

\displaystyle K^{e} = \frac{EA}{L}\left[\begin{array}{cc|cc} \cos^{2}\theta & \cos\theta\,\sin\theta & -\cos^{2}\theta & -\cos\theta\,\sin\theta \\[16pt] \cos\theta\,\sin\theta & \sin^{2}\theta & -\cos\theta\,\sin\theta & -\sin^{2}\theta \\[16pt]\hline -\cos^{2}\theta & -\cos\theta\,\sin\theta & \cos^{2}\theta & \cos\theta\,\sin\theta \\[16pt] -\cos\theta\,\sin\theta & -\sin^{2}\theta & \cos\theta\,\sin\theta & \sin^{2}\theta \end{array}\right]

The angle \theta measures the rotation that the bar needs to align with the global xy coordinate system.

There is an easy visual trick that I use to correctly determine \theta. Take the positive x axis—in particular, I think of a vector pointing to the right—and rotate it until it matches the vector that points from node i to node j. The amount of rotation needed is \theta. From this, we can compute the cosine and sine values.

For bar element 1, which goes from node A to B, we have \cos\theta = 0 and \sin\theta = 1.

\displaystyle K^{e} = \frac{EA}{L}\left[\begin{array}{cc|cc} 0 & 0 & 0 & 0 \\[16pt] 0 & 1 & 0 & -1 \\[16pt]\hline 0 & 0 & 0 & 0 \\[16pt] 0 & -1 & 0 & 1 \end{array}\right],\,\,\,\,\,\vec{u}^{e} = \left[\begin{array}{c} u_{A,x} \\[16pt] u_{A,y} \\[16pt]\hline u_{B,x} \\[16pt] u_{B,y} \end{array}\right]

For element 2, from node B to C, we have \cos\theta = 1 and \sin\theta = 0.

\displaystyle K^{e} = \frac{EA}{L}\left[\begin{array}{cc|cc} 1 & 0 & -1 & 0 \\[16pt] 0 & 0 & 0 & 0 \\[16pt]\hline -1 & 0 & 1 & 0 \\[16pt] 0 & 0 & 0 & 0 \end{array}\right],\,\,\,\,\,\vec{u}^{e} = \left[\begin{array}{c} u_{B,x} \\[16pt] u_{B,y} \\[16pt]\hline u_{C,x} \\[16pt] u_{C,y} \end{array}\right]

For element 3, from node C to D, \cos\theta = 0 and \sin\theta = -1.

\displaystyle K^{e} = \frac{EA}{L}\left[\begin{array}{cc|cc} 0 & 0 & 0 & 0 \\[16pt] 0 & 1 & 0 & -1 \\[16pt]\hline 0 & 0 & 0 & 0 \\[16pt] 0 & -1 & 0 & 1 \end{array}\right],\,\,\,\,\,\vec{u}^{e} = \left[\begin{array}{c} u_{C,x} \\[16pt] u_{C,y} \\[16pt]\hline u_{D,x} \\[16pt] u_{D,y} \end{array}\right]

For element 4, from node A to D, \cos\theta = 1 and \sin\theta = 0.

\displaystyle K^{e} = \frac{EA}{L}\left[\begin{array}{cc|cc} 1 & 0 & -1 & 0 \\[16pt] 0 & 0 & 0 & 0 \\[16pt]\hline -1 & 0 & 1 & 0 \\[16pt] 0 & 0 & 0 & 0 \end{array}\right],\,\,\,\,\,\vec{u}^{e} = \left[\begin{array}{c} u_{A,x} \\[16pt] u_{A,y} \\[16pt]\hline u_{D,x} \\[16pt] u_{D,y} \end{array}\right]

Finally, for element 5, from node A to C, \displaystyle\cos\theta = \frac{1}{\sqrt{2}} and \displaystyle\sin\theta = \frac{1}{\sqrt{2}}.

\displaystyle K^{e} = \frac{EA}{\sqrt{2}L}\left[\begin{array}{cc|cc} \displaystyle\frac{1}{2} & \displaystyle\frac{1}{2} & \displaystyle-\frac{1}{2} & \displaystyle-\frac{1}{2} \\[16pt] \displaystyle\frac{1}{2} & \displaystyle\frac{1}{2} & \displaystyle-\frac{1}{2} & \displaystyle-\frac{1}{2} \\[16pt]\hline \displaystyle-\frac{1}{2} & \displaystyle-\frac{1}{2} & \displaystyle\frac{1}{2} & \displaystyle\frac{1}{2} \\[16pt] \displaystyle-\frac{1}{2} & \displaystyle-\frac{1}{2} & \displaystyle\frac{1}{2} & \displaystyle\frac{1}{2} \end{array}\right],\,\,\,\,\,\vec{u}^{e} = \left[\begin{array}{c} u_{A,x} \\[16pt] u_{A,y} \\[16pt]\hline u_{C,x} \\[16pt] u_{C,y} \end{array}\right]

ii. Stiffness equation

Assemble the element matrices into a global matrix K \in \mathbb{R}^{8 \times 8}:

\displaystyle K = \frac{EA}{L}\left[\begin{array}{cc|cc|cc|cc} \displaystyle 1 + \frac{1}{2\sqrt{2}} & \displaystyle\frac{1}{2\sqrt{2}} & & & \displaystyle-\frac{1}{2\sqrt{2}} & \displaystyle-\frac{1}{2\sqrt{2}} & -1 & \\[16pt] \displaystyle\frac{1}{2\sqrt{2}} & \displaystyle 1 + \frac{1}{2\sqrt{2}} & & -1 & \displaystyle-\frac{1}{2\sqrt{2}} & \displaystyle-\frac{1}{2\sqrt{2}} & & \\[16pt]\hline & & 1 & & -1 & & & \\[16pt] & -1 & & 1 & & & & \\[16pt]\hline \displaystyle-\frac{1}{2\sqrt{2}} & \displaystyle-\frac{1}{2\sqrt{2}} & & & \displaystyle 1 + \frac{1}{2\sqrt{2}} & \displaystyle\frac{1}{2\sqrt{2}} & & \\[16pt] \displaystyle-\frac{1}{2\sqrt{2}} & \displaystyle-\frac{1}{2\sqrt{2}} & & & \displaystyle\frac{1}{2\sqrt{2}} & \displaystyle 1 + \frac{1}{2\sqrt{2}} & & -1 \\[16pt]\hline -1 & & & & & & 1 & \\[16pt] & & & & & -1 & & 1 \end{array}\right]

Note that \mbox{rank}(K) = 5, so we must have at least 3 displacement constraints in order to get a unique solution.

iii. Boundary conditions

From the drawing in the beginning, we can easily see that,

\boxed{\begin{array}{l} u_{C,x},\,u_{C,y} = 0 \\[16pt] f_{A,x},\,f_{A,y},\,f_{B,x},\,f_{B,y} = 0 \end{array}}

There is a duality between known displacements and unknown forces (and vice versa). We know the displacements of node C, but don’t know its reaction forces. We know the reaction forces of nodes A and B, but don’t know their displacements.

We have 2 displacement constraints and need 1 more to find a unique solution. What about node D, which moves along the inclined plane?

Consider a local coordinate system at node D. This means, one of the axes is parallel to the inclined plane. Since the node’s expected movement matches one of the axes, it will be easier to describe the displacement of node D.

Local coordinate system at node D
A local coordinate system at node D.

Without loss of generality, we can assume that the local coordinate system (\tilde{x},\,\tilde{y}) shares the same origin as the global coordinate system (x,\,y).

The local coordinate system shares the same origin as the global.
The local coordinate system shares the same origin as the global.

Our goal is to describe the displacement of node D in the local coordinate system. Keep in mind that the act of displacement is the same whether we consider global or local coordinates, but the numbers representing the displacement depend on which coordinate system we consider. More precisely, they depend on which basis we consider to represent the space of \mathbb{R}^{2}.

The basis for the local coordinate system, in terms of the global coordinates (the standard basis), is given by,

\tilde{e}_{1} = \left[\begin{array}{c} \cos(45^{\circ}) \\[16pt] \sin(45^{\circ}) \end{array}\right],\,\,\,\tilde{e}_{2} = \left[\begin{array}{c} -\sin(45^{\circ}) \\[16pt] \cos(45^{\circ}) \end{array}\right]

Since the act of displacement is the same,

\begin{aligned} \left[\begin{array}{c} u_{D,x} \\[16pt] u_{D,y} \end{array}\right] & =\, \tilde{u}_{D,x}\,\tilde{e}_{1} \,+\, \tilde{u}_{D,y}\,\tilde{e}_{2} \\[16pt] & =\, \tilde{u}_{D,x}\left[\begin{array}{c} \cos(45^{\circ}) \\[16pt] \sin(45^{\circ}) \end{array}\right] \,+\, \tilde{u}_{D,y}\left[\begin{array}{c} -\sin(45^{\circ}) \\[16pt] \cos(45^{\circ}) \end{array}\right] \\[16pt] & =\, \left[\begin{array}{cc} \cos(45^{\circ}) & -\sin(45^{\circ}) \\[16pt] \sin(45^{\circ}) & \cos(45^{\circ}) \end{array}\right] \left[\begin{array}{c} \tilde{u}_{D,x} \\[16pt] \tilde{u}_{D,y} \end{array}\right] \end{aligned}

We see that the global and local displacements are linearly related! (More information on this in the next blog post.)

Since \tilde{u}_{D,x} is free and \tilde{u}_{D,y} = 0, we have the following boundary conditions (BCs):

\begin{array}{l} \displaystyle u_{D,x} = \frac{1}{\sqrt{2}}\,\tilde{u}_{D,x} \\[16pt] \displaystyle u_{D,y} = \frac{1}{\sqrt{2}}\,\tilde{u}_{D,x} \end{array}

We see that the x and y displacements must equal (obvious in hindsight).

\boxed{\displaystyle u_{D,y} = u_{D,x}}

This equation gives us a starting point for the 3rd displacement constraint needed to solve the stiffness equation.

Similarly, the force at node D acts the same regardless of the coordinate system. The global and local numbers are related linearly:

\left[\begin{array}{c} f_{D,x} \\[16pt] f_{D,y} \end{array}\right] = \left[\begin{array}{cc} \cos(45^{\circ}) & -\sin(45^{\circ}) \\[16pt] \sin(45^{\circ}) & \cos(45^{\circ}) \end{array}\right] \left[\begin{array}{c} \tilde{f}_{D,x} \\[16pt] \tilde{f}_{D,y} \end{array}\right]

Since \tilde{f}_{D,x} = 0 and \tilde{f}_{D,y} is free, we must have,

\begin{array}{l} \displaystyle f_{D,x} = -\frac{1}{\sqrt{2}}\,\tilde{f}_{D,y} \\[16pt] \displaystyle f_{D,y} = \frac{1}{\sqrt{2}}\,\tilde{f}_{D,y} \end{array}

Hence,

\boxed{\displaystyle f_{D,y} = -f_{D,x}}

Recall the 8th line in the stiffness equation:

\displaystyle\frac{EA}{L}\Bigl[(-1)\,u_{C,y} + (1)\,u_{D,y}\Bigr] = f_{D,y}

Since u_{C,y} = 0, u_{D,y} = u_{D,x}, and f_{D,y} = -f_{D,x}, we end up with the following constraint for node D’s displacement in x direction:

\boxed{\displaystyle f_{D,x} = -\frac{EA}{L}\,u_{D,x}}

We now have 3 displacement constraints needed to solve the stiffness equation!

iv. Reduced system

When we apply the 3 constraints to our 8 \times 8 equation, the system reduces to a 5 \times 5 equation:

\displaystyle\frac{EA}{L}\left[\begin{array}{cc|cc|c} \displaystyle 1 + \frac{1}{2\sqrt{2}} & \displaystyle\frac{1}{2\sqrt{2}} & & & -1 \\[16pt] \displaystyle\frac{1}{2\sqrt{2}} & \displaystyle 1 + \frac{1}{2\sqrt{2}} & & -1 & \\[16pt]\hline & & 1 & & \\[16pt] & -1 & & 1 & \\[16pt]\hline -1 & & & & 2 \end{array}\right] \left[\begin{array}{c} u_{A,x} \\[16pt] u_{A,y} \\[16pt] u_{B,x} \\[16pt] u_{B,y} \\[16pt] u_{D,x} \end{array}\right] = \left[\begin{array}{c} -P \\[16pt] 0 \\[16pt] 0 \\[16pt] 0 \\[16pt] 0 \end{array}\right]

Note that the 5th entry of \vec{f}_{\tiny\mbox{reduced}}, which had been f_{D,x}, was changed to -(EA/L)\,u_{D,x}. Since the variable u_{D,x} is unknown, the coefficient -(EA/L) moved to the left-hand side of the equation. This is how we get the coefficient of 2 on row 5, column 5 of the matrix K_{\tiny\mbox{reduced}}.

Solve for the vector \vec{u}_{\tiny\mbox{reduced}}.

\displaystyle\left[\begin{array}{c} u_{A,x} \\[16pt] u_{A,y} \\[16pt] u_{B,x} \\[16pt] u_{B,y} \\[16pt] u_{D,x} \end{array}\right] = \frac{PL}{EA}\left[\begin{array}{c} -2 \\[16pt] 2 \\[16pt] 0 \\[16pt] 2 \\[16pt] -1 \end{array}\right]

Since we know the displacements, we can visualize the deformed structure:

Our structure before and after the deformation.
Before and after the deformation.

v. Postprocessing

Once we have the displacements, we can calculate the strains (if desired, stresses too) in each bar.

\begin{aligned} &\displaystyle\varepsilon_{1} = \frac{1}{L}\,(u_{B,y} - u_{A,y}) = \frac{P}{EA}\,(2 - 2) = 0 \\[16pt] &\displaystyle\varepsilon_{2} = \frac{1}{L}\,(u_{C,x} - u_{B,x}) = \frac{P}{EA}\,(0 - 0) = 0 \\[16pt] &\displaystyle\varepsilon_{3} = -\frac{1}{L}\,(u_{D,y} - u_{C,y}) = -\frac{P}{EA}\,(-1 - 0) = \frac{P}{EA}\,\,\mbox{(in tension)} \\[16pt] &\displaystyle\varepsilon_{4} = \frac{1}{L}\,(u_{D,x} - u_{A,x}) = \frac{P}{EA}\,(-1 + 2) = \frac{P}{EA}\,\,\mbox{(in tension)} \\[16pt] &\displaystyle\varepsilon_{5} = \frac{1}{2L}\,(u_{C,x} - u_{A,x}) + \frac{1}{2L}\,(u_{C,y} - u_{A,y}) = \frac{P}{2EA}\,(0 + 2) + \frac{P}{2EA}(0 - 2) = 0 \end{aligned}

From these strains, we can feel the deformed geometry. Since bar 5 (the diagonal one) did not change in length while bars 3 and 4 increased in length, the angle between bars 3 and 4 must have decreased. Bars 1 and 2 are simply lifted up.

vi. Remarks

There is an easier way to solve this problem. We recall that the difficulty of solving it came from determining the boundary condition at node D.

What if, instead, we consider the following coordinate system (tilt our head by 45^{\circ})?

The same problem in a different coordinate system
The same problem in a different coordinate system.

The roller BC becomes much easier to derive now.

\boxed{\begin{array}{l} \tilde{u}_{D,y} = 0 \\[16pt] \tilde{f}_{D,x} = 0 \end{array}}

c. Effect of symmetry

Last but not least, let’s look at how to simplify a problem when there is symmetry in the structure and boundary conditions.

A 8-bar truss that exhibits symmetry
A 8-bar truss that exhibits symmetry.

As usual, all bars are linearly elastic and isotropic. To make our hand calculations easy, we assign the material parameters E and \sqrt{2}A to bars 1, 2, 7, and 8, and parameters of E and A to bars 3, 4, 5, and 6.

With this setup, the truss exhibits symmetry about the line \overline{\mbox{BC}}, i.e. identical geometry, material, loading, and displacements occur on each side about this line.

When there is symmetry, we can reduce the problem by following these steps:

  • For loads that occur in the plane of symmetry, apply half of the total load to the reduced structure.
  • For elements in the plane of symmetry, take half of the cross-sectional area.
  • For nodes in the plane of symmetry, the displacement components that are normal (perpendicular) to the plane of symmetry must be set to 0.

Now, we can solve a simpler problem:

Use symmetry to reduce the number of unknowns (degrees of freedom)
Use symmetry to reduce the number of unknowns (degrees of freedom).

Note that bars 4 and 5 now have the (effective) cross-sectional area 0.5A.

i. Element stiffness matrices

For element 1 (node A to B), \displaystyle\cos\theta = \frac{1}{\sqrt{2}} and \displaystyle\sin\theta = \frac{1}{\sqrt{2}}.

\displaystyle K^{e} = \frac{E(\sqrt{2}A)}{\sqrt{2}L}\left[\begin{array}{cc|cc} 0.5 & 0.5 & -0.5 & -0.5 \\[16pt] 0.5 & 0.5 & -0.5 & -0.5 \\[16pt]\hline -0.5 & -0.5 & 0.5 & 0.5 \\[16pt] -0.5 & -0.5 & 0.5 & 0.5 \end{array}\right],\,\,\,\,\,\vec{u}^{e} = \left[\begin{array}{c} u_{A,x} \\[16pt] u_{A,y} \\[16pt]\hline u_{B,x} \\[16pt] u_{B,y} \end{array}\right]

For element 2 (node A to C), \displaystyle\cos\theta = \frac{1}{\sqrt{2}} and \displaystyle\sin\theta = -\frac{1}{\sqrt{2}}.

\displaystyle K^{e} = \frac{E(\sqrt{2}A)}{\sqrt{2}L}\left[\begin{array}{cc|cc} 0.5 & -0.5 & -0.5 & 0.5 \\[16pt] -0.5 & 0.5 & 0.5 & -0.5 \\[16pt]\hline -0.5 & 0.5 & 0.5 & -0.5 \\[16pt] 0.5 & -0.5 & -0.5 & 0.5 \end{array}\right],\,\,\,\,\,\vec{u}^{e} = \left[\begin{array}{c} u_{A,x} \\[16pt] u_{A,y} \\[16pt]\hline u_{C,x} \\[16pt] u_{C,y} \end{array}\right]

For element 3 (node A to D), \cos\theta = 1 and \sin\theta = 0.

\displaystyle K^{e} = \frac{EA}{L}\left[\begin{array}{cc|cc} 1 & 0 & -1 & 0 \\[16pt] 0 & 0 & 0 & 0 \\[16pt]\hline -1 & 0 & 1 & 0 \\[16pt] 0 & 0 & 0 & 0 \end{array}\right],\,\,\,\,\,\vec{u}^{e} = \left[\begin{array}{c} u_{A,x} \\[16pt] u_{A,y} \\[16pt]\hline u_{D,x} \\[16pt] u_{D,y} \end{array}\right]

For element 4 (node B to D), \cos\theta = 0 and \sin\theta = -1.

\displaystyle K^{e} = \frac{E(0.5A)}{L}\left[\begin{array}{cc|cc} 0 & 0 & 0 & 0 \\[16pt] 0 & 1 & 0 & -1 \\[16pt]\hline 0 & 0 & 0 & 0 \\[16pt] 0 & -1 & 0 & 1 \end{array}\right],\,\,\,\,\,\vec{u}^{e} = \left[\begin{array}{c} u_{B,x} \\[16pt] u_{B,y} \\[16pt]\hline u_{D,x} \\[16pt] u_{D,y} \end{array}\right]

For element 5 (node C to D), \cos\theta = 0 and \sin\theta = 1.

\displaystyle K^{e} = \frac{E(0.5A)}{L}\left[\begin{array}{cc|cc} 0 & 0 & 0 & 0 \\[16pt] 0 & 1 & 0 & -1 \\[16pt]\hline 0 & 0 & 0 & 0 \\[16pt] 0 & -1 & 0 & 1 \end{array}\right],\,\,\,\,\,\vec{u}^{e} = \left[\begin{array}{c} u_{C,x} \\[16pt] u_{C,y} \\[16pt]\hline u_{D,x} \\[16pt] u_{D,y} \end{array}\right]

ii. Stiffness equation

Assemble the element matrices to get K \in \mathbb{R}^{8 \times 8}:

\displaystyle K = \frac{EA}{L}\left[\begin{array}{cc|cc|cc|cc} 2 & & -0.5 & -0.5 & -0.5 & 0.5 & -1 & \\[16pt] & 1 & -0.5 & -0.5 & 0.5 & -0.5 & & \\[16pt]\hline -0.5 & -0.5 & 0.5 & 0.5 & & & & \\[16pt] -0.5 & -0.5 & 0.5 & 1 & & & & -0.5 \\[16pt]\hline -0.5 & 0.5 & & & 0.5 & -0.5 & & \\[16pt] 0.5 & -0.5 & & & -0.5 & 1 & & -0.5 \\[16pt]\hline -1 & & & & & & 1 & \\[16pt] & & & -0.5 & & -0.5 & & 1 \end{array}\right]

iii. Boundary conditions

Next, apply the boundary conditions:

\begin{array}{l} u_{A,x},\,u_{A,y},\,u_{B,x},\,u_{C,x},\,u_{D,x} = 0 \\[16pt] f_{B,y},\,f_{C,y},\,f_{D,y} = 0 \end{array}

Note that, had we considered the full problem instead, we would have 6 unknown displacements instead of 3.

iv. Reduced system

Solve the equation,

\displaystyle\frac{EA}{L}\left[\begin{array}{ccc} 1 & 0 & -0.5 \\[16pt] 0 & 1 & -0.5 \\[16pt] -0.5 & -0.5 & 1 \end{array}\right] \left[\begin{array}{c} u_{B,y} \\[16pt] u_{C,y} \\[16pt] u_{D,y} \end{array}\right] = \left[\begin{array}{c} 0 \\[16pt] 0 \\[16pt] -0.5P \end{array}\right]

to get the unknown displacements,

\displaystyle\left[\begin{array}{c} u_{B,y} \\[16pt] u_{C,y} \\[16pt] u_{D,y} \end{array}\right] = \frac{PL}{EA}\left[\begin{array}{c} -0.5 \\[16pt] -0.5 \\[16pt] -1 \end{array}\right]

Hence, the structure deforms in the following manner:

Our structure before and after the deformation.
Before and after the deformation.

v. Remarks

As we saw from above, symmetry helps us simplify a problem. For structural problems in 2D or 3D, reducing the number of unknowns (DOFs) can have a significant positive impact on speed and accuracy.

A classic example is the infinite plate with a circular hole that goes under uniaxial or biaxial tension. For this problem, we just need to model a quarter of the geometry.

Stress distribution for a finite plate with a circular hole
Stress distribution for a finite plate with a circular hole.