KineticForces — replace `ximag` offset with principal-value + residue for `nutype="zero"` pole handling

[logan-nc]

Motivation

The energy integrand in src/KineticForces/EnergyIntegration.jl has a real
pole wherever the resonance denominator vanishes:

denom = i*(l_eff*wb*sqrt(x) + n*(we + wd*x)) - nu

When nutype = "zero" (collisionless) there is no nu to regularize the
integral, so the current code relies on a contour offset ximag > 0:
integrate along (1e-15, x + i*ximag) instead of the real axis. This is
inherited from Fortran PENTRC and has two defects:

1. Accuracy: the offset introduces an O(ximag²) error in the integrated
   torque. There is no principled way to pick ximag — users currently set
   it by hand, and the default 0.0 silently disables pole handling.
2. Performance: when ximag > 0 the entire [0, xmax] range is
   integrated with a ComplexF64 kernel (_energy_integrand_complex),
   paying full complex sqrt/exp cost across the whole domain even
   though the pole is typically localized to a small x-window.

A proper principal-value + residue decomposition would be both more
accurate (exact in the ε → 0 limit) and faster (the PV pieces run in
the fast Float64 kernel landed in PR #224; only the small ε-semicircle
around each pole pays the complex kernel cost).

Performance context

PR #224 landed a Float64 real-axis energy kernel that's ~2× faster than
the old complex kernel on the ximag == 0 path. That makes the ximag == 0
fast path cheap, but users who currently need pole handling for
nutype=\"zero\" are stuck paying the slow complex kernel everywhere.
PV+residue would let them stay in the Float64 kernel for most of the
range, paying the complex kernel only on the ε-arc.

Proposed implementation

Replace the ximag-based contour offset with:

∫₀^xmax f(x) dx  =  PV ∫₀^xmax f(x) dx  +  i·π·Σ Res[f, xₚ]

where the principal value is computed as PV = lim_{ε→0} [∫₀^{xₚ-ε} + ∫_{xₚ+ε}^xmax]
plus a small-ε semicircular arc contribution of i·π·Res[f, xₚ] (upper
half-plane convention consistent with Landau causality).

Pole locations xₚ are the real roots of the denominator and can be
found analytically for the quadratic-in-sqrt(x) form. Residues follow
from f(xₚ) / denom'(xₚ).

TOML changes (follow-up PR)

Replace ximag = 0.0 with x_residue_radius = 0.1 in [PENT_CONTROL].
Semantics: residue-arc radius in x = E/T units. Default 0.1.

Regression risk (important)

A prior attempt at this rewrite broke the Solovev kinetic "calculated"
fullruns test: et[1] dropped from 34.176 → 20.475. Root cause was
never identified; the source edits were lost before bisection
completed. That test has since been disabled in PR #224 (the
kinetic_source = \"calculated\" path is out of scope for the current
perturbative FFS→PE→KF work and its baseline was never physics-validated).

Before shipping the PV+residue rewrite, the follow-up PR must:

1. Re-enable the kinetic_source = \"calculated\" fullruns test.
2. Establish a physics-based baseline for that test — either against
   Fortran PENTRC on the same Solovev input, or against an analytical
   NTV limit. Do not re-pin to a number captured from the code.
3. Validate that the PV+residue result matches the ximag → 0 limit of
   the old code on a smooth-test-case regression harness case with
   nutype = \"zero\".

Cross-references

• PR KineticForces — Complete PE→KF pipeline and validate DIIID-PENTRC benchmark #224: landed the fast Float64 real-axis kernel and deferred this work
• Plan file (local): /Users/nlogan/.claude/plans/validated-exploring-clarke.md (as of the commit landing KineticForces — Complete PE→KF pipeline and validate DIIID-PENTRC benchmark #224)
• Relevant files:
  • src/KineticForces/EnergyIntegration.jl (_energy_integrand_real, _energy_integrand_complex, integrate_energy_ode)
  • src/KineticForces/KineticForcesStructs.jl (control struct with ximag field)
  • TOML: every pentrc.toml with [PENT_CONTROL]