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.
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 and the other . The bars are joined by a roller at node A, while their other ends are fixed, as shown below:
A force of is applied to node A in the negative direction, while node A moves by an amount of in the positive direction. Both bars are assumed to be linearly elastic and isotropic, with Young’s modulus and cross-sectional area .
Determine the displacement of node A and the axial force in each bar. For a model problem, we can use the following parameters:
First, let’s agree on our notations. We will write the displacements of the nodes (also called degrees of freedom) as,
the bars’ internal forces as,
and the reaction forces as,
ii. Equilibrium equations
The net forces along and directions must be 0. To help us derive the equilibrium equations, we will draw free-body diagrams for each node.
For node A, we get the following system of equations:
For node B,
Finally, for node C,
We can combine these equations into a single matrix equation:
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:
Note, and denote the axial strains.
iii. Strain-displacement equations
For bar 1,
For bar 2,
iv. Stiffness equation
We combine the equilibrium, stress-strain, and strain-displacement equations to arrive at a single matrix equation:
The matrix 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 are linearly independent. On the other hand, each of the remaining four columns can be written as a linear combination of these two. Hence, .
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,
We can find the displacement of node A from the second equation:
For the model parameters, this simplifies to,
As we might expect, node A moves downwards by a small amount.
The strains are given by,
The axial (internal) forces are then,
b. Effect of coordinate system
Consider a truss with 5 bars. All bars are linearly elastic and isotropic, with Young’s modulus and cross-sectional area .
A force of is applied to node A in the negative direction. Notice that node D is on a roller that is inclined at an angle of .
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 for each bar and assemble the element matrices (stitch them together, so to speak) to create the stiffness matrix .
For a bar with nodes and , its element stiffness matrix is given by,
The angle measures the rotation that the bar needs to align with the global – coordinate system.
There is an easy visual trick that I use to correctly determine . Take the positive axis—in particular, I think of a vector pointing to the right—and rotate it until it matches the vector that points from node to node . The amount of rotation needed is . From this, we can compute the cosine and sine values.
For bar element 1, which goes from node A to B, we have and .
For element 2, from node B to C, we have and .
For element 3, from node C to D, and .
For element 4, from node A to D, and .
Finally, for element 5, from node A to C, and .
ii. Stiffness equation
Assemble the element matrices into a global matrix :
Note that , 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,
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.
Without loss of generality, we can assume that the local coordinate system shares the same origin as the global coordinate system .
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 .
The basis for the local coordinate system, in terms of the global coordinates (the standard basis), is given by,
Since the act of displacement is the same,
We see that the global and local displacements are linearly related! (More information on this in the next blog post.)
Since is free and , we have the following boundary conditions (BCs):
We see that the and displacements must equal (obvious in hindsight).
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:
Since and is free, we must have,
Recall the 8th line in the stiffness equation:
Since , , and , we end up with the following constraint for node D’s displacement in direction:
We now have 3 displacement constraints needed to solve the stiffness equation!
iv. Reduced system
When we apply the 3 constraints to our equation, the system reduces to a equation:
Note that the 5th entry of , which had been , was changed to . Since the variable is unknown, the coefficient 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 .
Solve for the vector .
Since we know the displacements, we can visualize the deformed structure:
Once we have the displacements, we can calculate the strains (if desired, stresses too) in each bar.
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.
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 )?
The roller BC becomes much easier to derive now.
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.
As usual, all bars are linearly elastic and isotropic. To make our hand calculations easy, we assign the material parameters and to bars 1, 2, 7, and 8, and parameters of and to bars 3, 4, 5, and 6.
With this setup, the truss exhibits symmetry about the line , 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:
Note that bars 4 and 5 now have the (effective) cross-sectional area .
i. Element stiffness matrices
For element 1 (node A to B), and .
For element 2 (node A to C), and .
For element 3 (node A to D), and .
For element 4 (node B to D), and .
For element 5 (node C to D), and .
ii. Stiffness equation
Assemble the element matrices to get :
iii. Boundary conditions
Next, apply the boundary conditions:
Note that, had we considered the full problem instead, we would have 6 unknown displacements instead of 3.
iv. Reduced system
Solve the equation,
to get the unknown displacements,
Hence, the structure deforms in the following manner:
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.