Function kompaneets_step_coupled_inplace

pub fn kompaneets_step_coupled_inplace(
    grid: &FrequencyGrid,
    delta_n: &mut [f64],
    theta_e: f64,
    theta_z: f64,
    dtau: f64,
    dcbr: Option<&DcbrCoupling<'_>>,
    rho_coupling: Option<&RhoECoupling>,
    ws: &mut KompaneetsWorkspace,
    max_dn_abs: f64,
    max_newton_iter: usize,
) -> (bool, f64, f64)

In-place Kompaneets + DC/BR step using pre-allocated workspace.

Modifies delta_n from old values to new values. Identical physics to kompaneets_step_nonlinear_coupled but avoids per-step heap allocations.

DC/BR is handled via backward Euler within the Newton iteration: the DC/BR residual dτ × em × (neq - Δn_new) and Jacobian dτ × em are added to the Kompaneets system. Backward Euler is unconditionally stable (amplification → 0 for stiff rates), avoiding the oscillation that Crank-Nicolson would produce for the stiff DC/BR absorption at low x.

max_dn_abs is the current max|Δn| used for adaptive Newton tolerance. Pass 0.0 for the tightest tolerance (equivalent to old fixed 1e-14).

Returns (converged, rho_e, last_correction). Note: the third field is the size of the last Newton correction |δx|, not the residual |F(x)|. At convergence |δx| < tol by construction; if the Newton loop exits via max_newton_iter, last_correction is the final step the solver attempted, which only upper-bounds the residual for contractive iterations. Treat it as a diagnostic, not a proof of small residual.