Scalar-type abstraction in generated code (enable float / fixed-point / autodiff) #84

Summary

The generated code hard-codes double as the scalar type everywhere. Introduce a scalar-type abstraction so the emitted code's number type is configurable — enabling float (memory/throughput), fixed-point, and most importantly automatic-differentiation types (dual numbers / operator-overloading AD) for parameter sensitivity and inverse problems (the plan's adjoint / model-identification follow-up).

This is the sibling of the LinearAlgebraEmitter abstraction (PR #83): same "library/representation woven into emitted text" smell, one rung lower (the element type itself).

Motivation / use cases

• Autodiff (the real prize). Re-instantiating the generated _compute with a dual-number / forward-AD scalar yields ∂(outputs)/∂(parameters) for free — calibration, sensitivity, inverse problems — without a separate adjoint emit path. Forward-mode symbolic diff already covers the algorithmic tangent (∂σ/∂ε); a scalar-type seam covers parameter derivatives cheaply.
• float for throughput/memory at per-QP scale; fixed-point for embedded targets (Phase 6 USERMAT-on-MCU-style scenarios).
• Verification: run the same recipe at double and at an AD type and cross-check the tangent against the AD-propagated derivative.

Current state — where double is baked in

• render_compute_function (recipe.h): every signature param + local — double const <param>, T& <out> (tensor) but double& <out> (scalar), the Newton iterates / residual / Jacobian locals.
• scalar_code_emit.h: emits std::pow/sin/exp/... — these must resolve for the chosen scalar (AD types rely on ADL / overloaded math; qualified std:: calls can defeat that — a real consideration).
• code_emit/linear_algebra_emitter.h: Eigen::Matrix<double,N,N> / arma::mat::fixed — the dense-solve types carry double.
• MOOSE/standalone tensor boundary: tmech::adaptor<double, …> (tmech is already templated on element type, so this side is more amenable).

Proposed design (two options to evaluate)

Option A — template the generated function on the scalar type. Emit template <class Scalar> inline void <Model>_compute(Scalar const …); the caller instantiates with double, float, or an AD type. Most flexible for AD (caller picks, even per-call). Tradeoffs: Eigen supports Matrix<Scalar,…> cleanly; Armadillo does not template on arbitrary scalars (so the Armadillo linalg backend would constrain this); std:: math must become AD-friendly (unqualified + using std::sin; or a math-shim). Interacts with the LinearAlgebraEmitter interface (its emit_newton_step hard-codes double — would need a scalar-type parameter).

Option B — a configurable concrete scalar type (per-target injection, exactly mirroring LinearAlgebraEmitter). A ScalarType / emitter advertises the type token ("double", "float", "MyADType") + any include; the target injects it; render/scalar_code_emit/LinearAlgebraEmitter consume it. Simpler, but a fixed choice per emit (can't instantiate the same output at multiple types).

Recommendation: explore A (templating) as the AD-aligned target, with B as the fallback where a concrete type is required (Armadillo backend, fixed-point). Reuse the established pattern: interface + per-target injection + advertised includes; the round-3 review's "usage-marker / per-target, not build-wide macro" lessons apply.

Scope / open questions

• Does templating force a header-only generated function (no separate .C)? MOOSE embeds the function inline already, so likely fine.
• The LinearAlgebraEmitter seam needs a scalar-type parameter (its doc already flags it as "dense double solve"); Armadillo's fixed double types bound Option A.
• AD-compatible math emission in scalar_code_emit (ADL vs qualified std::).
• Does tmech propagate an arbitrary scalar through the adaptor/expression templates end-to-end? (It's templated on element type, but worth verifying for AD types.)

Acceptance

• A recipe emits a _compute parameterized on the scalar type (Option A) or with a swapped concrete type (Option B).
• A test instantiates the generated function with a forward-AD scalar and recovers ∂output/∂parameter matching a finite-difference reference.
• The double default path stays byte-identical (no regression to existing emit).

Related

• PR Phase 3b-2b: coupled vector Newton (#35) #83 — LinearAlgebraEmitter (the pattern this mirrors; its emit_newton_step is the first consumer needing a scalar-type parameter).
• Plan follow-up: "Adjoint / reverse-mode differentiation" and "Symbolic regression / model identification" — a scalar-type seam is the cheap forward-mode enabler for the former's calibration use cases.