[BZOJ 3158]千钧一发
好有意思的题呢qwq
首先观察到一点，我们可以把所有装置按照\(a_i\)的奇偶性进行分类（也可以说是黑白染色？）。容易发现任意偶数的最大公约数都至少是2，所以任意偶数间不会互相冲突。然后任意两个奇数的平方和一定不是完全平方数，可以这么证一下（感谢金爷！）：
那两个奇数可以分别设成\(2x + 1\)和\(2y + 1\)，然后他们的平方和就是\(4(x^2 + y^2 + x + y) + 2\)。然后思考一点，就是一个奇数的平方一定是奇数，所以说那个平方和是一个偶数。但是，如果一个完全平方数是偶数，那么它的算术平方根一定是偶数。然而，一个偶数的平方一定是4的倍数。但我们求出的平方和膜4得2，所以不行。
所以说冲突只存在于奇偶性不同的数中。然后用无限大的边来表示冲突关系，最小割高一波即可。
代码（常数卡卡卡……只能用pb_ds的蛤希表了）：
/**************************************************************
Problem: 3158
User: danihao123
Language: C++
Result: Accepted
Time:9676 ms
Memory:130460 kb
****************************************************************/
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cctype>
#include <algorithm>
#include <utility>
#include <queue>
#include <set>
#include <vector>
#include <ext/pb_ds/assoc_container.hpp>
#include <ext/pb_ds/hash_policy.hpp>
typedef long long ll;
ll gcd(ll a, ll b) {
if(!b) {
return a;
} else {
return gcd(b, a % b);
}
}
__gnu_pbds::gp_hash_table<ll, bool> S2;
void process_S2() {
for(ll i = 1LL; i * i <= 2000000000000LL; i ++) {
S2[i * i] = true;
}
}
const int maxno = 1050;
const int maxm = 2000500;
int first[maxno];
int next[maxm], to[maxm], cap[maxm], flow[maxm];
int gcnt = 0;
void add_edge(int u, int v, int f) {
gcnt ++;
next[gcnt] = first[u];
first[u] = gcnt;
to[gcnt] = v;
cap[gcnt] = f;
flow[gcnt] = 0;
}
int rev(int i) {
if(1 & i) {
return i + 1;
} else {
return i - 1;
}
}
void ins_edge(int u, int v, int f) {
add_edge(u, v, f); add_edge(v, u, 0);
}
int n;
int s, t;
int d[maxno];
bool bfs() {
static bool vis[maxno];
memset(vis, 0, sizeof(vis));
d[s] = 1; vis[s] = true;
std::queue<int> Q; Q.push(s);
while(!Q.empty()) {
int u = Q.front(); Q.pop();
for(int i = first[u]; i; i = next[i]) {
int v = to[i];
if(cap[i] > flow[i] && !vis[v]) {
d[v] = d[u] + 1; vis[v] = true;
Q.push(v);
}
}
}
return vis[t];
}
int cur[maxno];
int dfs(int u, int a) {
if(a == 0 || u == t) return a;
int fl = 0;
for(int &i = cur[u]; i; i = next[i]) {
int v = to[i]; int f;
if(d[v] == d[u] + 1 && (f = dfs(v, std::min(cap[i] - flow[i], a))) > 0) {
fl += f; a -= f;
flow[i] += f; flow[rev(i)] -= f;
if(a == 0) break;
}
}
if(a > 0) d[u] = -1;
return fl;
}
int tot;
int dinic() {
int ans = 0;
while(bfs()) {
for(int i = 0; i <= tot; i ++) cur[i] = first[i];
ans += dfs(s, 0x7fffffff);
}
return ans;
}
const int maxn = 1005;
ll A[maxn]; int B[maxn];
int main() {
process_S2();
scanf("%d", &n); tot = n + 1;
s = 0, t = tot;
std::vector<int> odd, even;
for(int i = 1; i <= n; i ++) {
scanf("%lld", &A[i]);
if(1LL & A[i]) {
odd.push_back(i);
} else {
even.push_back(i);
}
}
int ans = 0;
for(int i = 1; i <= n; i ++) {
scanf("%d", &B[i]); ans += B[i];
if(1LL & A[i]) {
ins_edge(i, t, B[i]);
} else {
ins_edge(s, i, B[i]);
}
}
for(int i = 0; i < even.size(); i ++) {
int p1 = even[i];
for(int j = 0; j < odd.size(); j ++) {
int p2 = odd[j];
if(gcd(A[p1], A[p2]) == 1LL && S2.find(A[p1] * A[p1] + A[p2] * A[p2]) != S2.end()) {
ins_edge(p1, p2, 0x7f7f7f7f);
}
}
}
ans -= dinic();
printf("%d\n", ans);
return 0;
}
[BZOJ 1061][NOI2008]志愿者招募
膜了一发BYVoid的题解……竟然搞懂了
这首先是个非常简单的线性规划，就是有若干限制（每天要求的志愿者），形如：
\[P_i = X_a + X_b +\ldots + X_c\geq A_i\]
（这里用\(X_i\)表示第\(i\)类志愿者雇佣了多少个，\(P_i\)表示第\(i\)天的志愿者总数，\(A_i\)同原题面）
且最小化总费用。
既然我们我知道\(P_i\geq A_i\)，那么说明\(P_i\)一定可以表示为\(A_i + Y_i\)（\(Y_i\geq 0\)）。然后这样限制就变成了：
\[P_i = X_a + X_b +\ldots + X_c + Y_i = A_i\]
这个长得很像可以流量平衡的样子，但是流量的变动是表示不了的……
然后假设\(P_0 = 0\)且\(P_{n + 1} = 0\)，这样就可以对限制差分一下，我们就有了\(n + 1\)个限制，然后这个式子就可以流量平衡做了（因为不考虑常数项的话，同一变量只可能在两个限制中出现，并且一正一负，这样就可以通过一个限制送给另一个限制流量来表示了。至于常数项，通过源或者汇连边即可表达）。
然后由于志愿者无限，所以这个东西也一定有解……
我局的我这么渣各位看官看懂的可能性基本上是零……还是推荐BYVoid神犇的题解，比我透彻多了。
代码：
/**************************************************************
Problem: 1061
User: danihao123
Language: C++
Result: Accepted
Time:3164 ms
Memory:6824 kb
****************************************************************/
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cctype>
#include <utility>
#include <algorithm>
#include <queue>
#include <cassert>
typedef long long ll;
int m, n;
const int maxno = 100005;
const int maxe = 100005;
int first[maxno];
int next[maxe], from[maxe], to[maxe]; ll flow[maxe], cap[maxe];
ll cost[maxe];
int gcnt = 0;
void add_edge(int u, int v, ll f, ll c) {
gcnt ++;
next[gcnt] = first[u];
first[u] = gcnt;
from[gcnt] = u; to[gcnt] = v;
cap[gcnt] = f; flow[gcnt] = 0;
cost[gcnt] = c;
}
int rev(int i) {
return ((i - 1) ^ 1) + 1;
}
int ins_edge(int u, int v, ll f, ll c = 0LL) {
#ifdef LOCAL
printf("Inserting Edge (%d, %d, %lld, %lld)\n", u, v, f, c);
#endif
add_edge(u, v, f, c);
add_edge(v, u, 0, -c);
return gcnt - 1;
}
const ll LINF = 0x7f7f7f7f7f7f7f7fLL;
int tot;
bool spfa(int s, int t, ll &f, ll &c) {
static ll d[maxno];
static bool inq[maxno];
static ll a[maxno]; static int p[maxno];
std::fill(d, d + tot + 1, LINF);
std::fill(inq, inq + tot + 1, false);
std::fill(a, a + tot + 1, 0LL);
std::fill(p, p + tot + 1, 0);
d[s] = 0;
std::queue<int> Q; Q.push(s); inq[s] = true;
a[s] = LINF;
while(!Q.empty()) {
int u = Q.front(); Q.pop();
inq[u] = false;
for(int i = first[u]; i; i = next[i]) {
if(cap[i] > flow[i]) {
int v = to[i];
if(d[v] > d[u] + cost[i]) {
d[v] = d[u] + cost[i];
a[v] = std::min(a[u], cap[i] - flow[i]); p[v] = i;
if(!inq[v]) {
Q.push(v); inq[v] = true;
}
}
}
}
}
if(a[t] == 0LL) return false;
f += a[t];
c += a[t] * d[t];
for(int e = p[t]; e; e = p[from[e]]) {
flow[e] += a[t];
flow[rev(e)] -= a[t];
}
return true;
}
void mcmf(int s, int t, ll &f, ll &c) {
while(spfa(s, t, f, c));
}
const int maxn = 1005, maxm = 10005;
int E[maxn];
int A[maxn], D[maxn];
int main() {
scanf("%d%d", &n, &m);
int s = 0, t = n + 2; tot = n + 2;
for(int i = 1; i <= n; i ++) scanf("%d", &A[i]);
for(int i = 1; i <= n + 1; i ++) {
D[i] = A[i] - A[i - 1];
if(D[i] >= 0) {
E[i] = ins_edge(s, i, D[i]);
} else {
E[i] = ins_edge(i, t, -D[i]);
}
}
for(int i = 1; i <= n; i ++) {
ins_edge(i + 1, i, LINF);
}
while(m --) {
int S, T, C; scanf("%d%d%d", &S, &T, &C);
ins_edge(S, T + 1, LINF, C);
}
ll f = 0, c = 0; mcmf(s, t, f, c);
printf("%lld\n", c);
return 0;
}
[BZOJ 3171][TJOI2013]循环格
流量平衡入门中……
我竟然想民白了……
这个题就是要求每个点在且仅在一个环中，这样每个点的入、出度都是1。出度肯定是1，接下来我们想想怎么保证入度为1。
我们建两排点（也可以说是拆点？），一排表示流出，一排表示接受。然后流出点向相邻的接收点连边，费用的话考虑是否和原来这一格的方向不一样就行了。
这个不需要判断是否满流啥的……因为一定有解（比如说每个行构成一个环啥的）。
代码：
/**************************************************************
Problem: 3171
User: danihao123
Language: C++
Result: Accepted
Time:28 ms
Memory:13844 kb
****************************************************************/
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cctype>
#include <algorithm>
#include <utility>
#include <queue>
const int maxn = 20;
typedef long long ll;
char S[maxn][maxn];
int tot = 0; int num[maxn][maxn][2];
int R, C;
char D[4] = {'L', 'R', 'U', 'D'};
int tow(int i, int j, char dir, int flag) {
int dx = 0, dy = 0;
if(dir == 'L') {
dx = -1;
} else if(dir == 'R') {
dx = 1;
} else if(dir == 'U') {
dy = -1;
} else if(dir == 'D') {
dy = 1;
}
int nx = (i + dy + R) % R;
int ny = (j + dx + C) % C;
return num[nx][ny][flag];
}
const int maxno = 100005;
const int maxe = 400005;
int first[maxno];
int next[maxe], from[maxe], to[maxe], flow[maxe], cap[maxe];
ll cost[maxe];
int gcnt = 0;
void add_edge(int u, int v, int f, ll c) {
gcnt ++;
next[gcnt] = first[u];
first[u] = gcnt;
from[gcnt] = u; to[gcnt] = v;
cap[gcnt] = f; flow[gcnt] = 0;
cost[gcnt] = c;
}
int rev(int i) {
return ((i - 1) ^ 1) + 1;
}
void ins_edge(int u, int v, int f, ll c = 0LL) {
#ifdef LOCAL
printf("Inserting Edge (%d, %d, %d, %lld)\n", u, v, f, c);
#endif
add_edge(u, v, f, c);
add_edge(v, u, 0, -c);
}
const ll LINF = 0x7f7f7f7f7f7f7f7fLL;
bool spfa(int s, int t, int &f, ll &c) {
static ll d[maxno];
static bool inq[maxno];
static int a[maxno], p[maxno];
std::fill(d, d + tot + 1, LINF);
std::fill(inq, inq + tot + 1, false);
std::fill(a, a + tot + 1, 0);
std::fill(p, p + tot + 1, 0);
d[s] = 0;
std::queue<int> Q; Q.push(s); inq[s] = true;
a[s] = 0x7fffffff;
while(!Q.empty()) {
int u = Q.front(); Q.pop();
inq[u] = false;
for(int i = first[u]; i; i = next[i]) {
if(cap[i] > flow[i]) {
int v = to[i];
if(d[v] > d[u] + cost[i]) {
d[v] = d[u] + cost[i];
a[v] = std::min(a[u], cap[i] - flow[i]); p[v] = i;
if(!inq[v]) {
Q.push(v); inq[v] = true;
}
}
}
}
}
if(a[t] == 0) return false;
f += a[t];
c += (ll(a[t])) * d[t];
for(int e = p[t]; e; e = p[from[e]]) {
flow[e] += a[t];
flow[rev(e)] -= a[t];
}
return true;
}
void mcmf(int s, int t, int &f, ll &c) {
while(spfa(s, t, f, c));
}
int main() {
tot = 1; int s = 0, t = 1;
scanf("%d%d", &R, &C);
for(int i = 0; i < R; i ++) {
scanf("%s", S[i]);
for(int j = 0; j < C; j ++) {
num[i][j][0] = ++ tot;
num[i][j][1] = ++ tot;
}
}
for(int i = 0; i < R; i ++) {
for(int j = 0; j < C; j ++) {
ins_edge(s, num[i][j][0], 1);
ins_edge(num[i][j][1], t, 1);
for(int it = 0; it < 4; it ++) {
int u = num[i][j][0];
int v = tow(i, j, D[it], 1);
ins_edge(u, v, 1, (D[it] == S[i][j]) ? 0LL : 1LL);
}
}
}
int f = 0; ll c = 0;
mcmf(s, t, f, c);
printf("%lld\n", c);
return 0;
}
[BZOJ 1070][SCOI2007]修车
啊啊啊我为什么要去颓费用流啊……
这个题可以保证每个人的车最后都会被修，所以说也就是让你最小化总等待时间辣！
直接考虑给一个技术人员建时间轴是麻烦而不可取的。我们考虑一点，对于每一个技术人员，若总共修了\(t\)辆车，那么第\(i\)个来找他修车的人会影响后面人的等待时间，换言之我们可以认为第\(i\)个来修车的人给答案做出了\(c\cdot (t - i + 1)\)（这里用\(c\)表示他自己的修车时间，注意他也会影响他自己）的贡献，似乎可以苟最小费用最大流呢。
但问题在于我们并不能提前知道\(t\)是多少，并且这样的话费用不是递增的，一个人可能会优先考虑到后面的位置去修车。那么我们考虑倒过来，去考虑倒数第\(i\)个来修车的人的贡献，显然这个贡献是\(c\cdot i\)。然后接下来就肥肠简单了，，，
代码：
/**************************************************************
Problem: 1070
User: danihao123
Language: C++
Result: Accepted
Time:576 ms
Memory:13928 kb
****************************************************************/
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cctype>
#include <utility>
#include <algorithm>
#include <queue>
#include <cassert>
typedef long long ll;
int m, n;
const int maxno = 100005;
const int maxe = 400005;
int first[maxno];
int next[maxe], from[maxe], to[maxe], flow[maxe], cap[maxe];
ll cost[maxe];
int gcnt = 0;
void add_edge(int u, int v, int f, ll c) {
gcnt ++;
next[gcnt] = first[u];
first[u] = gcnt;
from[gcnt] = u; to[gcnt] = v;
cap[gcnt] = f; flow[gcnt] = 0;
cost[gcnt] = c;
}
int rev(int i) {
return ((i - 1) ^ 1) + 1;
}
void ins_edge(int u, int v, int f, ll c = 0LL) {
#ifdef LOCAL
printf("Inserting Edge (%d, %d, %d, %lld)\n", u, v, f, c);
#endif
add_edge(u, v, f, c);
add_edge(v, u, 0, -c);
}
const ll LINF = 0x7f7f7f7f7f7f7f7fLL;
int tot;
bool spfa(int s, int t, int &f, ll &c) {
static ll d[maxno];
static bool inq[maxno];
static int a[maxno], p[maxno];
std::fill(d, d + tot + 1, LINF);
std::fill(inq, inq + tot + 1, false);
std::fill(a, a + tot + 1, 0);
std::fill(p, p + tot + 1, 0);
d[s] = 0;
std::queue<int> Q; Q.push(s); inq[s] = true;
a[s] = 0x7fffffff;
while(!Q.empty()) {
int u = Q.front(); Q.pop();
inq[u] = false;
for(int i = first[u]; i; i = next[i]) {
if(cap[i] > flow[i]) {
int v = to[i];
if(d[v] > d[u] + cost[i]) {
d[v] = d[u] + cost[i];
a[v] = std::min(a[u], cap[i] - flow[i]); p[v] = i;
if(!inq[v]) {
Q.push(v); inq[v] = true;
}
}
}
}
}
if(a[t] == 0) return false;
f += a[t];
c += (ll(a[t])) * d[t];
for(int e = p[t]; e; e = p[from[e]]) {
flow[e] += a[t];
flow[rev(e)] -= a[t];
}
return true;
}
void mcmf(int s, int t, int &f, ll &c) {
while(spfa(s, t, f, c));
}
int pos[105][105][2];
int main() {
scanf("%d%d", &m, &n);
int s = 0, t = 1; tot = 1;
for(int i = 1; i <= n; i ++) {
ins_edge(s, ++ tot, 1);
}
for(int i = 1; i <= m; i ++) {
for(int j = 1; j <= n; j ++) {
pos[i][j][0] = ++ tot;
pos[i][j][1] = ++ tot;
ins_edge(tot - 1, tot, 1);
ins_edge(tot, t, 1);
}
}
for(int i = 1; i <= n; i ++) {
for(int j = 1; j <= m; j ++) {
int T; scanf("%d", &T);
for(int k = 1; k <= n; k ++) {
int id = pos[j][k][0]; int id2 = pos[j][k][1];
ins_edge(i + 1, id, 1, k * T);
}
}
}
int f = 0; ll c = 0;
mcmf(s, t, f, c);
#ifdef LOCAL
printf("%d\n", f);
#endif
printf("%.2lf\n", (double(c)) / (double(n)));
return 0;
}
[BZOJ 1449][JSOI2009]球队收益
二次函数费用流的入门题……我真太弱了……
可以给比赛、队伍都建点，然后\(S\)向每个比赛连容量为1的边，每个比赛向队伍连容量为\(1\)的边，来表示赢家是谁。这样一来一次比赛只有一个队伍赢了。
考虑怎么处理那个二次函数费用。费用函数为\(f(x, y) = C x^2 + D y^2\)，由于一个队伍的总比赛次数是已知的，所以\(y\)不需要，不妨假设一个队伍有\(t\)场比赛，则费用函数为\(f(x) = C x^2 + D(t - x)^2\)。
对这个函数做差分：\(\Delta f(x) = f(x + 1) - f(x) = 2Cx - 2D(t - x) + C + D\)，然后这个差分是单调不降的。因此我们对所有差分都从队伍向\(T\)连边（费用为这一点的差分），这样一来的话会优先选择费用小的边，保证这个队伍的费用一直合法。
代码：
/**************************************************************
Problem: 1449
User: danihao123
Language: C++
Result: Accepted
Time:732 ms
Memory:1448 kb
****************************************************************/
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cctype>
#include <algorithm>
#include <utility>
#include <queue>
#include <cassert>
typedef long long ll;
const int maxn = 5005;
const int maxm = 1005;
const int maxno = maxn + maxm + 5;
const int maxe = (1000 * 2 + 1000 * 3) * 2 + 50;
int first[maxno];
int next[maxe], from[maxe], to[maxe], flow[maxe], cap[maxe];
ll cost[maxe];
int gcnt = 0;
void add_edge(int u, int v, int f, ll c) {
gcnt ++;
next[gcnt] = first[u];
first[u] = gcnt;
from[gcnt] = u; to[gcnt] = v;
cap[gcnt] = f; flow[gcnt] = 0;
cost[gcnt] = c;
}
int rev(int i) {
return ((i - 1) ^ 1) + 1;
}
void ins_edge(int u, int v, int f, ll c) {
add_edge(u, v, f, c);
add_edge(v, u, 0, -c);
}
int n, m;
const ll LINF = 0x7f7f7f7f7f7f7f7fLL;
bool spfa(int s, int t, int &f, ll &c) {
static ll d[maxno];
static bool inq[maxno];
static int a[maxno], p[maxno];
std::fill(d, d + n + m + 2, LINF);
std::fill(inq, inq + n + m + 2, false);
std::fill(a, a + n + m + 2, 0);
std::fill(p, p + n + m + 2, 0);
d[s] = 0;
std::queue<int> Q; Q.push(s); inq[s] = true;
a[s] = 0x7fffffff;
while(!Q.empty()) {
int u = Q.front(); Q.pop();
inq[u] = false;
for(int i = first[u]; i; i = next[i]) {
if(cap[i] > flow[i]) {
int v = to[i];
if(d[v] > d[u] + cost[i]) {
d[v] = d[u] + cost[i];
a[v] = std::min(a[u], cap[i] - flow[i]); p[v] = i;
if(!inq[v]) {
Q.push(v); inq[v] = true;
}
}
}
}
}
if(a[t] == 0) return false;
f += a[t];
c += (ll(a[t])) * d[t];
for(int e = p[t]; e; e = p[from[e]]) {
flow[e] += a[t];
flow[rev(e)] -= a[t];
}
return true;
}
void mcmf(int s, int t, int &f, ll &c) {
while(spfa(s, t, f, c));
}
ll win[maxn], lose[maxn], C[maxn], D[maxn];
ll tot[maxn];
int main() {
scanf("%d%d", &n, &m);
ll ans = 0;
for(int i = 1; i <= n; i ++) {
scanf("%lld%lld%lld%lld", &win[i], &lose[i], &C[i], &D[i]);
tot[i] = win[i] + lose[i];
}
int s = 0, t = n + m + 1;
for(int i = 1; i <= m; i ++) {
int a, b; scanf("%d%d", &a, &b);
ins_edge(0, i + n, 1, 0);
ins_edge(i + n, a, 1, 0);
ins_edge(i + n, b, 1, 0);
tot[a] ++; tot[b] ++;
}
for(int i = 1; i <= n; i ++) {
ans += C[i] * win[i] * win[i] + D[i] * (tot[i] - win[i]) * (tot[i] - win[i]);
for(ll j = win[i]; j <= (tot[i] - lose[i] - 1LL); j ++) {
ins_edge(i, t, 1, 2LL * C[i] * j - 2LL * D[i] * (tot[i] - j) + C[i] + D[i]);
}
}
int f = 0; mcmf(s, t, f, ans);
printf("%lld\n", ans);
return 0;
}
[BZOJ 1857][SCOI2010]传送带
三分套三分入门题……
策略肯定是从A走到AB上一点，然后再走到CD上的一个点，再向D走。
很显然答案函数是一个关于那两个点下凸的东西（不会证？GeoGebra之类的东西画一下就好啦！还不如像我这样口胡），所以我们可以先对第一维三分，然后套上对第二维的三分……
代码：
/**************************************************************
Problem: 1857
User: danihao123
Language: C++
Result: Accepted
Time:244 ms
Memory:820 kb
****************************************************************/
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cctype>
#include <algorithm>
#include <utility>
#include <cmath>
typedef double R;
const R eps = 1e-6;
int sign(R x) {
if(fabs(x) < eps) {
return 0;
} else {
if(x < 0) return -1;
else return 1;
}
}
struct Point {
R x, y;
Point(R qx = 0, R qy = 0) {
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 *(R x, const Vector &v) {
return Point(v.x * x, v.y * x);
}
Vector operator *(const Vector &v, R x) {
return Point(v.x * x, v.y * x);
}
R dot(const Vector &a, const Vector &b) {
return a.x * b.x + a.y * b.y;
}
R times(const Vector &a, const Vector &b) {
return a.x * b.y - a.y * b.x;
}
R dist(const Point &a, const Point &b) {
return sqrt(dot(a - b, a - b));
}
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;
}
}
Point A, B, C, D;
R p, q, r;
Vector D_AB, D_DC;
R f(const Point &AB, const Point &CD) {
return (dist(AB, A) / p + dist(CD, D) / q + dist(AB, CD) / r);
}
R F(Point AB) {
R L = 0, R = 1;
int T = 200;
while(T --) {
double M1 = L + (R - L) / 3;
double M2 = R - (R - L) / 3;
Point P1 = D + M1 * D_DC;
Point P2 = D + M2 * D_DC;
double f1 = f(AB, P1), f2 = f(AB, P2);
if(f1 < f2) {
R = M2;
} else {
L = M1;
}
}
return f(AB, D + L * D_DC);
}
R solve() {
R L = 0, R = 1;
int T = 200;
while(T --) {
double M1 = L + (R - L) / 3;
double M2 = R - (R - L) / 3;
Point P1 = A + M1 * D_AB;
Point P2 = A + M2 * D_AB;
double F1 = F(P1), F2 = F(P2);
if(F1 < F2) {
R = M2;
} else {
L = M1;
}
}
return F(A + L * D_AB);
}
int main() {
scanf("%lf%lf%lf%lf", &A.x, &A.y, &B.x, &B.y);
scanf("%lf%lf%lf%lf", &C.x, &C.y, &D.x, &D.y);
scanf("%lf%lf%lf", &p, &q, &r);
D_AB = B - A; D_DC = C - D;
printf("%.2lf\n", solve());
return 0;
}
[BZOJ 3944]Sum
杜教筛模板题……（然而常数卡卡卡）
我还是讲一讲基础的东西吧……
先考虑求欧拉函数的前缀和，我们先找一个它的卷积（这里用\(f=\varphi\ast 1\)，且设\(g = 1\)），发现有（用\(S\)表示原函数的前缀和）：
\[\sum_{i = 1}^n f(i) = \sum_{i = 1}^n g(i)\times S(\lfloor\frac{n}{i}\rfloor) = \frac{n(n + 1)}{2}\]
第二个和式的第一项是\(g(1)\times S(n) = S(n)\)，也就是我们需要的答案。这样的话我们把第一项分割出来，后面的项大力求，然后用\(\frac{n(n + 1)}{2}\)减去那些项就得到了第一项。并且注意到\(\lfloor\frac{n}{i}\rfloor\)的取值种类是\(\sqrt{n}\)级别的，又可以大力优化一波。然后考虑预处理\(O(n^{\frac{2}{3}})\)级别的前缀和（时间和空间都允许），再优化一波。然后考虑把没预处理的状态扔一个哈希表里（反正看起来不会很多），又优化一波……
这样复杂度就是\(O(n^{\frac{2}{3}})\)的了……（又不会证了）
求莫比乌斯函数的前缀和的方法类似（并且更简单一点）……
代码（又被卡成苟了）：
/**************************************************************
Problem: 3944
User: danihao123
Language: C++
Result: Accepted
Time:11292 ms
Memory:127740 kb
****************************************************************/
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cctype>
#include <algorithm>
#include <utility>
#include <ext/pb_ds/assoc_container.hpp>
#include <ext/pb_ds/hash_policy.hpp>
typedef long long ll;
const int N = 3500000;
const int maxn = N + 5;
ll phi[maxn], phi_S[maxn];
ll mu[maxn], mu_S[maxn];
void sieve() {
static bool vis[maxn];
static int prm[maxn]; int cnt = 0;
phi[1] = 1; mu[1] = 1;
for(int i = 2; i <= N; i ++) {
if(!vis[i]) {
prm[cnt ++] = i;
mu[i] = -1; phi[i] = i - 1;
}
for(int j = 0; j < cnt && prm[j] * (ll(i)) <= (ll(N)); j ++) {
int v = prm[j] * i;
vis[v] = true;
if(i % prm[j] == 0) {
mu[v] = 0; phi[v] = phi[i] * prm[j];
break;
} else {
mu[v] = mu[i] * -1;
phi[v] = phi[i] * phi[prm[j]];
}
}
}
for(int i = 1; i <= N; i ++) {
phi_S[i] = phi_S[i - 1] + phi[i];
mu_S[i] = mu_S[i - 1] + mu[i];
}
}
__gnu_pbds::gp_hash_table<ll, ll> mu_d;
ll mu_sum(ll n) {
if(n <= (ll(N))) return mu_S[n];
if(mu_d.find(n) != mu_d.end()) {
return mu_d[n];
}
ll ans = 1;
for(ll i = 2; i <= n;) {
ll nx = n / (n / i);
ans -= (nx - i + 1LL) * mu_sum(n / i);
i = nx + 1LL;
}
mu_d[n] = ans;
return ans;
}
ll pre_sum(ll n) {
return (n * (n + 1LL) / 2LL);
}
__gnu_pbds::gp_hash_table<ll, ll> phi_d;
ll phi_sum(ll n) {
if(n <= (ll(N))) return phi_S[n];
if(phi_d.find(n) != phi_d.end()) {
return phi_d[n];
}
ll ans = pre_sum(n);
for(ll i = 2; i <= n;) {
ll nx = n / (n / i);
ans -= (nx - i + 1LL) * phi_sum(n / i);
i = nx + 1LL;
}
phi_d[n] = ans;
return ans;
}
int main() {
sieve();
int T; scanf("%d", &T);
while(T --) {
int n; scanf("%d", &n);
printf("%lld %lld\n", phi_sum(n), mu_sum(n));
}
return 0;
}
[BZOJ 4816][SDOI2017]数字表格
当年在考场上一点反演都不会啊啊啊啊啊啊
这题虽然牵扯到了对积的反演，但其实也不是很难。先让我们大力推导一波（逃
考虑枚举柿子中的\(f[(i, j)]\)中的最大公约数（设为\(k\)），然后式子变成了（这里不妨偷个懒，钦定\(n\leq m\)）：
\[\prod_{k = 1}^n f[k]^{\sum_{i = 1}^n \sum_{j = 1}^m\:[(i, j) = k]}\]
然后发现上面那个指数就是Problem b的式子，显然可化为\(\sum_{d = 1}^n \mu(d) \lfloor\frac{n}{kd}\rfloor \lfloor\frac{m}{kd}\rfloor\)。然后似乎化简没有余地了，但是根据直觉，我们可以钦定\(T = kd\)，然后枚举\(T\)。
然后柿子变成了：
\[\prod_{T = 1}^n (\prod_{k\mid T} f[k]^{\mu(\tfrac{T}{k})})^{\lfloor\tfrac{n}{T}\rfloor \lfloor\tfrac{m}{T}\rfloor}\]
柿子中的\(\prod_{k\mid T} f[k]^{\mu(\tfrac{T}{k})}\)是一个类似狄利克雷卷积的东西（取完对数就是了），可以枚举约数预处理。并且这个东西很不错的一点就是上面的指数是莫比乌斯函数，可以取到的值只有那三种，不需要快速幂。
然后剩下的部分就和通常的反演题一般无二了……
代码（卡线过的蒟蒻……求轻喷……但我在LOJ上跑的还挺快的）：
/**************************************************************
Problem: 4816
User: danihao123
Language: C++
Result: Accepted
Time:50840 ms
Memory:33048 kb
****************************************************************/
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <utility>
#include <cctype>
#include <cassert>
const int N = 1000000;
const int maxn = N + 5;
typedef long long ll;
const ll ha = 1000000007LL;
bool p_flag = false;
ll pow_mod(ll a, ll b) {
if(!b) return 1LL;
ll ans = pow_mod(a, b / 2LL);
ans = (ans * ans) % ha;
if(1LL & b) ans = (ans * a) % ha;
#ifdef LOCAL
if(p_flag) printf("%lld^%lld : %lld\n", a, b, ans);
if(b == ha - 2LL) p_flag = false;
#endif
return ans;
}
ll inv(ll x) {
// if(x == 1LL) p_flag = true;
return pow_mod(x, ha - 2LL);
}
int miu[maxn];
void sieve() {
static int prm[maxn]; int cnt = 0;
static bool vis[maxn];
miu[1] = 1;
for(int i = 2; i <= N; i ++) {
if(!vis[i]) {
miu[i] = -1;
prm[cnt ++] = i;
}
for(int j = 0; j < cnt && prm[j] * i <= N; j ++) {
int v = prm[j] * i;
vis[v] = true;
if(i % prm[j] == 0) {
miu[v] = 0;
break;
} else {
miu[v] = miu[i] * -1;
}
}
}
}
ll f[maxn], inv_d[maxn], F[maxn];
void process() {
sieve();
f[1] = 1LL; inv_d[1] = 1LL;
for(int i = 2; i <= N; i ++) {
f[i] = (f[i - 1] + f[i - 2]) % ha;
inv_d[i] = inv(f[i]);
assert(f[i] != 0LL); assert(inv_d[i] != 0LL);
}
for(int i = 1; i <= N; i ++) F[i] = 1LL;
for(int i = 1; i <= N; i ++) {
for(int j = i; j <= N; j += i) {
int res = j / i;
if(miu[res] == 1) {
F[j] = (F[j] * f[i]) % ha;
#ifdef LOCAL
if(j <= 3) printf("f[%d](%lld) -> F[%d]\n", i, f[i], j);
#endif
continue;
}
if(miu[res] == -1) {
F[j] = (F[j] * inv_d[i]) % ha;
#ifdef LOCAL
if(j <= 3) printf("inv_f[%d](%lld) -> F[%d]\n", i, inv_d[i], j);
#endif
}
}
}
F[0] = 1LL;
bool flag = true;
for(int i = 1; i <= N; i ++) {
F[i] = (F[i] * F[i - 1]) % ha;
if(flag && F[i] == 0LL) {
printf("SF[%d] has been 0.\n", i);
flag = false;
}
}
}
ll calc(ll n, ll m) {
if(n > m) std::swap(n, m);
ll ans = 1LL;
for(ll i = 1; i <= n;) {
ll nx = std::min(n / (n / i), m / (m / i));
#ifdef LOCAL
// printf("nx of %lld : %lld\n", i, nx);
#endif
ll mulv = pow_mod((F[nx] * inv(F[i - 1])) % ha, (n / i) * (m / i));
ans = (ans * mulv) % ha;
i = nx + 1LL;
}
return ans;
}
int main() {
process();
int T; scanf("%d", &T);
while(T --) {
int n, m; scanf("%d%d", &n, &m);
printf("%lld\n", calc(n, m));
}
return 0;
}
[BZOJ 5059]前鬼后鬼的守护
这不是可并堆裸题吗……我这种傻逼只会用线段树合并做
显然一个集合的中位数可以用平衡树或者权值线段树求出，然后为了方便合并我用了权值线段树。
然后考虑对于一个递减的序列，选择中位数一定更优，一个递增的序列当然xjb选即可。我们可以用一个栈来维护当前的所有答案区间，然后考虑一个新的点加进来。这个点的答案如果比栈首答案小，那么就违反了单调不降的要求，我们需要把两段合并。然后完了……
代码：
/**************************************************************
Problem: 5059
User: danihao123
Language: C++
Result: Accepted
Time:2348 ms
Memory:250520 kb
****************************************************************/
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cctype>
#include <cassert>
#include <algorithm>
#include <utility>
#include <stack>
const int BUFSIZE = 20 * 1024 * 1024;
int n;
struct Node {
Node *lc, *rc;
int sumv;
};
Node pool[BUFSIZE];
Node *nil, *cur;
void init_pool() {
cur = nil = pool;
nil -> lc = nil -> rc = nil;
nil -> sumv = 0;
}
Node *alloc_node(int v = 0, Node *lc = nil, Node *rc = nil) {
cur ++;
cur -> lc = lc; cur -> rc = rc;
cur -> sumv = v;
return cur;
}
void modify(Node* &o, int L, int R, int p, int v) {
if(o == nil) {
o = alloc_node(0);
}
o -> sumv += v;
if(L < R) {
int M = (L + R) / 2;
if(p <= M) {
modify(o -> lc, L, M, p, v);
} else {
modify(o -> rc, M + 1, R, p, v);
}
}
}
Node *merge(Node* &tag, Node* &src) {
if(src == nil) return tag;
if(tag == nil) {
return src;
}
tag -> sumv += src -> sumv;
tag -> lc = merge(tag -> lc, src -> lc);
tag -> rc = merge(tag -> rc, src -> rc);
return tag;
}
int kth(Node *o, int L, int R, int k) {
if(L == R) {
return L;
} else {
int M = (L + R) / 2;
if(k <= o -> lc -> sumv) {
return kth(o -> lc, L, M, k);
} else {
return kth(o -> rc, M + 1, R, k - o -> lc -> sumv);
}
}
}
void gou_tree(Node *o, int L, int R) {
printf("Tree (%d, %d, %d)\n", L, R, o -> sumv);
int M = (L + R) / 2;
if(L < R) {
gou_tree(o -> lc, L, M);
gou_tree(o -> rc, M + 1, R);
}
}
int lsiz;
const int maxn = 500005;
int A[maxn], A2[maxn];
void break_up() {
std::sort(A2 + 1, A2 + 1 + n);
lsiz = std::unique(A2 + 1, A2 + 1 + n) - A2 - 1;
for(int i = 1; i <= n; i ++) {
A[i] = std::lower_bound(A2 + 1, A2 + 1 + lsiz, A[i]) - A2;
}
}
struct Interval {
Node *rt;
int l, r; int siz, ans;
Interval() {
rt = nil;
siz = ans = 0;
}
Interval(int p) {
rt = nil; modify(rt, 1, lsiz, p, 1);
siz = 1; ans = p;
}
};
typedef long long ll;
ll solve() {
break_up(); init_pool();
std::stack<Interval> S;
for(int i = 1; i <= n; i ++) {
Interval nv(A[i]);
nv.l = nv.r = i;
while(!S.empty() && S.top().ans > nv.ans) {
Interval g = S.top(); S.pop();
#ifdef LOCAL
printf("ans of [%d, %d] : %d\n", g.l, g.r, A2[g.ans]);
#endif
nv.l = g.l; nv.siz += g.siz;
nv.rt = merge(nv.rt, g.rt);
#ifdef LOCAL
gou_tree(nv.rt, 1, lsiz);
#endif
nv.ans = kth(nv.rt, 1, lsiz, (nv.siz + 1) / 2);
}
S.push(nv);
}
ll ans = 0;
while(!S.empty()) {
Interval g = S.top(); S.pop();
#ifdef LOCAL
printf("ans of [%d, %d] : %d\n", g.l, g.r, A2[g.ans]);
#endif
for(int i = g.l; i <= g.r; i ++) {
ans += abs(A2[A[i]] - A2[g.ans]);
}
}
return ans;
}
int main() {
scanf("%d", &n);
for(int i = 1; i <= n; i ++) {
scanf("%d", &A[i]);
A2[i] = A[i];
}
printf("%lld\n", solve());
return 0;
}
[BZOJ 3527][ZJOI2014]力
现在才明白卷积真的是the deep, dark fantasies……（逃
首先约去\(q_i\)，得到：
\[E_j = \sum_{i < j}\frac{q_i}{(j - i)^2} - \sum_{j < i}\frac{q_i}{(j - i)^2}\]
注意到如果很容易得到等式前半部分的高效求法，后半部分把序列反过来就能做了。
那么我们会发现，设\(g_x = \frac{1}{x^2}\)，然后式子前半部分（姑且称作\(p_j\)）可表示为：
\[p_j = \sum_{i < j} q_i g_{j - i}\]
这不就是个卷积吗？FFT一波即可。
代码：
/**************************************************************
Problem: 3527
User: danihao123
Language: C++
Result: Accepted
Time:9984 ms
Memory:28952 kb
****************************************************************/
#include <cstdio>
#include <cstring>
#include <cctype>
#include <cstdlib>
#include <cmath>
#include <algorithm>
#include <utility>
#include <queue>
#include <cassert>
#include <complex>
typedef long double R;
typedef std::complex<R> C;
const R eps = 1e-3;
int sign(R x) {
if(fabs(x) < eps) {
return 0;
} else {
if(x < 0) {
return -1;
} else {
return 1;
}
}
}
const int maxn = 400005;
const double pi = acos(-1);
int bi, len;
int flip(int x) {
int ans = 0;
for(int i = 0; i < bi; i ++) {
if(x & (1 << i)) {
ans += 1 << (bi - i - 1);
}
}
return ans;
}
void fft(C *A, double g = 1) {
for(int i = 0; i < len; i ++) {
int R = flip(i);
#ifdef LOCAL
// printf("The flipping of %d is %d\n", i, R);
#endif
if(i < R) {
std::swap(A[R], A[i]);
}
}
for(int L = 1; L < len; L <<= 1) {
C xi_n(cos(pi / (double(L))), sin(g * pi / (double(L))));
for(int i = 0; i < len; i += L << 1) {
C xi(1, 0);
for(int j = i; j < i + L; j ++) {
C a = A[j], b = A[j + L];
A[j] = a + xi * b;
A[j + L] = a - xi * b;
xi = xi * xi_n;
}
}
}
}
int main() {
int n;
static C A[maxn], rd[maxn]; static R ans[maxn];
static R q[maxn];
scanf("%d", &n);
bi = 0; len = 1;
while(len <= 2 * n) {
len <<= 1;
bi ++;
}
for(int i = 0; i < n; i ++) {
scanf("%Lf", &q[i]); A[i] = q[i];
}
assert(sign(rd[0].real()) == 0);
for(int i = 1; i < n; i ++) {
rd[i] = 1.00 / ((R(i)) * (R(i)));
#ifdef LOCAL
printf("rd[%d] : %.5Lf\n", i, rd[i].real());
#endif
}
fft(A); fft(rd);
for(int i = 0; i < len; i ++) A[i] *= rd[i];
fft(A, -1);
for(int i = 0; i < n; i ++) {
ans[i] += A[i].real() / len;
#ifdef LOCAL
printf("delta_v of %d : %.5Lf\n", i, A[i].real() / len);
#endif
}
std::reverse(q, q + n);
for(int i = 0; i < len; i ++) A[i] = 0;
for(int i = 0; i < n; i ++) {
A[i] = q[i];
}
fft(A);
for(int i = 0; i < len; i ++) A[i] *= rd[i];
fft(A, -1);
for(int i = 0; i < n; i ++) {
ans[i] -= A[n - 1 - i].real() / len;
printf("%.3Lf\n", ans[i]);
}
return 0;
}
