Peak-wavelength / group-velocity blow-up at zero-net-momentum points (TC eye): parametric (e,m) reconstruction diverges

[mochell]

Summary

The grid→particle reconstruction c̄ = m·e / (2|m|²) (core_2D.GetVariablesAtVertex) — and the analysis-harness diagnostic, which uses the same formula — diverges where the net momentum |m| → 0 while energy e stays finite. This happens at the wind-free TC eye (and any crossing-sea / convergence zone): waves arriving from opposing directions deposit additively, so their momentum vectors cancel on the grid while energy adds. Peak wavelength λp ∝ |c̄|² then blows up.

Evidence

Stationary Holland CAT4, Rm=75 km, dx=10 km, dt=5 min, B-spline order P=2. The peak wavelength at the eye cell (x≈0, y≈0) reaches ~1e30 m (t=0 transient) and sits around 1e12–1e16 m for the whole run, vs. physical peak wavelengths of a few hundred metres.

Figure: lambdap_eye_timeseries_fix.png (eye λp vs time, log axis; "before" = unguarded, "after" = the partial mitigation below). Please attach from plots/hurricanes/comparison/.

Why it is not simply a missing clamp

• Energy is bounded both ways (ODESettings.log_energy_minimum/maximum); momentum has no floor in the reconstruction. minimal_state[2] (minimal momentum²) exists but is only a remesh keep/reset threshold, not a clamp on the reconstructed value.
• Note the inconsistency: the ODE physics already guards the same denominator — particle_waves_v6.jl αₚ uses max(speed(cx,cy), 1e-4). The reconstruction lacks the equivalent.
• Deeper issue: the parametric single-vector (e, m) representation cannot describe crossing seas / zero net direction at the eye — there is no single peak direction or wavelength there. A clamp bounds the symptom; the eye remains physically ill-defined.

Two unguarded sites: core_2D.GetVariablesAtVertex (remesh → simulation) and the harness diagnostic state_to_diags (output).

Potential path (partial mitigation, not a cure)

Cap |c̄| at a physical c_max by scaling the momentum vector up to m_min = e/(2·c_max) in its own direction (NOT flooring only the denominator — that sends c̄ → 0, which then blows up the ODE source αₚ ~ 1/|c̄|² and gives "step size underflow at t=0").

Prototype result (c_max = 25 m/s, applied in GetVariablesAtVertex only, diagnostics left honest):

• typical eye λp drops ~10 orders of magnitude (~1e15 → ~1e3 m),
• but periodic spikes to ~1e10 m remain (they line up with the over-remesh / pulsing resets — see Stationary forcing produces a non-steady radial 'pulsing'; instability worsens as dt shrinks (over-remeshing); Hs amplitude not dt-converged #64), and the t=0 transient (~1e30) is unchanged,
• max Hs and bulk λp are unchanged (20.75 → 20.76 m; median λp 444 → 445 m), so the cap is benign outside the eye.

So this cap is a partial internal mitigation, not a fix, and there is more to it: the residual spikes are coupled to the remesh instability, and the root limitation is the parametric model at crossing seas.

Directions to explore

1. A directional / two-component (or full-spectral) representation at convergence zones (proper fix).
2. Detect zero-net-momentum cells and fall back to a defined minimal/representative state.
3. The c_max cap above as a stopgap for simulation robustness.

Diagnostics should stay honest — do not mask this in state_to_diags; the diagnostic correctly exposes the singularity.

Related: #64 (over-remesh "pulsing"); surfaced after #60 removed the grid-axis anisotropy that previously masked it.