|
| 1 | +// src/translated_functions.c |
| 2 | +// Build with: python setup.py build_ext --inplace |
| 3 | +// Exposes: |
| 4 | +// generate_increasing_tuples(m, k) -> ndarray[intp] |
| 5 | +// custom_bit_length(x) -> int |
| 6 | +// x_is_a_superset_of_any_in_list(x, L) -> bool |
| 7 | +// lattice_slice2(m, k, minimal_non_faces) -> ndarray[uint64] |
| 8 | + |
| 9 | +#define PY_SSIZE_T_CLEAN |
| 10 | +#include <Python.h> |
| 11 | +#define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION |
| 12 | +#include <numpy/arrayobject.h> |
| 13 | +#include <stdint.h> |
| 14 | + |
| 15 | +/*---------------------------------------------------------------------- |
| 16 | + * Helpers |
| 17 | + *----------------------------------------------------------------------*/ |
| 18 | +static int64_t binom_int64(int64_t n, int64_t k) { |
| 19 | + if (k < 0 || k > n) return 0; |
| 20 | + if (k > n - k) k = n - k; |
| 21 | + int64_t c = 1; |
| 22 | + for (int64_t i = 1; i <= k; ++i) |
| 23 | + c = (c * (n - (k - i))) / i; |
| 24 | + return c; |
| 25 | +} |
| 26 | + |
| 27 | +/* Return 1 if (b & ~x) == 0 for any b in seq, 0 otherwise, -1 on error */ |
| 28 | +static int is_superset(unsigned long long x, PyObject *seq_fast) { |
| 29 | + Py_ssize_t n = PySequence_Fast_GET_SIZE(seq_fast); |
| 30 | + PyObject **items = PySequence_Fast_ITEMS(seq_fast); |
| 31 | + for (Py_ssize_t i = 0; i < n; ++i) { |
| 32 | + unsigned long long b = PyLong_AsUnsignedLongLong(items[i]); |
| 33 | + if (PyErr_Occurred()) return -1; |
| 34 | + if ((b & ~x) == 0ULL) return 1; |
| 35 | + } |
| 36 | + return 0; |
| 37 | +} |
| 38 | + |
| 39 | +/* -------------------------------------------------------------------------- */ |
| 40 | +/* intersection_of_codewords_from_bits(x, a) */ |
| 41 | +/* -------------------------------------------------------------------------- */ |
| 42 | +static uint64_t intersect_mask_with_words(unsigned long long mask, |
| 43 | + const uint64_t *words, |
| 44 | + npy_uintp len, |
| 45 | + int *err) |
| 46 | +{ |
| 47 | + if (mask == 0ULL) { |
| 48 | + *err = 1; return 0ULL; |
| 49 | + } |
| 50 | + uint64_t res = UINT64_MAX; int init = 0; npy_uintp idx = 0; |
| 51 | + unsigned long long tmp = mask; |
| 52 | + while (tmp) { |
| 53 | + if (tmp & 1ULL) { |
| 54 | + if (idx >= len) { *err = 1; return 0ULL; } |
| 55 | + res = init ? (res & words[idx]) : words[idx]; |
| 56 | + init = 1; |
| 57 | + if (res == 0ULL) return 0ULL; /* early exit */ |
| 58 | + } |
| 59 | + tmp >>= 1ULL; ++idx; |
| 60 | + } |
| 61 | + *err = 0; return res; |
| 62 | +} |
| 63 | + |
| 64 | + |
| 65 | + |
| 66 | +/*---------------------------------------------------------------------- |
| 67 | + * generate_increasing_tuples |
| 68 | + *----------------------------------------------------------------------*/ |
| 69 | +static PyObject *py_generate_increasing_tuples(PyObject *self, PyObject *args) { |
| 70 | + long long m_in, k_in; |
| 71 | + if (!PyArg_ParseTuple(args, "LL", &m_in, &k_in)) return NULL; |
| 72 | + if (m_in < 0 || k_in < 0 || k_in > m_in) { |
| 73 | + PyErr_SetString(PyExc_ValueError, "Require 0 <= k <= m and m,k >= 0"); |
| 74 | + return NULL; |
| 75 | + } |
| 76 | + int64_t m = (int64_t)m_in, k = (int64_t)k_in; |
| 77 | + int64_t ncomb = binom_int64(m, k); |
| 78 | + npy_intp dims[2] = {(npy_intp)ncomb, (npy_intp)k}; |
| 79 | + PyObject *arr = PyArray_SimpleNew(2, dims, NPY_INTP); |
| 80 | + if (!arr) return NULL; |
| 81 | + npy_intp *data = (npy_intp *)PyArray_DATA((PyArrayObject *)arr); |
| 82 | + |
| 83 | + if (k > 64) { |
| 84 | + Py_DECREF(arr); |
| 85 | + PyErr_SetString(PyExc_ValueError, "k too large; adjust buffer size"); |
| 86 | + return NULL; |
| 87 | + } |
| 88 | + npy_intp current[64]; |
| 89 | + for (npy_intp i = 0; i < k; ++i) current[i] = i; |
| 90 | + npy_intp row = 0; |
| 91 | + int done = 0; |
| 92 | + while (!done) { |
| 93 | + for (npy_intp j = 0; j < k; ++j) data[row * k + j] = current[j]; |
| 94 | + ++row; |
| 95 | + npy_intp idx = k - 1; |
| 96 | + while (idx >= 0 && current[idx] == m - k + idx) --idx; |
| 97 | + if (idx < 0) done = 1; |
| 98 | + else { |
| 99 | + ++current[idx]; |
| 100 | + for (npy_intp j = idx + 1; j < k; ++j) current[j] = current[j - 1] + 1; |
| 101 | + } |
| 102 | + } |
| 103 | + return arr; |
| 104 | +} |
| 105 | + |
| 106 | +/*---------------------------------------------------------------------- |
| 107 | + * custom_bit_length |
| 108 | + *----------------------------------------------------------------------*/ |
| 109 | +static PyObject *py_custom_bit_length(PyObject *self, PyObject *args) { |
| 110 | + unsigned long long x; |
| 111 | + if (!PyArg_ParseTuple(args, "K", &x)) return NULL; |
| 112 | + if (x == 0) return PyLong_FromLong(0); |
| 113 | + int len = 0; |
| 114 | + while (x) { ++len; x >>= 1; } |
| 115 | + return PyLong_FromLong(len); |
| 116 | +} |
| 117 | + |
| 118 | +/*---------------------------------------------------------------------- |
| 119 | + * x_is_a_superset_of_any_in_list |
| 120 | + *----------------------------------------------------------------------*/ |
| 121 | +static PyObject *py_x_is_superset_of_any(PyObject *self, PyObject *args) { |
| 122 | + PyObject *x_obj, *seq_obj; |
| 123 | + if (!PyArg_ParseTuple(args, "OO", &x_obj, &seq_obj)) return NULL; |
| 124 | + unsigned long long x = PyLong_AsUnsignedLongLong(x_obj); |
| 125 | + if (PyErr_Occurred()) return NULL; |
| 126 | + PyObject *seq_fast = PySequence_Fast(seq_obj, "L must be a sequence"); |
| 127 | + if (!seq_fast) return NULL; |
| 128 | + int res = is_superset(x, seq_fast); |
| 129 | + Py_DECREF(seq_fast); |
| 130 | + if (res < 0) return NULL; |
| 131 | + return PyBool_FromLong(res); |
| 132 | +} |
| 133 | + |
| 134 | +/*---------------------------------------------------------------------- |
| 135 | + * lattice_slice2 |
| 136 | + *----------------------------------------------------------------------*/ |
| 137 | +static PyObject *py_lattice_slice2(PyObject *self, PyObject *args) { |
| 138 | + long long m_in, k_in; PyObject *minfaces_obj; |
| 139 | + if (!PyArg_ParseTuple(args, "LLO", &m_in, &k_in, &minfaces_obj)) return NULL; |
| 140 | + long long m = m_in, k = k_in; |
| 141 | + if (k <= 1) { PyErr_Format(PyExc_ValueError, "The size k=%lld is not larger than 1", k); return NULL; } |
| 142 | + if (k > m) { PyErr_Format(PyExc_ValueError, "The size k=%lld is larger than m=%lld", k, m); return NULL; } |
| 143 | + if (m > 64) { PyErr_SetString(PyExc_ValueError, "m > 64 not supported"); return NULL; } |
| 144 | + |
| 145 | + PyObject *minfaces_fast = PySequence_Fast(minfaces_obj, "minimal_non_faces must be a sequence"); |
| 146 | + if (!minfaces_fast) return NULL; |
| 147 | + |
| 148 | + /* collect uint64 results first */ |
| 149 | + PyObject *tmp_list = PyList_New(0); |
| 150 | + if (!tmp_list) { Py_DECREF(minfaces_fast); return NULL; } |
| 151 | + |
| 152 | + unsigned char current[64]; |
| 153 | + for (long long i = 0; i < k; ++i) current[i] = (unsigned char)i; |
| 154 | + int done = 0; |
| 155 | + while (!done) { |
| 156 | + unsigned long long x = 0ULL; |
| 157 | + for (long long j = 0; j < k; ++j) x |= (1ULL << current[j]); |
| 158 | + int super = 0; |
| 159 | + if (PySequence_Fast_GET_SIZE(minfaces_fast) > 0) { |
| 160 | + super = is_superset(x, minfaces_fast); |
| 161 | + if (super < 0) { Py_DECREF(minfaces_fast); Py_DECREF(tmp_list); return NULL; } |
| 162 | + } |
| 163 | + if (!super) { |
| 164 | + PyObject *x_py = PyLong_FromUnsignedLongLong(x); |
| 165 | + if (!x_py || PyList_Append(tmp_list, x_py) < 0) { Py_XDECREF(x_py); Py_DECREF(minfaces_fast); Py_DECREF(tmp_list); return NULL; } |
| 166 | + Py_DECREF(x_py); |
| 167 | + } |
| 168 | + long long idx = k - 1; |
| 169 | + while (idx >= 0 && current[idx] == m - k + idx) --idx; |
| 170 | + if (idx < 0) done = 1; |
| 171 | + else { ++current[idx]; for (long long j = idx + 1; j < k; ++j) current[j] = current[j - 1] + 1; } |
| 172 | + } |
| 173 | + |
| 174 | + /* convert tmp_list -> ndarray[uint64] */ |
| 175 | + Py_ssize_t n_out = PyList_GET_SIZE(tmp_list); |
| 176 | + npy_intp dims[1] = {(npy_intp)n_out}; |
| 177 | + PyObject *arr = PyArray_SimpleNew(1, dims, NPY_UINT64); |
| 178 | + if (!arr) { Py_DECREF(minfaces_fast); Py_DECREF(tmp_list); return NULL; } |
| 179 | + uint64_t *arr_data = (uint64_t *)PyArray_DATA((PyArrayObject *)arr); |
| 180 | + for (Py_ssize_t i = 0; i < n_out; ++i) { |
| 181 | + PyObject *item = PyList_GET_ITEM(tmp_list, i); |
| 182 | + arr_data[i] = (uint64_t)PyLong_AsUnsignedLongLong(item); |
| 183 | + if (PyErr_Occurred()) { Py_DECREF(minfaces_fast); Py_DECREF(tmp_list); Py_DECREF(arr); return NULL; } |
| 184 | + } |
| 185 | + |
| 186 | + Py_DECREF(minfaces_fast); |
| 187 | + Py_DECREF(tmp_list); |
| 188 | + return arr; |
| 189 | +} |
| 190 | + |
| 191 | + |
| 192 | +/* ------------------------------------------------------------------------- */ |
| 193 | +/* intersection_of_codewords_from_bits(x, a) */ |
| 194 | +/* ------------------------------------------------------------------------- */ |
| 195 | +static PyObject *py_intersection_of_codewords(PyObject *self, PyObject *args) |
| 196 | +{ |
| 197 | + unsigned long long x; PyObject *arr_obj; |
| 198 | + if (!PyArg_ParseTuple(args, "KO", &x, &arr_obj)) return NULL; |
| 199 | + if (x == 0ULL) { |
| 200 | + PyErr_SetString(PyExc_ValueError, "x must be non-zero bitmask"); |
| 201 | + return NULL; |
| 202 | + } |
| 203 | + PyArrayObject *arr = (PyArrayObject *)PyArray_FROM_OTF(arr_obj, NPY_UINT64, NPY_ARRAY_IN_ARRAY); |
| 204 | + if (!arr) return NULL; |
| 205 | + if (PyArray_NDIM(arr) != 1) { |
| 206 | + PyErr_SetString(PyExc_ValueError, "a must be 1-D uint64 array"); |
| 207 | + Py_DECREF(arr); return NULL; |
| 208 | + } |
| 209 | + npy_uintp len = PyArray_SIZE(arr); |
| 210 | + const uint64_t *data = (const uint64_t *)PyArray_DATA(arr); |
| 211 | + |
| 212 | + uint64_t res = UINT64_MAX; int init = 0; npy_uintp idx = 0; |
| 213 | + unsigned long long tmp = x; |
| 214 | + while (tmp) { |
| 215 | + if (tmp & 1ULL) { |
| 216 | + if (idx >= len) { |
| 217 | + PyErr_SetString(PyExc_ValueError, "bit in x exceeds array length"); |
| 218 | + Py_DECREF(arr); return NULL; |
| 219 | + } |
| 220 | + res = init ? (res & data[idx]) : data[idx]; |
| 221 | + init = 1; |
| 222 | + if (res == 0ULL) { Py_DECREF(arr); return PyLong_FromUnsignedLongLong(0ULL); } |
| 223 | + } |
| 224 | + tmp >>= 1ULL; ++idx; |
| 225 | + } |
| 226 | + Py_DECREF(arr); |
| 227 | + return PyLong_FromUnsignedLongLong(res); |
| 228 | +} |
| 229 | + |
| 230 | + |
| 231 | + |
| 232 | +/*---------------------------------------------------------------------- |
| 233 | + * Module definition |
| 234 | + *----------------------------------------------------------------------*/ |
| 235 | +static PyMethodDef Methods[] = { |
| 236 | + {"generate_increasing_tuples", py_generate_increasing_tuples, METH_VARARGS, "Generate all increasing k‑tuples"}, |
| 237 | + {"custom_bit_length", py_custom_bit_length, METH_VARARGS, "Bit length of unsigned int"}, |
| 238 | + {"x_is_a_superset_of_any_in_list", py_x_is_superset_of_any, METH_VARARGS, "Return True if x supersets any b in list"}, |
| 239 | + {"lattice_slice2", py_lattice_slice2, METH_VARARGS, "Return uint64 ndarray of valid k-sets"}, |
| 240 | + {"intersection_of_codewords_from_bits", py_intersection_of_codewords, METH_VARARGS, "AND‑intersect selected codewords"}, |
| 241 | + {NULL, NULL, 0, NULL} |
| 242 | +}; |
| 243 | + |
| 244 | + |
| 245 | + |
| 246 | +static struct PyModuleDef moduledef = { |
| 247 | + PyModuleDef_HEAD_INIT, |
| 248 | + "translated_functions", |
| 249 | + "C translations of performance‑critical Numba functions", |
| 250 | + -1, |
| 251 | + Methods |
| 252 | +}; |
| 253 | + |
| 254 | +PyMODINIT_FUNC PyInit_translated_functions(void) { |
| 255 | + PyObject *m = PyModule_Create(&moduledef); |
| 256 | + if (!m) return NULL; |
| 257 | + import_array(); |
| 258 | + return m; |
| 259 | +} |
0 commit comments