3. Local sensitivity methods in CAE

Local sensitivity analysis answers a precise question: if a single input parameter p changes by an infinitesimal amount, how does a scalar response J (stress, compliance, natural frequency, peak temperature) change? The answer is the gradient ∂J/∂p evaluated at the current design point. In CAE practice this gradient is almost never available in closed form, because J is the output of a numerical solver. The three dominant strategies for recovering it are finite difference (FD), the complex-step derivative (CSD), and the adjoint method. Each trades accuracy against cost in a different way, and choosing among them is one of the most consequential decisions in a sensitivity study.

Finite difference: the baseline everyone starts from

The forward finite-difference approximation perturbs the parameter by a step h and differences the response:

∂J/∂p ≈ ( J(p + h) − J(p) ) / h        forward difference
∂J/∂p ≈ ( J(p + h) − J(p − h) ) / (2h)   central difference

The appeal is universality: any solver that accepts a parameter change can be wrapped in an FD loop, with no source-code access required. The cost is one full solve per perturbed parameter (two for central differences), so a model with 40 design variables needs 40–80 re-analyses per gradient evaluation. Accuracy is governed by two competing errors. Truncation error scales with h (forward) or h² (central), while subtractive cancellation error scales with machine epsilon divided by h. For double precision, a forward step around 10⁻⁶–10⁻⁷ of the parameter magnitude is a reasonable compromise; central differences tolerate larger steps and deliver roughly twice the accuracy at twice the cost.

The practical hazard is step-size sensitivity in noisy responses. Contact-status switching, adaptive remeshing, and iterative solver tolerances all introduce discontinuities that can make an FD gradient jump erratically as h shrinks. A robust workflow is to sweep h over several decades (e.g., 10⁻² to 10⁻⁸ relative) and plot the gradient; a stable plateau indicates a trustworthy value, while monotonic drift signals cancellation error and scatter signals solver noise. Always normalize the step to the parameter scale — a fixed absolute step of 1 mm is meaningless for a 0.02 mm fillet radius.

Complex-step differentiation: near-exact derivatives without cancellation

The complex-step method evaluates the response at a complex-perturbed parameter and extracts the imaginary part:

∂J/∂p ≈ Im[ J(p + i·h) ] / h

Because the subtraction of two nearly equal real numbers is eliminated, there is no subtractive cancellation. The step h can be made extremely small (10⁻²⁰ or smaller), and the result is accurate to machine precision — typically 10–15 significant digits, far beyond what FD can deliver. This makes CSD attractive for verifying analytic gradients or for problems where FD noise is unacceptable.

The catch is that the entire analysis chain must remain analytic in the complex plane. Every operation the solver performs on the perturbed parameter — constitutive laws, element formulations, contact logic, material tables — must support complex arithmetic. Branch cuts, absolute values, min/max comparisons, and table lookups with real-only interpolation all break the method or silently corrupt the imaginary part. Commercial solvers rarely expose a complex-step mode, so CSD is most often implemented in custom scripts, reduced-order models, or research codes where the full call stack is under control. When it is available, it is the cheapest high-accuracy option: one complex solve yields the derivative, versus two real solves for central FD.

Adjoint formulations: gradients for many parameters at the cost of one solve

When the number of design variables is large — topology optimization with thousands of element densities, or shape optimization with a dense surface parameterization — FD and CSD become prohibitive because their cost scales with the number of parameters. The adjoint method reverses the economics. For a response J(u, p) subject to the discretized equilibrium R(u, p) = 0, the total derivative is

dJ/dp = ∂J/∂p + λᵀ · ∂R/∂p

where the adjoint vector λ solves the linear system

(∂R/∂u)ᵀ · λ = −(∂J/∂u)ᵀ

The key observation: the adjoint system depends on the response J but not on the parameter p. One adjoint solve therefore serves every parameter in the gradient, and the per-parameter cost collapses to a sparse matrix–vector product ∂R/∂p · λ. For a model with 5,000 design variables, the adjoint gradient costs roughly one extra linear solve plus 5,000 cheap products, whereas FD would demand 5,000 full nonlinear re-analyses. The trade-off is implementation complexity: the adjoint requires the transposed tangent stiffness matrix and the partial derivatives ∂R/∂p and ∂J/∂u, which usually means either intrusive solver access or automatic differentiation of the residual assembly.

Putting the methods to work: a verification workflow

A defensible sensitivity study does not simply pick one method and trust it. The standard practice is to establish a reference gradient with the most accurate method available, then validate the production method against it. Consider a bracket model with a fillet radius r = 4.0 mm and a peak von Mises stress response J ≈ 212 MPa. A central FD sweep with relative steps from 10⁻² to 10⁻⁸ might yield gradients of 18.4, 17.9, 17.9, 17.8, 16.1 MPa/mm — the plateau at 17.8–17.9 identifies the converged value, while the collapse at 10⁻⁸ exposes cancellation error. A complex-step evaluation of the same model, if the solver supports it, should reproduce 17.85 MPa/mm to many digits. If an adjoint implementation is available, its result should match the CSD value to discretization accuracy; a mismatch of more than a few percent usually points to a missing term in ∂R/∂p — commonly a parameter that enters through a boundary condition or a material table rather than the element stiffness.

Two failure modes deserve explicit checks. First, inconsistent perturbation: in FD, the parameter must be changed everywhere it appears — geometry, mesh, loads, and material — or the gradient reflects only part of the physics. Second, non-smooth response: if the response is dominated by a max() operation (peak stress location switching between nodes), all three methods return a one-sided or undefined derivative at the switch point. The remedy is to smooth the response (p-norm or Kreisselmeier–Steinhauser aggregation) before differentiating, not to shrink the step.

In summary, finite difference remains the universal fallback and the natural first check; complex step delivers near-exact derivatives when the analysis chain is complex-safe; and the adjoint method is the only scalable choice when the parameter count is large. Matching the method to the parameter dimension, the available solver access, and the required accuracy is what separates a sensitivity study that informs design from one that merely produces numbers.