Z3py约束规划模型返回unsat问题排查求助

问题背景 #

我正在完成大学约束规划作业，任务是用Z3py结合SAT求解器解决多快递员配送问题：m名快递员配送n件物品，目标是最小化所有快递员的最大行驶距离。

相关数据：

• l：各快递员的最大载重数组
• s：每件物品的重量数组
• D：点间距离矩阵，D[i][j]表示物品i到物品j的距离，原点对应矩阵第n行/列，D[n][i]是原点到物品i的距离

核心约束：

1. 所有物品必须被配送
2. 每个快递员需遵守最大载重限制
3. 每个快递员必须从原点出发并返回原点

模型设计 #

模型使用两个布尔变量矩阵：

• p[j][k]：为True表示物品j位于路径的第k个位置
• b[i][j]：为True表示快递员i配送物品j

实现代码：

from z3 import *

def symMax(l):
  m = l[0]
  for v in l[1:]:
    m = If(v > m, v, m)
  return m

def precedes(a1, a2):
    if not a1:
        return True
    if not a2:
        return False
    return Or(And(a1[0] == True, a2[0] == False), And(a1[0] == a2[0], precedes(a1[1:], a2[1:])))

def multiple_couriers_problem_sat(m, n, l, s, D):
    # definition of a 3d array for the path of each courier
    p = [[Bool(f"p_{j}_{k}") for k in range(n//m+2)] for j in range(n)]
    # definition of a matrix that binds a courier with its items
    b = [[Bool(f"b_{i}_{j}") for j in range(n)] for i in range(m)]

    # solver instance
    solver = Optimize()
    solver.set("timeout", 300000)

    # channeling constraint between b and p
    for i in range(m):
      for j in range(n):
        solver.add(Implies(b[i][j], And(exactly_one(p[j]))))
    # channeling constraint between p and b
    for j in range(n):
      solver.add(Implies(And(exactly_one(p[j])), And(exactly_one([b[i][j] for i in range(m)]))))

    for j in range(n):
      solver.add(exactly_one([b[i][j] for i in range(m)]))

    # each courier load capacity is respected
    for i in range(m):
      solver.add(Sum([s[j] * b[i][j] for j in range(n)]) <= l[i])

    # each courier has to distribute at least one item
    for i in range(m):
      solver.add(at_least_one(b[i]))

    # ensures that the items that the courier delivers are in the first positions of the path
    for i in range(m):
     for k in range(n//m+1):
       solver.add(Implies(Sum([And(p[j][k], b[i][j])  for j in range(n)]) == 0, Sum([And(p[j][k+1], b[i][j])  for j in range(n)]) == 0))


    solver.add(PbEq([(x,1) for x in [p[j][0] for j in range(n)]], m)) 

    # symmetry breaking
    for i1 in range(m):
      for i2 in range(i1+1, m):
        if ((i1 != i2) and (l[i1] == l[i2])):
          solver.add(precedes(b[i1], b[i2]))

    # computing the distance
    d = [0 for i in range(m)]
    for i in range(m):
        total_distance = 0
        for k in range(n//m+1):
            for j1 in range(n):
                for j2 in range(n):
                    if j1 != j2:
                        total_distance += If(And(And(b[i][j1], p[j1][k]), And(b[i][j2], p[j2][k+1])), D[j1][j2], 0)
                total_distance += If(And(And(b[i][j1], p[j1][k]), Sum([b[i][t] * p[t][k+1] for t in range(n)]) == 0), D[j1][n] , 0)

        # Add the distance from the origin to the first point
        first_path = Sum([If(And(b[i][j], p[j][0]), D[n][j], 0) for j in range(n)])
        total_distance += first_path
        d[i] = total_distance

    # ascendent order of the courier travelled distances
    for i in range(m-1):
      solver.add(d[i] >= d[i+1])

    max_d = Int('max_d')
    max_d = symMax(d)

    solver.minimize(max_d)

    if solver.check() != unsat:
      return solver.model()
    else:
      return 'unsat'

问题实例（返回unsat） #

该模型在部分实例中可行，但针对以下输入返回unsat：

2
3
18 30
20 17 6
0 21 86 99 
21 0 71 80 
92 71 0 61 
59 80 61 0 

输入说明：

• 第一行m=2（快递员数量）
• 第二行n=3（物品数量）
• 第三行l=[18,30]（快递员最大载重）
• 第四行s=[20,17,6]（物品重量）
• 后续4行是4x4的距离矩阵D（n+1=4个点，含原点）

可行实例参考 #

以下实例可正常求解，解决方案如下：

输入实例 #

2
6
15 10
3 2 6 5 4 4
0 3 4 5 6 6 2 
3 0 1 4 5 7 3 
4 1 0 5 6 6 4 
4 4 5 0 3 3 2 
6 7 8 3 0 2 4 
6 7 8 3 2 0 4 
2 3 4 3 4 4 0 

解决方案 #

[['p_3_0', 'p_2_1', 'p_0_2'], ['p_1_0', 'p_4_1', 'p_5_2']]
14

注：每个子数组对应一名快递员的配送路径，隐含从原点出发、最后返回原点的逻辑。

求助点 #

无法理解上述问题实例返回unsat的原因，希望得到分析和解决建议。

内容的提问来源于stack exchange，提问作者Tech Rufy