Affected Repositories

• PyAutoArray (primary)

Branch Survey

No active.md claim on PyAutoArray (worktree_check_conflict exit 0). Note: an unregistered sparse-operator-oversampling-cache PyAutoArray worktree has 11 uncommitted files and may touch factory.py / sparse code — check before phase 2.

Suggested branch: feature/delaunay-nn-numba · worktree ~/Code/PyAutoLabs-wt/delaunay-nn-numba/

Implementation Steps

Phase 1 — autoarray/inversion/mesh/interpolator/sibson.py

1. New _sibson_from_tables_numba(queries, seeds, fallbacks, points, simplices_padded, simplex_neighbors, centres, r2, valid, C, M) built lazily with numba.njit(cache=True, parallel=True, nogil=True) + prange over queries (pattern: interferometer_numba/inversion_interferometer_numba_util.py:109-135). Per query, fixed-size scratch: BFS cavity walk over simplex_neighbors with the _contains_query tolerance 64*eps*max(r²,1) (cap C = SIBSON_MAX_CAVITY_TRIANGLES = 32, overflow flag); Watson circumcircle contributions (≤3C); group/compact by vertex id via insertion sort (replacing argsort + np.add.at); normalise.
2. Keep every edge case of _sibson_single_from_tables (L343-556): outside-mesh fallback (seed<0 → [outside_fallback_index], w=1, size 1), on-vertex fallback, on-edge barycentric row (128*eps, size 3), degenerate (invalid inserted circle / non-finite or |sum|≤tiny), neighbour overflow (unique_count > M, M = SIBSON_MAX_NEIGHBORS = 32) → NaN weights; all six outputs (mappings (Q,M) int32, sizes, weights (Q,M), cavity_sizes, overflow, degenerate) unchanged.
3. NumPy branch of sibson_mappings_weights_from_tables (L629-636) dispatches to the kernel when numba is importable and general.yaml numba.use_numba is true; else keeps the list comprehension, retained as the reference. JAX branch untouched.
4. Tests: new test_autoarray/inversion/pixelization/interpolator/test_sibson_numba.py — kernel vs reference on random meshes and constructed on-vertex / on-edge / outside-hull / degenerate / cavity-overflow cases; mappings, sizes and diagnostics identical, weights to ~1e-12. Existing test_sibson.py (incl. scipy_delaunay_nn fixed-seed regression L231), test_delaunay_nn.py, test_numerics_audit.py stay green.
5. Re-time with the Aris scripts (local scratchpad, collaborator data — never committed); numbers in the PR body.

Phase 2 — routing (autoarray/inversion/inversion/factory.py)
6. Interferometer: measure nnz_per_source_column_from(mapper) for DelaunayNN on the Aris mask (1405, 5526 vertices). If > 60, time direct_conv vs FFT for DelaunayNN to locate the real crossover; retune interferometer_numba_nnz_per_source_max (settings.py:377-412, general.yaml:21) or make the gate mesh-aware based on data. kernel_inputs_from already handles variable widths (CSR). Add a NaN-weight check in _use_interferometer_numba (L242-321) or the sparse.py:97-145 preconditions.
7. Imaging: confirm InversionImagingSparseNumba (imaging_numba/sparse.py:55, chosen at factory.py:89-163) is taken with a DelaunayNN mapper and that mapper_numba_util.data_slim_to_pixelization_unique_from handles widths up to 32; time Delaunay vs DelaunayNN on HST-scale imaging.
8. Tests: DelaunayNN cases in interferometer_numba/test_interferometer_numba.py routing tests (L309-459) and imaging in test_factory.py; likelihood parity numba vs FFT/mapping route.
9. Report table: inversion class per mesh × dataset, before/after, with Delaunay vs DelaunayNN timings.

Key Files

• autoarray/inversion/mesh/interpolator/sibson.py — Sibson tables, NumPy/JAX branches, scipy_delaunay_nn
• autoarray/numba_util.py, autoarray/config/general.yaml — numba conventions/config
• autoarray/inversion/inversion/factory.py, autoarray/settings.py — inversion routing + nnz gate
• autoarray/inversion/inversion/interferometer_numba/ — direct_conv kernel + kernel_inputs_from
• autoarray/inversion/inversion/imaging_numba/sparse.py — imaging numba inversion

Risks

• Kernel summation order differs from argsort grouping → weights compared to 1e-12, not bit-for-bit.
• numba compile time on first call (cached).

Make DelaunayNN fast on the numba/NumPy path

Type: feature
Target: PyAutoArray
Repos:

• PyAutoArray
  Difficulty: large
  Autonomy: supervised
  Priority: normal
  Memory: wiki/galaxies/sources/massive-ellipticals.md; wiki/lensing/sources/dark-matter-substructure.md; reading-queue.md
  Status: formalised
  Consequence: glance
  Witness: on the Aris interferometer setup (1405 vertices, use_jax=False) DelaunayNN steady log_likelihood_function time <= 2x Delaunay, with numba-vs-reference Sibson parity tests green
  Review-minutes: 3
  Unattended: ready

DelaunayNN is ~40x slower than Delaunay on the NumPy (use_jax=False) path because the
Sibson natural-neighbour interpolation is a pure-Python per-query loop. Make it fast with
a numba kernel so DelaunayNN costs about the same as Delaunay on CPU.

Evidence (measured 2026-10-15, laptop, PyAutoArray ccddfba)

Collaborator (Aris) SMA-scale interferometer data: 5,516,930 visibilities, 5,496 masked
pixels, sparse operator applied, AdaptSplit, use_positive_only_solver=False, threads
pinned to 1. Bare analysis.log_likelihood_function steady seconds per call:

Profile: almost all DelaunayNN time is in sibson_mappings_weights_from_tables ->
single() (autoarray/inversion/mesh/interpolator/sibson.py ~L613-644), reached via
scipy_delaunay_nn -> mappings_weights_for (~11k calls per likelihood: data points +
split points). The NumPy branch is a Python list comprehension over queries; each
_sibson_single_from_tables call runs ~50-100 tiny NumPy ops (argsort, np.add.at,
broadcast_to, where, stack, pad...) for at most ~32 cavity triangles / ~96 contributions,
~1.5 ms per query. The JAX branch is vectorised (vmap + lax.map) and is not affected.

Scripts/logs: scratchpad .../50c8c6cd-6ac2-47c2-9290-7ee73f9c4315/scratchpad/aris/step6
and PROGRESS2.md (temporary; numbers are reproduced above).

Fix

• Add an @njit (parallel=True / prange over queries) kernel for the NumPy path of
  sibson_mappings_weights_from_tables: cavity walk, Watson circumcircle contributions,
  compaction/dedup of vertex ids into fixed-size max_neighbors rows, normalisation.
• Reproduce every edge case of _sibson_single_from_tables exactly: outside-mesh
  fallback (seed_simplex < 0), on-vertex fallback, on-edge barycentric row, degenerate
  Watson circles and cavity/neighbour overflow -> NaN weights, plus the returned
  diagnostics (cavity_sizes, overflow, degenerate).
• Follow the existing numba conventions in interpolator/delaunay.py /
  rectangular.py (config numba: block, caching). JAX path untouched.

Numba likelihood paths must use it (imaging + interferometer)

The faster interpolation alone is not enough: the numba likelihood paths must actually
take DelaunayNN mappers.

• Interferometer: factory.py routes a NumPy inversion to the numba direct_conv
  curvature path (InversionInterferometerSparseNumba) only when mean non-zeros per
  source column <= interferometer_numba_nnz_per_source_max (60, tuned on Delaunay).
  DelaunayNN stencils carry more neighbours per data point (up to
  SIBSON_MAX_NEIGHBORS = 32 vs 3), so DelaunayNN may fall through to the slower FFT
  route. Measure nnz/column for DelaunayNN on the Aris mask; confirm which route it takes;
  if it misses direct_conv, make it take it (or re-measure the crossover) — and check the
  kernel handles variable-size Sibson rows correctly.
• Imaging: do the same audit for the numba imaging likelihood path — confirm DelaunayNN
  mappers are routed to it and handled correctly, and time imaging Delaunay vs DelaunayNN.
• Report, for both datasets, which inversion class each mesh is routed to before/after.

Acceptance

• Parity unit tests: numba kernel vs the current NumPy reference on random meshes and on
  constructed edge cases (query on vertex, on edge, outside hull, degenerate) — mappings
  and sizes identical, weights to ~1e-12.
• Re-time on the Aris setup: Sibson step from ~18 s to tens of ms; DelaunayNN per-call
  cost close to Delaunay (report both vertex counts).
• Imaging and interferometer numba paths are taken by DelaunayNN mappers (asserted in a test
  of the routing), with likelihood parity vs the FFT/reference route.
• Existing DelaunayNN tests and the workspace_test parity scripts still pass.