MV-12: The Grand Finale – The Finite Element Method (FEM)

Summary
The conclusion of our Mechanical Vibrations series. We bridge the gap between theory and the real world using FEM. Topics: Discretization, Shape Functions, Matrix Assembly, and a modern Julia implementation.

From a Single Spring to the Whole World.

We started our journey with a single mass on a spring. We evolved to two masses, then to Multi-DOF systems, and finally to continuous systems like strings and beams. We learned that while analytical methods are elegant, they hit a wall: continuous systems have infinite degrees of freedom and require complex partial differential equations.

But what if you need to find the natural frequency of an engine block, a Tesla chassis, or a jet turbine blade? These aren’t simple beams. There are no textbook formulas for them.

Enter the Finite Element Method (FEM).

FEM is the bridge between theory and reality. It allows us to shatter complex continuous structures into finite simple “elements,” turning an impossible calculus nightmare into a solvable matrix algebra problem. In this Grand Finale, based on Chapter 12 of Rao’s Mechanical Vibrations, we will explore the industry standard for structural analysis.

From Infinite to Finite: The Core Concept

Analytical solutions (like the wave equation) treat a structure as a continuous body with infinite degrees of freedom. This is mathematically elegant but practically impossible for anything other than a simple rod or beam.

FEM takes a pragmatic approach: “Divide and Conquer.”

Discretization: From Continuous to Discrete
Left: A continuous body with infinite complexity. Right: The same body 'meshed' into finite triangular elements. We solve for the displacement at the blue nodes.

The method involves three intuitive steps:

  1. Discretization (Meshing): Cut the complex structure into small, simple shapes called Elements (lines, triangles, cubes).
  2. Interpolation: Since we can’t calculate displacement everywhere, we only calculate it at specific points called Nodes. Inside each element, we assume the displacement follows a simple rule (like a straight line).
  3. Assembly: We connect these simple elements together to rebuild the complex structure mathematically, forming a massive system of linear equations ($[K]{u} = {F}$).

The Element Equation: Stiffness and Mass

Let’s look at the simplest Lego block in our set: the Bar Element (a rod subjected to axial loads). It has two nodes (1 at the start, 2 at the end) and moves only along its axis.

Bar Element: Stiffness and Mass
Top: A bar element interpolates displacement linearly between nodes $u_1$ and $u_2$. Bottom: Mass can be modeled as 'Lumps' at nodes or as a 'Consistent' distributed density.

The Stiffness Matrix: $[k]$ Derived

Where does this matrix come from? It’s not magic; it’s just Hooke’s Law disguised in algebra.

  1. Strain ($\epsilon$): If the displacement is linear, the strain is constant: $\epsilon = \frac{u_2 - u_1}{L}$.
  2. Stress ($\sigma$): From Hooke’s Law, stress is $\sigma = E\epsilon = E \frac{u_2 - u_1}{L}$.
  3. Force ($F$): Force is Stress $\times$ Area ($F = \sigma A$).
    • Force at Node 2 ($f_2$) = $+F = \frac{EA}{L}(u_2 - u_1)$
    • Force at Node 1 ($f_1$) = $-F = \frac{EA}{L}(u_1 - u_2)$

Writing this as a matrix system $\{f\} = [k]\{u\}$: $$ \begin{Bmatrix} f_1 \ f_2 \end{Bmatrix} = \frac{EA}{L} \begin{bmatrix} 1 & -1 \ -1 & 1 \end{bmatrix} \begin{Bmatrix} u_1 \ u_2 \end{Bmatrix} $$ This matrix relates the nodal forces to nodal displacements.

The Mass Matrix: $[m]$

How do we handle the mass? Since the element is continuous, mass is everywhere. We have two choices:

  1. Lumped Mass (The Lazy Way): We simply chop the mass in half and stick $m/2$ at Node 1 and $m/2$ at Node 2. This creates a diagonal matrix (zeros off-diagonal). It’s computationally cheap but less accurate.

  2. Consistent Mass (The Accurate Way): We respect the physics. We use the same linear shape function to calculate kinetic energy. This results in a full matrix, showing that the motion of Node 1 is dynamically coupled to Node 2. $$ [m] = \frac{\rho Al}{6} \begin{bmatrix} 2 & 1 \ 1 & 2 \end{bmatrix} $$

The Beam Element: Bending and Rotation

Things get more interesting when we look at a Beam, which bends instead of just stretching.

Beam Element Degrees of Freedom
Top: A beam element has 4 Degrees of Freedom: Vertical displacement ($v$) and Rotation ($\theta$) at each node. Bottom: Why we need Cubic functions—Linear interpolation creates 'kinks', while Cubic ensures smooth bending.

Why 4 Degrees of Freedom?

In a bar, we only cared about where the node is ($u$). In a beam, we also care about which way it points ($\theta$).

  • Total DOFs = 2 Nodes $\times$ 2 DOFs/Node = 4.
  • State Vector: ${q} = [v_1, \theta_1, v_2, \theta_2]^T$

Cubic Interpolation (Hermite Polynomials)

A simple straight line (Linear Interpolation) won’t work for beams. As visualised above, connecting straight lines creates “kinks” at the nodes. A kink means the slope changes instantly, which implies infinite curvature and infinite bending moment—physically impossible!

To insure a smooth slope (continuity), we use a Cubic Polynomial. This leads to a larger $4 \times 4$ stiffness matrix:

$$ [k] = \frac{EI}{l^3} \begin{bmatrix} 12 & 6l & -12 & 6l \\ 6l & 4l^2 & -6l & 2l^2 \\ -12 & -6l & 12 & -6l \\ 6l & 2l^2 & -6l & 4l^2 \end{bmatrix} $$

Transformation and Assembly: The Big Picture

Real structures aren’t just horizontal beams. They are 3D skeletons. A strut might be angled at 45°. To handle this, we perform coordinate transformations and then assemble everything into a “Super Matrix”.

FEM Assembly Process
Left: Transforming an element from its local axis ($x'$) to the global axis ($X$). Right: Assembling the Global Matrix. Notice how Element 1 (Blue) and Element 2 (Red) overlap at Node 2, summing their stiffnesses.
  1. Transformation ($ \vec{q} = [\lambda] \vec{Q} $): Each element calculates its stiffness in its own private “local” world. We then multiply by a rotation matrix $[\lambda]$ to translate this into the “global” language ($X, Y$) shared by the whole structure. $$ [K]_{global} = [\lambda]^T [k]_{local} [\lambda] $$

  2. Assembly (The “Lego” Step): This is the magic logic of FEM, but it comes from a simple physical principle: Node Equilibrium.

    • Element 1 pulls on Node 2 with stiffness $k_{22}^{(1)}$.
    • Element 2 pulls on Node 2 with stiffness $k_{22}^{(2)}$.
    • Node 2 is just a point; it must satisfy Newton’s Third Law. The total resistance force at Node 2 is simply the sum of the resistance from both elements.

    Therefore, at the shared index (2,2), the global stiffness is the sum of contributions: $$ K_{22} = k_{22}^{(1)} + k_{22}^{(2)} $$

The result is a massive system equation that looks exactly like our familiar MDOF equation, just bigger: $$ [M]\ddot{\vec{Q}} + [K]\vec{Q} = \vec{F} $$

Boundary Conditions: Stopping the Drift

When we first assemble the global matrices, the structure is technically “floating” in space. If you push it, it creates infinite acceleration (Singular Matrix). To solve this, we must apply Boundary Conditions (fixing it to the ground).

Applying Boundary Conditions
Left: A floating structure has rigid body motion (Singular). Fixing one node makes it solvable. Right: Mathematically, setting $u_1=0$ allows us to 'strike out' the corresponding rows and columns, reducing the matrix size.

The Strategy: If Node 1 is fixed to a wall ($u_1 = 0$), we don’t need to calculate its motion—we already know it!

  1. Delete the Row & Column corresponding to the constrained Degree of Freedom ($u_1$).
  2. Reduce the Matrix from size $N \times N$ to $(N-1) \times (N-1)$.
  3. The remaining matrix is now non-singular and ready to be solved.

Julia Example: Stepped Bar

Let’s put this into practice. We’ll find the natural frequencies of a Stepped Bar fixed at one end ($u_1=0$).

  • Element 1: Area $2A$, Length $L$.
  • Element 2: Area $A$, Length $L$.

We will use Julia for this calculation, leveraging its native LinearAlgebra library to solve the generalized eigenvalue problem.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
# Finite Element Analysis of a Stepped Bar in Julia
using LinearAlgebra

# 1. Define Parameters (Steel)
E = 200e9; rho = 7800; L_elem = 1.0; A0 = 1e-4

# 2. Define Element Matrices (Function)
function bar_element(E, A, L, rho)
    k = (E*A/L) * [1 -1; -1 1]
    m = (rho*A*L/6) * [2 1; 1 2] # Consistent Mass
    return k, m
end

# 3. Create Elements
k1, m1 = bar_element(E, 2*A0, L_elem, rho) # Element 1 (Larger)
k2, m2 = bar_element(E, A0, L_elem, rho)   # Element 2 (Smaller)

# 4. Global Assembly (3 Nodes)
K = zeros(3,3); M = zeros(3,3)
# Add Element 1 (Nodes 1-2)
K[1:2, 1:2] += k1; M[1:2, 1:2] += m1
# Add Element 2 (Nodes 2-3)
K[2:3, 2:3] += k2; M[2:3, 2:3] += m2

# 5. Apply Boundary Condition (u1 = 0)
# Remove Row 1 and Col 1
K_reduced = K[2:3, 2:3]
M_reduced = M[2:3, 2:3]

# 6. Solve Eigenvalue Problem (K*v = λ*M*v)
evals, evecs = eigen(K_reduced, M_reduced)
omega = sqrt.(evals)

println("Natural Frequencies: $(round.(omega, digits=2)) rad/s")
Stepped Bar Mode Shapes
The resulting Mode Shapes calculated by the FEM code. Mode 1 is the fundamental axial vibration, while Mode 2 shows a node (zero displacement) in the middle.

The code automatically handles the complex assembly and reduction, spitting out the exact frequencies we need. For a real car chassis, the logic is identical—just with a million elements instead of two.

This code illustrates the essence of FEM: generate element matrices, assemble them based on connectivity, apply constraints, and solve the eigenvalue problem.

Series Conclusion: The Circle of Vibration

We have come a long way. From a simple mass on a spring to the matrix assembly of complex structures, the journey of Mechanical Vibrations is compliant with one truth: Complexity is just a superposition of simple things.

Here is the complete roadmap of our 12-chapter journey:

graph TD
    subgraph "Phase 1: The Physics of One"
        MV01[MV-01: Modeling] --> MV02[MV-02: Free Vib]
        MV02 --> MV03[MV-03: Harmonic]
        MV03 --> MV04[MV-04: General]
    end

    subgraph "Phase 2: The Dance of Many"
        MV04 -.-> MV05
        MV05[MV-05: 2-DOF] --> MV06[MV-06: MDOF Matrix]
        MV06 --> MV07[MV-07: Eigenvalues]
    end

    subgraph "Phase 3: The Infinite Continuum"
        MV06 -.-> MV08
        MV08[MV-08: Continuous Strings/Beams] --> MV12[MV-12: Finite Element Method]
        MV12 --> MV_REAL[Real World Structures]
    end

    subgraph "The Toolbox"
        MV09[MV-09: Control]
        MV10[MV-10: Measurement]
        MV11[MV-11: Numerical Integration]
    end
    
    MV07 -.-> MV11
    MV11 -.-> MV12
    MV_REAL -.-> MV09
    MV_REAL -.-> MV10

    style MV01 fill:#3498db,color:#fff,stroke:#2980b9
    style MV12 fill:#e74c3c,color:#fff,stroke:#c0392b,stroke-width:4px
    style MV_REAL fill:#2ecc71,color:#fff,stroke:#27ae60

The Four Pillars

No matter if you are simulating a MEMS gyroscope or an earthquake-proof skyscraper, the physics always boils down to four immutable pillars:

ConceptThe Question it AsksThe Mathematical Role
Inertia ($m$)“How hard is it to move?”The Kinetic Energy Storage $\frac{1}{2}m\dot{x}^2$
Stiffness ($k$)“How hard is it to deform?”The Potential Energy Storage $\frac{1}{2}kx^2$
Damping ($c$)“Where does the energy go?”The Energy Dissipator (Heat)
Frequency ($\omega$)“How fast does it want to shake?”The System’s DNA $\sqrt{k/m}$

Final Thought: The Finite Element Method (MV-12) is not just a calculation tool; it is the bridge that connects our simplified textbook models to the chaotic reality of the physical world. By turning calculus into algebra, it empowers us to predict the future of any structure.

Thank you for joining me on this exploration of Mechanical Vibrations. Now, go simulate something beautiful.

References

Rao, S. S. (2018). Mechanical Vibrations (6th ed.). Pearson.