/*
=============================================================================================
SCE_CPYTHON_EXTENSION: SOVEREIGN REASONING & INT4 MoE CPYTHON NATIVE EXTENSION
=============================================================================================
100% Native CPython C-API implementation (Python.h):
1. INT4 Symmetric Group Quantization & Dequantization (Bitwise SIMD-speed unpack)
2. Symplectic Router Dynamics Integration (2nd-Order Critically Damped ODE)
3. LaSalle-Lyapunov Stability Manifold Evaluation (V(x) = x^T P x)
4. Fast POSIX Patch Sanitizer (CRLF -> LF and EOF Trailing Newline Enforcement)
5. Bradley-Terry (BT) / ELO Pairwise Judging Match Matrix
=============================================================================================
*/

#define PY_SSIZE_T_CLEAN
#include <Python.h>
#include <stdint.h>
#include <math.h>
#include <string.h>

/* --- 1. INT4 QUANTIZE (Float Array -> Packed Bytes & Scales) --- */
static PyObject* py_sce_int4_quantize(PyObject* self, PyObject* args) {
    PyObject* float_list_obj;
    int out_features, in_features, group_size;

    if (!PyArg_ParseTuple(args, "Oiii", &float_list_obj, &out_features, &in_features, &group_size)) {
        return NULL;
    }

    if (!PyList_Check(float_list_obj)) {
        PyErr_SetString(PyExc_TypeError, "Expected list of floats");
        return NULL;
    }

    Py_ssize_t total_elements = (Py_ssize_t)out_features * in_features;
    if (PyList_Size(float_list_obj) != total_elements) {
        PyErr_SetString(PyExc_ValueError, "List size does not match out_features * in_features");
        return NULL;
    }

    int groups_per_row = in_features / group_size;
    int packed_row_bytes = in_features / 2;

    PyObject* packed_bytes = PyBytes_FromStringAndSize(NULL, out_features * packed_row_bytes);
    if (!packed_bytes) return NULL;
    uint8_t* packed_ptr = (uint8_t*)PyBytes_AsString(packed_bytes);

    PyObject* scales_list = PyList_New(out_features * groups_per_row);
    if (!scales_list) {
        Py_DECREF(packed_bytes);
        return NULL;
    }

    for (int r = 0; r < out_features; ++r) {
        for (int g = 0; g < groups_per_row; ++g) {
            int start_idx = r * in_features + g * group_size;
            float max_val = 1e-8f;

            for (int k = 0; k < group_size; ++k) {
                PyObject* item = PyList_GET_ITEM(float_list_obj, start_idx + k);
                float val = (float)fabs(PyFloat_AsDouble(item));
                if (val > max_val) max_val = val;
            }

            float scale = max_val / 7.0f;
            PyList_SET_ITEM(scales_list, r * groups_per_row + g, PyFloat_FromDouble((double)scale));

            for (int k = 0; k < group_size; k += 2) {
                float v_even = (float)PyFloat_AsDouble(PyList_GET_ITEM(float_list_obj, start_idx + k));
                float v_odd  = (float)PyFloat_AsDouble(PyList_GET_ITEM(float_list_obj, start_idx + k + 1));

                int q_even = (int)roundf(v_even / scale);
                int q_odd  = (int)roundf(v_odd / scale);

                if (q_even < -8) q_even = -8; if (q_even > 7) q_even = 7;
                if (q_odd  < -8) q_odd  = -8; if (q_odd  > 7) q_odd  = 7;

                uint8_t u_even = (uint8_t)(q_even + 8);
                uint8_t u_odd  = (uint8_t)(q_odd + 8);

                int byte_idx = r * packed_row_bytes + (g * group_size + k) / 2;
                packed_ptr[byte_idx] = (u_odd << 4) | (u_even & 0x0F);
            }
        }
    }

    return Py_BuildValue("(NN)", packed_bytes, scales_list);
}

/* --- 2. INT4 DEQUANTIZE (Packed Bytes + Scales -> Float List) --- */
static PyObject* py_sce_int4_dequantize(PyObject* self, PyObject* args) {
    Py_buffer view;
    PyObject* scales_obj;
    int out_features, in_features, group_size;

    if (!PyArg_ParseTuple(args, "y*Oiii", &view, &scales_obj, &out_features, &in_features, &group_size)) {
        return NULL;
    }

    if (!PyList_Check(scales_obj)) {
        PyBuffer_Release(&view);
        PyErr_SetString(PyExc_TypeError, "Scales must be a list of floats");
        return NULL;
    }

    Py_ssize_t total_elements = (Py_ssize_t)out_features * in_features;
    PyObject* result_list = PyList_New(total_elements);
    if (!result_list) {
        PyBuffer_Release(&view);
        return NULL;
    }

    int groups_per_row = in_features / group_size;
    int packed_row_bytes = in_features / 2;
    const uint8_t* packed_ptr = (const uint8_t*)view.buf;

    for (int r = 0; r < out_features; ++r) {
        for (int g = 0; g < groups_per_row; ++g) {
            float scale = (float)PyFloat_AsDouble(PyList_GET_ITEM(scales_obj, r * groups_per_row + g));
            int group_start = g * group_size;

            for (int k = 0; k < group_size; k += 2) {
                int byte_idx = r * packed_row_bytes + (group_start + k) / 2;
                uint8_t b = packed_ptr[byte_idx];

                int8_t even_val = (int8_t)((b & 0x0F) - 8);
                int8_t odd_val  = (int8_t)(((b >> 4) & 0x0F) - 8);

                Py_ssize_t idx_even = (Py_ssize_t)r * in_features + group_start + k;
                Py_ssize_t idx_odd  = idx_even + 1;

                PyList_SET_ITEM(result_list, idx_even, PyFloat_FromDouble((double)(even_val * scale)));
                PyList_SET_ITEM(result_list, idx_odd,  PyFloat_FromDouble((double)(odd_val * scale)));
            }
        }
    }

    PyBuffer_Release(&view);
    return result_list;
}

/* --- 3. SYMPLECTIC ROUTER CRITICALLY DAMPED STEP --- */
static PyObject* py_sce_symplectic_step(PyObject* self, PyObject* args) {
    PyObject *z_list, *z_dot_list, *u_list;
    double omega, dt;

    if (!PyArg_ParseTuple(args, "OOOdd", &z_list, &z_dot_list, &u_list, &omega, &dt)) {
        return NULL;
    }

    Py_ssize_t n = PyList_Size(z_list);
    double omega_sq = omega * omega;
    double two_omega = 2.0 * omega;

    for (Py_ssize_t i = 0; i < n; ++i) {
        double z_val     = PyFloat_AsDouble(PyList_GET_ITEM(z_list, i));
        double z_dot_val = PyFloat_AsDouble(PyList_GET_ITEM(z_dot_list, i));
        double u_val     = PyFloat_AsDouble(PyList_GET_ITEM(u_list, i));

        double acc = omega_sq * (u_val - z_val) - two_omega * z_dot_val;
        z_dot_val += acc * dt;
        z_val     += z_dot_val * dt;

        PyList_SET_ITEM(z_dot_list, i, PyFloat_FromDouble(z_dot_val));
        PyList_SET_ITEM(z_list, i, PyFloat_FromDouble(z_val));
    }

    Py_RETURN_NONE;
}

/* --- 4. LASALLE-LYAPUNOV STABILITY VALUE EVALUATION --- */
static PyObject* py_sce_lyapunov_eval(PyObject* self, PyObject* args) {
    PyObject *x_list, *p_matrix_list;
    int dim;

    if (!PyArg_ParseTuple(args, "OOi", &x_list, &p_matrix_list, &dim)) {
        return NULL;
    }

    double v = 0.0;
    for (int i = 0; i < dim; ++i) {
        double row_sum = 0.0;
        for (int j = 0; j < dim; ++j) {
            double p_val = PyFloat_AsDouble(PyList_GET_ITEM(p_matrix_list, i * dim + j));
            double x_j   = PyFloat_AsDouble(PyList_GET_ITEM(x_list, j));
            row_sum += p_val * x_j;
        }
        double x_i = PyFloat_AsDouble(PyList_GET_ITEM(x_list, i));
        v += x_i * row_sum;
    }

    return PyFloat_FromDouble(v);
}

/* --- 5. FAST POSIX PATCH SANITIZER (Eliminates CRLF & guarantees EOF \n) --- */
static PyObject* py_sce_sanitize_patch(PyObject* self, PyObject* args) {
    const char* raw_patch;
    Py_ssize_t len;

    if (!PyArg_ParseTuple(args, "s#", &raw_patch, &len)) {
        return NULL;
    }

    if (len == 0) {
        return PyUnicode_FromString("");
    }

    // Allocate buffer for cleaned patch (max size: len + 2 for trailing \n)
    char* buf = (char*)malloc(len + 4);
    if (!buf) {
        return PyErr_NoMemory();
    }

    Py_ssize_t out_idx = 0;
    for (Py_ssize_t i = 0; i < len; ++i) {
        if (raw_patch[i] == '\r') {
            // Strip Carriage Return
            continue;
        }
        buf[out_idx++] = raw_patch[i];
    }

    // Enforce POSIX Trailing Newline if non-empty
    if (out_idx > 0 && buf[out_idx - 1] != '\n') {
        buf[out_idx++] = '\n';
    }
    buf[out_idx] = '\0';

    PyObject* res = PyUnicode_FromStringAndSize(buf, out_idx);
    free(buf);
    return res;
}

/* --- 6. BRADLEY-TERRY / ELO PAIRWISE UPDATE --- */
static PyObject* py_sce_bradley_terry_step(PyObject* self, PyObject* args) {
    double rating_a, rating_b;
    int a_won;
    double k_factor;

    if (!PyArg_ParseTuple(args, "ddid", &rating_a, &rating_b, &a_won, &k_factor)) {
        return NULL;
    }

    double exp_a = 1.0 / (1.0 + pow(10.0, (rating_b - rating_a) / 400.0));
    double exp_b = 1.0 - exp_a;

    double actual_a = a_won ? 1.0 : 0.0;
    double actual_b = a_won ? 0.0 : 1.0;

    double new_a = rating_a + k_factor * (actual_a - exp_a);
    double new_b = rating_b + k_factor * (actual_b - exp_b);

    return Py_BuildValue("(dd)", new_a, new_b);
}

/* Method Table */
static PyMethodDef SceCPythonMethods[] = {
    {"int4_quantize", py_sce_int4_quantize, METH_VARARGS, "Quantizes float list to INT4 bytes & scales."},
    {"int4_dequantize", py_sce_int4_dequantize, METH_VARARGS, "Dequantizes INT4 bytes back to float list."},
    {"symplectic_step", py_sce_symplectic_step, METH_VARARGS, "Critically damped symplectic router step."},
    {"lyapunov_eval", py_sce_lyapunov_eval, METH_VARARGS, "Calculates LaSalle-Lyapunov V(x)."},
    {"sanitize_patch", py_sce_sanitize_patch, METH_VARARGS, "Sanitizes patch to pure POSIX LF with guaranteed trailing newline."},
    {"bradley_terry_step", py_sce_bradley_terry_step, METH_VARARGS, "Updates ELO ratings via Bradley-Terry pairwise match."},
    {NULL, NULL, 0, NULL}
};

/* Module Definition */
static struct PyModuleDef sce_cpython_module = {
    PyModuleDef_HEAD_INIT,
    "sce_cpython",
    "Sovereign Cognitive Engine CPython Native Extension",
    -1,
    SceCPythonMethods
};

/* Module Initialization */
PyMODINIT_FUNC PyInit_sce_cpython(void) {
    return PyModule_Create(&sce_cpython_module);
}