materials: split J2 from Drucker-Prager, and give plasticity its own elastic tangent

[petlenz]

Three plasticity classes served by one templated small_strain_plasticity, all of them depending on an external elastic material for a stiffness they can only use if it is isotropic. Split into dedicated classes, each owning its own tangent.

Every step is verified against the previous implementation — bit-identical except where noted — and all pre-existing tests pass unchanged.

Performance

Two binaries built from the actual commits, run interleaved so each pair sees the same
machine state, 15 pairs, speedup computed per pair:

	old (median)	new (median)	paired speedup	range
J2	1466.2 ns/step	211.5	7.00x	5.54-8.84x
Drucker-Prager	1176.3 ns/step	652.4	1.79x	1.58-1.90x

The range is the honest figure. The same binary measured 1171-2253 ns across runs, so any
single-run comparison on this machine is worth about one significant digit; earlier drafts of
this description quoted 7.5x and 1.9x from sequential runs, which sat at the optimistic end of
these intervals.

Four independent changes, in order:

1. J2 split out. It used none of the yield-function generality and paid for all of it — a yield normal distinct from the flow normal (identical when associative), a "modified" equivalent stress carrying pressure coupling it does not have, an apex branch a cylinder cannot reach.

2. The rank-4 collapse. compute_tangent formed C_e : (dN/dsigma) : C_e explicitly — two rank-4 x rank-4 contractions per plastic step. Every flow normal here is deviatoric plus a constant volumetric term, so dN/dsigma is deviatoric in both index pairs, and for isotropic C_e:

C_e : (dN/dsigma) : C_e  ==  4 G^2 (dN/dsigma)

Verified standalone at 2.5e-16 against the explicit form, for the J2 and DP derivatives independently, before any material was touched.

3. The redundant stress contraction. The return map contracted twice — once for the trial stress, once for the returned stress. The second is sig_trial - dl (C_e : N), no contraction needed.

4. Plasticity owns its elastic tangent. See below.

The elastic dependency was dishonest, and harmful

All three classes read the stiffness from an elastic_source's tangent. That looked like flexibility and was not: the closed forms require isotropic C_e — it is what makes C_e : N = 2G N and N : C_e : N = 3G true. rk_plasticity assumed it twice, also through effective_modulus(G) = 3G. A caller wiring an anisotropic tangent would have got silently wrong answers.

It also pulled a linear_elasticity into every plasticity graph, whose own stress output is not merely unused there but wrong — it is C : eps with eps_p ignored, so once yielding starts it is not the stress of anything:

step	alpha	j2::stress	elastic::stress	error
40	0.0094	106.25	107.69	+1.4%
120	0.1094	306.25	323.08	+5.5%
200	0.2094	506.25	538.46	+6.4%

A postprocessor logging elastic::stress from a plasticity graph gets that, under a name that reads as authoritative.

All three now build C_e from moduli they hold. j2_plasticity and rk_plasticity gain a K parameter; drucker_prager_plasticity needs nothing new, since K_bulk and G were already required for the cone. A plasticity graph no longer contains an elastic material at all.

linear_elasticity is untouched and keeps its users — for a genuinely elastic model C : eps is the answer, and isotropic_damage consumes both its stress and its tangent legitimately.

This also unblocks #33

small_strain_plasticity took its yield function as a C++ object through a yield_function parameter, which the JSON reader cannot convert. That is why Drucker-Prager had to stay unregistered in #36 — a document naming it would silently get a default-constructed cone (eta = beta = k = 0) that builds, runs, never yields, and looks like elasticity.

eta, beta, K_bulk are now ordinary scalars. A DP model built entirely from JSON reproduces the C++ reference exactly:

DP from JSON: alpha=0.042067012632805115
  (C++ reference alpha = 0.042067012632805115)

Verification

• 50/50 tests, including the six-path DP tangent check and the apex branch
• J2 equivalence over 60 steps: max|dstress| = 0, max|dalpha| = 0
• DP and RK: bit-identical to all 17 digits against their pre-change output
• The one exception: dropping the redundant contraction is algebraically equal but ordered differently, so it agrees to ~1 ULP rather than bit-for-bit. Stated in that commit.
• Mutation-tested: a 0.5% error in either closed form fails the suite

A test that could not fail, fixed on the way. J2TangentTest bounded the whole run at 0.1, because the elastic->plastic transition step is genuinely inexact — a central difference straddling the yield surface averages two tangents. That made it blind: a 0.5% error in the tangent landed at 2.4e-4 and passed, while real plastic steps sit at 4.4e-10. Now bounded by regime, and the mutation that passed before now fails.

Measurement caveat

This machine's load moves absolute timings by 30%+, and I caught myself reporting a 300 ns "regression" that reversed sign under interleaving. The J2 figures are same-process A/B comparisons; the DP ones are interleaved across alternating runs. Ratios are meaningful; absolutes across commits are not.

Not done

drucker_prager_yield_function survives as an internal member rather than a template parameter — it holds the verified apex algebra, and rewriting that to save a file would trade real risk for a cosmetic gain.

Note on ordering

This targets feature/drucker-prager (#10), the root of the stack, so merging it means re-merging forward through nine descendants. Mechanical, but worth deciding alongside #10 rather than separately.

[petlenz]

Superseded — split into four PRs, because this had grown to seven commits spanning four independent arguments and could not be reviewed as one thing.

PR	argument	size
#41	one algebraic identity collapses two rank-4 contractions	+16/-4
#42	J2 does not need the yield-function generality	+213/-12
#43	with J2 gone, the template has one instantiation left — unblocks #33	+55/-74
#44	the elastic dependency is dishonest and harmful	+80/-41

Each stacks on the previous, and each was verified independently at 50/50 before pushing. #40 (the scalar-Newton split) is rebased onto #44.

The rebuilt stack is tree-identical to this branch — git diff between them is empty — so the split reorganised history without changing a byte of the result.

Closing this rather than repointing it, since its branch is the fully-combined form and keeping it open would mean two paths to the same tree.