1.   引　言

20世纪80年代以来, 随着人们对复杂网络的深入研究, 逐步认识到现实生活中的网络系统并非是孤立存在的单一网络系统^[1–3], 而是与其他网络之间存在或多或少的依赖关系, 从而形成更大规模的多网络依赖系统. 如电力网和通信网^[4]、交通运输网和供应链网^[5]、航空网和贸易网^[6]等存在相互依赖关系的网络. 2003年发生在意大利的大停电事故^[7], 其根源就是由于电力网和通信网之间的高度依赖关系. 为了研究网络之间相互作用对网络鲁棒性的影响, Buldyrev等^[8]通过研究电力网和通信网之间的相互依赖关系, 提出了完全一对一的节点相互依赖级联失效模型, 并且通过生成函数法^[9–11]开发了该相依网络理论数学分析框架.

此后学者在经典完全相依网络模型基础上进一步提出了各种相依网络模型. 最早Parshani等^[12]观察到现实网络并非一对一相互依赖关系, 建立了部分点相依网络模型, 发现耦合强度的变化可以使网络渗流行为表现一阶相变和二阶相变. 在上述部分相依模型基础上, 发现减弱耦合强度可以使网络从一阶相变转换为二阶相变^[13]. Gao等^[14]扩展了双层点相依网络模型为$n$个网络组成的NON网络, 研究发现对于ER组成的环状网络, 网络的最大联通片与$n$无关. Jiang等^[15]分析了多层网络中耦合强度对网络的跨层级联和层内级联的渗流行为, 即在初始节点故障时可以通过一阶或二阶渗流转变而崩溃. 考虑到耦合模式对相依网络的影响, Zhang等^[16]研究了多重依赖关系的网络, 即一层网络节点依赖另外一层网络多个节点, 另外一层节点只依赖该层一个节点, 具有多层依赖关系的层表现出一阶相变和混合相变, 但Dong等^[17]则研究两层网络都具有多重依赖关系(即节点多对多对应关系), 发现通过加强耦合强度和连接密度, 系统抵抗级联失效能力更强. 韩伟涛和伊鹏^[18]在具有多重依赖关系的网络层中添加节点失效冗余度$r$, 通过增大$r$明显提高了系统鲁棒性. Wang等^[19]考虑到现实网络中多个节点相互协同集结成组(超级节点), 同一组内一起失效或者生存, 结果显示这种群体渗流总是一阶相变, 且与超级节点大小分布无关.

不同于上述基于节点的网络模型, Gao等^[20]提出了边完全相互依赖网络模型. 此后, 在此模型基础上研究了部分边相依情况, 即随着耦合比例$q$的变化网络出现一阶相变到二阶相变的转换, 进一步研究了有向依赖链接的边耦合相互依赖网络模型, 网络A中$\beta$比例链路依赖于网络B, 随着依赖比例的变化, 两层网络都存在不同类型的相变过 程^[21]. Zhao等^[22]提出一个关于凝固网络视角的边相依网络模型, 通过开发基于凝固网络框架的理论来进一步揭示边相依网络相对于点相依网络更鲁棒的原因. Xie等^[23]在边相依网络模型基础上, 提出了加权边耦合的相依网络模型.

上述已有研究中网络相互依赖关系都体现在单一的网络功能元素之间, 即节点和节点或者边和边的相依关系, 称为一元耦合关系. 本文考虑到现实网络中, 两层网络之间的作用关系并非完全是节点或者连边这种单纯的依赖关系, 存在一层网络节点和另外一层网络连边之间的依赖关系. 现实网络例子如供应链系统中的生产商和物流网络, 生产商的生产和供应能力依赖于物流网络的运输和分配能力, 生产商的节点(如工厂)依赖于物流网络中的边(如运输网线路)来保证产品的及时性运输和分发, 而物流网络边的存在和强度取决于生产商的需求和供应链设计. 两者共同构成了二元相互依赖网络, 即节点和边的相互依存关系. 由此本文提出了节点和连边相互依赖的二元耦合网络(BINNEC). 基于渗流理论和自平衡方程, 建立详细的数学理论分析框架来对BINNEC渗流行为进行分析. 为考虑一般性, 文中分析了三种不同的网络结构, 验证了理论框架的数值解和计算机仿真结果一致, 研究了二元耦合网络渗流行为和鲁棒性.

2.   二元耦合网络模型

传统的相依网络渗流模型如图1(a), (b)所示, 为一对一相互依赖关系, 如节点和节点相互依赖、边和边的相互依赖模型, 即都为单一元素之间的耦合关系. 这种同元一对一依赖可使得两层网络的级联失效满足无反馈条件, 元素之间级联失效效应往往是对等的. 本文中相依网络模型简化如图1(c)所示, 网络A中失效的节点导致网络B中的多条连边失效, 网络B中连边失效, 反过来导致网络A中节点失效, 节点失效和边失效导致各自层中的破坏程度不完全一样. 下面将具体研究这种不对等级联失效效应的点-边相互依赖模型.

首先建立BINNEC模型, 第1层网络A的节点大小为${N_{\text{A}}}$, 度分布为${P_{\text{A}}}(k)$; 第2层网络B的节点大小为${N_{\text{B}}}$, 度分布为${P_{\text{B}}}(k)$. 根据网络的性质和一般性, 网络中边的数量大小远远大于该网络中节点的数量, 同时也考虑到网络A的节点单纯与网络B的边相互依赖. 为保证网络B的边全部和网络A的节点都有相应的依赖关系, 那么必然网络A节点相互依赖于网络B的多条边, 即存在网络A节点$j$随机依赖于网络B中$m$条边, 定义节点依赖群规模为$m$. 本文考虑单一依赖的情况, 网络平均度取偶数度值, 使得网络A节点$j$全部只依赖B网络$m$条边, 网络B的一条边随机依赖网络A的一节点. 考虑到节点$j$依赖多条边, 本文定义$\mu$($0 \leqslant \mu < 1$)为网络A节点$j$的边依赖群的失效容忍度, 即网络A某一节点$j$依赖的$m$条边失效比例超过$\mu$时该节点失效, 否则该节点继续存活. 为保证在不同失效容忍度$\mu$下, 网络A节点的依赖群失效边数$s$为整数. 因此, 定义节点规模为$m$的依赖群的存活条件为至多可以容忍$n$条连边($n = \left\lfloor {m\mu } \right\rfloor$)失效. 假设两层网络A和B的节点大小相同, 平均度$\left\langle {{k_{\text{A}}}} \right\rangle$和$\left\langle {{k_{\text{B}}}} \right\rangle$均为4, 意味着它们两者的边数量相同. 在网络A中随机选取一节点与网络B中任意两条边尚未建立相互依赖关系的边建立依赖关系, 直到网络B中全部边都与相依网络A建立起依赖关系, 假设网络A节点失效容忍度$\mu = 0$. 由于网络模型的功能分别被定义在网络A节点和网络B边上, 因此基于网络A的点渗流和网络B的边渗流进行分析. 假设随机攻击使A网络$1 - p$比例的节点失效, 那么同一层网络中与之相连的边和相互依赖的另外一层网络中的边也失效, 同时两层网络中脱离各自网络最大联通片(giant component, GC)的节点也随着失效. 当网络B中脱离最大联通片的节点失效导致的同层边失效时, 与失效的边相互依赖的节点也会失效. 这种级联失效过程不断反复直至网络达到稳定状态, 即没有其他的节点或边失效.

图2给出了点-边二元相互依赖网络的级联失效过程, 网络A的全部节点依赖于网络B的所有边. 图中黑色实线表示连接各层中属于GC节点的边, 蓝色虚线表示连接两层的相互依赖边, 黑色虚线表示网络内失效的边. 假设网络A随机节点${a_5}$和${a_6}$遭到攻击而失效, 导致与其相连的连边失效. ${a_5}$和${a_6}$失效导致网络A中节点${a_7}$因脱离GC而失效, 继而使得与网络A中失效节点有依赖关系的网络B中的相依边失效, 这导致网络B中脱离GC的边${b_5} \leftrightarrow {b_7}$失效, 同时边${b_5} \leftrightarrow {b_7}$失效导致网络A相依节点${a_4}$失效, ${a_4}$失效使网络B的边${b_3} \leftrightarrow {b_6}$失效, 进一步导致网络B的级联失效进程. 通过这个过程, 两个相互依赖网络中级联故障的交互传播, 使得故障从A网络传播到B网络, B网络又进一步将故障传播到A网络, 不断重复两者级联失效过程, 网络中没有其他节点和边的失效, 剩余节点和边均位于网络的最大相互连接巨型组件(mutual connected giant component, MCGC)中.

3.   理论分析框架

根据上述点-边二元相依耦合网络模型, 对于两层的二元耦合网络, 采用自平衡概率理论方法, 设置一个辅助变量参数$L$, 定义为随机选择的一条连边连接到最大互连通片(MCGC)的概率. 由图3参数定义可知, 对于本文二元网络模型, 可以假设${L_{\text{A}}}$为A网络中随机选择的一条连边连接到MCGC的概率, 那么$1 - {L_{\text{A}}}$可表示为随机选择的连边没有连接到MCGC的概率. 同理, ${L_{\text{B}}}$为B网络中随机选择的连边连接到MCGC的概率, $1- {L_{\text{B}}}$表示为随机的连边没有连接到MCGC的概率.

对于A网络的随机选择的节点, 不同于节点耦合网络(一层网络的节点与另外一层网络的节点相互关联), 同样也区别于边耦合的网络(两层网络之间连边存在相互依赖关系), 节点和边相互耦合的二元耦合网络中网络A的节点与网络B中节点无直接相互作用, 但是与另外一层网络的边存在耦合关联. 对于本文相依网络模型, 网络A的节点存活的概率即连接到该网络GC的概率, 不仅要求该节点在网络A中存活, 同时也要求该节点的边依赖群失效比例不超过失效容忍度$\mu$.

因此, 重点在于找出A网络与B网络边依赖群之间的联系. 首先求出在网络A中随机选择一个节点$u$, 该节点在A网络GC中的概率为

$${D_{\text{A}}} = \sum\nolimits_k{P_{\text{A}}}(k)[1 - {{(1 - {L_{\text{A}}})}^k}],$$

(1)

其中${P_{\text{A}}}(k)$是随机选择的节点$u$度为$k$的概率, ${(1 - {L_{\text{A}}})^k}$表示该节点$k$条向外延伸的边都没有连接到GC的概率, 那么$1 - {(1 - {L_{\text{A}}})^k}$则表示该节点$k$条向外延神的边中至少有一条连接到GC的概率.

对于求出网络B中随机选择的连边存在于该网络GC的概率, 即存活的概率(需要区分的是这里与随机选择一条边连接到GC的概率不同). 本文中求出随机的连边不在巨型组件的概率, 表示该边完全失效, 即该连边的两头连接的节点均不在巨型组件当中. 如图3参数定义示意图, 当顺着该连边${D'}$找到一个端点, 该端点具有度值为k的概率${f_{\text{B}}}(k) = {{{P_{\text{B}}}(k)k}}/{{ \langle k \rangle }}$, ${P_{\text{B}}}(k)$表示网络B中随机选择的节点其度值为$k$的概率, 那么网络B内巨型组件不存在该节点的概率为

$$\tilde D_{\text{B}}^{1} = \sum\nolimits_{k'} {{f_{\text{B}}}} \left( {k'} \right){\left( {1 - {L_{\text{B}}}} \right)^{k' - 1}},$$

(2)

其中, ${(1 - {L_{\text{B}}})^{{k'} - 1}}$表示该节点除了该边以外其他${k'} - 1$条向外延伸的边都没有连接到GC的概率.

若该条边的两端点都连接一个节点, 沿着这条边的两端点都没有连接到巨型组件的概率, 即该边失效的概率:

$$\tilde D_{\text{B}}^{2} = {\left( {\tilde D_{\text{B}}^{1}} \right)^2} = {\left[ {\sum\nolimits_{k'} {{f_{\text{B}}}} \left( {k'} \right){{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}} \right]^2}.$$

(3)

反之, 该条边在巨型组件内而没有失效的概率为

$$D_{\text{B}}^{2} = 1 - {\left[ {\sum\nolimits_{k'} {{f_{\text{B}}}} \left( {k'} \right){{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}} \right]^2}.$$

(4)

对于BINNEC, 使用类似于序参数$S$来表示MCGC的大小. 将${S_{\text{A}}}$(${S_{\text{B}}}$)定义为网络A(B)中随机选择的节点属于MCGC的概率, 也是网络中最大联通片的相对大小. 网络A任一节点依赖于网络B的$m$条边, 要保证BINNEC中网络A节点存活, 首先需要考虑网络A随机选择的节点在该层网络巨型组件内, 概率为$\displaystyle\sum\nolimits_k{P_{\text{A}}}(k)[1 - {{(1 - {L_{\text{A}}})}^k} ]$. 同时, 由于网络A节点与另外一层网络连边存在直接联系, 因此需要考虑两层网络之间的约束条件, 该节点在网络A中存活的同时其对应相依关系的$m$条边也应满足存活的条件, 即网络A中节点对应于网络B中的边依赖群应满足失效比例不超过$\mu$, 该概率为

$$\begin{split} \;& M = \sum\limits_{s = 0}^{\left\lfloor {m\mu } \right\rfloor } {\bigg( {\begin{array}{*{20}{c}} m \\ s \end{array}} \bigg)} \left[ {{{\left( {\widetilde {D_{\text{B}}^{2}} } \right)}^s} \cdot {{\left( {D_{\text{B}}^{2}} \right)}^{m - s}}} \right] \\ =\;& \sum\limits_{s = 0}^{\left\lfloor {m\mu } \right\rfloor } {\bigg( {\begin{array}{*{20}{c}} m \\ s \end{array}} \bigg)} {\bigg[ {\sum\limits_{k'} {{f_{\text{B}}}} \left( {k'} \right){{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}} \bigg]^{2s}} \\ & \times \bigg\{ {1 - {{\bigg[ {\sum\limits_{k'} {{f_{\text{B}}}} \left( {k'} \right){{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}} \bigg]}^2}} \bigg\}^{m - s}, \end{split}$$

(5)

其中, ${ \Big[ {\displaystyle\sum\nolimits_{k'} {{f_{\text{B}}}} \left( {k'} \right){{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}} \Big]^{2 s}}$表示网络B中随机选择失效$s$条边的概率; $\Big\{ 1 - \Big[ \displaystyle\sum\nolimits_{k'} {{f_{\text{B}}}} \left( {k'} \right)\times {{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}} \Big]^2 \Big\}^{m - s}$为网络B中除去失效边$s$后依赖群存活的概率; $s$表示网络A中节点$j$依赖群的失效边数; $\left\lfloor {m\mu } \right\rfloor$表示小于等于$m\mu$的最大整数.

因此, 网络A随机选择的节点在MCGC存活的概率为

$$\begin{split} {S_{\text{A}}}\;& = {D_{\text{A}}} \cdot M = \sum _k {P_{\text{A}}}(k) [1 - {{(1 - {L_{\text{A}}})}^k}] \\ &\times \sum\limits_{s = 0}^{\left\lfloor {m\mu } \right\rfloor } {\bigg( {\begin{array}{*{20}{c}} m \\ s \end{array}} \bigg)} {\bigg[ {\sum\limits_{k'} {{f_{\text{B}}}} \left( {k'} \right){{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}} \bigg]^{2s}}\\ &\times \bigg\{ {1 - {{\bigg[ {\sum\limits_{k'} {{f_{\text{B}}}} \left( {k'} \right){{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}} \bigg]}^2}} \bigg\}^{m - s}. \end{split}$$

(6)

对于网络B, 随机选择的节点属于MCGC的概率与网络A不同, 由于网络B节点与网络A节点或边无直接相互作用关系. 因此该概率为

$${S_{\text{B}}} = \sum\nolimits_k {{P_{\text{B}}}({k'})} \left[ {1 - {{(1 - {L_{\text{B}}})}^{{k'}}}} \right].$$

(7)

如随机攻击网络A中$1 - p$比例的节点使其失效, 那么网络中只剩下$p$比例的节点保留在网络中. 因此随机选择的节点存活的概率变为原来的$p$倍, 则有:

$$\begin{split} {S_{\text{A}}} \;&= p \cdot {D_{\text{A}}} \cdot M = p \cdot \sum _k{P_{\text{A}}}(k)[1 - {{(1 - {L_{\text{A}}})}^k }] \\ &\times\sum\limits_{s = 0}^{\left\lfloor {m\mu } \right\rfloor } {\left( {\begin{array}{*{20}{c}} m \\ s \end{array}} \right)} {\bigg[ {\sum\limits_{k'} {{f_{\text{B}}}} \left( {k'} \right){{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}} \bigg]^{2s}}\\ &\times{\bigg\{ {1 - {{\bigg[ {\sum\limits_{k'} {{f_{\text{B}}}} \left( {k'} \right){{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}} \bigg]}^2}} \bigg\}^{m - s}}. \\[-1pt] \end{split}$$

(8)

随机攻击网络A中$1 - p$比例的节点使其失效, 由于网络B节点与网络A节点或者边无直接相互作用关系. 因此该概率不受影响:

$${S_{\text{B}}} = \sum\nolimits_k {{P_{\text{B}}}({k'})} \left[ {1 - {{(1 - {L_{\text{B}}})}^{{k'}}}} \right].$$

(9)

接下来, 推导关于${L_{\text{A}}}$和${L_{\text{B}}}$的自平衡方程. 图3中网络A随机选择的一条边${L_{\text{A}}}$, 其中一头端点$u$需连接到A网络巨型组件中, 概率为$\displaystyle\sum\nolimits_k {{f_{\text{A}}}} (k)[1 - {(1 - {L_{\text{A}}})^{k - 1}}]$. 既保证了该边在层内存活, 同时也需要保证该边的端点依赖B网络的边依赖群在失效容忍度$\mu$内存活. 这要求随机选择的边既需要满足A网络的层内条件, 同时也应该满足两层网络的层间条件. 因此, ${L_{\text{A}}}$自平衡方程为

$$\begin{split} {L_{\text{A}}} \;&= \sum\limits_k {{f_{\text{A}}}} (k)[1 - {(1 - {L_{\text{A}}})^{k - 1}}] \cdot \sum\limits_{s = 0}^{\left\lfloor {m\mu } \right\rfloor } {\left( {\begin{array}{*{20}{c}} m \\ s \end{array}} \right)} \\ &\times{\bigg[ {\sum\limits_{k'} {{f_{\text{B}}}} \left( {k'} \right){{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}} \bigg]^{2s}}\\ &\times \bigg\{ {1 - {{\bigg[ {\sum\limits_{k'} {{f_{\text{B}}}} \left( {k'} \right){{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}} \bigg]}^2}} \bigg\}^{m - s},\\[-1pt] \end{split}$$

(10)

式中, ${f_{\text{A}}}\left( k \right)$和${f_{\text{B}}}\left( {k'} \right)$表示网络A和网络B随机选择的边找到GC的节点, 其度值分布为$k$和$k'$的概率.

如果通过在网络A中随机移除$1 - p$比例的节点, 使该网络中随机保留$p$比例的节点, 那么随机选择的一条边连接MCGC的概率变为

$$\begin{split} {L_{\text{A}}} \;&= p \cdot \sum\limits_k {{f_{\text{A}}}} (k)[1 - {(1 - {L_{\text{A}}})^{k - 1}}] \cdot \sum\limits_{s = 0}^{\left\lfloor {m\mu } \right\rfloor } {\left( {\begin{array}{*{20}{c}} m \\ s \end{array}} \right)} \\ &\times{\bigg[ {\sum\limits_{k'} {{f_{\text{B}}}} \left( {k'} \right){{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}} \bigg]^{2s}}\\ &\times{\bigg\{ {1 - {{\bigg[ {\sum\limits_{k'} {{f_{\text{B}}}} \left( {k'} \right){{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}} \bigg]}^2}} \bigg\}^{m - s}}. \end{split}$$

(11)

网络B中${L_{\text{B}}}$的自平衡方程与网络A中${L_{\text{A}}}$的自平衡方程类似, 同样需要满足层内条件和层间条件. 首先, 从网络B中随机选择一条边连接到GC的概率为

$$D_{\text{B}}^{1} = 1 - \tilde D_{\text{B}}^{1} = \sum\limits_k {{f_{\text{B}}}\left( {k'} \right)\left[ {1 - {{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}} \right]} .$$

(12)

考虑到网络B中这条随机选择的边与网络A节点存在直接依赖关系, 为保证该边不失效, 那么其依赖的网络A节点应该同时存活. 从网络A中随机选择一个节点在巨型组件中的概率为$\displaystyle \sum\nolimits_k {P_{\text{A}}}(k)[1 - {{(1 - {L_{\text{A}}})}^k}]$, 该网络A节点同时依赖于网络B中除了该条依赖边的其他$m - 1$条随机依赖的边, 要求网络B这$m - 1$条依赖边在失效容忍度$\mu$内存活, 网络A中该节点才能存活. 网络B中一条随机选择的边存活的概率为$D_{\text{B}}^{2}$. 由此, 可以写出${L_{\text{B}}}$的自平衡方程为

$$\begin{split} &{L_{\text{B}}} = \sum\limits_k {{f_{\text{B}}}\left( {k'} \right)[1 - {{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}]} \\ &\times\sum {_k{P_{\text{A}}}(k)[1 - {{(1 - {L_{\text{A}}})}^k}} ] \cdot \\ &\times\sum\limits_{s = 0}^{\left\lfloor {m\mu } \right\rfloor } {\left( {\begin{array}{*{20}{c}} {m - 1} \\ s \end{array}} \right)} \cdot {\bigg[ {\sum\limits_{k'} {{f_{\text{B}}}} \left( {k'} \right){{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}} \bigg]^{2s}}\\ &\times{\bigg\{ {1 - {{ \bigg[ {\sum\limits_{k'} {{f_{\text{B}}}} \left( {k'} \right){{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}} \bigg]}^2}} \bigg\}^{m - 1 - s}}. \\[-1pt] \end{split}$$

(13)

如果随机攻击网络A中$1 - p$比例的节点使之失效, 那么网络A中仅剩余$p$比例的节点, 此时:

$$\begin{split} &{L_{\text{B}}} = \sum\limits_k {{f_{\text{B}}}\left( {k'} \right)[1 - {{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}]} \cdot p \\ &\times \sum {_k{P_{\text{A}}}(k)[1 - {{(1 - {L_{\text{A}}})}^k}} ] \\ &\times\sum\limits_{s = 0}^{\left\lfloor {m\mu } \right\rfloor } {\left( {\begin{array}{*{20}{c}} {m - 1} \\ s \end{array}} \right)} \cdot {\bigg[ {\sum\limits_{k'} {{f_{\text{B}}}} \left( {k'} \right){{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}} \bigg]^{2s}}\\ &\times { \bigg\{ {1 - {{\bigg[ {\sum\limits_{k'} {{f_{\text{B}}}} \left( {k'} \right){{\left( {1 - {L_{\text{B}}}} \right)}^{k' - 1}}} \bigg]}^2}} \bigg\}^{m - 1 - s}}. \\[-1pt] \end{split}$$

(14)

在改变两层网络的平均度使$m$取得不同的值, 即网络A节点全部依赖于B网络$m$条边, 同时给定保留比例参数$p$和失效容忍度$\mu$, 可以联立(11)式和(14)式求解方程, 得到相应${L_{\text{A}}}$和${L_{\text{B}}}$的值, 代入(8)式和(9)式分别求出网络中${S_{\text{A}}}$和${S_{\text{B}}}$的值.

通过分析${S_{\text{A}}}$和${S_{\text{B}}}$与$p$的关系, 得到BINNEC网络的渗流行为, 即其相变类型和相变阈值${p_{\text{c}}}$. 为了具体得到网络发生渗流行为时的关键阈值${p_{\text{c}}}$, 将方程转换为

$$\left\{ {\begin{aligned} & {{L_{\text{A}}} = {F_1}(p,{\text{ }}{L_{\text{B}}},{\text{ }}m)}, \\ & {{L_{\text{B}}} = {F_2}(p,{\text{ }}{L_{\text{A}}},{\text{ }}m)}. \end{aligned}} \right.$$

(15)

如果系统发生不连续相变, 即为一阶相变阈值点$p_{\text{c}}^{1}$, 则在${p_{\text{c}}} = p_{\text{c}}^{1}$处, ${L_{\text{A}}} = {F_1}(p, {\text{ }}{L_{\text{B}}}, {\text{ }}m)$和${L_{\text{B}}} = {F_2}(p, {L_{\text{A}}}, {\text{ }}m)$所表示的曲线正好相切, 即:

$$\frac{{\partial {F_1}(p_{\text{c}}^{1},{\text{ }}{L_{\text{B}}},{\text{ }}m)}}{{\partial {L_{\text{B}}}}} \cdot \frac{{\partial {F_2}(p_{\text{c}}^{1},{\text{ }}{L_{\text{A}}},{\text{ }}m)}}{{\partial {L_{\text{A}}}}} = 1.$$

(16)

通过(16)式可以得到系统非连续相变关键阈值$p_{\text{c}}^{1}$.

如果系统发生的是连续相变, 即为二阶相变阈值点$p_{\text{c}}^{2}$. 为求得网络的相变阈值点, 将${L_{\text{A}}} \to 0$代入(11)式求得$L_{\text{B}}^{0}$, 继续将求得的$L_{\text{B}}^{0}$和${L_{\text{A}}} \to 0$代入(14)式即可得到关键阈值$p_{\text{c}}^{2}$的解.

4.   理论和仿真分析

第3节推导出了二元耦合网络的理论公式, 下面假设二元耦合网络中子网络A和子网络B的节点数目均为N = 10⁵, 平均度相同. 将该理论分析框架应用于3种经典的网络结构, 即随机规则网络(RR)、随机网络(ER)和无标度网络(SF), 验证上述理论的正确性并且进一步分析系统的相变行为及其鲁棒性.

4.1   RR BINNEC结果分析

假设网络A和网络B平均度相同, 即$\left\langle {{k_{\text{A}}}} \right\rangle = \left\langle {{k_{\text{B}}}} \right\rangle = \left\langle k \right\rangle$. 由网络性质可知${P_{\text{A}}}(k) = 1$, ${P_{\text{B}}}(k) = 1$, 将这些公式分别代入(11)式和(14)式, 可得到辅助参量${L_{\text{A}}}$和${L_{\text{B}}}$的自平衡方程式:

$$\begin{split} {L_{\text{A}}} =\;& p \cdot \left[ {1 - {{(1 - {L_{\text{A}}})}^{\left\langle k \right\rangle - 1}}} \right] \cdot \sum\limits_{s = 0}^{\left\lfloor {m\mu } \right\rfloor } {\left( {\begin{array}{*{20}{c}} m \\ s \end{array}} \right)}\\ &\times {\left[ {{{(1 - {L_{\text{B}}})}^{2(\left\langle k \right\rangle - 1)}}} \right]^s} \\ &\times {\left[ {1 - {{(1 - {L_{\text{B}}})}^{2(\left\langle k \right\rangle - 1)}}} \right]^{m - s}}, \end{split}$$

(17)

$$\begin{split} {L_{\text{B}}} =\;& \left[ {1 - {{(1 - {L_{\text{B}}})}^{\left\langle k \right\rangle - 1}}} \right] \cdot p \cdot \left[ {1 - {{(1 - {L_{\text{A}}})}^{\left\langle k \right\rangle }}} \right] \\ & \times\sum\limits_{s = 0}^{\left\lfloor {m\mu } \right\rfloor } {\left( {\begin{array}{*{20}{c}} {m - 1} \\ s \end{array}} \right)} \cdot {\left[ {{{(1 - {L_{\text{B}}})}^{2(\left\langle k \right\rangle - 1)}}} \right]^s}\\ &\times{\left[ {1 - {{(1 - {L_{\text{B}}})}^{2(\left\langle k \right\rangle - 1)}}} \right]^{m - 1 - s}}. \\[-1pt] \end{split}$$

(18)

根据(8)式和(9)式, 可以得到网络A和B的最大巨型组件${S_{\text{A}}}$和${S_{\text{B}}}$的方程式:

$$\begin{split} {S_{\text{A}}} =\;& p \cdot \left[ {1 - {{(1 - {L_{\text{A}}})}^k}} \right] \cdot \sum\limits_{s = 0}^{\left\lfloor {m\mu } \right\rfloor } {\bigg( {\begin{array}{*{20}{c}} m \\ s \end{array}} \bigg)} \\ &\times {\left[ {{{(1 - {L_{\text{B}}})}^{2(\left\langle k \right\rangle - 1)}}} \right]^s} \\ &\times {\left[ {1 - {{(1 - {L_{\text{B}}})}^{2(\left\langle k \right\rangle - 1)}}} \right]^{m - s}}, \end{split}$$

(19)

$${S_{\text{B}}} = \left[ {1 - {{(1 - {L_{\text{B}}})}^{\left\langle k \right\rangle - 1}}} \right].$$

(20)

在随机攻击A网络$1 - p$比例的节点下, 给定参数$p$, $m$和不同的失效容忍度$\mu$代入(17)式和(18)式, 得到辅助参量${L_{\text{A}}}$和${L_{\text{B}}}$的值, 接着代入(19)式和(20)式得到网络A和B存活的最大巨型组件S_A和S_B与$p$的关系如图4所示. 图中不同颜色实线代表不同边依赖群$m$下的理论计算结果, 不同颜色符号代表不同$m$值下计算机仿真结果. 图中可以很明显看出理论与仿真完全吻合, 说明了理论分析的正确性.

图4(a), (b)为RR BINNEC在$\mu \in [0, {\text{ }}{1 \mathord{\left/ {\vphantom {1 m}} \right. } m})$(表示网络A节点的边依赖群只要失效一条边, 则该节点对应失效的情况)范围内S_A和S_B的大小图. 通过比较, 在相同的失效容忍度$\mu$前提下, 网络A和网络B的渗流阈值与依赖群规模$m$成反比, 即网络A节点依赖群规模$m$越大, 网络鲁棒性越强.

通过纵向比较图4(a), (c)或图4(b), (d), 面对失效容忍度$\mu$在两种极端的情况下(即网络A节点的边依赖群只要失效一边, 则对应依赖的节点失效; 网络A节点的边依赖群全部失效, 则对应依赖节点才失效), 失效容忍度$\mu$越大, 其对应的渗流阈值越小, 即网络鲁棒性越强.

图5(a)中图例m2-0表示依赖群规模$m = 2$时, 网络A节点的失效容忍能力不足以容忍任意一条边失效, 即边依赖群中任何一条边失效时, 则该节点失效; 同样类似的m2-1表示依赖群规模$m = 2$时, 网络A节点的失效容忍能力可以至多容忍1条边失效, 即边依赖群中任意2条及以上边失效时, 则该节点失效. 类似的图5(b)中图例定义也是如此.

由图5清晰可知, 对于边依赖群规模$m = 2$或者$m = 3$的网络, 改变依赖群的失效容忍度$\mu$, 两个子网络的一阶相变点始终一样, 表明网络A和网络B渗流阈值一致性不受失效容忍度$\mu$的影响. 不同的是随着初始攻击比例增大, 网络A的最大联通片几乎随着攻击节点比例$1 - p$成线性相关减小; 而网络B最大联通片在一个比较大的值减小, $p$接近一阶相变点$p_{\text{c}}^{1}$时才突变为零.

4.2   ER BINNEC结果分析

假设网络A和网络B的平均度分别为$\left\langle {{k_{\text{A}}}} \right\rangle$和$\left\langle {{k_{\text{B}}}} \right\rangle$, 则ER网络的度分布分别为

$$\left\{ {\begin{aligned} & {{P_{\text{A}}}(k) = \frac{{{{\left\langle {{k_{\text{A}}}} \right\rangle }^k}}}{{k!}}{{\text{e}}^{ - \left\langle {{k_{\text{A}}}} \right\rangle }}}, \\ &{{P_{\text{B}}}(k) = \frac{{{{\left\langle {{k_{\text{B}}}} \right\rangle }^k}}}{{k!}}{{\text{e}}^{ - \left\langle {{k_{\text{B}}}} \right\rangle }}}. \end{aligned}} \right.$$

(21)

将上述代入(11)式和(14)式, 可得到辅助参量${L_{\text{A}}}$和${L_{\text{B}}}$的自平衡方程式:

$$\begin{split} {L_{\text{A}}} =\;& p \cdot \left( {1 - {{\text{e}}^{\left\langle {{k_{\text{A}}}} \right\rangle {L_{\text{A}}}}}} \right) \cdot \sum\limits_{s = 0}^{\left\lfloor {m\mu } \right\rfloor } {\left( {\begin{array}{*{20}{c}} m \\ s \end{array}} \right)} \\ &\times {\left( {{{\text{e}}^{ - 2\left\langle {{k_{\text{A}}}} \right\rangle {L_{\text{B}}}}}} \right)^s}{\left( {1 - {{\text{e}}^{ - 2\left\langle {{k_{\text{A}}}} \right\rangle {L_{\text{B}}}}}} \right)^{m - s}},\end{split}$$

(22)

$$\begin{split}{L_{\text{B}}} = \;&\left( {1 - {{\text{e}}^{ - \left\langle {{k_{\text{B}}}} \right\rangle {L_{\text{B}}}}}} \right) \cdot p \cdot \left( {1 - {{\text{e}}^{ - \left\langle {{k_{\text{A}}}} \right\rangle {L_{\text{A}}}}}} \right) \\ &\times \sum\limits_{s = 0}^{\left\lfloor {m\mu } \right\rfloor } {\left( {\begin{array}{*{20}{c}} {m - 1} \\ s \end{array}} \right)} \cdot {\left( {{{\text{e}}^{ - 2\left\langle {{k_{\text{B}}}} \right\rangle {L_{\text{B}}}}}} \right)^s}\\ &\times{\left( {1 - {{\text{e}}^{ - 2\left\langle {{k_{\text{B}}}} \right\rangle {L_{\text{B}}}}}} \right)^{m - 1 - s}}.\end{split}$$

(23)

根据(8)式和(9)式, 可以得到网络A和B的最大巨型组件${S_{\text{A}}}$和${S_{\text{B}}}$的方程式:

$$\begin{split}{S_{\text{A}}} = \;&p \cdot \left( {1 - {{\text{e}}^{\left\langle {{k_{\text{A}}}} \right\rangle {L_{\text{A}}}}}} \right) \cdot \sum\limits_{s = 0}^{\left\lfloor {m\mu } \right\rfloor } {\left( {\begin{array}{*{20}{c}} m \\ s \end{array}} \right)}\\ &\times {\left( {{{\text{e}}^{ - 2\left\langle {{k_{\text{B}}}} \right\rangle {L_{\text{B}}}}}} \right)^s}{\left( {1 - {{\text{e}}^{ - 2\left\langle {{k_{\text{B}}}} \right\rangle {L_{\text{B}}}}}} \right)^{m - s}},\end{split}$$

(24)

$${S_{\text{B}}} = \left( {1 - {{\text{e}}^{ - \left\langle {{k_{\text{B}}}} \right\rangle {L_{\text{B}}}}}} \right).$$

(25)

为了具体得到网络发生渗流行为时的关键阈值${p_{\text{c}}}$, (22)式和(23)式可以表示为

$$\left\{ {\begin{aligned} &{{L_{\text{A}}} = {F_1}(p,{\text{ }}{L_{\text{B}}},{\text{ }}m,{\text{ }}\mu )} ,\\ & {{L_{\text{B}}} = {F_2}(p,{\text{ }}{L_{\text{A}}},{\text{ }}m,{\text{ }}\mu )} . \end{aligned}} \right.$$

(26)

对于给定的依赖群规模$m$和失效容忍度$\mu$, 根据$m$和$\left\langle k \right\rangle$之间的关系. 如果系统发生不连续相变, 即为一阶相变阈值点$p_{\text{c}}^{1}$, 则在${p_{\text{c}}} = p_{\text{c}}^{1}$处, ${L_{\text{A}}} = {F_1}(p, {\text{ }}{L_{\text{B}}})$和${L_{\text{B}}} = {F_2}(p, {\text{ }}{L_{\text{A}}})$所表示的曲线正好相切, 即:

$${\left. {1 = \frac{{\partial {F_1}(p,{\text{ }}{L_{\text{B}}})}}{{\partial {L_{\text{B}}}}} \cdot \frac{{\partial {F_2}(p,{\text{ }}{L_{\text{A}}})}}{{\partial {L_{\text{A}}}}}} \right|_{p = {p_{\text{c}}^{1}}}}.$$

(27)

通过(27)式可以求得一组解$p_{\text{c}}^{1}$, 即对应网络一阶相变的情况.

图6为ER BINNEC在极端值$\mu$和不同m时, S_A和S_B与p的关系图, 其中网络节点规模$N = {10^5}$, 结果为重复30次的平均值, 符号代表仿真分析结果, 实线表示理论分析结果, 可以看出理论分析结果与仿真分析结果完全吻合, 再次验证了本文理论分析框架的正确性. 此外, ER BINNEC在$\mu \in [0, {\text{ }}{1 \mathord{\left/ {\vphantom {1 m}} \right. } m})$和$\mu \in [{{m - 1} \mathord{\left/ {\vphantom {{m - 1} m}} \right. } m}, {\text{ }}1)$两种极端情况下, 即网络A节点对应边依赖群失效一边, 则对应依赖的节点失效; 网络A节点对应边依赖群全部失效, 则对应依赖节点才失效. 通过对比分析, 无论失效容忍度$\mu$在哪种情况下, 网络节点对应的依赖群规模越大, 两个子网网络A和B的渗流阈值越小, 即网络鲁棒性越强.

图7给出了依赖群规模$m = 2$和$m = 3$时, 失效容忍度在$\mu \in [0, {\text{ }}1)$范围内, 网络A和网络B的渗流阈值对比. 由图7(a), (b)可知依赖群规模$m = 2, {\text{ }}3$时, 节点的边依赖群能够容忍不同的边失效, 随着可容忍失效边增加, 网络渗流阈值均逐渐减小, 意味着网络抵抗级联失效的能力增强, 这点结论可以很好的从现实角度被理解.

4.3   SF BINNEC结果分析

根据文献[20]无标度网络的度分布可以被 描述为

$$\begin{split} & \\[0mm] & {P_{\text{A}}}(k) = {P_{\text{B}}}(k) \\ =\;& \frac{{\displaystyle\sum\nolimits_{{k_{\min }}}^{{k_{\max }}} {{{\left( {k + 1} \right)}^{1 - \lambda }} - {k^{1 - \lambda }}} }}{{{{\left( {{k_{\max }} + 1} \right)}^{1 - \lambda }} - {{\left( {{k_{\min }} + 1} \right)}^{1 - \lambda }}}}. \end{split}$$

(28)

假设式中${k_{\min }} = 2$, ${k_{\max }} = 70$, $\lambda = 2.7$或$\lambda = 2.1$, 将这些参数代入(8)式—(14)式和(28)式, 采用图解法或数值分析法求得无标度网络一阶相变的关键阈值$p_{\text{c}}^{1}$.

由图8可以看出, 在SF BINNEC对于不同的依赖群规模$m = 2$和$m = 3$, 随着依赖群规模由2到3, 网络A和B的渗流阈值均变小, 即网络的鲁棒性增强.

图9(a)表示在$\lambda = 2.716$、平均度$\left\langle k \right\rangle = 4$、依赖群规模$m = 2$, 对应不同$\mu$值下, 两个子网络的最大巨型组件$S$与$p$的关系图; 图9(b)表示在$\lambda = 2.1521$、平均度$\left\langle k \right\rangle = 6$、依赖群规模$m = 3$, 对应不同$\mu$值, 两个子网络的最大巨型组件$S$与$p$的关系图. 由图9可知$m = 2, {\text{ }}3$时, 随着可容忍失效边增加, 网络渗流阈值明显减小, 网络鲁棒性显著增强, 这点结论同ER BINNEC一样.

通过上述RR, ER, SF网络的分析. 为了进一步直观分析依赖群规模$m$和失效容忍度对网络鲁棒性的影响. 从RR, ER, SF三种网络中取出ER-ER网络进行分析, 图10为网络依赖群规模m分别取2, 3, 4和5时网络A的一阶相变阈值$p_{\text{c}}^{1}$与$\mu$的关系图. 图中不同颜色的实线分别代表不同的依赖群规模m. 考虑到网络B中失效的连边必须为整数, 因此对$m\mu$进行向下取整操作, 从而网络的一阶相变阈值在$m\mu$非整数范围内为直线, $m\mu$整数时呈现阶梯式下降现象.

图10中网络在不同的依赖群下, 均发生一阶相变. 当依赖群规模m一定时, 随着参数$\mu$的逐渐增大, 网络A和网络B的一阶相变阈值呈现出阶梯状的降低趋势, 即系统的鲁棒性随着参数$\mu$的增大而增大. 这一结论可以很好的理解, 当网络A节点的边依赖群规模一定时, 网络A节点的失效容忍度$\mu$越大意味着网络该节点越不容易失效, 系统鲁棒性也就越强.

由图10明显可以看出$m2$, $m3$, $m4$和$m5$面对不同的参数$\mu$时, 它们之间的渗流阈值始终是逐级降低, 所以当失效容忍度$\mu$一定时, 网络A和B的一阶相变阈值随着依赖群m增大而逐渐降低, 即系统鲁棒性随着依赖规模的增大而增强. 这一结论可以理解为, 网络A的节点失效容忍度$\mu$一定时, 该节点的边依赖群规模m越大, 一方面导致$m\mu$变大即容忍该节点失效的能力变强, 另外一方面该节点的边依赖群规模的增大, 可以使失效组合增多, 同样可以降低该节点相依失效的概率.

5.   度分布的影响

网络的度分布在一定程度上反映了其自身拓扑结构, 度分布越宽意味着节点的度数多样性, 即大多数节点的度数较低但同时存在少数度数非常高的节点. 三种经典网络结构(RR网络、ER网络和无标度网络)因存在不同的大度节点, 从而使得度分布逐渐变宽. 为进一步分析两层网络的度分布对BINNEC网络的鲁棒性影响, 假设不同BINNEC网络的平均度都为6, 即RR网络、ER网络、$\lambda = 2.15$的SF网络的平均度均相同, 网络的依赖群规模都为$m = 3$, 通过理论分析和仿真可以得到图11. 由图11可知, 在失效容忍度$\mu \in [0, {\text{ }}1/3)$时, 即网络A节点不允许边依赖群有边失效, ${p_{{\text{c}} \cdot {\text{RR}}}} = 0.42$$<$${p_{{\text{c}} \cdot {\text{ER}}}} = 0.43$$<$${p_{{\text{c}} \cdot {\text{SF}}}} = 0.45$, 一阶相变阈值随着度分布越宽变得越大, 意味着鲁棒性变弱. 相反, 对于失效容忍度在$\mu \in [1/3, {\text{ }}1)$时, 即网络A节点存活条件为至多可以容忍${n_1} = 1$和${n_2} = 2$条连边失效下, 网络的一阶相变阈值都为$p{}_{{\text{c}} \cdot {\text{RR}}}$$>$${p_{{\text{c}} \cdot {\text{ER}}}}$$>$${p_{{\text{c}} \cdot {\text{SF}}}}$, 意味着随着度分布变宽网络的鲁棒性反而变强.

该结论可以这样理解, 当网络A节点的失效容忍度不足以容忍边依赖群一条边失效时, 意味着边依赖群有一条边失效, 则该节点便失效. 因此网络度分布越宽, 那么网络A存在更多度大的节点依赖网络B连边, 边依赖群中一条边失效便可导致节点失效, 网络越来越脆弱. 而在网络A中节点的失效容忍度可以容忍边依赖群至少失效一条边时, 随着度分布变宽, 网络存在越来越多度大的节点不容易失效, 网络B存活的边也越多, 意味着依赖群数量也越多, 网络节点越不容易失效, 网络的鲁棒性变强.

6.   结　论

由于现实世界中网络之间的相互依赖关系比较复杂, 不只是节点和节点之间、边和边之间单一的相互依赖关系, 也存在节点和边之间相互依赖的情况. 本文提出节点和边耦合的二元相互依赖网络(BINNEC), 研究考虑节点依赖于数条边的情况. 通过数值仿真研究和理论分析发现, RR网络、ER网络和无标度网络构成的BINNEC在面对随机失效下, 均表现出一阶相变行为; 随着网络依赖群规模m增大, 网络的一阶相变阈值减小, 即网络鲁棒性变强; 在网络边依赖群规模m一定时, 失效容忍度$\mu$增大意味着容忍边失效的能力增强, 抵抗级联失效的能力也就越强; 在网络$\left\langle k \right\rangle$, m和$\mu$确定的条件下, 网络A节点的失效容忍度不足以容忍边依赖群一条边失效时, 随着网络度分布变宽, 网络的鲁棒性变弱; 当网络A节点的失效容忍度可以容忍边依赖群至少一条边失效时, 随着网络度分布变宽, 网络的鲁棒性反而增强.