J2 Plasticity: Implicit Return Mapping

Rate-independent metal plasticity integrated the way finite-element codes do it at every integration point: a trial elastic step, a check of the yield function, and an implicit return to the yield surface. Two models, both in plain NumPy:

• 1D with combined linear isotropic + linear kinematic hardening, under stress or strain control, with a time-step convergence study.
• 3D J2 (von Mises) with Armstrong-Frederick nonlinear kinematic hardening and power-law isotropic hardening, under strain control, applied to a thin-walled tube in combined axial-torsion loading — the classic test for nonproportional cyclic plasticity.

pip install numpy scipy matplotlib
python run_examples.py          # every figure below, ~1 minute
python j2_plasticity.py         # built-in physics checks

Yield function $\phi = |\sigma - q| - (S_y + K\alpha)$, with backstress $q$ evolving as $\dot q = H,\dot\varepsilon^p$ and $\alpha$ the cumulative plastic strain. Each step: freeze the plastic variables, compute the trial state, evaluate $\phi_{trial}$. If it is positive the plastic increment that returns exactly to the yield surface has a closed form:

Stress-controlled cycling	Time-step convergence

The middle panel shows the discipline of the scheme: $\phi = 0$ to machine precision on every plastic step and $\phi < 0$ during elastic unloading. The convergence plot shows first-order accuracy of backward Euler.

The Armstrong-Frederick "dynamic recovery" term $-\gamma\mathbf q\dot\lambda$ makes the backstress saturate (at $2C/3\gamma$ in uniaxial tension), which is what produces realistic hysteresis loops and ratcheting.

Implicit update. Backward Euler on the backstress gives $\mathbf q_{n+1} = [\mathbf q_n + \tfrac{2}{3}C\Delta\lambda,\mathbf N]/(1+\gamma\Delta\lambda)$, and a short calculation shows the updated relative stress $\mathbf S_{n+1} - \mathbf q_{n+1}$ is parallel to

so the flow direction is known once $\Delta\lambda$ is, and the entire tensor update collapses to one scalar equation:

A tube under axial strain $\varepsilon_{zz}$ and torsional shear $\varepsilon_{z\theta}$ has stress-free lateral surfaces, $\sigma_{rr} = \sigma_{\theta\theta} = 0$. The tube examples therefore use mixed control: the axial and shear strains are imposed while the two lateral strains are solved each step so the lateral stresses vanish (integrate(..., free_lateral=True)). Prescribing all other strains as zero instead — a common simplification in coursework — adds a large spurious hydrostatic stress ($K,\varepsilon_{zz}$, about 2.7 GPa here) to the axial component.

The headline figure shows an Inconel-like material (E = 210 GPa, ν = 0.26, $S_y$ = 980 MPa, n = 1/5.1, h = 300 MPa, C = 26 000 MPa, γ = 50) under $\varepsilon_{zz} = A\sin 2\pi t$, $\varepsilon_{z\theta} = (A/2)\sin\pi t$, $A = 4S_y/E$, for 2 s. The 2:1 frequency ratio makes the loading nonproportional: the butterfly-shaped shear loop and the phase plot wandering outside the initial yield ellipse are the signatures.

Ramp to the same final strains three ways — axial first then shear, shear first then axial, or both together:

Same final strain state, stresses that differ by a factor of two and more: plasticity is path dependent, and the component loaded last carries most of the stress because the backstress has rotated toward it. This is why incremental integration along the actual history is necessary — a deformation-theory shortcut from final strains cannot get these numbers.

plasticity1d.py     1D model (stress and strain control)
j2_plasticity.py    3D J2 model with AF + power-law hardening, tube loading programs
run_examples.py     regenerates every figure
figures/
matlab/             the original 2014 MATLAB scripts, kept for reference

• J. C. Simo and T. J. R. Hughes, Computational Inelasticity, Springer, 1998 — the return-mapping framework.
• P. J. Armstrong and C. O. Frederick, "A mathematical representation of the multiaxial Bauschinger effect," CEGB Report RD/B/N731, 1966.
• J. L. Chaboche, "Constitutive equations for cyclic plasticity and cyclic viscoplasticity," Int. J. Plasticity 5, 247–302, 1989.

Aiden Azarnoush — originally developed in MATLAB for a graduate plasticity course (2014); rewritten in Python.

MIT — see LICENSE.