Skip to content

Changelog

Notable changes in fastcore-rs. The version number is shared across all three surfaces — the Rust crate fastcore, Python's navis-fastcore and R's nat.fastcore — so a release often touches only some of them; where that matters it is called out.

Tags, source archives and the original announcements are on GitHub.

0.13.0 (2026-08-13)

A mesh can be stripped of the surface it keeps on its inside. An invagination is a piece of membrane that bulges into the cell, usually the boundary of a mitochondrion or a vesicle touching it from within. Segmented meshes are full of them and they are ruinous for anything that walks the surface: each one is a tunnel a skeletonisation's wave front can take a shortcut through, or split on. On the 3.5 M-face mesh this was built against they account for 3,076 of its handles, and removing them takes the skeleton from 3,227 leafs to ~500 without losing a real branch.

The property that defines one — the cell encloses it — is also what makes it invisible from outside, so drop_internals asks that directly:

vertices, faces, keep, passes = fastcore.drop_internals(vertices, faces)

One pass scores every face by the fraction of a ray spray that escapes, diffuses that field over the faces, cuts below a threshold, drops the components left enclosing no volume — the inside-out shreds of pocket wall — and caps the holes. Then it repeats, because capping a pocket mouth turns a partially open neighbour into a fully buried one; the passes converge fast, burying 18.7% of faces, then 0.5%, then 0.2%. Only vertices are ever removed, since caps re-use the ones already on the boundary, so keep is all a caller needs to carry a vertex map or any other per-vertex annotation across the repair. Note that this also closes openings the mesh arrived with — a neurite truncated at the edge of the dataset, say — which for skeletonisation is what you want: a stump left open splits the wave front just as a pocket mouth does.

openness is that score on its own, for callers who want the field rather than the repair — to colour a mesh by it, or to pick a threshold. Outer membrane lands at 0.5-1.0 and the wall of an invagination at exactly 0, with the ground between them thinly populated, which is why threshold is not really a tuning parameter: anything in 0.05-0.10 does the same thing, and what says so is the boundary-edge count, flat across that range and then exploding as the cut outruns the capping.

The ray casting is ours rather than a collision library's, because the question is much smaller than the one those answer: not what a ray hit, nor where, nor which hit came first, only whether it got out. The BVH is a binned-SAH tree whose only query returns on the first triangle it touches — no nearest-hit ordering, no shrinking far bound, no hit record. Against the ncollpyde prototype on that mesh, 118.8s -> 3.4s, for a repaired mesh differing by 527 faces in 2.85 M and a skeleton well inside the noise floor (531 leafs against 523).

Two things the prototype could not promise. The spray is drawn by hashing (face, ray) rather than from a stream, so the answer does not depend on how the faces were split across threads — which is what lets a pass re-cast only near the previous pass's caps and still provably agree with a full sweep. And faces must be wound outward: rays follow the face normals, so a consistently inward mesh reads as entirely buried and comes back empty. Inconsistent winding is worse, because it fails quietly — the faces that disagree read as buried and are cut out of otherwise healthy membrane.

Python and Rust (fastcore::internals) only; there is no R wrapper yet.

Caps no longer span the opening they are closing. Two changes, one in each half of the capping path, both aimed at the same artefact: a hole coming back closed by a single fan of triangles reaching from one vertex to every other.

trace_loops now returns rings that are simple — no vertex twice. The walk is greedy, so at a pinch, where several boundary edges meet at one point, it could leave and re-enter the same vertex; what it traced was then a figure of eight rather than a polygon, and no triangulator is defined on one of those. Such a walk is now cut where it crosses itself and the pieces handed on separately. Same half-edges, grouped the way the caller can use.

triangulate_rings has a third attempt before it gives up. A ring gets past the first only if its flattening self-intersects, and since both of the existing attempts are projections — the Newell normal, then the best-fit plane — they tend to fail together, leaving the fan. Worth being clear that this is not the same as the ring being un-planar: a gently curved ring can cast a crossed shadow and a folded one need not, and on one real mesh the rings that fell through were on average more planar than the ones that did not. The new attempt ear-clips in three dimensions without flattening at all, taking the cheapest ear by area + 0.05 * perimeter². That is a heuristic about the shape of the cap rather than a guarantee about it, which is why it sits after the two that can promise more.

Over the 933 openings of that mesh the three rungs take 94.0%, 1.1% and 4.9%, so the best-fit plane does still earn its place — it rescues about one in five of the rings the Newell normal cannot flatten. Its other entry, a ring whose signed areas cancel so exactly that there is no area-weighted normal at all, did not come up once.

On a 3.5 M-face neuron mesh with heavy invagination, capping 933 openings: the longest edge any cap triangle carries drops from 4.2 µm to 1.4 µm, and the number over a micron from 1003 to 5. The worst case, a 1223-vertex opening around a soma, goes from a median cap edge of 2.7 µm to 65 nm. Rings small enough that a fan was acceptable are unaffected in any way that shows.

mesh_connected_components can read connectivity across faces, not only across vertices. It now takes a connectivity argument with three readings, each strictly finer than the one before it, and each step drops a kind of junction.

"vertex" (the default, and the old behaviour) joins two vertices whenever a face names them both. "face" joins two faces wherever they share an edge, which drops the pinch points:

# Two triangles meeting at vertex 2 and nowhere else
faces = np.array([[0, 1, 2], [2, 3, 4]], dtype=np.uint32)

fastcore.mesh_connected_components(faces, n_vertices=5)
# array([0, 0, 0, 0, 0], dtype=uint32)   -- one component: you can walk through the pinch

fastcore.mesh_connected_components(faces, connectivity="face")
# array([0, 1], dtype=uint32)            -- two: you cannot step across it

"manifold" joins two faces only across an edge carrying exactly two of them, which drops the seams — an edge three or more faces deep belongs to no single surface:

# Three fins meeting along the spine (1, 2)
faces = np.array([[1, 2, 3], [1, 2, 4], [1, 2, 5]], dtype=np.uint32)

fastcore.mesh_connected_components(faces, connectivity="face")
# array([0, 0, 0], dtype=uint32)   -- the spine is a shared edge like any other

fastcore.mesh_connected_components(faces, connectivity="manifold")
# array([0, 1, 2], dtype=uint32)   -- three faces on it, so it joins nothing

So "face" splits a mesh into the pieces you could walk across — trimesh's split(only_watertight=False) — and "manifold" splits it into pieces that are surfaces, each with a well-defined inside, which is what you want before asking a piece for its volume or its winding. That last one is trimesh.graph.face_adjacency exactly, down to its group_rows(edges_sorted, require_count=2), and there is a test pinning the two together on a mesh where the distinction bites.

Both face readings label per face, which is the only place that answer can live: a pinch vertex belongs to several face components at once, so there is no per-vertex form of it. A boundary edge — one face — joins nothing under any of the three.

n_vertices is now optional, since the face readings do not need it, and threads is new: those two have to group the 3F edges the faces name, which they do with a parallel sort. They then differ by a single test on how many faces each edge came back with, so "manifold" costs no more than "face". The vertex pass is unchanged — still one serial sweep over a single integer array, still the default.

Each of those arguments belongs to particular connectivities, and passing one where it does not apply is an error rather than something quietly dropped — the rule smooth_mesh already follows for its per-method parameters. mesh_connected_components(faces, n_vertices, connectivity="face") would otherwise hand back an array of a different length than the caller sized it for, which is the one failure that looks like success.

In Rust the face readings are a second function, mesh::mesh_face_components, taking a manifold_only flag rather than a name; in R it is connectivity = "manifold", as in Python, with the usual partial matching.

segment_coords can hand back one flat array instead of a list. Pass flat=True (it is keyword-only) and the segments come back as a single array with a row of NaNs after each one — the form every plotting backend actually draws.

# One (M, 3) array, segments separated by NaN rows, rather than a list of (n, 3) arrays.
coords = fastcore.segment_coords(ids, parents, xyz, flat=True)

# `node_colors` is padded to match, so it still lines up row for row.
coords, colors = fastcore.segment_coords(ids, parents, xyz, node_colors=rgba, flat=True)

Callers that draw a whole skeleton as one line collection were splitting into per-segment arrays only to concatenate them straight back together; on a few hundred neurons that round-trip cost more than the segmentation did. Building the padded array directly is ~3x faster than the list form when colours are along for the ride.

The separators have to be NaN, so integer coordinates are promoted to float. Integer colours are not: those rows are never drawn, so they are padded with zeros and a uint8 RGB array stays uint8.

The segmentation now crosses into Python as CSR, which sped up the list forms too. The pyo3 generate_segments used to hand back a Vec<Vec<i32>> — a PyLong per node index, which the caller converted straight back into numpy arrays. Together that cost roughly twice what generating the segments did. It now returns (nodes, offsets, lengths) in the same CSR form as trace_loops and matches_above: segment i is nodes[offsets[i]:offsets[i + 1]].

Both Python wrappers take that path, so this is a straight speed-up and not an API change — navis_fastcore.generate_segments still returns its list of arrays. On a few hundred neurons that is ~2x, and segment_coords(flat=False) ~3.5x: it also stopped using np.split, which swaps axes twice per segment.

Only the extension module changed. fastcore::dag::generate_segments, which the R bindings call, is untouched — hence Rust and R see nothing of this release.

0.12.0 (2026-08-10)

The skeleton smoothers take any per-node column, not just coordinates. smooth_skeleton now accepts an (N, K) field and smooth_skeleton_gaussian a separate values array, so a radius, a confidence or anything else numeric smooths by the same code as an x — matching what navis.smooth_skeleton's to_smooth can ask for. A (N, ) field comes back as (N, ).

# The window is a node count, so there is nothing to measure and the field is the only array.
xyzr = fastcore.smooth_skeleton(ids, parents, np.column_stack([coords, radius]), window=5)

# The kernel is a distance along the neurite, so the geometry stays a separate argument.
r = fastcore.smooth_skeleton_gaussian(ids, parents, coords, sigma=2000, values=radius)

That asymmetry is deliberate and is the whole reason the Gaussian gained an argument rather than a wider one. Its weights come from distance along the neurite; hand it a radius column as though it were geometry and "distance" becomes the cumulative absolute change in radius — a plausible-looking number and a meaningless kernel, and nothing downstream could tell. Columns are independent either way, so stacking them is exactly equivalent to a call each and one pass cheaper. Smoothing the coordinates costs exactly what it did before: with no separate geometry to gather, that path still reads its arc lengths out of the very buffer it is averaging.

R gets the same values argument on both functions, taking a numeric vector or an (N, K) matrix and handing back whichever shape it was given.

Note that this does not make the smoothers numerically interchangeable with navis, which they have never been: navis smooths with a trailing rolling(window, min_periods=1).mean() — half a window of lag towards each segment's distal end — and lets branch points move. Here the window is centred, shrinks symmetrically at segment ends, and endpoints do not move at all.

Two fixes to the smoothers, both of which change results from 0.11.0. All three surfaces, since the R crate bundles the core.

smooth_skeleton_gaussian reflected each segment about its ends once. A kernel wider than the segment ran off the end of the mirrored copy as well, and was then summing more neurite on one side of a node than the other — which bends a perfectly straight twig, the one thing the reflection exists to prevent. Under sigma=2000, a straight 10-node twig at 200 nm spacing bowed by 353 nm, more than the node spacing itself; most of a real arbour is segments that short. Worse, raising truncate did not help, because the error is not truncation: the sum converges to the one-reflection answer, not to the right one. The walk now turns around as many times as the kernel needs, so a straight line comes back exactly — and so does any linear ramp, the extension being odd. Nothing changes where truncate * sigma fits inside the segment; there a single reflection already covered the kernel.

smooth_skeleton rounded an even window up: window / 2 gave 4 the same two-node reach as 5. The Rust, Python and R docs all said it rounds down to the odd value below, and down is the better reading — window is a budget, and averaging six nodes when five were asked for is the surprising direction to round — so the code now matches what they promise. 2 therefore lands on 1, joining 0 and 1 as a no-op.

Packing shapes onto a page. Three new primitives for laying neurons out as a collage — rasterize_segments turns line work into binary masks, pack_masks arranges those so that no two share a pixel, and pack_rectangles is the cheaper bounding-box alternative (MaxRects).

masks = fastcore.rasterize_segments(coords, edges, scale=40, pad=1)
positions, variant, grid = fastcore.pack_masks(masks, (1300, 980))

Packing the shapes rather than their boxes is the point: a bounding box is mostly empty for anything branching, and on masks a neuron may reach into another's empty space — even sit inside a loop of it — as long as no cable meets.

The way to write this in numpy is a cross correlation, which scores every position at once so the free ones are the zeros. It is the wrong computation: an exact overlap count everywhere, in floating point, when the question is boolean and the answer is wanted at one position. Both cost models put a total order on positions that does not depend on the shape, so the best free position is the first free one — scanned in that order over bit-packed rows, 64 pixels per instruction, stopping there. On 200 arbors and a 1300x980 page: rasterising 32 ms → 1.8 ms, packing 1.55 s → 25 ms bottom-up and 4.39 s → 69 ms under a cost surface. It grows as O(N² res²), so the gap widens with resolution.

The search for a position runs on every core, and so do the variants of a shape against each other; a shape with enough edges splits its own walk, so one outsized mesh in a neuron list does not set the pace for the whole batch. Placement itself is sequential and stays that way — each shape goes down against the page the one before it left behind.

Both packers order items largest-first with a stable sort where np.argsort defaults to an unstable one, so ties go in input order; otherwise the results are identical to the correlation, bit for bit.

0.11.0 (2026-08-07)

Drawing a mesh flat, in one pass instead of six. project_mesh_2d takes a mesh and a view and hands back the polygons a 2-D renderer draws: projected onto the view plane, back faces dropped, sorted furthest-first so that painting them gives correct occlusion, and laid out as closed rings — plus the bounding box, reduced on the way past.

rings, bbox, ix, depth, normals = fastcore.project_mesh_2d(
    vertices, faces, xy_ix=(0, 1), depth_ix=2, front=1
)

Written the obvious vectorised way in numpy those are six passes over hundreds of megabytes, four of which exist only to feed the next one — 1.27 s and ~900 MB of intermediates on an 8.4M-vertex, 16.9M-face neuron. Fused and parallel: 133 ms, same faces, same geometry, same box.

The cull is the interesting part. Whether a face points at the viewer is the sign of its normal's depth component, and that component is a 2x2 determinant of the two columns being projected onto — so it never forms the other two components of the cross product and never reads the depth column, while applying exactly the test the full cross product would.

order=False skips the sort and the depths, for a caller filling the whole mesh as one path in one colour: a nonzero-winding fill is blind to the order its subpaths arrive in. normals=False skips the face normals, for a caller that is not shading. Smooth shading is the one thing this does not do — averaged vertex normals need every face, back-facing ones included, so they cannot come from the survivors alone.

In navis.plot2d this is most of a 17M-face mesh's plotting time: 3.32 s to 0.19 s, and the peak memory from 4.75 GB to 1.76 GB, for a byte-identical PNG. Python only for now; nat.fastcore does not have it.

See Meshes.

You can stop fastcore from taking every core. set_num_threads sizes the thread pool for the whole process, and get_num_threads reads it back — in Python and in R. The default is still every core the process can see, which is right for one call in one process and wrong the moment the caller is the one running things in parallel: nothing tells a worker process that it is one of twenty, so every worker claims the whole machine.

That is not a theoretical cost. On a 224-core node, navis.heal_skeleton(nl, parallel=True, n_cores=20) runs 20 workers x 224 threads = 4480 threads over 224 cores, and healing 40 skeletons of 200k nodes measured slower than the same work on a single worker (6.71 s vs 5.10 s) while burning 2.3x the CPU. Capping each worker to one thread: 3.60 s, at a sixth of the CPU.

import functools, navis, navis_fastcore as fastcore

navis.compute.worker_init_hooks.append(
    functools.partial(fastcore.set_num_threads, 1)
)

The pool is built once per process, so call it before anything else; calling it again with the same value is a no-op (which is what makes it safe in a hook that fires per chunk), and with a different value raises. RAYON_NUM_THREADS remains the no-code-change equivalent.

heal_skeleton, stitch_fragments and geodesic_pairs also grew the per-call threads argument the rest of the library already had — they were the last parallel entry points without one. See Threads, and scripts/profile-heal-parallel.py for measuring where your own workload sits.

Capping the holes a cut leaves in a mesh. Four functions, in Python and in R: boundary_halfedges finds every edge of a mesh with only one face on it, exposed_halfedges finds only the ones a subset is about to open (given the faces before the cut), trace_loops walks either into closed rings, and triangulate_rings ear-clips those shut. Only faces are ever added, never vertices, so every vertex index a caller already holds still means what it meant — which is what lets the cap be applied after a subset rather than during it.

halfedges = fastcore.boundary_halfedges(faces)
rings, offsets = fastcore.trace_loops(halfedges)
faces = np.vstack((faces, fastcore.triangulate_rings(rings, offsets, vertices)))

Almost all of the win is in the first one. Grouping the 3F edges a face array names is the whole cost of finding a boundary, and the obvious numpy spelling — np.unique(keys, return_inverse=True, return_counts=True) — is a stable argsort: 75 ms of an 84 ms call on a 578k-face mesh. That is not a formulation problem, which is worth saying because it is the kind that usually is: the bare np.sort of the same keys is already 51 ms, so no rearrangement inside numpy can win. Sorting bare u64 keys in parallel and taking a second pass over the faces to recover each boundary edge's direction brings the call to 8 ms. On the same mesh with ~23k holes punched into it, end to end: 224 ms of numpy becomes 11 ms. On the subset path — 400 twig cuts exposing 4.3k half-edges — 10.7 ms becomes 0.87 ms.

There is one caveat worth stating rather than burying: boundary_halfedges, exposed_halfedges and trace_loops reproduce the numpy implementation's output half-edge for half-edge and ring for ring, but triangulate_rings does not. It agrees on about 93% of rings and returns an equally valid alternative on the rest — same triangle count, same total oriented area, same winding. Depend on the hole being closed the right way round, not on the exact triangles.

trace_loops is greedy rather than a cycle basis, and that is deliberate. At a non-manifold boundary vertex several half-edges leave at once; taking whichever is still free puts every half-edge in exactly one ring, so the whole boundary is covered. networkx.cycle_basis — what trimesh.repair.fill_holes uses — quietly drops the edges that are not part of a simple cycle, and those holes stay open.

The ear-clipping is a Rust port of mapbox's earcut rather than a binding to it, which means the triangulation no longer needs the mapbox_earcut extension module at all — and that turns out to fix a hang as well as remove a dependency. Greedy tracing can walk back through a non-manifold boundary vertex, leaving a ring that names the same vertex twice; on a punched neuron mesh roughly 10% of rings are like that, and mapbox_earcut loops forever on the best-fit-plane retry one of them provokes. This implementation falls back to a fan, which is the right answer for a polygon that touches itself.

Mesh smoothing, without the shrinkage and without the translation. smooth_mesh moves a triangle mesh's vertices to take the noise out of its surface, leaving the face array, the vertex count and the vertex order untouched — so anything indexed by vertex is still attached to the vertex it was attached to.

smoothed = fastcore.smooth_mesh(faces, vertices)

That default is Taubin's λ|μ filter rather than the plain Laplacian, and it is worth saying why. The Laplacian step removes high frequencies quickly and low ones slowly, and a closed surface's enclosed volume is a low frequency: at lamb=0.5 and five iterations — what navis.smooth_mesh ships today — a neuron mesh comes out having lost 88% of its volume. Taubin alternates a shrinking λ pass with an inflating μ pass tuned so the two cancel below a cut-off frequency; on the same fixture it holds its volume to within 5%. method="laplacian" and method="humphrey" (the HC filter of Vollmer et al.) are there when you want them.

weights="cotangent" is the discrete Laplace–Beltrami operator, and the reason to reach for it is that the uniform umbrella cannot tell "this vertex is off the surface" from "this vertex has closer neighbours on one side", so on uneven tessellation it slides vertices along the surface. Cotangent weights are a function of the shape rather than of the triangulation, so they move vertices along the normal instead — on a UV sphere they drift less than half as far for the same amount of smoothing. trimesh builds its operator once from the input geometry and reuses it; here the geometry-dependent weightings are recomputed from the current positions every pass, which is the flow they are supposed to discretise rather than a snapshot taken before the first step. That is free, because the weights are never materialised at all.

The volume correction scales about the mesh's centroid, not about the origin. This is the one place the result deliberately differs from trimesh.smoothing.filter_laplacian, which is what navis.smooth_mesh calls, and the difference is not cosmetic. Upstream rescales by (vol_before / vol_after) ** (1/3) about the origin, which is not a shape operation: on the 722817260 test neuron at navis' own defaults it displaces the result by 41 µm, and the mesh is 19–26 µm across. It is also not translation invariant — the same mesh smoothed at two different offsets comes out two different shapes, and far enough from the origin the volume ratio goes negative and the cube root returns NaN — and it divides by the smoothed volume, so a mesh with a hole big enough to make that zero is a ZeroDivisionError rather than a diagnostic. Scaling about the centroid is the same size change with none of that, and where the ratio is genuinely undecidable (a flat sheet, both volumes exactly zero) the mesh comes back unscaled with a warning rather than silently wrong.

The correction also runs once, at the end, which is not an approximation of running it every iteration but exactly equal to it: every filter here is an affine combination of a vertex and a normalised average of its neighbours, and those commute with a uniform scaling.

preserve_border pins the rim of an open mesh — a boundary vertex's one-ring lies entirely to one side of it, so without this the rim rolls inwards under any of these filters — and lock freezes an arbitrary set on top of that, bitwise, while still letting those vertices pull on their neighbours. Same name and same meaning as simplify_mesh's lock, since it is the same concept.

Ten iterations with the volume correction on a 421k-vertex / 881k-face mesh: 5.42 s with trimesh, 0.03 s here (0.06 s on cotangent weights). The arithmetic was never the cost — the sparse matrix–vector product is 42 ms of upstream's 5.4 s. What it spends is 57% building the operator, where vertex_neighbors is a list of 421k Python lists costing 636 MB of heap for 10 MB of vertices, and another 40% in the volume constraint's per-iteration vertices[faces] gather. Against trimesh on the uniform umbrella all three filters agree to ~1e-11 on coordinates spanning 25,880 units, which is what the test suite checks; the deliberate divergences are the volume correction, the per-pass weight recomputation, and that one Taubin iteration here is a full λ/μ pair rather than a half-step.

Available in Python and in R. See Meshes.

Mesh simplification that remembers where every vertex went. simplify_mesh decimates a triangle mesh by quadric-error edge collapse and returns, alongside the smaller mesh, a vertex_map: for each vertex of the original, the index of the vertex of the simplified mesh it ended up in, or -1 if it did not survive.

That map is the point. Every other simplifier hands back a mesh and nothing else, so per-vertex data — synapses, radii, compartment labels — is orphaned by the operation, and the usual workaround of re-attaching by nearest neighbour afterwards is both slower and wrong: a collapse moves its survivor to the quadric-optimal point, which is often nearer some other vertex than the one that actually merged into it. bincount over the map replaces the spatial query entirely.

lock pins a set of vertices: a locked vertex is never merged into another and never moved, so it comes back at bitwise the same coordinates. It may still absorb its neighbours, which is what keeps a face target reachable when the pinned set is large. simplify_mesh_lossless is the other mode — collapse only what costs nothing, run to a fixed point — for shedding over-tessellation rather than hitting a budget.

This is a port of Sven Forstmann's Simplify.h (MIT), the algorithm pyfqmr wraps, rather than a wrapper around an existing crate. Wrapping was the preferred route and none of the candidates survived two constraints. meshopt has exactly the right semantics but vendors C++ and builds it through cc, which would cost the pyodide wheel and the R source tarball — the same reason flate2 and kodama are pinned to pure-Rust backends. alum and baby_shark are pure Rust but built on a halfedge and a corner table respectively, so they need manifold input: the first returns Err(ComplexVertex), the second silently drops the offending faces and returns an empty mesh on any build error. Meshes out of EM segmentation have edges shared by three faces as a matter of course. And none of the three expose a collapse map, so each would have needed a fork anyway. The algorithm as written is flat index arrays with no adjacency invariants to violate, which is what makes non-manifold input merely data.

Because it is the same algorithm, the port is checked against pyfqmr directly: on a clean mesh the two produce identical face arrays and positions agreeing to ~1e-12, across face-count ratios, aggressiveness settings and both border modes. Three deliberate divergences, all about degenerate geometry: upstream normalises vectors unconditionally, and since every comparison against the resulting NaN is false, NaN silently defeats the two guards that exist to reject a bad collapse. Zero-area faces are dropped, absent normals are represented rather than faked, and a collapse landing on top of a neighbour is rejected.

Speed is a side effect rather than the aim, but it is not worse: measured end-to-end against pyfqmr on UV spheres at ratio 0.1, 19.5 ms vs 30.5 ms at 80k faces, 88 ms vs 129 ms at 319k, and 259 ms vs 395 ms at 979k — about 1.5x, scaling linearly, at roughly 4M input faces/second on one core. Single-threaded by nature — each collapse invalidates its own neighbourhood — but the GIL is released for the duration, so simplifying several meshes from a thread pool does scale. Deterministic run to run; the result does depend on face order, as it does for every implementation of this family.

No new dependencies. Available on all three surfaces as simplify_mesh and simplify_mesh_lossless.

Six ways to change how densely a skeleton is sampled, in a new downsample module on all three surfaces. Until now fastcore could only reduce a skeleton topologically, with simplify_skeleton dropping every slab node; there was nothing that reduced it geometrically. See the new Downsampling page.

  • downsample_skeleton — keep every Nth node of each segment. This is navis.downsample_neuron.
  • simplify_rdp / simplify_vw — the same job done by shape rather than by counting. Ramer-Douglas-Peucker drops a node when removing it would move the path by less than epsilon; Visvalingam-Whyatt removes whichever node contributes least area, which sheds detail more evenly when the simplification is aggressive.
  • resample_skeleton — the inverse: interpolated nodes at a fixed spacing, so a skeleton whose node density varies tenfold between neurites comes out even. It reports the input edge and fraction each new node came from, so radii (or anything else per-node) interpolate in one expression.
  • smooth_skeleton / smooth_skeleton_gaussian — move the coordinates and nothing else, over a window of nodes or a distance along the neurite.

All six work on the linear segments between roots, branch points and leafs and never touch the nodes at their ends, so the topology is unchanged: same leafs, same branch points, same tree. The three that drop nodes return the same tuple as simplify_skeleton, all take a preserve list of nodes that must survive whatever the rule decides, and all carry each dropped chain's length into the edge that replaces it, so total cable length and geodesic distances are preserved exactly.

Every one of them that changes the node table also says where your data went. A skeleton rarely travels alone — synapses, soma tags and manual annotations hang off particular nodes — and renumbering the nodes strands them. The four node-dropping functions and resample_skeleton now return a node_map: one entry per input node, naming the output node its data belongs to now. It is the same direction as simplify_mesh's vertex_map, and total, so unlike that one it needs no -1 masking: every input node names exactly one output node, the nearest along the neurite, ties going towards the root.

For resample_skeleton this is the half that could not be derived. source/alpha carries a per-node column forward, but it does not invert — an input node that fell between two output nodes has no output row of its own — so there was no way to answer "which new node does this synapse belong to" short of a spatial query, which is wrong across a hairpin. The two smoothers have no node_map and need none: they move coordinates only, so every node keeps its ID and its parent.

The map falls out of walks these functions already do, so it is cheap rather than free. On a 1M-node arbor: resample_skeleton 177.5 → 184.0 ms, simplify_rdp 48.7 → 52.3 ms, downsample_skeleton 39.2 → 44.4 ms, simplify_skeleton 8.3 → 11.1 ms. The share is largest where the function does least — simplify_skeleton is a bare pointer chase, so a 4-byte-per-node array and one store per dropped node are a real fraction of it, and 3 ms per million nodes.

Breaking: simplify_skeleton returns a 4-tuple rather than a 3-tuple in Python and Rust — it shares its rewiring with the three new droppers, so it gained the same node_map. R is unaffected: the return is a list, and node_map is a new element.

The three linestring algorithms are ports of Chris L. Barnes' simples (MIT), adapted for skeletons: iterative rather than recursive (a 100k-node unbranched neurite would overflow the stack), triangle areas from the Gram determinant rather than Heron's formula (which loses most of its digits on exactly the sliver triangles Visvalingam-Whyatt searches for), and a smoothing kernel over distance along the neurite rather than between the points. No new dependencies.

0.10.1 (2026-08-03)

heal_skeleton is reproducible again. Healing the same fragmented neuron twice could return two different skeletons. The bridge search runs in parallel and prunes with a per-fragment bound shared across threads; when two nodes tied for their fragment's shortest bridge, which one got to report it came down to which thread read that bound first. Ties are routine here, because skeleton coordinates come off a lattice. The healed skeleton was never wrong — the total added cable was identical every time, since the bound cannot drop below a fragment's true minimum — but which of several equally short bridges it used varied from run to run.

Each query now searches just past the shared bound, so every node achieving its fragment's minimum reports it however the threads interleave, and equal-length candidates are settled on their endpoints. Same bridges, same total length, same answer every run. Affects heal_skeleton and stitch_fragments on all three surfaces.

0.10.0 (2026-08-01)

Six new tree primitives, filling the gap between what navis asks igraph for internally and what fastcore could answer. Each looks like a general graph algorithm but is a linear pass over the parent vector on a rooted forest, so building a graph object to answer it costs more than the answer does. See the new Topology page.

  • descendants / paths_to_root — the two directions of the same walk: everything below a node, and everything above it. descendants replaces igraph's subcomponent(v, mode="IN") and is what makes "cut the skeleton here" a masking operation rather than a graph rebuild.
  • reroot — re-orient a forest at given nodes, reversing only the edges between each new root and the old one. Components nobody names are left byte-identical. Generalises topo::reroot_rewire, which takes one preferred root plus a set of new edges.
  • contract_nodes — collapse groups of nodes onto a representative and rewire what is left. Edges internal to a group are dropped rather than turned into self-loops; a mapping that would close a cycle is refused rather than silently returning a non-forest.
  • simplify_skeleton — keep only roots, leafs and branch points, with each replacement edge carrying the total length of the chain it stands in for, so cable length is preserved exactly. On the example skeleton that is 4332 nodes down to 1290.
  • adjacency — the CSR triple (indptr, indices, data) of the skeleton's adjacency matrix, with column indices sorted within each row. Handing back the raw arrays rather than a matrix keeps this package free of a scipy dependency.

Four more, completing the set navis needs:

  • longest_path / longest_paths — the longest path from a node to its root, and the n longest taken in turn with each peeled off before the next is sought. Not the NP-hard general problem: in a rooted forest every maximal path is fixed by its start node, so this is a distances-to-root question. Ties break towards the lowest index, matching numpy.argmax, which is what navis' implementation relies on for a stable answer.
  • betweenness — betweenness centrality in O(N) rather than Brandes' O(V·E). Shortest paths in a tree are unique, so the count through a node is a closed form: descendants × ancestors when directed, and a sum of products over the parts it separates when not. Counts are int64, because an undirected 100k-node skeleton reaches ~5e9.
  • descendant_counts — how many nodes, or how many of a given set, lie strictly below each node. See the note below on why this exists.

betweenness is not navis' betweeness_centrality(from_=...)

navis' from_ branch does not compute betweenness at all. It walks root→source paths and tallies every node except the source, which counts, for each node, how many of from_ lie below it — a descendant count. That is why descendant_counts is a separate function rather than a sources argument here: naming it betweenness would have made the two behaviours indistinguishable at the call site. navis' find_main_branchpoint(method="betweenness") — which is that function's default — is the one caller, and wants descendant_counts.

Three graph primitives, the ones skeletor still needs a graph library for. Where the tree primitives above serve navis, these serve mesh skeletonization — see the Meshes page.

  • parents_from_edges — orient an edge list into a rooted forest: one parent per node, -1 at the roots, cycles broken. minimum_spanning_tree picks which edges survive; this picks which way they point, which is what turns a bag of undirected edges into something you can walk, root or write out as SWC. It also returns the order the nodes settled in, which is free (the search visits them in it) and is exactly the relabelling that makes parents precede their children. One search covers the whole graph rather than one per component: the shortest-path-tree-per-component construction costs O(components x n_nodes) in output alone, which on a skeleton shattered into four thousand fragments is a 2 GB array for an answer that is one column. At 100k nodes, 2.9 ms against igraph's 14 ms for one arbor — and 2.7 ms against 4370 ms once it fragments, because igraph pays per component and this does not.
  • bridges — which edges may not be dropped without disconnecting their component (Tarjan, on an explicit stack, so a 200k-node strip does not overflow). The counterpart to minimum_spanning_tree rather than a variant of it: the MST asks what to keep to stay connected, this asks what may not be removed. Parallel edges are honoured — two nodes joined twice are joined by a cycle, so neither edge is a bridge, which is why this does not share the deduplicated adjacency the geodesic searches use. 2.6 ms against igraph's 13.5 ms at 100k nodes, 2.2 ms against 207 ms fragmented.
  • geodesic_mst_mesh / geodesic_mst_graph — the minimum spanning tree over a subset of nodes, weighted by geodesic distance through the graph they were carved out of, without materialising the k x k distance matrix. That matrix is k**2 distances to use k - 1 of them — 400 MB at k = 10_000 before the O(k^2) MST itself, and k searches to fill it. Following Mehlhorn's distance-network construction, one multi-source sweep partitions every node by nearest subset member and each edge straddling two cells offers one candidate; an MST over those is an MST of the full distance network. The cost is flat in k because it is one sweep whatever k is: at 100k nodes, 12.7 ms at k=250 and 8.3 ms at k=4000, against 187 ms and 7820 ms for the matrix route. Returned weights are exactly the geodesic distances between the pairs they join, so they are usable as lengths. limit bounds how far apart two nodes may be and still be joined, and prunes the sweep rather than merely discarding results.

The R bindings caught up. 22 new functions in nat.fastcore — every tree and graph primitive above, plus the ones that had accumulated unbound before it — taking R from 39 documented capabilities to 58, and to 66 of 77 with the clustering pair below. The signatures are the ones the Python side settled on, translated to R conventions: 0-based node indices throughout (as the rest of the DAG family already used), roots and "no such node" as -1, and multi-value results as a named list rather than a tuple.

New optional arguments have R defaults, so adjacency(parents) and parents_from_edges(edges, n) work without spelling out a NULL per argument. Note this is not true of the bindings that predate this release — those still require every argument positionally, which is worth fixing before 1.0.

Two things are worth doing before 1.0 alongside that. Argument errors currently reach R as Error: User function panicked: <name>, because extendr discards the panic payload and the R layer does no validation of its own — where Python raises a message naming the offending value. The fix belongs in the three shared converters, not per function. And what remains unbound is GeodesicGraph and its methods (a stateful pointer class, so a different kind of job) plus five of the NBLAST/matches helpers.

Every one of these is pinned against igraph in the parity suite and against brute-force references under hypothesis, across a fixture matrix that includes 100k-deep chains — the traversals are iterative precisely because the recursive versions segfault there.

The geodesic searches run in float64 on request. Dijkstra sums one weight per hop, so a path of k hops carries up to k roundings; at float32 and k in the tens of thousands — a densely sampled arbor, a fine mesh — the drift becomes visible against an exact answer, and comparisons against scipy.sparse.csgraph, which works in float64 unconditionally, stop agreeing to the last bits. It also matters when weights span a wide dynamic range, since fl(du + w) loses w entirely once du exceeds it by 2^24.

The width is now a type parameter on the core rather than baked into it: Adjacency, the search scratch, both kernels and every driver in mesh are generic over a new mesh::Weight trait, implemented for f32 and f64. The heap key stays an integer compare on the raw IEEE bits — see Weight::Bits for why that is an exact ordering and not an approximation — so the float32 path is unchanged, in both results and speed.

In Python the rule is your dtype in, your dtype out, the one linkage already follows for score matrices: float64 weights give float64 distances, anything else gives float32.

fastcore.geodesic_matrix_graph(edges, n, weights=w.astype(np.float64))  # -> float64
fastcore.geodesic_matrix_graph(edges, n, weights=w, dtype=np.float64)   # -> float64

A new dtype argument overrides that in either direction, on geodesic_matrix_graph, geodesic_matrix_mesh, geodesic_nearest_mesh, geodesic_farthest_mesh, geodesic_predecessors and the two geodesic_mst_*. Only something carrying a float64 dtype counts as having asked: a list of Python floats does not, because np.asarray([1.0, 2.0]) is float64 by numpy's default rather than by anyone's intent, and honouring it would quietly double the output of every caller who passes one.

This changes the return dtype for existing callers

If you already pass a float64 weights array, you were getting it cast down to float32 and a float32 result; you now get float64 — twice the output memory, and about 10% slower. Pass dtype=np.float32 to keep the old behaviour. Callers passing lists, int arrays or no weights at all are unaffected.

The mesh functions default to float32 and take dtype alone, with no input dtype read off vertices. Those are coordinates, taken at float64 either way — each edge length is computed from them at that width and rounded once on the way in — so reading the distances' width off them would flip nearly every existing call to float64 and double the largest thing this library allocates. A full V x V matrix is already 107 GB at V = 164k.

geodesic_path, geodesic_clusters, parents_from_edges and minimum_spanning_tree have no dtype argument, because none of them returns a distance — but all four honour the weights' own width, since it decides which route or which tie wins.

GeodesicGraph stays float32. It is the type for "large graph, many small queries", which is exactly the case where float32 is the right width and where doubling the several node-sized arrays it keeps resident across a whole run would be felt.

R gets a precision argument (32 or 64, default 32) on the eleven corresponding functions, matching nblast(precision = ). R has no float32 type, so unlike Python there is nothing to read the choice off — weights arrive as doubles whatever the caller meant by them, and the result goes back as doubles either way — so this buys accuracy, not a different return type.

Interface polish, ahead of 1.0. Small things that cost nothing to change now and get expensive once the API is frozen.

  • Integer returns follow one rule, written down in the Python overview: a node id is uint32; a node id needing a -1 sentinel, or a dense label such as a cluster id, is int32; and a position in an array you passed in is int64. The point is the last — int64 tells you the values index your array rather than the graph, so nodes[out] is the node-id form and out alone is not. Four returns moved to fit (see Breaking). The rule governs the index-space API; the tree functions work in ID space and hand back values in the dtype of the node_ids you gave them, which the overview now states explicitly.
  • from navis_fastcore import * exports functions only. The package had no __all__, so import * took the default and dragged in seven submodule objects while dropping __version__ for starting with an underscore. It is now composed from the submodules' own __all__, so a new function is exported by listing it in one place. parent_dist is exported for the first time along with it: it was public and documented but in no __all__, so fastcore.parent_dist did not resolve.
  • has_cycles is callable from Python. The core has had it all along and R binds it, but on the Python side it lived in the extension module only, used internally by the scipy interop shim — so the one function that tells you whether the input to everything else is well-formed was the one you could not call. It is now fastcore.has_cycles(node_ids, parent_ids), in ID space like the rest of the tree family, with a parent ID that is not a node treated as a root rather than as a cycle.
  • Four more that the Rust had and Python did not. An audit of the crate against the bindings, prompted by the one above; what it turned up is now bound, tested and documented. Only the first needed new Rust — the rest were already compiled in and merely private.
    • leaf_order — SciPy's leaves_list: the order to place the leaves in so a dendrogram draws without crossing branches. It was the one core function with no pyo3 wrapper at all, which meant the one step of the clustering story that still required scipy was drawing it. Iterative, so a 200k-observation chain does not need 200k stack frames, and it rejects a linkage matrix naming a cluster that does not exist yet rather than walking off the end of it.
    • nblast_pairs — NBLAST of an explicit (query, target) list, one score per pair rather than a matrix. k pairs cost k comparisons instead of n_query x n_target, which is the point when a cheaper filter has already told you which cells you care about. Smart NBLAST has used this internally for its full-resolution pass since it shipped; it is the same primitive, so a target whose whole column you request reproduces that column of nblast exactly.
    • reroot_rewire — turn an edited edge set back into a rooted forest. Step 2 of heal_skeleton, for callers who choose their own edges rather than taking the minimal bridges from stitch_fragments. Distinct from reroot, which re-orients an unchanged forest and leaves untouched components byte-identical: once the edge set moves there is no unchanged to preserve.
    • symmetrize — the in-place symmetrise, on its own. linkage and condensed_distances already fold it into their fused pass, so this is for when something else has to read the matrix. It is the case numpy cannot do cheaply: (M + M.T) / 2 builds two full n x n temporaries and even np.add(M, M.T, out=M) still builds one, where this allocates nothing.
  • CmtkRegistration.domain, the spline warp's domain box, which R has had as cmtk_domain. Points outside [0, domain] have no spline value — CMTK prints FAILED and xform returns NaN — so this is how you predict a NaN instead of reconstructing the box from .spacing and .dims yourself.
  • The same two, back the other way, in R. symmetrize and leaf_order are now exported from nat.fastcore as well, so the clustering family reads the same on both surfaces. Two differences are forced by R rather than chosen: symmetrize returns a copy, because R's value semantics forbid writing to the caller's matrix (it is still one n x n against (m + t(m)) / 2's two), and leaf_order takes an hclust or its merge matrix rather than a SciPy linkage matrix, returning a 1-based ordering in the same form as hclust$order. It agrees with stats::hclust's own order element on the trees that package builds, which is what "same child order" has to mean.
  • The capability tables gained a Clustering section. linkage, condensed_distances, symmetrize and leaf_order and their R counterparts had never appeared in them, so the one family where all three surfaces differ in the object they hand back — linkage matrix, dist, hclust — was the one you could not look up. Every public Python name now appears in some row.
  • GeodesicGraph.subset validated its nodes argument in three places — the Python wrapper, the binding layer and the core. The binding-layer copy is gone; the wrapper now uses the same unique=True check every other node subset in the package goes through.

Breaking

  • spanning_forest is now parents_from_edges. It sat one word away from minimum_spanning_tree while answering a different question — that one picks which edges survive, this picks which way they point — and the new name says what it hands back, which is the parent vector the rest of the package consumes. minimum_spanning_tree keeps its name: it is what scipy, igraph and networkx all call this, including scipy's behaviour of returning a forest when the input is disconnected. No alias, since spanning_forest was never in a release. Affects all three surfaces.
  • Node ids come back as uint32 where they used to be int64, under the rule above: unique_edges' edges, contract_vertices, and parents_from_edges' order (the last of these unreleased). Their index / inverse companions stay int64 — those are positions in the 3F edge list, not node ids. Python callers who fed these straight back in were already coercing to uint32; those coercions are now no-ops rather than copies. R is unaffected — the bindings already narrowed to R integers.
  • matches_above' indices is now int64, not uint32. It is a position along the scanned axis of the scores matrix you passed in — the same quantity top_matches returns, which was already int64. The two disagreed on the dtype of one thing. The ragged array is the larger of the two, so this does cost memory; consistency at the call site is worth more than the width. MatchError::AxisTooLong goes with it — it existed only to refuse a scanned axis longer than u32::MAX, which is no longer a limit.
  • generate_segments now measures a segment's length from its first node to its last. It previously summed the weight vector over every node in the segment, including the terminal one — but a segment stops at a branch point, whose own child→parent edge continues into the parent segment. Every segment ending at a branch point was therefore one edge too long; segments ending at a root were already correct, because a root's weight slot is 0. Unweighted lengths change with it, from a node count to an edge count, so that weights=None stays equivalent to weights=ones — the same correction dist_to_root had in 0.6.0. Segments themselves are unchanged; only lengths moves (and, where lengths tie differently, the order they are sorted in). Affects all three surfaces.

Fixes

  • results carrying a "no such node" sentinel no longer wrap around on uint64 node IDs. geodesic_nearest, geodesic_farthest and heal_skeleton built their output with np.full(..., -1, dtype=node_ids.dtype) and then wrote -1 into it, where it wraps — so an unreachable source or a root came back as 18446744073709551615, on exactly the uint64 IDs segmentation backends hand out. All of them now go through the helper that already handled this for reroot and friends, which promotes to int64 when the ID dtype cannot represent the sentinel. The sentinel marks a source with no reachable target (so, any skeleton with more than one component) and the healed skeleton's root, which is what made this easy to hit.
  • _ids_to_indices no longer raises on an empty ID array when the node and target dtypes differ — it took max() of both unconditionally. Reachable from any function taking an optional set of node IDs (descendant_counts(targets=[]) and friends) whenever node IDs are uint64 and the target array is int64, which is navis' normal convention.
  • geodesic_matrix(directed=True) no longer leaks across zero-weight edges when sources or targets are given. The partial backend used depth as a proxy for ancestry, which holds only while every edge weight is strictly positive: a zero-weight edge gives an ancestor the same depth, so it slipped through the guard and was written at distance 0 — reporting a parent as reachable from its child's direction. Coincident nodes are routine in traced and resampled skeletons. The all-by-all backend, geodesic_nearest, geodesic_farthest and geodesic_pairs were unaffected.

0.9.0 (2026-07-28)

mesh.GeodesicGraph — build the adjacency index once, query it many times. The geodesic free functions each build an index, answer one question and throw it away: the right trade for a single sweep, the wrong one for algorithms asking many small questions of one graph. Every existing query is available as a method (distances, nearest, farthest, predecessors, path, clusters, components), so migrating is mechanical, and subset carves an induced subgraph out of the built CSR instead of sending you back to numpy to mask and renumber an edge list. 500 short-path queries on a 40k-vertex mesh run ~100x faster as methods. It buys nothing measurable when one query already sweeps the graph.

Four operations are new, because they only make sense against a graph you keep:

  • grow — a connected region of a fixed number of nodes (or of attached cloud points), plus each one's distance to the seed. 33x per patch on a 160k-vertex mesh; 37x over a 53k-point cloud.
  • farthest_seed — the next farthest-point seed, for spreading patches evenly. An incrementally folded and self-pruning distance field plus a lazily-corrected max-heap: 2560 seeds on a 160k-vertex mesh, 92.5 s → 0.35 s (265x).
  • ball — everything within a radius of a set of nodes, how far, and which source is nearest. Returns the ball itself rather than three node-sized arrays with the ball buried in them.
  • set_weights — re-weight edges in place, O(edits log valence) against the O(E) of a rebuild. TEASAR zeroing each path it extracts is the motivating case.

Pinned against navis' own pure-Python implementations (navis.ml.chunk) and scipy. One caveat: the search is float32 where those references are float64, so on graphs with edge lengths that tie in float32 the settle order of equally-distant nodes can differ. Regions are still exact balls — verified against a float64 oracle.

Pyodide / JupyterLite wheels (experimental). A wasm32-unknown-emscripten wheel is now built, run against the full test suite under Pyodide in CI, and published to PyPI alongside the native wheels under PEP 783. Emscripten cannot spawn threads, so NBLAST and the transforms run serially there and are not interruptible — a documented degradation rather than a runtime panic.

Python-facing release; R is unchanged.

0.8.0 (2026-07-28)

NBLAST peak memory cut ~45% — 2.22 GB → 1.26 GB on 6.9M points. Two independent reductions: f32 coordinate storage, and aann 0.3.0's u32 neighbourhood graph. Scoring is unchanged — the accumulator, tangent dot products and every score-matrix lookup still run in f64; only storage and the descent's own distance comparisons narrow. Python selects the width from the input dtype, so float32 points/vect now reach the index with no copy (previously they were silently upcast, which cost more than passing f64). Expect ~1e-4 relative movement on scores: on 40k real neurons, 98.7% of k-NN rows are bit-identical and top-1 agreement is 99.98%.

dotprops() — tangent vectors and alpha from a bare point cloud. An exact k-d tree k-NN plus a parallel Jacobi symmetric-3x3 eigensolve, replacing cKDTree.query + N SVDs: 10.5x at N=10k, 15.5x at N=100k. This was the last thing pulling in scipy, so navis-fastcore no longer imports it anywhere. Also available as Dotprop.from_points().

Graph primitives on mesh, straight off an edge list — no graph object to build first:

  • level_set_components — connected components of every label's induced subgraph in one DSU sweep, replacing a per-level loop (12.2 ms → 0.3 ms on a 41k-vertex mesh)
  • geodesic_path / geodesic_predecessors — the route, not just its length (11.7x including build; 2.0x against a cached igraph object)
  • geodesic_clusters — greedy partition into connected clusters of bounded geodesic radius
  • connected_components_graph, contract_vertices, minimum_spanning_tree

New APIs are Rust + Python. R gets the synced core and the NBLAST memory work; wrappers for the new mesh/points functions are still to come.

0.7.3 (2026-07-20)

New

  • thin plate spline and moving least squares landmark transforms (TpsTransform, MlsTransform) — for when there is no registration file, only matched landmarks
  • nblast_knn: k nearest neighbours without a score matrix

Fixes

  • MlsTransform.xform gained a reverse parameter
  • custom scoring matrices are now checked for conformance
  • picked up the upstream tie-breaker fix in aann
  • Python: Dotprops is exposed at top-level

0.7.2 (2026-07-20)

New

  • NBLAST support functions across Rust/Python/R: linkage and condensed_distances in Python, fast_hclust, nblast_dist and nblast_hclust in R — hierarchical clustering straight off a score matrix, without ever casting it (at 100k a side that would quietly materialise tens of GB)

0.7.1 (2026-07-19)

New

  • unique_edges: the unique undirected edges of a triangle mesh

0.7.0 (2026-07-15)

New

  • Rust re-implementations of CMTK and Elastix transforms — both several orders of magnitude faster than the original binaries (streamxform / transformix), and neither tool needs to be installed. Bonus: fastcore can invert Elastix transforms, which Elastix itself cannot.
  • functions to extract top matches (N, threshold, percentile) from NBLAST matrices: indices, values = fastcore.top_matches(scores, 5, skip_self=True)

Improvements

  • another ~2x speed-up for NBLAST, from upstream improvements in aann

0.6.1 (2026-07-12)

New

  • subtree_height() returns, for each node, the geodesic distance to the farthest leaf
  • dist_to_root() returns, for each node, the distance to its root

0.6.0 (2026-07-12)

New

  • geodesic_farthest(parents, sources, targets) implements the opposite of geodesic_nearest()
  • geodesic_matrix_mesh and geodesic_farthest/nearest_mesh compute geodesic distances on meshes:
    • they operate directly on faces/vertices; in navis this needs a costly round trip meshigraphcsgraph.dijkstra
    • they use threads
    • they allow defining sources and targets (csgraph only does sources), which cuts the memory footprint drastically and can be faster

Improvements

  • better performance (speed / scaling / memory) for classify_nodes, strahler_order, geodesic_matrix, geodesic_nearest, all_dists_to_root, has_cycles, geodesic_distances, geodesic_pairs, node_indices_* and extract_parent_child
  • guarded geodesic_distances against segfaults on very large (100ks of nodes) neurons
  • geodesic_distances now computes in f64 internally but still writes f32, improving relative error from 3.9e-06 to 7.0e-08

Fixes

  • generate_segments now returns segment lengths in the correct order

Breaking

  • geodesic_pairs now correctly returns -1 for unreachable pairs (previously an incorrect 1)
  • dist_to_root now counts edges instead of nodes (previously the root was reported at distance 1)

0.5.1 (2026-07-11)

Fixes

  • refactored heal_skeleton to deal with some pathological cases

0.5.0 (2026-07-11)

New

  • skeleton healing functions (Rust + Python/R bindings)

0.4.0 (2026-07-09)

Improvements

  • NBLAST is rebased on aann for all-nearest-neighbour lookup plus shull for Delaunay triangulation — together a ~5x speed-up over navis' built-in implementation
  • precompiled binaries for the R bindings are now on R-universe

0.3.0 (2026-06-28)

New

  • geodesic_nearest returns, for each query node, the closest node among a set of targets and the distance to it. The full distance matrix is never materialised, so it scales.

0.2.0 (2026-06-15)

New

  • connected_components_mesh(): fast connected components on meshes, much faster than e.g. conversion to igraph and subsequent graph-based CC

Improvements

  • large performance increase for connected_components(), especially when the skeleton has many roots
  • speed-up for segment_coords()
  • technically supporting Python's free-threading (3.14+), though this is largely untested

0.1.0 (2025-03-13)

Improvements

  • support for 16-bit node IDs
  • better handling of node and parent IDs with different dtypes
  • csgraph drop-ins: added documentation and proper checks for whether the input is a rooted tree

Fixes

  • fixed an incorrect parameter name in the csgraph drop-ins

0.0.9 (2025-02-06)

New

  • generate_segments accepts a lengths parameter, sorts segments by total length when given, and returns the segment lengths

0.0.8 (2025-01-05)

New

  • prune_twigs() gained a mask parameter to restrict pruning to parts of the neuron

0.0.7 (2024-09-16)

Fixes

  • geodesic_matrix() no longer returns the distance matrix in the wrong order when sources and/or targets are given

0.0.6 (2024-09-16)

New

  • geodesic_pairs: geodesic distances between given node pairs

0.0.5 (2024-09-07)

Improvements

  • more flexible about 32-bit vs 64-bit node/parent IDs
  • more checks before calling into Rust, so exceptions are more helpful

0.0.4 (2024-07-25)

New

  • reorganised into a monorepo holding the fastcore Rust crate, the navis-fastcore Python bindings and the nat.fastcore R bindings
  • navis-fastcore: drop-in replacements for some csgraph functions

Fixes

  • fastcore: break_segments no longer drops the last node of a sequence

0.0.3a (2024-07-09)

New

  • new functions: prune_twigs(), strahler_index(), break_segments() and classify_nodes()

Improvements

  • synapse_flow_centrality(): added a mode parameter (centrifugal, centripetal or sum)

Fixes

  • geodesic_matrix() now actually respects directed=True

Breaking

  • the module was renamed from fastcore to navis-fastcore to avoid name clashes:

    pip install navis-fastcore
    
    import navis_fastcore as fastcore
    

0.0.2 (2024-06-27)

New

Improvements

  • a Rust implementation of geodesic_matrix for specific sources and/or targets, which is faster and much more memory efficient than the all-by-all as long as the number of sources/targets is small-ish. fastcore.geodesic_matrix picks the right one for you.

Breaking

  • removed the modifier parameter from fastcore.segment_coords

0.0.1 (2024-06-22)

Mostly a release to test publishing via CI.