5D Library › 复杂系统相变、临界慢化与生态韧性的非线性动力学 · English

第五章 级联失效:冗余容量与相互依赖

一个部件,一大片停电

2003 年 8 月 14 日,北美东北部与加拿大安大略省大面积停电,五千多万人受到影响,约 61,800 兆瓦的负荷中断。美加联合调查组的最终报告把原因归为四类:对电力系统缺乏足够的理解、缺乏足够的态势感知、树木修剪不足,以及可靠性协调机构的诊断支持不足 (U.S.–Canada Power System Outage Task Force 2004)。起点是俄亥俄州的几条输电线路:一条线路断开后,电力转移到其他线路上,使它们过载、跳闸,负荷再次转移。这种「失效—负荷转移—再失效」的过程叫做级联失效。它与折叠分岔一样,系统长期看起来正常,然后在很短的时间里大面积失效;不同之处在于,它由许多部件之间的连接造成。

容量—负荷模型

证据

莫特与赖 2002。 网络中每一对节点之间沿最短路径传送一个单位的「流」,节点的负荷 LiL_i 是经过它的最短路径数(介数);容量 Ci=(1+α)LiC_i = (1 + \alpha) L_i,α≥0\alpha \ge 0 是冗余参数。移除一个节点后负荷重新分配,过载的节点失效,负荷再次分配,直到不再有节点过载;损害用级联后最大连通分量与原来之比 GG 衡量。在约 5,000 个节点的无标度网络中,随机移除一个节点几乎不引起级联;移除负荷最高的节点则损害很大:α=0.2\alpha = 0.2 时超过 60% 的节点受到影响,即使 α=1\alpha = 1,GG 仍减少 20% 以上。所有节点度都为 3 的均匀网络中,即使 α=0.05\alpha = 0.05 也不引起级联。在美国西部电网的拓扑中,度分布接近指数分布但负荷分布偏斜,按负荷攻击一个节点,即使 α=1\alpha = 1,最大连通分量也缩小到不到一半。作者的结论是:脆弱性来自负荷的异质性 (Motter and Lai 2002)。

推导

一个 300 个节点的无标度网络(巴拉巴西—阿尔伯特模型,每个新节点连 2 个老节点),与一个节点数相同、连接数几乎相同(600 条对 597 条)、所有节点度都为 4 的均匀网络。负荷用布兰德斯算法计算。对一系列 α\alpha,分别移除负荷最高的节点与随机 5 个节点之一(取平均),模拟级联。网络规模是示意的,远小于莫特与赖的网络。

推导

# Overload cascades (Motter & Lai 2002). Load L_i = number of shortest paths between other
# pairs that pass through i (Brandes's algorithm); capacity C_i = (1 + alpha) L_i0. Remove one
# node, recompute loads, remove every overloaded node, repeat. G = largest component / N.
# Scale-free (Barabasi-Albert, N = 300, m = 2, 597 links) vs homogeneous (all degrees 4, 600 links).
import random
from collections import deque
rng = random.Random(10)
N = 300

def ba(m=2):
    adj = {i: set() for i in range(N)}; targets = []
    for i in range(m + 1):
        for j in range(i): adj[i].add(j); adj[j].add(i); targets += [i, j]
    for v in range(m + 1, N):
        chosen = set()
        while len(chosen) < m: chosen.add(rng.choice(targets))
        for t in chosen: adj[v].add(t); adj[t].add(v); targets += [v, t]
    return adj

def regular(k=4):
    while True:                                   # configuration model, retry until simple
        stubs = [v for v in range(N) for _ in range(k)]; rng.shuffle(stubs)
        adj = {i: set() for i in range(N)}; ok = True
        for a, b in zip(stubs[::2], stubs[1::2]):
            if a == b or b in adj[a]: ok = False; break
            adj[a].add(b); adj[b].add(a)
        if ok: return adj

def loads(adj, alive):
    L = dict.fromkeys(alive, 0.0)
    for s in alive:                               # Brandes (2001), unweighted
        pred = {v: [] for v in alive}; sig = dict.fromkeys(alive, 0); sig[s] = 1
        dist, order, q = {s: 0}, [], deque([s])
        while q:
            v = q.popleft(); order.append(v)
            for w in adj[v]:
                if w not in alive: continue
                if w not in dist: dist[w] = dist[v] + 1; q.append(w)
                if dist[w] == dist[v] + 1: sig[w] += sig[v]; pred[w].append(v)
        delta = dict.fromkeys(alive, 0.0)
        for w in reversed(order):
            for v in pred[w]: delta[v] += sig[v] / sig[w] * (1 + delta[w])
            if w != s: L[w] += delta[w]
    return {v: x / 2 for v, x in L.items()}       # each unordered pair counted twice

def giant(adj, alive):
    seen, best = set(), 0
    for s in alive:
        if s in seen: continue
        seen.add(s); st, n = [s], 0
        while st:
            v = st.pop(); n += 1
            for w in adj[v]:
                if w in alive and w not in seen: seen.add(w); st.append(w)
        best = max(best, n)
    return best

def cascade(adj, L0, alpha, first):
    alive, rounds = set(adj) - {first}, 0
    while True:
        L = loads(adj, alive)
        over = {v for v in alive if L[v] > (1 + alpha) * L0[v] + 1e-9}
        if not over: return giant(adj, alive) / N, rounds
        alive -= over; rounds += 1

if __name__ == "__main__":
    sf, h = ba(), regular()
    L0, Lh = loads(sf, set(sf)), loads(h, set(h))
    top, toph = max(L0, key=L0.get), max(Lh, key=Lh.get)
    print(f"scale-free: max degree {max(len(v) for v in sf.values())}, top load {L0[top]:.0f}"
          f" (mean {sum(L0.values()) / N:.0f}); homogeneous: top load {Lh[toph]:.0f}")
    randoms, res = rng.sample(range(N), 5), {}
    print(" alpha   G after removing: top-load node (rounds)  a random node (mean of 5)"
          "  homogeneous, top-load")
    for alpha in (0.0, 0.1, 0.2, 0.3, 0.5, 0.8, 1.2, 1.4, 2.0):
        g1, r1 = cascade(sf, L0, alpha, top); res[alpha] = g1
        g2 = sum(cascade(sf, L0, alpha, v)[0] for v in randoms) / 5
        print(f"  {alpha:3.1f}            {g1:5.2f} ({r1})                      {g2:5.2f}"
              f"                   {cascade(h, Lh, alpha, toph)[0]:5.2f}")
    grid = [i / 10 for i in range(31)]                # alpha_c on a 0.1 grid up to 3.0
    gs = {a: res[a] if a in res else cascade(sf, L0, a, top)[0] for a in grid}
    ac = min(a for a in grid if all(gs[b] >= 0.95 for b in grid if b >= a))
    print(f"alpha_c (G >= 0.95 at every alpha from here to 3.0, step 0.1): {ac}"
          f"  (G = {gs[ac]:.2f}; at {ac - 0.1:.1f}: {gs[round(ac - 0.1, 1)]:.2f})")
# scale-free: max degree 28, top load 7510 (mean 395); homogeneous: top load 720
#  alpha   G after removing: top-load node (rounds)  a random node (mean of 5)  homogeneous, top-load
#   0.0             0.02 (1)                       0.46                    0.02
#   0.1             0.17 (3)                       0.99                    1.00
#   0.2             0.44 (4)                       1.00                    1.00
#   0.3             0.63 (3)                       1.00                    1.00
#   0.5             0.79 (2)                       1.00                    1.00
#   0.8             0.86 (2)                       1.00                    1.00
#   1.2             0.90 (3)                       1.00                    1.00
#   1.4             0.96 (2)                       1.00                    1.00
#   2.0             0.97 (2)                       1.00                    1.00
# alpha_c (G >= 0.95 at every alpha from here to 3.0, step 0.1): 1.4  (G = 0.96; at 1.3: 0.91)
左:容量—负荷模型中级联后的最大连通分量随冗余容量 \alpha 的变化(300 个节点)。右:两个相互依赖的随机网络(平均度 4)中仍在运转的节点比例随保留比例 p 的变化:单一网络在 p = 0.25 处连续地消失;一对一完全依赖时(q = 1)在 p \approx 0.614 处突然从约三分之一跳到零;只有 20% 的节点有依赖关系时(q = 0.2),崩溃又变为连续的。

无标度网络的负荷极不均匀:负荷最高的节点承担 7,510 条最短路径,平均每个节点 395 条;均匀网络中最高只有 720。α=0\alpha = 0 时任何变化都会使某些节点过载,随机移除一个节点 GG 也降到 0.46。只要有 10% 的冗余,随机移除与均匀网络中的定向移除就几乎不造成损害(0.99 与 1.00);但移除无标度网络中负荷最高的节点,GG 只剩 0.17,α=0.5\alpha = 0.5 时为 0.79,α=1.2\alpha = 1.2 时为 0.90。以「从这个 α\alpha 起直到 3.0,GG 都不低于 0.95」为操作性定义,按 0.1 的步长扫描,这个网络的临界冗余 αc=1.4\alpha_c = 1.4(G=0.96G = 0.96;α=1.3\alpha = 1.3 时为 0.91):每个节点都要有比日常负荷多约 140% 的容量,才能挡住对一个关键节点的攻击。这是一个网络、一个被攻击节点的结果,不是普适常数;换一个随机生成的网络,这个数字会变。

冗余的多少之外,冗余的分配同样重要。均匀的 10% 余量足以应付随机故障;关键节点的失效需要大得多的余量,而这些余量平时几乎永远用不上,在只看日常效率的预算中最容易被砍掉。第四章所说的「抽干冗余」,在网络中就是把 α\alpha 压到恰好够应付随机故障的水平。

注

模型的限度。 模型假设负荷沿最短路径分布、失效节点永久移除、没有人干预;真实电网的潮流由物理定律决定,操作员会切负荷、解列。模型说明的是负荷异质性与冗余的关系,不预测任何真实电网的失效。同样的逻辑常被用作金融传染的类比;本书不就任何金融系统作判断,也不构成投资建议。

相互依赖的网络

证据

布尔德列夫等人 2010。 他们以 2003 年 9 月 28 日意大利大停电为例:发电站停机导致互联网通信节点失效,后者又使更多发电站失效。模型的假设是:两个网络 A 与 B 节点数相同,各自内部随机连接、彼此不相关;A 中每个节点恰好依赖 B 中的一个节点,反之亦然(一对一、双向);一个节点只有在属于本网络的最大连通分量、并且它依赖的节点也在运转时才能运转。随机移除 A 中比例为 1−p1 – p 的节点后,失效在两个网络之间来回传递。对平均度同为 ⟨k⟩\langle k\rangle 的两个随机网络,相互连通的巨分量存在的临界值为 pc=2.4554/⟨k⟩p_c = 2.4554/\langle k\rangle,在 pcp_c 处巨分量突然跳到零(一阶相变)。对单个网络而言,度分布越宽越能抵抗随机失效;对相互依赖的网络则相反 (Buldyrev et al. 2010)。

推导

两个平均度为 4 的随机网络,各 4,000 个节点。理论:单个网络保留比例 pp 的节点后,巨分量满足 S=p(1−e−⟨k⟩S)S = p\,(1 – e^{-\langle k\rangle S});一对一相互依赖时 P=p(1−f)2P = p\,(1 – f)^2,f=exp⁡[−⟨k⟩p(1−f)2]f = \exp[-\langle k\rangle p (1 – f)^2]。模拟:只有比例为 qq 的节点有依赖关系(帕尔沙尼等人的「部分依赖」),迭代失效直到稳定,报告 A 中仍在运转的比例。

推导

# Interdependent networks (Buldyrev et al. 2010): two random networks A and B, mean degree 4.
# Assumptions: a fraction 1-p of A fails at random; a node works only if it is in the giant
# component of its own network AND its partner in the other network works. One-to-one
# dependence: P = p (1-f)^2, f = exp(-k p (1-f)^2), jump at p_c = 2.4554/k; a single network:
# S = p (1 - exp(-k S)), continuous at 1/k. Simulation: N = 4,000; a share q of nodes has a
# partner (partial dependence, Parshani et al. 2010).
import math, random
rng = random.Random(12)
K, N = 4.0, 4000

def theory_pair(p):
    f = 0.0
    for _ in range(3000): f = math.exp(-K * p * (1 - f) ** 2)
    return p * (1 - f) ** 2

def theory_single(p):
    s = 1.0
    for _ in range(3000): s = p * (1 - math.exp(-K * s))
    return s

def er():
    adj = [[] for _ in range(N)]
    for _ in range(int(K * N / 2)):
        a, b = rng.randrange(N), rng.randrange(N)
        if a != b: adj[a].append(b); adj[b].append(a)
    return adj

def giant(adj, alive):                       # alive flags of the largest component
    seen, best = [False] * N, []
    for s in range(N):
        if not alive[s] or seen[s]: continue
        comp, st = [s], [s]; seen[s] = True
        while st:
            for w in adj[st.pop()]:
                if alive[w] and not seen[w]: seen[w] = True; st.append(w); comp.append(w)
        if len(comp) > len(best): best = comp
    out = [False] * N
    for v in best: out[v] = True
    return out

def mutual(A, B, partner, p):
    aliveA, aliveB = [rng.random() < p for _ in range(N)], [True] * N
    while True:
        aliveA = giant(A, aliveA)
        for v in range(N):
            if partner[v] is not None and not aliveA[v]: aliveB[partner[v]] = False
        aliveB = giant(B, aliveB); changed = False
        for v in range(N):
            if partner[v] is not None and aliveA[v] and not aliveB[partner[v]]:
                aliveA[v] = False; changed = True
        if not changed: return sum(aliveA) / N

if __name__ == "__main__":
    print(f"theory: single network p_c = {1 / K:.3f}; coupled pair p_c = {2.4554 / K:.4f}")
    A, B = er(), er()
    perm = list(range(N)); rng.shuffle(perm)
    print("   p    theory: single  pair    simulated: q = 1   q = 0.5   q = 0.2")
    for p in (0.90, 0.80, 0.70, 0.65, 0.62, 0.60, 0.50, 0.40, 0.30):
        sims = []
        for q in (1.0, 0.5, 0.2):
            partner = [perm[v] if rng.random() < q else None for v in range(N)]
            sims.append(mutual(A, B, partner, p))
        print(f"  {p:.2f}       {theory_single(p):.3f}     {theory_pair(p):.3f}"
              f"            {sims[0]:.3f}     {sims[1]:.3f}     {sims[2]:.3f}")
# theory: single network p_c = 0.250; coupled pair p_c = 0.6139
#    p    theory: single  pair    simulated: q = 1   q = 0.5   q = 0.2
#   0.90       0.873     0.838            0.837     0.857     0.867
#   0.80       0.762     0.709            0.689     0.733     0.752
#   0.70       0.647     0.558            0.560     0.626     0.648
#   0.65       0.588     0.460            0.482     0.584     0.591
#   0.62       0.552     0.369            0.359     0.536     0.540
#   0.60       0.527     0.000            0.000     0.498     0.513
#   0.50       0.398     0.000            0.000     0.362     0.399
#   0.40       0.257     0.000            0.000     0.212     0.249
#   0.30       0.094     0.000            0.000     0.017     0.066

单个网络的巨分量在 p=0.25p = 0.25 处连续地消失:保留 60% 时还有 52.7% 连在一起。一对一完全依赖时,崩溃来得早得多也突然得多:pc=0.6139p_c = 0.6139;保留 62% 时还有 36.9% 在运转(模拟 35.9%),保留 60% 时只剩 0。多失去 2% 的节点,系统就从三分之一在运转变为全部瘫痪。依赖减弱后图景又变回来:q=0.5q = 0.5 时保留 60% 仍有 49.8% 在运转,q=0.2q = 0.2 时与单个网络几乎一样。帕尔沙尼、布尔德列夫与哈夫林证明,降低耦合强度会使相变从一阶变为二阶 (Parshani et al. 2010)。

证据

互联的利与弊。 布鲁米特、德苏扎与莱希特在两个稀疏互联的网络上(包括由美国东南部两个电网得到的网络)研究沙堆式负荷级联:少量互联使过剩负荷可以被另一个网络吸收,大级联的概率明显降低;但互联超过一个最优水平(在他们的理想化随机网络中约为 7.5%)后,大级联又变得更可能;为了压制本网络的大级联而互联,还可能增大整个系统中最大的级联 (Brummitt et al. 2012)。

所以,「相互依赖使系统变脆弱」与「互联使系统更稳健」都只在一定条件下成立:前者依赖强的一对一依赖与「不在巨分量中就失效」的规则,后者依赖适度的、能分担负荷的连接。真实基础设施的依赖往往是部分的、有空间结构的、有备用电源与人工干预的。帕尔沙尼等人的结果还说明,崩溃是否「突然」本身是可以通过设计改变的:减少强依赖,崩溃变为渐进的,人们就有时间在半途察觉并干预。

本章小结

在容量—负荷模型中,负荷越不均匀,越容易被一个关键节点的失效拖垮:300 个节点的无标度网络中,10% 的冗余足以应付随机故障,但挡住对负荷最高节点的攻击需要约 140%。最关键的冗余平时最像浪费。两个一对一相互依赖的随机网络在保留约 61% 的节点时突然全部瘫痪,单个网络要到 25% 才连续地消失;依赖减弱后崩溃又变为连续的。这些结论都依赖明确的模型假设。

参考文献

Brummitt, Charles D., Raissa M. D’Souza, and E. A. Leicht. 2012. “Suppressing Cascades of Load in Interdependent Networks.” Proceedings of the National Academy of Sciences 109 (12): E680–89. https://doi.org/10.1073/pnas.1110586109.
Buldyrev, Sergey V., Roni Parshani, Gerald Paul, H. Eugene Stanley, and Shlomo Havlin. 2010. “Catastrophic Cascade of Failures in Interdependent Networks.” Nature 464 (7291): 1025–28. https://doi.org/10.1038/nature08932.
Motter, Adilson E., and Ying-Cheng Lai. 2002. “Cascade-Based Attacks on Complex Networks.” Physical Review E 66 (6): 065102. https://doi.org/10.1103/PhysRevE.66.065102.
Parshani, Roni, Sergey V. Buldyrev, and Shlomo Havlin. 2010. “Interdependent Networks: Reducing the Coupling Strength Leads to a Change from a First to Second Order Percolation Transition.” Physical Review Letters 105 (4): 048701. https://doi.org/10.1103/PhysRevLett.105.048701.
U.S.–Canada Power System Outage Task Force. 2004. Final Report on the August 14, 2003 Blackout in the United States and Canada: Causes and Recommendations. U.S. Department of Energy; Natural Resources Canada.

Fan, Y. (2026). Phase Transitions, Critical Slowing Down and Resilience in Complex Dynamical Systems: From Saddle-Node Bifurcations and Early-Warning Signals to Cascading Failures, Read Beside Classical Chinese Thought on Incipient Change [复杂系统相变、临界慢化与生态韧性的非线性动力学] (Version 1.0.0). 5D Library. https://doi.org/10.5281/zenodo.23770082 · Licensed under Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International (CC BY-NC-ND 4.0).