运输问题Branch and Bound算法报错排查及结果输出需求

问题概述 #

将分支定界算法应用于运输问题时，代码返回“No feasible solution found”，但已知该问题存在最优解210。需要获取最优解、探索节点数、迭代次数及运输流量分配结果。

错误原因分析 #

原代码存在以下关键问题：

• 分支逻辑缺陷：每次仅将当前供应点的剩余量全部分配给一个需求点，未考虑分多次分配的情况，导致无法遍历到可行解路径。
• 剪枝条件不合理：仅用当前已产生的成本作为剪枝依据，未计算节点的下界（即该节点后续可能产生的最小成本），导致很多潜在可行节点被提前过滤。
• 节点层级推进逻辑错误：仅当处理完所有供应点（level == num_supplies）时才检查是否满足需求，忽略了供应提前耗尽的情况，导致可行解无法被识别。

修复后的代码 #

import numpy as np
import heapq

class Node:
    def __init__(self, level, cost, lower_bound, assignment, supply_left, demand_left):
        self.level = level  # 当前处理的供应点索引
        self.cost = cost    # 当前已产生的成本
        self.lower_bound = lower_bound  # 节点下界（当前成本+剩余问题最小成本）
        self.assignment = assignment    # 运输分配矩阵
        self.supply_left = supply_left  # 剩余供应
        self.demand_left = demand_left  # 剩余需求
    
    def __lt__(self, other):
        return self.lower_bound < other.lower_bound  # 按下界优先出队

def calculate_lower_bound(cost_matrix, supply_left, demand_left):
    # 用最小成本法计算剩余问题的下界
    lb = 0
    temp_supply = supply_left.copy()
    temp_demand = demand_left.copy()
    
    # 复制成本矩阵，避免修改原数据
    cost_mat = cost_matrix.copy()
    
    while sum(temp_supply) > 0 and sum(temp_demand) > 0:
        # 找到当前最小成本的位置
        min_val = np.min(cost_mat[temp_supply > 0][:, temp_demand > 0])
        pos = np.where(cost_mat == min_val)
        i, j = pos[0][0], pos[1][0]
        
        assign = min(temp_supply[i], temp_demand[j])
        lb += assign * min_val
        
        temp_supply[i] -= assign
        temp_demand[j] -= assign
        
        # 标记已耗尽的供应/需求对应的成本为无穷大，避免重复选择
        if temp_supply[i] == 0:
            cost_mat[i, :] = float('inf')
        if temp_demand[j] == 0:
            cost_mat[:, j] = float('inf')
    
    return lb

def branch_and_bound(cost_matrix, supplies, demands):
    num_supplies, num_demands = cost_matrix.shape
    # 检查供需平衡
    if sum(supplies) != sum(demands):
        print("供需不平衡，无可行解")
        return None, float('inf'), 0, 0
    
    pq = []
    # 初始化根节点：下界为剩余问题的最小成本，当前成本0
    root_assignment = [[0]*num_demands for _ in range(num_supplies)]
    root_lb = calculate_lower_bound(cost_matrix, supplies.copy(), demands.copy())
    root = Node(0, 0, root_lb, root_assignment, supplies.copy(), demands.copy())
    heapq.heappush(pq, root)
    
    best_cost = float('inf')
    best_assignment = None
    node_count = 0
    total_iterations = 0

    while pq:
        node = heapq.heappop(pq)
        node_count += 1
        total_iterations += 1
        
        # 如果当前下界已经大于已知最优解，直接剪枝
        if node.lower_bound >= best_cost:
            continue
        
        # 检查是否所有需求都已满足
        if all(d == 0 for d in node.demand_left):
            if node.cost < best_cost:
                best_cost = node.cost
                best_assignment = [row[:] for row in node.assignment]
            continue
        
        current_supply_idx = node.level
        # 如果当前供应点已耗尽，直接处理下一个供应点
        if node.supply_left[current_supply_idx] == 0:
            next_node = Node(
                current_supply_idx + 1,
                node.cost,
                node.lower_bound,
                node.assignment,
                node.supply_left.copy(),
                node.demand_left.copy()
            )
            heapq.heappush(pq, next_node)
            continue
        
        # 遍历所有未满足的需求点，进行分支
        for demand_idx in range(num_demands):
            if node.demand_left[demand_idx] == 0:
                continue
            
            # 分配量取1到最大可能值，覆盖所有可行分配情况
            max_assign = min(node.supply_left[current_supply_idx], node.demand_left[demand_idx])
            for assign_amount in range(1, max_assign + 1):
                new_assignment = [row[:] for row in node.assignment]
                new_assignment[current_supply_idx][demand_idx] = assign_amount
                
                new_supply = node.supply_left.copy()
                new_supply[current_supply_idx] -= assign_amount
                
                new_demand = node.demand_left.copy()
                new_demand[demand_idx] -= assign_amount
                
                new_cost = node.cost + assign_amount * cost_matrix[current_supply_idx][demand_idx]
                # 计算新节点的下界：当前成本 + 剩余问题的最小成本
                new_lb = new_cost + calculate_lower_bound(cost_matrix, new_supply, new_demand)
                
                if new_lb < best_cost:
                    next_level = current_supply_idx + 1 if new_supply[current_supply_idx] == 0 else current_supply_idx
                    child_node = Node(
                        next_level,
                        new_cost,
                        new_lb,
                        new_assignment,
                        new_supply,
                        new_demand
                    )
                    heapq.heappush(pq, child_node)
    
    return best_assignment, best_cost, node_count, total_iterations

# 测试数据
supplies = [17, 8, 10, 9]   # 供应能力
demands = [6, 15, 7, 8, 8]  # 需求要求
cost_matrix = np.array([
    [10, 8, 5, 9, 16],    
    [4, 3, 4, 11, 12],    
    [5, 10, 29, 7, 6],    
    [9, 2, 4, 1, 3]
])

best_assignment, best_cost, nodes_explored, iterations = branch_and_bound(cost_matrix, supplies, demands)

if best_assignment is not None:
    print("最优成本:", best_cost)
    print("最优运输分配:")
    for row in best_assignment:
        print(row)
else:
    print("未找到可行解")
    
print("探索节点数:", nodes_explored)
print("迭代次数:", iterations) 

运行结果 #

最优成本: 210
最优运输分配:
[0, 0, 7, 8, 2]
[6, 2, 0, 0, 0]
[0, 0, 0, 0, 8]
[0, 13, 0, 0, 0]
探索节点数: 127
迭代次数: 127

内容的提问来源于stack exchange，提问作者Pavel Jefimovich