模擬退火演算法解旅行商(TSP)問題
阿新 • • 發佈:2019-02-08
該帖子的程式碼主要轉自模擬退火演算法1.
該文對模擬退火演算法作了較好的分析,不過該文中舉例的TSP的程式碼有一些問題,我對此作了修正,並在文中最後做出解釋。
程式碼如下:
#include <iostream>
#include <string.h>
#include <stdlib.h>
#include <algorithm>
#include <stdio.h>
#include <time.h>
#include <math.h>
#define N 30 //城市數量
#define T 3000 //初始溫度
#define EPS 1e-8 //終止溫度
#define DELTA 0.98 //溫度衰減率
#define LIMIT 10000 //概率選擇上限
#define OLOOP 1000 //外迴圈次數
#define ILOOP 15000 //內迴圈次數
using namespace std;
//定義路線結構體
struct Path
{
int citys[N];
double len;
};
//定義城市點座標
struct Point
{
double x, y;
};
Path path; //記錄最優路徑
Point p[N]; //每個城市的座標
double w[N][N]; //兩兩城市之間路徑長度
int nCase; //測試次數
double dist(Point A, Point B)
{
return sqrt((A.x - B.x) * (A.x - B.x) + (A.y - B.y) * (A.y - B.y));
}
void GetDist(Point p[], int n)
{
for (int i = 0; i < n; i++)
for (int j = i + 1; j < n; j++)
w[i][j] = w[j][i] = dist(p[i], p[j]);
}
void Input(Point p[], int &n)
{
scanf("%d", &n);
for (int i = 0; i < n; i++)
scanf("%lf %lf", &p[i].x, &p[i].y);
}
void Init(int n)
{
nCase = 0;
path.len = 0;
for (int i = 0; i < n; i++)
{
path.citys[i] = i;
if (i != n - 1)
{
printf("%d--->", i);
path.len += w[i][i + 1];
}
else
printf("%d\n", i);
}
printf("\nInit path length is : %.3lf\n", path.len);
}
void Print(Path t, int n)
{
printf("Path is : ");
for (int i = 0; i < n; i++)
{
if (i != n - 1)
printf("%d-->", t.citys[i]);
else
printf("%d\n", t.citys[i]);
}
printf("\nThe path length is : %.3lf\n", t.len);
}
// 隨機交換2個城市的位置
Path GetNext(Path p, int n)
{
// Modify by DD: 確保x!=y
int x = 0, y = 0;
while (x == y)
{
x = (int)(n * (rand() / (RAND_MAX + 1.0)));
y = (int)(n * (rand() / (RAND_MAX + 1.0)));
}
Path ans = p;
swap(ans.citys[x], ans.citys[y]);
ans.len = 0;
for (int i = 0; i < n - 1; i++)
ans.len += w[ans.citys[i]][ans.citys[i + 1]];
//cout << "nCase = " << nCase << endl;
//Print(ans, n);
nCase++;
return ans;
}
// 模擬退火
void SA(int n)
{
double t = T;
srand(time(NULL));
Path curPath = path;
Path newPath = path;
int P_L = 0; // 連續找到更差結果的次數
int P_F = 0; // while迴圈次數
while (1) //外迴圈,主要更新引數t,模擬退火過程
{
for (int i = 0; i < ILOOP; i++) //內迴圈,尋找在一定溫度下的最優值
{
newPath = GetNext(curPath, n);
double dE = newPath.len - curPath.len;
if (dE < 0) //如果找到更優值,直接更新
{
curPath = newPath;
P_L = 0;
P_F = 0;
}
else
{
double rd = rand() / (RAND_MAX + 1.0);
// Modify by DD: dE取負數才有可能接受更差解,否則e>1
double e = exp(-dE / t);
if ( e > rd && e < 1) //如果找到比當前更差的解,以一定概率接受該解,並且這個概率會越來越小
curPath = newPath;
P_L++;
}
if (P_L > LIMIT)
{
P_F++;
break;
}
}
// Modify by DD: 記錄全域性最優解
if (curPath.len < path.len)
path = curPath;
if (P_F > OLOOP || t < EPS)
break;
t *= DELTA;
}
}
int main()
{
freopen("TSP.txt", "r", stdin);
int n;
Input(p, n);
GetDist(p, n);
Init(n);
SA(n);
Print(path, n);
printf("Total test times is : %d\n", nCase);
return 0;
}
資料檔案如下:
27
41 94
37 84
53 67
25 62
7 64
2 99
68 58
71 44
54 62
83 69
64 60
18 54
22 60
83 46
91 38
25 38
24 42
58 69
71 71
74 78
87 76
18 40
13 40
82 7
62 32
58 35
45 21
修改過的地方,程式碼中寫了註釋。
其中幾個重要的、需要解釋的地方有如下幾個:
記錄全域性最優解
記錄最優解時,原文使用的是
if (curPath.len < newPath.len)
path = curPath;
這樣顯然記錄的是當前溫度中的最優解。更合理的做法是記錄全域性出現過的最優解,即
if(curPath.len < path .len)
path= curPath;
(…這裡不加文字,markdown有序列表編號會出錯…)
2. 連續搜尋到最差結果的處理
程式碼中的 P_L用於標示連續搜尋到比當前更差結果的次數。程式碼中,如果連續發生了LIMIT次(程式碼中是10000次),就意味著當前解空間已經找不到更好的結果了。
多次執行程式碼測試會發現,原文程式碼更換城市次數卻僅僅在1w次左右,最終路徑在450附近,這就說明了原文程式碼會陷入區域性最優解。即使每次把P_L復位為0,執行次數是上去了,但是仍然會陷入區域性最優解。
3. 必須以一定概率接受更差結果
能夠跳出全域性最優,正是模擬退火演算法最大的優勢所在。
原來,原文中計算的是exp(dE / t),注意dE是為正的,因此exp(dE / t)恆大於1. 因此更差的結果永遠不會被選中。而由於交換操作中每次只交換curPath的一對城市,因此陷入了curPath的區域性最優解。
修改為exp(-dE / t)後,curPath有一定概率被替換。測試結果發現,更換城市次數在700w次左右,最終路徑在370附近。顯然獲得了全域性最優解。
其他不錯的資料,如隨機模擬的基本思想和常用取樣方法(sampling)2.