4. Verify Solver Settings Against Analytical and Benchmark Solutions

Why Solver Verification Comes Before Model Validation

Before any injury prediction result can be trusted, the solver itself must be shown to reproduce solutions whose answers are already known. This is verification, and it is distinct from validation: verification asks "are we solving the equations correctly?" while validation asks "are we solving the right equations?" For an impact and injury mechanics workflow, three solver subsystems dominate the error budget — explicit time integration, contact treatment, and hourglass control for reduced-integration elements. Each can be exercised against an analytical or benchmark problem with a closed-form or published reference answer, and each has a characteristic failure signature when misconfigured.

The practical objective is a small, fast-running verification suite that you can re-execute whenever the deck, solver version, or hardware changes, producing pass/fail evidence against stated tolerances. The suite does not need to be large; it needs to be diagnostic — each case must isolate one mechanism so that a failure points to a specific setting.

Selecting Verification Problems That Map to Injury Mechanics

Choose problems whose physics overlaps the phenomena you will later predict: wave propagation through tissue-like solids, large-deformation contact between compliant bodies, and bending/shear of thin structures where reduced integration is attractive. The table below gives a compact suite with reference solutions and the solver setting each case interrogates.

Cases 1–3 have closed-form answers and run in seconds; case 4 is a published benchmark that exercises all three subsystems together and acts as the integration check. Keep the analytical cases as the fast gate and the benchmark as the confirmation.

Configuring the Time-Integration Verification Case

The 1D bar case is the cleanest probe of the explicit time step. Model a bar of length L = 0.2 m, cross-section 10 mm × 10 mm, with E = 2.0 GPa and ρ = 1000 kg/m³ (representative of a soft tissue surrogate). The analytical wave speed is c = √(2.0×10⁹ / 1000) ≈ 1414 m/s, so a stress step applied at one end should arrive at the far end after t = L/c ≈ 141 µs. Apply a step pressure at the loaded end, fix the far end, and record the axial stress at a gauge point near the far end.

The critical stable time step for a lumped-mass explicit scheme scales as Δt_cr ≈ L_e/c, where L_e is the smallest element characteristic length. With 100 elements along the bar, L_e = 2 mm, giving Δt_cr ≈ 1.4 µs. Run at a safety factor of 0.9 and confirm the solver-reported stable step matches your hand calculation within a few percent. If the solver's reported step is far smaller than expected, a single distorted element or an over-stiff contact penalty is usually the cause — check the element quality report and the contact stiffness scale before touching the time-step controls.

Two failure signatures matter here. First, if the step arrives early or the amplitude overshoots, bulk viscosity (the linear and quadratic damping terms) is too aggressive and is smearing or amplifying the front; reduce the quadratic term and re-run. Second, if the solution oscillates behind the front without settling, the step is at or beyond the stability limit — the run may complete but the result is numerical noise, not physics.

Configuring the Contact Verification Case

The Hertzian impact case isolates contact stiffness and penalty scaling. Impact a rigid sphere of radius R = 10 mm onto an elastic half-space with E = 5 MPa and ν = 0.45 (soft tissue range). The effective modulus is E* = E/(1−ν²) ≈ 6.29 MPa. For a prescribed peak indentation δ = 0.5 mm, the analytical peak force is F = (4/3)E*√R·δ^3/2 ≈ (4/3)(6.29×10⁶)(0.1)^0.5(0.5×10⁻³)^1.5 ≈ 0.94 N. Compare the simulated peak contact force against this value.

Penalty contact is sensitive to the penalty stiffness scale: too soft and the sphere penetrates visibly, under-predicting force; too stiff and the stable time step collapses or the contact chatters. Start with a penalty scale that keeps maximum penetration below 1% of the sphere radius, verify the peak force, then confirm the time step has not dropped disproportionately. If penetration exceeds the target while force is low, increase the penalty scale; if the run becomes unstable, the penalty is over-driven and the time step is being throttled by contact — reduce the scale and refine the contact surface instead.

Configuring the Hourglass-Control Verification Case

Reduced-integration solid and shell elements are prone to zero-energy hourglass modes. The simply supported beam case exposes this directly: a beam under mid-span load modeled with reduced-integration elements and no hourglass control will show a deflection that grows without bound or a visibly checkered deformation pattern, while the closed-form small-strain deflection is a fixed number. Run the beam with hourglass control enabled at a stiffness scale of 0.05–0.1 of the element's elastic stiffness, and compare mid-span deflection to the analytical value.

The diagnostic value is in the comparison: if deflection matches the analytical value with control on but diverges with control off, the control is doing its job and the scale is adequate. If deflection is stiffened relative to analytical (too low), the hourglass stiffness is over-applied and is adding artificial rigidity — reduce the scale. Record both runs so the verification evidence shows the control is necessary and correctly tuned, not merely present.

Running the Suite and Reading the Results

Execute the analytical cases first, in order, because each isolates one subsystem and a failure is unambiguous. Record for each case: the solver-reported stable time step, the peak or arrival observable, the analytical reference value, the percentage error, and pass/fail against the stated tolerance. The Taylor bar benchmark runs last; if it fails while the three analytical cases pass, the interaction between subsystems (typically contact penalty throttling the time step, or hourglass control interacting with contact) is the suspect, and you should re-examine the combined settings rather than any single one.

A useful discipline is to keep the verification deck under version control alongside the solver version string and hardware description, so that a later change in any of these can be traced to a specific verification result. When a case fails, the failure signature — early wave arrival, low contact force with visible penetration, or a checkered deformation pattern — tells you which setting to adjust before you re-run, rather than prompting a blind sweep of parameters.

Verification Record Template

The following record captures the minimum evidence needed to demonstrate that solver settings reproduce known solutions. Fill one row per case and retain it with the model version.

Once all four cases pass, the solver settings for time integration, contact, and hourglass control are verified against known solutions, and the model is ready to move from verification into material and coupon-level validation. The verification record becomes the baseline against which any future solver upgrade or deck change is re-checked.