spn_reachprops.h

Go to the documentation of this file.

/*
 * Copyright (c) 2012-2026, QORE Lab, Imperial College London
 * All rights reserved.
 */
#ifndef LINE_API_SPN_SPN_REACHPROPS_H
#define LINE_API_SPN_SPN_REACHPROPS_H
 
/**
 * @file
 * @ingroup api_spn
 * Behavioural properties of a stochastic Petri net from its reachability
 * graph (Murata 1989, sec. IV; Haverkort 1998, sec. 14.2), under the GSPN
 * firing rule of spn_reachgraph, so inhibitor arcs, priorities and immediate
 * events are honoured. Verdicts are tri-state: 1 true, 0 false, -1 unknown.
 * On a complete graph:
 *   bounded = 1 and bound(l) is the largest token count of level l;
 *   deadlocks are the markings where nothing may fire;
 *   live(e) = 1 (L4-live) iff e labels an edge leaving a marking of every
 *     bottom strongly connected component;
 *   dead(e) = 1 iff e never fires;
 *   reversible = 1 iff there is one bottom SCC and it holds m0.
 * On an unbounded or capped graph, bounded is 0 or -1, bound is NaN (Inf on
 * the proven unbounded levels), a found deadlock makes deadlockfree, every
 * live(e) and (when it is not m0) reversible 0, and the other verdicts are
 * -1. MATLAB twin: matlab/src/api/spn/spn_reachprops.m.
 */
 
#include <algorithm>
#include <cstddef>
#include <limits>
#include <string>
#include <vector>
 
#include "line/api/spn/spn_reachgraph.h"
 
namespace line {
namespace spn {
 
struct SpnReachProps {
    std::string status;                       ///< the graph status
    std::size_t nstates = 0;
    int bounded = -1;
    std::vector<double> bound;                ///< per level; NaN unknown, Inf unbounded
    int deadlockfree = -1;
    std::vector<std::size_t> deadlocks;       ///< 0-based marking indices
    std::vector<std::vector<long long>> deadmarkings;
    std::vector<int> live, dead;              ///< per event, tri-state
    int reversible = -1;
    int nbscc = -1;                           ///< number of bottom SCCs (-1 unknown)
    SpnReachGraph graph;
};
 
namespace detail {
 
/** Strongly connected components (iterative Tarjan); comp[i] is a 0-based label. */
inline std::vector<std::size_t> spn_scc(std::size_t ns, const std::vector<std::vector<std::size_t>>& Ed,
                                        std::size_t& ncomp) {
    std::vector<std::vector<std::size_t>> adj(ns);
    for (std::size_t r = 0; r < Ed.size(); ++r) adj[Ed[r][0]].push_back(Ed[r][1]);
    const std::size_t none = std::numeric_limits<std::size_t>::max();
    std::vector<std::size_t> index(ns, none), low(ns, 0), comp(ns, none), stack, it(ns, 0);
    std::vector<bool> on(ns, false);
    std::size_t idx = 0;
    ncomp = 0;
    for (std::size_t s = 0; s < ns; ++s) {
        if (index[s] != none) continue;
        std::vector<std::size_t> call(1, s);
        index[s] = low[s] = idx++;
        stack.push_back(s);
        on[s] = true;
        while (!call.empty()) {
            const std::size_t v = call.back();
            if (it[v] < adj[v].size()) {
                const std::size_t w = adj[v][it[v]++];
                if (index[w] == none) {
                    index[w] = low[w] = idx++;
                    stack.push_back(w);
                    on[w] = true;
                    call.push_back(w);
                } else if (on[w]) {
                    low[v] = std::min(low[v], index[w]);
                }
                continue;
            }
            call.pop_back();
            if (!call.empty()) low[call.back()] = std::min(low[call.back()], low[v]);
            if (low[v] == index[v]) {
                for (;;) {
                    const std::size_t w = stack.back();
                    stack.pop_back();
                    on[w] = false;
                    comp[w] = ncomp;
                    if (w == v) break;
                }
                ++ncomp;
            }
        }
    }
    return comp;
}
 
}  // namespace detail
 
template <class T>
SpnReachProps spn_reachprops(const qn::NetworkStruct<T>& sn,
                             const std::vector<double>& init = std::vector<double>(),
                             std::size_t maxstates = 0) {
    SpnReachProps rp;
    rp.graph = spn_reachgraph(sn, init, maxstates);
    const SpnReachGraph& g = rp.graph;
    const detail::SpnNet net = detail::spn_build_net(sn, "spn_reachprops");
    const std::size_t ns = g.M.size(), E = net.E, n = net.n;
    rp.status = g.status;
    rp.nstates = ns;
    for (std::size_t i = 0; i < ns; ++i) {
        const std::vector<bool> f = net.fireable(g.M[i]);
        if (std::find(f.begin(), f.end(), true) == f.end()) {
            rp.deadlocks.push_back(i);
            rp.deadmarkings.push_back(g.M[i]);
        }
    }
    std::vector<bool> fired(E, false);
    for (std::size_t r = 0; r < g.E.size(); ++r) fired[g.E[r][2]] = true;
    if (g.status == "complete") {
        rp.bounded = 1;
        rp.bound.assign(n, 0.0);
        for (std::size_t i = 0; i < ns; ++i)
            for (std::size_t l = 0; l < n; ++l)
                rp.bound[l] = std::max(rp.bound[l], static_cast<double>(g.M[i][l]));
        rp.deadlockfree = rp.deadlocks.empty() ? 1 : 0;
        std::size_t nc = 0;
        const std::vector<std::size_t> comp = detail::spn_scc(ns, g.E, nc);
        std::vector<bool> bottom(nc, true);
        for (std::size_t r = 0; r < g.E.size(); ++r)
            if (comp[g.E[r][0]] != comp[g.E[r][1]]) bottom[comp[g.E[r][0]]] = false;
        rp.live.assign(E, 1);
        std::size_t nb = 0;
        for (std::size_t b = 0; b < nc; ++b) {
            if (!bottom[b]) continue;
            ++nb;
            std::vector<bool> lab(E, false);
            for (std::size_t r = 0; r < g.E.size(); ++r)
                if (comp[g.E[r][0]] == b) lab[g.E[r][2]] = true;
            for (std::size_t e = 0; e < E; ++e)
                if (!lab[e]) rp.live[e] = 0;
        }
        rp.dead.assign(E, 0);
        for (std::size_t e = 0; e < E; ++e) rp.dead[e] = fired[e] ? 0 : 1;
        rp.reversible = (nb == 1 && bottom[comp[0]]) ? 1 : 0;
        rp.nbscc = static_cast<int>(nb);
    } else {
        rp.bounded = g.status == "unbounded" ? 0 : -1;
        rp.bound.assign(n, std::numeric_limits<double>::quiet_NaN());
        for (std::size_t l = 0; l < n; ++l)
            if (g.unboundedlevels[l]) rp.bound[l] = std::numeric_limits<double>::infinity();
        if (rp.deadlocks.empty()) {
            rp.deadlockfree = -1;
            rp.live.assign(E, -1);
            rp.reversible = -1;
        } else {
            rp.deadlockfree = 0;
            rp.live.assign(E, 0);
            bool other = false;
            for (std::size_t i = 0; i < rp.deadlocks.size(); ++i) other = other || rp.deadlocks[i] != 0;
            rp.reversible = other ? 0 : -1;
        }
        rp.dead.assign(E, -1);
        for (std::size_t e = 0; e < E; ++e)
            if (fired[e]) rp.dead[e] = 0;
        rp.nbscc = -1;
    }
    return rp;
}
 
}  // namespace spn
}  // namespace line
 
#endif  // LINE_API_SPN_SPN_REACHPROPS_H