Source code for telemetry_anomdet.models.deep.codegen

# src/telemetry_anomdet/models/deep/codegen.py

"""
C code generation for a distilled KANGDN.

Takes the torch-free specification produced by
:func:`~telemetry_anomdet.models.deep.distill.extract_kan_gdn` and emits a
self-contained C translation unit: every learned parameter becomes a
``static const`` array and the scoring path becomes ordinary loops. The result
has no dependencies beyond ``math.h``, allocates nothing, and exposes two entry
points, so it ports to flight hardware without modification.

The generated code mirrors
:class:`~telemetry_anomdet.models.deep.distill.KANGDNNumpy` statement for
statement, which is what makes it checkable: the same window scored by the torch
detector, the NumPy evaluator and the compiled C must agree.

Three properties of the distilled form are exploited:

* **The graph is frozen**, so neighbours are emitted as explicit index lists and
  the attention softmax runs over ``topk + 1`` terms per node rather than over
  all nodes.
* **Attention is separable**, so each node needs two dot products
  (``a_src . g_i`` and ``a_dst . g_i``) instead of a pairwise concatenation.
* **The spline knots are fixed**, so every Cox-de Boor denominator is known at
  generation time. They are tabulated per recursion level and index, which turns
  each division in the inner loop into a multiply and keeps the result exact for
  any knot vector rather than only for evenly spaced ones.

Coding standard
---------------

The output follows the NASA/JPL Power of Ten rules for safety critical code:

1. No ``goto``, ``setjmp``, or recursion.
2. Every loop bound is a compile-time constant.
3. No allocation. Scratch buffers are function-scope ``static``.
4. Functions stay short; the scoring path is split into single-purpose steps.
5. At least two assertions per function, routed through ``KANGDN_ASSERT``.
6. Data is declared at the smallest workable scope; no file-scope mutable state.
7. Every function returns a status code and validates its parameters, and every
   call site checks the returned status.
8. The preprocessor is limited to include guards, the interface dimensions, and
   the assertion hook. Knot vectors are ``static const`` data rather than
   function-like macros, and internal dimensions are enumerations.
9. Pointer use is limited to a single dereference. No function pointers.
10. The output is intended to compile clean under
    ``-Wall -Wextra -Werror -pedantic``.

Integrators should define ``KANGDN_ASSERT`` to route failures into the mission's
fault handler. The default delegates to ``<assert.h>``, which aborts, and is
appropriate for ground testing rather than flight.

Numeric type is selectable. ``float`` matches the single-precision FPUs on the
relevant targets and is the deployment mode; ``double`` reproduces the NumPy
evaluator to near machine precision and exists to separate porting errors from
precision loss.
"""

from __future__ import annotations

from string import Template

import numpy as np

from ...version import __version__

# Banner shared by every generated file, substituted in as ${provenance}.
#
# Deliberately no generation timestamp. It records the generator version and
# nothing else, so regenerating from the same fitted detector reproduces the
# sources byte for byte and a diff shows only real changes.
_PROVENANCE = f"""\
 * Generated by telemetry_anomdet {__version__}. Do not edit; regenerate from
 * the fitted detector instead.
 *
 * The generated interface is unstable until telemetry_anomdet 1.0.0: symbol
 * names, macros and the entry point signatures may change in any 0.x release.
 * Record the version above with the artifact, and regenerate rather than
 * hand-patching when upgrading."""

_HEADER = Template("""\
/* ${prefix}.h : distilled KANGDN detector.
 *
${provenance}
 *
 * Entry points take one window of raw, unscaled telemetry: ${window_size}
 * timesteps by ${n_features} channels, row major, so element (t, f) lives at
 * index t * ${up}_N_FEATURES + f. Scaling, forecasting, the deviation
 * normalisation and the threshold comparison all happen internally.
 *
 *   ${prefix}_real_t win[${up}_WINDOW_SIZE][${up}_N_FEATURES];
 *   ${prefix}_real_t score;
 *   if (${prefix}_score(&win[0][0], &score) == ${up}_OK) { ... }
 *
 * The window is a flat pointer rather than a two-dimensional array parameter on
 * purpose. ISO C before C23 forbids passing a pointer to a non-const array
 * where a pointer to a const array is expected, so a caller filling a mutable
 * buffer from sensor readings could not pass it without a cast.
 *
 * Both functions return a status code and deliver their result through an out
 * parameter. Callers must check the status: a non-zero value means the result
 * was not written and must not be acted upon.
 */
#ifndef ${guard}
#define ${guard}

/* Generator version, so a deployed artifact can report its own provenance. */
#define ${up}_VERSION "${version}"

#define ${up}_WINDOW_SIZE ${window_size}
#define ${up}_N_FEATURES  ${n_features}
#define ${up}_THRESHOLD   ${threshold}

/* Status codes. */
#define ${up}_OK        0
#define ${up}_ERR_NULL  1
#define ${up}_ERR_RANGE 2

typedef ${real} ${prefix}_real_t;

/* The detector is C, but embedded projects are frequently C++ (Arduino among
   them), so give it C linkage when included from a C++ translation unit. */
#ifdef __cplusplus
extern "C" {
#endif

/* Graph deviation score for one window; higher means more anomalous. */
int ${prefix}_score(const ${prefix}_real_t *window, ${prefix}_real_t *out_score);

/* Sets *out_flag to 1 when the score exceeds the distilled threshold. */
int ${prefix}_is_anomaly(const ${prefix}_real_t *window, int *out_flag);

#ifdef __cplusplus
}
#endif

#endif /* ${guard} */
""")

_PRELUDE = Template("""\
/* ${prefix}.c : distilled KANGDN detector.
 *
${provenance}
 *
 * Follows the NASA/JPL Power of Ten rules. No allocation, no recursion, fixed
 * loop bounds, status codes on every function, and no file-scope mutable state.
 * Learned parameters are static const and therefore flash resident.
 */
#include "${prefix}.h"
#include <math.h>
#include <stddef.h>

/* Assertion hook. Define ${up}_ASSERT before building to route failures into a
   mission fault handler; the default aborts and suits ground testing only. */
#ifndef ${up}_ASSERT
#include <assert.h>
#define ${up}_ASSERT(cond) assert(cond)
#endif

enum {
    N_NODES = ${n_nodes},
    EMBED   = ${embed_dim},
    WINDOW  = ${window},
    N_NBR   = ${n_nbr},
    N_SCORE = ${n_score}
};

/* Channels allowed to raise an alarm. Every channel still feeds the graph and
   the forecast; only these contribute to the deviation score. */
${score_channels}

static const ${real} LEAKY = ${leaky};
/* Any score beyond this is treated as corruption rather than an anomaly. */
static const ${real} SCORE_LIMIT = ${score_limit};

/* exp(-x), guarded on both sides. */
static ${real} kan_exp_neg(${real} x)
{
    ${real} y;
    ${up}_ASSERT(isfinite(x) != 0);
    y = ${expf}(-x);
    ${up}_ASSERT(y >= (${real})0);
    return y;
}
""")

_LAYER = Template("""
/* ---- KAN layer ${name}: ${in_f} -> ${out_f} ---- */
enum {
    ${up}_IN     = ${in_f},
    ${up}_OUT    = ${out_f},
    ${up}_ORDER  = ${order},
    ${up}_NKNOTS = ${n_knots},
    ${up}_NCOEFF = ${n_coeff}
};

${knots}${inv_steps}${base_w}${spline_w}
/* B-spline bases at x via the Cox-de Boor recursion. Every denominator is a
   precomputed reciprocal, so the inner loop multiplies and never divides. */
static int ${name}_bsplines(${real} x, ${real} out[${up}_NCOEFF])
{
    ${real} b[${up}_NKNOTS - 1];
    int i;
    int k;
    int n;

    if (out == NULL) {
        return ${up_prefix}_ERR_NULL;
    }
    ${up_prefix}_ASSERT(isfinite(x) != 0);
    ${up_prefix}_ASSERT(${up}_NCOEFF == (${up}_NKNOTS - 1 - ${up}_ORDER));

    for (i = 0; i < (${up}_NKNOTS - 1); ++i) {
        b[i] = ((x >= ${name}_knots[i]) && (x < ${name}_knots[i + 1]))
             ? (${real})1 : (${real})0;
    }

    n = ${up}_NKNOTS - 1;
    for (k = 1; k <= ${up}_ORDER; ++k) {
        const ${real} *inv = &${name}_inv_span[(k - 1) * (${up}_NKNOTS - 1)];
        ${up_prefix}_ASSERT(n >= 2);
        /* b[i] is written before b[i + 1] is read on the following step, so the
           update is safe in place. The right-hand span at index i is the
           left-hand span at index i + 1, so one table serves both. */
        for (i = 0; i < (n - 1); ++i) {
            b[i] = ((x - ${name}_knots[i]) * inv[i] * b[i])
                 + ((${name}_knots[i + k + 1] - x) * inv[i + 1] * b[i + 1]);
        }
        --n;
    }

    for (i = 0; i < ${up}_NCOEFF; ++i) {
        out[i] = b[i];
    }
    return ${up_prefix}_OK;
}

/* y[o] = sum_i ( base_w[o][i] * silu(x[i]) + sum_c spline_w[o][i][c] * B_c(x[i]) ) */
static int ${name}_forward(const ${real} *x, ${real} *y)
{
    ${real} basis[${up}_NCOEFF];
    int i;
    int o;
    int c;
    int status;

    if ((x == NULL) || (y == NULL)) {
        return ${up_prefix}_ERR_NULL;
    }
    ${up_prefix}_ASSERT(${up}_IN > 0);
    ${up_prefix}_ASSERT(${up}_OUT > 0);

    for (o = 0; o < ${up}_OUT; ++o) {
        y[o] = (${real})0;
    }
    for (i = 0; i < ${up}_IN; ++i) {
        const ${real} xi = x[i];
        const ${real} s = xi / ((${real})1 + kan_exp_neg(xi));
        status = ${name}_bsplines(xi, basis);
        if (status != ${up_prefix}_OK) {
            return status;
        }
        for (o = 0; o < ${up}_OUT; ++o) {
            const ${real} *sw = &${name}_spline_w[(o * ${up}_IN + i) * ${up}_NCOEFF];
            ${real} acc = ${name}_base_w[o * ${up}_IN + i] * s;
            for (c = 0; c < ${up}_NCOEFF; ++c) {
                acc += sw[c] * basis[c];
            }
            y[o] += acc;
        }
    }
    return ${up_prefix}_OK;
}
""")

_SCORING = Template("""
/* ---- scoring path ---- */

/* Buffers are flat with explicit index arithmetic. ISO C before C23 rejects
   passing a pointer to a non-const array where a pointer to a const array is
   expected, so two-dimensional parameters cannot be const qualified. Flat
   arrays keep const correctness and match the layout of the parameter tables. */

/* Scale the window and split it into per-node context and the final target. */
static int kan_prepare(const ${prefix}_real_t *window, ${real} *context, ${real} *target)
{
    int t;
    int f;

    if ((window == NULL) || (context == NULL) || (target == NULL)) {
        return ${up}_ERR_NULL;
    }
    ${up}_ASSERT(N_NODES == ${up}_N_FEATURES);
    ${up}_ASSERT(WINDOW == (${up}_WINDOW_SIZE - 1));

    for (f = 0; f < N_NODES; ++f) {
        for (t = 0; t < WINDOW; ++t) {
            context[(f * WINDOW) + t] = ${scale_ctx};
        }
        target[f] = ${scale_tgt};
    }
    return ${up}_OK;
}

/* h_i = W x_i + b */
static int kan_features(const ${real} *context, ${real} *node_h)
{
    int i;
    int d;
    int t;

    if ((context == NULL) || (node_h == NULL)) {
        return ${up}_ERR_NULL;
    }
    ${up}_ASSERT(EMBED > 0);
    ${up}_ASSERT(WINDOW > 0);

    for (i = 0; i < N_NODES; ++i) {
        for (d = 0; d < EMBED; ++d) {
            ${real} acc = feat_b[d];
            for (t = 0; t < WINDOW; ++t) {
                acc += feat_w[(d * WINDOW) + t] * context[(i * WINDOW) + t];
            }
            node_h[(i * EMBED) + d] = acc;
        }
    }
    return ${up}_OK;
}

/* The attention score a^T [g_i | g_j] is separable, so each node needs only the
   two projections of g_i = [v_i | h_i] rather than a pairwise concatenation. */
static int kan_attention(const ${real} *node_h, ${real} *att_src, ${real} *att_dst)
{
    int i;
    int d;

    if ((node_h == NULL) || (att_src == NULL) || (att_dst == NULL)) {
        return ${up}_ERR_NULL;
    }
    ${up}_ASSERT(N_NODES > 0);
    ${up}_ASSERT(EMBED > 0);

    for (i = 0; i < N_NODES; ++i) {
        ${real} s = (${real})0;
        ${real} u = (${real})0;
        for (d = 0; d < EMBED; ++d) {
            const ${real} v = embedding[(i * EMBED) + d];
            const ${real} h = node_h[(i * EMBED) + d];
            s += (attn_src[d] * v) + (attn_src[EMBED + d] * h);
            u += (attn_dst[d] * v) + (attn_dst[EMBED + d] * h);
        }
        att_src[i] = s;
        att_dst[i] = u;
    }
    return ${up}_OK;
}

/* Softmax over one node's frozen neighbour list. */
static int kan_weights(int node, const ${real} *att_src, const ${real} *att_dst,
                       ${real} *weights)
{
    ${real} maxlog;
    ${real} total = (${real})0;
    int n;

    if ((att_src == NULL) || (att_dst == NULL) || (weights == NULL)) {
        return ${up}_ERR_NULL;
    }
    ${up}_ASSERT(node >= 0);
    ${up}_ASSERT(node < N_NODES);

    for (n = 0; n < N_NBR; ++n) {
        const ${real} v = att_src[node] + att_dst[neighbours[(node * N_NBR) + n]];
        weights[n] = (v >= (${real})0) ? v : (LEAKY * v);
    }
    maxlog = weights[0];
    for (n = 1; n < N_NBR; ++n) {
        if (weights[n] > maxlog) {
            maxlog = weights[n];
        }
    }
    for (n = 0; n < N_NBR; ++n) {
        weights[n] = ${expf}(weights[n] - maxlog);
        total += weights[n];
    }
    ${up}_ASSERT(total > (${real})0);
    for (n = 0; n < N_NBR; ++n) {
        weights[n] /= total;
    }
    return ${up}_OK;
}

/* z_i = activation( sum_j alpha_ij h_j ) * v_i, over the frozen neighbour list. */
static int kan_aggregate(const ${real} *node_h, const ${real} *att_src,
                         const ${real} *att_dst, ${real} *node_z)
{
    ${real} agg[EMBED];
    ${real} weights[N_NBR];
    int i;
    int d;
    int n;
    int status;

    if ((node_h == NULL) || (att_src == NULL) || (att_dst == NULL) || (node_z == NULL)) {
        return ${up}_ERR_NULL;
    }
    ${up}_ASSERT(N_NBR > 0);
    ${up}_ASSERT(N_NBR <= N_NODES);

    for (i = 0; i < N_NODES; ++i) {
        status = kan_weights(i, att_src, att_dst, weights);
        if (status != ${up}_OK) {
            return status;
        }
        for (d = 0; d < EMBED; ++d) {
            agg[d] = (${real})0;
        }
        for (n = 0; n < N_NBR; ++n) {
            const int j = (int)neighbours[(i * N_NBR) + n];
            for (d = 0; d < EMBED; ++d) {
                agg[d] += weights[n] * node_h[(j * EMBED) + d];
            }
        }
        status = act_forward(agg, &node_z[i * EMBED]);
        if (status != ${up}_OK) {
            return status;
        }
        for (d = 0; d < EMBED; ++d) {
            node_z[(i * EMBED) + d] *= embedding[(i * EMBED) + d];
        }
    }
    return ${up}_OK;
}

/* Forecast each node, normalise its deviation, and take the maximum. */
static int kan_deviation(const ${real} *node_z, const ${real} *target, ${real} *out_score)
{
    ${real} best = (${real})0;
    int k;
    int status;

    if ((node_z == NULL) || (target == NULL) || (out_score == NULL)) {
        return ${up}_ERR_NULL;
    }
    ${up}_ASSERT(N_SCORE > 0);
    ${up}_ASSERT(SCORE_LIMIT > (${real})0);

    for (k = 0; k < N_SCORE; ++k) {
        const int i = (int)score_channels[k];
        ${real} pred = (${real})0;
        ${real} err;
        ${real} normed;

        status = head_forward(&node_z[i * EMBED], &pred);
        if (status != ${up}_OK) {
            return status;
        }
        err = ${fabsf}(pred - target[i]);
        normed = ${fabsf}(err - err_median[i]) / (err_iqr[i] + (${real})1e-9);
        if (normed > best) {
            best = normed;
        }
    }

    /* Corruption checks: a non-finite or implausible score must not be treated
       as a detection. */
    if (isfinite(best) == 0) {
        return ${up}_ERR_RANGE;
    }
    if (best > SCORE_LIMIT) {
        return ${up}_ERR_RANGE;
    }
    *out_score = best;
    return ${up}_OK;
}

int ${prefix}_score(const ${prefix}_real_t *window, ${prefix}_real_t *out_score)
{
    /* Scratch lives at function scope, static to keep it off the stack. */
    static ${real} context[N_NODES * WINDOW];
    static ${real} target[N_NODES];
    static ${real} node_h[N_NODES * EMBED];
    static ${real} node_z[N_NODES * EMBED];
    static ${real} att_src[N_NODES];
    static ${real} att_dst[N_NODES];
    int status;

    if ((window == NULL) || (out_score == NULL)) {
        return ${up}_ERR_NULL;
    }
    ${up}_ASSERT(N_NODES > 0);
    ${up}_ASSERT(EMBED > 0);

    status = kan_prepare(window, context, target);
    if (status != ${up}_OK) {
        return status;
    }
    status = kan_features(context, node_h);
    if (status != ${up}_OK) {
        return status;
    }
    status = kan_attention(node_h, att_src, att_dst);
    if (status != ${up}_OK) {
        return status;
    }
    status = kan_aggregate(node_h, att_src, att_dst, node_z);
    if (status != ${up}_OK) {
        return status;
    }
    return kan_deviation(node_z, target, out_score);
}

int ${prefix}_is_anomaly(const ${prefix}_real_t *window, int *out_flag)
{
    ${prefix}_real_t score = (${prefix}_real_t)0;
    int status;

    if ((window == NULL) || (out_flag == NULL)) {
        return ${up}_ERR_NULL;
    }
    ${up}_ASSERT(SCORE_LIMIT > (${real})0);

    status = ${prefix}_score(window, &score);
    if (status != ${up}_OK) {
        return status;
    }
    ${up}_ASSERT(isfinite(score) != 0);
    *out_flag = (score > (${real})${up}_THRESHOLD) ? 1 : 0;
    return ${up}_OK;
}
""")


_VECTORS = Template("""\
/* ${prefix}_vectors.h : golden vectors for the distilled detector.
 *
${provenance}
 *
 * Each row is one window, flattened row major to match the layout
 * ${prefix}_score expects. The expected scores are the reference evaluator's
 * output in double precision, so a port is validated by scoring each window and
 * comparing, with no host tooling required on the target.
 *
 * Compare as |got - want| <= ATOL + RTOL * |want|. The absolute term matters:
 * the score is a difference of differences, so a score far below the threshold
 * carries large relative error while its absolute error stays negligible. The
 * expected flags are the decisive check, since they are what a port must
 * reproduce for the detector to behave identically.
 */
#ifndef ${guard}
#define ${guard}

#include "${prefix}.h"

enum { ${up}_N_VECTORS = ${n_vectors} };

/* Tolerances appropriate to the generated numeric type. */
static const double ${up}_TEST_RTOL = ${rtol};
static const double ${up}_TEST_ATOL = ${atol};

static const ${prefix}_real_t ${up}_TEST_WINDOWS[${up}_N_VECTORS]
                                                [${up}_WINDOW_SIZE * ${up}_N_FEATURES] = {
${windows}
};

static const double ${up}_TEST_SCORES[${up}_N_VECTORS] = {
${scores}
};

/* 1 where the reference evaluator flags the window as anomalous. */
static const int ${up}_TEST_FLAGS[${up}_N_VECTORS] = {
${flags}
};

#endif /* ${guard} */
""")


def _fmt(value: float, digits: int, suffix: str) -> str:
    """Format one scalar as a C floating point literal."""
    if not np.isfinite(value):
        raise ValueError(f"Non-finite parameter in the distilled model: {value}")
    text = f"{float(value):.{digits}g}"
    if not any(c in text for c in ".eE"):
        text += ".0"
    return text + suffix


def _c_array(name: str, values, real: str, digits: int, suffix: str, per_line: int = 6) -> str:
    """Emit a flat ``static const`` array declaration for any-rank ``values``."""
    flat = np.asarray(values, dtype=float).reshape(-1)
    literals = [_fmt(v, digits, suffix) for v in flat]
    lines = [
        "    " + ", ".join(literals[start : start + per_line])
        for start in range(0, len(literals), per_line)
    ]
    body = ",\n".join(lines) if lines else "    0"
    return f"static const {real} {name}[{flat.size}] = {{\n{body}\n}};\n"


def _layer_knots(grid) -> np.ndarray:
    """
    Return the single knot vector a KAN layer uses for every input feature.

    One vector is emitted per layer, so the knots must agree across features.
    They are exported verbatim rather than reconstructed from a start and a
    spacing: the grid is built in float32, so a reconstruction differs from the
    stored values by roughly 1e-7 relative, which would show up as a spurious
    mismatch against the reference evaluator.

    Raises
    ------
    NotImplementedError
        If the knot vector differs between input features.
    """
    grid = np.asarray(grid, dtype=float)
    if not np.allclose(grid, grid[0], rtol=1e-9, atol=1e-12):
        raise NotImplementedError(
            "C generation requires one knot vector shared across input features."
        )
    return grid[0]


def _inverse_spans(knots: np.ndarray, order: int) -> np.ndarray:
    """
    Reciprocal Cox-de Boor denominators, one row per recursion level.

    Entry ``[k - 1, i]`` is ``1 / (knots[i + k] - knots[i])``, the divisor the
    recursion applies at level ``k`` and index ``i``. Precomputing these turns
    every division in the inner loop into a multiply, and tabulating per index
    rather than per level keeps the result exact for any knot vector instead of
    only for evenly spaced ones. Entries with no valid span are zero and are
    never read.
    """
    width = len(knots) - 1
    table = np.zeros((order, width))
    for k in range(1, order + 1):
        for i in range(width):
            if i + k < len(knots):
                span = knots[i + k] - knots[i]
                if span != 0.0:
                    table[k - 1, i] = 1.0 / span
    return table


def _emit_layer(layer: dict, name: str, prefix: str, real: str, digits: int, suffix: str) -> str:
    """Emit the constants and forward functions for one distilled KAN layer."""
    order = int(layer["spline_order"])
    knots = _layer_knots(layer["grid"])
    n_knots = int(layer["grid"].shape[1])
    inv = _inverse_spans(knots, order)

    return _LAYER.substitute(
        name=name,
        up=name.upper(),
        up_prefix=prefix.upper(),
        real=real,
        in_f=int(layer["in_features"]),
        out_f=int(layer["out_features"]),
        order=order,
        n_knots=n_knots,
        n_coeff=n_knots - 1 - order,
        knots=_c_array(f"{name}_knots", knots, real, digits, suffix),
        inv_steps=_c_array(f"{name}_inv_span", inv, real, digits, suffix),
        base_w=_c_array(f"{name}_base_w", layer["base_weight"], real, digits, suffix),
        spline_w=_c_array(f"{name}_spline_w", layer["spline_weight"], real, digits, suffix),
    )


def _emit_vectors(extracted, test_windows, prefix, dtype, digits, suffix) -> str:
    """
    Emit golden vectors: input windows beside the reference evaluator's scores.

    Scores are always written in double precision, whatever type the detector is
    generated in, so the same expectations serve a float build and a double one.
    The tolerance is chosen to suit the generated type, leaving each target's
    harness to do nothing but compare.
    """
    from .distill import KANGDNNumpy

    windows = np.asarray(test_windows, dtype=float)
    if windows.ndim != 3:
        raise ValueError(
            f"test_windows must be 3-D (n, window_size, n_features), got {windows.shape}"
        )

    expected = KANGDNNumpy(extracted).decision_function(windows)
    flat = windows.reshape(windows.shape[0], -1)

    rows = []
    for row in flat:
        literals = [_fmt(v, digits, suffix) for v in row]
        chunks = ["        " + ", ".join(literals[i : i + 6]) for i in range(0, len(literals), 6)]
        rows.append("    {\n" + ",\n".join(chunks) + "\n    }")

    return _VECTORS.substitute(
        prefix=prefix,
        up=prefix.upper(),
        provenance=_PROVENANCE,
        guard=f"{prefix.upper()}_VECTORS_H",
        n_vectors=windows.shape[0],
        # Single precision loses relative accuracy on scores that are small
        # compared with the intermediate magnitudes, so an absolute term is
        # needed alongside the relative one.
        rtol="1e-4" if dtype == "float" else "1e-9",
        atol="1e-5" if dtype == "float" else "1e-12",
        windows=",\n".join(rows),
        scores=",\n".join(f"    {v:.17g}" for v in expected),
        flags=",\n".join(f"    {int(v > extracted['threshold'])}" for v in expected),
    )


[docs] def generate_c( extracted: dict, prefix: str = "kangdn", dtype: str = "float", test_windows=None, ) -> dict[str, str]: """ Generate a self-contained C implementation of a distilled KANGDN. Parameters ---------- extracted : dict Output of :func:`~telemetry_anomdet.models.deep.distill.extract_kan_gdn`. prefix : str, default="kangdn" Identifier prefix for the emitted symbols, macros and file names. dtype : {"float", "double"}, default="float" C type for parameters and arithmetic. ``float`` is the deployment mode; ``double`` reproduces the NumPy evaluator to near machine precision. test_windows : np.ndarray, optional Windows of shape ``(n, window_size, n_features)``. When given, a third file of golden vectors is emitted: each window flattened as the entry point expects, beside the reference evaluator's score for it. A port is then validated on the target itself, by scoring each window and comparing, with no host tooling in the loop. Returns ------- dict ``{f"{prefix}.h": ..., f"{prefix}.c": ...}``, source text ready to write, plus ``{prefix}_vectors.h`` when ``test_windows`` is given. Raises ------ NotImplementedError If a KAN layer's knot vector differs between input features, if neighbour counts vary between nodes, or if the model has more nodes than a byte index can address. ValueError If ``dtype`` is unrecognised or a parameter is not finite. """ if dtype not in ("float", "double"): raise ValueError(f"dtype must be 'float' or 'double', got {dtype!r}") real = dtype score_channels = extracted.get("score_channels") if score_channels is None: score_channels = list(range(int(extracted["net"]["n_nodes"]))) score_channels = [int(c) for c in score_channels] if extracted.get("smoothing") is not None: # EWMA carries one value of state per node between windows. The emitted # entry points are pure functions of a single window, so supporting it # requires a state block and a reset call in the interface. raise NotImplementedError( "C generation does not yet support a detector fitted with smoothing; " "the EWMA recursion is stateful across windows." ) # 9 significant digits round-trips float32, 17 round-trips float64. digits = 9 if dtype == "float" else 17 suffix = "f" if dtype == "float" else "" expf = "expf" if dtype == "float" else "exp" fabsf = "fabsf" if dtype == "float" else "fabs" net = extracted["net"] n_nodes = int(net["n_nodes"]) embed_dim = int(net["embed_dim"]) window = int(net["window"]) scaled = extracted["scaler_mean"] is not None if n_nodes > 255: raise NotImplementedError( f"Neighbour indices are emitted as bytes; {n_nodes} nodes exceeds that." ) # The frozen graph becomes explicit neighbour lists, so the softmax runs over # a fixed small count instead of every node. adj = np.asarray(net["adj"], dtype=bool) neighbours = [np.flatnonzero(row) for row in adj] counts = {len(n) for n in neighbours} if len(counts) != 1: raise NotImplementedError( f"Expected a fixed neighbour count per node, got {sorted(counts)}" ) n_nbr = counts.pop() up = prefix.upper() parts = [ _PRELUDE.substitute( prefix=prefix, up=up, provenance=_PROVENANCE, real=real, expf=expf, n_nodes=n_nodes, embed_dim=embed_dim, window=window, n_nbr=n_nbr, n_score=len(score_channels), score_channels=( f"static const unsigned char score_channels[{len(score_channels)}] = {{\n" + ",\n".join( " " + ", ".join(str(c) for c in score_channels[i : i + 12]) for i in range(0, len(score_channels), 12) ) + "\n};\n" ), leaky=_fmt(net["leaky_slope"], digits, suffix), score_limit=_fmt(1.0e9, digits, suffix), ), _emit_layer(net["activation"], "act", prefix, real, digits, suffix), _emit_layer(net["out"], "head", prefix, real, digits, suffix), "\n/* ---- encoder parameters ---- */\n", _c_array("embedding", net["embedding"], real, digits, suffix), _c_array("feat_w", net["feat_weight"], real, digits, suffix), _c_array("feat_b", net["feat_bias"], real, digits, suffix), _c_array("attn_src", net["attn_src"], real, digits, suffix), _c_array("attn_dst", net["attn_dst"], real, digits, suffix), "static const unsigned char neighbours[" + f"{n_nodes} * {n_nbr}] = {{\n" + ",\n".join(" " + ", ".join(str(int(j)) for j in row) for row in neighbours) + "\n};\n", "\n/* ---- scoring parameters ---- */\n", ] if scaled: parts.append(_c_array("scaler_mean", extracted["scaler_mean"], real, digits, suffix)) parts.append(_c_array("scaler_scale", extracted["scaler_scale"], real, digits, suffix)) parts.append(_c_array("err_median", extracted["err_median"], real, digits, suffix)) parts.append(_c_array("err_iqr", extracted["err_iqr"], real, digits, suffix)) if scaled: scale_ctx = f"(window[(t * {up}_N_FEATURES) + f] - scaler_mean[f]) / scaler_scale[f]" scale_tgt = ( f"(window[(({up}_WINDOW_SIZE - 1) * {up}_N_FEATURES) + f]" f" - scaler_mean[f]) / scaler_scale[f]" ) else: scale_ctx = f"window[(t * {up}_N_FEATURES) + f]" scale_tgt = f"window[(({up}_WINDOW_SIZE - 1) * {up}_N_FEATURES) + f]" parts.append( _SCORING.substitute( prefix=prefix, up=up, real=real, expf=expf, fabsf=fabsf, scale_ctx=scale_ctx, scale_tgt=scale_tgt, ) ) header = _HEADER.substitute( prefix=prefix, up=up, provenance=_PROVENANCE, version=__version__, guard=f"{up}_H", real=real, window_size=window + 1, n_features=n_nodes, threshold=_fmt(extracted["threshold"], digits, suffix), ) files = {f"{prefix}.h": header, f"{prefix}.c": "".join(parts)} if test_windows is not None: files[f"{prefix}_vectors.h"] = _emit_vectors( extracted, test_windows, prefix, dtype, digits, suffix ) return files
[docs] def write_c( extracted: dict, out_dir, prefix: str = "kangdn", dtype: str = "float", test_windows=None, ) -> list: """ Generate the C sources and write them into ``out_dir``. Returns the list of written paths. See :func:`generate_c` for the parameters. """ from pathlib import Path out_dir = Path(out_dir) out_dir.mkdir(parents=True, exist_ok=True) written = [] generated = generate_c(extracted, prefix=prefix, dtype=dtype, test_windows=test_windows) for filename, source in generated.items(): path = out_dir / filename path.write_text(source, encoding="utf-8") written.append(path) return written