算法设计与分析课作业【week15-16】Capacitated Facility Location Problem

本文探讨了设施选址与顾客分配的问题,旨在最小化设施启动成本和分配成本总和,同时确保设施容量不超过需求。介绍了使用贪心算法和模拟退火法解决此问题的方法,包括算法实现细节和对71个数据样例的实验结果对比。

摘要生成于 C知道 ,由 DeepSeek-R1 满血版支持, 前往体验 >

题目

Suppose there are n facilities and m customers. We wish to choose:
(1) which of the n facilities to open
(2) the assignment of customers to facilities
The objective is to minimize the sum of the opening cost and the assignment cost.
The total demand assigned to a facility must not exceed its capacity.

数据样例

第一行为设备的数量n和顾客的数量m。

接下来n行是每个设备的容量和启动的费用。

接下来的 m / 10 行 是每个顾客需要处理的数量需求。

接下来的 n * m / 10 行 是每个设备对于每个顾客所选择需要花费的费用。

我选择使用贪心算法和模拟退火法来寻找解。

贪心算法

贪心算法在该题目的使用方法是:在顾客能够选择的工厂中,选择顾客在这些工厂中需要花的最少费用的工厂。即在容量足够的工厂中,找顾客的assignment cost最少。这是一次贪心,原本我选择贪心的是选择assignment cost + open cost 最少的工厂,但发现结果却并不是很好,直接贪心 assignment cost而不管open cost 结果更好。所以便选择了这个方式,并且为了让贪心的结果能更好,我贪心了十次,每次通过逆序法来调整顾客选择工厂的顺序,然后重新贪心计算一次,在这十次中找出最好的一个来输出。

逆序法的使用:随机产生两个数p1, p2(p1 < p2), 代表两个被挑选的顾客,将从p1位置到p2的位置的顾客顺序反转,如:在 1 2 3 4 5 6 7 调到了2 和 6,逆序后变成 1 6 5 4 3 2 7。

贪心算法python代码如下:

import numpy as np
import time
import pandas as pd

# 定义顾客的类,记录顾客的需求和对于每个设备需要的assignment cost
class customer:
    def __init__(self, demand):
        self.demand = demand
    
    def setAssignmentCost(self, assignment_cost):
        self.assignment_cost = np.array(assignment_cost).copy()

class findMinCost:
    
    def __init__(self, filepath):
        self.readData(filepath)
    
    # 读取数据
    def readData(self, filepath):
        fr = open(filepath)
        lines = fr.readlines()
        firstline = lines[0].strip().split()
        self.facility_num = int(firstline[0])
        self.customer_num = int(firstline[1])
        line_num = 1
        # get data of facilitys
        self.facilitys = []
        while line_num <= self.facility_num:
            line = lines[line_num].strip().split()
            tempFacility = [int(line[0]), int(line[1])]
            self.facilitys.append(tempFacility)
            line_num += 1
        self.facilitys = np.array(self.facilitys)
        # get demand of customers
        hadread = line_num
        self.customers = []
        while line_num - hadread < self.customer_num / 10:
            # print(line_num)
            line = lines[line_num].strip().split()
            for demand in line:
                tempCustomer = customer(int(float(demand)))
                self.customers.append(tempCustomer)
            line_num += 1
        hadread = line_num
        # get assignment cost of customers
        everyFcostforCustomer = self.customer_num / 10
        readCount = 0
        assignment_cost = []
        fac_assignment_cost = []
        while line_num - hadread < self.customer_num * self.facility_num / 10:
            line = lines[line_num].strip().split()
            for as_cost in line:
                assignment_cost.append(int(float(as_cost)))
            readCount += 1
            if readCount == everyFcostforCustomer:
                fac_assignment_cost.append(assignment_cost)
                readCount = 0
                assignment_cost = []
            line_num += 1
        fac_assignment_cost = np.array(fac_assignment_cost)
        for index in range(self.customer_num):
            self.customers[index].setAssignmentCost(fac_assignment_cost[:, index])
    
    # 贪心算法
    def greedy(self, times):
        start_time = time.time()
        self.min_opencost = float('inf')
        self.min_assigncost = float('inf')
        self.fa_select = np.zeros(self.customer_num, np.int32)
        self.open_status = np.zeros(self.facility_num, np.int32)
        for t in range(times):
            
            skip = False
            open_cost = 0
            assign_cost = 0
            fa_capacity = self.facilitys[:, 0].copy()
            fa_select = np.zeros(self.customer_num, np.int32)
            open_status = np.zeros(self.facility_num, np.int32)
            # 遍历每一个顾客
            for cusIndex in range(len(self.customers)):
                # 找出能够选择的设备
                fa_canselect = []            
                for index in range(self.facility_num):
                    if fa_capacity[index] > self.customers[cusIndex].demand:
                        fa_canselect.append(index)

                # 若没有可以选择的设备,则本次没有解
                if len(fa_canselect) == 0:
                    skip = True
                    break
                self.cost_sum = self.customers[cusIndex].assignment_cost.copy()
                # 这次原本想利用 assignment cost 和open cost 的和进行贪心的代码
                # 效果不好,于是弃用
                # for index in range(self.facility_num):
                #     if open_status[index] == 0:
                #         self.cost_sum[index] += self.facilitys[index][1]
                # self.cost_sum = self.customers[cusIndex].assignment_cost + self.facilitys[:, 1]
                
                # print(self.cost_sum)
                # 找到费用最少的设备
                select = fa_canselect[0]
                min_cost = self.cost_sum[select]
                for index in fa_canselect:
                    if min_cost > self.cost_sum[index]:
                        select = index
                        min_cost = self.cost_sum[index]
                # 若设备没开,则计算open cost
                if open_status[select] == 0:
                    open_status[select] = 1
                    open_cost += self.facilitys[select][1]
                fa_capacity[select] -= self.customers[cusIndex].demand
                fa_select[cusIndex] = select
                assign_cost += self.customers[cusIndex].assignment_cost[select]
            if times != 1:
                self.changeCustomerOrder()
            if skip:
                continue
            # print("open cost: {}, assignment cost: {}, sum :{}".format(open_cost, assign_cost, open_cost + assign_cost))
            # print(open_status)
            # print(fa_select)
            # 记录花费最少的一次
            if open_cost + assign_cost < self.min_opencost + self.min_assigncost:
                self.min_opencost = open_cost
                self.min_assigncost = assign_cost
                self.open_status = open_status.copy()
                self.fa_select = fa_select.copy()
            
        print("the best one: open cost: {}, assignment cost: {}, sum :{}".format(self.min_opencost, self.min_assigncost, self.min_opencost + self.min_assigncost))
        print(self.open_status)
        print(self.fa_select)
        end_time = time.time()
        return self.min_opencost + self.min_assigncost, end_time - start_time
    
    # 逆序法获得新的顾客数组
    def changeCustomerOrder(self):
        pos1 = 0
        pos2 = 0
        while pos1 == pos2 or (pos1 == 0 and pos2 == self.customer_num - 1):
            pos1 = int(np.random.uniform(0, self.customer_num))
            pos2 = int(np.random.uniform(0, self.customer_num))
            if pos1 > pos2:
                pos1, pos2 = pos2, pos1
        while pos1 < pos2:
            self.customers[pos1], self.customers[pos2] = self.customers[pos2], self.customers[pos1]
            pos1 += 1
            pos2 -= 1


if __name__ == "__main__":
    file_names = []
    costs = []
    times = []
    for i in range(71):
        file_names.append("p" + str(i + 1))
        test = findMinCost("Instances/p" + str(i + 1))
        print("p" + str(i + 1) + ": ")
        # 贪心10次
        cost, t = test.greedy(10)
        costs.append(cost)
        times.append(t)
    dataframe = pd.DataFrame({'file': file_names, 'result': costs, 'time': times})
    dataframe.to_csv("greedy_result.csv", index=False, sep=',')

模拟退火法

模拟退火法是通过设定初始温度,每一个温度内多次寻找新解,如果找到了更好的解,便更新当前解,如果得到的解没有比当前好,则通过概率判断是否选择该坏解,以此跳出局部最优,达到全局最优解。

模拟退火法中,寻找新解的方式刚开始是随机找到一个顾客,改变该顾客选择的工厂,计算更换后的open cost + assignment cost是否变小来获得新解,但发现这样的方式并没有找到很好的结果。所以换了另一种找新解的方式,随机选择两个顾客,交换他们的选择的工厂,这样只需计算assignment cost的差别,通过这种方式来获得新解的结果较好,跟贪心算法的结果比较不相上下。

为了优化结果,我在使用前首先用贪心算法来获得初解,以此来争取获得更低的费用。

模拟退火法的python代码如下:

import numpy as np
import time
import pandas as pd

class customer:
    def __init__(self, demand):
        self.demand = demand
    
    def setAssignmentCost(self, assignment_cost):
        self.assignment_cost = np.array(assignment_cost).copy()

class findMinCost:
    
    def __init__(self, filepath):
        self.readData(filepath)

    def readData(self, filepath):
        fr = open(filepath)
        lines = fr.readlines()
        # print(len(lines))
        firstline = lines[0].strip().split()
        self.facility_num = int(firstline[0])
        self.customer_num = int(firstline[1])
        # print(self.facility_num)
        # print(self.customer_num)
        line_num = 1
        # get data of facilitys
        self.facilitys = []
        while line_num <= self.facility_num:
            line = lines[line_num].strip().split()
            tempFacility = [int(line[0]), int(line[1])]
            self.facilitys.append(tempFacility)
            line_num += 1
        self.facilitys = np.array(self.facilitys)
        # get demand of customers
        hadread = line_num
        self.customers = []
        while line_num - hadread < self.customer_num / 10:
            # print(line_num)
            line = lines[line_num].strip().split()
            for demand in line:
                tempCustomer = customer(int(float(demand)))
                self.customers.append(tempCustomer)
            line_num += 1
        hadread = line_num
        # get assignment cost of customers
        everyFcostforCustomer = self.customer_num / 10
        readCount = 0
        assignment_cost = []
        fac_assignment_cost = []
        while line_num - hadread < self.customer_num * self.facility_num / 10:
            # print(line_num - hadread)
            line = lines[line_num].strip().split()
            for as_cost in line:
                assignment_cost.append(int(float(as_cost)))
            readCount += 1
            if readCount == everyFcostforCustomer:
                fac_assignment_cost.append(assignment_cost)
                readCount = 0
                assignment_cost = []
            line_num += 1
        fac_assignment_cost = np.array(fac_assignment_cost)
        for index in range(self.customer_num):
            self.customers[index].setAssignmentCost(fac_assignment_cost[:, index])
        # print(self.facilitys[1].openging_cost)
        # print(self.customers[1].assignment_cost)

    def greedy(self):
        self.open_cost = 0
        self.assign_cost = 0
        self.fa_capacity = self.facilitys[:, 0].copy()
        self.fa_select = np.zeros(self.customer_num, np.int32)
        self.open_status = np.zeros(self.facility_num, np.int32)
        for cusIndex in range(len(self.customers)):

            fa_canselect = []            
            for index in range(self.facility_num):
                if self.fa_capacity[index] > self.customers[cusIndex].demand:
                    fa_canselect.append(index)

            self.cost_sum = self.customers[cusIndex].assignment_cost.copy()
            
            select = fa_canselect[0]
            min_cost = self.cost_sum[select]
            for index in fa_canselect:
                if min_cost > self.cost_sum[index]:
                    select = index
                    min_cost = self.cost_sum[index]
            
            if self.open_status[select] == 0:
                self.open_status[select] = 1
                self.open_cost += self.facilitys[select][1]
            self.fa_capacity[select] -= self.customers[cusIndex].demand
            self.fa_select[cusIndex] = select
            self.assign_cost += self.customers[cusIndex].assignment_cost[select]

        # self.printData()
    
    # 输出结果
    def printData(self):
        print(" open cost: {}, assignment cost: {}, sum :{}".format(self.open_cost, self.assign_cost, self.open_cost + self.assign_cost))
        print(self.open_status)
        print(self.fa_select)
    
    # 通过选择一个客户来更换设备的SA,结果不好所以弃用
    def SA1(self, temperature_limit, temperature, times, raduceRate):
        self.initSAsalution()
        while temperature > temperature_limit:
            times *= 1.01
            temperature_times = times
            while temperature_times > 0:
                cus_choose = int(np.random.uniform(0, self.customer_num))
                fac_choose = int(np.random.uniform(0, self.facility_num))
                while fac_choose == self.fa_select[cus_choose] or self.fa_capacity[fac_choose] < self.customers[cus_choose].demand:
                    fac_choose = int(np.random.uniform(0, self.facility_num))
                origin_fac_choose = int(self.fa_select[cus_choose])
                opencost_offset = 0
                assigncost_offset = 0
                if self.open_status[fac_choose] == 0:
                    opencost_offset += self.facilitys[fac_choose][1]
                if self.fa_capacity[origin_fac_choose] + self.customers[cus_choose].demand == self.facilitys[origin_fac_choose][0]:
                    opencost_offset -= self.facilitys[origin_fac_choose][1]
                assigncost_offset = self.customers[cus_choose].assignment_cost[fac_choose] - self.customers[cus_choose].assignment_cost[origin_fac_choose]
                offset = opencost_offset + assigncost_offset
                print("opencost-of: {}, assigncost-of: {}, sum-of: {}".format(opencost_offset, assigncost_offset, offset))
                if opencost_offset + assigncost_offset < 0 or np.exp(-offset / temperature) >= np.random.rand():
                    if self.open_status[fac_choose] == 0:
                        self.open_status[fac_choose] = 1
                    if self.fa_capacity[origin_fac_choose] + self.customers[cus_choose].demand == self.facilitys[origin_fac_choose][0]:
                        self.open_status[origin_fac_choose] = 0
                    self.fa_select[cus_choose] = fac_choose
                    self.open_cost += opencost_offset
                    self.assign_cost += assigncost_offset
                temperature_times -= 1
                print("T: {}, opencost: {}, assigncost: {}, sum: {}".format(temperature, self.open_cost, self.assign_cost, self.open_cost+self.assign_cost))
            temperature *= raduceRate
        self.printData()
    
    # 通过交换2个客户选择的设备SA
    def SA2(self, temperature_limit, temperature, times, raduceRate):
        start_time = time.time()
        self.initSAsalution()
        while temperature > temperature_limit:
            times *= 1.01
            temperature_times = times
            while temperature_times > 0:
                cus1 = 0
                cus2 = 0
                while cus1 == cus2:
                    cus1 = int(np.random.uniform(0, self.customer_num))
                    cus2 = int(np.random.uniform(0, self.customer_num))
                if (self.fa_capacity[self.fa_select[cus1]] + self.customers[cus1].demand < self.customers[cus2].demand or 
                    self.fa_capacity[self.fa_select[cus2]] + self.customers[cus2].demand < self.customers[cus1].demand):
                    continue
                # 计算assigncost 的差值
                assigncost_offset = (self.customers[cus1].assignment_cost[self.fa_select[cus2]]
                                    + self.customers[cus2].assignment_cost[self.fa_select[cus1]]
                                    - self.customers[cus1].assignment_cost[self.fa_select[cus1]]
                                    - self.customers[cus2].assignment_cost[self.fa_select[cus2]])
                # print("assigncost_offset: {}".format(assigncost_offset))
                
                #判断是否可以接受新解
                if assigncost_offset < 0 or np.exp(-assigncost_offset / temperature) >= np.random.rand():
                    self.assign_cost += assigncost_offset
                    self.fa_capacity[self.fa_select[cus1]] += self.customers[cus1].demand - self.customers[cus2].demand
                    self.fa_capacity[self.fa_select[cus2]] += self.customers[cus2].demand - self.customers[cus1].demand
                    self.fa_select[cus1], self.fa_select[cus2] = self.fa_select[cus2], self.fa_select[cus1]
                temperature_times -= 1
                # print("T: {}, opencost: {}, assigncost: {}, sum: {}".format(temperature, self.open_cost, self.assign_cost, self.open_cost+self.assign_cost))
            temperature *= raduceRate
        self.printData()
        end_time = time.time()
        return self.open_cost + self.assign_cost, end_time - start_time
    
    # 贪心获得初解                
    def initSAsalution(self):
        self.greedy()


if __name__ == "__main__":
    file_names = []
    costs = []
    times = []
    for i in range(71):
        file_names.append("p" + str(i + 1))
        test = findMinCost("Instances/p" + str(i + 1))
        print("p" + str(i + 1) + ": ")
        cost, t = test.SA2(5, 100, 3000, 0.97)
        costs.append(cost)
        times.append(t)
    dataframe = pd.DataFrame({'file': file_names, 'result': costs, 'time': times})
    dataframe.to_csv("sa_result2.csv", index=False, sep=',')

最终71个样例的结果:

样例贪心算法结果模拟退火法结果贪心算法时间模拟退火法时间
p1932293390.0087.268
p2812680390.0087.428
p310126100250.0097.04
p411942120250.0097.124
p5937591700.0217.378
p6806179180.0067.727
p71006198560.0057.697
p812061118560.0068.165
p9904090400.0056.665
p10772678900.0066.571
p11972697260.0066.837
p1211726117260.0086.768
p1312032120320.0096.495
p14918091800.016.652
p1513180131800.0097.55
p1617180171800.0096.918
p1712032120320.0096.787
p18918091800.0087.162
p1913180131800.0097.179
p2017180171800.0097.07
p2112032120320.0096.673
p22918091800.0096.92
p2313180131800.0096.58
p2417180171800.0097.014
p2519207191950.0356.656
p2616106161490.0346.712
p2721602215290.0346.808
p2826972269620.0356.651
p2919305191640.0347.219
p3016239161180.0346.929
p3121639215330.0357.001
p3227039269190.0356.737
p3319055191500.0456.985
p3415989161690.0427.014
p3521389215170.0436.809
p3626789268930.0476.865
p3719055192290.0526.763
p3815989162020.0436.801
p3921389214940.0436.623
p4026789268820.0426.608
p41711771430.016.937
p42995799720.0156.512
p4312448125380.0166.575
p44726872180.0347.118
p45984898480.0456.548
p4612639126390.0416.585
p47660264120.0387.226
p48904490720.046.729
p4912420124580.046.879
p50996196030.0417.49
p5111351113920.0547.033
p5210364103750.047.374
p5312470123940.0487.08
p54987299160.047.247
p5511934119460.0487.129
p5623882241320.1166.877
p5732882329430.1226.454
p5853882540000.1026.64
p5939121392010.1146.6
p6023882240100.1146.522
p6132882330310.1076.696
p6253882539100.1036.467
p6339121391210.1016.452
p6423882239540.1156.613
p6532882328820.1126.54
p6653882540710.1166.774
p6739671397540.126.643
p6823882238870.1016.602
p6932882330000.1116.844
p7053882539430.1146.585
p7139121391830.1136.571

输出样式

其中每一个样例的选择都会在控制台输出,由于样例过多,这里就不贴了,代码一跑就可以看到了。

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值