Connectivity

Connectivity detection for multi-block structured (Plot3D) meshes.

The algorithm runs in three phases:

  1. Phase 1 – Candidate pairing: Block pairs whose axis-aligned bounding boxes (AABBs) overlap are identified by candidate_neighbor_pairs(). Only these pairs proceed to expensive point matching.

  2. Phase 2 – Face matching: For each candidate pair, find_matching_blocks() compares every outer face of block i against every outer face of block j. Each comparison (get_face_intersection()) tries three strategies in order:

    • Step 1 (same-size fast path) – if both faces have the same number of points, try all 8 permutations via _try_permutations_with_transpose().

    • Step 2 (sub-region) – for different-size faces, locate the smaller face’s diagonal corners on the larger face to extract a sub-region, then try the same 8 permutations.

    • Step 3 (fallback) – per-point geometric matching for any remaining cases.

  3. Phase 3 – Fresh-face validation: Outer faces left unmatched after Phase 2 (because a partner block’s face pool was consumed by an earlier pair) are re-checked against fresh (un-consumed) outer faces of neighbouring blocks.

Orientation system

When two faces match, the module records how their local (u, v) axes relate. Because each of the two varying axes can independently be reversed and the two axes can be swapped, there are 2 x 2 x 2 = 8 distinct orientations, encoded as 2x2 signed permutation matrices in PERMUTATION_MATRICES.

The index into that array is computed with the bit formula:

index = u_reversed | (v_reversed << 1) | (swapped << 2)

where u_reversed / v_reversed are booleans indicating whether the traversal direction along that axis is flipped between the two faces, and swapped indicates whether the u and v axes of face 1 map to v and u of face 2 (a cross-plane match). Indices 0–3 are the four direct (non-swapped) permutations; indices 4–7 are the four transposed (swapped) permutations.

JSON export convention

The exported permutation_index follows the diagonal (lb/ub) convention:

  • In-plane matches (perm 0–3): permutation_index is -1 because the traversal direction is fully encoded in block2’s lb/ub ordering.

  • Cross-plane matches (perm 4–7): permutation_index is the actual index because lb/ub alone cannot represent an axis swap.

The permutation_matrix field always contains the actual 2x2 matrix.

plot3d.connectivity.CELL_CORNER_OFFSETS: List[Tuple[int, int, int]] = [(1, 0, 0), (0, 1, 0), (0, 0, 1), (1, 1, 0), (1, -1, 0), (1, 0, 1), (1, 0, -1), (0, 1, 1), (0, 1, -1), (1, 1, 1), (1, 1, -1), (1, -1, 1), (1, -1, -1)]

The 13 index offsets that enumerate every unordered pair of corners of a hexahedral cell exactly once: 3 edges, 6 face diagonals, 4 body diagonals. Each pair of the 8 corners differs by some (di, dj, dk) in {-1, 0, 1}^3 \ {0}; taking the sign in which the first non-zero component is positive picks one of the two orderings.

plot3d.connectivity.TOL_FLOOR = 1e-06

Absolute floor for the matching tolerance – the value used unconditionally before the tolerance became adaptive. adaptive_tolerance() never returns anything smaller, so no interface detected before this change can be lost by it.

plot3d.connectivity.adaptive_tolerance(blocks: List[Block]) → float[source]

Derive the node-matching tolerance for blocks from the mesh itself.

Mirrors plot3d-rs connectivity::adaptive_tolerance (commit c7e9cf6, src/connectivity.rs).

A fixed absolute tolerance silently encodes an assumption that coordinates are of order 1. Coordinate storage loses precision in proportion to magnitude, so two copies of the same physical interface node – one written by each of the two blocks that share it – can differ by an amount that grows with the coordinate magnitude, and a fixed 1e-6 stops bridging that gap well before “exotic” magnitudes.

Formula:

scale = max |coordinate|
noise = _TOL_RELATIVE_NOISE * scale
h     = min non-zero distance between two corners of one cell
        (edges, face diagonals and body diagonals)
tol   = clamp(noise, TOL_FLOOR, max(TOL_FLOOR, _TOL_SPACING_FRACTION * h))

Floor (TOL_FLOOR, 1e-6): the result is never below the historical fixed constant, so this cannot remove a match that is found today. Meshes with coordinates of magnitude <= 1 get exactly 1e-6, bit-for-bit, and skip the (more expensive) spacing pass entirely.

Ceiling (a quarter of the finest cell corner spacing): bounds the false-positive risk of a node finding the wrong partner once storage noise pushes its true partner just outside tolerance. h is the smallest distance between any two corners of one cell (not just edges, since a sheared cell’s short face diagonal can be closer than an edge), so bounding the tolerance at 0.25 * h keeps every wrong corner at least three tolerances from the true partner.

plot3d.connectivity.candidate_neighbor_pairs(blocks: List[Block], tol: float = 1e-06)[source]

Returns candidate block pairs whose AABBs overlap or nearly touch.

This replaces the former centroid-distance approach which only considered the N nearest blocks and could miss neighbours for L-shaped or elongated geometries. AABB overlap is both more robust and more correct.

Parameters:

  • blocks (List[Block]) – list of all your blocks

  • tol (float) – AABB expansion tolerance

Returns:

candidate (i, j) pairs with i < j

Return type:

List[Tuple[int,int]]

plot3d.connectivity.combinations_of_nearest_blocks(blocks: List[Block], nearest_nblocks: int = 4)[source]

Returns the indices of the nearest N blocks based on their centroid.

Deprecated since version Use: candidate_neighbor_pairs() instead for AABB-based pairing.

Parameters:

  • blocks (List[Block]) – list of all your blocks

  • nearest_nblocks (int) – number of nearest blocks to consider

Returns:

combinations of nearest blocks

Return type:

List[Tuple[int,int]]

plot3d.connectivity.connectivity(blocks: List[Block], tol: float | None = None)[source]

Returns a dictionary outlining the connectivity of the blocks along with any exterior surfaces.

Each face match dict includes an orientation sub-dict with:

  • permutation_index: -1 for in-plane matches (direction encoded in lb/ub), or the actual index (4-7) for cross-plane matches.

  • plane: 'in-plane' or 'cross-plane'.

  • permutation_matrix: the actual 2x2 signed permutation matrix.

Parameters:

  • blocks (List[Block]) – List of all blocks in multi-block plot3d mesh

  • tol (float, Optional) – Euclidean node-matching tolerance. Defaults to None, which derives the tolerance from blocks via adaptive_tolerance(). Passing TOL_FLOOR reproduces the fixed-tolerance behaviour this module had before the tolerance became adaptive.

Returns:

All matching faces formatted as a list of { ‘block1’: {‘block_index’, ‘lb’, ‘ub’} } (List[Dict]): All exterior surfaces formatted as a list of { ‘block_index’, ‘lb’, ‘ub’, ‘id’ }

Return type:

(List[Dict])

plot3d.connectivity.connectivity_fast(blocks: List[Block], use_minmax: bool = False, tol: float | None = None)[source]

Find connectivity by GCD-reducing blocks first for speed.

Computes the minimum GCD across all block dimensions, reduces all blocks uniformly, runs connectivity(), then scales bounds back up.

Face match dicts follow the diagonal (lb/ub) convention: permutation_index is -1 for in-plane, actual index for cross-plane.

Parameters:

  • blocks (List[Block]) – List of blocks to find connectivity for.

  • use_minmax (bool) – If True, normalise lb/ub to strict min/max order (IMIN,JMIN,KMIN → IMAX,JMAX,KMAX) and recompute the permutation matrix accordingly. Default is False (traversal order).

  • tol (float, Optional) – Euclidean node-matching tolerance. Defaults to None, which derives the tolerance from blocks via adaptive_tolerance(). The tolerance is always derived from the full-resolution blocks passed in here, before GCD reduction – reduction leaves storage noise unchanged but multiplies cell size by the GCD, so deriving from the reduced grid would give a ceiling up to gcd times looser. This matches plot3d-rs’s connectivity_fast invariant that it and connectivity() use one tolerance for a given mesh.

Returns:

Face matches with orientation info. (List[Dict]): Outer (non-connected) faces.

Return type:

(List[Dict])

plot3d.connectivity.demote_to_outer(outer_faces: List[dict], rejected: List[dict]) → None[source]

Append both faces of each rejected proposal to outer_faces IN PLACE.

Mirrors plot3d-rs connectivity::demote_to_outer (commit 0e1b1a1). Deduplicates by (block_index, lb, ub) against faces already present in outer_faces and against faces from other rejected proposals in this same batch. Newly-appended faces get fresh sequential id values continuing from the highest id already present in outer_faces, following the same convention connectivity() uses when it first builds outer_faces_formatted (sequential id starting at 1).

plot3d.connectivity.face_matches_to_dict(face1: Face, face2: Face, block1: Block, block2: Block)[source]

Makes sure the diagonal of face 1 match the diagonal of face 2

Parameters:

  • face1 (Face) – Face 1 with block index

  • face2 (Face) – Face 2 with block index

  • block1 (Block) – Block 1

  • block2 (Block) – Block 2

Returns:

dictionary describing the corner matches

Return type:

(dict)

plot3d.connectivity.find_matching_blocks(block1: Block, block2: Block, block1_outer: List[Face], block2_outer: List[Face], tol: float = 1e-06)[source]

Takes two blocks and finds all matching pairs

Parameters:

  • block1 (Block) – Any plot3d Block that is not the same as block2

  • block2 (Block) – Any plot3d Block that is not the same as block1

  • block1_outer (List[Face]) – outer faces for block 1.

  • block2_outer (List[Face]) – Outer faces for block 2

  • tol (float, Optional) – tolerance to use. Defaults to 1E-6

Note

This function was changed to be given an input of outer faces for block 1 and block 2. Outer faces can change and we should use the updated value

Returns:

containing

  • df (pandas.DataFrame): corners of matching pair as block1_corners,block2_corners ([imin,jmin,kmin],[imax,jmax,kmax]), ([imin,jmin,kmin],[imax,jmax,kmax])

  • block1_outer (List[Face]):

  • block2_outer (List[Face]):

Return type:

(tuple)

plot3d.connectivity.get_face_intersection(face1: Face, face2: Face, block1: Block, block2: Block, tol: float = 1e-06)[source]

Get the index of the intersection between two faces located on two different blocks.

Three-step approach:

Step 1 (fast): Same-size faces — try 8 permutations (4 direction + 4 transposed). Step 2 (subregion): Different-size faces — find smaller face’s corners on larger face to identify subregion, then try 8 permutations on subregion. Step 3 (fallback): Per-point geometric matching for any remaining cases.

Parameters:

  • face1 (Face) – An exterior face

  • face2 (Face) – An exterior face from a different block

  • block1 (Block) – block containing face1

  • block2 (Block) – block containing face2

  • tol (float) – matching tolerance

Returns:

containing

  • (pandas.DataFrame): dataframe with matches. Columns = i1, j1, k1, i2, j2, k2

  • (List[Face]): any split faces from block 1

  • (List[Face]): any split faces from block 2

Return type:

(Tuple)

plot3d.connectivity.normalize_face_matches(face_matches: list) → list[source]

Convert face match lb/ub to strict min/max order.

By default the library encodes traversal direction in lb/ub ordering (lb is not necessarily < ub). This function normalises every face match so that lb = [IMIN, JMIN, KMIN] and ub = [IMAX, JMAX, KMAX] on both sides, and recomputes the orientation / permutation matrix accordingly.

The resulting PM satisfies:

Face_A(IMIN,JMIN,KMIN -> IMAX,JMAX,KMAX) * PM
    = Face_B(IMIN,JMIN,KMIN -> IMAX,JMAX,KMAX)
Parameters:

face_matches (list of dict) – Face match dicts as returned by connectivity_fast() or connectivity(). Each dict must have block1 and block2 sub-dicts with lb and ub keys. If a match key (point-match DataFrame) is present it is used to recompute the orientation; otherwise the orientation is recomputed from the new min/max bounds.

Returns:

A new list with the same structure, but lb/ub in min/max order and orientation updated.

Return type:

list of dict

plot3d.connectivity.revalidate_full_resolution(blocks: List[Block], proposed: List[dict], transforms: List[Callable[[ndarray], ndarray]], tol: float, stage: str) → Tuple[List[dict], List[dict]][source]

Re-certify every reduced-grid face-match proposal against the ORIGINAL full-resolution blocks.

Mirrors plot3d-rs connectivity::revalidate_full_resolution (commit 0e1b1a1). GCD reduction can occasionally produce a proposal that looks valid on the coarse grid but does not actually hold node-for-node at full resolution (e.g. an interior node perturbed independently of its neighbours during mesh generation, invisible at reduced resolution). blocks must be the untouched, full-resolution mesh – never the GCD-reduced copy – since the whole point is to check against data the reduced grid could not see.

Each transform in transforms is tried in order (first one that certifies wins); a proposal that carries a declared orientation (a populated orientation.permutation_matrix) must certify AS DECLARED via correspondence.certify_permutation() – no search, so a wrong declared orientation is a definite failure, never silently replaced by a different one found via search. An undeclared proposal (no orientation key, e.g. a self-match record) is searched via correspondence.certify_correspondence().

A single aggregated RuntimeWarning is raised (not one per rejected proposal) naming the count of demotions and the worst-case discrepancy among them, pulled from whichever correspondence.MappingFailure subclass each rejection raised.

Returns:

the proposals that certified, and the ones that did not (in original relative order within proposed).

Return type:

(kept, rejected)

plot3d.connectivity.select_multi_dimensional(T: ndarray, dim1: tuple, dim2: tuple, dim3: tuple)[source]
Takes a block (T) and selects X,Y,Z from the block given a face’s dimensions

theres really no good way to do this in python

Parameters:

  • T (np.ndarray) – arbitrary array so say a full matrix containing X

  • dim1 (tuple) – 20,50 this selects X in the i direction from i=20 to 50

  • dim2 (tuple) – 40,60 this selects X in the j direction from j=40 to 60

  • dim3 (tuple) – 10,20 this selects X in the k direction from k=10 to 20

Returns:

returns X or Y or Z given some range of I,J,K

Return type:

np.ndarray