// crates/ml/src/cuda_pipeline/apply_pearls_kernel.cu
//
// SP4 GPU-only Pearls A+D applicator (2026-05-01).
//
// Replaces the old host-side `apply_pearls_to_slot` helper per
// `feedback_no_cpu_compute_strict.md` and `pearl_cold_path_no_exception_to_gpu_drives.md`:
// CPU compute path is strictly forbidden — any formula (EMA / Pearl A
// sentinel detect / adaptive α blend / Wiener variance update) belongs on
// GPU regardless of frequency.
//
// One launch applies Pearls A+D to N consecutive ISV slots, each with its
// own `[sample_var, diff_var, x_lag]` Wiener-state triple. The kernel
// runs a single-thread device-side loop because each slot has its own
// state — no parallelism across slots is meaningful (N is at most 33 for
// the param-group oracle and typically 1..9).
//
// Stream-ordered with the producer kernel that wrote `scratch_buf` (no
// host sync between producer and applicator), so this is graph-capture
// compatible by construction. Pre-refactor host-side `pearls_ad_update`
// pattern broke graph capture at `aux_label_scale_ema` mid-step
// (CUDA_ERROR_STREAM_CAPTURE_UNSUPPORTED — smoke `smoke-test-v9kjv`).
//
// `ALPHA_META` is now passed as a kernel parameter so Pearl 4 can use
// a halved rate (5e-4) while all SP4 and other SP5 producers continue
// to use the SP4 default (1e-3). `EPS_DIV` remains a compile-time
// constant — it is numerical-stability-only and must never vary.
//
// Constants (`EPS_DIV`) MUST match
// `cuda_pipeline::sp4_wiener_ema::EPS_DIV` bit-for-bit.
// Single-source-of-truth lives in `sp4_wiener_ema.rs` for documentation;
// this duplicate is required because CUDA C++ can't include the Rust
// const directly. The new GPU oracle test asserts agreement within
// fp32 ULP.

extern "C" __global__ void apply_pearls_ad_kernel(
    const float* __restrict__ scratch_buf,    /* [scratch_idx_base + n_slots) device-side step observations */
    int          scratch_idx_base,
    float*       __restrict__ isv_buf,        /* [isv_idx_base + n_slots) device-side ISV */
    int          isv_idx_base,
    float*       __restrict__ wiener_state_buf, /* [wiener_offset_base + 3*n_slots) packed [sv, dv, xl] triples */
    int          wiener_offset_base,
    int          n_slots,
    float        alpha_meta               /* meta-EMA rate for Wiener variance tracking;
                                             SP4 default = 1e-3; Pearl 4 uses 5e-4 */
) {
    // Single-thread device-side loop. Every slot has its own state; no
    // cross-slot parallelism is meaningful. Guard avoids redundant work
    // if the launcher accidentally requests >1 thread.
    if (blockIdx.x != 0 || threadIdx.x != 0) return;

    constexpr float EPS_DIV    = 1.0e-8f;

    for (int i = 0; i < n_slots; ++i) {
        const int scratch_i = scratch_idx_base   + i;
        const int isv_i     = isv_idx_base       + i;
        const int wiener_i  = wiener_offset_base + 3 * i;

        const float step_obs    = scratch_buf[scratch_i];
        const float prev_x_mean = isv_buf[isv_i];
        const float sample_var  = wiener_state_buf[wiener_i + 0];
        const float diff_var    = wiener_state_buf[wiener_i + 1];
        const float x_lag       = wiener_state_buf[wiener_i + 2];

        // Degenerate-zero short-circuit: producer kernels write 0.0 when
        // their input is empty (kernel hasn't run, or every sample is
        // exactly 0). Updating Wiener state with a zero observation would
        // advance `x_lag` and suppress the next genuine first-observation
        // (Pearl A sentinel). Skip writeback entirely. Mirrors the
        // pre-refactor host-side `if step_obs == 0.0 { continue; }` guard
        // that lived on every launcher's Phase 2 loop.
        if (step_obs == 0.0f) continue;

        // Pearl A: first-observation bootstrap (Pearl A sentinel branch).
        // Both `prev_x_mean == 0` and `x_lag == 0` indicate either step 0
        // post-fold-reset or a never-touched slot. Pearl D's blend formula
        // would yield 0 at t=0 (`(1 − 0) · 0 + 0 · obs = 0`) — Pearl A's
        // direct replacement is mathematically required.
        if (prev_x_mean == 0.0f && x_lag == 0.0f) {
            isv_buf[isv_i] = step_obs;
            wiener_state_buf[wiener_i + 0] = 0.0f;
            wiener_state_buf[wiener_i + 1] = 0.0f;
            wiener_state_buf[wiener_i + 2] = step_obs;
            continue;
        }

        // Pearl D: Wiener-optimal adaptive α blend. Variances themselves
        // EMA at uniform meta-α. α* tracks the residual-vs-step variance
        // ratio, naturally landing near 0 on stationary signals (track
        // anchor) and near 1 on rapid step-changes (track new value).
        const float sample_diff = step_obs - prev_x_mean;
        const float x_lag_diff  = step_obs - x_lag;

        const float new_sample_var =
            (1.0f - alpha_meta) * sample_var + alpha_meta * sample_diff * sample_diff;
        const float new_diff_var =
            (1.0f - alpha_meta) * diff_var   + alpha_meta * x_lag_diff  * x_lag_diff;

        // EPS_DIV protects 0/0 at t=1 when both variances start near 0.
        const float alpha_opt =
            new_diff_var / (new_diff_var + new_sample_var + EPS_DIV);
        const float new_x_mean =
            (1.0f - alpha_opt) * prev_x_mean + alpha_opt * step_obs;

        isv_buf[isv_i] = new_x_mean;
        wiener_state_buf[wiener_i + 0] = new_sample_var;
        wiener_state_buf[wiener_i + 1] = new_diff_var;
        wiener_state_buf[wiener_i + 2] = step_obs;
    }

    // Make all writes globally visible before the captured graph hands
    // control to the next producer / consumer kernel. Same fence the
    // producer kernels issue after writing `scratch_buf`.
    __threadfence_system();
}
