SPDX-License-Identifier: AGPL-3.0-or-later

Math

Mathematical foundations: category theory, topological observables, and differential geometry on coupling graphs.

Category Theory

sc_neurocore.math.category_theory

Category-theoretic functors mapping stochastic, quantum, and bio domains.

CategoryObject dataclass

Bases: Generic[T]

Domain-tagged value transported between computational categories.

Source code in src/sc_neurocore/math/category_theory.py

@dataclass
class CategoryObject(Generic[T]):
    """Domain-tagged value transported between computational categories."""

    data: T
    domain: str

Morphism

Named structure-preserving map between two :class:CategoryObject domains.

Source code in src/sc_neurocore/math/category_theory.py

class Morphism:
    """Named structure-preserving map between two :class:`CategoryObject` domains."""

    def __init__(self, func: Callable[[Any], Any], name: str):
        self.func = func
        self.name = name

    def __call__(self, obj: CategoryObject[Any]) -> CategoryObject[Any]:
        """Apply the morphism, tagging the result with the morphism name."""
        return CategoryObject(data=self.func(obj.data), domain=self.name)

__call__(obj)

Apply the morphism, tagging the result with the morphism name.

Source code in src/sc_neurocore/math/category_theory.py

def __call__(self, obj: CategoryObject[Any]) -> CategoryObject[Any]:
    """Apply the morphism, tagging the result with the morphism name."""
    return CategoryObject(data=self.func(obj.data), domain=self.name)

CategoryTheoryBridge

Functors mapping between the stochastic, quantum, and bio domains.

Source code in src/sc_neurocore/math/category_theory.py

class CategoryTheoryBridge:
    """Functors mapping between the stochastic, quantum, and bio domains."""

    @staticmethod
    def stochastic_to_quantum(bitstream: np.ndarray[Any, Any]) -> np.ndarray[Any, Any]:
        """Map a bitstream probability ``p`` to the quantum amplitude pair."""
        p = np.mean(bitstream)
        # Quantum state |psi> = sqrt(p)|1> + sqrt(1-p)|0>
        alpha = np.sqrt(1 - p)
        beta = np.sqrt(p)
        return np.array([alpha, beta])

    @staticmethod
    def quantum_to_bio(state_vector: np.ndarray[Any, Any]) -> float:
        """Map quantum probability ``|beta|^2`` to a concentration in ``[0, 10]`` uM."""
        prob_1 = np.abs(state_vector[1]) ** 2
        concentration = prob_1 * 10.0
        return float(concentration)

    @staticmethod
    def bio_to_stochastic(concentration: float, length: int = 100) -> np.ndarray[Any, Any]:
        """Map a concentration to a Bernoulli bitstream of the given length."""
        p = np.clip(concentration / 10.0, 0, 1)
        rands = np.random.random(length)
        bitstream: np.ndarray[Any, Any] = (rands < p).astype(np.uint8)
        return bitstream

    def get_functor(self, source: str, target: str) -> Morphism:
        """Return the morphism mapping the ``source`` domain to ``target``."""
        if source == "Stochastic" and target == "Quantum":
            return Morphism(self.stochastic_to_quantum, "Functor: Sto->Quant")
        if source == "Quantum" and target == "Bio":
            return Morphism(self.quantum_to_bio, "Functor: Quant->Bio")
        if source == "Bio" and target == "Stochastic":
            return Morphism(self.bio_to_stochastic, "Functor: Bio->Sto")
        raise ValueError(f"No morphism from {source} to {target}")

stochastic_to_quantum(bitstream) staticmethod

Map a bitstream probability p to the quantum amplitude pair.

Source code in src/sc_neurocore/math/category_theory.py

@staticmethod
def stochastic_to_quantum(bitstream: np.ndarray[Any, Any]) -> np.ndarray[Any, Any]:
    """Map a bitstream probability ``p`` to the quantum amplitude pair."""
    p = np.mean(bitstream)
    # Quantum state |psi> = sqrt(p)|1> + sqrt(1-p)|0>
    alpha = np.sqrt(1 - p)
    beta = np.sqrt(p)
    return np.array([alpha, beta])

quantum_to_bio(state_vector) staticmethod

Map quantum probability |beta|^2 to a concentration in [0, 10] uM.

Source code in src/sc_neurocore/math/category_theory.py

@staticmethod
def quantum_to_bio(state_vector: np.ndarray[Any, Any]) -> float:
    """Map quantum probability ``|beta|^2`` to a concentration in ``[0, 10]`` uM."""
    prob_1 = np.abs(state_vector[1]) ** 2
    concentration = prob_1 * 10.0
    return float(concentration)

bio_to_stochastic(concentration, length=100) staticmethod

Map a concentration to a Bernoulli bitstream of the given length.

Source code in src/sc_neurocore/math/category_theory.py

@staticmethod
def bio_to_stochastic(concentration: float, length: int = 100) -> np.ndarray[Any, Any]:
    """Map a concentration to a Bernoulli bitstream of the given length."""
    p = np.clip(concentration / 10.0, 0, 1)
    rands = np.random.random(length)
    bitstream: np.ndarray[Any, Any] = (rands < p).astype(np.uint8)
    return bitstream

get_functor(source, target)

Return the morphism mapping the source domain to target.

Source code in src/sc_neurocore/math/category_theory.py

def get_functor(self, source: str, target: str) -> Morphism:
    """Return the morphism mapping the ``source`` domain to ``target``."""
    if source == "Stochastic" and target == "Quantum":
        return Morphism(self.stochastic_to_quantum, "Functor: Sto->Quant")
    if source == "Quantum" and target == "Bio":
        return Morphism(self.quantum_to_bio, "Functor: Quant->Bio")
    if source == "Bio" and target == "Stochastic":
        return Morphism(self.bio_to_stochastic, "Functor: Bio->Sto")
    raise ValueError(f"No morphism from {source} to {target}")

Topological Observables

Geometric invariants computed on SCPN coupling graphs: winding number from phase dynamics, Ollivier-Ricci curvature on edges, sheaf consistency defect, and connection curvature from parallel transport.

ollivier_ricci_curvature evaluates the graph-metric definition on the non-negative coupling support: lazy random-walk measures are transported by Wasserstein-1 distance over shortest-hop graph distances, and invalid matrices or node indices fail closed before transport.

Ollivier-Ricci curvature backends

The curvature solve is the one compute-bound observable in this module — each node pair runs an exact successive-shortest-path min-cost flow to obtain the Wasserstein-1 (earth-mover) distance between the two lazy random-walk measures. It therefore carries a polyglot accelerator chain selected through the backend argument ("auto", "rust", "julia", "go", "mojo", "python"). "auto" prefers the Rust path (shipped in the sc_neurocore_engine wheel) and falls back to the pure-NumPy reference. Every accelerator reproduces the NumPy reference to machine epsilon because all five implementations share the same deterministic Bellman-Ford iteration order, so the chosen augmenting paths — and the floating-point accumulation of the transport cost — coincide.

Python

from sc_neurocore.math.topology import ollivier_ricci_curvature

# auto picks the fastest available compiled backend, else pure NumPy
kappa = ollivier_ricci_curvature(coupling_matrix, i=0, j=7)

# force a specific backend (raises if that backend is not built)
kappa_rust = ollivier_ricci_curvature(coupling_matrix, 0, 7, backend="rust")

The lighter observables (winding number, sheaf defect, connection curvature) are single vectorised NumPy expressions for which NumPy is already the fastest path; they are intentionally not accelerated.

Measured performance

Reproduce with python benchmarks/bench_topology.py --json benchmarks/results/bench_topology.json. Workload: 15 curvature solves per sweep over weighted random graphs of N = 20, 50, 100; median of 7 repeats. These figures are non-isolated (loaded developer workstation, Python 3.12 / NumPy 2.3) and are functional/regression evidence, not isolated-core release numbers.

backend	median (ms)	speedup vs NumPy	parity Δ vs NumPy
python (NumPy)	3385.19	1.00×	0
mojo	66.80	50.68×	7.8e-16
rust	68.15	49.67×	7.8e-16
go	90.40	37.45×	7.8e-16
julia	93.16	36.34×	3.3e-16

Mojo and Rust are within measurement noise of each other on this loaded host; "auto" selects Rust because it ships in the wheel and needs no local build.

sc_neurocore.math.topology.winding_number(phases)

Compute the winding number of a phase trajectory around S^1.

The winding number counts how many times the phase wraps around the circle [0, 2*pi). It is a topological invariant — continuous deformations of the trajectory cannot change it.

Parameters

phases : np.ndarray, shape (T,) Time series of phase values (radians).

Returns

int Number of complete windings (positive = counterclockwise).

Source code in src/sc_neurocore/math/topology.py

def winding_number(phases: np.ndarray[Any, Any]) -> int:
    """Compute the winding number of a phase trajectory around S^1.

    The winding number counts how many times the phase wraps around
    the circle [0, 2*pi). It is a topological invariant — continuous
    deformations of the trajectory cannot change it.

    Parameters
    ----------
    phases : np.ndarray, shape (T,)
        Time series of phase values (radians).

    Returns
    -------
    int
        Number of complete windings (positive = counterclockwise).
    """
    diffs = np.diff(phases)
    # Unwrap: large jumps indicate wrapping
    diffs = np.where(diffs > np.pi, diffs - 2 * np.pi, diffs)
    diffs = np.where(diffs < -np.pi, diffs + 2 * np.pi, diffs)
    return int(np.round(np.sum(diffs) / (2 * np.pi)))

sc_neurocore.math.topology.ollivier_ricci_curvature(knm, i, j, backend='auto')

Compute Ollivier-Ricci curvature between nodes i and j on the coupling graph.

Ollivier (2009), "Ricci curvature of Markov chains on metric spaces." The curvature kappa(i,j) measures how much the neighborhoods of i and j overlap. Positive curvature = neighborhoods converge (community structure). Negative curvature = neighborhoods diverge (bottleneck).

kappa(i,j) = 1 - W1(mu_i, mu_j) / d(i,j) where mu_i is the lazy random walk distribution from node i, and W1 is the Wasserstein-1 distance on the unweighted support graph (an exact successive-shortest-path min-cost flow).

Parameters

knm : np.ndarray, shape (N, N) Coupling matrix (non-negative, not necessarily symmetric). i, j : int Node indices. backend : {"auto", "rust", "julia", "go", "mojo", "python"} Acceleration backend selector. auto prefers Rust when the sc_neurocore_engine wheel is built, else the pure-NumPy path. The named backends force a specific path and raise RuntimeError when that backend is unavailable. Every backend reproduces the NumPy reference to float64 round-off.

Returns

float Ollivier-Ricci curvature. Returns 0.0 for self or disconnected pairs.

Source code in src/sc_neurocore/math/topology.py

def ollivier_ricci_curvature(
    knm: np.ndarray[Any, Any], i: int, j: int, backend: str = "auto"
) -> float:
    """Compute Ollivier-Ricci curvature between nodes i and j on the coupling graph.

    Ollivier (2009), "Ricci curvature of Markov chains on metric spaces."
    The curvature kappa(i,j) measures how much the neighborhoods of i and j
    overlap. Positive curvature = neighborhoods converge (community structure).
    Negative curvature = neighborhoods diverge (bottleneck).

    kappa(i,j) = 1 - W1(mu_i, mu_j) / d(i,j)
    where mu_i is the lazy random walk distribution from node i,
    and W1 is the Wasserstein-1 distance on the unweighted support graph
    (an exact successive-shortest-path min-cost flow).

    Parameters
    ----------
    knm : np.ndarray, shape (N, N)
        Coupling matrix (non-negative, not necessarily symmetric).
    i, j : int
        Node indices.
    backend : {"auto", "rust", "julia", "go", "mojo", "python"}
        Acceleration backend selector. ``auto`` prefers Rust when the
        ``sc_neurocore_engine`` wheel is built, else the pure-NumPy path.
        The named backends force a specific path and raise ``RuntimeError``
        when that backend is unavailable. Every backend reproduces the
        NumPy reference to float64 round-off.

    Returns
    -------
    float
        Ollivier-Ricci curvature. Returns 0.0 for self or disconnected pairs.
    """
    if backend not in ("auto", "rust", "julia", "go", "mojo", "python"):
        raise ValueError(f"backend must be auto/rust/julia/go/mojo/python, got {backend!r}")

    graph = _validate_coupling_graph(knm)
    n_nodes = graph.shape[0]
    i = _validate_node_index("i", i, n_nodes)
    j = _validate_node_index("j", j, n_nodes)
    if i == j:
        return 0.0

    if backend == "rust" and not _HAS_RUST_TOPOLOGY:
        raise RuntimeError(
            "Rust topology backend requested but py_ollivier_ricci_curvature "
            "is not available; install the sc_neurocore_engine wheel."
        )
    if backend == "julia" and not _ensure_julia_loaded():
        raise RuntimeError(
            "Julia topology backend requested but juliacall + the "
            "accel/julia/math/topology.jl module is not available."
        )
    if backend == "go" and not _ensure_go_loaded():
        raise RuntimeError(
            "Go topology backend requested but libtopology.so is not built; "
            "run `cd src/sc_neurocore/accel/go/topology && "
            "go build -buildmode=c-shared -o libtopology.so topology.go`."
        )
    if backend == "mojo" and not _ensure_mojo_loaded():
        raise RuntimeError(
            "Mojo topology backend requested but libtopology.so is not built; "
            "run `cd src/sc_neurocore/accel/mojo/math && "
            "mojo build --emit shared-lib -o libtopology.so topology.mojo`."
        )

    if backend == "rust" or (backend == "auto" and _HAS_RUST_TOPOLOGY):
        return _ollivier_ricci_rust(graph, i, j)
    if backend == "julia":
        return _ollivier_ricci_julia(graph, i, j)
    if backend == "go":
        return _ollivier_ricci_go(graph, i, j)
    if backend == "mojo":
        return _ollivier_ricci_mojo(graph, i, j)
    return _ollivier_ricci_python(graph, i, j)

sc_neurocore.math.topology.sheaf_consistency_defect(phases, knm)

Compute the sheaf consistency defect for the SCPN phase state.

In sheaf theory, a global section exists iff the gluing conditions are satisfied on all overlaps. For the SCPN, the coupling matrix defines the overlaps, and the phase differences weighted by coupling measure the failure to glue.

defect = (1/N^2) * sum_{i,j} |K_ij| * |1 - cos(theta_i - theta_j)|

When phases are synchronized (all equal), defect = 0. When phases are maximally incoherent, defect approaches max(|K|).

This is equivalent to (1 - Kuramoto_R) weighted by coupling.

Parameters

phases : np.ndarray, shape (N,) Phase values (radians) for each layer/oscillator. knm : np.ndarray, shape (N, N) Coupling matrix.

Returns

float Sheaf consistency defect >= 0. Zero means globally coherent.

Source code in src/sc_neurocore/math/topology.py

def sheaf_consistency_defect(phases: np.ndarray[Any, Any], knm: np.ndarray[Any, Any]) -> float:
    """Compute the sheaf consistency defect for the SCPN phase state.

    In sheaf theory, a global section exists iff the gluing conditions
    are satisfied on all overlaps. For the SCPN, the coupling matrix
    defines the overlaps, and the phase differences weighted by coupling
    measure the failure to glue.

    defect = (1/N^2) * sum_{i,j} |K_ij| * |1 - cos(theta_i - theta_j)|

    When phases are synchronized (all equal), defect = 0.
    When phases are maximally incoherent, defect approaches max(|K|).

    This is equivalent to (1 - Kuramoto_R) weighted by coupling.

    Parameters
    ----------
    phases : np.ndarray, shape (N,)
        Phase values (radians) for each layer/oscillator.
    knm : np.ndarray, shape (N, N)
        Coupling matrix.

    Returns
    -------
    float
        Sheaf consistency defect >= 0. Zero means globally coherent.
    """
    N = len(phases)
    diffs = phases[np.newaxis, :] - phases[:, np.newaxis]
    cost = np.abs(knm) * (1.0 - np.cos(diffs))
    return float(cost.sum() / (N * N))

sc_neurocore.math.topology.connection_curvature(phases, knm)

Compute the connection curvature from PGBO phase dynamics.

The PGBO covariant derivative u_mu = dphi_mu - alpha * A_mu defines a U(1) connection. The curvature F_{ij} = K_{ij} * cos(theta_i - theta_j) measures the obstruction to parallel transport between layers i and j.

Parameters

phases : np.ndarray, shape (N,) Phase values. knm : np.ndarray, shape (N, N) Coupling matrix.

Returns

np.ndarray, shape (N, N) Connection curvature matrix. Diagonal is zero.

Source code in src/sc_neurocore/math/topology.py

def connection_curvature(
    phases: np.ndarray[Any, Any], knm: np.ndarray[Any, Any]
) -> np.ndarray[Any, Any]:
    """Compute the connection curvature from PGBO phase dynamics.

    The PGBO covariant derivative u_mu = dphi_mu - alpha * A_mu
    defines a U(1) connection. The curvature F_{ij} = K_{ij} * cos(theta_i - theta_j)
    measures the obstruction to parallel transport between layers i and j.

    Parameters
    ----------
    phases : np.ndarray, shape (N,)
        Phase values.
    knm : np.ndarray, shape (N, N)
        Coupling matrix.

    Returns
    -------
    np.ndarray, shape (N, N)
        Connection curvature matrix. Diagonal is zero.
    """
    diffs = phases[np.newaxis, :] - phases[:, np.newaxis]
    curvature: np.ndarray[Any, Any] = knm * np.cos(diffs)
    return curvature