保度换边优化 R 指标，得到"洋葱"结构（Schneider 等 2011）

# 随机选两条边 (i,j),(k,l) 改为 (i,k),(j,l)（或 (i,l),(j,k)），保持每个节点的度不变；只有当自适应度攻击下的 R 增大时才接受。
# 这里用小网络（N=200）和有限次尝试，仅演示思想；度同配系数 r 用来观察"洋葱化"往往伴随同配性升高。
# 解读：R 明显增大而度分布不变；优化后网络的度同配系数通常上升，高度节点彼此相连形成核心。
import numpy as np
import networkx as nx

def R_adaptive_degree(G):
    H = G.copy(); N = G.number_of_nodes(); total = 0.0
    for _ in range(N):
        v = max(H.degree, key=lambda x: x[1])[0]; H.remove_node(v)
        total += max((len(c) for c in nx.connected_components(H)), default=0) / N
    return total / N

G = nx.barabasi_albert_graph(200, 2, seed=5)
G_init = G
rng = np.random.default_rng(0)
R = R_adaptive_degree(G)
print("初始 R = %.3f，同配系数 r = %.3f" % (R, nx.degree_assortativity_coefficient(G)))
edges = list(G.edges()); accepted = 0
for t in range(1500):
    (i, j), (a, b) = [edges[x] for x in rng.choice(len(edges), 2, replace=False)]
    if len({i, j, a, b}) < 4:
        continue
    new = [(i, a), (j, b)] if rng.random() < 0.5 else [(i, b), (j, a)]
    if any(G.has_edge(*e) for e in new):
        continue
    H = G.copy(); H.remove_edges_from([(i, j), (a, b)]); H.add_edges_from(new)
    if not nx.is_connected(H):
        continue
    Rn = R_adaptive_degree(H)
    if Rn > R:
        G, R = H, Rn; edges = list(G.edges()); accepted += 1
print("接受换边 %d 次，优化后 R = %.3f，同配系数 r = %.3f" % (accepted, R, nx.degree_assortativity_coefficient(G)))

import matplotlib.pyplot as plt
def attack_s(G):
    H = G.copy(); N = G.number_of_nodes(); s = [1.0]
    for _ in range(N):
        H.remove_node(max(H.degree, key=lambda x: x[1])[0])
        s.append(max((len(c) for c in nx.connected_components(H)), default=0) / N)
    return np.array(s)
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(10, 4))
for H, lab, ls in [(G_init, "original BA", "--"), (G, "optimized", "-")]:
    s = attack_s(H)
    ax1.plot(np.linspace(0, 1, len(s)), s, ls=ls, lw=2, label="%s (R = %.3f)" % (lab, s[1:].mean()))
ax1.set_xlim(0, 0.5); ax1.set_xlabel("fraction removed (adaptive degree attack)"); ax1.set_ylabel("largest component / N")
ax1.legend(); ax1.set_title("Degree-preserving swaps raise R")
dg = dict(G.degree()); pos = nx.kamada_kawai_layout(G)
nx.draw_networkx_edges(G, pos, ax=ax2, edge_color="0.6", alpha=0.4)
nx.draw_networkx_nodes(G, pos, ax=ax2, node_size=[min(120, 40 + 4 * dg[v]) for v in G], node_color=[dg[v] for v in G], cmap="viridis")
ax2.set_title("Optimized network (size/color = degree)"); ax2.axis("off")
plt.tight_layout(); plt.show()