4. Extracting injury criteria from simulation data

Raw simulation output — nodal displacements, element stresses, rigid-body accelerations — is not an injury criterion. A crashworthiness or occupant-safety model produces time histories for thousands of nodes and dozens of parts, but the metrics that regulators, rating programs, and internal design targets actually reference are a small set of derived scalars: HIC15, HIC36, NIC, chest deflection, chest compression rate, femur force, and a handful of others. The extraction step is where those scalars are computed from the time histories, and it is also where most reporting errors originate. This chapter covers the computation of the principal head, neck, and chest criteria, the sign and filtering conventions that govern them, and a reproducible Python workflow that turns solver output into a validated metric table.

Head injury: HIC and its variants

The Head Injury Criterion is defined over a resultant translational acceleration at the head center of gravity:

HIC = max over (t1, t2) of [ (t2 - t1) * ( (1/(t2-t1)) * ∫ a(t) dt )^2.5 ]

with the integral taken over the resultant acceleration a(t) in g, and the window t2 - t1 bounded by 15 ms (HIC15) or 36 ms (HIC36). Two implementation details dominate accuracy. First, the acceleration must be the resultant magnitude, not a single axis; computing HIC on the X component alone is a common and silent error. Second, the time history must be filtered per the applicable standard — SAE J211 channel frequency class 1000 for head acceleration — before the integral is evaluated. Filtering after computing HIC changes the answer and is not equivalent.

A direct double loop over all sample pairs is O(n²) and becomes slow on 10 kHz data. The standard efficient approach uses a cumulative integral and a two-pointer sweep, since for a fixed t1 the HIC integrand is monotonic in t2 up to the window limit:

import numpy as np

def hic(a_g, dt, window=0.015):
    a = np.asarray(a_g, dtype=float)
    n = len(a)
    # cumulative integral of acceleration (velocity-like term)
    cum = np.concatenate(([0.0], np.cumsum(a) * dt))
    best = 0.0
    t1_idx = 0
    for t2 in range(1, n):
        # advance t1 so that window constraint holds
        while (t2 - t1_idx) * dt > window:
            t1_idx += 1
        for t1 in range(t1_idx, t2):
            span = (t2 - t1) * dt
            if span <= 0:
                continue
            mean_a = (cum[t2] - cum[t1]) / span
            val = span * mean_a ** 2.5
            if val > best:
                best = val
    return best

In practice the inner loop is still quadratic; production code restricts t1 to a small neighborhood or uses a vectorized sliding-window formulation. The point for this course is that the window constraint and the mean-acceleration term must both be enforced explicitly — a frequent bug is to omit the span multiplier, which inflates HIC by orders of magnitude.

Neck injury: NIC and Nij

The Neck Injury Criterion (NIC) is used primarily in rear-impact and whiplash assessment and is defined on the relative acceleration and velocity between the head and torso at the occipital condyle level:

NIC = 0.2 * a_rel + v_rel^2

where a_rel is the relative horizontal acceleration (m/s²) and v_rel the relative horizontal velocity (m/s), both computed from the head and T1 accelerations integrated and differentiated consistently. The 0.2 coefficient carries units and must not be dropped. NIC is evaluated over the impact interval and reported as a peak; the time of the peak is as important as the value, because it indicates whether the load path is compression-dominated or extension-dominated.

For frontal and side impact, the more relevant neck metric is Nij, the normalized neck injury criterion combining axial force and bending moment:

Nij = Fz / Fzc + My / Myc

where Fzc and Myc are the critical intercept values from the applicable regulation (for a 50th-percentile male, approximately 6806 N in tension, 6160 N in compression, and 310 N·m in extension). Nij is computed separately for tension-extension, tension-flexion, compression-extension, and compression-flexion, and the reported value is the maximum of the four. Sign conventions for Fz and My must match the regulation's coordinate system; a flipped sign silently swaps tension for compression and produces a plausible but wrong Nij.

Chest deflection and compression rate

Chest deflection is the maximum inward displacement of the sternum relative to the spine, measured along the chest normal. In a finite element model this is not a single node value — it is the relative displacement between a sternum reference node and a spinal reference node, projected onto the chest normal direction. The extraction must therefore define both reference nodes and the projection axis explicitly:

def chest_deflection(sternum_xyz, spine_xyz, normal):
    rel = sternum_xyz - spine_xyz          # shape (n, 3)
    return rel @ np.asarray(normal)        # signed projection, mm

The sign convention matters: positive deflection is inward (compression), negative is outward. Report the peak compression and the time of peak, and also the rate of compression, since viscous criteria and some rating protocols use deflection rate. A deflection of 40 mm at 3 m/s is a different injury risk than 40 mm at 8 m/s, and the scalar alone hides that.

For the chest, the Viscous Criterion (VC) combines deflection and rate:

VC = max over t of [ (dD/dt) * (D(t) / D_max) ]

where D(t) is the instantaneous deflection and D_max the maximum deflection over the event. VC is expressed in m/s and is typically evaluated over the interval where deflection is between 0 and 40% of its peak, because the product is ill-conditioned near zero deflection.

Other metrics and the extraction pipeline

Femur force is the axial force along the femur axis, reported as peak compression; tibia index combines axial force and bending moment normalized by critical values. Pelvis force, abdominal pressure, and lower-extremity metrics follow the same pattern: identify the load path, project onto the anatomical axis, apply the standard's filter, and take the peak within the specified window.

The pipeline that ties these together has four stages, and each stage has a failure mode worth checking:

• Channel identification. Map solver output channels to metric inputs by name and part ID. A renamed node set or a re-meshed part silently breaks the mapping; assert that every expected channel exists before computing.
• Filtering. Apply the correct SAE J211 class per channel (CFC 1000 for head acceleration, CFC 600 for chest, CFC 180 for femur). Filtering must be applied to the raw time history, not to the derived scalar.
• Computation. Evaluate each criterion over its specified window with explicit sign and unit handling. Keep the time of peak alongside the value.
• Validation. Compare against a reference run or a hand-checked single-channel calculation. A metric table that has never been cross-checked is a liability, not a deliverable.

The figure below shows the data flow from solver channels through filtering and criterion computation to the validated metric table, with the sign and window conventions annotated at each stage.

Worked extraction example

Consider a frontal impact model sampled at 10 kHz over 120 ms. The head acceleration channel is a resultant in g, already filtered at CFC 1000. Running the HIC routine with a 15 ms window returns a value of 612 at t_peak = 0.058 s. The chest sternum and spine nodes are projected onto the chest normal, giving a peak inward deflection of 34.2 mm at 0.071 s and a VC of 0.41 m/s. The femur axial force peaks at 4.8 kN compression. These are illustrative values for a teaching model, not measured results; the point is the shape of the output — each metric carries a value, a time, a unit, and a sign convention, and each is traceable to a named channel and filter class.

When a value looks wrong, the diagnostic order is fixed: confirm the channel exists and is the intended one, confirm the filter class and that it was applied to the raw history, confirm the unit and sign, then confirm the window. In practice the first two checks resolve most discrepancies, and the third resolves nearly all the rest. A metric table that passes the validation gate is the correct input to the automated Excel reporting stage that follows.