Irons constant-strain patch test

The canonical FE consistency test: a 2×2×2 hex patch with an interior node displaced from its regular position, subjected to boundary displacements \(u = \varepsilon_0\, x\) consistent with a uniform axial strain field, must recover \(\varepsilon = \varepsilon_0\) everywhere to machine precision.

Reference: Irons, B. M. and Razzaque, A. (1972) “Experience with the patch test for convergence of finite elements,” The Mathematical Foundations of the Finite Element Method, pp. 557-587.

A discretisation that fails this test cannot pass the inf-sup / consistency conditions and therefore cannot converge at the optimal rate. A pass at machine precision says the element kernel is wiring \(B\), \(D\), and the isoparametric mapping correctly under non-trivial geometry.

from __future__ import annotations

import numpy as np
import pyvista as pv

from femorph_solver.validation.problems import IronsPatchTest

Problem statement

problem = IronsPatchTest()
pv_ = problem.published_values[0]
print(f"{problem.name}: {problem.description}")
print(f"\n  source : {pv_.source}")
print(f"  formula: {pv_.formula}")
print(f"  target : {pv_.value} (tolerance {pv_.tolerance:.0e})")

patch_test: Irons constant-strain patch test: prescribed linear boundary displacement must reproduce a uniform strain to machine precision.

  source : Irons & Razzaque 1972
  formula: max_e |eps_xx^h(x_e) - eps0| == 0 (machine precision)
  target : 0.0 (tolerance 1e-09)

Run the check — single-shot validate()

results = problem.validate()
r = results[0]
print(f"\n  computed max_strain_error: {r.computed:.3e}")
print(f"  passes tolerance          : {r.passed}")

computed max_strain_error: 8.132e-20
passes tolerance          : True

Figure — the distorted patch

The validation’s build_model offsets the single interior node by 10 % of the edge length so the isoparametric mapping is non-trivial; rendering the resulting mesh is a nice sanity check that the geometry actually has the promised distortion.

m = problem.build_model()
plotter = pv.Plotter(off_screen=True)
plotter.add_mesh(m.grid, show_edges=True, opacity=0.7, color="#b8cde3")
plotter.add_points(
    np.asarray(m.grid.points),
    render_points_as_spheres=True,
    point_size=12,
    color="red",
)
plotter.view_isometric()
plotter.camera.zoom(1.1)
plotter.show()

Acceptance

The 1e-9 tolerance is conservative — we measure ~1e-19 on this build. A regression above 1e-12 would indicate a precision-dropping change in the assembly / solve pipeline.

assert r.passed, f"patch test failed: {r.computed:.3e} > {pv_.tolerance:.1e}"