很久没有来这个Blog写东西了。之前因为决定去腾讯，所以竞赛的事情放下了很长一段时间，现在重新拾起。
zhk大神要准备为WF拼一次，我和Talker与他并肩奋斗，于是zhk给我们推荐了一些题目去做，于是，就有了这个做到泪奔的Cross the Wall。
闲话少说，该题目为DP中的战斗机——斜率优化，在最简单的斜率优化基础上添加了一维K，总的来说思想还是没太大变化。之前上网找了很多资料，也去各位大神的Blog取过经，花了一周的时间才断断续续地搞清楚。下面就来简单总结一下自己这段时间学习斜率优化的心得。
1.DP方程
斜率优化的典型DP方程为dp(i) = min{dp(j) + f(i,j)}，对于这道题，f(i,j) = h(j+1)*w(i)。加上挖洞的个数k，完整的DP方程为dp(i, k) = min{dp(j, k-1) + h(j+1) * w(i)};由于dp(i, k)只与dp(j, k-1)有关，则可以使用滚动数组，在此将DP方程简化为dp(i) = min{dp(j) + h(j+1) * w(i)};
若dp(i) == dp(j) + h(j+1) * w(i) （即当j时，能使dp(i, k)最小），则对于j'(0 <= j' < i && j' != j)，有dp(j) + h(j+1) * w(i) < dp(j') + h(j'+1) * w(i);移项可得dp(j) - dp(j') < w(i) * (h(j'+1) - h(j+1));把h作为x，dp作为y，则w为斜率。至此，斜率优化的雏形初现。
2.决策点维护
由dp(j) - dp(j') < w(i) * (h(j'+1) - h(j+1))分析以及题目需要求最小值可知，决策点必然是所有可能作为决策点的点集中最靠近下方的，换句话说就是要求这个点集的左链即凸包的左半部分。至于这个决策点按照什么顺序产生，如何维护，都是大同小异的。要注意的是，h和w是有单调性的。在这道题中，我是保证h单减，w单增。
用堆栈p来维护左链，这样，若dp(i) = dp(j) + h(j+1) * w(i), (j < i)，则在决策点的堆栈中，p(idx)作为决策点，有p(idx) = (h(j+1), dp(j));根据斜率优化的表达式以及凸包的性质可知，在堆顶元素之前，线段的斜率都是小于w(i)的。由于w是递增的，那么在计算dp(i+1)时，决策点不会退化到当前决策点p(idx)之前，这样就可以更快地找到决策点，降低复杂度。
3.实例
在分析算法的时候，自己写了个简单的数据来帮助理解。
4 100
99 100
100 99
200 50
201 49
则h[] = {100, 99, 50, 49}, w[] = {99, 100, 200, 201};
具体的就不细讲了。写这篇东西一来是自己总结，二来是希望能够给需要的人一些帮助。这道题dp值可能很大，要用int64，在输出的时候要用%I64d，否则会出错。就是这里我用%lld，让我查了一个下午...
太久不写代码始终是有些生疏了，很多问题需要注意，加油吧~
HDU 3669 Cross the Wall
2 #include <string.h>
3 #include <algorithm>
4
5 #define DEBUG 0
6
7 using namespace std;
8
9 const int MaxN = 50010, MaxK = 110;
10 typedef long long i64;
11 int N, K;
12
13 struct Rect {
14 int W, H;
15 }rect[MaxN+10];
16
17 struct Point {
18 i64 x, y;
19 Point() {}
20 Point(const i64 X, const i64 Y): x(X), y(Y) {}
21 Point operator - (const Point P) const {
22 return Point(x-P.x, y-P.y);
23 }
24 Point operator + (const Point P) const {
25 return Point(x+P.x, y+P.y);
26 }
27 };
28
29 i64 mul(Point a, Point b) {
30 return a.x*b.y - a.y*b.x;
31 }
32
33 struct Line {
34 int sp, i;
35 Point point[MaxN+10];
36 Line(){}
37 void init() {
38 sp = -1;
39 i = 0;
40 }
41 bool insert(Point P) {
42 while(sp > 0 && mul(P-point[sp], P-point[sp-1]) < 0) //维护凸包
43 sp--;
44 point[++sp] = P;
45 return true;
46 }
47 i64 optimize(int k) {
48 while(i < sp && point[i].y + point[i].x*k >= point[i+1].y + point[i+1].x*k) //斜率优化
49 i++;
50 return point[i].y + point[i].x*k;
51 }
52 }line;
53
54 int x[MaxN], y[MaxN], yn;
55
56 int ID(int num) {
57 int s = -1, e = yn, mid = (s + e) >> 1;
58 while(s + 1 < e) {
59 if(y[mid] == num)
60 return mid;
61 if(y[mid] > num)
62 e = mid;
63 else
64 s = mid;
65 mid = (s + e) >> 1;
66 }
67 }
68
69 void optimize() {
70 for(int i = 1; i <= N; ++i)
71 y[i-1] = rect[i].W;
72 sort(y, y + N);
73 yn = unique(y, y + N) - y;
74 memset(x, 0, sizeof(x));
75 for(int i = 1; i <= N; ++i) {
76 int idx = ID(rect[i].W);
77 x[idx] = max(x[idx], rect[i].H);
78 }
79 N = yn;
80 for(int i = 0; i < yn; ++i) {
81 rect[i+1].H = x[i];
82 rect[i+1].W = y[i];
83 }
84 int sp = 1;
85 for(int i = 1; i <= N; ++i) {
86 if(sp == 0) {
87 rect[++sp] = rect[i];
88 continue;
89 }
90 if(rect[i].W == rect[sp].W)
91 continue;
92 if(rect[i].H >= rect[sp].H) {
93 sp--; i--;
94 }
95 else if(rect[i].H < rect[sp].H)
96 rect[++sp] = rect[i];
97 }
98 N = sp;
99 #if DEBUG
100 for(int i = 1; i <= N; ++i)
101 printf("%d %d\n", rect[i].W, rect[i].H);
102 #endif
103 }
104
105 i64 dp[2][MaxN+10];
106
107 void work() {
108 optimize(); //预处理所有矩形
109 int cur = 1, pre = 0;
110 for(int i = 1; i <= N; ++i) { //初始化dp[pre]为挖一个洞的花费
111 dp[pre][i] = 1LL*rect[i].W * rect[1].H;
112 }
113 int cnt = 0;
114 for(int k = 2; k <= K; ++k) { //从挖两个洞开始计算
115 line.init(); //左链初始化
116 line.insert(Point(rect[1].H, 0)); //预先插入决策点(h[1], dp[x][0])。
117 #if DEBUG
118 for(int i = 1; i <= N; ++i)
119 printf("%15lld", dp[pre][i]);
120 printf("\n");
121 #endif
122 for(int i = 1; i <= N; ++i) {
123 dp[cur][i] = line.optimize(rect[i].W); //计算斜率优化出来的dp值
124 if(i < N) {
125 line.insert(Point(rect[i+1].H, dp[pre][i])); //将之前计算的dp值作为可能的决策点插入
126 }
127 }
128 cur = (cur+1)%2;
129 pre = 1 - cur;
130 if(dp[cur][N] == dp[pre][N])
131 cnt++;
132 if(cnt > 5)
133 break;
134 }
135 printf("%I64d\n", dp[pre][N]);
136 }
137
138 int main() {
139 #if DEBUG
140 freopen("test.in", "r", stdin);
141 freopen("test.out", "w", stdout);
142 #endif
143 while(scanf("%d%d", &N, &K) == 2) {
144 for(int i = 1; i <= N; ++i)
145 scanf("%d%d", &rect[i].W, &rect[i].H);
146 work();
147 }
148 return 0;
149 }
150
