"""Original ideal crush-force curves; Python 3.14.4, declared NumPy/Matplotlib."""
from pathlib import Path
import numpy as np
import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt

stroke = 0.20
x = np.linspace(0, stroke, 401)
flat = np.full_like(x, 4000.0)
ramp = 8000.0 * x / stroke
assert np.isclose(np.trapezoid(flat, x), 800.0)
assert np.isclose(np.trapezoid(ramp, x), 800.0)
fig, axes = plt.subplots(1, 2, figsize=(7.4, 3.2), sharey=True, layout="constrained", facecolor="white")
for ax, force, label, color in zip(axes, [flat, ramp], ["Constant load: peak 4 kN", "Rising load: peak 8 kN"], ["#27638f", "#a34b23"]):
    ax.fill_between(x, force / 1000, color=color, alpha=.22)
    ax.plot(x, force / 1000, color=color, lw=2.2)
    ax.text(.025, 1.4, "Work = 800 J", fontsize=11)
    ax.set(xlim=(0, .2), ylim=(0, 9), xlabel="Crush displacement (m)", title=label)
    ax.set_xticks([0, .1, .2])
    ax.grid(alpha=.2)
axes[0].set_ylabel("Resisting force (kN)")
fig.suptitle("Equal absorbed energy can require different peak loads", fontsize=12)
fig.savefig(Path("crashworthiness-equal-work.png"), dpi=140, facecolor="white")