Visualizing Bending Moment Diagrams with Matplotlib for Structural Engineers
title: "Visualizing Bending Moment Diagrams with Matplotlib for Structural Engineers"
published: true
description: "A practical Python tutorial to generate accurate bending moment and shear force diagrams for simple and continuous beams — no FEA software needed."
tags: python, engineering, matplotlib, tutorial
If you've ever worked on a structural project, you know that bending moment diagrams (BMDs) are the backbone of beam design. Every reinforced concrete beam, every steel section check starts here.
The problem? Most structural engineers either rely on black-box software that gives them the diagram without the intuition, or they draw it by hand — which works but doesn't scale.
In this tutorial, we'll build a fully functional BMD and SFD (Shear Force Diagram) generator in pure Python using NumPy and Matplotlib. No FEA software, no license fee, no black box.
By the end, you'll be able to generate diagrams for:
- Simply supported beams with point loads
- Beams with uniformly distributed loads (UDL)
- Continuous beams with multiple spans
💡 This is inspired by the calculation workflows we use at Structalis, a structural engineering consultancy based in Paris, where we deal with reinforced concrete, steel, and timber structures daily under Eurocode standards.
1. Prerequisites
pip install numpy matplotlib
You should be comfortable with basic Python. Structural mechanics knowledge helps but isn't strictly required — we'll cover the theory as we go.
2. Theory in 60 Seconds
For a beam in static equilibrium, at any cross-section at position x:
- Shear Force V(x) = algebraic sum of all vertical forces to the left of x
- Bending Moment M(x) = algebraic sum of moments of all forces to the left of x
The key relationships:
dV/dx = -q(x) (distributed load intensity)
dM/dx = V(x) (shear is the derivative of moment)
Sign convention (Eurocode / French standard):
- Positive moment → sagging (tension at the bottom)
- Positive shear → left side up, right side down
3. Simply Supported Beam — Point Load
Let's start with the classic case: a simply supported beam of span L, with a single point load P at position a from the left support.
P
↓
A───────────────B
|←──a──→|←─b──→|
←──────L────────→
Ra Rb
Reactions
Ra = P × b / L
Rb = P × a / L
Shear Force and Bending Moment expressions
For 0 ≤ x < a:
V(x) = Ra
M(x) = Ra × x
For a ≤ x ≤ L:
V(x) = Ra - P
M(x) = Ra × x - P × (x - a)
Python implementation
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.patches as mpatches
def simply_supported_point_load(L, P, a):
"""
Simply supported beam with a single point load.
Parameters:
L (float): Beam span in meters
P (float): Point load in kN (positive = downward)
a (float): Distance from left support to load in meters
Returns:
x, V, M arrays
"""
b = L - a
# Support reactions
Ra = P * b / L
Rb = P * a / L
# Discretize the beam
x = np.linspace(0, L, 1000)
V = np.where(x < a, Ra, Ra - P)
M = np.where(x < a, Ra * x, Ra * x - P * (x - a))
print(f"Ra = {Ra:.2f} kN | Rb = {Rb:.2f} kN")
print(f"Mmax = {max(M):.2f} kN·m at x = {a:.2f} m")
return x, V, M
def plot_diagrams(x, V, M, title="Beam Diagrams"):
"""Plot SFD and BMD side by side."""
fig, axes = plt.subplots(3, 1, figsize=(12, 9), sharex=True)
fig.suptitle(title, fontsize=14, fontweight='bold')
L = x[-1]
# --- Beam schema ---
ax0 = axes[0]
ax0.plot([0, L], [0, 0], 'k-', linewidth=5, solid_capstyle='round')
ax0.plot(0, 0, '^', markersize=15, color='steelblue')
ax0.plot(L, 0, 'o', markersize=12, color='steelblue')
ax0.set_xlim(-0.1, L + 0.1)
ax0.set_ylim(-1, 1)
ax0.axis('off')
ax0.set_title('Beam Schema', loc='left', fontsize=10)
# --- Shear Force Diagram ---
ax1 = axes[1]
ax1.fill_between(x, V, 0,
where=(V >= 0), color='steelblue', alpha=0.4, label='V > 0')
ax1.fill_between(x, V, 0,
where=(V < 0), color='tomato', alpha=0.4, label='V < 0')
ax1.plot(x, V, color='steelblue', linewidth=1.5)
ax1.axhline(0, color='black', linewidth=0.8, linestyle='--')
ax1.set_ylabel('V (kN)')
ax1.set_title('Shear Force Diagram', loc='left', fontsize=10)
ax1.legend(fontsize=8)
ax1.grid(True, alpha=0.3)
# --- Bending Moment Diagram ---
ax2 = axes[2]
ax2.fill_between(x, M, 0,
where=(M >= 0), color='seagreen', alpha=0.4, label='M > 0 (sagging)')
ax2.fill_between(x, M, 0,
where=(M < 0), color='orange', alpha=0.4, label='M < 0 (hogging)')
ax2.plot(x, M, color='seagreen', linewidth=1.5)
ax2.axhline(0, color='black', linewidth=0.8, linestyle='--')
ax2.set_ylabel('M (kN·m)')
ax2.set_xlabel('x (m)')
ax2.set_title('Bending Moment Diagram', loc='left', fontsize=10)
ax2.legend(fontsize=8)
ax2.grid(True, alpha=0.3)
# Annotate max moment
idx_max = np.argmax(np.abs(M))
ax2.annotate(f'Mmax = {M[idx_max]:.1f} kN·m',
xy=(x[idx_max], M[idx_max]),
xytext=(x[idx_max] + 0.3, M[idx_max] * 0.8),
arrowprops=dict(arrowstyle='->', color='black'),
fontsize=9, color='darkgreen')
plt.tight_layout()
plt.savefig('bmd_point_load.png', dpi=150, bbox_inches='tight')
plt.show()
# --- Run it ---
L = 6.0 # meters
P = 50.0 # kN
a = 2.0 # meters from left support
x, V, M = simply_supported_point_load(L, P, a)
plot_diagrams(x, V, M, title=f"Simply Supported Beam — Point Load P={P} kN at a={a} m")
Output:
Ra = 33.33 kN | Rb = 16.67 kN
Mmax = 66.67 kN·m at x = 2.00 m
4. Simply Supported Beam — Uniformly Distributed Load (UDL)
Now the more common case in practice: a uniformly distributed load q (kN/m) over the full span.
q (kN/m)
↓↓↓↓↓↓↓↓↓↓↓↓↓↓↓↓↓↓
A──────────────────B
Ra Rb
Reactions and equations:
Ra = Rb = q × L / 2
V(x) = Ra - q × x
M(x) = Ra × x - q × x² / 2
Mmax = q × L² / 8 (at midspan x = L/2)
def simply_supported_udl(L, q):
"""
Simply supported beam with uniformly distributed load.
Parameters:
L (float): Beam span in meters
q (float): Distributed load in kN/m (positive = downward)
Returns:
x, V, M arrays
"""
Ra = Rb = q * L / 2
x = np.linspace(0, L, 1000)
V = Ra - q * x
M = Ra * x - q * x**2 / 2
Mmax = q * L**2 / 8
print(f"Ra = Rb = {Ra:.2f} kN")
print(f"Mmax = {Mmax:.2f} kN·m at midspan (x = {L/2:.2f} m)")
return x, V, M
# --- Run it ---
L = 8.0 # meters
q = 15.0 # kN/m
x, V, M = simply_supported_udl(L, q)
plot_diagrams(x, V, M, title=f"Simply Supported Beam — UDL q={q} kN/m, L={L} m")
Output:
Ra = Rb = 60.00 kN
Mmax = 120.00 kN·m at midspan (x = 4.00 m)
This parabolic moment diagram is what you'll see in virtually every floor beam design. At Structalis, this is the baseline we verify before moving to any reinforced concrete or steel section check under Eurocode 2 or 3.
5. Combined Loads — Point Load + UDL
Real structures rarely have a single load type. Let's combine both.
def combined_loads(L, q, P, a):
"""
Simply supported beam with UDL + point load.
By superposition (linear elasticity):
M_total(x) = M_udl(x) + M_point(x)
"""
b = L - a
# Reactions
Ra_udl = q * L / 2
Ra_P = P * b / L
Ra = Ra_udl + Ra_P
Rb_udl = q * L / 2
Rb_P = P * a / L
Rb = Rb_udl + Rb_P
x = np.linspace(0, L, 1000)
# Shear (by superposition)
V_udl = Ra_udl - q * x
V_P = np.where(x < a, Ra_P, Ra_P - P)
V = V_udl + V_P
# Moment (by superposition)
M_udl = Ra_udl * x - q * x**2 / 2
M_P = np.where(x < a, Ra_P * x, Ra_P * x - P * (x - a))
M = M_udl + M_P
print(f"Ra = {Ra:.2f} kN | Rb = {Rb:.2f} kN")
print(f"Mmax = {max(M):.2f} kN·m at x ≈ {x[np.argmax(M)]:.2f} m")
return x, V, M
# --- Run it ---
L = 7.0 # meters
q = 10.0 # kN/m
P = 30.0 # kN
a = 3.0 # meters from left support
x, V, M = combined_loads(L, q, P, a)
plot_diagrams(x, V, M, title=f"Combined Loads — UDL {q} kN/m + Point Load {P} kN at {a} m")
6. Continuous Beam — Two Spans
A continuous beam over three supports is the bread-and-butter of slab and floor beam design. We'll use the Three-Moment Equation (Clapeyron's theorem) for the analytical solution.
For two equal spans with UDL:
M0 = M2 = 0 (simply supported ends)
M1 = -q × L² / 8 (hogging moment at intermediate support)
def two_span_continuous_udl(L1, L2, q):
"""
Two-span continuous beam under UDL using Three-Moment Equation.
Parameters:
L1, L2 (float): Span lengths in meters
q (float): UDL in kN/m (same on both spans)
Returns:
x, V, M arrays (concatenated for both spans)
"""
# Three-Moment Equation: 2*M1*(L1+L2) = -q/4*(L1^3 + L2^3)
M1 = -q * (L1**3 + L2**3) / (8 * (L1 + L2))
# Reactions
Ra = q * L1 / 2 + M1 / L1
Rc = q * L2 / 2 + M1 / L2
Rb = q * (L1 + L2) - Ra - Rc
print(f"M1 (support B) = {M1:.2f} kN·m (hogging)")
print(f"Ra = {Ra:.2f} kN | Rb = {Rb:.2f} kN | Rc = {Rc:.2f} kN")
# Span 1: A → B
x1 = np.linspace(0, L1, 500)
V1 = Ra - q * x1
M_b_left = 0 # Ma = 0
M1_val = M1
# Linear correction for moment at B
M1_arr = Ra * x1 - q * x1**2 / 2 + (M1_val / L1) * x1
# Span 2: B → C
x2 = np.linspace(0, L2, 500)
V2 = (q * L2 / 2 - M1 / L2) - q * x2
M2_arr = (q * L2 / 2 - M1 / L2) * x2 - q * x2**2 / 2 + M1 * (1 - x2 / L2)
# Concatenate
x_full = np.concatenate([x1, x2 + L1])
V_full = np.concatenate([V1, V2])
M_full = np.concatenate([M1_arr, M2_arr])
return x_full, V_full, M_full
# --- Run it ---
L1 = L2 = 5.0 # meters
q = 20.0 # kN/m
x, V, M = two_span_continuous_udl(L1, L2, q)
plot_diagrams(x, V, M, title=f"Two-Span Continuous Beam — UDL {q} kN/m, L1=L2={L1} m")
Output:
M1 (support B) = -62.50 kN·m (hogging)
Ra = 37.50 kN | Rb = 125.00 kN | Rc = 37.50 kN
Notice the hogging moment at the intermediate support — this is what drives the top reinforcement in a continuous RC slab.
7. Full Reusable Module
Here's everything packaged as a clean module you can drop into any project:
# beam_diagrams.py
import numpy as np
import matplotlib.pyplot as plt
class SimpleBeam:
def __init__(self, L: float):
self.L = L
self.loads = []
self._x = np.linspace(0, L, 2000)
def add_point_load(self, P: float, a: float):
"""Add a point load P (kN) at position a (m) from left support."""
self.loads.append(('point', P, a))
return self
def add_udl(self, q: float, x_start: float = 0, x_end: float = None):
"""Add a UDL q (kN/m) from x_start to x_end."""
if x_end is None:
x_end = self.L
self.loads.append(('udl', q, x_start, x_end))
return self
def solve(self):
"""Compute reactions, V(x) and M(x) by superposition."""
x = self._x
L = self.L
V = np.zeros_like(x)
M = np.zeros_like(x)
Ra_total = 0
Rb_total = 0
for load in self.loads:
if load[0] == 'point':
_, P, a = load
b = L - a
Ra = P * b / L
Rb = P * a / L
Ra_total += Ra
Rb_total += Rb
V_i = np.where(x < a, Ra, Ra - P)
M_i = np.where(x < a, Ra * x, Ra * x - P * (x - a))
elif load[0] == 'udl':
_, q, xs, xe = load
length = xe - xs
resultant = q * length
centroid = xs + length / 2
Ra = resultant * (L - centroid) / L
Rb = resultant * centroid / L
Ra_total += Ra
Rb_total += Rb
def udl_contribution(xi, q=q, xs=xs, xe=xe, Ra=Ra):
loaded = np.clip(xi - xs, 0, xe - xs)
V_i = Ra - q * loaded
M_i = Ra * (xi - xs) * (xi >= xs) - q * loaded**2 / 2
return V_i, M_i
V_i, M_i = udl_contribution(x)
V += V_i
M += M_i
self.Ra = Ra_total
self.Rb = Rb_total
self.x = x
self.V = V
self.M = M
print(f"Ra = {Ra_total:.2f} kN | Rb = {Rb_total:.2f} kN")
print(f"Mmax = {max(M):.2f} kN·m at x = {x[np.argmax(M)]:.2f} m")
return self
def plot(self, title="Beam Analysis"):
plot_diagrams(self.x, self.V, self.M, title=title)
return self
# Usage example
if __name__ == "__main__":
beam = SimpleBeam(L=8.0)
beam.add_udl(q=12.0) # 12 kN/m full span
beam.add_point_load(P=40.0, a=3.0) # 40 kN at 3 m
beam.solve()
beam.plot(title="Custom Beam — UDL 12 kN/m + Point Load 40 kN at 3 m")
8. Going Further
This tutorial covers linear elastic static analysis — the foundation of any structural check. From here, you can:
-
Plug into Eurocode checks : use Mmax to verify the required steel area
As = Med / (fyd × z)for reinforced concrete under EC2 - Connect to OpenSeesPy for non-linear analysis and seismic loads
- Export to PDF with ReportLab to generate automated calculation notes
- Read from IFC with IfcOpenShell to extract beam geometry directly from your Revit model
If you're working on a real structural project and need expert review of your beam design — from preliminary sizing to full Eurocode justification — Structalis provides structural engineering services across France, covering reinforced concrete, steel, and timber structures.
Conclusion
In just ~100 lines of Python, we've built a tool that:
✅ Computes support reactions by superposition
✅ Generates accurate V(x) and M(x) diagrams
✅ Handles point loads, UDL, and combined loads
✅ Solves continuous beams with the Three-Moment Equation
✅ Packages everything into a reusable SimpleBeam class
The code is intentionally transparent — no black box, every formula is visible and traceable. That's exactly how structural calculation notes should work.
Found this useful? Drop a ❤️ and share with an engineer friend. Questions or improvements? Open a discussion below — I read every comment.
References:
- EN 1990: Eurocode — Basis of structural design
- EN 1992-1-1: Eurocode 2 — Design of concrete structures
- Timoshenko, S. & Gere, J. — Mechanics of Materials
- Structalis — Bureau d'études structure
Top comments (0)