3. Surrogate modeling techniques for crash CAE

Crash CAE is expensive: a single full-vehicle explicit dynamic run with a detailed finite element model (often 5–15 million elements, 100–200 ms simulated time) can take 6–24 wall-clock hours on a 64-core node. Design of Experiments (DOE) studies that require hundreds of runs are therefore impractical if every sample is a full simulation. Surrogate models — cheap mathematical approximations trained on a limited set of high-fidelity runs — are the enabling technology for DOE-driven optimization and uncertainty quantification (UQ). This chapter compares the three dominant surrogate families used in crash engineering: Gaussian process regression (GPR), polynomial chaos expansion (PCE), and neural network (NN) surrogates, with particular attention to the strongly nonlinear, discontinuous, and heteroscedastic responses that characterize crashworthiness metrics.

Why crash responses are hard to surrogate

Crash outputs differ fundamentally from the smooth responses typical of structural statics. Three properties dominate surrogate selection:

  • Non-smoothness and discontinuities. Buckling transitions, spot-weld failure, and contact switching produce step-like changes in intrusion and acceleration responses. A small change in a gauge or material parameter can flip a deformation mode.
  • High-frequency content. Acceleration pulses contain energy up to several kHz; peak values are sensitive to time-step and mesh details, making them noisy targets.
  • Heteroscedasticity and non-Gaussianity. Response scatter often grows with severity, and outputs such as intrusion are bounded below by zero, violating Gaussian assumptions.

These properties mean that no single surrogate is universally best; the choice must be matched to response type, sample budget, and the intended downstream use (optimization, sensitivity analysis, or probabilistic UQ).

Gaussian process regression

A GPR surrogate models the response as a realization of a Gaussian random field: f(x) ~ GP(m(x), k(x, x′)), where the mean function m is typically constant or linear and the covariance kernel k encodes smoothness assumptions. The most common choice is the anisotropic squared-exponential (RBF) kernel with automatic relevance determination (ARD):

k(x, x') = σ_f² · exp( −½ · Σ_i ((x_i − x'_i) / ℓ_i)² ) + σ_n² · δ(x, x')

The length-scales ℓ_i are fitted by maximizing the marginal likelihood, which simultaneously performs a form of feature relevance analysis: a large ℓ_i indicates the response is insensitive to that variable. Prediction at a new point returns a posterior mean and variance, giving GPR its key advantage for UQ — an intrinsic, calibrated uncertainty estimate that grows away from training data.

For crash applications, the standard RBF kernel is often too smooth. Practitioners use Matérn kernels (ν = 3/2 or 5/2) to relax smoothness assumptions, and consider non-stationary or deep kernels when response variance changes across the design space. A practical workflow for a frontal-intrusion response with five design variables (rail gauge, material yield, thickness, trigger position, and bumper stiffness) might use 60–120 Latin hypercube samples, ARD length-scales, and a nugget term σ_n² to absorb solver noise. Typical achievable accuracy on held-out runs is R² ≈ 0.90–0.97 for intrusion-type responses, but only 0.6–0.8 for peak acceleration — a gap that reflects the high-frequency nature of the latter.

GPR scales as O(n³) in training and O(n²) per prediction due to kernel matrix inversion. Beyond roughly 1,000–2,000 samples, use sparse approximations (inducing points), or switch to a NN surrogate.

Polynomial chaos expansion

PCE represents the response as a spectral expansion in orthogonal polynomials of the input random variables:

Y(x) ≈ Σ_{α ∈ A} c_α · Ψ_α(ξ(x))

where ξ are standardized inputs, Ψ_α are multivariate orthonormal polynomials (Hermite for Gaussian inputs, Legendre for uniform, Laguerre for exponential), and the multi-index set A is truncated by total degree p. The coefficients c_α can be computed non-intrusively by projection (quadrature or regression) from the same simulation samples used for any DOE. Once fitted, PCE yields analytical moments — mean and variance come directly from the coefficients — and Sobol sensitivity indices are obtained by grouping coefficient variances by variable. This makes PCE exceptionally efficient for UQ: a single fit provides mean, standard deviation, and global sensitivity at negligible cost.

The limitation is the "curse of dimensionality": the number of terms grows as (p + d)! / (p! d!) for d inputs and degree p. For d = 10 and p = 3, that is 286 terms; for d = 20, it exceeds 1,700. Sparse PCE (using ℓ₁ regularization or adaptive basis selection) mitigates this when only a few variables and low-order interactions matter — common in crash problems where a handful of gauges dominate. However, PCE struggles with discontinuous responses: polynomial bases exhibit Gibbs-type oscillations near mode-switching boundaries, and convergence stalls. For such responses, either increase the sample count substantially, use piecewise or multi-element PCE, or prefer GPR/NN.

Neural network surrogates

Feed-forward multilayer perceptrons (MLPs) with 2–4 hidden layers and 64–256 units per layer, ReLU or tanh activations, and L2 regularization are the workhorse NN surrogate. They scale to large sample sets (10³–10⁵), handle high-dimensional inputs, and can approximate discontinuous functions given enough capacity. Training uses mini-batch Adam with early stopping on a validation split; standardization of inputs and outputs is essential because crash responses span orders of magnitude (intrusion in mm, acceleration in g).

For time-series outputs such as acceleration pulses, two architectures are common: (1) a single MLP predicting a fixed-length discretized pulse (e.g., 200 time points), and (2) a sequence model (LSTM or temporal CNN) mapping design variables to the full time history. The latter captures temporal correlation but requires more data. A practical compromise is to predict scalar crash metrics (peak acceleration, intrusion, HIC) with an MLP and reserve sequence models for pulse-shape studies.

NN surrogates provide no native uncertainty estimate. For UQ, use ensembles (train 5–10 networks with different initializations and average predictions; ensemble spread approximates epistemic uncertainty) or Monte Carlo dropout at inference. Calibration must be checked: ensemble spread frequently underestimates true error on extrapolation, so validate against held-out simulations before trusting probabilistic outputs.

Selecting and validating a surrogate

The figure below summarizes the decision logic and the characteristic behavior of each family on nonlinear crash responses.

Surrogate Selection for Nonlinear Crash Responses Decision criteria Sample budget n n < 200 → GPR / PCE n > 1000 → NN Response smoothness smooth → PCE discontinuous → GPR / NN Input dimension d d ≤ 10 → PCE feasible d > 15 → sparse PCE / NN Uncertainty needed? yes → GPR (native) NN → ensemble Validate on held-out runs: R², RMSE, PIT / coverage GPR Matérn 5/2 + ARD kernel Posterior mean ± variance Handles discontinuities O(n³) training cost Best for n < 500 Intrusion R² ≈ 0.90–0.97 PCE Hermite / Legendre basis Analytical mean & variance Sobol indices from coefficients Gibbs oscillation at jumps Best for smooth, d ≤ 10 Sparse PCE for d > 10 Neural Network MLP 2–4 layers, 64–256 units Scales to n > 1000 Uncertainty via ensemble High-dimensional inputs Response behavior design variable intrusion true response (mode switch) GPR: tracks jump, smooths slightly PCE: Gibbs oscillation near jump NN: sharp jump, needs more data Validation metrics R² ≥ 0.90 for intrusion PIT histogram uniform 95% CI coverage ≈ 0.95 Peak acceleration is harder: R² often 0.6–0.8 → filter pulse or predict integrated metrics instead

Validation must go beyond a single R². Use k-fold cross-validation on the DOE samples, and reserve an independent test set of 10–20% of runs. For UQ-oriented surrogates, check probabilistic calibration: the probability integral transform (PIT) histogram should be approximately uniform, and the empirical coverage of 90% and 95% intervals should match nominal levels. A surrogate that predicts means well but has overconfident intervals will corrupt any downstream reliability estimate.

Practical guidance and common failure modes

Three failure patterns recur in crash surrogate projects. First, extrapolation beyond the DOE box: all three families degrade sharply outside the sampled range, and GPR variance is the only reliable warning — always constrain optimization to the convex hull of the training design. Second, noise-dominated targets: peak acceleration and similar high-frequency metrics contain solver noise that no surrogate can learn; either filter the signal, use integrated metrics (HIC, 3 ms clip), or average over multiple time-step settings. Third, insufficient sampling near discontinuities: mode-switching boundaries need adaptive sampling (e.g., GPR-variance-driven infill) rather than uniform Latin hypercube alone.

Do not report surrogate accuracy from training data. Crash responses are frequently overfitted by flexible NN models; always hold out runs and, where possible, validate the final surrogate against one or two full simulations at design points selected by the optimizer.

A defensible workflow for a typical crash UQ study is: (1) define 5–10 design variables and their distributions; (2) generate 80–150 Latin hypercube samples; (3) run full simulations and extract scalar metrics; (4) fit GPR with Matérn kernel as the baseline, and PCE or NN as alternatives; (5) compare cross-validated R², RMSE, and interval coverage; (6) select the surrogate whose accuracy and calibration meet the UQ tolerance; (7) propagate input uncertainty through the surrogate via Monte Carlo (10⁵–10⁶ draws, seconds of compute) to obtain response distributions and Sobol indices. The surrogate is not the deliverable — the calibrated probabilistic response is.