python调用SCIP求解p-center和p-median选址问题

1. p-中心选址问题

1.1 问题定义

p-中心问题(p-Center Problem) 是一类经典的离散选址问题，已被证明是NP-难问题。p-中心问题是研究其他选址问题的基础，在物流设施规划、通讯系统设计、军队、医院、紧急情况和有服务标准承诺的服务行业等诸多领域具有十分广阔的应用背景。

• p-中心问题是指选定 p 个位置建造服务设施，使得所有客户到最近设施距离中的最大值最小，也就是说，使得最坏情况最优、最大损失最小。

In such problems, we require coverage of all demands, but we seek to locate a given number of facilities in such a way that minimizes coverage distance

• 顶点p-中心问题（vertex center problem）是指待建服务设施的位置都位于网络图的顶点上。

If facility locations are restricted to the nodes of the network, the problem is a vertex center problem.

1.2 数学模型

数学模型见参考文献[1]

• 决策变量
  

• 参数

  • $i$: 需求点的索引
  • $j$: 设施点的索引
  • $d_{ij}$ : 需求点 $i$ 到设施点 $j$ 的距离
  • $P$ : 设施的总数量

2. p-中位选址问题

2.1 问题定义

p-中位问题(p-Median Problem) 是指在备选集里选定$p$个点建造设施，使得需求点到离它最近设施的加权距离总和最小。p-中位问题对模拟现实世界比如公共场所、仓库等的选址非常有用，目标是寻求最小费用。也属于NP-难问题。

The P-median problem uses this measure of effectiveness, and is stated as follows: Find the location of P facilities so as to minimize the total demand-weighted travel distance between demands and facilities.

2.2 数学模型

数学模型见参考文献[1]

• 决策变量
  

• 参数

  • $i$: 需求点的索引
  • $j$: 设施点的索引
  • $d_{ij}$ : 需求点 $i$ 到设施点 $j$ 的距离
  • $P$ : 设施的总数量
  • $h_i$ : 需求点$i$的需求量（理解为权重）

大白话：

• p-中心问题与p-中位问题都是在备选的设施点集合里选p个点。
• 最主要的差别是目标函数，这也决定了他们的应用场景不同。 p-中心问题是希望最小化每个需求点到设施的距离的最大距离，这样不会出现有些距离很大有些距离很小，相对平衡；p-中位问题是希望带权重的距离总和最小。

---

3. python代码实现

3.1 调用SCIP求解p-中心问题与p-中位问题

完整代码如下：

import pandas as pd
import numpy as np
import pyscipopt as opt
import random

def p_center_problem():
    # ==========测试数据==========
    random.seed(1)
    num_nodes, num_facilities, P = 100, 50, 10
    distance_matrix = np.random.randint(1, 20, size=(num_nodes, num_facilities))

    # ==========建立模型==========
    model = opt.Model('p-Center Problem')

    # ==========定义变量==========
    facility = {}
    for j in range(num_facilities):
        facility[j] = model.addVar(vtype='B', name='facility_' + str(j))
        
    assign_serve = {}
    for i in range(num_nodes):
        for j in range(num_facilities):
            assign_serve[i, j] = model.addVar(vtype='B', name='node_' + str(i) + '_facility_' + str(j))
            
    max_distance = model.addVar(vtype='C', name='max_distance')

    # ==========定义约束==========
    # 约束1: 选择的设施的总数量 = P
    model.addCons(opt.quicksum(facility[j] for j in range(num_facilities)) == P)

    # 约束2: 每个需求点都必须被唯一一个设施点满足
    for i in range(num_nodes):
        model.addCons(opt.quicksum(assign_serve[i, j] for j in range(num_facilities)) == 1)

    # 约束3: 如果需求点i被设施点j满足，则代表设施点j被选择。反之亦然
    for i in range(num_nodes):
        for j in range(num_facilities):
            model.addCons(facility[j] >= assign_serve[i, j])
            
    # 约束4: 计算每个需求点最大距离     
    for i in range(num_nodes):
        model.addCons(opt.quicksum(distance_matrix[i][j] * assign_serve[i, j] for j in range(num_facilities)) <= max_distance)
    
    # ==========定义目标==========
    # 目标: 最小化每个需求点的最大距离
    model.setObjective(max_distance)
    model.setMinimize()    
    
    # ==========求解模型========== 
    model.optimize()
    
    select_facility = []
    for i in range(num_facilities):
        if model.getVal(facility[i]) > 0:
            select_facility.append(i)

    assign_nodel_facility = []
    for i in range(num_nodes):
        for j in range(num_facilities):
            if model.getVal(assign_serve[i, j]) > 0:
                assign_nodel_facility.append((i,j))

    print(f'Objective value = {model.getObjVal()}\nselect facility = {select_facility}\nassign_nodel_facility={assign_nodel_facility}')

def p_median_problem():
    # ==========测试数据==========
    random.seed(1)
    num_nodes, num_facilities, P = 100, 50, 10
    distance_matrix = np.random.randint(1, 20, size=(num_nodes, num_facilities))
    weight_demand = np.random.randint(1, 5, size=num_nodes)

    # ==========建立模型==========
    model = opt.Model('p-Median Problem')

    # ==========定义变量==========
    facility = {}
    for j in range(num_facilities):
        facility[j] = model.addVar(vtype='B', name='facility_' + str(j))
        
    assign_serve = {}
    for i in range(num_nodes):
        for j in range(num_facilities):
            assign_serve[i, j] = model.addVar(vtype='B', name='node_' + str(i) + '_facility_' + str(j))
            
    # ==========定义约束==========
    # 约束1: 选择的设施的总数量 = P
    model.addCons(opt.quicksum(facility[j] for j in range(num_facilities)) == P)

    # 约束2: 每个需求点都必须被唯一一个设施点满足
    for i in range(num_nodes):
        model.addCons(opt.quicksum(assign_serve[i, j] for j in range(num_facilities)) == 1)

    # 约束3: 如果需求点i被设施点j满足，则代表设施点j被选择。反之亦然
    for i in range(num_nodes):
        for j in range(num_facilities):
            model.addCons(facility[j] >= assign_serve[i, j])
    
    # ==========定义目标==========
    # 目标: 最小化每个带权重需求点的总和
    model.setObjective(opt.quicksum(weight_demand[i] * distance_matrix[i,j] * facility[j] for i in range(num_nodes) for j in range(num_facilities)))
    model.setMinimize()    
    
    # ==========求解模型========== 
    model.optimize()
    
    select_facility = []
    for i in range(num_facilities):
        if model.getVal(facility[i]) > 0:
            select_facility.append(i)

    assign_nodel_facility = []
    for i in range(num_nodes):
        for j in range(num_facilities):
            if model.getVal(assign_serve[i, j]) > 0:
                assign_nodel_facility.append((i,j))

    print(f'Objective value = {model.getObjVal()}\nselect facility = {select_facility}\nassign_nodel_facility={assign_nodel_facility}')

if __name__ == '__main__':
    # python调用SCIP求解p-中心选址问题
    p_center_problem()
    # python调用SCIP求解p-中位选址问题
    p_median_problem()

3.2 算法结果

• p-中心问题求解结果
  

• p-中位问题求解结果
  

参考文献

[1] Owen S H , Daskin M S . Strategic facility location: A review[J]. Euro.j.oper.res, 1998, 111(3):423-447.