题目
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个样例的结果:
样例 | 贪心算法结果 | 模拟退火法结果 | 贪心算法时间 | 模拟退火法时间 |
p1 | 9322 | 9339 | 0.008 | 7.268 |
p2 | 8126 | 8039 | 0.008 | 7.428 |
p3 | 10126 | 10025 | 0.009 | 7.04 |
p4 | 11942 | 12025 | 0.009 | 7.124 |
p5 | 9375 | 9170 | 0.021 | 7.378 |
p6 | 8061 | 7918 | 0.006 | 7.727 |
p7 | 10061 | 9856 | 0.005 | 7.697 |
p8 | 12061 | 11856 | 0.006 | 8.165 |
p9 | 9040 | 9040 | 0.005 | 6.665 |
p10 | 7726 | 7890 | 0.006 | 6.571 |
p11 | 9726 | 9726 | 0.006 | 6.837 |
p12 | 11726 | 11726 | 0.008 | 6.768 |
p13 | 12032 | 12032 | 0.009 | 6.495 |
p14 | 9180 | 9180 | 0.01 | 6.652 |
p15 | 13180 | 13180 | 0.009 | 7.55 |
p16 | 17180 | 17180 | 0.009 | 6.918 |
p17 | 12032 | 12032 | 0.009 | 6.787 |
p18 | 9180 | 9180 | 0.008 | 7.162 |
p19 | 13180 | 13180 | 0.009 | 7.179 |
p20 | 17180 | 17180 | 0.009 | 7.07 |
p21 | 12032 | 12032 | 0.009 | 6.673 |
p22 | 9180 | 9180 | 0.009 | 6.92 |
p23 | 13180 | 13180 | 0.009 | 6.58 |
p24 | 17180 | 17180 | 0.009 | 7.014 |
p25 | 19207 | 19195 | 0.035 | 6.656 |
p26 | 16106 | 16149 | 0.034 | 6.712 |
p27 | 21602 | 21529 | 0.034 | 6.808 |
p28 | 26972 | 26962 | 0.035 | 6.651 |
p29 | 19305 | 19164 | 0.034 | 7.219 |
p30 | 16239 | 16118 | 0.034 | 6.929 |
p31 | 21639 | 21533 | 0.035 | 7.001 |
p32 | 27039 | 26919 | 0.035 | 6.737 |
p33 | 19055 | 19150 | 0.045 | 6.985 |
p34 | 15989 | 16169 | 0.042 | 7.014 |
p35 | 21389 | 21517 | 0.043 | 6.809 |
p36 | 26789 | 26893 | 0.047 | 6.865 |
p37 | 19055 | 19229 | 0.052 | 6.763 |
p38 | 15989 | 16202 | 0.043 | 6.801 |
p39 | 21389 | 21494 | 0.043 | 6.623 |
p40 | 26789 | 26882 | 0.042 | 6.608 |
p41 | 7117 | 7143 | 0.01 | 6.937 |
p42 | 9957 | 9972 | 0.015 | 6.512 |
p43 | 12448 | 12538 | 0.016 | 6.575 |
p44 | 7268 | 7218 | 0.034 | 7.118 |
p45 | 9848 | 9848 | 0.045 | 6.548 |
p46 | 12639 | 12639 | 0.041 | 6.585 |
p47 | 6602 | 6412 | 0.038 | 7.226 |
p48 | 9044 | 9072 | 0.04 | 6.729 |
p49 | 12420 | 12458 | 0.04 | 6.879 |
p50 | 9961 | 9603 | 0.041 | 7.49 |
p51 | 11351 | 11392 | 0.054 | 7.033 |
p52 | 10364 | 10375 | 0.04 | 7.374 |
p53 | 12470 | 12394 | 0.048 | 7.08 |
p54 | 9872 | 9916 | 0.04 | 7.247 |
p55 | 11934 | 11946 | 0.048 | 7.129 |
p56 | 23882 | 24132 | 0.116 | 6.877 |
p57 | 32882 | 32943 | 0.122 | 6.454 |
p58 | 53882 | 54000 | 0.102 | 6.64 |
p59 | 39121 | 39201 | 0.114 | 6.6 |
p60 | 23882 | 24010 | 0.114 | 6.522 |
p61 | 32882 | 33031 | 0.107 | 6.696 |
p62 | 53882 | 53910 | 0.103 | 6.467 |
p63 | 39121 | 39121 | 0.101 | 6.452 |
p64 | 23882 | 23954 | 0.115 | 6.613 |
p65 | 32882 | 32882 | 0.112 | 6.54 |
p66 | 53882 | 54071 | 0.116 | 6.774 |
p67 | 39671 | 39754 | 0.12 | 6.643 |
p68 | 23882 | 23887 | 0.101 | 6.602 |
p69 | 32882 | 33000 | 0.111 | 6.844 |
p70 | 53882 | 53943 | 0.114 | 6.585 |
p71 | 39121 | 39183 | 0.113 | 6.571 |
输出样式
其中每一个样例的选择都会在控制台输出,由于样例过多,这里就不贴了,代码一跑就可以看到了。