SSCFLP（Single Source Capacity Facility Location Problem)-CSDN博客

本文链接：https://blog.csdn.net/qq_36347365/article/details/84964053

文章目录

以下问题中客户与工厂的矩阵读反了！！！
问题：SSCFLP
分析
题解
四种算法结果比较
结论
参考资料

以下问题中客户与工厂的矩阵读反了！！！

问题：SSCFLP

SSCFLP，也就是单源设施选址问题。问题简述如下：

有m个工厂，n个用户（以下下标中i表示的是工厂，j表示的是用户）
每一个工厂的费用为 $f_{i}$ ，同时每一个工厂也是有容量的，为 $b_{i}$
每一个用户的需求为 $a_{j}$
同时每一个工厂服务对应的不同的用户的费用是不一样的，需要用一个二维数组表示为 $c_{ji}$
条件：
每一个用户只能够选择一个工厂服务，并且每一个工厂所服务的用户的总需求不能够超过自己的容量
求解：
最小的总的费用：有开的工厂的费用加上工厂服务每一个用户的费用合成。

分析

这一道题主要的难点在于添加了容量限制与工产费用限制，其中每个工厂开设需要一定的费用，解决工厂开设的费用与用户需求、客户选择服务工厂以及工厂容量问题四者之间的平衡问题成为了解题的关键之处。
比较精确的解法是分支定界法，也就是整数线性规划问题，因为这一道题是可以表示称最小化问题加上一系列的约束条件。表示为如下：
$x_{ji}$ 表示的是二维数组 $c_{ji}$ 前面的系数，只有当用户j选择工厂i时该系数才为1，否则为0
最小化问题如下： $Minimize(\sum_{i = 1}^{m}f_{i}y_{i} + \sum_{i = 1}^{m}\sum_{j = 1}^{n}x_{ji}c_{ji})$
约束条件：
$\begin{cases} \sum_{j = 1}^{n}a_{j}x_{ji} \le b_{i}, i=1, ...,m;j=1, ...,n \\ \sum_{i = 1}^{m}x_{ji} = 1, j = 1, ..., n \\ x_{ji} \le y_{i}, i=1, ...,m;j=1, ...,n \\ x_{ji} \in \{ 0, 1\}, i=1, ...,m;j=1, ...,n \\ y_i \in \{0, 1\}, i=1, ...,m; \\ \end{cases}$
而精确求解大规模问题会造成巨大的时间消耗。因为线性规划大多使用对偶式等，在超高维问题中容易造成时间爆炸性问题，因此实现一种求的较优解（10%）以及缩短求解时间是有必要的。
最简单但是不能达到最优解的方法是贪心算法，贪心的规则是：
1. 维持两类工厂（已经打开的与仍未打开的）
2. 在已经打开的工厂中寻找符合条件的（不会超出工厂总的服务量）的最小的用户需求 $f e e 1$ ，以及在未打开的工厂中寻找最小的用户需求 $f e e 2$ ，将后者加上工厂的费用 $cost_{factory}$ （可以除以 $size_{customer}$ ）然后两者进行对比，选取两者之中的较小者
  $min_{chose} = min\{fee1, fee2 + \frac {cost_{factory}}{size_{customer}}\} (size_{customer}可替换为其他值)$
3. 如果选择未开的工厂则将未开的工厂加入已开的工厂中
4. 重复以上的过程直到所有的客户都已经选择服务工厂重复以上的过程直到所有的客户都已经选择服务工厂重复以上的过程直到所有的客户都已经选择服务工厂重复以上的过程直到所有的客户都已经选择服务工厂
5. 贪心算法的模型是在原来的四元关系的基础上的简化：
上面的贪心算法的规则可以自定义：可以将前提假设为所有的工厂都打开等，但是都有各自的缺点，因为没有处理好平衡的问题，在工厂开设费用较大或者工厂的数目较多的时候会出现极端的现象。
在贪心算法的基础上我们可以采用模拟退火算法，将贪心算法产生的解作为模拟退火的初始状态，在新的状态产生的过程中采用局部贪心的策略，并行计算四种新状态产生的方法：
1. 两点交换
2. 两点之间逆序
3. 单点插入
4. 变异（重要）

在上面的四种方法中选择结果最优的结果，其中变异是比较重要的，因为是对每个工厂的选择数目进行改变。而增加变异率是有必要的，所以在进行变异的过程中如果变异的效果比较好的话就继续进行变异的操作，从而使得解空间更加丰富。

在模拟退火的过程中如果不接受差解则转化为局部搜索算法（LS），也就是爬山法。
因此上面的四种方法都有各自的优缺点。可以列表个如下：因此上面的四种方法都有各自的优缺点。可以列表个如下：

算法	时间复杂度	求解效果
分支定界（整数规划问题）	复杂（未证明是多项式时间）	精确
贪心算法	简单（线性时间）	粗糙
局部搜索法	较简单（线性时间）	较精确，但容易陷入局部最优解
基于贪心的模拟退火算法	较简单（线性时间）	较精确

题解

工具：python3.7 （numpy、re、pulp、random、time）

方法1：分支定界法（整数规划）

代码实现（请慎重运行，个别测例以小时为单位）

from numpy import *
import re
from pulp import *
import time


class Problem:
    # 读取文件初始化数据
    def __init__(self, filename):
        self.data = open(filename, "r")
        [self.factory_size, self.customer_size] = self.list_str_to_int(self.data.readline().split())
        # 字典， 用于制定方程式
        self.factorys_cost = dict()
        self.factorys_capacity = dict()
        self.customers_need = dict()
        self.c2f_fee = dict()

        # 列表，用于存储索引
        self.factorys_cost_list = []
        self.factorys_capacity_list = []
        self.customers_need_list = []
        self.c2f_fee_list = []
        for i in range(self.factory_size):
            [a, b] = self.list_str_to_int(self.data.readline().split())
            self.factorys_capacity[i] = a
            self.factorys_capacity_list.append(i)
            self.factorys_cost[i] = b
            self.factorys_cost_list.append(i)
        count = 0
        while count < self.customer_size:
            temp = self.list_str_to_int(re.split(r'[\s\.]+', self.data.readline()))
            l = len(temp)
            for i in range(l):
                self.customers_need[count + i] = temp[i]
                self.customers_need_list.append(count + i)
            count = len(self.customers_need.keys())
        for i in range(self.customer_size):
            count = 0
            while count < self.factory_size:
                temp = self.list_str_to_int(re.split(r'[\s\.]+', self.data.readline()))
                l = len(temp)
                for j in range(l):
                    self.c2f_fee[i, j+count] = temp[j]
                    self.c2f_fee_list.append((i, j+count))
                count += l
        # print(self.customers_need)
        # print("_" *100)
        # print(self.factorys_capacity)
        # print("_"*100)
        # print(self.c2f_fee)
        # print("_" * 100)
        # print(self.factorys_cost)
        # self.data.close()

    # 用于字符串列表转化为int型
    def list_str_to_int(self, str_list):
        int_list = []
        for i in str_list:
            if i == "":
                continue
            int_list.append(int(i))
        return int_list


def remove_none_str(list1):
    while "" in list1:
        list1.remove("")


if __name__ == '__main__':
    open_file = open("Linear.md", "a")
    for p in range(71):
        num = str(p+1)
        file = "./Instances/p" + num
        problem = Problem(filename=file)
        # 声明参数变量
        factory_open = LpVariable.dicts("factory_oepn", problem.factorys_cost_list, 0, 1, LpInteger)
        # print(factory_open)
        c2f_which = LpVariable.dicts("c2f_which", problem.c2f_fee_list, 0, 1, LpInteger)
        # print(c2f_which)
        # 创建极小化问题
        prob = LpProblem("SSCFLP", LpMinimize)
        a =  lpSum([factory_open[i] * problem.factorys_cost[i] for i in problem.factorys_capacity_list]) + lpSum([c2f_which[i] * problem.c2f_fee[i] for i in problem.c2f_fee_list])
        prob += a, "Total cost"
        # 客户的约束，要求每一个客户只能够有一个工厂
        for i in problem.customers_need_list:
            b = (lpSum(c2f_which[i, j] for j in problem.factorys_cost_list) == 1)
            prob += b
        # 对工厂的约束，服务的客户的数目不能够超过工厂的服务度量
        for i in problem.factorys_cost_list:
            b = lpSum(problem.customers_need[j] * c2f_which[j, i] for j in problem.customers_need_list) <= problem.factorys_capacity[i]
            prob += b
        # 对矩阵中所有的元素的约束
        for i in problem.factorys_cost_list:
            for j in problem.customers_need_list:
                b = c2f_which[j, i] <= factory_open[i]
                prob += b
        # 求解线性方程组
        start_time = time.time()
        prob.solve(PULP_CBC_CMD(fracGap = 0.001))
        calculate_time = time.time() - start_time
        # 打印状态
        # print("Status:", LpStatus[prob.status])
        # 打印出所有的数据结果，处理数据结果
        factory_result = zeros([2, problem.factory_size], dtype=int)
        customer_result = zeros([3, problem.customer_size], dtype=int)
        for i in range(problem.factory_size):
            factory_result[0, i] = i
        for i in range(problem.customer_size):
            customer_result[0, i] = i
        factory_result_list = []
        # 读取结果数据
        for v in prob.variables():
            if "factory_oepn_" in v.name:
                b = re.split(r"[\s\_]+", v.name)
                remove_none_str(b)
                factory_result[1, int(b[2])] = v.varValue
                if v.varValue == 1:
                    factory_result_list.append(int(b[2]))
            elif "c2f_which_" in v.name and v.varValue == 1:
                b = re.split(r"[\s\_\(\)\,]+", v.name)
                remove_none_str(b)
                customer_result[1, int(b[2])] = int(b[3])
                customer_result[2, int(b[2])] = problem.c2f_fee[int(b[2]), int(b[3])]
        # print(customer_result)
        # print(sort(factory_result_list))
        # print(value(prob.objective))
        # print(calculate_time)
        # 写入文件
        open_file.write("- " + str("test: " + str(num)) + "\n")
        open_file.write("- total fee: " + str(value(prob.objective)) + "\n")
        open_file.write("- open_factory_list: " + str(sort(factory_result_list)) + "\n")
        open_file.write("- time: " + str(calculate_time) + "\n")
        open_file.write("\n")
        open_file.write("| customer")
        for j in range(len(customer_result[0, :])):
            open_file.write(" | " + str(customer_result[0, j]))
        open_file.write(" |\n")
        open_file.write("| ------ ")
        for j in range(len(customer_result[0, :])):
            open_file.write("| ------ ")
        open_file.write("|\n")
        open_file.write("| factory")
        for j in range(len(customer_result[1, :])):
            open_file.write("| " + str(customer_result[1, j]))
        open_file.write(" |\n")
        open_file.write("| fee")
        for j in range(len(customer_result[2, :])):
            open_file.write("| " + str(customer_result[2, j]))
        open_file.write(" |\n")
        open_file.write("\n\n\n")
    open_file.close()
    # print("Total Cost of Transportation = ", value(prob.objective))
    # 打印结果
    # prob.writeLP("1.txt")
    # 检验和
    # result = 0
    # for i in range(problem.customer_size):
    #     result += (problem.c2f_fee[i, customer_result[1, i]])
    # for j in range(problem.factory_size):
    #     if factory_result[1][j] == 1:
    #         result += (problem.factorys_cost[j])
    # print(result)

结果数据

见：SSCFLP Benchmark Data 线性规划算法结果

方法2：贪心算法

代码实现

from numpy import *
import re, time


class Factory:
    def __init__(self, capacity, open_cost):
        self.capacity = capacity
        self.open_cost = open_cost

    def whether(self, data):
        if self.capacity >= data:
            return True
        else:
            return False

    def minus(self, data):
        self.capacity -= data


class Customer:
    def __init__(self, need):
        self.need = need


class Problem:
    # 读取文件初始化数据
    def __init__(self, filename):
        self.data = open(filename, "r")
        [self.factory_size, self.customer_size] = self.list_str_to_int(self.data.readline().split())
        self.factorys = dict()
        self.customers = dict()
        self.c2f_fee = zeros([self.customer_size, self.factory_size], dtype=int)
        for i in range(self.factory_size):
            [a, b] = self.list_str_to_int(self.data.readline().split())
            self.factorys[i] = Factory(a, b)
        count = 0
        while count < self.customer_size:
            temp = self.list_str_to_int(re.split(r'[\s\.]+', self.data.readline()))
            l = len(temp)
            for i in range(l):
                self.customers[count + i] = Customer(temp[i])
            count = len(self.customers.keys())
        self.max_num = 0
        for i in range(self.customer_size):
            count = 0
            while count < self.factory_size:
                temp = self.list_str_to_int(re.split(r'[\s\.]+', self.data.readline()))
                self.max_num = max(self.max_num, array(temp).max())
                l = len(temp)
                for j in range(l):
                    self.c2f_fee[i][j+count] = temp[j]
                count += l
        # print("max", self.max_num)
        self.data.close()

    # 用于每一个获取被满足的用户的最小值之后的置为最大值，避免重取
    # 用于每一个达到最大的服务数量的工厂的重置，以避免重取
    def reset(self, customer_id=-1, factory_id=-1):
        if customer_id != -1:
            for i in range(self.factory_size):
                self.c2f_fee[customer_id][i] = self.max_num + 1
        else:
            for i in range(self.customer_size):
                self.c2f_fee[i][factory_id] = self.max_num + 1

    # 获取指定的工厂的最小值
    # 获取全图中的最小值
    # 返回最小值与坐标
    # 复杂度：o(m*n)
    def find_min(self, factory_id=[]):
        r = 0
        c = 0
        if factory_id != []:
            min_fee = self.c2f_fee[0][factory_id[0]]
            c = factory_id[0]
            for j in factory_id:
                for i in range(self.customer_size):
                    if min_fee > self.c2f_fee[i][j]:
                        r = i
                        c = j
                        min_fee = self.c2f_fee[i][j]
            return min_fee, [r, c]
        else:
            min_fee = self.c2f_fee[0][0]
            for i in range(self.customer_size):
                for j in range(self.factory_size):
                    if min_fee > self.c2f_fee[i][j]:
                        r = i
                        c = j
                        min_fee = self.c2f_fee[i][j]
            return min_fee, [r, c]

    # 用于字符串列表转化为int型
    def list_str_to_int(self, str_list):
        int_list = []
        for i in str_list:
            if i == "":
                continue
            int_list.append(int(i))
        return int_list


    # 贪心算法
    # 第一次是寻找最小的数值，因为所有的工厂都还没有开放
    # 第二次是寻找剩余的数组中的最小的数值 以及已开放工厂中的最小的数值进行比较，选择两者之中最小的数值
if __name__ == '__main__':
    # 打开写入的文件
    open_file = open("Greedy.md", "a")
    for i in range(71):
        num = str(i + 1)
        # 读取测例
        file = "./Instances/p" + num
        problem = Problem(filename=file)
        # 记录打开的工厂
        open_list = []
        # 记录用户列表
        customer_list = []
        # 记录未打开的工厂
        not_open_list = [i for i in range(problem.factory_size)]
        pairs = zeros([3, problem.customer_size], dtype=int)
        for j in range(problem.customer_size):
            pairs[0, j] = j
        # 寻找最小值
        min_fee, [r, c] = problem.find_min()
        fee = 0
        start_time = time.time()
        # 满足条件插入
        if problem.factorys[c].whether(problem.customers[r].need):
            open_list.append(c)
            not_open_list.remove(c)
            customer_list.append(r)
            pairs[1, r] = c
            pairs[2, r] = min_fee
            problem.factorys[c].minus(problem.customers[r].need)
            fee += (problem.factorys[c].open_cost + min_fee)
            problem.reset(customer_id=r)
        # 寻找未打开的工厂中的最小值
        min_fee1, [r1, c1] = problem.find_min(factory_id=not_open_list)
        while len(customer_list) < problem.customer_size:
            if len(not_open_list) > 0:
                min_fee1, [r1, c1] = problem.find_min(factory_id=not_open_list)
            min_fee2, [r2, c2] = problem.find_min(factory_id=open_list)
            # 分两种情况插入
            if len(not_open_list) > 0 and problem.factorys[c1].whether(problem.customers[r1].need) and (
                    min_fee2 != problem.max_num + 1 and problem.factorys[
                c1].open_cost + min_fee1 < min_fee2 or min_fee2 == problem.max_num + 1):
                open_list.append(c1)
                problem.factorys[c1].minus(problem.customers[r1].need)
                fee += problem.factorys[c1].open_cost + min_fee1
                customer_list.append(r1)
                pairs[1, r1] = c1
                pairs[2, r1] = min_fee1
                not_open_list.remove(c1)
                problem.reset(customer_id=r1)
            elif problem.factorys[c2].whether(problem.customers[r2].need):
                customer_list.append(r2)
                pairs[1, r2] = c2
                pairs[2, r2] = min_fee2
                problem.factorys[c2].minus(problem.customers[r2].need)
                fee += min_fee2
                problem.reset(customer_id=r2)
            # 已经打开的某个工厂满了，进行重置，置为最大值
            elif not problem.factorys[c2].whether(problem.customers[r2].need):
                problem.reset(factory_id=c2)
        # 总的计算时间
        calculate_time = time.time() - start_time
        # 文件写入
        open_file.write("- " + str("test: " + str(i)) + "\n")
        open_file.write("- fee:" + str(fee) + "\n")
        # 写入打开的工厂
        open_file.write("- " + str(open_list) + "\n")
        open_file.write("- time" + str(calculate_time) + "\n")
        open_file.write("\n")
        # 写入每一个用户，从0-customer_size
        open_file.write("| customer" )
        for j in range(len(pairs[0, :])):
            open_file.write(" | " + str(pairs[0, j]))
        open_file.write(" |\n")
        open_file.write("| ------ ")
        for j in range(len(pairs[0, :])):
            open_file.write("| ------ ")
        open_file.write("|\n")
        # 写入每一个用户的选择
        open_file.write("| factory")
        for j in range(len(pairs[1, :])):
            open_file.write("| " + str(pairs[1, j]))
        open_file.write(" |\n")
        # 写入每一个选择的费用
        open_file.write("| fee")
        for j in range(len(pairs[2, :])):
            open_file.write("| " + str(pairs[2, j]))
        open_file.write(" |\n")
        open_file.write("\n\n\n")
    open_file.close()

结果数据

见 SSCFLP Benchmark Data 贪心算法结果

方法3：基于贪心的模拟退火算法（以及局部搜索法）

代码实现

from numpy import *
import re
import time
import random
import copy

global x_index
global y_index
global problem
global distance


class Factory:
    def __init__(self, capacity, open_cost):
        self.capacity = capacity
        self.open_cost = open_cost

    def whether(self, data):
        if self.capacity >= data:
            return True
        else:
            return False

    def minus(self, data):
        self.capacity -= data


class Customer:
    def __init__(self, need):
        self.need = need


class Problem:
    # 读取文件初始化数据
    def __init__(self, filename):
        self.data = open(filename, "r")
        [self.factory_size, self.customer_size] = self.list_str_to_int(self.data.readline().split())
        self.factory_capacity = zeros([self.factory_size], dtype=int)
        self.factory_cost = zeros([self.factory_size], dtype=int)
        self.customer = zeros([self.customer_size], dtype=int)
        self.c2f_fee = zeros([self.customer_size, self.factory_size], dtype=int)
        self.factorys = dict()
        self.customers = dict()
        for i in range(self.factory_size):
            [a, b] = self.list_str_to_int(self.data.readline().split())
            self.factory_capacity[i] = a
            self.factory_cost[i] = b
            self.factorys[i] = Factory(a, b)
        count = 0
        while count < self.customer_size:
            temp = self.list_str_to_int(re.split(r'[\s\.]+', self.data.readline()))
            l = len(temp)
            for i in range(l):
                self.customer[count] = temp[i]
                self.customers[count] = Customer(temp[i])
                count += 1
        self.max_num = 0
        for i in range(self.customer_size):
            count = 0
            while count < self.factory_size:
                temp = self.list_str_to_int(re.split(r'[\s\.]+', self.data.readline()))
                self.max_num = max(self.max_num, array(temp).max())
                l = len(temp)
                for j in range(l):
                    self.c2f_fee[i][j+count] = temp[j]
                count += l
        # print(self.factory_capacity)
        # print(self.factory_cost)
        # print(self.c2f_fee)
        # print(self.customer)
        self.c2f_fee2 = copy.deepcopy(self.c2f_fee)
        self.data.close()

    # 用于每一个获取被满足的用户的最小值之后的置为最大值，避免重取
    # 用于每一个达到最大的服务数量的工厂的重置，以避免重取
    def reset(self, customer_id=-1, factory_id=-1):
        if customer_id != -1:
            for i in range(self.factory_size):
                self.c2f_fee[customer_id][i] = self.max_num + 1
        else:
            for i in range(self.customer_size):
                self.c2f_fee[i][factory_id] = self.max_num + 1

    # 获取指定的工厂的最小值
    # 获取全图中的最小值
    # 返回最小值与坐标
    # 复杂度：o(m*n)
    def find_min(self, factory_id=[]):
        r = 0
        c = 0
        if factory_id != []:
            min_fee = self.c2f_fee[0][factory_id[0]]
            c = factory_id[0]
            for j in factory_id:
                for i in range(self.customer_size):
                    if min_fee > self.c2f_fee[i][j]:
                        r = i
                        c = j
                        min_fee = self.c2f_fee[i][j]
            return min_fee, [r, c]
        else:
            min_fee = self.c2f_fee[0][0]
            for i in range(self.customer_size):
                for j in range(self.factory_size):
                    if min_fee > self.c2f_fee[i][j]:
                        r = i
                        c = j
                        min_fee = self.c2f_fee[i][j]
            return min_fee, [r, c]

    # 用于字符串列表转化为int型
    def list_str_to_int(self, str_list):
        int_list = []
        for i in str_list:
            if i == "":
                continue
            int_list.append(int(i))
        return int_list

    def valid(self, arr):
        factorys = zeros([self.factory_size], dtype=int)
        sum_fee = 0
        judge = True
        for i in range(len(arr)):
            factorys[arr[i]] += self.customer[i]
            sum_fee += self.c2f_fee2[i][arr[i]]
        for i in range(self.factory_size):
            if factorys[i] > 0:
                sum_fee += self.factory_cost[i]
            if factorys[i] > self.factory_capacity[i]:
                judge = False
        return judge, sum_fee


# 第一种交换方式：两个点进行交换
def method1(arr):
    array1 = copy.deepcopy(arr)
    array1[x_index], array1[y_index] = arr[y_index], arr[x_index]
    tag, fee_sum = problem.valid(array1)
    return tag, fee_sum, array1


# 第二种交换方式，两个点之间的点进行逆序
def method2(arr):
    array2 = copy.deepcopy(arr)
    num_min = min(x_index, y_index)
    num_max = max(x_index, y_index)
    for i in range(num_max - num_min + 1):
        array2[num_min + i] = arr[num_max - i]
    tag, fee_sum = problem.valid(array2)
    return tag, fee_sum, array2


# 第三种交换方式，一个点的前移
def method3(arr):
    array3 = copy.deepcopy(arr)
    num_min = min(x_index, y_index)
    num_max = max(x_index, y_index)
    temp = arr[num_max]
    for i in range(num_max - num_min):
        array3[num_max - i] = arr[num_max - i - 1]
    array3[num_min] = temp
    tag, fee_sum = problem.valid(array3)
    return tag, fee_sum, array3


# 第四种交换方式：变异
def method4(arr):
    tag, fee_sum = problem.valid(arr)
    array4 = copy.deepcopy(arr)
    num1 = random.randint(0, problem.factory_size - 1)
    num2 = random.randint(0, problem.factory_size - 1)
    array4[x_index] = num1
    array4[y_index] = num2
    tag2, fee_sum2 = problem.valid(array4)
    dimension = len(arr)
    # 如果不合法
    if not tag2:
        return tag2, fee_sum2, array4
    # 如果合法
    else:
        # 如果比原来的值小则继续进行变异
        if fee_sum2 <= fee_sum:
            array5 = copy.deepcopy(array4)
            while tag2 and fee_sum2 < fee_sum:
                tag, fee_sum = tag2, fee_sum2
                array5 = copy.deepcopy(array4)
                num1 = random.randint(0, problem.factory_size - 1)
                num2 = random.randint(0, problem.factory_size - 1)
                index1 = random.randint(0, dimension - 1)
                index2 = random.randint(0, dimension - 1)
                while index1 == index2:
                    index2 = random.randint(0, dimension - 1)
                array4[index1] = num1
                array4[index2] = num2
                tag2, fee_sum2 = problem.valid(array4)
            return tag, fee_sum, array5
        # 否则返回大的值
        else:
            return tag2, fee_sum2, array4


# 总的交换处理函数，四种交换方式取最优（局部贪心）
def swap(arr, dimension):
    global x_index
    global y_index
    # 产生两个随机值（位置）
    x_index = random.randint(0, dimension - 1)
    y_index = random.randint(0, dimension - 1)
    # 如果相同则继续产生
    while x_index == y_index:
        y_index = random.randint(0, dimension - 1)
    result1 = method1(arr)
    result2 = method2(arr)
    result3 = method3(arr)
    result4 = method4(arr)
    tag, fee_sum  = problem.valid(arr)
    array_result = copy.deepcopy(arr)
    if not (result1[0] or result2[0] or result3[0] or result4[0]):
        return fee_sum, array_result
    elif result1[0]:
        fee_sum = result1[1]
        array_result = copy.deepcopy(result1[2])
    elif result2[0]:
        fee_sum = result2[1]
        array_result = copy.deepcopy(result2[2])
    elif result3[0]:
        fee_sum = result3[1]
        array_result = copy.deepcopy(result3[2])
    elif result4[0]:
        fee_sum = result4[1]
        array_result = copy.deepcopy(result4[2])


    if result1[0] and result1[1] < fee_sum:
        fee_sum = result1[1]
        array_result = copy.deepcopy(result1[2])
    if result2[0] and result2[1] < fee_sum:
        fee_sum = result2[1]
        array_result = copy.deepcopy(result2[2])
    if result3[0] and result3[1] < fee_sum:
        fee_sum = result3[1]
        array_result = copy.deepcopy(result3[2])
    if result4[0] and result4[1] < fee_sum:
        fee_sum = result4[1]
        array_result = copy.deepcopy(result4[2])
    return fee_sum, array_result


    # 模拟退火算法核心(注释掉if条件中的or部分就是局部搜索算法
def sscflp_sa(initial_temp, rate, arr, dimension):
    global distance
    while initial_temp > 0.1:
        for i in range(500):
            new_result = swap(arr, dimension)
            diff = abs(distance - new_result[0])
            # 注视掉or后面的就是局部搜索算法
            if new_result[0] < distance or random.random() < math.exp(-1 * (diff / initial_temp)):
                distance = copy.deepcopy(new_result[0])
                arr = copy.deepcopy(new_result[1])
        initial_temp = initial_temp * rate
    return arr, distance


if __name__ == '__main__':
    global problem
    global distance
    # for i in range(71):
    #     num = str(i + 1)
    #     file = "./Instances/p" + num
    open_file = open("SA1.md", "a")
    for k in range(71):
        num = str(k + 1)
        file = "./Instances/p" + num
        problem = Problem(filename=file)
        # 以下部分是贪心算法
        open_list = []
        customer_list = []
        not_open_list = [i for i in range(problem.factory_size)]
        pairs = zeros([3, problem.customer_size], dtype=int)
        for j in range(problem.customer_size):
            pairs[0, j] = j
        min_fee, [r, c] = problem.find_min()
        fee = 0
        if problem.factorys[c].whether(problem.customers[r].need):
            open_list.append(c)
            not_open_list.remove(c)
            customer_list.append(r)
            pairs[1, r] = c
            pairs[2, r] = min_fee
            problem.factorys[c].minus(problem.customers[r].need)
            fee += (problem.factorys[c].open_cost + min_fee)
            problem.reset(customer_id=r)
        min_fee1, [r1, c1] = problem.find_min(factory_id=not_open_list)
        while len(customer_list) < problem.customer_size:
            if len(not_open_list) > 0:
                min_fee1, [r1, c1] = problem.find_min(factory_id=not_open_list)
            min_fee2, [r2, c2] = problem.find_min(factory_id=open_list)
            if len(not_open_list) > 0 and problem.factorys[c1].whether(problem.customers[r1].need) and (
                    min_fee2 != problem.max_num + 1 and problem.factorys[
                c1].open_cost + min_fee1 < min_fee2 or min_fee2 == problem.max_num + 1):
                open_list.append(c1)
                problem.factorys[c1].minus(problem.customers[r1].need)
                fee += problem.factorys[c1].open_cost + min_fee1
                customer_list.append(r1)
                pairs[1, r1] = c1
                pairs[2, r1] = min_fee1
                not_open_list.remove(c1)
                problem.reset(customer_id=r1)
            elif problem.factorys[c2].whether(problem.customers[r2].need):
                customer_list.append(r2)
                pairs[1, r2] = c2
                pairs[2, r2] = min_fee2
                problem.factorys[c2].minus(problem.customers[r2].need)
                fee += min_fee2
                problem.reset(customer_id=r2)
            elif not problem.factorys[c2].whether(problem.customers[r2].need):
                problem.reset(factory_id=c2)



        # 此处开始模拟退火算法
        initial_result = array(pairs[1, :])
        tag, distance = problem.valid(initial_result)
        start_time = time.time()
        arr, distance = sscflp_sa(100, 0.95, initial_result, problem.customer_size)
        # 总的计算时间
        calculate_time = time.time() - start_time
        # print(distance)
        # print(arr)
        fee = distance
        open_factory_list = []
        for i in range(len(arr)):
            pairs[1, i] = arr[i]
            pairs[2, i] = problem.c2f_fee2[i, arr[i]]
            if arr[i] not in open_factory_list:
                open_factory_list.append(arr[i])
        # print(calculate_time)
        # print(pairs)
        # print(open_factory_list)
        # 文件写入
        open_file.write("- " + str("test: " + str(k)) + "\n")
        open_file.write("- fee:" + str(fee) + "\n")
        # 写入打开的工厂
        open_file.write("- open-factory-list: " + str(open_factory_list) + "\n")
        open_file.write("- time" + str(calculate_time) + "\n")
        open_file.write("\n")
        # 写入每一个用户，从0-customer_size
        open_file.write("| customer")
        for j in range(len(pairs[0, :])):
            open_file.write(" | " + str(pairs[0, j]))
        open_file.write(" |\n")
        open_file.write("| ------ ")
        for j in range(len(pairs[0, :])):
            open_file.write("| ------ ")
        open_file.write("|\n")
        # 写入每一个用户的选择
        open_file.write("| factory")
        for j in range(len(pairs[1, :])):
            open_file.write("| " + str(pairs[1, j]))
        open_file.write(" |\n")
        # 写入每一个选择的费用
        open_file.write("| fee")
        for j in range(len(pairs[2, :])):
            open_file.write("| " + str(pairs[2, j]))
        open_file.write(" |\n")
        open_file.write("\n\n\n")
    open_file.close()

    # print(problem.valid(arr))

结果数据

模拟退火：SSCFLP Benchmark Data 基于贪心策略的模拟退火算法结果
LS局部搜索：SSCFLP Benchmark Data 基于贪心策略的局部搜索算法结果

四种算法结果比较

下面有空缺的是因为时间超长运行不出来的（以小时为单位的）

	结果						时间消耗(s）
	分支定界法	贪心算法	LS局部搜索	SA模拟退火	目前的最优解	目前的最差解	分支定界法	贪心算法	LS局部搜索	SA模拟退火
0	9907	12458	10017	10075	9907	12458	0.27	0.01	31.4	29.91
1	8908	10733	9123	8917	8908	10733	0.24	0.01	30	30.04
2	10374	13877	10547	10574	10374	13877	0.22	0.01	29.04	29.78
3	11774	15701	12263	12169	11774	15701	0.87	0.01	29.4	29.04
4	9921	11925	10154	10238	9921	11925	0.53	0.01	28.92	28.86
5	8994	10905	9170	8994	8994	10905	0.54	0.01	28.82	28.88
6	10625	12505	10768	10825	10625	12505	1.06	0.01	28.72	28.8
7	12225	14105	12346	12280	12225	14105	1.47	0.01	28.67	28.74
8	9530	11609	9677	9677	9530	11609	0.23	0.01	28.56	28.64
9	8778	10890	8958	8917	8778	10890	0.16	0.01	28.59	28.89
10	10146	12467	10354	10392	10146	12467	0.41	0.01	28.51	28.6
11	11312	14687	12272	12357	11312	14687	1.2	0.01	28.59	28.58
12	8919	11114	10123	9671	8919	11114	1.74	0.02	31.1	31.15
13	7908	9947	8515	8372	7908	9947	1.27	0.02	31.17	31.2
14	9508	11979	10476	10529	9508	11979	2.24	0.02	31.07	31.08
15	11108	14111	12623	12603	11108	14111	4.79	0.02	30.93	30.84
16	8815	11688	9941	9673	8815	11688	0.97	0.02	31.09	31.15
17	7800	9831	8241	8202	7800	9831	0.46	0.02	31.17	31.33
18	9441	12069	10397	10322	9441	12069	2.48	0.02	31.05	31.12
19	11041	14327	12436	12396	11041	14327	9.35	0.02	30.9	30.86
20	8706	12453	9867	9399	8706	12453	0.65	0.02	31.03	31.17
21	7780	10815	8400	8232	7780	10815	0.45	0.02	31.1	31.34
22	9308	13013	10591	9942	9308	13013	2.01	0.02	30.99	31.09
23	10551	15130	12257	12257	10551	15130	2.17	0.02	30.85	30.78
24	11008	22673	12096	12925	11008	22673	3.71	0.26	77.38	77.72
25	10198	17957	10725	11762	10198	17957	2.79	0.26	77.06	77.71
26	11795	23274	13027	13760	11795	23274	8.78	0.25	77.2	77.49
27	13195	25613	16069	16004	13195	25613	33.12	0.26	77.03	77.17
28	11844	17548	13184	13774	11844	17548	962.82	1.08	77.82	77.83
29	10757	17665	11742	11960	10757	17665	127.59	0.27	78.31	77.69
30	12757	18652	14019	14579	12757	18652	672.33	0.26	80.11	77.69
31	14757	21052	16398	17134	14757	21052	593.92	0.25	77.35	77.44
32	11076	23118	12103	13222	11076	23118	2.39	0.24	77.22	77.5
33	10254	17298	10843	11684	10254	17298	2.05	0.25	77.11	77.9
34	11850	23280	12807	13538	11850	23280	4.4	0.28	77.19	77.49
35	13450	25619	16107	16016	13450	25619	6.91	0.33	76.99	77.33
36	10889	31174	12591	13710	10889	31174	0.98	0.37	77.35	77.67
37	10076	19329	11162	11375	10076	19329	0.88	1.5	77.35	77.63
38	11476	31810	13636	13954	11476	31810	1.9	0.3	77.27	77.62
39	12685	38647	17010	16827	12685	38647	2.89	0.32	77.02	77.43
40	7891	9448	8295	8498	7891	9448	0.9	0.04	46.29	46.75
41	6451	8017	7160	7157	6451	8017	3.35	0.06	43.57	43.78
42	5599	7476	6638	6573	5599	7476	3.01	0.39	41.5	41.99
43	11322	13707	11628	11911	11322	13707	1.39	0.07	46.32	46.65
44	8364	10305	8734	9384	8364	10305	2.76	0.23	44.08	44.36
45	6934	9848	7960	7932	6934	9848	4	0.28	41.45	41.69
46	12460	14039	12822	12924	12460	14039	4.73	0.05	46.31	46.95
47	9605	11767	9990	10368	9605	11767	12.96	0.08	44.08	44.34
48	7367	10238	8412	8582	7367	10238	4.1	0.06	41.93	42.04
49	8721	11645	8989	8812	8721	11645	2.2	0.04	50.63	50.98
50	7390	10339	8821	8749	7390	10339	5.68	0.07	52.58	52.99
51	14189	17136	14332	14496	14189	17136	3.1	0.04	50.77	51.02
52	11368	13923	11573	12355	11368	13923	2	0.08	52.61	53.02
53	13617	17019	13910	13924	13617	17019	4.73	0.04	50.76	51.04
54	10493	12011	10793	11129	10493	12011	8.19	0.07	52.63	52.87
55	24575	33496	26877	26514	24575	33496	145.72	0.46	99.49	99.88
56	29374	38596	31856	31960	29374	38596	2785.24	0.46	99.64	99.83
57	\	50496	43637	43819	43637	50496	\	0.48	99.41	100.26
58	31886	42202	35246	35401	31886	42202	4809.23	0.45	99.5	100
59	24279	32807	27084	25733	24279	32807	1.84	0.48	99.25	99.53
60	28406	36483	30393	30304	28406	36483	28.17	0.6	99.12	99.3
61	36112	44183	38083	38105	36112	44183	5336	1.91	99.01	99.28
62	30293	39811	34203	33626	30293	39811	27.8	0.61	99.32	99.25
63	24277	40039	29894	27874	24277	40039	1.56	0.65	99.16	99.26
64	28347	44270	33377	33126	28347	44270	30.39	0.57	98.83	99.19
65	34999	49870	38683	38701	34999	49870	154.29	0.48	98.93	99.85
66	30099	40479	33626	33650	30099	40479	33.36	1.06	99.3	99.46
67	24289	36278	28431	27946	24289	36278	2.42	1.66	99.24	99.19
68	\	38278	31700	31556	31556	38278	\	0.46	99.28	99.19
69	\	46278	38659	38705	38659	46278	\	0.46	99.156	99.13
70	30348	42755	35193	34936	30348	42755	96.77	0.7	99.12	99.16

结论

从上面的结果中我们可以看到线性规划的结果是最好的，但是在一些情况下会出现计算时间爆炸的现象，因为数据产生的图形太复杂的原因吧，而且很多情况下运行时间比较小的原因是pulp库优化的结果。
模拟退火与局部搜索的结果都差不多，而且时间随规模的变化不会很大，结果大多不超过最优解的10%（以线性规划为标准）
贪心算法的结果则比较差，但是运行的时间比较短。
在模拟退火的算法过程可以有很多的花样，比如可以用工厂是否打开作为标准（打开为1，关闭为0），然后使用贪心算法进行计算出每种状态下的费用。以及可以融合遗传算法进行在模拟退火的算法过程可以有很多的花样，比如可以用工厂是否打开作为标准（打开为1，关闭为0），然后使用贪心算法进行计算出每种状态下的费用，但是这种方法得到的结果一般会比上面的结果差，因为工厂形成的解空间规模是： $2^{m}$ ，用户的是： $m^{n}$ 。另外还可以融合遗传算法进行求解。