Connectivity

FlowGeometries.Connectivity.AbstractImageConventionType
AbstractImageConvention

How a ball query treats a periodic direction: NearestImage or AllImages.

Singleton types rather than a Bool, for the reason given above about stencils: the traversal branches on this per candidate, so a runtime value leaves the coordinate expression unresolved and the whole walk allocates — measured at 14 KB for one query on a 32² grid, against nothing when it is a type.

source
FlowGeometries.Connectivity.CSRConnectivityType
CSRConnectivity{VN,VP}

Sparse neighbor list: node i owns nbrs[ptr[i]:ptr[i+1]-1]. The two buffers are typed independently, and their element type is any Integer, so a large mesh can hold Int32 indices (half the memory and bandwidth of Int64, and the width GPU kernels want).

source
FlowGeometries.Connectivity.ConnectedType
Connected()
Connected(stencil)

The cells within ball reachable from the seed through cells that are themselves within ball and active: the connected component of the ball that contains the seed. A graph query — it expands along adjacency and prunes at the ball's edge — so it is a subset of Unrestricted. The two agree whenever the ball is connected under the adjacency, which a maskless Cartesian StructuredGrid always is: each index step moves monotonically in one coordinate, so a staircase path to a cell never leaves that cell's own distance. Connected is strictly smaller wherever a mask or a boundary separates two parts of the ball.

Adjacency has to be named, since only a node set carries its own: with no argument, direction-1 adjacency (Stencils.Axial(1)) on the index-space architectures and the stored neighbour lists on an UnstructuredGrid. Connected(stencil) sets it explicitly for the former, and is an error for the latter, which has no index space for a stencil to mean anything in.

This is not a cheaper way to compute Unrestricted. Take cells P (the seed), Q and R, adjacent only as P–Q–R, with d(P,Q) = 1.2r and d(P,R) = 0.8r. A walk outward from P that drops any cell farther than r stops at Q and never reaches R, though R is inside the ball. That walk is not a broken Unrestricted — it is exactly Connected, which is why it cannot produce the other one.

source
FlowGeometries.Connectivity.IndexTopologyType
IndexTopology(size, periodic, mask)
IndexTopology(grid)

Extent, wrapping and activity per dimension — the whole of what a neighbor computation reads. Coordinates, cell measure and geometry never enter one, so a sampling can hand this over directly rather than materializing a grid (axes, dense measure, full mask) to be read once and dropped. A curvilinear grid is the N = 2 case of the same algorithm, not a separate one.

mask === nothing means every cell is active, and costs no storage and no load.

source
FlowGeometries.Connectivity.MetricTopologyType
MetricTopology(grid; index = nothing)

Everything a distance query reads that depends on the grid alone — the counterpart of IndexTopology for the metric path.

A stencil query reads (size, periodic, mask) and IndexTopology carries it. A ball query additionally needs the tightest per-direction step bound, which sizes the search window, and — where there are no separable axes to bound with — a spatial index.

Constructing one is O(1) and allocates nothing: the per-axis reductions live on the grid already, as Grids.AxisStats, computed once when it was built. So the default topology on every query costs nothing and there is no hoisting to remember. What is still worth hoisting is the index, which is not built by default because a k-d tree inside a single query would cost more than the scan it replaces — foreach_within and mapreduce_within do that hoisting for a sweep, and indexed does it explicitly.

It is not cached on the grid: grid types are immutable and the Adapt extension reconstructs them field-by-field for a device, so a mutable cache field would break both that and thread safety.

index holds a spatial index when one is available, for the architectures with no separable axes to bound.

source
FlowGeometries.Connectivity.StencilNeighborsType
StencilNeighbors{G,N,S}

Lazy neighbor sequence of one cell of an index-topology grid: iterating it walks the stencil offsets and yields the linear index of each in-range, active neighbor.

Nothing is stored, so a traversal that visits every cell allocates nothing at all — where returning a freshly built Vector per cell would cost two heap allocations per cell. Use neighbors! to write into a caller-supplied buffer, or collect this to materialize it.

source
FlowGeometries.Connectivity.UnrestrictedType
Unrestricted()

Every cell whose centre lies within ball. The default, and a purely spatial query: with a spatial index it costs O(log n + m) and never consults adjacency.

Note that the result is a ball, not a connected patch — with a mask, or a concave domain, it can contain cells that are near the seed in space but reachable from it only by leaving the ball.

source
FlowGeometries.Connectivity._csr_from_candidatesMethod
_csr_from_candidates(emit!, n, maxdeg) -> CSRConnectivity

Build CSR for n nodes whose degree is bounded by maxdeg. emit!(buf, lo, i) writes node i's candidate neighbors into buf[lo+1 : lo+maxdeg] and returns how many it wrote; duplicates, self and out-of-range entries are removed here. emit! must touch only that slice and carry no state between calls — nodes are emitted concurrently under a threaded backend.

Candidates land in one n*maxdeg block that is then compacted in place down to the exact CSR — one allocation for the neighbor list, one for the offsets.

source
FlowGeometries.Connectivity._csr_from_undirected_edgesMethod
_csr_from_undirected_edges(nnodes, edges) -> CSRConnectivity

Build CSR from an undirected edge list, by counting degrees first and then filling each node's slot range directly — no intermediate per-node vectors, and no reallocation while filling.

source
FlowGeometries.Connectivity._cubed_neighborMethod

Exact face-neighbor under offset (di,dj) for the gnomonic cubed sphere matching SphericalSampling._cubed_face_to_xyz / cubed_sphere_points!.

Panel interiors stay on-face. Crossing an edge uses cube face adjacency with index maps derived by matching cube XYZ along shared edges. Diagonal (corner) exits return (0,0,0) — no unique adjacent face.

source
FlowGeometries.Connectivity._default_node_areasMethod
_default_node_areas(sampling, geometry, λ, φ) -> Vector

Cell areas to use when the caller supplies none.

Dispatched on whether the sampling is equal-area, because a uniform 4πR²/N is exact for one family and simply wrong for the others: measured, an icosahedral geodesic's dual cells span a min/max ratio of 0.69, so a uniform default would silently corrupt every area-weighted integral on it. Non-equal-area samplings therefore get their true Voronoi dual areas, or — if the tessellation extension is not loaded — the clear error _voronoi_areas already raises, rather than a plausible wrong number.

source
FlowGeometries.Connectivity._default_node_areasMethod
_default_node_areas(::CubedSphereSampling, geometry, λ, φ)

Exact cell areas in closed form. A cubed-sphere cell is the spherical quadrilateral cut by its own panel coordinates, so its area is the spherical excess of the two triangles through its four corner directions — no tessellation, no convex hull, no optional dependency. The (n+1)² corner directions per panel are built once and shared by the four cells that meet at each, rather than re-derived per cell.

source
FlowGeometries.Connectivity._icosahedral_dual_areasMethod
_icosahedral_dual_areas(geometry, verts, triangles, nvert) -> Vector

Exact spherical-Voronoi dual-cell areas from the mesh's own triangulation — no convex hull, no optional dependency.

Each triangle is divided among its three vertices by the three arcs from its circumcenter O to its edge midpoints. Those arcs are the perpendicular bisectors of the edges — O is equidistant from all three vertices and each midpoint is equidistant from its two — so vertex a's share is the spherical quadrilateral (a, M_ab, O, M_ca), exactly its Voronoi cell restricted to that triangle. The three shares tile the triangle, so accumulating over all 20ν² triangles tiles the sphere and the areas sum to 4πR² identically. Ordering the cells' corners is never needed, so there is no per-vertex sortperm and no incident-triangle list.

source
FlowGeometries.Connectivity._sort_unique_filter!Method
_sort_unique_filter!(buf, lo, m, self, n) -> Int

Sort buf[lo+1 : lo+m] in place, drop duplicates, self-references, and out-of-range indices, and return the surviving count (left packed at the front of the slice).

source
FlowGeometries.Connectivity._yin_yang_areasMethod
_yin_yang_areas(geometry, nlon, nlat) -> Vector

Exact cell areas in closed form. A Yin–Yang cell is a lat–lon patch in its own panel frame, so it integrates to R² Δλ (sin(φ+Δφ/2) - sin(φ-Δφ/2)) = R² Δλ 2sin(Δφ/2) cos φ — independent of λ, and identical on the two panels because yang is a rigid rotation of yin.

Each panel's cells then sum to exactly √2 (3π/2) R², its full [-3π/4, 3π/4] × [-π/4, π/4] box. The two panels overlap by construction, so the areas sum to 3√2 π R² — 6.07% more than the sphere, at every resolution. That excess is the grid's real geometry, not a discretisation error: integrating over both panels needs a partition-of-unity weight for the shared region, which is a modelling choice the consumer makes on top of these areas.

source
FlowGeometries.Connectivity.adjacency_matrixMethod
adjacency_matrix(grid_or_conn; kwargs...) -> Matrix{Bool}

Dense n × n adjacency over the n nodes.

This allocates n² bytes, which is quadratic in the node count and therefore quartic in the side of a 2-D grid: a 1000×1000 grid has 10⁶ nodes and so needs ~10¹² bytes. Dense adjacency is for small node counts and for testing. For anything else use sparse_adjacency_matrix, which stores nedges entries instead of , or neighbors!, which answers neighbour queries from the index stencil with no graph storage at all.

source
FlowGeometries.Connectivity.ball_scratchMethod
ball_scratch() -> Vector{Int}

A candidate buffer to hand to repeated ball queries through their scratch argument, so an indexed query reuses one allocation instead of making one per call. One buffer per task.

source
FlowGeometries.Connectivity.build_connectivityMethod
build_connectivity(sampling, nlat; nlon, mask, periodic, stencil, active_only)

Sampling topology straight to CSR. A tensor-product sampling's neighbor graph is fixed by its axis LENGTHS and longitude wrapping alone, so the axes themselves are never evaluated — for Gauss–Legendre that is an O(n²) root solve.

source
FlowGeometries.Connectivity.build_connectivityMethod
build_connectivity(::YinYangSampling, nlon, nlat; stencil=Axial(1)) -> CSRConnectivity

Panel-local face/vertex stencils on yin then yang. Global ordering matches SphericalSampling.spherical_points! for Yin–Yang: yin (nlon×nlat, lon fastest), then yang with the same panel indexing. Overlap is not cross-linked — that is the standard Yin–Yang discrete topology (panels couple through interpolation, not shared mesh edges).

source
FlowGeometries.Connectivity.build_connectivity_withinFunction
build_connectivity_within(grid; ball, active_only=true) -> CSRConnectivity

Materialize the CSR adjacency of every pair of cells within ball of each other — the bulk form of neighbors_within, row k holding exactly what the per-cell query returns for cell k.

Symmetric by construction, since the metric is.

On the architectures with no separable axes to bound a window with — curvilinear and node grids — the default topology is indexed when NearestNeighbors is loaded, making the build O(n log n) rather than O(n²); pass topology = MetricTopology(grid) for the scanning build.

Rows are balls, i.e. Unrestricted, and there is no Connected form: reachability within one cell's ball is not a symmetric relation — a bridge cell can lie in one ball and not the other — so such a graph would not be an adjacency.

source
FlowGeometries.Connectivity.connected_componentsMethod
connected_components(grid; stencil = Stencils.Axial(1), active = true) -> (labels, ncomponents)

Label the connected components of the active region (or of the inactive region with active = false), by flood fill honouring the grid's own wrapping. labels is 0 off the region and 1:ncomponents on it.

source
FlowGeometries.Connectivity.count_holesMethod
count_holes(grid; stencil = Stencils.Axial(1)) -> Int

How many connected inactive regions are fully enclosed by active cells — the number of holes in the active region, and so an estimate of its first Betti number.

A region that reaches a non-wrapping edge is outside rather than enclosed. Along a wrapping direction there is no edge to reach, so enclosure there is decided by the fill alone.

source
FlowGeometries.Connectivity.fold_atFunction
fold_at(f, init, grid, p; ball, active_only=true, topology, scratch) -> acc

Fold acc = f(acc, J, d) over every cell J within ball of the point p, d being its distance. p is a coordinate in the grid's own coordinates, written any way a point is accepted elsewhere; fold_within is this at a cell centre.

There is no cell to exclude, so unlike the cell-seeded form every cell within ball is visited.

source
FlowGeometries.Connectivity.fold_withinFunction
fold_within(f, init, grid, I...; ball, images=NearestImage(), self=false, active_only=true)

Fold acc = f(acc, J, d) over every cell J within ball of cell I, d being its distance. The traversal every distance query here is built on; the accumulator is threaded through as a value rather than mutated, so nothing is captured and nothing is boxed.

images selects how a periodic direction is treated — NearestImage (the default, each cell once, the convention neighbors_within! exposes) or AllImages, which visits every image of a cell that lands inside the ball, each carrying its own displacement rather than a reduced one.

AllImages is what a periodic convolution needs: on a torus of period L, f̄(x) = Σₖ ∫ K(x − y − kL) f(y) dy, so where the kernel support exceeds L/2 one cell contributes through several images at different displacements, and keeping only the nearest drops the rest. Below L/2 the two conventions coincide exactly. It also widens the search — the window becomes the uncapped ceil(r/s) per periodic direction rather than metric_window's one-turn cap — and it is refused where a periodic direction is angular rather than a translation.

self = true also folds the centre cell, at distance zero. A neighbour set excludes it, which is the default; a convolution needs it, and it carries the kernel's largest weight, so omitting it is not a small error.

reach selects the ball (Unrestricted) or the part of it reachable from the seed without leaving it (Connected); topology and scratch are as in neighbors_within!.

source
FlowGeometries.Connectivity.foreach_withinMethod
foreach_within(f, grid; ball, …) -> nothing

Call f(I, J, d) for every cell I of grid and every cell J within ball of it. The same hoisting as mapreduce_within; use this one when f writes rather than reduces.

Under a threaded backend, f runs on disjoint spans of cells concurrently, so what it writes has to be determined by I — the same contract the connectivity builders keep.

source
FlowGeometries.Connectivity.gradient_planFunction
gradient_plan(grid; stencil=Stencils.Axial(1), active_only=true, conn=nothing) -> GradientPlan

Build the least-squares gradient of grid — the geometry of it, with no field involved. See Discretization.GradientPlan for what it is and why it is that; apply it with Discretization.gradient!.

This is the counterpart of apply_stencil! for the two architectures that have no separable axis to difference along: a CurvilinearGrid, whose neighbours come from its index topology, and an UnstructuredGrid, whose come from its stored adjacency. conn overrides the neighbour set; otherwise one is built, from stencil where the architecture takes one.

Surface fields only — the tangent plane is two-dimensional, so the grid's coordinates must be a (λ, φ) or (x, y) pair.

A masked cell gets no coefficients at all and reads zero gradient, and an inactive neighbour is not offered to the fit, on the same rule as everywhere else: not determined by the active data, so not invented.

source
FlowGeometries.Connectivity.healpix_neighbors!Method
healpix_neighbors!(out, nside, ipix0) -> n_written

RING-scheme topological neighbors of 0-based pixel ipix0. Writes up to 8 0-based neighbor indices into out, skipping the neighbors that do not exist at the eight singular pixels. Order: SW, W, NW, N, NE, E, SE, S.

source
FlowGeometries.Connectivity.interiorMethod
interior(grid; stencil = Stencils.Axial(1)) -> Array{Bool}

Which active cells have their whole stencil active and in range. false at a domain edge that does not wrap, and beside any masked-out cell.

source
FlowGeometries.Connectivity.is_symmetric_adjacencyMethod
is_symmetric_adjacency(conn) -> Bool

Whether j ∈ N(i) implies i ∈ N(j) throughout. O(nedges + nnodes), by comparing the graph with its transpose rather than searching a row per edge; guards the shortcut that reads a CSR as a CSC.

source
FlowGeometries.Connectivity.k_nearest!Function
k_nearest!(idx, dist, grid, I...; k, active_only=true, …) -> n

Write the k cells nearest to cell I into idx, and their distances into dist, nearest first. Returns how many were written, which is fewer than k only when the grid holds fewer candidates. The cell itself is excluded, as in neighbors_within!.

Exact under the geometry's own metric, on every architecture. It searches a ball, widens it until k cells have been seen, and keeps the k smallest in a bounded heap — so the answer never depends on the starting radius, and no candidate list is materialized. topology, scratch and reach behave as they do for neighbors_within!; an indexed topology makes each round a range query.

Ties at equal distance are broken by linear index, so the result is reproducible.

source
FlowGeometries.Connectivity.mapreduce_withinMethod
mapreduce_within(f, op, init, grid; ball, …) -> value

Reduce f(I, J, d) with op over every cell I of grid and every cell J within ball of it, d being the distance. The bulk counterpart of fold_within.

Everything that depends on the grid rather than the query is built once and reused across all n cells — above all the spatial index, which is what makes the sweep O(n log n) instead of O(n²) on a curvilinear or node grid. Writing the loop by hand gets the topology for free, since that is O(1), but not the index; measured at 9× on a 9 216-cell curvilinear grid.

op must be associative; chunks are reduced in index order, so a threaded backend gives the same answer as the serial default rather than one that depends on scheduling.

source
FlowGeometries.Connectivity.metric_bandFunction
metric_band(grid, dim, coord_t, coord_n, ball) -> T

The exact half-width along direction dim of the part of the row at coord_n that lies within ball of a point at coord_t, in that direction's own coordinate units. coord_t and coord_n are coordinates on the other direction of a two-direction grid.

metric_window returns a bounding box, which is the right answer for a query that then filters on distance. A separable sweep — a prefix sum along a row, a row-by-row convolution — cannot filter, and using the box instead of the exact extent costs it the exactness that made it worth doing. This is the same geodesic solve, resolved per row rather than maximised into a box.

Returns a negative number where the row is out of reach entirely, so band < 0 is the empty test. A row the ball covers completely gives the half-width of the whole direction (π in longitude).

On a sphere, for dim = 1, this inverts the spherical law of cosines:

\[|Δλ| ≤ \arccos\left(\frac{\cos(r/R) - \sin φ_t \sin φ_n}{\cos φ_t \cos φ_n}\right)\]

with the empty band, the whole circle and a pole at either end all falling out of the same expression — the pole case being where the denominator vanishes and the separation stops depending on λ at all, handled here once rather than by each caller.

source
FlowGeometries.Connectivity.metric_windowFunction
metric_window(grid, I, ball) -> NTuple{N,Int}

Per-direction index half-width guaranteed to contain every cell within ball of cell I.

Each direction is bounded through its smallest gap — Grids.minimum_spacing, and the seam gap where the direction wraps — which is one number on a uniform axis and an O(N) scan on a stretched one. MetricTopology holds those gaps, so the four-argument form does no scanning at all; the three-argument form builds a topology per call. On a spherical or ellipsoidal grid the longitude cut additionally walks the latitude window for its smallest cosφ, so it costs the window it returns.

This is a bound, not the answer: it is the window neighbors_within! scans before filtering on the geometry's own distance. It never under-covers, which is why it is geometry-specific. On a spherical grid one longitude step spans R·cosφ·Δλ, so the λ half-width is taken at the latitude in the window nearest a pole rather than at the cell's own latitude — at a polar cell every longitude is in range, and the window says so.

source
FlowGeometries.Connectivity.metric_windowMethod
metric_window(grid, ball) -> NTuple{N,Int}
metric_window(grid, ball, topology) -> NTuple{N,Int}

The window valid for every cell of grid, rather than for one of them: the per-cell form maximised over the grid, which is what sizing a cache or a footprint table needs.

O(1). Taking maximum of the per-cell form would be O(N) for something the cached Grids.AxisStats already determines: the smallest gap per axis is stored, and the extreme |cos φ| over a latitude axis is at one of its two ends, since |cos| on [-π/2, π/2] is largest in the middle. No cos per row, and no scan.

Conservative by construction — it is the per-cell window at the worst cell — so it never under-covers.

source
FlowGeometries.Connectivity.neighbors_within!Function
neighbors_within!(out, grid, I...; ball, active_only=true, topology, scratch, reach) -> n_written

Write the linear indices of every cell whose centre lies within ball of cell I. ball is a Stencils.MetricBall or a bare radius in the geometry's length units.

Distance is the geometry's own — great-circle on a sphere, Vincenty on a spheroid, the chord where a third direction is present — so the neighbourhood is a genuine metric ball, not a box. The cell itself is excluded, matching stencil semantics where the zero offset is not a neighbour. A periodic direction wraps, each cell appears at most once, and its coordinate is taken by minimum image, so the seam neither shortens nor lengthens a distance.

Cost is metric_windowO(1) per direction on any separable axis, uniform or stretched, given the per-direction minimum steps that topology carries — times one distance evaluation per candidate. Size the buffer with nneighbors_within; there is no fixed count, since how many cells fall within a fixed distance varies from cell to cell on any non-uniform or curved grid.

Three arguments matter for anything beyond a single query:

  • topology — a MetricTopology, the grid invariants a ball query reads. The default is O(1) and allocation-free, so leaving it out costs nothing. On a curvilinear or node grid, pass indexed to make each query O(log n + m) rather than a scan of every cell — or use foreach_within, which builds the index once for a whole sweep.
  • scratch — a candidate buffer from ball_scratch, one per task. Accepted on every grid type and used where a query goes through a spatial index, i.e. on a curvilinear or node grid; a separable window has no candidate list to buffer. With one an indexed query allocates nothing; without one it allocates its candidate list per call, 480 bytes whatever the grid size. The time difference is within noise — the allocation is the reason to pass one.
  • reachUnrestricted (the ball, and the default) or Connected (the part of it reachable from I without leaving it). Note that a ball is not a connected patch: with a mask or a concave domain it can contain cells reachable from the seed only by going outside ball.
source
FlowGeometries.Connectivity.neighbors_withinMethod
neighbors_within(grid, p; ball, …) -> Vector{Int}
nneighbors_within(grid, p; ball, …) -> Int

Cells within ball of the point p, and how many. A coordinate rather than the cell indices the other methods take, which is what distinguishes them; it may be a Tuple, NamedTuple, AbstractVector or SVector, as everywhere else a point is accepted.

source
FlowGeometries.Connectivity.sparse_adjacency_csc!Method
sparse_adjacency_csc!(colptr, rowval, conn) -> nedges

Fill caller-owned CSC structure arrays (length(colptr) ≥ nnodes+1, length(rowval) ≥ nedges) for the adjacency of conn, so that entry (i, j) is set iff j is a neighbor of i.

This is the direct route to a sparse matrix: CSR and CSC are the same layout transposed, so the structure is obtained by one counting pass and one placement pass over the existing neighbor list — no coordinate triples are materialized, nothing is sorted, and no permutation vector is built. Row indices come out ascending within each column for free, because the placement pass walks nodes in order. The running cursors live in colptr itself and are shifted back at the end, so this needs no scratch beyond the two output arrays.

source
FlowGeometries.Connectivity.sparse_adjacency_matrix!Function
sparse_adjacency_matrix!(colptr, rowval, nzval, conn) -> SparseMatrixCSC

Assemble the adjacency matrix into caller-owned buffers and wrap them without copying, so a repeated build reuses storage instead of allocating a new matrix each time. Requires using SparseArrays.

The returned matrix ALIASES the three buffers, so reusing them for a later call invalidates any matrix built from them earlier. Buffers longer than needed are trimmed to fit (nnodes+1 and nedges), since a matrix must own arrays of exactly the right length.

source
FlowGeometries.Connectivity.structured_gridMethod
structured_grid([T], sampling, nlat; geometry, nlon, mask, periodic) -> StructuredGrid

Build a spherical StructuredGrid from a tensor-product sampling (Clenshaw–Curtis, Gauss–Legendre, Driscoll–Healy, McEwen–Wiaux, lat–lon, …). Longitude periodicity is auto-detected unless periodic is set.

T is the element type to build in, and defaults to the geometry's own — a geometry fixes the width of every coordinate and metric factor computed against it. Naming T carries the geometry to that width rather than letting it promote the grid back.

source
FlowGeometries.Connectivity.unstructured_gridMethod
unstructured_grid(::AbstractScatteredSphericalSampling, λ, φ; geometry, k, radius, areas, mask)

Build an UnstructuredGrid on an arbitrary (λ, φ) point set. The points are the caller's, so there is no resolution parameter. Adjacency comes from a k-d-tree query (k nearest, or everything within a physical radius) and cell areas default to the spherical Voronoi dual; see Grids.UnstructuredGrid for the extension each needs.

source
FlowGeometries.Connectivity.unstructured_gridMethod
unstructured_grid([T], sampling, args...; geometry, areas, mask) -> UnstructuredGrid

Points from spherical_points plus exact sampling topology from build_connectivity. Default cell areas are uniform (4π R² / N on a sphere, 1 on Cartesian).

T is the element type to build in, and defaults to the geometry's own, as for structured_grid.

source