rainboy的解析

§ 解析

和诗人小G很像啊

先写出状态转移方程,设f(i)f(i)表示前i个玩具可以达到的最小值

f(i)=min{f(j)+val(j+1,i)} f(i) = \min\{f(j) + val(j+1,i)\}

其中val(j,i)val(j,i)表示[j,i][j,i]区间内的所有的玩具放到一个容器内的时候得到的值

val(j,i)=(ij+k=jiCkL)2 val(j,i) = ( i-j + \sum_{k=j}^{i}C_k - L)^2

轻松写出一下朴素的n2n^2的算法如下

点击
cpp
copy
        
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
//Author by [Rainboy](https://github.com/rainboylvx) //date: 2024-04-18 15:43:43 #include <bits/stdc++.h> using namespace std; const int maxn = 1e6+5; typedef long long ll; int n,L; int a[maxn]; //前缀和 ll pre_sum[maxn]; ll f[maxn]; ll p[maxn]; ll val(int j,int i) { ll len = i-(j+1)+pre_sum[i] - pre_sum[j]; len -=L; return len * len; } void init() { std::cin >> n >> L; for(int i = 1;i <= n ;++i ) // i: 1->n { std::cin >> a[i]; pre_sum[i] = pre_sum[i-1]; pre_sum[i] += a[i]; } } void debug_f() { cout << "f[i]: "; for(int i = 1;i <= n ;++i ) // i: 1->n { cout << f[i] << " "; } cout << endl; cout << "p[i]: "; for(int i = 1;i <= n ;++i ) // i: 1->n { cout << p[i] << " "; } cout << endl; } int main (int argc, char *argv[]) { init(); f[0] = 0; for(int i = 1;i <= n ;++i ) // 前i个元素 { f[i] = 0x7f7f7f7f7f7f7f7f; for(int j = 0 ;j < i;j++) { ll t = f[j] + val(j,i); if( f[i] >= t) { f[i] = t; p[i] = j; } } } cout << f[n] << endl; #ifdef DEBUG debug_f(); #endif return 0; }

显然要优化的:

§ 决策单调性

这个状态转移方程有决策单调性.下面证明

我们先把题目转化成一个数学问题.

设在一个数轴上有这几个点a<b<x<ya < b < x < y

-------a--------b----------x----------y---------

于是得到g(a,x)g(a,x)

g(a,x)=f(a)+val(a+1,x)=f(a)+(k=a+1xCk+x(a+1)L)2=f(a)+(sum(x)sum(a)+xa1L)2=f(a)+(sum(x)+x(sum(a)+a)(L+1))2=f(a)+(SxSaL1)2(1) \begin{aligned} g(a,x) &= f(a) + val(a+1,x) \\ &= f(a) + (\sum_{k=a+1}^{x}C_k + x-(a+1) - L)^2 \\ &= f(a) + (sum(x) - sum(a) + x-a-1 - L)^2 \\ &= f(a) + (sum(x)+x - (sum(a) +a ) -(L+1))^2 \\ &= f(a) + (S_x - S_a -L_1)^2 \end{aligned} \tag 1

同理,可以得到g(b,x)g(b,x)

g(b,x)=f(b)+(SxSbL1)2(2) g(b,x) = f(b) + (S_x - S_b -L_1)^2 \tag 2

那么先假定g(b,x)<g(a,x)g(b,x) < g(a,x)成立,那么问题变成证明

g(b,x)<g(a,x)g(b,y)<g(a,y) g(b,x) < g(a,x) \Rightarrow g(b,y) < g(a,y)

上式恒成立.

同样得到g(a,y),g(b,y)g(a,y),g(b,y)

g(a,y)=f(a)+(SySaL1)2(3) g(a,y) = f(a) + (S_y - S_a -L_1)^2 \tag 3
g(b,y)=f(b)+(SySbL1)2(4) g(b,y) = f(b) + (S_y - S_b -L_1)^2 \tag 4

观察发现(3),(4)(3),(4)式对于(1),(2)(1),(2)式来说,只有Sy,SyS_y,S_y不同.根据题意,显然Sy>Sx>0S_y > S_x > 0,那么就设Sy=Sx+v,v>0S_y = S_x+v, v> 0,那么得到新的式子

g(a,y)=f(a)+(Sx+vSaL1)2g(b,y)=f(b)+(Sx+vSbL1)2 \begin{align} g(a,y) = f(a) + (S_x+v - S_a -L_1)^2 \tag 5 \\ g(b,y) = f(b) + (S_x+v - S_b -L_1)^2 \tag 6 \\ \end{align}

根据数学直觉,若想知道(5),(6)(5),(6)式的在小关系,那么只需要知道对于(1),(2)(1),(2)式之间的增量之间的大小关系.

g(a,y)=f(a)+(Sx+vSaL1)2=f(a)+(v+(SxSaL1))2=f(a)+v2+2v(SxSaL1)+(SxSaL1)2 \begin{aligned} g(a,y) &= f(a) + (S_x+v - S_a -L_1)^2 \\ &= f(a) + (v+(S_x - S_a -L_1))^2 \\ &= f(a) + v^2+2v(S_x - S_a -L_1) +(S_x - S_a -L_1)^2 \\ \end{aligned}

那么增加的量就是v2+2v(SxSaL1)v^2+2v(S_x - S_a -L_1).同理(6)(6)式的增量就是v2+2v(SxSbL1)v^2+2v(S_x - S_b -L_1).比较两者之间的大小

v2+2v(SxSaL1)(v2+2v(SxSbL1))=2v(SxSaL1Sx+SbL1)=2v(SbSa)>0 \begin{aligned} &v^2+2v(S_x - S_a -L_1) - ( v^2+2v(S_x - S_b -L_1)) \\ &= 2v( S_x - S_a -L_1 -S_x + S_b -L_1) \\ &= 2v(S_b - S_a) > 0 \end{aligned}

因为Sb>SaS_b > S_a一定成立.

所以题目的状态转移方程是满足决策单调性

然后根据[ Rbook: 决策单调性]可以写出来下的代码

点击
cpp
copy
        
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
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
// 使用 决策单调性 , 实现的 //Author by [Rainboy](https://github.com/rainboylvx) //date: 2024-04-19 14:49:04 #include <bits/stdc++.h> using namespace std; const int maxn = 1e6+5; typedef long long ll; int n,L; int a[maxn]; //前缀和 ll pre_sum[maxn]; ll f[maxn]; //存入栈中的点 struct node { int p; //最优决策点 int l,r; // [l,r]的决策点都是p }; // 求[j+1,i]放在一起的值 ll val(int j,int i) { ll len = i-(j+1)+pre_sum[i] - pre_sum[j]; len -=L; return len * len; } void init() { std::cin >> n >> L; for(int i = 1;i <= n ;++i ) // i: 1->n { std::cin >> a[i]; pre_sum[i] = pre_sum[i-1]; pre_sum[i] += a[i]; } } //栈模板 template<typename T = int,int siz = maxn> struct mystack{ T sta[siz+5]; int head = 0; void clear() { head = 0;} void push(T a) { sta[head++] = a;} void pop(){head--;} T top() { return sta[head-1];} bool empty() { return head == 0;} int size() { return head;} }; mystack<node> sta; //二分的代码 int mid(int l,int r) { return (l+r) >> 1; //这是最快的写法 } //bs_find = binary search find template <typename F> int bs_find(int l,int r,F check) { while( l < r) { int m = mid(l,r); if( check(m)) //成立 r = m; else //不成立,抛弃左半边 l = m+1; } return l ; } //输出栈内的内容 void debug_sta() { cout << "p[i]: \n"; // for(int i = 0 ;i< sta.size();i++) { // node t = sta.sta[i]; // for(int j = t.l; j<=t.r ;j++) // cout << t.p << " "; // } for(int i = 0 ;i< sta.size();i++) { cout << "p: " << sta.sta[i].p << " "; cout << "l: " << sta.sta[i].l << " "; cout << "r: " << sta.sta[i].r << "\n"; } cout << endl; } void debug_f() { cout << "f[i]: "; for(int i = 1;i <= n ;++i ) // i: 1->n { cout << f[i] << " "; } cout << endl; } int main (int argc, char *argv[]) { init(); f[0] = 0; sta.push({0,1,n}); for(int i = 1;i <= n ;++i ) // 前i个元素 { #ifdef DEBUG cout << "before calc i= " << i << endl; debug_f(); debug_sta(); cout << endl; #endif f[i] = 0x7f7f7f7f7f7f7f7f; // 先找到i的决策点p,然后计算出f[i]的值 //找到包含i的区间 int pos = bs_find(0,sta.size(),[i](int m){ //判断这个m所在的区间在i的的左边 //判断方法: 当前node的r大于等于i, return sta.sta[m].r >= i; }); pos = sta.sta[pos].p; // 决策点 f[i] = f[pos] + val(pos,i); // 单调栈 // 从后面向前找,如果点i比当前节点的p(决策点)优秀,就出栈 while( !sta.empty()) { node t = sta.top(); //是否更优秀 ll org = f[t.p] + val(t.p,t.l); ll now = f[i] +val(i,t.l); if( now <= org) { //更加优秀 sta.pop(); } else { //不是更优秀 //在区间 [t.l,t.r+1] 找到开头的那个点 sta.pop(); int pos = bs_find(t.l,t.r+1,[p=t.p,i](int m){ ll org = f[p] + val(p,m); ll now = f[i] +val(i,m); return now <= org; }); if( pos > n){ // sta.push(t); } else { sta.push({t.p,t.l,pos-1}); sta.push({i,pos,n}); } break; } } #ifdef DEBUG cout << "after calc i= " << i << endl; debug_f(); debug_sta(); cout << endl; #endif } cout << f[n] << endl; return 0; }

§ slope

cpp
copy
        
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
115
/* author: Rainboy email: rainboylvx@qq.com time: 2023年 06月 06日 星期二 21:20:13 CST */ #include <bits/stdc++.h> using namespace std; typedef long long ll; const int maxn = 1e6+5,maxe = 1e6+5; //点与边的数量 int n,L; /* 定义全局变量 */ using db = long double; //前缀和 long long sum[maxn]; ll dp[maxn]; //斜率优化 struct slope_optimize { ll k(int i) { return -2 * B(i);} ll B(int i) {return sum[i] +i -L -1;} // A(j) 就是xj ll A(int j) {return sum[j] +j;} //dj 就是 yj ll d(int j) {return dp[j] + A(j)*A(j);} //两个点的斜率 double slope(int i,int j) { return (d(j) - d(i)) *1.0 / ( A(j) - A(i)); } ll calc_dp(int i,int j) { return k(i)*A(j) +d(j)+B(i)*B(i); } } slope_opt; //队列 struct myque { int q[maxn]; int head,tail; int size() const { return tail - head; } bool empty() const { return size() == 0; } void push(int a) { q[tail++] = a; } void del(int i) { while( size() > 1 && slope_opt.slope(q[tail-2],i) < slope_opt.slope(q[tail-2], q[tail-1])) tail--; } //单调队列不停的删除头部 void update(double v) { while( size() > 1 && slope_opt.slope(q[head],q[head+1]) < v) ++head; } int front() const { return q[head]; } int back() const { return q[tail-1]; } } que; void init() { cin >> n >> L; for(int i=1;i<=n;++i){ cin >> sum[i]; sum[i] += sum[i-1]; } } int main(int argc,char * argv[]){ cin.sync_with_stdio(false); cin.tie(0); cout.tie(0); init(); //读取数据 que.push(0); //尝试每一个点 for(int i=1;i<=n;i++) { que.update(- slope_opt.k(i) ); dp[i] = slope_opt.calc_dp(i, que.front()); que.del(i); que.push(i); } cout << dp[n]; return 0; }