2. Full-scale finite element impact simulation

Before any reduced-order model can be trusted, it must be validated against the high-fidelity reference it is meant to replace. In automotive crashworthiness that reference is a full-scale, explicit-dynamics finite element (FE) simulation of the complete body-in-white plus closures, powertrain, and occupant-restraint subsystems. This chapter explains what such a model actually is, how the explicit solver advances it in time, why contact and material nonlinearity dominate the cost, and where the computational bottlenecks concentrate. Understanding these mechanics is a prerequisite for deciding which physics a reduced-order model may legitimately simplify.

Why Explicit Dynamics Governs Crash Simulation

A frontal offset impact event lasts roughly 100–150 ms, but the stiffest structural modes and contact interactions demand stable time steps on the order of 0.5–2 µs. The semi-discrete equation of motion solved at each step is:

M · a(t) = F_ext(t) − F_int(t)
F_int = Σ_elements ∫ Bᵀ σ dV
a(t) = M⁻¹ · (F_ext − F_int)

With a lumped mass matrix, M is diagonal, so no system of equations is solved — the acceleration is obtained by simple division per degree of freedom. This is the defining advantage of explicit integration: cost per step scales linearly with element count, and no global stiffness matrix is assembled or factorized. The trade-off is conditional stability. The central-difference scheme is stable only if the time step satisfies the Courant–Friedrichs–Lewy (CFL) condition:

Δt ≤ Δt_crit = 2 / ω_max ≈ L_e / c
c = sqrt(E / ρ)   (acoustic wave speed)

For steel, c ≈ 5100 m/s. A 5 mm shell element therefore yields Δt_crit ≈ 1 µs. A 120 ms event at 1 µs requires about 120,000 time steps. This step count, multiplied by the per-step cost of a model with 5–15 million elements, is the fundamental reason full-scale crash simulation is expensive.

Contact: The Dominant Nonlinearity

Structural crashworthiness is fundamentally a contact problem. The vehicle must absorb energy through controlled folding, and that folding is driven by self-contact between sheet-metal panels, plus contact with the rigid barrier, the ground, and the occupant. LS-DYNA and comparable explicit codes use penalty-based contact, where a penetration δ generates a restoring force:

F_contact = k_penalty · δ
k_penalty = f · (K · A²) / V   (per segment)

Here K is the bulk modulus, A the segment area, V the element volume, and f a scale factor (typically 0.1–1.0). The penalty stiffness directly affects the stable time step: a very stiff penalty contact can reduce Δt_crit below the element-based value, forcing a smaller global step. Contact search itself is a bottleneck. Each step requires a broad-phase bucket sort to find candidate segment pairs, followed by a narrow-phase geometric test. For a model with 8 million shell elements, contact search can consume 30–50% of total runtime.

Material Nonlinearity and Strain-Rate Effects

Automotive steels, aluminum alloys, and polymers all exhibit rate-dependent plasticity. The Johnson–Cook constitutive model is widely used for metals:

σ = (A + B·ε_p^n) · (1 + C·ln(ε̇/ε̇₀)) · (1 − T*^m)

where A is yield stress, B and n describe hardening, C captures strain-rate sensitivity, and the thermal term T* accounts for adiabatic heating. For a typical DP780 dual-phase steel, representative values are A = 450 MPa, B = 800 MPa, n = 0.5, C = 0.015. At a strain rate of 1000 s⁻¹, the rate term alone raises flow stress by roughly 10–15%. Neglecting rate sensitivity in a reduced-order model can under-predict energy absorption by a similar margin.

Each integration point must update stress from the strain increment, return-map to the yield surface, and update internal state variables. With shell elements using 5 through-thickness integration points and 3–5 in-plane points, a 10-million-element model performs 150–250 million constitutive evaluations per step. This is the second major bottleneck after contact.

Where the Wall-Clock Time Goes

The following table summarizes the typical runtime distribution for a full-scale frontal impact model on a modern HPC cluster (example values, representative of a 10-million-element model on 256 cores):

Two structural bottlenecks deserve emphasis. First, mass scaling is often applied to raise Δt_crit for small elements, but excessive scaling adds artificial inertia and corrupts the energy balance. A practical guideline is to keep added mass below 1–2% of total model mass and to verify that kinetic energy remains a small fraction of internal energy throughout the event. Second, parallel scalability degrades beyond a few hundred cores because contact search requires global communication of segment data. Domain decomposition must balance both element count and contact surface area across ranks.

Reading the Results: Energy Balance as the Primary Health Metric

Before extracting any deformation or intrusion result, verify the global energy balance. In LS-DYNA this is reported in glstat and matsum. The total energy (internal + kinetic + contact + hourglass + sliding interface) should remain approximately constant after the initial kinetic energy is imposed. A monotonic drift indicates numerical instability. Specific thresholds used in practice:

• Hourglass energy should stay below 5% of peak internal energy. Higher values mean the shell formulation is under-integrated and the response is non-physical.
• Added mass from mass scaling should remain below 2% of total model mass.
• Contact energy should be positive and bounded; negative contact energy signals penetration or a bad penalty stiffness.
• Kinetic-to-internal energy ratio should decay smoothly as the vehicle crushes; a sudden spike means a part has been ejected or a contact has failed.

Implications for Reduced-Order Modeling

The bottlenecks described above define what a reduced-order model must preserve and what it may approximate. Contact topology and the spatial distribution of plastic hinges cannot be averaged away without losing the folding mechanism. Strain-rate sensitivity must be retained in the constitutive surrogate because it shifts flow stress by 10–15%. The time-step constraint, however, is a solver artifact — a reduced-order model based on modal decomposition or proper orthogonal decomposition (POD) can use much larger time steps because it does not resolve the acoustic wave. This is precisely the opportunity that the remainder of this course exploits: replace the expensive explicit time integration with a projection onto a low-dimensional subspace, while keeping the contact and material physics that govern the crash response.

In the next stage of the course, these full-scale simulations serve as the training and validation data source. The quality of the reduced-order model is bounded by the quality of this reference — which is why energy balance verification and contact definition review are not optional steps but the foundation of the entire modeling chain.