用最大匹配找驱动节点，并用 Kalman 判据验证

# 把有向图拆成二分图：每个节点 u 有出端 ('o', u) 和入端 ('i', u)，边 u->v 连接 ('o',u)-('i',v)。
# 用 Hopcroft–Karp 求最大匹配后，入端未被匹配的节点就是驱动节点；若全部匹配，则任选 1 个节点。
# 再给随机有向网络的边赋随机权重，只在这些驱动节点上施加输入，检查 Kalman 秩是否等于 N。
# 解读：秩 = N 说明最大匹配给出的驱动节点集足以控制（几乎所有权重下）；N_D 是最小驱动节点数，集合本身不唯一。
import numpy as np
import networkx as nx

def driver_nodes(G):
    Bg = nx.Graph()
    top = [("o", u) for u in G]
    Bg.add_nodes_from(top)
    Bg.add_nodes_from([("i", v) for v in G])
    Bg.add_edges_from((("o", u), ("i", v)) for u, v in G.edges())
    M = nx.bipartite.hopcroft_karp_matching(Bg, top_nodes=top)
    D = [v for v in G if ("i", v) not in M]
    return D if D else [next(iter(G))]

G = nx.gnm_random_graph(12, 18, directed=True, seed=4)
D = driver_nodes(G)
print("驱动节点:", D, " N_D =", len(D))

n = G.number_of_nodes()
rng = np.random.default_rng(0)
A = np.zeros((n, n))
for u, v in G.edges():
    A[v, u] = rng.uniform(0.5, 1.5)
A /= max(1.0, np.abs(np.linalg.eigvals(A)).max())   # 缩放，避免高次幂数值溢出
B = np.zeros((n, len(D)))
for j, d in enumerate(D):
    B[d, j] = 1.0
C = np.hstack([np.linalg.matrix_power(A, k) @ B for k in range(n)])
print("Kalman 秩 =", np.linalg.matrix_rank(C), "/", n)

import matplotlib.pyplot as plt
fig, ax = plt.subplots(figsize=(6.4, 4.2))
pos = nx.kamada_kawai_layout(G.to_undirected())
nx.draw_networkx_edges(G, pos, ax=ax, edge_color="0.5", alpha=0.6, arrows=True, arrowsize=10, node_size=250)
nx.draw_networkx_nodes(G, pos, ax=ax, node_size=250, node_color=["tab:red" if v in D else "0.85" for v in G])
nx.draw_networkx_labels(G, pos, ax=ax, font_size=8)
ax.set_title("Driver nodes from maximum matching (red): N_D = %d, rank = %d / %d"
             % (len(D), np.linalg.matrix_rank(C), n), fontsize=10)
ax.axis("off")
plt.tight_layout(); plt.show()