LOSER
利用局部搜索方法以及模拟退火算法来求解带容量设备选址问题,对比两种方式的优劣程度,及其讲述使用的技巧。对于CFLP的71个测例来进行测试,写出实验结果的数据。实验中采用多种邻域操作的局部搜索local search策略尝试解决相同规模的CFLP问题,并与相同局部搜索的模拟退火算法进行对比。从算法结果来看,模拟退火算法求得的最优解要比局部搜索得出的最优解要好。通过该实验得出结论,局部搜索容易陷入局部最优解,而模拟退火利用其一定概率接受差解的特性跳出局部最优,从而能逐渐收敛到全局最优。
给定n个工厂和m个顾客,开工厂需要一定的费用,一个工厂有一定的容量限制,每个顾客也有一定的需求,而每个顾客要选取某个工厂也需要一定的分配费用,现在要求找出一个分配方案,把顾客分配给不同的工厂,然后在可以满足所有顾客需求的前提下让所有的花费(开工厂的花费和分配的花费)最小。
这显然是一个NP-hard问题,因为情况数非常的多,而且也没有什么固定的好的策略,所以现在我们就用贪心算法和模拟退火法两种方法来解决这个问题。
局部搜索(Local Search Algorithm)是指能够无穷接近最优解的能力,而全局收敛能力是指找到全局最优解所在大致位置的能力。
局部搜索是解决最优化问题的一种启发式算法。对于某些计算起来非常复杂的最优化问题,比如各种NP完全问题,要找到最优解需要的时间随问题规模呈指数增长,因此诞生了各种启发式算法来退而求其次寻找次优解,是一种近似算法(Approximate algorithms),以时间换精度的思想。局部搜索就是其中的一种方法。
算法的关键问题是:
局部最优问题
步长问题
起始点问题
模拟退火算法(Simulated Annealing Algorithm)是元启发式搜索算法的其中一种。相比起效率低下的随即搜索算法,元启发式搜索算法借助了演化思想和集群智能思想,对解的分布有规律的复杂问题有良好的效果。
模拟退火法(Simulated Annealing)是克服爬山法缺点的有效方法,所谓退火是冶金专家为了达到某些特种晶体结构重复将金属加热或冷却的过程,该过程的控制参数为温度T。模拟退火法的基本思想是,在系统朝着能量减小的趋势这样一个变化过程中,偶尔允许系统跳到能量较高的状态,以避开局部极小点,最终稳定到全局最小点。
算法的关键问题是:
初温与接受温度的设置
局部搜索的方法
退火系数
某一温度下的内循环次数
随机挑选出一单元,并给他一个随机的位移k,求出系统因此而产生的能量变化∆Ek
若∆Ek小于等于0,该位移可采纳,而变化后的系统状态可作为下次变化的起点
若∆Ek大于0,位移后的状态可采纳的概率如下公式所示:式中T为温度,然后从(0,1)区间均匀分布的随机数中挑选一个数R。若R小于Pk,则将变化后的状态作为下次的起点;否则,将变化前的状态作为下次的起点
根据文件中的设备数量,用户数量读取数据,并将其放在vector中,这里我利用了5个vector。
vector<int> facilityResult; // 安排哪些工厂开放(01表示)
vector<int> customerResult; // 安排每个顾客去哪个设备
vector<int> facilitySequence; // 设备开放的顺序
vector<int> customerSequence; // 客户加入的顺序
vector<spare> spareCapacity; // 设备的剩余空间
除此以外,为了搜索操作的需要,要对于剩余容量进行排序,或者对于用户花费进行排序,所以这里需要一些结构体。结构体根据value或cost进行排序,id保持不变,指示唯一的用户或设备。
// 用户对于不同设备的花费
struct customer_cost
{
float value;
int id;
};
// 剩余空间
struct spare
{
int id;
int capacity;
};
// 顾客结构体
struct customer
{
float demand;
std::vector<customer_cost> cost;
};
// 设备结构体
struct facility
{
int capacity;
int cost;
};
一些参数,前四个用于保存文件中的数据,只会在读取文件的时候改变,后面为只读。
后四个参数为模拟退火算法使用的参数,通过调整这些参数能获取更好的解。
int nFacility; // 工厂数量
int nCustomer; // 顾客数量
std::vector<customer> customerVector; // 顾客的信息
std::vector<facility> facilityVector; // 设备的信息
const double coolSpeed = 0.99; // 退火速度
const double acceptT = 0.00001; // 接受温度
const double initialT = 2000; // 初始温度
const int turn = 1000; // 内循环次数
统计一个解的长度
评价函数是寻找最优解的一个关键,CFLP问题可以利用解的总花费代价作为评价该解的优秀程度。
计算花费代价,包括两部分,其一是开启设备的花费,其二是用户选择该设备的花费,两者求和即是该解的总花费
// 统计一个解的总花费
void CalculateTotalCost(){
totalCost = 0;
// 设备花费
for (int i = 0; i < nFacility; ++i)
{
if (facilityResult[i])
{
totalCost += facilityVector[i].cost;
}
}
// 用户花费
for (int i = 0; i < nCustomer; ++i)
{
for (int j = 0; j < nFacility; ++j)
{
if (customerResult[i] == customerVector[i].cost[j].id)
{
totalCost += customerVector[i].cost[j].value;
break;
}
}
}
}
判断一个解的合法性
一个解是否合法需要检查,每一个设备的容量是否小于被选中该设备的用户需求之和,否则认为该解是不合法的,超出设备的容量。
// 检查该解的合法性
bool CheckSolution(){
int facilityCapacity[50];
for(int i = 0; i < nFacility; i++){
facilityCapacity[i] = 0;
}
// 检查容量
for (int i = 0; i < nCustomer; ++i)
{
facilityCapacity[customerResult[i]] += customerVector[i].demand;
}
for (int i = 0; i < nFacility; ++i)
{
if (facilityCapacity[i] > facilityVector[i].capacity)
{
return false;
}
}
return true;
}
随机生成初始解
这里先随机生成设备开放顺序以及顾客加入顺序
然后根据这两个顺序来遍历用户与设备,这里填充的法则是容量First Fit
当该设备可以装下时,继续放这个设备;直到放不下才放到下一个设备。
// 随机生成一个解,根据开放的顺序,顺序填充顾客
void RandomGenerate(){
srand(time(0));
int nChange1 = rand() % nFacility;
int nChange2 = rand() % nCustomer;
for (int i = 0; i < nFacility; ++i)
{
facilitySequence[i] = i;
}
// 随机设备开放顺序
for (int i = 0; i < nChange1; ++i)
{
int pos1 = rand() % nFacility;
int pos2 = rand() % nFacility;
int temp = facilitySequence[pos1];
facilitySequence[pos1] = facilitySequence[pos2];
facilitySequence[pos2] = temp;
}
// 随机顾客加入顺序
for (int i = 0; i < nCustomer; ++i){
customerSequence[i] = i;
}
for (int i = 0; i < nChange2; ++i)
{
int pos1 = rand() % nCustomer;
int pos2 = rand() % nCustomer;
int temp = customerSequence[pos1];
customerSequence[pos1] = customerSequence[pos2];
customerSequence[pos2] = temp;
}
// 按照前面的顺序来填充一个解
int j = 0; // 统计设备的序号
facilityResult[facilitySequence[0]] = 1;
for (int i = 0; i < nCustomer; ++i)
{
// 当当前设备装不下用户的需求,查找一个容量能装下的设备
while(customerVector[customerSequence[i]].demand > spareCapacity[facilitySequence[j]].capacity){
j++;
j = j % 10;
// 新的
if (facilityResult[facilitySequence[j]] == 0)
{
facilityResult[facilitySequence[j]] = 1;
spareCapacity[facilitySequence[j]].capacity = facilityVector[facilitySequence[j]].capacity;
}
}
// 将用户装入设备序号j
spareCapacity[facilitySequence[j]].capacity -= customerVector[i].demand;
customerResult[i] = facilitySequence[j];
}
}
四种局部搜索的策略
这里,我们设置了四种局部搜索的策略来进行找邻域中的解。
(1) 交换两设备的开关情况
随机找到两个设备,要求这两个设备是一开一关的,将两个设备的开关情况交换。
需要的操作:
将a设备的用户提取出来
将a设备的容量还原
将提取出来的用户放入设备b中,当b放不下后,再往后遍历
更改设备b的容量
重新统计该解的花费
count = 40;
while(count > 0){
int i = rand() % nFacility;
int j = rand() % nFacility;
if (facilityResult[i] == 1 && facilityResult[j] == 0)
{
facilityResult[i] = 0;
facilityResult[j] = 1;
// 需要改变的用户
std::vector<int> changeCustomer;
for (int k = 0; k < nCustomer; ++k)
{
if (customerResult[k] == i)
{
changeCustomer.push_back(k);
}
}
// 返还空间给设备
for (int x = 0; x < nFacility; ++x)
{
if(spareCapacity[x].id == i){
spareCapacity[x].capacity = facilityVector[x].capacity;
break;
}
}
for (int k = 0; k < changeCustomer.size(); ++k)
{
// 采用容量best fit,根据容量的大小来安排
sort(spareCapacity.begin(),spareCapacity.end());
for(int x = 0; x < nFacility; x++){
if(customerVector[changeCustomer[k]].demand < spareCapacity[x].capacity){
spareCapacity[x].capacity -= customerVector[changeCustomer[k]].demand;
customerResult[changeCustomer[k]] = spareCapacity[x].id;
facilityResult[spareCapacity[x].id] = 1;
break;
}
}
}
break;
}
count--;
(2) 提取两个设备的所有用户出来,重新分配
随机找到两个设备,要求这两个设备都具有用户的,将提取出来的用户重新分配。
需要的操作:
将a, b设备的用户提取出来
将a, b设备的容量还原
将提取出来的用户按照容量Best Fit重新放置
更改设备的容量
重新统计该解的花费
int i1 = rand() % nFacility;
int i2 = rand() % nFacility;
while(i1 == i2) i2 = rand() % nFacility;
// 随机搜索一些用户,拿出来重新放置
int count = nCustomer;
// 需要改变的用户
std::vector<int> changeCustomer;
for (int k = 0; k < count; ++k)
{
if (customerResult[k] == i1)
{
changeCustomer.push_back(k);
// 返还空间给设备
for (int x = 0; x < nFacility; ++x)
{
if(spareCapacity[x].id == i1){
spareCapacity[x].capacity += customerVector[k].demand;
break;
}
}
}
}
for (int k = 0; k < count; ++k)
{
if (customerResult[k] == i2)
{
changeCustomer.push_back(k);
// 返还空间给设备
for (int x = 0; x < nFacility; ++x)
{
if(spareCapacity[x].id == i2){
spareCapacity[x].capacity += customerVector[k].demand;
break;
}
}
}
}
facilityResult[i1] = 0;
facilityResult[i2] = 0;
// 从下一个开始装填
int j = (i1+i2) % nFacility;
for (int k = 0; k < changeCustomer.size(); ++k)
{
// 采用best fit
sort(spareCapacity.begin(),spareCapacity.end());
for(int x = 0; x < nFacility; x++){
if(customerVector[changeCustomer[k]].demand < spareCapacity[x].capacity){
spareCapacity[x].capacity -= customerVector[changeCustomer[k]].demand;
customerResult[changeCustomer[k]] = spareCapacity[x].id;
//cout << "spareCapacity[x].id " << spareCapacity[x].id << endl;
facilityResult[spareCapacity[x].id] = 1;
break;
}
}
}
break;
(3) 提取所有设备的最后两个用户出来,重新分配
要求用户数大于等于2的设备,提取出来的用户重新分配。
需要的操作:
将所有设备的用户提取出来
将所有设备的容量填补
将提取出来的用户按照价格Best Fit重新放置
更改设备的容量
重新统计该解的花费
// 将每个设备的后面4个客户拿出来,重新分配
// 统计每个设备的顾客
vector<vector<int> > v(nFacility);
// 需要改变的用户
std::vector<int> changeCustomer;
for (int i = 0; i < nCustomer; ++i)
{
v[customerResult[i]].push_back(i);
}
for (int i = 0; i < v.size(); ++i)
{
if(v[i].size() > 2){
int temp1 = v[i].back();
v[i].pop_back();
// 返还空间给设备
for (int x = 0; x < nFacility; ++x)
{
if(spareCapacity[x].id == i){
spareCapacity[x].capacity += customerVector[temp1].demand;
break;
}
}
int temp2 = v[i].back();
v[i].pop_back();
// 返还空间给设备
for (int x = 0; x < nFacility; ++x)
{
if(spareCapacity[x].id == i){
spareCapacity[x].capacity += customerVector[temp2].demand;
break;
}
}
changeCustomer.push_back(temp1);
changeCustomer.push_back(temp2);
}
}
for (int k = 0; k < changeCustomer.size(); ++k)
{
// 采用价格best fit
std::vector<customer_cost> vc;
vc = customerVector[k].cost;
sort(vc.begin(),vc.end());
for(int x = 0; x < nFacility; x++){
int fa_subscript = 0;
// 找下标
for(int q = 0; q < nFacility; q++){
if(vc[x].id == spareCapacity[q].id){
fa_subscript = q;
break;
}
}
// 减容量
if(spareCapacity[fa_subscript].capacity > customerVector[k].demand){
spareCapacity[fa_subscript].capacity -= customerVector[k].demand;
customerResult[changeCustomer[k]] = vc[x].id;
facilityResult[vc[x].id] = 1;
break;
}
}
}
break;
设置循环次数,局部搜索
局部搜索仅接受好解,当遇到差解时,保留旧解。循环多次直到解不再改变,这时说明局部搜索进入了局部最优解,而且无法跳出来。
Solution SA() {
//之后用到很多随机数,设定种子
srand(time(0));
int i = 0;
Solution s;
// 查看初解
Solution best = s;
//设定一个固定的循环次数
for (i = 0; i < turn; i++) {
Solution next = s;
getNewSolution(next);//获取邻域解
if (Accept(next, s, t)) {
s = next;//接受新解为当前最优
if(next < best){//记录全局最优
best = next;
}
}
}
return best;
}
其中Accept函数仅介绍比当前解要好的解,而getNewSolution就是通过之前定义的领域查找操作获取新解。
在局部搜索的基础上,不改变局部搜索的策略,仅仅添加上温度的控制,以此来让算法在一定的概率下接受差解。
参数设定
定义初始温度为20000,接受温度为0.0001,退火系数为0.99(即$T_{i+1} = 0.99 * T_i$)
每个温度内循环为1000次
将局部搜索加入到每一个温度下的循环中,并且找到一个新解时候,判断其与旧解的差值。利用这个插值以及温度来设定概率,用这个概率判断是否接受这个差解,当然好解仍然是全部接受。
执行完一次温度的1000次内循环,进行降温操作。
为了提高精度,在降温到一定程度,我还使用了二次升温的操作,能继续跳出局部最优。
Solution best = s;
int addT = 0;
while (t > acceptT) {
//设定一个固定的循环次数
for (i = 0; i < turn; i++) {
Solution next = s;
getNewSolution(next);//获取邻域解
//cout << next.CheckSolution() << " " << i << endl;
if (Accept(next, s, t)) {
s = next;//接受新解为当前最优
if(next < best){//记录全局最优
best = next;
}
}
}
t *= r;//降温
//升温
if(addT < 3 && t < 20){
t *= 20;
addT++;
}
}
实验环境: Windows10
Result(SA) | Time(s)(SA) | Result(LS) | Time(s)(LS) | |
---|---|---|---|---|
p1 | 15054 | 3.599 | 18436 | 1.320 |
p2 | 15730 | 4.235 | 19856 | 1.217 |
p3 | 18795 | 4.689 | 24013 | 1.356 |
p4 | 22208 | 4.761 | 29014 | 1.157 |
p5 | 18231 | 4.16 | 20165 | 1.323 |
p6 | 19997 | 4.404 | 34103 | 1.029 |
p7 | 21155 | 4.075 | 25612 | 1.594 |
p8 | 20788 | 4.259 | 24658 | 1.358 |
p9 | 18834 | 4.019 | 26874 | 1.364 |
p10 | 15319 | 4.707 | 21684 | 1.277 |
p11 | 19323 | 3.407 | 25674 | 1.471 |
p12 | 19152 | 3.876 | 26453 | 1.351 |
p13 | 16688 | 4.881 | 24351 | 0.987 |
p14 | 18437 | 5.584 | 21368 | 1.375 |
p15 | 15328 | 4.741 | 19874 | 1.312 |
p16 | 16377 | 5.135 | 18795 | 1.489 |
p17 | 14281 | 5.038 | 19634 | 1.687 |
p18 | 18375 | 4.976 | 19874 | 1.547 |
p19 | 18699 | 4.794 | 20120 | 1.214 |
p20 | 17147 | 5.239 | 24310 | 1.200 |
p21 | 15569 | 5.272 | 21587 | 1.219 |
p22 | 14727 | 5.402 | 19574 | 1.354 |
p23 | 18181 | 4.823 | 20148 | 1.114 |
p24 | 16539 | 5.032 | 24687 | 1.365 |
p25 | 54957 | 7.822 | 84123 | 1.987 |
p26 | 66687 | 8.808 | 84712 | 2.359 |
p27 | 62994 | 7.691 | 78998 | 2.126 |
p28 | 83972 | 8.832 | 98564 | 1.784 |
p29 | 71849 | 8.144 | 95678 | 1.859 |
p30 | 42215 | 7.794 | 65987 | 1.879 |
p31 | 61148 | 7.893 | 66874 | 2.003 |
p32 | 55206 | 7.586 | 65239 | 2.410 |
p33 | 47642 | 8.265 | 75462 | 1.983 |
p34 | 88701 | 8.368 | 95611 | 2.153 |
p35 | 52371 | 7.268 | 66782 | 1.984 |
p36 | 55141 | 7.025 | 79621 | 2.264 |
p37 | 69995 | 8.635 | 98561 | 2.361 |
p38 | 57220 | 8.546 | 77856 | 2.223 |
p39 | 62602 | 8.587 | 85463 | 1.925 |
p40 | 92181 | 8.286 | 98711 | 2.314 |
p41 | 12846 | 4.085 | 15687 | 1.256 |
p42 | 12260 | 6.124 | 15655 | 1.597 |
p43 | 11387 | 6.221 | 13200 | 1.548 |
p4 | 14150 | 3.99 | 14895 | 1.268 |
p45 | 12393 | 5.669 | 13547 | 1.344 |
p46 | 11054 | 5.971 | 15314 | 1.562 |
p47 | 12982 | 3.371 | 14355 | 1.435 |
p48 | 11069 | 6.018 | 16423 | 1.998 |
p49 | 10414 | 7.082 | 13251 | 2.013 |
p50 | 15268 | 3.790 | 15468 | 1.246 |
p51 | 18369 | 7.856 | 20114 | 1.985 |
p52 | 18539 | 4.471 | 26478 | 2.014 |
p53 | 15378 | 7.315 | 20431 | 2.687 |
p54 | 21268 | 3.911 | 25143 | 1.654 |
p55 | 17370 | 6.965 | 19562 | 1.664 |
p56 | 44390 | 12.102 | 56899 | 3.012 |
p57 | 54296 | 11.502 | 98543 | 2.987 |
p58 | 48486 | 12.382 | 56214 | 2.887 |
p59 | 59698 | 10.41 | 85356 | 3.124 |
p60 | 48793 | 11.909 | 63512 | 3.210 |
p61 | 46821 | 12.259 | 59823 | 3.214 |
p62 | 52271 | 11.585 | 64985 | 2.963 |
p63 | 43909 | 11.592 | 47351 | 2.799 |
p64 | 44307 | 12.712 | 53684 | 3.569 |
p65 | 48465 | 11.768 | 58164 | 2.689 |
p66 | 47556 | 12.533 | 60012 | 3.045 |
p67 | 47062 | 12.556 | 56841 | 3.512 |
p68 | 46467 | 12.680 | 46295 | 3.246 |
p69 | 43157 | 12.321 | 48652 | 2.869 |
p70 | 53930 | 11.740 | 85612 | 2.694 |
p71 | 49072 | 11.782 | 62548 | 3.001 |
使用单纯的局部搜索方法来解决CFLP问题,非常容易陷入局部最优,这样就会导致解的准确度很差,但是我们可以先利用局部搜索的解来作为模拟退火的初解,这样会加快模拟退火的准确度,所以局部搜索在对于一些局部峰值较少的问题还是很好的,因为其收敛速度很快,其解的性质也不会太差。
对于这个CFLP问题的领域搜索:
对于容量
对于花费
经过我设置搜索策略的观察,往往对于花费的BEST FIT操作能找到更好的解,但是如果仅仅使用花费的贪心搜索,这跟贪心算法没什么区别,且达不到最优解,因为会收到容量的限制。故我在搜索策略中加入对于容量的FIRST FIT以及BEST FIT,然而这与花费的结果是有矛盾冲突的,所以也需要再三权衡。
关于局部搜索与模拟退火的对比:局部搜索收敛速度较快,但易陷入局部最优,对于初值不太敏感;而模拟退火算法收敛速度较慢,算法的准确度对于参数的设置以及初解的好坏敏感。对于不同的测试例子也有不同的反馈结果,需要进行调参。