[LibreOJ 2249][NOI2014]购票
在紧张刺激的等待之后终于肝掉了这道题……
本题的DP方程长成这样(其中\(a\)指\(v\)的某个满足距离限制的祖先,\(d_v\)指\(v\)到根的路径长):
\[f_v = min(f_a + p_v(d_v - d_a) + q_v)\]
化简之后发现:
\[f_v = q_v + p_v d_v + min(f_a - p_v d_a)\]
利用\(min\)中那一块很容易发现是截距式……但是问题在于,我们的转移来源是树上的连续一段祖先,怎样维护他们的凸包?
答案很狂暴啊……用树链剖分套上向量集那题的线段树套凸包,然后完了……
(注意一点细节:本题因为数据范围过大,故可能存在两个向量叉乘爆long long,所以在求凸包时如果直接用叉积判断是否需要删点会炸掉,建议用斜率判断)
代码:
#include <cstdio> #include <cstring> #include <cstdlib> #include <cctype> #include <algorithm> #include <utility> #include <vector> #include <queue> #include <deque> #include <cmath> #include <set> #include <climits> using ll = long long; using T = ll; using R = long double; const R eps = 1e-8; int sign(R x) { if(fabsl(x) < eps) { return 0; } else { if(x > (R(0.00))) { return 1; } else { return -1; } } } struct Point { T x, y; Point(T qx = 0LL, T qy = 0LL) { x = qx; y = qy; } }; using Vector = Point; Vector operator +(const Vector &a, const Vector &b) { return Vector(a.x + b.x, a.y + b.y); } Vector operator -(const Point &a, const Point &b) { return Vector(a.x - b.x, a.y - b.y); } Vector operator *(const Vector &a, T lam) { return Vector(a.x * lam, a.y * lam); } Vector operator *(T lam, const Vector &a) { return Vector(a.x * lam, a.y * lam); } inline T dot(const Vector &a, const Vector &b) { return (a.x * b.x + a.y * b.y); } inline T times(const Vector &a, const Vector &b) { return (a.x * b.y - a.y * b.x); } inline bool cmp(const Point &a, const Point &b) { if(a.x == b.x) { return a.y < b.y; } else { return a.x < b.x; } } inline R slope(const Vector &a) { R dx = a.x, dy = a.y; return (dy / dx); } inline void andrew(Point *P, int L, std::vector<Point> &bot, std::vector<Point> &top) { std::sort(P + 1, P + 1 + L, cmp); for(int i = 1; i <= L; i ++) { if(i != 1 && (P[i].x == P[i - 1].x)) continue; while(bot.size() > 1 && sign(slope(P[i] - bot.back()) - slope(bot.back() - bot[bot.size() - 2])) <= 0) { bot.pop_back(); } bot.push_back(P[i]); } for(int i = L; i >= 1; i --) { if(i != L && (P[i].x == P[i + 1].x)) continue; while(top.size() > 1 && sign(slope(P[i] - top.back()) - slope(top.back() - top[top.size() - 2])) <= 0) { top.pop_back(); } top.push_back(P[i]); } std::reverse(top.begin(), top.end()); } const int maxn = 200005; const int maxno = maxn << 2; const int N = 200000; bool zen[maxno]; std::vector<Point> bot[maxno], top[maxno]; Point P[maxn]; inline void maintain(int o, int L, int R) { static Point tmp[maxn]; const int lc = o << 1, rc = o << 1 | 1; const bool used = zen[o]; zen[o] = (zen[lc] && zen[rc]); if(zen[o] != used) { std::copy(P + L, P + R + 1, tmp + 1); int len = R - L + 1; andrew(tmp, len, bot[o], top[o]); } } void modify(int o, int L, int R, const int &p, const Point &v) { if(L == R) { zen[o] = true; P[L] = v; bot[o].push_back(v); top[o].push_back(v); } else { const int M = (L + R) / 2; if(p <= M) { modify(o << 1, L, M, p, v); } else { modify(o << 1 | 1, M + 1, R, p, v); } maintain(o, L, R); } } inline T calc_ans(T k, const Point &v) { return v.y - k * v.x; } inline T search(const std::vector<Point> &vec, const T &k) { int l = 0, r = vec.size() - 1; while(r - l > 2) { int lm = (l * 2 + r) / 3, rm = (2 * r + l) / 3; if((calc_ans(k, vec[lm]) > calc_ans(k, vec[rm]))) { l = lm; } else { r = rm; } } T ans = LLONG_MAX; for(int i = l; i <= r; i ++) { ans = std::min(ans, calc_ans(k, vec[i])); } return ans; } T query(int o, int L, int R, const int &ql, const int &qr, const T &k) { if(ql <= L && R <= qr) { return search(bot[o], k); } else { int M = (L + R) / 2; T ans = LLONG_MAX; if(ql <= M) { ans = std::min(ans, query(o << 1, L, M, ql, qr, k)); } if(qr > M) { ans = std::min(ans, query(o << 1 | 1, M + 1, R, ql, qr, k)); } return ans; } } int first[maxn]; int next[maxn << 1], to[maxn << 1]; ll dist[maxn << 1]; void add_edge(int u, int v, ll d) { static int cnt = 0; cnt ++; next[cnt] = first[u]; first[u] = cnt; to[cnt] = v; dist[cnt] = d; } int fa[maxn], dep[maxn], hson[maxn]; ll d[maxn]; int siz[maxn]; int bs[maxn][18]; void dfs_1(int x, int f = -1, int depth = 1) { fa[x] = bs[x][0] = f; dep[x] = depth; siz[x] = 1; int max_siz = 0; for(int i = first[x]; i; i = next[i]) { int v = to[i]; if(v != f) { d[v] = d[x] + dist[i]; dfs_1(v, x, depth + 1); siz[x] += siz[v]; if(siz[v] > max_siz) { hson[x] = v; max_siz = siz[v]; } } } } int dfn[maxn], tid[maxn], up[maxn]; void dfs_2(int x, int a = 1, int f = 0) { static int cnt = 0; dfn[x] = ++ cnt; tid[cnt] = x; up[x] = a; if(hson[x]) { dfs_2(hson[x], a, x); } else { return; } for(int i = first[x]; i; i = next[i]) { int v = to[i]; if(v != f && v != hson[x]) { dfs_2(v, v, x); } } } int k_anc(int x, ll k) { int yx = x; for(int j = 17; j >= 0; j --) { int a = bs[x][j]; if(a != -1 && d[yx] - d[a] <= k) { x = a; } } #ifdef LOCAL printf("%d's %lld-th anc : %d\n", yx, k, x); #endif return x; } int n; ll get_up(int x, int anc, ll k) { ll ans = LLONG_MAX; while(up[x] != up[anc]) { ans = std::min(ans, query(1, 1, n, dfn[up[x]], dfn[x], k)); x = fa[up[x]]; } return std::min(ans, query(1, 1, n, dfn[anc], dfn[x], k)); } ll p[maxn], q[maxn], l[maxn]; ll f[maxn]; void dp(int x) { #ifdef LOCAL printf("processing %d...\n", x); printf("d : %lld\n", d[x]); #endif if(x != 1) { #ifdef LOCAL printf("b : %lld\n", get_up(fa[x], k_anc(x, l[x]), p[x])); #endif f[x] = get_up(fa[x], k_anc(x, l[x]), p[x]) + d[x] * p[x] + q[x]; } else { f[x] = 0; } #ifdef LOCAL printf("ans : %lld\n", f[x]); #endif modify(1, 1, n, dfn[x], Point(d[x], f[x])); for(int i = first[x]; i; i = next[i]) { int v = to[i]; dp(v); } } int main() { int t; scanf("%d%d", &n, &t); for(int i = 2; i <= n; i ++) { int father; T s; scanf("%d%lld%lld%lld%lld", &father, &s, &p[i], &q[i], &l[i]); add_edge(father, i, s); } memset(bs, -1, sizeof(bs)); dfs_1(1); dfs_2(1); for(int j = 1; (1 << j) < n; j ++) { for(int i = 1; i <= n; i ++) { int a = bs[i][j - 1]; if(a != -1) { bs[i][j] = bs[a][j - 1]; } } } dp(1); for(int i = 2; i <= n; i ++) { printf("%lld\n", f[i]); } return 0; }
[LibreOJ 2353][NOI2007]货币兑换
emmm做了一下这道神题……(我可能是少有的用动态凸包苟的?)
首先DP方程长这样:
\[f_i = max(f_{i - 1}, f_j\cdot\frac{A_iR_j+B_i}{A_jR_j+B_j})\]
然后这个方程炒鸡复杂……首先\(f_{i - 1}\)不要管了,然后设\(a_i = \frac{f_i}{A_iR_i + B_i}\)。在xjb推了一番之后我们终于得到了截距式……
\[-a_j R_j \frac{A_i}{B_i} + \frac{f_i}{B_i} = a_j\]
但是这玩意太毒瘤了……斜率不可能单调的,这还好,在凸壳上二分/三分一下即可。但问题在于,横坐标也不单调……
这个时候就需要动态维护凸包了(其实是我不会CDQ),我直接把我向量集那题的二进制分组线段树搬了过来……(逃
代码:
#include <cstdio> #include <cstring> #include <cstdlib> #include <cctype> #include <algorithm> #include <utility> #include <vector> #include <cmath> #include <climits> #include <deque> #include <cassert> using R = double; const R eps = 1e-8; int sign(R x) { if(fabs(x) < eps) { return 0; } else { if(x > 0.00) { return 1; } else { return -1; } } } struct Point { R x, y; Point(R qx = 0, R qy = 0) { x = qx; y = qy; } }; using Vector = Point; Vector operator +(const Vector &a, const Vector &b) { return Vector(a.x + b.x, a.y + b.y); } Vector operator -(const Point &a, const Point &b) { return Vector(b.x - a.x, b.y - a.y); } Vector operator *(const Vector &a, R lam) { return Vector(a.x * lam, a.y * lam); } Vector operator *(R lam, const Vector &a) { return Vector(a.x * lam, a.y * lam); } inline R dot(const Vector &a, const Vector &b) { return (a.x * b.x + a.y * b.y); } inline R times(const Vector &a, const Vector &b) { return (a.x * b.y - a.y * b.x); } inline bool cmp(const Point &a, const Point &b) { if(sign(a.x - b.x) == 0) { return a.y < b.y; } else { return a.x < b.x; } } inline void andrew(Point *P, int L, std::vector<Point> &bot, std::vector<Point> &top) { std::sort(P + 1, P + 1 + L, cmp); for(int i = 1; i <= L; i ++) { if(i != 1 && sign(P[i].x - P[i - 1].x) == 0) continue; while(bot.size() > 1 && sign(times(P[i] - bot.back(), bot.back() - bot[bot.size() - 2])) >= 0) { bot.pop_back(); } bot.push_back(P[i]); } for(int i = L; i >= 1; i --) { if(i != L && sign(P[i].x - P[i + 1].x) == 0) continue; while(top.size() > 1 && sign(times(P[i] - top.back(), top.back() - top[top.size() - 2])) >= 0) { top.pop_back(); } top.push_back(P[i]); } std::reverse(top.begin(), top.end()); } const int maxn = 1000005; const int N = 1000000; const int maxno = maxn << 2; bool zen[maxno]; std::vector<Point> bot[maxno], top[maxno]; Point P[maxn]; inline void maintain(int o, int L, int R) { static Point tmp[maxn]; const int lc = o << 1, rc = o << 1 | 1; const bool used = zen[o]; zen[o] = (zen[lc] && zen[rc]); if(zen[o] != used) { std::copy(P + L, P + R + 1, tmp + 1); int len = R - L + 1; andrew(tmp, len, bot[o], top[o]); } } void modify(int o, int L, int R, const int &p, const Point &v) { if(L == R) { zen[o] = true; P[L] = v; bot[o].push_back(v); top[o].push_back(v); } else { const int M = (L + R) / 2; if(p <= M) { modify(o << 1, L, M, p, v); } else { modify(o << 1 | 1, M + 1, R, p, v); } maintain(o, L, R); } } inline R calc_ans(R k, const Point &v) { return v.y - k * v.x; } inline R search(const std::vector<Point> &vec, const R &k) { int l = 0, r = vec.size() - 1; while(r - l > 2) { int lm = (l * 2 + r) / 3, rm = (2 * r + l) / 3; if(sign(calc_ans(k, vec[lm]) - calc_ans(k, vec[rm])) == 1) { r = rm; } else { l = lm; } } R ans = INT_MIN; for(int i = l; i <= r; i ++) { ans = std::max(ans, calc_ans(k, vec[i])); } return ans; } R query(int o, int L, int R, const int &ql, const int &qr, const double &k) { if(ql <= L && R <= qr) { return search(top[o], k); } else { int M = (L + R) / 2; double ans = INT_MIN; if(ql <= M) { ans = std::max(ans, query(o << 1, L, M, ql, qr, k)); } if(qr > M) { ans = std::max(ans, query(o << 1 | 1, M + 1, R, ql, qr, k)); } return ans; } } int n, s; R A[maxn], B[maxn], Rate[maxn]; R f[maxn]; R dp() { static double a[maxn]; f[0] = s; f[1] = s; a[1] = f[1] / (A[1] * Rate[1] + B[1]); modify(1, 1, n, 1, Point(a[1] * Rate[1], a[1])); for(int i = 2; i <= n; i ++) { f[i] = query(1, 1, n, 1, i - 1, -A[i] / B[i]) * B[i]; f[i] = std::max(f[i], f[i - 1]); a[i] = f[i] / (A[i] * Rate[i] + B[i]); if(i < n) modify(1, 1, n, i, Point(a[i] * Rate[i], a[i])); } return f[n]; } int main() { scanf("%d%d", &n, &s); for(int i = 1; i <= n; i ++) { scanf("%lf%lf%lf", &A[i], &B[i], &Rate[i]); } printf("%.3lf\n", dp()); return 0; }
[LibreOJ 2035][SDOI2016]征途
又做了一道适合我这种沙茶做的题qwq
考虑方差,它满足这个式子:
\[Var(X) = E[X^2] - (E[X])^2\]
对于这个题展开,发现后半部分是一个常量(\((\frac{s_n}{m})^2\))。最小化的就是前半部分,前半部分相当于要求你求一个序列划分,使得各部分平方和的平均值最小。这个值乘上\(m^2\)之后就变成了各部分平方和乘上\(m\)。
然后就很简单了……DP方程化简之后是这样的:
\[f_{t, i} = ms_i^2 + max(f_{t - 1, j} + ms_j^2 - 2ms_i s_j)\]
求截距式形式,得:
\[2ms_i s_j + d_i = f_j + ms_j^2\]
然后用类似于上上题的方法搞就行。还有我想想啥时候补一下斜率优化的tutorial……
代码:
#include <cstdio> #include <cstring> #include <cstdlib> #include <cctype> #include <algorithm> #include <utility> #include <deque> #include <cmath> #include <climits> typedef long long ll; typedef ll T; struct Point { T x, y; Point(T qx = 0LL, T qy = 0LL) { x = qx; y = qy; } }; typedef Point Vector; Vector operator +(const Vector &a, const Vector &b) { return Vector(a.x + b.x, a.y + b.y); } Vector operator -(const Point &a, const Point &b) { return Vector(a.x - b.x, a.y - b.y); } Vector operator *(const Vector &a, T lam) { return Vector(a.x * lam, a.y * lam); } Vector operator *(T lam, const Vector &a) { return Vector(a.x * lam, a.y * lam); } inline T dot(const Vector &a, const Vector &b) { return (a.x * b.x + a.y * b.y); } inline T times(const Vector &a, const Vector &b) { return (a.x * b.y - a.y * b.x); } const int maxn = 3005; T f[maxn][maxn]; T S[maxn]; int n; T m; void dp() { for(int t = 1; t <= m; t ++) { std::deque<Point> Q; Q.push_back(Point(0LL, 0LL)); for(int i = 1; i <= n; i ++) { ll k = 2 * m * S[i]; Vector st(1LL, k); while(Q.size() > 1 && times(Q[1] - Q[0], st) > 0LL) { Q.pop_front(); } f[t][i] = m * S[i] * S[i]; f[t][i] += Q.front().y - Q.front().x * k; #ifdef LOCAL printf("f[%d][%d] : %lld\n", t, i, f[t][i]); #endif if(t == 1) continue; Vector ins(S[i], f[t - 1][i] + m * S[i] * S[i]); while(Q.size() > 1 && times(ins - Q.back(), Q.back() - Q[Q.size() - 2]) > 0LL) { Q.pop_back(); } Q.push_back(ins); } } } int main() { scanf("%d%lld", &n, &m); for(int i = 1; i <= n; i ++) { scanf("%lld", &S[i]); } for(int i = 1; i <= n; i ++) { S[i] += S[i - 1]; } dp(); printf("%lld\n", f[m][n] - S[n] * S[n]); return 0; }
[BZOJ 3156]防御准备
又做了一个简单的斜率优化题TAT
首先,设\(f_i\)表示最后一个放塔的点事\(i\)时的最优解,那么将原方程化简得:
\[f_i = a_i + \frac{i^2 - i}{2} + max(-ij + \frac{j^2 + j}{2} + f_j)\]
然后求直线形式,得到:
\[ij + d_i = f_j + \frac{j^2 + j}{2}\]
用类似于玩具装箱(上一篇题解)的方式搞一搞即可。
代码:
/************************************************************** Problem: 3156 User: danihao123 Language: C++ Result: Accepted Time:2496 ms Memory:16480 kb ****************************************************************/ #include <cstdio> #include <cstring> #include <cstdlib> #include <cctype> #include <algorithm> #include <utility> #include <deque> #include <cmath> #include <climits> typedef long long ll; typedef ll T; struct Point { T x, y; Point(T qx = 0LL, T qy = 0LL) { x = qx; y = qy; } }; typedef Point Vector; Vector operator +(const Vector &a, const Vector &b) { return Vector(a.x + b.x, a.y + b.y); } Vector operator -(const Point &a, const Point &b) { return Vector(a.x - b.x, a.y - b.y); } Vector operator *(const Vector &a, T lam) { return Vector(a.x * lam, a.y * lam); } Vector operator *(T lam, const Vector &a) { return Vector(a.x * lam, a.y * lam); } inline T dot(const Vector &a, const Vector &b) { return (a.x * b.x + a.y * b.y); } inline T times(const Vector &a, const Vector &b) { return (a.x * b.y - a.y * b.x); } const int maxn = 1000005; T f[maxn], a[maxn]; int n; void dp() { f[1] = a[1]; std::deque<Point> Q; Q.push_back(Point(1LL, f[1] + 1LL)); for(T i = 2; i <= n; i ++) { T k = i; Vector st(1LL, k); while(Q.size() > 1 && times(Q[1] - Q[0], st) > 0LL) { Q.pop_front(); } f[i] = a[i] + ((i * i - i) / 2LL); f[i] += Q.front().y - Q.front().x * i; Point ins(i, f[i] + ((i * i + i) / 2LL)); while(Q.size() > 1 && times(ins - Q.back(), Q.back() - Q[Q.size() - 2]) > 0LL) { Q.pop_back(); } Q.push_back(ins); } } int main() { scanf("%d", &n); for(int i = n; i >= 1; i --) { scanf("%lld", &a[i]); } dp(); ll ans = f[n]; #ifdef LOCAL printf("f[n] : %lld\n", ans); #endif for(T i = 1; i < n; i ++) { #ifdef LOCAL printf("f[%d] : %lld\n", i, f[i]); #endif ans = std::min(ans, f[i] + (n - i + 1LL) * (n - i) / 2LL); } printf("%lld\n", ans); return 0; }
[BZOJ 1010][HNOI2008]玩具装箱toy
很久之前是学过并写过斜率优化的……但是很快就忘了。现在感觉自己理解了,感觉是真的懂了……抽空写篇文章解释一下吧……
先单独说这一个题。将DP方程完全展开,并且设\(P_i = S_i + i\),\(c = L + 1\),可得:
\[f_i = c^2 + P_i^2 - 2P_i c + max(P_j^2 + 2P_j c + f_j - 2P_i P_j)\]
然后\(c^2 + P_i^2 - 2P_i c\)这部分是常数项不需要管了,我们就想想max里面那些(姑且设之为\(d_i\))咋整好了。
设\(d_i = P_j^2 + 2P_j c + f_j - 2P_i P_j\),稍作移项,得:
\[2P_i P_j + d_i = P_j^2 + 2P_j c + f_j\]
于是乎,\(d_i\)可以看做斜率为\(2P_i\)的直线过点\((P_j, P_j^2 + 2P_j c + f_j)\)得到的截距。而那些点我们之前都知道了,问题就变成了已知斜率,求过某点集中的点的最大截距。
想象一个固定斜率的直线从下往上扫,那么碰到的第一个点就是最优解。首先这个点一定在下凸壳上,其次下凸壳上这点两侧的线段的斜率肯定一个比\(2P_i\)大另一个比它小。并且最好的一点是这个斜率还是单调的,那么分界点一定是单调递增的。
代码:
/************************************************************** Problem: 1010 User: danihao123 Language: C++ Result: Accepted Time:132 ms Memory:2416 kb ****************************************************************/ #include <cstdio> #include <cstring> #include <cstdlib> #include <cctype> #include <algorithm> #include <utility> #include <deque> #include <cmath> typedef long long ll; typedef ll T; struct Point { T x, y; Point(T qx = 0LL, T qy = 0LL) { x = qx; y = qy; } }; typedef Point Vector; Vector operator +(const Vector &a, const Vector &b) { return Vector(a.x + b.x, a.y + b.y); } Vector operator -(const Point &a, const Point &b) { return Vector(a.x - b.x, a.y - b.y); } Vector operator *(const Vector &a, T lam) { return Vector(a.x * lam, a.y * lam); } Vector operator *(T lam, const Vector &a) { return Vector(a.x * lam, a.y * lam); } inline T dot(const Vector &a, const Vector &b) { return (a.x * b.x + a.y * b.y); } inline T times(const Vector &a, const Vector &b) { return (a.x * b.y - a.y * b.x); } const int maxn = 50005; T C[maxn], S[maxn], P[maxn]; T f[maxn]; int n; ll c; void process() { for(int i = 1; i <= n; i ++) { S[i] = S[i - 1] + C[i]; P[i] = S[i] + (ll(i)); } } void dp() { std::deque<Point> Q; Q.push_back(Point(0LL, 0LL)); for(int i = 1; i <= n; i ++) { ll k = 2 * P[i]; Vector st(1, k); while(Q.size() > 1 && times(Q[1] - Q[0], st) > 0LL) { Q.pop_front(); } f[i] = c * c + P[i] * P[i] - 2LL * P[i] * c; f[i] += Q.front().y - k * Q.front().x; #ifdef LOCAL printf("f[%d] : %lld\n", i, f[i]); #endif Vector ins(P[i], f[i] + P[i] * P[i] + 2LL * P[i] * c); while(Q.size() > 1 && times(ins - Q.back(), Q.back() - Q[Q.size() - 2]) > 0LL) { #ifdef LOCAL printf("Deleting (%lld, %lld)...\n", Q.back().x, Q.back().y); #endif Q.pop_back(); } Q.push_back(ins); #ifdef LOCAL printf("Inserting (%lld, %lld)...\n", ins.x, ins.y); #endif } } int main() { scanf("%d%lld", &n, &c); c ++; for(int i = 1; i <= n; i ++) { scanf("%lld", &C[i]); } process(); dp(); printf("%lld\n", f[n]); return 0; }