文章目录
- Capacitated Facility Location Problem
- [更多技术博客 http://vilins.top/](http://vilins.top/)
- 题目
- 分析
- 实验工具
- 解决算法
- 模拟退火法
- 禁忌搜索
- 算法局部设计细节
- 模拟退火法
- 禁忌搜索
- 结果整理
- 两种做法的对比
- 原始输出结果
- 程序运行说明
- 源码
- github
- [更多技术博客 http://vilins.top/](http://vilins.top/)
Capacitated Facility Location Problem
更多技术博客 http://vilins.top/
题目
- 要求
Suppose there are n facilities and m customers. We wish to choose:- which of the n facilities to open.
- 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.
- 输入格式
- 结果要求
分析
选址问题在生产生活,物流,甚至军事中都有这非常广泛的应用,以成为运筹学中经典的问题之一,其目的是确定一个或者多个设施的地址,使得判断标准获得优化的同时,及时向顾客提供其所需的产品和服务。工厂选址问题不仅影响到工厂的利润和市场竞争力,甚至决定了工厂的生死存亡。所以,工厂选址的研究具有重要意义。
这类问题在论文中出现得比较多,而网上现成的以实现的办法几乎没有,很多都是从理论层面上分析解决这类问题的可行性,而没有实操去证明其可行性。由于这类问题属于优化类问题,在线性时间很难求得其最优解,所以我们选择算法的时候一般选择优化类算法,动态规划之类的算法是行不通的,一般用于解决这类问题的办法有蚁群算法,遗传算法,贪婪算法,模拟退火法,禁忌搜索算法,近似算法等。
具体到这个题目,我们首先分析其限制的条件,这是一类有容量限制的工厂问题,具体限制如下:
- 每个工厂开放有一定的花费opening cost。
- 有一定的顾客数量,顾客对工厂具有需求,也就是demand,任何工厂都可以向顾客提供需求,但这里的限制是每个顾客的需求只需一个工厂来提供,这也是对问题的简化。
- 每个工厂具有一定的容量限制,工厂向顾客提供的需求不能超过该工厂的容量限制。
- 每个顾客到工厂有一定的路程,称为assignment cost。
- 花费要尽可能少
总花费(total cost) = 工厂开放花费(opening cost) + 顾客分配到某个工厂以后的路程花费(assginment cost)实验工具
语言:python
工具:Pycharm
解决算法
模拟退火法
- 算法框架
- 细节实现
模拟退火法的设计有几个关键的因素:初始解
这里模拟退火法目的就是为了跳出局部最优解的,所以初始解重要,但又没有那么重要,初始解设计得好的话,模拟退火法可以很快得到收敛,使用的时间相对较短一点;但初始解很差,模拟退火法求得的解也没什么差距,只不过需要的时间需要的长一些。这里我尝试过几种不同的初始化的方法,包括一下几种,但测试后发现,加入贪婪因素后求解的效率表现得很差,反而最后两种初始化的效果不错,而最后一种表现得最好,个人猜测是这个解的形式与理论上最优解类似。- 用全贪婪获得初始解
- 用部分贪婪获得初始解(百分比)
- 简单的填满前面几个工厂
- 完全随机化选择部分工厂开放,但不让所有的工厂开放
新解的产生
新解的产生设计到我们搜索范围的大小,这里有两个关键因素,一个是重复产生新解的次数,一个是新解产生的方法;有足够多的次数来产生邻居,那么我们找到邻居的最优解的概率就越高,这里设置为200次就可以了;而产生邻居的方式要根据具体问题而定,这个问题我采取随机产生一个顾客和工厂的下标,将该顾客分配给该工厂,有一下情况发生:- 工厂没有开放的情况并且工厂的容量满足顾客需求,将该顾客分配给该工厂。
- 如果工厂已经开放,并且其剩下的容量足够提供顾客的需求,那么将顾客分配给该工厂。
- 如果工厂容量不足以提供顾客需求,重复生成随机数操作,直到找到一个新解为止。
初温
初温的设置与模拟退火跳出局部有关系,在温度较高的时候,解跳出局部最优的可能性就越大,所以初温应该有一定的高度,让解跳出局部,这里设置的初温是120度。降温速率
降温速率与邻居的效率有关,降温幅度底,那么找到邻居最优解的邻居解的产生可能性就越大,参考降温速率为0.95-0.99。终止温度
终止温度与最终找到的最优解有关,温度足够低,那么当前解接受差解的概率就会很低,终止温度底,那么保证最终获得的解是一个不错的解。
禁忌搜索
- 算法框架
- 细节实现
同样,禁忌搜索也有几个关键的因素初始解
这里的初始解的产生很模拟退火采用的是同一种初始化的方法,效果表现得还不错,具体描述见上。重复次数
这里重复的次数是保证在足够多的解空间中对结果进行搜索,实验证明,这里我们需要对解空间的大小以及搜索需要的时间进行权衡;解空间大的时候,搜索需要的时间也必然会变得长一些,找出的解也可能会好一些;而解空间小,搜索的花费没有那么大,找出的最优解可能没有那么好,权衡之下,重复20000次是比较好的选择。新解的产生
新解产生根模拟退火新解产生一致,具体见上。禁忌表的设计
原理的示意图如下
禁忌表的引入可以有效避免环的产生,避免那些冗余得搜索,根据这样的原理,那么很自然要设计的一个参数就是在几步以内不可以搜索同样得解呢?也就是在几步以内这个解要保存在禁忌表中,当有同样的解在禁忌表发现的时候,我们直接跳过这个解,不去搜索它。当搜索了足够长的步数以后再将这个解从禁忌表中移除,那么这个禁忌表的长度如何设计呢?有下面两种策略:- 固定长度20
- 根据解空间而变化
可能是由于数据集的适应度比较好,所以固定长度20的时候,算法的表现已经足够好了。
当如果根据上述方式来设计的话,会存在一个明显的问题,那就是时间的花费,禁忌搜索本来需要的时间就比较长了,如果中间还要加上搜索禁忌表的花费的话,那么所用的时间复杂度将会变得异常高,对于这样的情况,我们可以采用这样的方法进行优化,
这样只需要O(1)的时间就可以查询到该解是否在禁忌表中,当当前维护的count>表中对应元素的值的时候,表示我们可以重新搜索这个解,否则这个解在禁忌的状态。
算法局部设计细节
模拟退火法
初始解的产生,采用完全随机化选择部分工厂开放,但不让所有的工厂开放,符合一般最优化表现形式,表现最好。
def init_solution(self): for i in range(self.customer_num): self.current_state_of_customer.append(-1) for i in self.demand: self.total_demand += i self.current_capacity = self.capacity.copy() total_cost = 0 temp = [] while total_cost < self.total_demand: r = random.randint(0, self.facility_num - 1) if r not in temp: temp.append(r) total_cost += self.capacity[r] index = 0 for i in range(self.customer_num): if self.current_capacity[temp[index]] - self.demand[i] >= 0: self.current_capacity[temp[index]] -= self.demand[i] self.current_state_of_customer[i] = temp[index] else: index += 1 i -= 1 # self.current_state_of_customer[i] = index # self.current_capacity[index] -= self.demand[i] self.current_cost = self.calculate_cost(self.current_state_of_customer) self.new_capacity = self.current_capacity.copy() self.new_state_of_customer = self.current_state_of_customer.copy()初始化解,利用局部贪婪,可以调节贪婪的比例,但效果不太好
def init_solution(self): for i in range(self.customer_num): self.current_state_of_customer.append(-1) for i in self.demand: self.total_demand += i self.current_capacity = self.capacity.copy() for i in range(self.customer_num): ran_select = random.random() if ran_select > 0.2: min = 100000 index = 0 for j in range(self.facility_num): if min > self.assignment[j][i]: min = self.assignment[j][i] index = j if self.current_capacity[index] >= self.demand[i]: self.current_capacity[index] -= self.demand[i] self.current_state_of_customer[i] = index else: tag = True while tag: ran = random.randint(0, self.facility_num - 1) if self.current_capacity[ran] >= self.demand[i]: self.current_capacity[ran] -= self.demand[i] self.current_state_of_customer[i] = ran tag = False else: tag_out = True while tag_out: ran = random.randint(0, self.facility_num - 1) if self.current_capacity[ran] >= self.demand[i]: self.current_capacity[ran] -= self.demand[i] self.current_state_of_customer[i] = ran tag_out = False self.current_cost = self.calculate_cost(self.current_state_of_customer) self.new_capacity = self.current_capacity.copy() self.new_state_of_customer = self.current_state_of_customer.copy() print('cost ', self.current_cost) print('current cap ', self.current_capacity) print('current cu ', self.current_state_of_customer)产生新解的方式
def gen_new_solution(self): tag = True self.new_capacity = self.current_capacity.copy() self.new_state_of_customer = self.current_state_of_customer.copy() while tag: ran_customer = random.randint(0, self.customer_num - 1) ran_facility = random.randint(0, self.facility_num - 1) # print("ran ", ran_facility, ran_customer) if self.new_capacity[ran_facility] >= self.demand[ran_customer]: tag = False self.new_capacity[ran_facility] -= self.demand[ran_customer] origin = self.new_state_of_customer[ran_customer] self.new_state_of_customer[ran_customer] = ran_facility if origin == -1: break self.new_capacity[origin] += self.demand[ran_customer] else: continue退火过程
def simulated_annealing(self): global fp_anneal time_start = time.time() while self.T > 0.01: i = 0 while i < self.repeat: i += 1 self.gen_new_solution() self.current_cost = self.calculate_cost(self.current_state_of_customer) new_cost = self.calculate_cost(self.new_state_of_customer) if new_cost < self.current_cost: self.current_capacity = self.new_capacity.copy() self.current_state_of_customer = self.new_state_of_customer.copy() # print("current cost ", new_cost) else: num = math.exp((self.current_cost - new_cost) / self.T) ran = random.random() # print("num ", num) # print("ran ", ran) if num >= ran: self.current_capacity = self.new_capacity.copy() self.current_state_of_customer = self.new_state_of_customer.copy() # print("current cost ", new_cost) self.T *= self.cooling_rate time_end = time.time() used_time = time_end - time_start fp_anneal.write('time ' + str(used_time) + '\n') fp_anneal.write('cost ' + str(self.current_cost) + '\n') temp = [] result_state = [] for i in range(self.customer_num): if self.current_state_of_customer[i] not in temp: temp.append(self.new_state_of_customer[i]) temp.sort() for i in range(self.facility_num): result_state.append(0) for i in temp: result_state[i] = 1 # print('cost ' + str(self.current_cost)) # print(result_state) # print(self.current_state_of_customer) # print(self.current_capacity) fp_anneal.write('facility state ' + str(result_state) + '\n') fp_anneal.write('current state of customer ' + str(self.current_state_of_customer) + '\n\n')禁忌搜索
新解产生与模拟退火类似
邻居的产生
def gen_nei_solution(self): tag = True self.nei_capacity = self.current_capacity.copy() self.nei_state_of_customer = self.current_state_of_customer.copy() ran_customer = 0 origin = 0 while tag: ran_customer = random.randint(0, self.customer_num - 1) ran_facility = random.randint(0, self.facility_num - 1) # print("ran ", ran_facility, ran_customer) if self.nei_capacity[ran_facility] >= self.demand[ran_customer]: tag = False self.nei_capacity[ran_facility] -= self.demand[ran_customer] origin = self.nei_state_of_customer[ran_customer] self.nei_state_of_customer[ran_customer] = ran_facility if origin == -1: break self.nei_capacity[origin] += self.demand[ran_customer] else: continue return ran_customer, origin禁忌搜索过程
def tabu_search(self): global fp_tabu self.nei_state_of_customer = self.current_state_of_customer.copy() self.nei_capacity = self.current_capacity.copy() # print(self.current_state_of_customer) # print(self.nei_capacity) self.tabu_list = [[0 for i in range(self.customer_num)]for i in range(self.facility_num)] # print(tabu_list) repeat = 30000 count = 0 nei_count = int(self.customer_num/4) best_solution_customer = self.current_state_of_customer.copy() best_solution_capacity = self.current_capacity.copy() best_cost = self.calculate_cost(best_solution_customer) time_start = time.time() while count < repeat: count += 1 i = nei_count self.gen_new_solution() new_cost = self.calculate_cost(self.new_state_of_customer) ran_customer = 0 origin = 0 while i > 1: i -= 1 ran_customer, origin = self.gen_nei_solution() nei_cost = self.calculate_cost(self.nei_state_of_customer) # print("nei ", nei_cost) # print("new ", new_cost) if nei_cost < new_cost: self.new_state_of_customer = self.nei_state_of_customer.copy() self.new_capacity = self.nei_capacity.copy() new_cost = nei_cost if new_cost < best_cost: # print("new ", new_cost) best_solution_customer = self.new_state_of_customer.copy() best_solution_capacity = self.new_capacity.copy() best_cost = new_cost if self.tabu_list[origin][ran_customer] < count: self.tabu_list[origin][ran_customer] = count + 20 self.current_state_of_customer = self.new_state_of_customer.copy() self.current_capacity = self.new_capacity.copy() time_end = time.time() used_time = time_end - time_start temp = [] result_state = [] for i in range(self.customer_num): if best_solution_customer[i] not in temp: temp.append(best_solution_customer[i]) temp.sort() for i in range(self.facility_num): result_state.append(0) for i in temp: result_state[i] = 1 fp_tabu.write('time ' + str(used_time) + '\n') fp_tabu.write('cost ' + str(best_cost) + '\n') fp_tabu.write('facility state ' + str(result_state) + '\n') fp_tabu.write('customer state ' + str(best_solution_customer) + '\n\n') print(best_cost) print(result_state) print(best_solution_customer) print(best_solution_capacity)结果整理
| file | result | time | result | time |
|---|---|---|---|---|
| 模拟退火法 | 模拟退火法 | 禁忌搜索 | 禁忌搜索 | |
| p1 | 9033 | 5.123326063156128 | 8985 | 9.848682165145874 |
| p2 | 8005 | 5.528823137283325 | 7940 | 9.718816757202148 |
| p3 | 9314 | 5.139236688613892 | 9678 | 9.466161012649536 |
| p4 | 11628 | 5.307790279388428 | 10954 | 9.69756555557251 |
| p5 | 9131 | 5.344730377197266 | 9137 | 10.06436276435852 |
| p6 | 7781 | 5.396562337875366 | 7907 | 9.97243070602417 |
| p7 | 9577 | 5.525223731994629 | 9625 | 10.066269397735596 |
| p8 | 11439 | 5.533200263977051 | 11433 | 10.113304615020752 |
| p9 | 8593 | 4.830053329467773 | 9031 | 8.910385608673096 |
| p10 | 7648 | 4.938756942749023 | 7724 | 9.172668695449829 |
| p11 | 9226 | 5.148226976394653 | 9062 | 8.851669073104858 |
| p12 | 11033 | 4.917876720428467 | 10277 | 8.693775177001953 |
| p13 | 9312 | 5.169209718704224 | 9604 | 10.038275718688965 |
| p14 | 7164 | 5.242977142333984 | 8204 | 10.189182996749878 |
| p15 | 8964 | 5.441417932510376 | 10335 | 9.89920711517334 |
| p16 | 11925 | 5.371662139892578 | 12496 | 9.891655445098877 |
| p17 | 9278 | 5.450403690338135 | 9795 | 9.87352442741394 |
| p18 | 7781 | 5.39055061340332 | 8265 | 10.18008804321289 |
| p19 | 9267 | 5.236989974975586 | 10408 | 9.878618240356445 |
| p20 | 11758 | 5.060497522354126 | 12593 | 9.863529920578003 |
| p21 | 9295 | 4.858996391296387 | 9103 | 9.683191061019897 |
| p22 | 7738 | 5.125286340713501 | 8088 | 9.918416261672974 |
| p23 | 9654 | 5.08836030960083 | 10070 | 9.597949266433716 |
| p24 | 11263 | 4.70444393157959 | 11800 | 9.354219436645508 |
| p25 | 12281 | 12.890488147735596 | 13589 | 67.44565677642822 |
| p26 | 11623 | 13.553741693496704 | 11254 | 70.20729756355286 |
| p27 | 13574 | 13.028119087219238 | 13377 | 66.27394366264343 |
| p28 | 15729 | 13.24556589126587 | 17572 | 67.49193954467773 |
| p29 | 13510 | 13.5068519115448 | 13683 | 72.6292462348938 |
| p30 | 11894 | 13.90679669380188 | 11797 | 72.82930493354797 |
| p31 | 14823 | 13.666414260864258 | 14412 | 72.44930481910706 |
| p32 | 17733 | 12.95035457611084 | 16708 | 72.29048991203308 |
| p33 | 12536 | 12.981273174285889 | 12343 | 67.8763542175293 |
| p34 | 11582 | 13.536755084991455 | 11142 | 69.67934989929199 |
| p35 | 14487 | 13.469933271408081 | 13099 | 68.77509903907776 |
| p36 | 15688 | 13.147795677185059 | 17320 | 67.9144344329834 |
| p37 | 12200 | 12.65013074874878 | 12248 | 67.00625324249268 |
| p38 | 11283 | 13.485923528671265 | 11170 | 67.53133893013 |
| p39 | 13480 | 13.172777891159058 | 12862 | 65.43386697769165 |
| p40 | 15682 | 13.089987754821777 | 15569 | 64.73360013961792 |
| p41 | 7008 | 8.108281373977661 | 7906 | 24.470109462738037 |
| p42 | 6932 | 7.624598264694214 | 6312 | 21.45367980003357 |
| p43 | 6529 | 6.635222673416138 | 6733 | 17.159975051879883 |
| p44 | 7206 | 8.201093912124634 | 8197 | 27.290011882781982 |
| p45 | 7384 | 7.414166450500488 | 7489 | 21.240010738372803 |
| p46 | 6166 | 6.826710224151611 | 7208 | 17.667112588882446 |
| p47 | 6434 | 8.080354928970337 | 7249 | 25.64135766029358 |
| p48 | 6299 | 7.8948752880096436 | 6438 | 22.5309157371521 |
| p49 | 6329 | 6.82571268081665 | 6764 | 17.12479257583618 |
| p50 | 9313 | 8.631909847259521 | 10880 | 29.890727758407593 |
| p51 | 7913 | 8.980969190597534 | 8939 | 31.077434062957764 |
| p52 | 9347 | 8.536164283752441 | 9726 | 33.53241181373596 |
| p53 | 9946 | 9.425786256790161 | 11141 | 32.7076473236084 |
| p54 | 9178 | 8.499296188354492 | 10560 | 33.5239634513855 |
| p55 | 8465 | 9.157528400421143 | 8318 | 32.04251146316528 |
| p56 | 22713 | 21.777714252471924 | 23903 | 144.6936798095703 |
| p57 | 33581 | 19.876826524734497 | 30914 | 142.13889384269714 |
| p58 | 56501 | 18.72487449645996 | 46049 | 138.9957413673401 |
| p59 | 32738 | 20.65970516204834 | 38781 | 140.89498019218445 |
| p60 | 23201 | 23.380425214767456 | 24435 | 128.45606589317322 |
| p61 | 31446 | 18.320961236953735 | 29426 | 122.40536570549011 |
| p62 | 63518 | 16.87701916694641 | 38882 | 126.1437463760376 |
| p63 | 32027 | 18.300015449523926 | 34958 | 120.8425030708313 |
| p64 | 23347 | 21.732863426208496 | 24081 | 118.14805221557617 |
| p65 | 28483 | 17.28874921798706 | 32030 | 115.24073910713196 |
| p66 | 40809 | 15.934371948242188 | 38021 | 116.89702320098877 |
| p67 | 29907 | 17.89909291267395 | 29351 | 125.08517384529114 |
| p6 |