有容量集合覆盖选址问题的降阶回溯算法
2020-04-11尚春剑宁爱兵彭大江张惠珍
尚春剑，宁爱兵，彭大江，张惠珍
(上海理工大学 管理学院，上海 200093)
1 引 言
集合覆盖问题(Location Set Covering Problem，简称LSCP)最早由Toregas[1]等人提出，是运筹学领域中一个经典的NP-Hard问题[2]，除非P=NP，否则该问题不存在多项式时间的精确算法.集合覆盖问题要求设施以最少的开设数量来覆盖所有需求点，在实际生活中有许多应用场合，如设施选址、车辆调度、资源分配、电路设计、网络安全等领域.该问题的重要性以及应用的广泛性引起了学者们广泛的关注，目前针对该类问题研究的算法主要分为三类：第一类是精确算法，主要基于branch-and-bound算法和branch-and-cut算法[3，4]，20世纪70年代，Toregas和ReVelle[5，6]提出了精确求解LSCP的行和列简化算法，这类精确算法虽然能够得到最优解，但是没有研究问题本身的数学性质，求解速度相对较慢.第二类是近似算法，HOCHBAUM[7]提出了一个求解LSCP时间复杂度为O(n3)的近似算法；GROSSMAN[8]等人对9种不同的LSCP近似算法进行了比较研究，包括几种贪婪变量、分数松弛、随机化算法和神经网络算法，近似算法虽然可以在多项式时间内得到一个常数近似比的解，但缺点是无法取得最优解.第三类是启发式算法，Vasko和Wilson[9]提出了求解LSCP的启发式算法；Alminana和Pastor[10]提出了求解LSCP改进形式的代理启发式算法，目前有许多学者致力于该问题的启发式算法的研究[11-14]，启发式算法的求解速度相对较快，但是缺点是容易陷入局部最优解.
本文将集合覆盖问题的模型应用到有容量设施选址问题中，约束条件中添加了设施容量限制和需求点需求量约束，其中需求点的需求量可以分割，即一个需求点可由多个设施提供服务，集合覆盖问题是有容量集合覆盖问题的一个特殊子问题，求解有容量集合覆盖问题的难度大于集合覆盖问题，Current和Storbeck[15]研究了此类有容量约束的LSCP.
由于每个设施的容量有限，同时要满足每个需求点的需求量，在精确算法领域中，有容量约束的LSCP问题研究难度较大，相关研究相对较少.本文针对上述现有算法的不足之处，首先研究有容量集合覆盖问题具有的数学性质，这些数学性质能缩小问题规模，降低算法时间复杂度；然后利用上下界子算法对问题剪枝，加快问题求解速度；之后设计的分配子算法解决了传统精确算法在有容量设施选址问题中的难点；最后结合数学性质、上界子算法、下界子算法、分配子算法，设计出一种能够降低问题规模从而加快问题求解速度并且能求出问题最优解的降阶回溯算法.
2 问题描述
2.1 数学符号及解释
以下为数学模型中涉及的数学符号：
ei：表示第i个需求点；
fj：表示第j个设施；
m：表示需求点的个数；
n：表示设施的个数；
E={ei|i=1，2，…，m}：表示所有需求点的集合；
F={fj|j=1，2，…，n}：表示设施的集合；
di：表示第i个需求点的需求量；
cj：表示第j个设施的容量；
xj∈{0，1}：若设施fj开设，则xj=1，否则xj=0；
yi，j：表示需求点ei的需求中被分配给设施fj的部分；
A(j)：表示第j个设施fj可服务的需求点集合；
B(i)：表示可以覆盖第i个需求点ei的设施集合。
以下为后文设计的算法中涉及的数学符号：
aj：若设施fj的容量有剩余，即rj> 0，则aj=1，否则aj=0；
gi：若需求点ei的需求量完全满足，即ki=0，则gi=1，否则gi=0；
deg(ei)：表示需求点ei的度，即需求点ei可被deg(ei)个设施服务；
deg(fj)：表示设施fj的度，即设施fj可服务deg(fj)个需求点；
N(vi)：表示二分图G中顶点vi的邻点集；
b：算法在某一状态下的目标函数值下界，为局部变量；
u：算法在某一状态下的目标函数值上界，为全局变量；
F1：一定开设的设施集合，若F1中任一设施不开设则不能取得最优解，为全局变量；
F0：一定不开设的设施集合，若F0中任一设施开设则不能取得最优解，为全局变量；
F5=FF1F0：暂时未确定是否开设的设施集合，为全局变量；
Ftemp：后文算法中将某些设施临时放入集合Ftemp中，为局部变量；
FF1：F5中的设施在回溯算法分情况时假设开设的设施放入集合FF1中，为全局变量；
FF0：F5中的设施在回溯算法分情况时假设不开设的设施放入集合FF0中，为全局变量；
FF5=F5FF1FF0：回溯算法分情况时暂时未处理的设施放入集合FF5中，为全局变量；
Fbest：算法在当前状态下已知最好的目标函数值对应的开设设施集合，为全局变量；
best_q：回溯算法在当前状态下已知最好的目标函数值，best_q= |Fbest|，为全局变量；
cur_n：当前累计开设设施数，为局部变量；
cur_i：回溯算法的当前搜索层数，为局部变量。
2.2 数学模型
本文用二分图表示有容量集合覆盖问题：将E和F中的元素分别作为二分图的两个顶点子集E和F，若E中的ei可辐射F中的fj或F中的fj可辐射E中的ei，即ei与fj之间存在路径，则将ei与fj连线，否则不连线，将所有路径放入集合X中.处理后得到一个图G= (V，X)，其中V=E∪F，显然图G是一个二分图.
该有容量集合覆盖选址问题的具体模型如下[11]：
(1)
目标函数(1)表示在满足全部约束条件下，最小化开设设施的数量；约束(2)表示对于任意一个需求点ei的全部需求均被完全满足；约束(3)表示开设的设施所服务的需求点的需求量总和不能超过所开设设施本身的容量；约束(4)表示设施是否开设的决策变量xj取值为0或1，xj=0表示设施fj不开设，xj=1表示设施fj开设；约束(5)中yi，j表示需求点ei的需求中被分配给设施fj的部分，取值为di与cj的最小值.该问题已证明是NP-Hard问题[2]，除非P=NP，否则该问题不存在多项式时间的精确算法.
3 数学性质及子算法设计
3.1 数学性质
性质1.若图G是非连通图，则可对图G的连通分量分别看作单独的子问题，分别进行求解.
证明：由于子问题之间不存在通路，求解时相互独立，互不影响，因此可以分别求解.
性质2.对于每个设施的容量rj和需求点的需求量ki，若全体rj和ki之间存在公约数，则可将当前问题中的容量和需求量同时约去该公约数，降低了问题求解过程中计算的复杂程度.
证明：由于对设施的容量和需求点的需求量同时进行约分，对于某个设施服务某个需求点的能力不产生影响，问题转化为其等价问题，数值更小便于求解.
证明：因为fj在图G中所起到的覆盖作用完全可以由fh来代替，并且fh的剩余容量能满足N(fh)中所有需求点的需求量.
证明：反证法，由于满足该x个需求点的设施总容量恰好等于x个需求点的总需求量，若覆盖x个需求点的设施中有设施未开设，则这x个需求点的需求量不能完全被满足，因此fj∈B(i)的设施均开设，且服务于这x个需求点.
性质7.对于任意一个已确定开设的设施fj，且rj≥∑ei∈A(j)ki，则ei∈A(j)的需求点均由设施fj服务.
证明：由于目标函数是最小化开设设施数，fj开设且剩余容量大于等于其所覆盖的需求点的需求量之和，则ei∈A(j)的需求点均由设施fj服务使得该资源尽量更多的被利用，若其中有需求点不由设施fj服务，则浪费设施fj的容量而造成要占用其他设施的容量.
性质8.若需求点ei满足deg(ei)=1，则ei对应的B(i)中唯一的设施fj一定开设，放入F1中，且服务ei.
证明：由于问题要求全覆盖，因此需求点ei必须被服务，若需求点ei只能被一个设施fj服务，则设施fj一定开设且为需求点ei提供服务.
性质9.在问题以及在问题处理过程中出现的子问题中，若对于任意一个已确定开设的设施fj，满足aj=1且deg(fj)=1，则fj一定服务于其当前对应的A(j)中唯一的需求点ei.
证明：由于设施fj已开设且有剩余容量，目标函数要求开设的设施数量最少，则设施fj剩余的容量应当服务于当前对应的A(j)中唯一的需求点ei，否则会浪费设施fj剩余的容量.
性质10.在假设某设施fj开设的情况下，此时若下界b大于之前的上界u，则设施fj一定不开设，其中上下界求解算法详见3.3和4.4.
证明：新的下界大于之前的上界，说明设施fj开设则不可能获得最优解，因此应关闭设施fj.
性质11.在假设某设施fj不开设的情况下，此时若下界b大于之前的上界u，则设施fj一定开设.
证明：新的下界大于之前的上界，说明设施fj关闭不可能获得最优解，因此应开设设施fj.
性质12.在假设某设施fj不开设的情况下，若FF5中所有设施均开设，但分配子算法无解，则设施fj一定开设.
证明：若设施fj不开设，则该问题无解，因此设施fj一定开设.
性质13.算法执行过程中，若ki=0，则删除对应的需求点ei及其邻接边，更新图G；若rj=0，则删除对应的设施点fj及其邻接边，更新图G.
证明：显然，当需求点ei完全被满足就可以删除；同样，当设施点fj剩余容量为零时也可以删除.
3.2 分配子算法
在集合覆盖选址问题中，首先要从许多候选设施点中选取开设的设施，然后将需求点分配到所选取的设施上，经典的集合覆盖问题选址中，需求点的不同分配顺序不影响目标函数值，而本文研究的有容量集合覆盖选址问题，由于设施和需求点都有容量的限制，因此每个需求点的不同分配顺序会使得所得到的目标函数值不同.于是在算法设计中，有容量集合覆盖选址问题不仅有设施的选择过程，还有分配过程.
对于本文研究的问题，在开设设施已定的情况下，只要通过合理分配方法在多项式时间内找到一个可行解即可.本文通过下面的分配子算法对需求点进行分配，该子算法的具体内容如下：
Step 1.初始化所有需求点，令所有gi=0；
Step 2.根据性质4判断开设的设施集合是否不可行，若不可行，则该问题无解，分配子算法结束，否则执行下一步；
Step 3.判断问题是否满足数学性质5，此时得到一个解，不需要执行分配子算法的后面步骤，将所有需求点标记gi=1，则已全部得到满足的需求点个数dis=m，分配子算法结束；若不满足数学性质5，执行下一步；
Step 4.将需求点和设施分别按照其度的大小升序排列标号，即度越小的结点越先处理.将满足数学性质7的需求点标记为gi=1，若dis=m，分配子算法结束；若dis Step 5.下面的分配采用求最大流的算法进行分配，具体操作步骤如下：设立一个超级源点和超级汇点，将超级源点连向所有的需求点，每条边的流量为所连需求点的需求量；将所有设施连向超级汇点，每条边的流量为设施的容量.当ei∈A(j)时，将需求点ei连接到设施fj，并将该边的容量设为min{di，cj}.然后按照最大流算法求解：计算从超级源点到超级汇点的最大流[15]，将超级源点到需求点的弧是饱和弧的需求点标记为gi=1，若已全部得到满足的需求点个数dis=m，此时分配可行，分配子算法结束；若dis 证明：当问题所开设设施确定之后，问题的关键在于在多项式时间内怎样把已开设设施的容量合理地分配到所有需求点且使所有需求点的容量得到满足.分配子算法中用到的数学性质已经在前文证明过，若最大流中的超级源点到所有的需求点都是饱和弧，由最大流算法可知，每一个需求点的需求量都能得到满足；若最大流中的超级源点到某一个需求点不是饱和弧，由最大流算法可知，至少有一个需求点的需求量不能得到满足，因此此时不存在可行解.由文献[15-21]可知，最大流问题能在多项式时间内解决，其中文献[15]中的最大流算法时间复杂度为O(mnlog(n2/ m))，因此该分配子算法可以在多项式时间内结束. 本文设计的下界子算法的具体步骤如下： Step 1.初始化b=|F1|，Ftemp={ }，计算开设的设施容量之和∑fj∈F1∪Ftemprj； Step 3.在图G中，将设施按照其容量降序排列，开设容量最大的设施fj放入设施集合Ftemp中，跳到Step 2. 证明：该算法选出的设施集合Ftemp中任一设施的容量都大于等于未选中设施的容量，若选出的开设的设施个数小于|Ftemp|，那么此时开设设施的容量之和一定小于所有需求点的需求量之和，所以此时没有可行解，因此开设的设施数一定大于等于|F1∪Ftemp|. 事实上，任何一个可行解对应的目标函数值均能作为该问题的上界.本文先利用前面所介绍的数学性质，将问题进行降阶处理，之后执行如下的上界子算法求出上界： 将所有需求点按照其需求量降序排列，将设施按照其容量降序排列；依次针对每个需求点ei，检查ei的邻接结点N(ei)是否存在已开设的设施，若存在，则N(ei)中开设的设施为ei服务，若不存在，选取ei邻接结点中容量最大的设施开设并服务需求点ei，直到需求点ei的需求量完全被满足，每个需求点的需求量均被满足，算法结束. Step 1.初始化Ftemp={ }，i=1，将所有需求点按照其需求量降序排列，将设施按照其容量降序排列； Step 2.若ki>0，则分以下三种情况处理： 情况1：若N(ei)中不存在开设的设施，此时按照下面步骤处理： 1)若N(ei)中的设施剩余容量之和大于等于ei的剩余需求量，则选取N(ei)中剩余容量最大的设施fmax开设并为ei提供服务，Ftemp=Ftemp∪{fmax}；若N(ei)中的设施剩余容量之和小于ei的剩余需求量，则本次分配不可行，令上界u=+∞，结束上界子算法； 2)若fmax的剩余容量大于等于ei的剩余需求量，那么此时ei的所有需求都能得到满足，情况1的处理结束； 3)若fmax的剩余容量小于ei的剩余需求量，那么此时ei的部分需求不能得到满足，跳到情况1的步骤(1). 情况2：若N(ei)中存在开设的设施且N(ei)中开设的设施剩余容量之和大于等于ei的剩余需求量，此时按照下面步骤处理： 1)选取N(ei)中已开设且剩余容量最大的设施fmax为ei提供服务，Ftemp=Ftemp∪{fmax}； 2)若fmax的剩余容量大于等于ei的剩余需求量，那么此时ei的所有需求都能得到满足，情况2的处理结束； 3)若fmax的剩余容量小于ei的剩余需求量，那么此时ei的部分需求不能得到满足，跳到情况2的步骤(1)； 情况3：若N(ei)中存在开设的设施且N(ei)中开设的设施剩余容量之和小于ei的剩余需求量，此时N(ei)中已开设设施的全部剩余容量为ei提供服务，(Ftemp=Ftemp∪{fj}(fj∈N(ei)且xj=1)，把N(ei)中已开设的设施移除，然后跳到情况1处理. Setp 3.若ki= 0，i=i+1； Step 4.若i>m或F5=∅，输出开设设施数|F1∪Ftemp|作为上界u，本上界子算法结束；否则跳到Step 2. 证明：若本子算法求出的解是可行解，那么最优解肯定小于等于可行解对应的目标函数值；若本子算法全部设施开设都不能满足需求点的需求，此时u=+∞，可以作为目标函数的上界. 降阶回溯算法包括降阶算法和回溯算法两个部分，降阶算法通过结合前文介绍的数学性质对问题进行降阶，从而缩小问题的规模，降低求解的难度；回溯算法采用深度优先的方法搜索问题的解空间，从根结点出发，对每一个结点判断该情况下其理论下界b是否小于上界u，若满足则继续向下深度优先搜索，否则进行剪枝，逐级向上回溯. 降阶算法步骤如下： Step 1.初始化cur_i=1，cur_n=0，best_q=+∞，u=0，F1=F0={ }，F5=F={fj|j=1，2，…，n}； Step 2.根据3.1介绍的数学性质确定一定开设和一定不开设的设施，对问题进行降阶处理，若某设施fj一定开设，则F1=F1∪{fj}，F5=F5 Step 3.根据3.4的上界子算法计算问题上界u； Step 4.对F5中的每一个设施fj，若满足性质10，即当fj开设时有下界b>u，则fj一定不开设，令F0=F0∪{fj}，F5=F5 执行4.1介绍的降阶算法对问题规模进行降阶处理后，调用下面的回溯子算法backtrack(1).回溯算法对设施集合FF5=F5中的每一个设施分为以下2种情况进行处理： 1)若设施fj开设，并搜索对应的左子树； 2)若设施fj不开设，搜索对应的右子树. 回溯子算法backtrack(cur_i)的详细步骤如下： Step 1.初始化cur_i=1，cur_n=0，best_q=+∞，u=0，FF1=FF0={}，FF5=F5； Step 2.当cur_i>|FF5|或|FF5|=0，说明搜索到达叶子结点，此时调用分配子算法，对需求点进行分配，此时若dis=m，则得到一个可行解.若对应目标函数值Q Step 3.情况1：开设设施fcur_i，在此情况下计算下界b，若b≤best_q，表明此时可以开设设施fcur_i，令FF1=FF1∪{fcur_i}，FF5=FF5 Step 4.算法跳到搜索树的上一层之前，恢复FF1和FF5到Step 3的初始状态，FF1=FF1({fcur_i}∪FF1_temp1)，FF5=FF5∪({fcur_i}∪FF1_temp1∪FF1_temp0)，FF0=FF0FF1_temp0； Step 5.情况2：不开设设施fcur_i，在此情况下计算下界b，若b≤best_q，则说明不开设fcur_i可能取得最优解，令FF0=FF0∪{fcur_i}，FF5=FF5 Step 6.算法跳到搜索树的上一层之前，恢复FF0和FF5到Step 5的初始状态，FF0=FF0({fcur_i}∪FF0_temp0)，FF5=FF5∪({fcur_i}∪FF0_temp1∪FF0_temp0)，FF1=FF1FF0_temp1. 算法结束后，当前的最优目标函数值best_q和对应开设的设施Fbest就是整个问题的一个最优解. 有了前面的子算法，就可以执行下面的主算法： Step 1.先调用降阶子算法； Step 2.调用上下界子算法，计算问题的上界u和下界b； Step 3.调用回溯子算法backtrack(1). 本文研究的集合覆盖选址问题规模即设施个数为n，该算法在搜索空间中最多产生2n个叶子结点，在降阶算法被调用后，问题规模减少为|F5|，令k=|F5|，因此算法在最坏情况下的时间复杂度为O(2k)，k≤n. 许多学者提出不同算法针对集合覆盖选址问题进行求解[22]，这些算法主要分为精确算法、近似算法和启发式算法.近似算法与启发式算法的优点在于求解速度快，但是一般无法得到问题的最优解，不能保证解的质量.本文设计的精确算法，首先保证了所求的解为最优解；其次，相对于其他一般的精确算法而言，本文通过研究问题的数学性质，这些数学性质能对问题进行降阶，降低了问题规模，设计的上下界算法对搜索树进行合理剪枝，缩小了问题的解空间，只对部分解子树搜索，使得时间复杂度由O(2n)降低为O(2k)，其中k=|F5|且k≤n. 为了更清晰的说明本文算法，下面给出一个示例对算法执行的整个过程进行说明. 示例1:如图1所示，现有6个需求点ei，7个设施fj，若设施fj能为需求点ei提供服务，则用实线连接.现从中选出若干个设施服务于需求点，使得每个需求点的需求量均得到满足，求所选取设施个数的最小值Q. 对示例1的分析过程可描述如下： 图1 设施服务需求点服务二分图Fig.1 Bipartite graph of facility service demand point service 经过分析，图1为非连通图，根据性质1，可将示例1分成两个单独的子问题分别进行求解.{j1、j6、i1、i6}可构成子问题1，{j2、j3、j4、j5、j7、i2、i3、i4、i5}可构成子问题2，问题分解后如图2所示. 图2 问题分解后二分图Fig.2 Bipartite graph after problem decomposition 1)对于子问题1，n1=2：经判断，子问题1不满足性质4；根据性质3，j1、j6有N(v6)⊆N(v1)且c1≥(d1+d6)，则将设施j6及其关联边从图中删除，即设施j6关闭；根据性质8，设施j1开设，为需求点j1提供2个需求，为需求点j6提供3个需求.将设施j1标记a1=1，需求点i1、i6标记g1=1，g6=1，则子问题1中有dis=n1=2. 2)对于子问题2，n2=5：经判断，子问题2不满足性质4；子问题2中{j2、j3、j4、j5、j7、i2、i3、i4、i5}的容量或需求量存在公约数3，根据性质2，可将当前问题中的容量和需求量同时约去3；根据性质8，设施j5开设，标记a5=1，为需求点i5提供1个需求，将需求点i5标记g5=1；根据性质9，设施j5为需求点i4提供1个需求，需求点i4剩余需求量k4=2. 3)问题经过降阶处理后，更新问题规模如图3所示. 4)执行回溯算法，执行过程如图4所示的二叉树. ① 在搜索过程中，从根结点0出发，计算该结点下界b=5，上界u=6.假设j2开设，进入左子树，该结点下界b=5，小于上界u=6，因此假设j3开设，进入左子树，该结点下界b=5，上界u=6，因此假设j4开设，该结点下界b=5，上界u=6，因此假设j7开设，到达结点1.调用分配子算法，结点1处x2=1，x3=1，x4=1，x7=1构成可行解，此时目标函数值为6； 图3 问题降阶处理后二分图Fig.3 Bipartite graph after problem reduction ② 结点1搜索完成，回溯到上一层.假设设施j7不开设，进入右子树，到达结点2.根据数学性质4，结点2处x2=1，x3=1，x4=1，x7=0，不能构成可行解； 图4 解空间二叉树Fig.4 Solution space binary tree ③ 结点2搜索完成，回溯到上一层.假设设施j4不开设，进入右子树，假设设施j7开设，到达结点3.调用分配子算法，结点3处x2=1，x3=1，x4=0，x7=1构成可行解，此时目标函数值为5； ④ 结点3搜索完成，回溯到上一层.假设设施j7不开设，进入右子树，到达结点4，根据数学性质4，结点4处x2=1，x3=1，x4=0，x7=0，不能构成可行解； ⑤ 结点4搜索完成，回溯到上一层.假设设施j3不开设，进入右子树，到达结点5，根据数学性质4，j2不开设时不能构成可行解，于是剪枝，结点5以下的子树不再搜索； ⑥ 结点5搜索完成，回溯到上一层.假设设施j2不开设，进入右子树，到达结点6，根据数学性质4，j2不开设时不能构成可行解，于是剪枝，结点6以下的子树不再搜索. 通过上述算法得到该问题的最优解集Fbest={j1，j2，j3，j5，j7}，最优目标函数值Q=5，Fbest={j1，j2，j3，j5，j7}其中一个可行的分配如下：j1为需求点i1提供2个需求，为需求点i6提供3个需求；j2为需求点i2提供6个需求，为需求点i4提供3个需求；j3为需求点i3提供6个需求，为需求点i4提供3个需求；j5为需求点i4提供3个需求，为需求点i5提供3个需求；j7为需求点i3提供6个需求，算法结束.从这个例子可以清楚地看到，数学性质能降低问题的规模，上下界子算法能对问题解空间大量剪枝，加快了算法的执行速度. 在精确算法领域，有容量约束的LSCP问题的相关研究较少，本文通过研究有容量集合覆盖选址问题设计了一种精确算法，首先总结出该问题的数学性质并给出证明，这些数学性质能对问题进行降阶从而缩小问题的规模，同时这些数学性质不仅能用于本文的算法中，而且可以应用于其他各种算法；然后利用上界子算法、下界子算法和分配子算法设计出一种能够快速降低问题规模并且能求出最优解的降阶回溯算法.通过理论证明分析以及示例1中对算法执行过程的演示可以看出，本文设计的算法不仅能求出最优解，还能缩小问题的规模，降低算法时间复杂度，加快问题的求解速度.3.3 下界子算法
3.4 上界子算法
4 降阶回溯算法
4.1 降阶子算法
4.2 回溯子算法
4.3 主算法
4.4 算法时间复杂度与对比分析
5 示例分析
6 结 论
