2. Nonlinear finite element analysis essentials for optimization

Structural optimization with nonlinear CAE constraints is fundamentally different from linear-based optimization. In a linear workflow, the stiffness matrix is assembled once, and every design iteration reuses the same factorization pattern; response scales proportionally with load, and the optimizer can rely on smooth, predictable gradients. Once nonlinearity enters the picture, the relationship between design variables and structural response becomes path-dependent, history-dependent, and frequently non-smooth. Before any optimizer can be trusted, the underlying nonlinear finite element analysis (NFEA) must be understood as a solver with its own convergence behavior, limit states, and failure modes. This chapter establishes the three nonlinearity sources, the mechanics of the Newton–Raphson solution process, the meaning of convergence criteria, and the limit states that optimization constraints must capture.

2.1 The Three Sources of Nonlinearity

Nonlinearity in structural FEA arises from three distinct physical mechanisms, and in real optimization problems they usually appear in combination.

Geometric nonlinearity occurs when displacements are large enough that the strain–displacement relationship is no longer linear, or when the structure's stiffness changes due to its deformed shape. The classic indicator is the ratio of displacement to a characteristic structural dimension. As a rule of thumb, when maximum displacement exceeds roughly 5–10% of the smallest structural dimension, geometric effects become significant. Two formulations are common: total Lagrangian (TL), which references the original configuration, and updated Lagrangian (UL), which references the current configuration. TL is preferred for moderate rotations and hyperelastic materials; UL is standard for contact and large-strain metal forming. The practical consequence for optimization is that the tangent stiffness matrix K_T now depends on the displacement field u, so each Newton iteration must reassemble and refactorize.

Material nonlinearity appears whenever the stress–strain relation is not proportional. This includes plasticity (metal yielding), hyperelasticity (rubber, elastomers), viscoelasticity, and damage or softening. For metal plasticity, the key ingredients are the yield criterion (von Mises for ductile metals), the flow rule (associated or non-associated), and the hardening law (isotropic, kinematic, or combined). The return-mapping algorithm integrates the plastic strain increment at each Gauss point, and the consistent tangent modulus C_ep must be used to preserve quadratic convergence of Newton's method. Using the continuum tangent instead of the consistent tangent is a common cause of slow or stalled convergence.

Contact nonlinearity is the most disruptive for optimization because it introduces inequality constraints and non-smooth force–displacement relations. Contact status can switch between open, sticking, and sliding at any iteration, changing the active constraint set. Penalty methods, Lagrange multipliers, and augmented Lagrangian methods each trade off constraint satisfaction against conditioning. For optimization, the critical issue is that contact forces are non-differentiable at status transitions, which can break gradient-based sensitivity analysis.

Newton–Raphson Iteration with Three Nonlinearity Sources Geometric K_T(u) reassembly TL / UL formulation Material Return mapping Consistent tangent C_ep Contact Status: open/stick/slide Non-smooth forces Assemble K_T, F_ext, F_int Residual R = F_ext − F_int Solve K_T Δu = R Update u ← u + Δu Converged? ‖R‖/‖F‖ < tol Next load step or optimization iteration No: reassemble Limit States Monitored as Optimization Constraints Limit-point buckling det(K_T) → 0 Negative eigenvalue Plastic collapse Load factor λ_c Plastic strain limit Contact failure Penetration > tol Sliding instability Snap-through Displacement jump Arc-length control

2.2 The Newton–Raphson Solution Process

Nonlinear FEA solves the residual equation R(u) = F_ext − F_int(u) = 0 iteratively. At each iteration k, the tangent stiffness K_T = ∂F_int/∂u is assembled, the linear system K_T Δu = R is solved, and the displacement is updated. The process repeats until the residual norm falls below a tolerance. For optimization, three properties of this loop matter:

  • Quadratic convergence near the solution requires the consistent tangent. If convergence degrades to linear, the optimizer's finite-difference gradients become noisy and expensive.
  • Load stepping is essential when the structure exhibits limit points. A single load step to the full load will diverge; incremental loading with adaptive step size is required.
  • Line search and arc-length methods handle snap-through and snap-back by controlling the load parameter alongside the displacement increment.

A typical convergence criterion combines force and displacement norms. For example, in Abaqus the default residual tolerance is 5×10−3 relative to the time-averaged force, with a displacement correction tolerance of 10−2. In ANSYS, the force convergence criterion is 0.005 × ‖F‖ by default. These are starting points, not universal truths: for optimization, tolerances should be tightened to roughly 10−4–10−5 so that small design perturbations produce distinguishable responses.

Common pitfall: Loosening convergence tolerances to speed up optimization iterations introduces artificial noise into the design response. The optimizer may then exploit numerical artifacts rather than real structural behavior, producing designs that fail when re-analyzed with tight tolerances.

2.3 Limit States and Their Role in Optimization

A limit state is a condition beyond which the structure no longer satisfies its design requirements. In nonlinear optimization, the most relevant limit states are:

Limit State Physical Indicator Numerical Detection Optimization Constraint
Limit-point buckling Loss of stability under proportional loading Smallest eigenvalue of K_T crosses zero λmin(K_T) ≥ λcr
Plastic collapse Unbounded deformation under constant load Load–displacement curve becomes horizontal λcollapse ≥ SF × λservice
Snap-through Sudden displacement jump Negative slope in load–displacement path Arc-length solver; limit displacement
Contact penetration Physically inadmissible overlap Penetration exceeds penalty tolerance gpen ≤ gallow
Material failure Rupture, damage initiation Equivalent plastic strain exceeds limit εp,eq ≤ εallow

For gradient-based optimization, each limit state must be expressed as a differentiable constraint. Eigenvalue constraints are notoriously non-smooth when eigenvalues cross, and contact status changes introduce discontinuities. Practical strategies include smoothing the minimum eigenvalue via a p-norm aggregation, using a smoothed contact penalty, and adopting semi-analytical or adjoint sensitivity formulations that account for the active constraint set.

2.4 A Worked Example: Buckling-Constrained Column Optimization

Consider a slender steel column of length L = 2.0 m, hollow circular cross-section with outer diameter D and wall thickness t. The design variables are D and t, bounded by 40 mm ≤ D ≤ 120 mm and 2 mm ≤ t ≤ 10 mm. The objective is to minimize mass. The constraints are:

  1. Euler buckling: Pcr = π²EI / (K L)² ≥ 2.5 × Pservice, with Pservice = 80 kN and K = 1.0.
  2. Local buckling: D/t ≤ 80 for the chosen steel grade.
  3. Yield: σ = P/A ≤ σy / 1.5, with σy = 355 MPa.

In a linear eigenvalue analysis, the buckling constraint is smooth and the optimizer converges in a few iterations. However, if the column has an initial geometric imperfection of L/1000, the nonlinear load–displacement path softens, and the limit load is lower than the Euler prediction. The optimization must then use the nonlinear limit load as the constraint. A typical workflow:

# Pseudo-workflow for nonlinear buckling-constrained optimization
for each design iteration:
    update D, t
    build FE model with imperfection amplitude = L/1000
    apply incremental axial displacement
    run Newton-Raphson with arc-length control
    extract limit load P_limit from load-displacement curve
    compute constraint g = 2.5*P_service - P_limit
    compute sensitivity dg/dD, dg/dt via adjoint method
    pass g and gradients to optimizer (e.g., SQP or MMA)
    update design variables
until convergence

With D = 80 mm, t = 4 mm, the linear Euler load is approximately Pcr ≈ 210 kN, giving a safety factor of 2.6. But the nonlinear limit load with imperfection drops to roughly 175 kN, a safety factor of 2.19 — below the required 2.5. The optimizer must therefore increase t or D until the nonlinear limit load satisfies the constraint. This example illustrates why linear buckling constraints are insufficient: they can overestimate capacity by 15–25% for realistic imperfection levels.

Practical guidance: When setting up nonlinear constraints for optimization, always verify the constraint value at the initial design using both linear and nonlinear analyses. If the discrepancy exceeds 10%, the linear constraint is unreliable and the nonlinear path must be used throughout the optimization.

2.5 Convergence Troubleshooting for Optimization Loops

Convergence failures inside an optimization loop are expensive: each failed analysis wastes a full iteration. Common symptoms and remedies:

  • Residual oscillates without decreasing: Likely a contact status chattering or a non-smooth material tangent. Remedy: increase contact stabilization, use a smoothed friction law, or switch to a quasi-Newton solver (BFGS) for the global equilibrium.
  • Divergence after a few iterations: Load step too large near a limit point. Remedy: reduce initial increment, enable automatic time stepping with a minimum increment of 10−5 of the step.
  • Convergence stalls at a fixed residual level: Inconsistent tangent or insufficient penalty stiffness. Remedy: verify the consistent tangent implementation; increase penalty factor gradually.
  • Negative eigenvalues of K_T during assembly: Structure has passed a stability limit. Remedy: switch to arc-length or displacement control; do not continue with load control.

For optimization, the most robust approach is to wrap the nonlinear solver in a fallback strategy: if the primary solver fails, retry with smaller increments and a more tolerant convergence criterion, then verify the final design with a strict analysis. This two-tier approach balances speed during exploration with accuracy at the final design point.

2.6 Summary

Nonlinear FEA is not a black box that optimization can call blindly. Geometric, material, and contact nonlinearities each impose specific requirements on formulation, tangent stiffness, and solution strategy. Convergence criteria must be tight enough to produce smooth design responses, and limit states — buckling, plastic collapse, snap-through, contact penetration, and material failure — must be expressed as differentiable constraints. The next chapter builds on this foundation by examining how sensitivity analysis is extended to nonlinear problems, where the adjoint method must account for the converged nonlinear state and the active contact set.