最短路问题求解 | Python+Gurobi-COPT-SCIP

写在前面

沉寂了很久，接下来小编要朝着目标努力了。

目前初步规划了2期内容：用 Python + 三类主流求解器（Gurobi、COPT、SCIP），把网络流问题和铁路运输经典问题两条线一以贯之地讲透。第 1 期讲网络流家族，第 2 期沿着铁路运输组织规划流程逐个建模求解。

把路网抽象一下：

• 节点：家、若干路口/地标、公司

• 边：路口间的路段，带通行时间（分钟）

• 起点：家

• 终点：公司

• 目标：家到公司的总通行时间最小

这就是一个最典型的最短路问题：在带权有向图上，求从起点到终点权值和最小的路径。

为了让模型落地且便于读者手算验证，我们构造一个小型路网。

Gurobi 是学术界最常用的商业求解器，学术免费、性能强劲。

from gurobipy import GRB,Model,quicksum# ============================================================#1. 定义案例数据# ============================================================N = [1, 2, 3, 4, 5, 6]A = [(1, 2),(1, 3),(2, 3),(2, 4),(3, 4),(3, 5),(4, 6),(5, 4),(5, 6)]C = {(1, 2): 6,(1, 3): 4,(2, 3): 2,(2, 4): 2,(3, 4): 1,(3, 5): 2,(4, 6): 7,(5, 4): 1,(5, 6): 3}s = 1t = 6b = {1: 1,    2: 0,    3: 0,    4: 0,    5: 0,    6: -1}# ============================================================# 2. 建模与求解# ============================================================def build_and_solve(N: list, A: list, C: dict,b: dict, s: int, t: int):    # ---- 建模 ----        model = Model("ShortestPath_Gurobi")       # 决策变量 x_{ij} in {0,1}        x = model.addVars(A, vtype=GRB.BINARY, name="x")        # 目标: min sum c_{ij} * x_{ij}        model.setObjective( quicksum(C[i, j] * x[i, j] for (i, j) in A),GRB.MINIMIZE)        # 约束: 流平衡        for i in N:        out_expr = quicksum(x[i, j] for (ii, j) in A if ii == i)        in_expr = quicksum(x[j, i] for (j, ii) in A if ii == i)                model.addConstr(out_expr - in_expr == b[i], name=f"flow_balance_{i}")      # ---- 求解 ----        model.optimize()        # ---- 结果 ----        status = model.Status        if model.Status == GRB.OPTIMAL:                chosen = [(i, j) for (i, j) in A if x[i, j].X > 0.5]                total = model.ObjVal                path = reconstruct_path(chosen, s, t)                return { "status": "OPTIMAL", "total_cost": total, "path": path, "chosen_edges": chosen, }    return {"status": f"status={status}", "total_cost": None, "path": None, "chosen_edges": None}def reconstruct_path(chosen, s, t):        """由选中的边集合按出度链还原有序路径。"""        nxt = {i: j for (i, j) in chosen}        path = [s]        cur = s        while cur != t:              cur = nxt[cur]              path.append(cur)          return path# ============================================================# 3. 主流程# ============================================================def main():        res = build_and_solve(N, A, C,b, s, t)        if res["status"] == "OPTIMAL":                print(f"求解状态: OPTIMAL")                print(f"目标值：{res['total_cost']}")                print(f"选中边: {res['chosen_edges']}")                print(f"最短路径: {res['path']}")       else:                print(f"求解状态: {res['status']} (无可行径路)")if __name__ == "__main__":        main()

COPT（Cardinal Optimizer，杉数科技）是国产高性能商业求解器，对 MIP 求解性能突出。

from coptpy import COPT, Envr# ============================================================# 1. 定义案例数据# ============================================================N = [1, 2, 3, 4, 5, 6]A = [ (1, 2),(1, 3),(2, 3),(2, 4),(3, 4),(3, 5),(4, 6),(5, 4),(5, 6)]C = { (1, 2): 6, (1, 3): 4,(2, 3): 2,(2, 4): 2,(3, 4): 1,(3, 5): 2,(4, 6): 7,(5, 4): 1,(5, 6): 3}s = 1t = 6b = { 1: 1,2: 0, 3: 0, 4: 0,5: 0,6: -1}# ============================================================# 2. 建模与求解# ============================================================def build_and_solve(N: list, A: list, C: dict,b: dict, s: int, t: int):        # ---- 建模 ----        model = Envr().createModel("ShortestPath_COPT")        # 决策变量 x_{ij} in {0,1}        x = model.addVars(A,vtype=COPT.BINARY, nameprefix="x")        # 目标: min sum c_{ij} * x_{ij}        model.setObjective(sum(C[i, j] * x[i, j] for (i, j) in A), COPT.MINIMIZE)        # 约束: 流平衡        for i in N:                out_expr = sum(x[i, j] for (ii, j) in A if ii == i)                in_expr = sum(x[j, i] for (j, ii) in A if ii == i)                model.addConstr(out_expr - in_expr == b[i], name=f"flow_balance_{i}")        # ---- 求解 ----        model.solve()        # ---- 结果 ----        status = model.status        if status == COPT.OPTIMAL:                chosen = [(i, j) for (i, j) in A if x[i, j].X > 0.5]                total = float(model.getAttr(COPT.Attr.BestObj))                path = reconstruct_path(chosen, s, t)                return {"status": "OPTIMAL", "total_cost": total, "path": path, "chosen_edges": chosen}        return {"status": f"status={status}", "total_cost": None, "path": None, "chosen_edges": Nonedef reconstruct_path(chosen, s, t):        """由选中的边集合按出度链还原有序路径。"""        nxt = {i: j for (i, j) in chosen}        path = [s]        cur = s        while cur != t:                cur = nxt[cur]                path.append(cur)            return path# ============================================================#3. 主流程# ============================================================def main():        res = build_and_solve(N, A, C,b, s, t)        if res["status"] == "OPTIMAL":                print(f"求解状态: OPTIMAL")                print(f"目标值：{res['total_cost']}")                print(f"选中边: {res['chosen_edges']}")                print(f"最短路径: {res['path']}")        else:                print(f"求解状态: {res['status']} (无可行径路)")if __name__ == "__main__":        main()

SCIP 是开源求解器（Apache-2.0 协议），无需 license ，适合学习、原型验证和受限部署场景。

from pyscipopt import Model, quicksum# ============================================================# 1. 数据读取# ============================================================N = [ 1, 2, 3, 4, 5, 6]A = [ (1, 2),(1, 3), (2, 3), (2, 4),(3, 4),(3, 5),(4, 6),(5, 4),(5, 6)]C = {    (1, 2): 6,    (1, 3): 4,    (2, 3): 2,    (2, 4): 2,    (3, 4): 1,    (3, 5): 2,    (4, 6): 7,    (5, 4): 1,    (5, 6): 3}s = 1t = 6b = {    1: 1,    2: 0,    3: 0,    4: 0,    5: 0,    6: -1}# ============================================================# 2. 建模与求解# ============================================================def build_and_solve(N: list, A: list, C: dict,b: dict, s: int, t: int):        # ---- 建模 ----        model = Model("ShortestPath_SCIP")        # 决策变量 x_{ij} in {0,1}        x = {}        for (i, j) in A:            x[i, j] = model.addVar(lb=0, ub=1, vtype="BINARY", name=f"x_{i}_{j}")        # 目标: min sum c_{ij} * x_{ij}        model.setObjective(quicksum(C[i, j] * x[i, j] for (i, j) in A), "minimize")        # 约束: 流平衡        for i in N:                out_expr = quicksum(x[i, j] for (ii, j) in A if ii == i)                in_expr = quicksum(x[j, i] for (j, ii) in A if ii == i)                model.addCons(out_expr - in_expr == b[i], name=f"flow_balance_{i}")        # ---- 求解 ----        model.optimize()        # ---- 结果 ----        status = model.getStatus()        if status == "optimal":                chosen = [(i, j) for (i, j) in A if model.getVal(x[i, j]) > 0.5]                total = model.getObjVal()                path = reconstruct_path(chosen, s, t)                return {"status": "OPTIMAL", "total_cost": total, "path": path, "chosen_edges": chosen}        return {"status": f"status={status}", "total_cost": None, "path": None, "chosen_edges": None}def reconstruct_path(chosen, s, t):        """由选中的边集合按出度链还原有序路径。"""        nxt = {i: j for (i, j) in chosen}        path = [s]        cur = s        while cur != t:                cur = nxt[cur]                path.append(cur)        return path# ============================================================# 3. 主流程# ============================================================def main():        res = build_and_solve(N, A, C, b, s, t)        if res["status"] == "OPTIMAL":                print(f"求解状态: OPTIMAL")                print(f"目标值：{res['total_cost']}")                print(f"选中边: {res['chosen_edges']}")                print(f"最短路径: {res['path']}")        else:                print(f"求解状态: {res['status']} (无可行径路)")if__name__ == "__main__":        main()