约束可满足性问题(constraint satisfaction problem, 简称CSP)是人工智能中的重要研究领域^[1-3].CSP问题是三元组〈X, D, C〉:X表示变量的集合, 记为X={x₁, …, x_n}; D是关于变量x₁, …, x_n的取值域的集合, 记为D={D₁, …, D_n}; C是约束的集合, 每一个约束c∈C表示为c=〈σ, ρ〉, 其中, 约束范围σ是变量的列表, 约束关系ρ是关于σ中变量取值域的笛卡尔积的子集.可满足性判定问题(satisfiability, 简称SAT)是典型的CSP问题, 即:给定一个合取范式公式(conjunctive normal form, 简称CNF)F, SAT判定问题指是否存在一组指派使得F为真.

SAT问题是计算机科学中的核心问题, 其理论和应用研究是计算机与数理逻辑界学者共同关注的一个重大问题.随着SAT问题的深入研究, 各类具有特殊结构的SAT问题受到越来越多的学者重视.限制子句的长度或变量出现的次数, 实质上是在特殊结构下研究SAT问题, 目前研究较多的是k-SAT问题.当k≥3时, k-SAT问题是著名的NP完全问题.已经发现:在随机k-SAT实例中, 子句个数与变元个数的比值a是一个重要参数, 它不仅对公式的可满足性产生影响, 还对公式可满足性判定的难易程度产生影响.一方面, 在可满足性问题的相变现象的研究中, 随机统计现象表明:满足与不可满足之间出现一种临界现象, 对于随机k-SAT问题, 存在可满足的相变点a_d, 当a < a_d时, 实例高概率可满足; 当a > a_d, 实例高概率不可满足^{[4, 5]}.在随机3-SAT问题中, 尽管人们不知道a_d的确切值, 但研究结果表明:a_d至少为3.52^[6], 至多为4.506^[7].另一方面, 大多数优秀的可满足性判定算法能够解决的实例与a密切相关.如:基于纯文字规则的算法能够有效求解a < 1.6区域的随机3-SAT实例(如DPLL算法); 基于正文字或负文字出现次数最多而设计的算法可有效求解a < 3.52区域的随机3-SAT实例(如Zchaff算法); WALKSAT算法可有效求解a稍高于3.95区域的随机3-SAT实例^{[8, 9]}.

近年来, 使用统计物理中的腔域方法^{[10, 11]}给出了可满足相变点的近似值a_d≃4.267, 且存在阈值a_c≃3.921.称a < a_c为易解SAT区域, 称a_c < a < a_d为难解SAT区域.位于a < a_c区域的可满足实例的解空间形成一个较大的解集簇, 启发式局部搜索算法容易求解; 位于a_c < a < a_d区域的实例, 解空间被分解成许多很小的解集簇, 这些解集簇之间相互分开, 且相距很远, 不可能通过改变有限数目变量的赋值, 使得可满足指派从一个解集簇转到另一个解集簇, 使用局部搜索算法求解时收敛速度很慢, 且求解难度较大^[12].SAT问题的相变现象进一步催生了概率方法在SAT问题中的研究和应用.基于腔域方法设计的求解可满足性问题的信息传播算法, 在难解区域的随机3-SAT实例上非常有效.如:置信传播(belief propagation, 简称BP)算法可有效求解a < 3.95区域的随机3-SAT实例; 调查传播(survey propagation, 简称SP)算法能够有效求解a < 4.26区域的随机3-SAT实例.SP算法是目前求解随机3-SAT问题最为有效的算法, 能够在稍高于线性时间内求解位于难解SAT区域的具备10⁷变量规模的实例, 几乎能够有效求解接近相变点的可满足性实例.对SP算法而言, 难解SAT区域变窄^[12-15].

信息传播算法不仅对固定定义域的随机可满足性问题非常有效(如随机k-SAT问题), 还可用于求解具有增长定义域的随机可满足性问题.许可等人在NP完全问题相变现象的研究中, 提出了具有精确相变和难解实例的RB模型.该模型是典型的具有增长定义域的随机可满足实例产生模型, 被广泛应用于构造难解性问题^[16].随后, 文献[17, 18]分别从理论上和实验上说明:在相变的阈值附近, 求解RB模型随机实例的难度随着问题规模的增加以指数形式增长.文献[19]表明, 基于变量熵的BP算法在求解RB模型实例时非常有效.同样, 信息传播算法可用于求解量化布尔公式(quantified boolean formulae, 简称QBF)的可满足性.殷明浩等人设计了一种基于SP的启发式算法HSPQBF用于求解QBF问题, 实验结果表明, 该算法有效^[20].

近几年, 信息传播算法得到了计算机界学者的重视, 已取得了一些令人瞩目的研究成果.例如, Ravanbakhsh和Greiner通过引入参数, 提出了求解CSP问题的PBP (perturbed BP)算法和PSP (perturbed SP)算法.实验结果表明:在随机CSP实例上, PBP算法效果优于BP算法.同样, 对于随机SAT实例, PSP算法也优于SP算法^[21].特别是在组合优化问题方面, 信息传播算法的应用已取得丰硕的成果, 如MAX-SAT问题、Maximum weight matching问题、Minimum cost network flow问题、图的着色、最大独立集等^[22-27].

尽管信息传播算法求解组合优化问题非常有效, 但仍有两个问题有待在理论上深入研究:

1)信息传播算法是否收敛;

2)如果信息传播算法收敛, 得到的结果是否有效, 即, 近似度如何.如果信息传播不收敛, 这种情况下算法失效.

已经知道:因子图为树型结构的实例, 信息传播算法能够有效收敛; 但对于因子图含有多个回路的图型结构的实例, 信息传播算法不总有效, 常表现为不收敛^[13].

目前, 对信息传播算法的收敛性研究已取得了一些成果.如, 文献[14]中给出了信息传播算法的一般形式, 其本质是变量边缘概率分布的一种近似.文献[28]中分析了有限树型结构上信息传播算法的复杂性.文献[29]中分析了因子图含有一个环的图型结构的实例上信息传播算法的收敛性, 表明信息传播算法能够有效收敛, 且得到的结果是变量边缘分布的有效近似.文献[30]表明:对于任意结构的高斯图模型实例, 信息传播算法能够正确收敛.Tatikonda和Jordan分析了基于计算树的Gibbs测度序列的唯一性与BP算法收敛性之间的关系, 导出了BP算法收敛的一个充分条件^[31].Heskes利用Bethe自由能量最小值的唯一性, 给出了BP算法有唯一固定点的有效条件, 分析了信息传播算法与统计物理中自由能量函数之间的关系, 指出信息传播算法收敛后得到的固定点是自由能量函数的稳定点, 并给出了BP算法收敛到唯一固定点的充分条件, 但未对信息传播算法是否收敛做出理论分析^[32].最近几年, 一些学者根据信息更新函数的特征, 并利用函数的压缩映射原理, 给出了信息传播算法收敛的一些充分条件.如, Ihler直接从BP算法的信息更新函数出发, 获得了传递错误信息的界, 并利用该错误界给出了信息传播算法收敛的一个充分条件^[33].Mooij和Kappen利用向量空间中信息更新函数的压缩映射原理导出了BP算法收敛的充分条件, 并指出:如果更新函数系数矩阵的谱半径严格小于1, 则BP算法收敛^[34].Shi和Schonfeld等人得到了传递错误信息的一个更紧致的界, 利用这个界, 研究了信息传播算法的动态行为, 并给出了算法收敛的一个充分条件^[35].文献[23]表明:BP算法在立方级复杂度内能够有效收敛到最大权重匹配问题的最优解.文献[25]中分析了求解最小代价网络流(MCF)的BP算法的收敛性, 研究结果表明:在亚多项式时间内, BP算法收敛到MCF问题的最优解, 并且这个最优解具有唯一性.与此同时, 文献中还提供了简化了的BP算法, 能够对MCF问题给出多项式时间随机近似方案(FPRAS).文献[36]中给出了信息传播算法收敛的一个概率条件, 研究结果表明:BP算法能够给出最大权重匹配问题和最小代价网络流问题的最优解, 且算法迭代次数高概率的有多项式界.文献[37]中分析了BP算法的复杂性, 提出SBP (stochastic BP)算法, 该算法的消息更新复杂性较BP算法的复杂性低了一个指数级别.

警示传播(warning propagation, 简称WP)算法是最为基础的信息传播算法, 对WP算法的收敛性分析, 有助于其他信息传播算法的收敛性分析.文献[38]中, Feige和Mossel等人分析了WP算法的收敛性, 但主要局限于植入指派的随机可满足实例产生模型${\cal P}_{n,p}^{plant},$即:给定一个指派, 以概率p随机选取满足指派的子句之集构成实例.具体地讲, 在${\cal P}_{n,p}^{plant},$模型中取概率p充分大, 其本质是高概率地最多只有一个赋值满足该实例.文献[39]证明了WP算法在G(n, k, p)模型产生的3CNF公式上的收敛性, 但只证明了p取值很小的一部分实例.目前, 对警示传播算法的性能及收敛性的理论分析仍然不完善.

将WP算法的警示信息取值从{0, 1}松弛为[0, 1], 利用向量空间中压缩函数的性质, 并借鉴文献[40]中的方法, 给出了WP算法收敛的一个充分条件.具体地讲, 对于CNF公式F={a, b, …}, 含有n个变元x₁, x₂, …, x_n.公式F的因子图中边的序列集合记为D={a→i:(a, i)∈E}, 构造|D|×|D|的非负矩阵A, A中的对角元素为${A_{a \to i,a \to i}} = \sqrt 2 \delta $, 其他元素为${A_{a \to i,b \to j}} = {1_{V(j)\backslash a}}(b)\,{1_{V(a)\backslash i}}(j)\,\theta ( - J_j^aJ_j^b)$, δ为常数.如果非负矩阵A的谱半径ρ(A)严格小于1, 则WP算法收敛, 且与初始信息无关.其中, 1_X(x)是标识函数, 若x∈X, 则1_X(x)=1;若x∉X, 则1_X(x)=0.最后, 选取了两组规模不同的随机3-SAT实例进行实验模拟, 结果表明:当子句与变元的比值α < 1.8时, 该判定条件有效.

1 警示传播算法

设F={C₁, C₂, …, C_m}为CNF公式, 含有n个变元x₁, x₂, …, x_n, 用i代表变元x_i.公式F可以用二分图G=(C∪X, E)表示, 称G为因子图.其中, 变元结点集为X={1, 2, …, n}, 子句结点集为C={C₁, C₂, …, C_m}.图G中的边分为两类:实边和虚边.

实边:(C_i, j)∈E⇔子句C_i含正文字x_j;

虚边:(C_i, j)∈E⇔子句C_i含负文字¬x_j.

为了书写方法, 不妨用字母a, b, …分别代表子句C₁, C₂, ….

· V(a):表示出现在子句a中的变元集合, V(a)=V⁺(a)∪V^-(a).其中,

➢ V⁺(a):表示出现在子句a中的正文字对应的变元标识集合;

➢ V^-(a):表示出现在子句a中的负文字对应的变元标识集合;

➢ V(a)\i:V⁺(a)-{i}.

· V(j):表示含变元x_j的子句集合, V(j)=V⁺(j)∪V^-(j).其中,

➢ V⁺(j):表示变元x_j正出现的子句集合;

➢ V^-(j):表示变元x_j负出现的子句集合;

➢ V(j)\a:V⁺(j)-{a}.

$J_j^a$是一个标识参数, 若x_j出现在子句a中, 则$J_j^a = - 1$; 若¬x_j出现在子句a中, 则$J_j^a = 1$.在因子图的每条边(a, i)上, 定义WP算法中的消息传递u_a→i(常称为警示信息).u_a→i表示子句a的可满足性对变元x_i的取值倾向.WP算法的消息更新迭代方程如下:

${u_{a \to i}}(t) = \prod\limits_{j \in V(a)\backslash i} {\theta \left( { - J_j^a\left( {\sum\limits_{b \in V(j)\backslash a} {J_j^b{u_{b \to j}}(t - 1)} } \right)} \right)} $

(1)

t表示迭代次数, θ(x)是截尾函数.如果x≤0, 则θ(x)=0;否则, θ(x)=1.若a中仅包含变元x_i, 则置u_a→i=1.当WP算法收敛时, 根据警示信息固定变元x_i的赋值:

${H_i} = - \sum\limits_{b \in V(i)} {J_i^bu_{b \to i}^*} $

(2)

如果H_i > 0, 则x_i=1;如果H_i < 0, 则x_i=0;否则, x_i暂不赋值.一般地, 将公式(1)改写为如下形式:

${u_{a \to i}}(t) = \prod\limits_{j \in V(a)\backslash i} {\theta \left( {J_j^a\left( {\sum\limits_{b \in {V^ + }(j)\backslash a} {{u_{b \to j}}(t - 1) - } \sum\limits_{b \in {V^ - }(j)\backslash a} {{u_{b \to j}}(t - 1)} } \right)} \right)} $

(3)

设

${h_{j \to a}} = \sum\limits_{b \in {V^ + }(j)\backslash a} {{u_{b \to j}}(t - 1) - } \sum\limits_{b \in {V^ - }(j)\backslash a} {{u_{b \to j}}(t - 1)} $

(4)

称h_j→a为腔域.若变元x_j仅出现在a中, 则置h_j→a=0.求解CNF公式F的WP (warning propagation)算法如下:

Warning Propagation (CNF formula F).

1.构造相应的因子图G(F);

2.给因子图上的所有消息边u_a→i(t=0)随机赋值0或1;

3.重复如下过程, 直到算法收敛(可设置最大迭代步t_max强迫算法结束):

3.1   对G(F)中的边随机排列;

3.2   根据随机边序列, 利用公式(1)更新消息u_a→i;

4.根据H_i计算部分指派ψ, 对公式F进行简化;

5.返回ψ;

WP算法中的收敛, 是指算法某t次迭代得到的信息u_a→i(t)与t-1次的信息u_a→i(t-1)一致.当算法收敛时, 得到警示信息的固定点, 进而可固定部分变元的赋值, 并返回部分指派y; 否则, 算法不收敛, 返回失败.可以看出, 算法的收敛性对算法的性能起关键作用.有如下结论:

定理1^[13].如果公式对应的因子图为树型结构, 那么WP算法收敛.

然而, 在带有环的图型结构的因子图实例上, WP算法不总收敛.目前, 对算法收敛性的理论分析仍然不完善.因此, 有必要对WP算法的收敛性给予理论解释.

2 主要结论及证明

警示传播算法是一种迭代算法, 算法的收敛性影响算法的性能.如果警示传播算法能够有效收敛, 则可以用来固定部分变元的赋值, 从而对原问题进行简化.然而, 对于因子图含有回路的实例, 警示传播算法可能不收敛.下面我们给出警示传播算法收敛的充分条件.

定理2.如果非负方阵M的谱半径ρ(M)严格小于1, 则WP算法收敛, 且与初始信息无关.其中, 矩阵M中的元素为${M_{a \to i,a \to i}} = \sqrt 2 \delta ,{M_{a \to i,b \to j}} = {1_{V(j)\backslash a}}(b)\,{1_{V(a)\backslash i}}(j)\,\theta ( - J_j^aJ_j^b)$, 1_X(x)是标识函数.若x∈X, 则1_X(x)=1;若x∉X, 则1_X(x)=0, δ为非负常数.

在WP算法中, 信息u_a→i取值为0或1.为了证明定理2, 我们将WP算法的消息取值从{0, 1}松弛为[0, 1].具体地讲:对于因子图的每条边(a, i), 定义消息传递η_a→i∈[0, 1], 将警示传播算法的消息更新迭代方程(1)更改为

${\eta _{a \to i}}(t) = (\delta {({\eta _{a \to i}}(t - 1))^2} - \delta ({\eta _{a \to i}}(t - 1)) + 1)\prod\limits_{j \in V(a)\backslash i} {\theta \left( { - J_j^a\left( {\sum\limits_{b \in V(j)\backslash a} {J_j^b{\eta _{b \to j}}(t - 1)} } \right)} \right)} $

(5)

其中, δ为非负常数.

相对于公式(1), 在公式(5)中, 我们构造了项H(η_a→i(t-1))=δ(η_a→i(t-1))²-δ(η_a→i(t-1))+1.注意到:当η_a→i∈{0, 1}时, H(η_a→i(t-1))=1.进一步, 如果公式(5)中消息η_a→i(t=0)的初始取值仍然为0或1, 则公式(5)退化为公式(1).我们的想法是:通过分析公式(5)的收敛性来解释公式(1)的收敛性, 进而给出WP算法收敛的有效条件.

设实例的因子图中边的序列集合为D={a→i:(a, i)∈E}, 定义向量空间V=[0, 1]^|D|, (V, ||·||)是一个赋范空间.设向量x∈V, l₁-范数定义为

$||x|{|_1} = \sum\limits_{i = 1}^{|D|} {|{x_i}|} $

(6)

l_∞-范数定义为

$||x|{|_\infty } = \mathop {\max }\limits_{i \in \{ 1,...,|D|\} } |{x_i}|$

(7)

向量空间V上的范数可诱导一个V上的度量d(v, w)=||v-w||, v, w∈V.显然, 度量空间(V, d)是完备的^[40].设函数f:V→V, 如果存在0≤K < 1, 对于任意的v, w∈V, 使得d(f(v)-f(w))≤Kd(v, w), 则称函数f是压缩函数.如果f:V→V是完备度量空间(V, d)中的压缩函数, 那么对于任意的x∈V, f收敛到唯一固定点x_∞.也即, x, f(x), f²(x), …, x_∞.

引理1.设函数f:V→V, d是V上的度量.对于某个N∈ℕ, 假设f^N是关于d的压缩函数, 那么f收敛到唯一固定点x_∞, 且与初始信息无关.

证明:对于任意的x∈V, 考虑如下N-1个序列:

$\begin{array}{l} x,{f^N}(x),{f^{2N}}(x),...,\\ f(x),{f^{N + 1}}(x),{f^{2N + 1}}(x),...,\\ {\rm{ }} \vdots \\ {f^{N - 1}}(x),{f^{2N - 1}}(x),{f^{3N - 1}}(x),... \end{array}$

因为f^N是关于d的压缩函数, 所以上述每一个序列都收敛到x_∞.因此, 序列x, f(x), f²(x), …一定收敛到x_∞.所以函数f收敛到唯一固定点x_∞, 且与初始信息无关.

设映射A:V→V, (V, ||·||)是一个赋范空间, 诱导A的矩阵范数:

$||A|| = \mathop {\sup }\limits_{v \in V} ||Av||$

(8)

||A||₁-范数定义为

$||A|{|_1} = \mathop {\max }\limits_{j \in \{ 1,...,|D|\} } \sum\limits_{i = 1}^{|D|} {|{A_{ij}}|} $

(9)

||A||_∞-范数定义为

$||A|{|_\infty } = \mathop {\max }\limits_{i \in \{ 1,...,|D|\} } \sum\limits_{j = 1}^{|D|} {|{A_{ij}}|} $

(10)

证毕.

若函数f:V→V可微, 构造雅可比矩阵A.依据中值定理有如下引理.

引理2^[40].设(V, ||·||)是赋范空间, 映射f:V→V可微, 那么对于任意的x, y∈V, 有公式(11)成立:

$||f(x) - f(y)|| \le ||x - y|| \cdot \mathop {\sup }\limits_{z \in V} ||f'(z)||$

(11)

根据引理2, 如果$\mathop {\sup }\limits_{z \in V} ||f'(z)|| < 1$, 那么函数f收敛到固定点x_∞.假设算法中的消息更新并行执行, 对于任意的向量η∈V, 函数f的分量f(η)_a→i定义为f(η)_a→i=η_a→i.自然, 消息更新方程(5)的并行执行过程可以用函数f的迭代过程来刻画.下面从引理2出发导出WP算法收敛的一个充分条件.

为了求解需要, 我们对公式(5)进行参数化.通过观察, 取sinhλ_a→i=η_a→i, 则公式(5)变为如下形式:

$\sinh {\lambda _{a \to i}}(t) = (\delta {\sinh ^2}{\lambda _{a \to i}}(t - 1) - \delta (\sinh {\lambda _{a \to i}}(t - 1)) + 1)\\\prod\limits_{j \in V(a)\backslash i} {\theta \left( { - J_j^a\sum\limits_{b \in V(j)\backslash a} {J_j^b\sinh {\lambda _{b \to j}}(t - 1)} } \right)} $

(12)

由于η_a→i∈[0, 1], 所以${\lambda _{a \to i}} \in [0,\log (1 + \sqrt 2 )]$.定义向量空间$V' = {[0,\log (1 + \sqrt 2 )\;]^{|D|}}$, 设函数h:V'→V', 对于任意的向量λ∈V', 函数h的分量h(λ)_a→i定义为h(λ)_a→i=l_a→i.设i, j∈V(a), i≠j, b∈V(j)\a.令:

$G({\lambda _{a \to i}}(t - 1)) = \delta {\sinh ^2}{\lambda _{a \to i}}(t - 1) - \delta (\sinh {\lambda _{a \to i}}F({\lambda _{c \to k}}(t - 1)) =\\ \prod\limits_{k \in V(a)\backslash i,j} {\theta \left( { - J_k^a\sum\limits_{c \in V(k)\backslash a} {J_k^c\sinh {\lambda _{c \to k}}(t - 1)} } \right)} (t - 1) + 1$

对函数h求导, 根据公式(12), 我们有:

$(h'(\lambda ))_{b \to j}^{a \to i} = \frac{{\partial {\lambda _{a \to i}}(t)}}{{\partial {\lambda _{b \to j}}(t - 1)}} = \frac{{G({\lambda _{a \to i}}(t - 1)) \cdot F({\lambda _{c \to k}}(t - 1))}}{{\sqrt {1 + {{\sinh }^2}{\lambda _{a \to i}}(t)} }}\theta ( - J_j^aJ_j^b\cosh {\lambda _{b \to j}}(t - 1))$

(13)

注意, 上述是关于变量λ_b→j的求导, 对于所有的子句c∈V(j)\{a, b}, λ_c→j在公式(12)中均被看作常量.同样, 对于所有的k∈V(a)\{i, j}, d∈V(k)\a, λ_d→k也被看做常量处理.令$y({\lambda _{b \to j}}) = - J_j^a\sum\limits_{b \in V(j)\backslash a} {J_j^b\sinh {\lambda _{b \to j}}} $, 考虑到截尾函数θ(·)在公式(12)中的实际作用, 先通过对函数y(λ_b→j)进行求导, 再进行截尾函数θ(·)处理.因此, 公式(12)通过求导后, 截尾函数为$\theta ( - J_j^aJ_j^b\cosh {\lambda _{b \to j}}(t - 1)).$

又因为coshλ_b→j(t-1) > 0, 所以$\theta ( - J_j^aJ_j^b\cosh {\lambda _{b \to j}}(t - 1)) = \theta ( - J_j^aJ_j^b)$.注意到:

$\mathop {\sup }\limits_{\lambda \in {{[0,\log (1 + \sqrt 2 )]}^{|D|}}} \frac{{G({\lambda _{a \to i}}(t - 1)) \cdot F({\lambda _{c \to k}}(t - 1))}}{{\sqrt {1 + {{\sinh }^2}{\lambda _{a \to i}}(t)} }} = 1$

(14)

因此,

$|(h'(\lambda ))_{b \to j}^{a \to i}| = \left| {\frac{{\partial {\lambda _{a \to i}}(t)}}{{\partial {\lambda _{b \to j}}(t - 1)}}} \right| \le \;\theta ( - J_j^aJ_j^b)$

(15)

而

$\begin{array}{l} h'(\lambda )_{a \to i}^{a \to i} = \frac{{\partial {\lambda _{a \to i}}(t)}}{{\partial {\lambda _{a \to i}}(t - 1)}}\\ {\rm{ }} = \frac{{\delta \cosh {\lambda _{a \to i}}(t - 1)(2\sinh {\lambda _{a \to i}}(t - 1) - 1)\prod\limits_{j \in V(a)\backslash i} {\theta \left( { - J_j^a\sum\limits_{b \in V(j)\backslash a} {J_j^b\sinh {\lambda _{b \to j}}(t - 1)} } \right)} }}{{\sqrt {1 + {{\sinh }^2}{\lambda _{a \to i}}(t)} }}\\ {\rm{ }} \le \delta \cosh {\lambda _{a \to i}}(t - 1)(2\sinh {\lambda _{a \to i}}(t - 1) - 1) \end{array}$

(16)

由于sinhλ_a→i=η_a→i∈[0, 1], 因此$h(\lambda )_{a \to i}^{a \to i} = \frac{{\partial {\lambda _{a \to i}}(t)}}{{\partial {\lambda _{a \to i}}(t - 1)}} \le \delta \cosh {\lambda _{a \to i}}(t - 1).$

又因为$\cosh {\lambda _{a \to i}}(t - 1) = \sqrt {1 + {{\sinh }^2}{\lambda _{a \to i}}(t - 1)} \le \sqrt 2 ,$

所以,

$|(h'(\lambda ))_{a \to i}^{a \to i}| = \left| {\frac{{\partial {\lambda _{a \to i}}(t)}}{{\partial {\lambda _{a \to i}}(t - 1)}}} \right| \le \sqrt 2 \delta $

(17)

根据代数知识, 利用函数h的导数构造雅可比矩阵A.由压缩函数的性质及引理2可知:如果A的范数严格小于1, 则h为压缩函数, 进而收敛到唯一固定点, 且与初始信息无关.有如下结论:

定理3.对于||A||₁-范数, 如果

$\mathop {\max }\limits_b \mathop {\max }\limits_{j \in V(b)} \sum\limits_{a \in V(j)\backslash b} {\theta ( - J_j^aJ_j^b)} < 1 - \sqrt 2 \delta $

(18)

则公式(5)在||A||₁-范数下收敛, 且与初始信息无关.

证明:根据公式(9)、公式(15)、公式(17), 有:

$||A|{|_1}\; = \;||h'(\lambda )|{|_1} = \mathop {\max }\limits_{b \to j} \sum\limits_{a \to i} {\left| {\frac{{\partial {\lambda _{a \to i}}}}{{\partial {\lambda _{b \to j}}}}} \right|} \le \sqrt 2 \delta + \mathop {\max }\limits_b \mathop {\max }\limits_{j \in V(b)} \sum\limits_{a \in V(j)\backslash b} {\theta ( - J_j^aJ_j^b)} $

(19)

由压缩函数的性质和引理2可知:如果公式(20)成立, 则函数h收敛到唯一固定点, 且与初始信息无关:

$\mathop {\max }\limits_b \mathop {\max }\limits_{j \in V(b)} \sum\limits_{a \in V(j)\backslash b} {\theta ( - J_j^aJ_j^b)} < 1 - \sqrt 2 \delta $

(20)

根据代数知识, 函数f收敛到某一固定点, 且与初始信息无关.所以, 公式(5)在||A||₁-范数下收敛, 且与初始信息无关.

同理可证如下结论.

定理4.对于||A||_∞-范数, 如果:

$\mathop {\max }\limits_a \mathop {\max }\limits_{j \in V(a)} \sum\limits_{b \in V(j)\backslash a} {\theta ( - J_j^aJ_j^b)} < 1 - \sqrt 2 \delta $

(21)

则公式(5)在||A||_∞-范数下收敛, 且与初始信息无关.

证明:根据公式(10)、公式(15)、公式(17), 有:

$||A|{|_\infty }\; = \;||h'(\lambda )|{|_\infty } = \mathop {\max }\limits_{a \to i} \sum\limits_{b \to j} {\left| {\frac{{\partial {\lambda _{a \to i}}}}{{\partial {\lambda _{b \to j}}}}} \right|} \le \sqrt 2 \delta + \mathop {\max }\limits_a \mathop {\max }\limits_{j \in V(a)} \sum\limits_{b \in V(j)\backslash a} {\theta ( - J_j^aJ_j^b)} $

(22)

由压缩函数的性质和引理2可知:如果公式(23)成立, 则函数h收敛到唯一固定点, 且与初始信息无关:

$\mathop {\max }\limits_a \mathop {\max }\limits_{j \in V(a)} \sum\limits_{b \in V(j)\backslash a} {\theta ( - J_j^aJ_j^b)} < 1 - \sqrt 2 \delta $

(23)

根据代数知识, 函数f收敛到某一固定点, 且与初始信息无关.所以, 则公式(5)在||A||_∞-范数下收敛, 且与初始信息无关.

设M是一个方阵, γ₁, …, γ_n是M的特征值, 称γ(M)={γ₁, …, γ_n}为方阵M的谱, ρ(M)=sup{|γ_i|:γ_i∈γ(M)}为方阵M的谱半径.

定理5.设函数h:V'→V'可微, 进一步假设|h'(x)|≤M, M是不依赖于x的非负矩阵.如果ρ(M) < 1, 那么对于任意x∈V', h收敛到唯一固定点x_∞.也即, x, h(x), h²(x), …, x_∞(注意:|h'(x)|≤M表示由h'(x)构成的雅可比矩阵中的每个元素的绝对值不大于M中对应元素的值).

证明:对于任意的x∈V', 假设n=1, 2, …, 我们有:

$||({h^n})'(x)|| = \left\| {\prod\limits_{i = 1}^n {h'({h^{i - 1}}(x))} } \right\| = {\left\| {\left| {\prod\limits_{i = 1}^n {h'({h^{i - 1}}(x))} } \right|} \right\|_1} \le {\left\| {\prod\limits_{i = 1}^n {|h'({h^{i - 1}}(x))} |} \right\|_1} \le {\left\| {\prod\limits_{i = 1}^n M } \right\|_1} = ||{M^n}|{|_1}$

(24)

由Gelfand定理^[40]:$\mathop {\lim }\limits_{n \to \infty } ||{M^n}||_1^{1/n}\; = \;\rho (M)$.我们选取ε > 0, 使得ρ(M)+ε < 1.对于某个N∈ℕ, ||M^N||₁=(ρ(M)+ε)^N < 1.因此, ||(h^N)'(x)|| < 1.所以, ||(h^N)(x)||是压缩函数.

由引理1可知:函数h(x)收敛到唯一固定点x_∞, 且与初始信息无关.

对于任意的i, j∈V(a), i≠j, b∈V(j)\a, 定义|D|×|D|的非负方阵M.M中的元素:

${M_{a \to i,a \to i}} = \sqrt 2 \delta ,{M_{a \to i,b \to j}} = {1_{V(j)\backslash a}}(b)\,{1_{V(a)\backslash i}}(j)\,\theta ( - J_j^aJ_j^b).$

其中, 1_X(x)是标识函数, 若x∈X, 则1_X(x)=1;若x∉X, 则1_X(x)=0, δ为非负常数.根据公式(15)、公式(17)和定理5, 我们有如下结论:

定理6.如果方阵M的谱半径ρ(M)严格小于1, 则迭代方程(5)式收敛, 且与初始信息无关.

证明:根据公式(15)、公式(17)和定理5, 函数h收敛到唯一固定点, 且与初始信息无关.进一步地, 函数f收敛到唯一固定点, 且与初始信息无关.所以, 方程(5)收敛, 且与初始信息无关.

假设M的特征值γ, 根据代数知识, 我们有:|γ|·||x||=||γx||=||Mx||≤||M||·||x||.因此, 定理6的条件弱于定理3和定理4.

设向量空间V″={0, 1}^|D|⊂V, 定义函数g:V″→V″.假设WP算法中的信息更新是并行执行, 对于任意的向量u∈V″, 函数g的分量g(u)_a→i定义为g(u)_a→i=u_a→i.自然, WP算法的执行过程可用函数g的迭代过程来刻画.根据函数f和g的构造, 显然, 对于任意的x∈V″, 有f(x)=g(x).

引理3.如果函数f是一个压缩映射, 那么函数g也为压缩映射.

证明:由于f是一个压缩映射, 对于任意的x∈V, 存在唯一固定点x_∞, 有x, f(x), f²(x), …, x_∞, 也即:函数f收敛到x_∞, 且与初始信息无关.因此, 对于任意的y∈V″, 仍然有y, f(y), f²(y), …, x_∞.而当y∈V″时, 有f(y)=g(y).进而有y, g(y), g²(y), …, x_∞.因此, 函数g为压缩映射, 收敛到固定点x_∞, 且与初始信息无关.

根据定理6和引理3, 至此, 定理2得证.

注意到:定理2中非负矩阵的元素由公式的结构决定.因此, 参数δ的值取决于公式的结构.显然, 对于定理1中的结论, 容易使用定理2来证明.为了方便理解, 下面举例说明定理2的使用.

例1:设CNF公式F₁={(x₁∨x₂), (x₂∨x₃), …, (x_n-1∨x_n), (x_n∨x₁)}, 其因子图含有一个环.因此, 矩阵M的谱半径$\rho (M) = \sqrt 2 \delta $.若取$0 < \delta < \frac{1}{{\sqrt 2 }}$, 则ρ(M) < 1.根据定理2, 则WP算法在公式F₁上收敛.

同样, WP算法在公式F₂={(¬x₁∨¬x₂), (¬x₂∨¬x₃), …, (¬x_n-1∨¬x_n), (¬x_n∨¬x₁)}上也收敛.

3 数值实验及分析

在数值实验中, 使用模型G(n, 3, m)产生随机3-SAT实例, 其中, n表示变元个数, m表示子句个数, 子句大小为3.模型G(n, 3, m)按如下方式产生随机3-SAT实例:随机均匀地在可能的子句(这样的子句共有${2^3} \cdot \left( {\begin{array}{*{20}{c}} n\\ 3 \end{array}} \right)$个)中选择m=αn个子句, 被选中的子句之集构成实例.在算法中, 取最大迭代次数t_max=10³.

我们选取G(n, 3, m)随机实例进行模拟实验.图 1是WP算法分别在变元规模n=20, n=40, n=60上的运行情况.图 1(a)是算法的收敛性随着约束参数α的变化情况, 图 1(b)是算法的计算时间随着约束参数α的变化情况.图中每个数据点由G(n, 3, m)模型产生的100个随机实例构成.可以看出:当参数α < 3.5时, WP算法高概率收敛; 当α > 4.5时, 算法高概率不收敛.在α > ≈3.5点, WP算法的收敛性发生了突变.同样, 当参数α < 3.5时, WP算法的运行时间较小; 当α > 3.5时, WP算法的运行时间突然增大.

对于随机3-SAT实例, α是重要的参数, 它影响实例的可满足性和判定难度.存在可满足性的相变点a_d, 当α < a_d时, 实例高概率可满足; 当α > a_d时, 实例高概率不可满足.同时, 存在第2相变点a_c, a < a_c为易解SAT区域, a > a_c为难解SAT区域.大量的材料表明:a_d至少为3.52, 至多为4.506, a_c≃3.921.注意到, WP算法在3.52 < α < 4.506区域的收敛情况不能从概率上给出确切的回答.因此, 警示传播算法能够有效判定的随机3-SAT实例最多为α < 3.52区域.从图 1(a)中可以看出:随着问题规模n的增大, WP算法在G(n, 3, m)实例上收敛的概率趋向于一个稳定值, 即, 收敛的概率曲线随着n的增大几乎趋于不变.我们知道, 因子图的结构直接影响WP算法的收敛性.事实上, 当参数α > a_c时, G(n, 3, m)模型产生的随机实例的因子图结构比较复杂, 求解难度较大.从图 1(b)中我们也可以看出:大约在α > 3.5处, 算法的运行时间急剧上升.

为了验证定理2中的判定条件, 我们分别取了n=20和n=60的G(n, 3, m)实例数据集进行模拟实验.从图 2中可以看出:当α < 1.8时, 定理2能够判定WP算法的收敛性.事实上, 文献[39]已经表明:当约束参数α较小时, G(n, 3, m)实例的因子图的结构内部呈现树形状.因此, 根据定理2, 容易判定WP算法收敛.随着α的增大, 因子图结构中开始出现圈, 并且圈的数量随着α的进一步增大而急剧增加.由于定理2中的非负矩阵的元素构成与因子图的结构有关, 当因子图的结构变得复杂时, 非负矩阵的谱半径不再小于1.这时, 定理2中的判定条件失效.如图 2所示:当α > 1.8时, 定理2中的判定条件逐渐失效.

如果一个因子图包含环, 且只含正文字或只含负文字的实例F, 根据定理2容易判断, 则WP算法在F上收敛.我们分别选取3组不同变元规模的由G(n, 3, m)模型产生的随机3SAT可满足实例进行实验模拟, 这里取α=2.0.这些实例只含正文字或只含负文字, 每组数据有200个实例(其中, 有100个只包含正文字的实例和100个只包含负文字的实例).实验结果表明:在这种随机实例上, WP算法收敛.表 1是WP算法的收敛情况统计表.用符号“+”表示只含正文字的实例; 用符号“-”表示只含负文字的实例.Var_num表示实例的变元规模, Instance表示实例集, Convergence表示WP算法收敛的概率.

4 结束语

本文分析了WP算法的收敛性, 通过引入参数, 将WP算法的信息取值从{0, 1}松弛为[0, 1], 利用向量空间中压缩函数的性质, 给出了WP算法收敛的一个充分条件.大量的实验结果表明, 该方法有效.同时, 该方法提供了理解WP算法性能的理论依据, 并为分析其他信息传播算法的收敛性奠定了理论基础.进一步工作是:给出WP算法更弱的收敛条件.