hdu 5067 Harry And Dig Machine tsp问题

缘起

学习经典的 TSP 问题的N种姿势. 即售货商问题(又叫做中国邮政问题,或者货郎担问题).

TSP是著名的NP问题,并且所有其他NP问题可以(在多项式时间)内约化(reduce)到TSP问题上来

hdu 5067 Harry And Dig Machine

TSP是著名的NP完全问题. 已经被证明无法给出多项式时间的算法.但是竞赛中依旧可能出现数据范围较小的题目.

分析

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
有n*m棋盘上给定不超过10个城市,从左上角城市出发,每个城市只能走一次,再回到左上角城市,求最小代价。

【输入】
多样例. 对于每个样例,首先是两个正整数n,m (1<=n,m<=50)
然后n行,每行m个正整数a[i][j],只要它>0就表示此处是城市.

【输出】
移动一格用1时间,对每个样例输出最短完成时间.

【样例输入】
3 3
0 0 0
0 100 0
0 0 0
2 2
1 1
1 1

【样例输出】
4
4

【限制】
1s

tsp问题在小数据规模下有如下几种算法

  1. 状压dp
  2. dfs(贪心剪枝)
  3. 分支限界
  4. 记忆化搜索(或称为备忘录算法)

约莫四种姿势. 下面逐一使用这几种姿势艹它.

状压DP

其实首先读入棋盘,然后得到这不超过10个城市的坐标. 然后for循环建立这些城市的邻接矩阵(城市索引最小为0). 然后令

dp[S][u] 是已经经过的城市集合S(注意,从城市0出发, 因为最后还要回到0,所以不认为一开始0属于S, S一开始是空集),目前在u(所以u属于S)这种状态下完成剩下任务的所需要的最小代价. 则

显然有

1
2
3
4
dp[V][0] = 0; // 已经经历了所有顶点, 并且当前已经在0处了(即完成了TSP), 则完成剩下任务(其实已经不存
在剩下任务了)最小代价显然就是0
dp[S][v] = min{dp[S+{u}][u]+cos[u][v]}, 其中v的取值范围是在S中,而u的取值范围是u不在S中,这表示
此时售货商在v顶点已经经过了S之后, 根据下一步(到u去)所有选择取最小.

因为n<=10, 所以可以使用二进制状态压缩S.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
#//include "stdafx.h"

#include <stdio.h>
#include <algorithm>
using namespace std;
#include <string.h>
//#define LOCAL
#define ABS(x,y) ((x)>(y)?((x)-(y)):((y)-(x)))

int n,m, g[15][15], b[55][55],tot,dp[1<<15][15]; // g是各个城市之间的邻接矩阵, tot是城市的个数

struct City
{
int x,y;
}c[15];

int h(const City &a, const City &b) // 2个城市之间的距离
{
return ABS(a.x, b.x) + ABS(a.y, b.y);
}

int main()
{
#ifdef LOCAL
freopen("d:\\data.in", "r", stdin);
//freopen("d:\\my.out", "w", stdout);
#endif
while(~scanf("%d%d", &n,&m))
{
tot = 1; // 左上角城市首先算进去
for (int i = 0; i<n; i++)
{
for (int j = 0; j<m;j++)
{
scanf("%d", b[i]+j);
if (b[i][j]&&(i||j)) // 如果是非左上角的城市(因为左上角肯定要算进去的)
{
c[tot].x = i, c[tot].y = j;
tot++;
}
}
}
if (tot==1) // 如果只有左上角一个城市,则显然代价为0,一开始这里没注意,wa了一发
{
puts("0");
continue;
}
for (int i = 0; i<tot; i++)
{
for (int j = i+1; j<tot; j++)
{
g[i][j] = g[j][i] = h(c[i], c[j]);
}
} // 初始化邻接矩阵
memset(dp, 0x3f, sizeof(dp));
dp[(1<<tot)-1][0] = 0;
for (int s = (1<<tot)-2; s; s--)
{
for (int v = 0; v<tot; v++)
{
if (s>>v&1)
{
for (int u = 0; u<tot; u++)
{
if (!(s>>u&1))
{
dp[s][v] = min(dp[s][v], dp[s|(1<<u)][u]+g[v][u]);
}
}
}
}
}
int ans = 0x3f3f3f3f;
for (int i = 1; i<tot; i++)
{
ans = min(ans, dp[1<<i][i]+g[0][i]);
}
printf("%d\n", ans);
}
return 0;
}

ac情况

Status Accepted
Time 15ms
Memory 3148kB
Length 1403
Lang G++
Submitted 2019-09-05 06:28:31
Shared
RemoteRunId 30508668

记忆化搜索

首先是记忆化搜索(打表dfs), 这个就不多说什么了. 只是反向dp而已.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
//#include "stdafx.h"

#include <stdio.h>
#include <algorithm>
using namespace std;
#include <string.h>
//#define LOCAL
#define ABS(x,y) ((x)>(y)?((x)-(y)):((y)-(x)))

int n,m, g[15][15], b[55][55],tot,dp[1<<15][15]; // g是各个城市之间的邻接矩阵, tot是城市的个数,dp是打好的表

struct City
{
int x,y;
}c[15];

int h(const City &a, const City &b) // 2个城市之间的距离
{
return ABS(a.x, b.x) + ABS(a.y, b.y);
}

int rec(int s, int v)
{
if (dp[s][v]!=0x3f3f3f3f)
{
return dp[s][v];
}
for (int u = 0; u<tot; u++)
{
if (!(s>>u&1))
{
dp[s][v] = min(dp[s][v], rec(s|(1<<u),u)+g[v][u]);
}
}
return dp[s][v];
}

int main()
{
#ifdef LOCAL
freopen("d:\\data.in", "r", stdin);
//freopen("d:\\my.out", "w", stdout);
#endif
while(~scanf("%d%d", &n,&m))
{
tot = 1; // 左上角城市首先算进去
for (int i = 0; i<n; i++)
{
for (int j = 0; j<m;j++)
{
scanf("%d", b[i]+j);
if (b[i][j]&&(i||j)) // 如果是非左上角的城市(因为左上角肯定要算进去的)
{
c[tot].x = i, c[tot].y = j;
tot++;
}
}
}
if (tot==1) // 如果只有左上角一个城市,则显然代价为0
{
puts("0");
continue;
}
for (int i = 0; i<tot; i++)
{
for (int j = i+1; j<tot; j++)
{
g[i][j] = g[j][i] = h(c[i], c[j]);
}
} // 初始化邻接矩阵
memset(dp, 0x3f, sizeof(dp));
dp[(1<<tot)-1][0] = 0; // 初始化表
int ans = 0x3f3f3f3f;
for (int i = 1; i<tot; i++)
{
ans = min(ans, rec(1<<i,i)+g[0][i]);
}
printf("%d\n", ans);
}
return 0;
}

但是就像我曾经说过的,记忆化搜索其实很慢(经常被T).

Status Time Limit Exceeded
Length 1376
Lang G++
Submitted 2019-09-05 06:41:57
Shared
RemoteRunId 30508670

但是,其实上面的状压是有优化空间的——因为有些状态他喵的就是inf啊~ 所以只需要打上标记flag就可以了——该标记表明他喵的它就是inf,不需要继续深搜计算了! 这样就可以节省很多时间.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
//#include "stdafx.h"

#include <stdio.h>
#include <algorithm>
using namespace std;
#include <string.h>
//#define LOCAL
#define ABS(x,y) ((x)>(y)?((x)-(y)):((y)-(x)))

int n,m, g[15][15], b[55][55],tot,dp[1<<15][15]; // g是各个城市之间的邻接矩阵, tot是城市的个数,dp是打好的表
bool flag[1<<15][15]; // flag[i][j]=true表示dp[i][j]计算过了

struct City
{
int x,y;
}c[15];

int h(const City &a, const City &b) // 2个城市之间的距离
{
return ABS(a.x, b.x) + ABS(a.y, b.y);
}

int rec(int s, int v)
{
if (flag[s][v])
{
return dp[s][v];
}
for (int u = 0; u<tot; u++)
{
if (!(s>>u&1))
{
dp[s][v] = min(dp[s][v], rec(s|(1<<u),u)+g[v][u]);
}
}
flag[s][v] = true; // dp[s][v]已经算过了
return dp[s][v];
}

int main()
{
#ifdef LOCAL
freopen("d:\\data.in", "r", stdin);
//freopen("d:\\my.out", "w", stdout);
#endif
while(~scanf("%d%d", &n,&m))
{
tot = 1; // 左上角城市首先算进去
for (int i = 0; i<n; i++)
{
for (int j = 0; j<m;j++)
{
scanf("%d", b[i]+j);
if (b[i][j]&&(i||j)) // 如果是非左上角的城市(因为左上角肯定要算进去的)
{
c[tot].x = i, c[tot].y = j;
tot++;
}
}
}
if (tot==1) // 如果只有左上角一个城市,则显然代价为0
{
puts("0");
continue;
}
for (int i = 0; i<tot; i++)
{
for (int j = i+1; j<tot; j++)
{
g[i][j] = g[j][i] = h(c[i], c[j]);
}
} // 初始化邻接矩阵
memset(dp, 0x3f, sizeof(dp));
memset(flag, 0, sizeof(flag));
dp[(1<<tot)-1][0] = 0; // 初始化表
int ans = 0x3f3f3f3f;
flag[(1<<tot)-1][0] = true;
for (int i = 1; i<tot; i++)
{
ans = min(ans, rec(1<<i,i)+g[0][i]);
}
printf("%d\n", ans);
}
return 0;
}

哇哈哈哈~ 为记忆化搜索正名了~

Status Accepted
Time 15ms
Memory 3636kB
Length 1520
Lang G++
Submitted 2019-09-05 07:10:45
Shared
RemoteRunId 30508682

dfs+贪心剪枝

因为数据规模较小. 所以dfs亦是可取的. 当然,就算n<=10, n!种状态(最多10层的递归)也是可能让代码表现不佳的(毕竟本题要求<1s). 所以适当使用贪心进行剪枝也是必要的.

tsp如果暴搜的话,复杂度是O((N-1)!)(减1是因为最后的城市是确定了的)

而因为城市不超过10, 所以采用 next_permutation 暴搜试一波~

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
//#include "stdafx.h"

#include <stdio.h>
#include <algorithm>
using namespace std;
#include <string.h>
//#define LOCAL
#define ABS(x,y) ((x)>(y)?((x)-(y)):((y)-(x)))

int n,m, g[15][15], b[55][55],tot; // g是各个城市之间的邻接矩阵, tot是城市的个数

struct City
{
int x,y;
}c[15];

int h(const City &a, const City &b) // 2个城市之间的距离
{
return ABS(a.x, b.x) + ABS(a.y, b.y);
}

int getans(int a[], int tot) // 0->a[0]->...->a[tot-1]->0
{
int ans = 0;
for (int i =0;i<tot-1; i++)
{
ans+=g[a[i]][a[i+1]];
}
ans+=g[0][a[0]];
ans+=g[a[tot-1]][0];
return ans;
}

int main()
{
#ifdef LOCAL
freopen("d:\\data.in", "r", stdin);
//freopen("d:\\my.out", "w", stdout);
#endif
while(~scanf("%d%d", &n,&m))
{
tot = 1; // 左上角城市首先算进去
for (int i = 0; i<n; i++)
{
for (int j = 0; j<m;j++)
{
scanf("%d", b[i]+j);
if (b[i][j]&&(i||j)) // 如果是非左上角的城市(因为左上角肯定要算进去的)
{
c[tot].x = i, c[tot].y = j;
tot++;
}
}
}
if (tot==1) // 如果只有左上角一个城市,则显然代价为0
{
puts("0");
continue;
}
for (int i = 0; i<tot; i++)
{
for (int j = i+1; j<tot; j++)
{
g[i][j] = g[j][i] = h(c[i], c[j]);
}
}
int ans = 0x3f3f3f3f;
int a[15];
tot--; // 不算起点城市(因为暴搜其实是遍历其余城市的全排列. 因为最后要回到的城市是确定了的)
for (int i = 0; i<tot; i++)
{
a[i] = i+1;
}
do
{
ans = min(ans, getans(a,tot));
} while (next_permutation(a,a+tot));
printf("%d\n", ans);
}
return 0;
}

是不是有一种没天理的感觉? 哇哈哈哈~

Status Accepted
Time 327ms
Memory 1228kB
Length 1382
Lang G++
Submitted 2019-09-05 07:05:16
Shared
RemoteRunId 30508679

不说笑了,言归正传, 该剪枝还是要剪枝的. 不然很容易卡你(例如n开到20就嗝屁了). 毕竟327ms依旧很慢

剪枝的策略是:贪心(ps:本题如果用贪心来做是错误的~ 但是贪心可以用来快速生成一个近似解来剪枝)

每次选择距离距离当前城市最近的城市. 然后最后一步跳回0

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
//#include "stdafx.h"

#include <stdio.h>
#include <algorithm>
using namespace std;
#include <string.h>
//#define LOCAL
#define ABS(x,y) ((x)>(y)?((x)-(y)):((y)-(x)))

int n,m, g[15][15], b[55][55],tot,ans; // g是各个城市之间的邻接矩阵, tot是城市的个数,dp是打好的表
bool visited[15];

struct City
{
int x,y;
}c[15];

int h(const City &a, const City &b)
{
return ABS(a.x, b.x) + ABS(a.y, b.y);
}

int greedy()
{
memset(visited, 0, sizeof(visited));
int cur = 0,ans = 0,res = tot-1;
while(res)
{
int _min=0x3f3f3f3f,_minindex = -1;
for (int i = 1; i<tot; i++)
{
if (!visited[i] && _min>g[cur][i])
{
_min = g[cur][i];
_minindex = i;
}
}
ans += g[cur][_minindex]; // 此次选择了_minindex
cur = _minindex;
visited[cur] = true;
res--;
}
return ans+g[cur][0]; // 最后回到0
}

void dfs(int step, int curcity, int cur) // 当前决定的是第step个城市, 而迄今已经产生的消耗是cur,而且现在在curcity这个城市
{
if (step==tot) // 表明已经得到一组方案了
{
if (ans>cur+g[curcity][0]) // 从curcity回到0
{
ans = cur+g[curcity][0];
}
return;
}
if (cur>=ans) // 剪枝
{
return;
}
for (int i = 1; i<tot; i++) // 选择下一个城市
{
if (visited[i])
{
continue;
}
visited[i] = true;
dfs(step+1, i, cur+g[curcity][i]); // 选择城市i
visited[i] = false; // 因为是搜索, 所以要改回来
}
}

int main()
{
#ifdef LOCAL
freopen("d:\\data.in", "r", stdin);
//freopen("d:\\my.out", "w", stdout);
#endif
while(~scanf("%d%d", &n,&m))
{
tot = 1; // 左上角城市首先算进去
for (int i = 0; i<n; i++)
{
for (int j = 0; j<m;j++)
{
scanf("%d", b[i]+j);
if (b[i][j]&&(i||j)) // 如果是非左上角的城市(因为左上角肯定要算进去的)
{
c[tot].x = i, c[tot].y = j;
tot++;
}
}
}
if (tot==1) // 如果只有左上角一个城市,则显然代价为0
{
puts("0");
continue;
}
for (int i = 0; i<tot; i++)
{
for (int j = i+1; j<tot; j++)
{
g[i][j] = g[j][i] = h(c[i], c[j]);
}
}
ans = greedy(); // 先用贪心O(n^2)时间求一波近似解
memset(visited, 0, sizeof(visited));
dfs(1, 0, 0);
printf("%d\n", ans);
}
return 0;
}

ac情况(贪心剪枝优化到了93ms, 暴搜是300+ms)

Status Accepted
Time 93ms
Memory 1232kB
Length 1930
Lang G++
Submitted 2019-09-05 09:00:31
Shared
RemoteRunId 30508829

分支限界

tsp本质上是要找一条哈密顿回路:从图中的任意一点出发,路途中经过图中每一个结点当且仅当一次,则成为哈密顿回路.

采用优先队列式分支限界法. 关于分支限界,【1】将它的本质说透了. 强烈推荐看看 ^_^

本题的分支限界采用的真是代价函数显然就是到达当前路径最后一个节点P产生的代价. 而启发式函数是什么呢? 就是P以及剩下所有节点出发的最短路径之和(这部分其实我们可以预处理,即下面代码中的minOut数组).

很显然,这里的启发式函数是<=真实走完剩余路产生的代价函数. 所以是ok的. 这里采用活动节点是路径. 但是等到堆变空才结束算法. 从而要进行剪枝.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
//#include "stdafx.h"

#include <stdio.h>
#include <algorithm>
#include <queue>
using namespace std;
#include <string.h>
//#define LOCAL
#define ABS(x,y) ((x)>(y)?((x)-(y)):((y)-(x)))

int n,m, g[15][15], b[55][55],tot,ans,minout[15],minoutsum; // g是各个城市之间的邻接矩阵, tot是城市的个数,minout[i]表示顶点i出发的最短边长度, minoutsum是minout数组的和

struct Node
{
int g,h,res,last; // g是已经产生的代价, h是启发式代价, res是剩余未访问的节点, 使用状态压缩(因为tot较小),last是本节点表示的当前处在的城市
Node(int g, int h, int res, int last):g(g),h(h), res(res),last(last){}
bool operator <(const Node &o) const
{
return g+h==o.g+o.h?g<o.g:(g+h>o.g+o.h); // g+h越大, 优先级越低,如果相等的话, g小的(即确确实实攥在手里的少,虚的成分多)则优先级低
}
};

struct City
{
int x,y;
}c[15];

int h(const City &a, const City &b)
{
return ABS(a.x, b.x) + ABS(a.y, b.y);
}

void bb() // branches and bound
{
priority_queue<Node> pq;
pq.push(Node(0, minoutsum, (1<<tot)-2, 0)); // 例如tot=3, 则1<<tot-2得到6,即110, 因为0号城市不在考虑范围内(它一定是最后回到的城市)
while(!pq.empty())
{
Node top = pq.top();
pq.pop();
if (!top.res) // 已经不剩了,说明找到了一条哈密顿回路
{
if (ans>top.g+g[top.last][0])
{
ans = top.g+g[top.last][0];
}
}
for (int i = 1; i<tot; i++)
{
if (top.res>>i&1) // 如果i号城市尚未访问的话
{
Node nxt = Node(top.g+g[top.last][i], top.h-minout[i], top.res&~(1<<i), i); // 下一个准备开拓的节点
if (nxt.g+nxt.h>=ans) // 如果它根本不可能成为最优节点的话(因为它作弊了都>=当前ans)
{
continue; // 分支限界中的剪枝
}
pq.push(nxt); // 有希望就进入堆中
}
}
}
}

int main()
{
#ifdef LOCAL
freopen("d:\\data.in", "r", stdin);
//freopen("d:\\my.out", "w", stdout);
#endif
while(~scanf("%d%d", &n,&m))
{
tot = 1; // 左上角城市首先算进去
minoutsum = 0;
ans = 0x3f3f3f3f;
for (int i = 0; i<n; i++)
{
for (int j = 0; j<m;j++)
{
scanf("%d", b[i]+j);
if (b[i][j]&&(i||j)) // 如果是非左上角的城市(因为左上角肯定要算进去的)
{
c[tot].x = i, c[tot].y = j;
tot++;
}
}
}
if (tot==1) // 如果只有左上角一个城市,则显然代价为0
{
puts("0");
continue;
}
for (int i = 0; i<tot; i++)
{
for (int j = i+1; j<tot; j++)
{
g[i][j] = g[j][i] = h(c[i], c[j]);
}
}
for (int i = 1; i<tot; i++)
{
minout[i] = 0x3f3f3f3f;
for (int j = 1; j<tot; j++)
{
if (minout[i]>g[i][j])
{
minout[i] = g[i][j];
}
}
minoutsum+=minout[i];
} // 得到minout数组
bb();
printf("%d\n", ans);
}
return 0;
}

ac情况 (水过水过, 看来 分支限界剪枝还是不够强力呀~)

Status Accepted
Time 499ms
Memory 10464kB
Length 2271
Lang G++
Submitted 2019-09-05 14:54:34
Shared
RemoteRunId 30510777

参考

【1】https://yfsyfs.github.io/2019/09/05/%E5%88%86%E6%94%AF%E9%99%90%E7%95%8C-bfs-%E6%90%9C%E7%B4%A2-%E7%9A%84%E4%B8%80%E8%88%AC%E6%80%A7%E8%AE%BA%E8%BF%B0/