Matrix Exponentiation & Min-Plus (Tropical) Semiring Matrix Multiplication for Fixed-Length Paths

1. Formal Invariant — Semiring Associativity & Transition Operators

2. Binary Exponentiation Pipeline

3. State-Space Compilation Tricks

• Order-$K$ Linear Recurrences: If $f(n) = c_1 f(n-1) + c_2 f(n-2) + \dots + c_K f(n-K)$, the state vector is $[f(n-1), f(n-2), \dots, f(n-K)]^T$ and the companion matrix has $[c_1, c_2, \dots, c_K]$ in row $0$ and $1$s on the subdiagonal.
• Affine Constants & Polynomial Terms: If $f(n) = c_1 f(n-1) + d \cdot n + e$, augment the state vector with $[n, 1]^T$ using the binomial expansion $(n+1) = 1 \cdot n + 1 \cdot 1$.
• At Most $N$ Edges vs Exactly $N$ Edges: In the tropical $(\min, +)$ semiring, adding a zero-weight self-loop $A_{i,i} = 0$ at every vertex allows a walk to 'wait' at vertex $i$, converting 'shortest walk of exactly $N$ edges' into 'shortest walk of at most $N$ edges'.
• Multiple Queries on the Same Graph: Precompute powers of two $T^{2^0}, T^{2^1}, \dots, T^{2^{60}}$ in $O(K^3 \log N)$ time once. Then each query for walk length $N_q$ multiplies a $1 \times K$ vector by at most $\log_2 N_q$ matrices in $O(K^2 \log N_q)$ time instead of $O(K^3 \log N_q)$!

4. Hardware & Cache Optimization

5. Complexity Matrix

6. Production-Grade C++20 Implementation (Both Standard & Min-Plus Semirings)

#include <vector>
#include <cstdint>
#include <limits>
#include <algorithm>

// Tropical Min-Plus Matrix Multiplication for Shortest Walk of Exactly N Edges
struct MinPlusMatrix {
    static constexpr int64_t INF = std::numeric_limits<int64_t>::max() / 4;
    int n;
    std::vector<int64_t> mat; // Flat 1D row-major array for cache locality

    explicit MinPlusMatrix(int dim, int64_t fill_val = INF)
        : n(dim), mat(dim * dim, fill_val) {}

    static MinPlusMatrix identity(int dim) {
        MinPlusMatrix I(dim, INF);
        for (int i = 0; i < dim; ++i) I(i, i) = 0; // Tropical multiplicative identity = 0
        return I;
    }

    int64_t& operator()(int r, int c) { return mat[r * n + c]; }
    int64_t operator()(int r, int c) const { return mat[r * n + c]; }

    // Cache-friendly i-k-j loop order
    MinPlusMatrix operator*(const MinPlusMatrix& other) const {
        MinPlusMatrix res(n, INF);
        for (int i = 0; i < n; ++i) {
            for (int k = 0; k < n; ++k) {
                int64_t a_ik = (*this)(i, k);
                if (a_ik >= INF) continue;
                for (int j = 0; j < n; ++j) {
                    int64_t cand = a_ik + other(k, j);
                    if (cand < res(i, j)) {
                        res(i, j) = cand;
                    }
                }
            }
        }
        return res;
    }
};

MinPlusMatrix power(MinPlusMatrix base, int64_t exp) {
    MinPlusMatrix res = MinPlusMatrix::identity(base.n);
    while (exp > 0) {
        if (exp & 1) res = res * base;
        base = base * base;
        exp >>= 1;
    }
    return res;
}