[BZOJ 4919][Lydsy1706月赛]大根堆
ao神啊这题，，，
很显然的思路是设计状态，然后线段树合并或者线段树启发式合并转移来优化一下，但是真的不好写……
考虑求LIS的那个二分+单调数组来做，那个东西的本质其实就是维护了一堆各种长度的断点一样的东西（考虑答案为\(n\)，\(n - 1\)……的情况，在值的选取上各个情况之间会有断点）。那么我们就用一个multiset来维护断点，然后启发式合并即可。
考虑一棵子树，把根的子树进行合并，根加进来之后对断点会有何影响？那么考虑现在那个multiset里那个东西的后继（要lower_bound，因为如果有和根权值一样的点的话拐点不会发生变化。这里讨论有后继的情况），如果要取到那个点的长度的LIS的话，现在只需要小于等于根的权值的点就行了，取到这个点也不会使LIS变大（它不能和根共存），所以那个点就不是拐点了。并且无论是否出现这种情况，根那个点都会成为一个新的拐点。
然后DFS一遍就行啦。
代码：
/**************************************************************
Problem: 4919
User: danihao123
Language: C++
Result: Accepted
Time:992 ms
Memory:25648 kb
****************************************************************/
#include <cstdio>
#include <set>
#include <vector>
const int maxn = 200005;
std::vector<int> G[maxn];
std::multiset<int> *st[maxn];
void merge(int x, int y) {
if(st[x] -> size() < st[y] -> size()) std::swap(st[x], st[y]);
while(st[y] -> size() > 0) {
std::multiset<int>::iterator it = st[y] -> begin();
int v = *it; st[y] -> erase(it);
st[x] -> insert(v);
}
}
int d[maxn];
void dfs(int x) {
for(int i = 0; i < G[x].size(); i ++) {
int v = G[x][i];
dfs(v);
merge(x, v);
}
std::multiset<int>::iterator it = st[x] -> lower_bound(d[x]);
if(it != (st[x] -> end())) st[x] -> erase(it);
st[x] -> insert(d[x]);
}
int main() {
int n; scanf("%d", &n);
for(int i = 1; i <= n; i ++) {
int f; scanf("%d%d", &d[i], &f);
G[f].push_back(i);
}
for(int i = 1; i <= n; i ++) {
st[i] = new std::multiset<int>();
}
dfs(1);
printf("%d\n", st[1] -> size());
for(int i = 1; i <= n; i ++) {
delete st[i];
}
return 0;
}
[洛谷 P4177][CEOI2008]order
啊啊啊只会做水题了……
很显然是最小割吧……长得很像最大权闭合子图，但又不一样的。
那么就把最大权闭合子图中的无穷边（我称之为违约边）的费用改成租用的费用。这样一来，如果选了计划（不割计划）还不购买机器（割去机器），那就只能支付租金了（割违约边）。
这其实也算是一种套路了吧？觉得以前见过很多次但是有很多叫法……
代码：
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cctype>
#include <utility>
#include <algorithm>
#include <queue>
const int maxn = 200005;
const int maxm = 2000005;
const int maxe = maxm << 2;
int n, m;
int first[maxn];
int next[maxe], to[maxe], flow[maxe], cap[maxe];
int gcnt = 0;
inline void add_edge(int u, int v, int f) {
gcnt ++;
next[gcnt] = first[u];
first[u] = gcnt;
to[gcnt] = v;
flow[gcnt] = 0; cap[gcnt] = f;
}
inline int rev(int i) {
return (((i - 1) ^ 1) + 1);
}
inline void ins_edge(int u, int v, int f) {
add_edge(u, v, f);
add_edge(v, u, 0);
}
int d[maxn];
int s, t, tot;
inline bool bfs() {
static bool vis[maxn];
std::fill(vis, vis + 1 + tot, false);
std::fill(d, d + 1 + tot, 0);
std::queue<int> Q;
Q.push(s); vis[s] = true; d[s] = 1;
while(!Q.empty()) {
int u = Q.front(); Q.pop();
for(int i = first[u]; i; i = next[i]) {
int v = to[i];
if(!vis[v] && cap[i] > flow[i]) {
d[v] = d[u] + 1; vis[v] = true;
Q.push(v);
}
}
}
return vis[t];
}
int cur[maxn];
int dfs(int x, int a) {
if(a == 0 || x == t) return a;
int ans = 0;
for(int &i = cur[x]; i; i = next[i]) {
int v = to[i];
int f;
if(d[v] == d[x] + 1 && (f = dfs(v, std::min(a, cap[i] - flow[i]))) > 0) {
ans += f; a -= f;
flow[i] += f; flow[rev(i)] -= f;
if(a == 0) break;
}
}
if(a > 0) d[x] = -1;
return ans;
}
int dinic() {
int ans = 0;
while(bfs()) {
for(int i = 0; i <= tot; i ++) cur[i] = first[i];
ans += dfs(s, 0x7fffffff);
}
return ans;
}
int main() {
int n, m; scanf("%d%d", &n, &m);
s = 0; t = tot = n + m + 1;
int ans = 0;
for(int i = 1; i <= n; i ++) {
int w, p; scanf("%d%d", &w, &p);
ins_edge(s, i, w); ans += w;
for(int j = 1; j <= p; j ++) {
int id, co; scanf("%d%d", &id, &co);
ins_edge(i, id + n, co);
}
}
for(int i = n + 1; i <= n + m; i ++) {
int co; scanf("%d", &co);
ins_edge(i, t, co);
}
printf("%d\n", ans - dinic());
return 0;
}
[LibreOJ 2137][ZJOI2015]诸神眷顾的幻想乡
复习SAM了呢~
那个度数为1的点至多有20个的条件非常神奇……让我们想想怎么用。
我们发现，（钦定根后）在树上有一条路径是弯的这种情况非常不好统计，但如果是直着下来就很好说（把所有从根到一个点的路径扔到SAM里，然后是经典题）。但是，任何一条路一定会在某个度数为1的点为根的情况下是直的（可以意会一下吧（逃
然后我们从那20个点每个点当根DFS一遍，把搞出来的每一条从根开始的路径放前缀Trie里。然后对前缀Trie构一个SAM就行啦~
代码：
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cctype>
#include <algorithm>
#include <utility>
#include <queue>
#include <vector>
const int BUFSIZE = 128 * 1024 * 1024;
char buf[BUFSIZE];
void *alloc(size_t size) {
static char *cur = buf;
if(cur - buf + size > BUFSIZE) {
return malloc(size);
} else {
char *ret = cur; cur += size;
return ret;
}
}
const int maxn = 100005;
std::vector<int> G[maxn];
int deg[maxn];
int ssiz;
struct TNode {
TNode *ch[10];
};
TNode *alloc_tn() {
auto ret = (TNode*)alloc(sizeof(TNode));
memset(ret -> ch, 0, sizeof(ret -> ch));
return ret;
}
TNode *step(TNode *o, int c) {
if(!(o -> ch[c])) {
o -> ch[c] = alloc_tn();
}
return (o -> ch[c]);
}
TNode *trt;
int col[maxn];
void dfs(int x, int fa, TNode *last) {
TNode *np = step(last, col[x]);
for(auto v : G[x]) {
if(v != fa) {
dfs(v, x, np);
}
}
}
struct Node {
int len; Node *fa;
Node *ch[10];
};
std::vector<Node*> pool;
Node *alloc_node(int len = 0, Node *fa = NULL) {
Node *ret = (Node*)alloc(sizeof(Node));
ret -> len = len; ret -> fa = fa;
memset(ret -> ch, 0, sizeof(ret -> ch));
pool.push_back(ret);
return ret;
}
Node *rt;
Node *extend(int c, Node *last) {
Node *np = alloc_node(last -> len + 1);
Node *p = last;
while(p != NULL && p -> ch[c] == NULL) {
p -> ch[c] = np; p = p -> fa;
}
if(p == NULL) {
np -> fa = rt;
} else {
Node *q = p -> ch[c];
if(q -> len == p -> len + 1) {
np -> fa = q;
} else {
Node *nq = alloc_node(p -> len + 1, q -> fa);
memcpy(nq -> ch, q -> ch, sizeof(q -> ch));
q -> fa = np -> fa = nq;
while(p != NULL && p -> ch[c] == q) {
p -> ch[c] = nq; p = p -> fa;
}
}
}
return np;
}
void dfs_2(TNode *o, Node *last) {
for(int i = 0; i < ssiz; i ++) {
TNode *v = o -> ch[i];
if(!v) continue;
Node *np = extend(i, last);
dfs_2(v, np);
}
}
using ll = long long;
int main() {
int n; scanf("%d%d", &n, &ssiz);
for(int i = 1; i <= n; i ++) {
scanf("%d", &col[i]);
}
trt = alloc_tn();
for(int i = 1; i <= n - 1; i ++) {
int u, v; scanf("%d%d", &u, &v);
G[u].push_back(v); G[v].push_back(u);
deg[u] ++; deg[v] ++;
}
for(int i = 1; i <= n; i ++) {
if(deg[i] == 1) {
dfs(i, 0, trt);
}
}
rt = alloc_node();
dfs_2(trt, rt);
ll ans = 0;
for(auto p : pool) {
if(p -> fa) {
ans += p -> len - p -> fa -> len;
}
}
printf("%lld\n", ans);
return 0;
}
[LibreOJ 2034][SDOI2016]排列计数
我这种池沼只会做水题陶冶身心了……
考虑用\(f_i\)表示长为\(i\)的完全错位全排列的个数，那么有\(i\)个错位的长度为\(n\)的全排列的数量就是\(\binom{n}{i} f_i\)。从这一点出发，可以得到一个式子（枚举错位的有几个）：
\[n! = \sum_{i = 0}^n \binom{n}{i} f_i\]
考虑使用二项式反演，得到：
\[
\begin{aligned}
f_n &= \sum_{i = 0}^n (-1)^{n - i}\binom{n}{i} i!\\
&= \sum_{i = 0}^n (-1)^{n - i}\frac{n!}{(n - i)!}\\
&= n!\sum_{i = 0}^n \frac{(-1)^{n - i}}{(n - i)!}\\
&= n!\sum_{i = 0}^n \frac{(-1)^i}{i!}
\end{aligned}
\]
后面的和式可以预处理，然后就做完啦~
代码：
#include <cstdio>
#include <cstring>
using ll = long long;
const ll ha = 1000000007LL;
inline ll pow_mod(const ll &a, ll b) {
ll ans = 1LL, res = a;
while(b) {
if(1LL & b) ans = (ans * res) % ha;
res = (res * res) % ha;
b >>= 1;
}
return ans;
}
inline ll inv(ll x) {
return pow_mod(x, ha - 2LL);
}
const int N = 1000000;
const int maxn = N + 5;
ll fac[maxn], f[maxn];
inline void process() {
fac[0] = 1LL;
for(int i = 1; i <= N; i ++) {
fac[i] = (fac[i - 1] * (ll(i))) % ha;
}
for(int i = 0; i <= N; i ++) {
f[i] = (i % 2 == 1) ? (ha - 1LL) : 1LL;
f[i] = (f[i] * inv(fac[i])) % ha;
if(i > 0) f[i] = (f[i - 1] + f[i]) % ha;
}
}
int main() {
process();
int T; scanf("%d", &T);
while(T --) {
int n, m; scanf("%d%d", &n, &m);
int p = n - m;
ll bs = (fac[p] * f[p]) % ha;
ll C = (fac[n] * inv((fac[p] * fac[m]) % ha)) % ha;
printf("%lld\n", (bs * C) % ha);
}
return 0;
}
[BZOJ 5093][Lydsy1711月赛]图的价值
统计的时候，考虑每个店对每个图的贡献，可以看出答案是：
\[n2^{\tfrac{(n - 1)(n - 2)}{2}}\,\sum_{i = 0}^{n - 1} \binom{n - 1}{i} i^k\]
求和号前面还好说，主要看求和号后面的咋搞。
后面的那一部分是个经典题吧……但我还是来推导一番吧：
\[
\begin{aligned}
\sum_{i = 0}^{m} \binom{m}{i} i^k&= \sum_{i = 0}^{m} \binom{m}{i}\sum_{j = 0}^k S(k, j) j!\binom{i}{j}\\
&= \sum_{j = 0}^k S(k, j) j! \sum_{i = 0}^m \binom{m}{i}\binom{i}{j} \\
&= \sum_{j = 0}^k S(k, j) m^{\underline{j}} \sum_{i = 0}^{m} \binom{n - j}{i - j} \\
&= \sum_{j = 0}^k S(k, j) m^{\underline{j}} 2^{m - j}
\end{aligned}
\]
然后如果我们能高效的求出某一行的第二类斯特林数，那么问题就迎刃而解了。
考虑斯特林反演（这本质上是一个容斥，用二项式反演也可以推导）：
\[S(k, j) = \sum_{i = 0}^j \frac{(-1)^{j - i}}{(j - i)!} \cdot \frac{i^k}{i!}\]
这个东西很显然是一个卷积……而且模数还很友好（UOJ模数），用NTT求出一行的第二类斯特林数即可。
代码：
/**************************************************************
Problem: 5093
User: danihao123
Language: C++
Result: Accepted
Time:22632 ms
Memory:25824 kb
****************************************************************/
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cctype>
#include <algorithm>
#include <utility>
typedef long long ll;
const ll ha = 998244353LL;
const ll bs = 3LL;
inline ll pow_mod(const ll &a, ll b) {
ll ans = 1LL, res = a;
while(b) {
if(1LL & b) ans = (ans * res) % ha;
res = (res * res) % ha;
b >>= 1;
}
return ans;
}
inline ll inv(const ll &x) {
return pow_mod(x, ha - 2LL);
}
inline int flip(const int &bi, const int &x) {
int ans = 0;
for(int i = 0; i < bi; i ++) {
if((1 << i) & x) {
ans |= (1 << (bi - i - 1));
}
}
return ans;
}
inline void ntt(ll *A, const int &bi, const int &len, bool flag = false) {
for(int i = 0; i < len; i ++) {
int v = flip(bi, i);
if(v < i) std::swap(A[v], A[i]);
}
for(int L = 1; L < len; L <<= 1) {
ll xi = pow_mod(bs, (ha - 1LL) / (ll(L << 1)));
if(flag) xi = inv(xi);
for(int i = 0; i < len; i += (L << 1)) {
ll w = 1LL;
for(int j = i; j < i + L; j ++) {
ll x = A[j], y = A[j + L];
A[j] = (x + (w * y) % ha) % ha;
A[j + L] = (x - (w * y) % ha + ha) % ha;
w = (w * xi) % ha;
}
}
}
}
const int maxn = 800005;
ll A[maxn], B[maxn];
ll fac[maxn], down[maxn];
int main() {
ll n; int k; scanf("%lld%d", &n, &k);
if(n == 1) {
puts("0"); return 0;
}
int len = 1, bi = 0;
while(len <= (2 * k)) {
len <<= 1; bi ++;
}
fac[0] = 1LL;
for(int i = 1; i <= k; i ++) {
fac[i] = (fac[i - 1] * (ll(i))) % ha;
}
for(int i = 0; i <= k; i ++) {
A[i] = (inv(fac[i]) * pow_mod(ha - 1LL, i)) % ha;
}
for(int i = 0; i <= k; i ++) {
B[i] = inv(fac[i]);
B[i] = (B[i] * pow_mod(i, k)) % ha;
}
ntt(A, bi, len); ntt(B, bi, len);
for(int i = 0; i < len; i ++) {
A[i] = (A[i] * B[i]) % ha;
}
ntt(A, bi, len, true);
ll inv_l = inv(len);
for(int i = 0; i < len; i ++) {
A[i] = (A[i] * inv_l) % ha;
}
#ifdef LOCAL
for(int i = 0; i <= k; i ++) {
printf("%lld ", A[i]);
}
puts("");
#endif
down[0] = 1LL;
for(int i = 1; i <= std::min(ll(k), n - 1LL); i ++) {
down[i] = (down[i - 1] * (n - (ll(i)))) % ha;
}
ll ans = 0LL;
for(int i = 1; i <= std::min(ll(k), n - 1LL); i ++) {
ll delta = (A[i] * down[i]) % ha;
delta = (delta * pow_mod(2, n - 1LL - (ll(i)))) % ha;
ans = (ans + delta) % ha;
}
ans = (ans * pow_mod(2LL, (n - 1LL) * (n - 2LL) / 2LL)) % ha;
ans = (ans * n) % ha;
printf("%lld\n", ans);
return 0;
}
[HDU 5780]gcd
有个很好康的结论：
\[\gcd(x^a - 1, x^b - 1) = x^{\gcd(a, b)} - 1\]
然后尝试去用常见套路（枚举gcd）化简柿子，得到：
\[\sum_{k = 1}^n (x^k - 1)\sum_{1\leq a,b\leq n} [\gcd(a, b) = k]\]
推到这里了，看到那个\(\sum_{1\leq a,b\leq n} [\gcd(a, b) = k]\)可能很多同学准备直接上莫比乌斯函数了……其实并不需要，其实那个柿子就是：
\[2\sum_{i = 1}^{\lfloor \tfrac{n}{k}\rfloor} \varphi (i) - 1\]
这个意义还是很显然的……更好的一点是这个柿子可以预处理出来，然后整除分块直接搞一波即可。
代码：
#include <cstdio>
typedef long long ll;
const ll ha = 1000000007LL;
const int maxn = 1000005;
const int N = 1000000;
ll phi[maxn], phi_S[maxn];
void sieve() {
static bool vis[maxn];
static int prm[maxn]; int cnt = 0;
phi[1] = 1LL;
for(int i = 2; i <= N; i ++) {
if(!vis[i]) {
prm[cnt ++] = i;
phi[i] = i - 1;
}
int v;
for(int j = 0; j < cnt && (v = i * prm[j]) <= N; j ++) {
vis[v] = true;
if(i % prm[j] == 0) {
phi[v] = (phi[i] * (ll(prm[j]))) % ha;
break;
} else {
phi[v] = (phi[i] * phi[prm[j]]) % ha;
}
}
}
for(int i = 1; i <= N; i ++) {
phi_S[i] = (phi_S[i - 1] + phi[i]) % ha;
}
for(int i = 1; i <= N; i ++) {
phi_S[i] = ((2LL * phi_S[i]) % ha - 1LL + ha) % ha;
}
}
ll pow_mod(ll a, ll b) {
ll ans = 1LL, res = a;
while(b) {
if(1LL & b) ans = (ans * res) % ha;
res = (res * res) % ha;
b >>= 1;
}
return ans;
}
ll inv(ll x) {
return pow_mod(x, ha - 2LL);
}
ll pre_sum(ll x, ll r) {
ll p = pow_mod(x, r);
p = (p - 1LL + ha) % ha;
p = (p * inv(x - 1LL)) % ha;
return p;
}
ll seg_sum(ll x, int a, int b) {
if(x == 1LL) return 0;
ll len = b - a + 1;
ll ret = (pre_sum(x, b + 1) - pre_sum(x, a) + ha) % ha;
ret = (ret - len + ha) % ha;
return ret;
}
ll calc(ll x, int n) {
ll ans = 0;
for(int i = 1; i <= n;) {
int v = n / i;
int nx = n / v;
ll delta = (seg_sum(x, i, nx) * phi_S[v]) % ha;
ans = (ans + delta) % ha;
i = nx + 1;
}
return ans;
}
int main() {
sieve();
int T; scanf("%d", &T);
while(T --) {
int x, n; scanf("%d%d", &x, &n);
printf("%lld\n", calc(x, n));
}
return 0;
}
[51Nod 1355]斐波那契的最小公倍数
MIN-MAX容斥第一题……
众所周知斐波那契数有个非常好的性质：
\[\gcd(\{f_i | i\in T\}) = f_{\gcd\{T\}}\]
但很不幸的事情是，lcm没有类似的性质。那我们就考虑怎么把lcm改成gcd。
根据MIN-MAX容斥，可以得到这个式子（可以理解为对指数的容斥。这里设\(S\)为全集）：
\[\prod_{T\subseteq S, T\neq\emptyset} f_{\gcd\{T\}}^{(-1)^{|T| + 1}}\]
gcd套在上面非常烦人，考虑用狄利克雷卷积去化解它。考虑有个函数\(g\)满足\(f_n = \prod_{d|n} g(d)\)，然后xjb化简一下之后柿子变成：
\[\prod_{d} g(d)^{\sum_{T\subseteq S, T\neq\emptyset, \forall x\in T, d | x}(-1)^{|T| + 1}}\]
再去考虑那个指数。当且仅当\(d\in S\)时，\(g(d)\)才会对答案做贡献。并且无论\(S\)中\(d\)的约数有多少个（只要大于），最后指数肯定都是1。因为其实大于零的情况那个指数也是一个所有最小值都是1的MIN-MAX容斥……这样求出来最大值肯定是1了，那个指数自然就是1。
于是乎柿子变成了：
\[\prod_{\exists x\in S, d | x} g(d)\]
然后搞出来\(S\)中所有数的约数很容易。当务之急就是求\(g\)了（这里其实就是莫比乌斯反演了）。不过球莫比乌斯函数反而让问题复杂化了……移下项可以发现：
\[g(n) = f_n\cdot (\prod_{d | n, d\neq n} g(d))^{-1}\]
然后就很愉快的搞出来了（逃
代码：
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cmath>
#include <cctype>
#include <algorithm>
#include <utility>
#include <queue>
const int maxn = 1000005;
using ll = long long;
const ll ha = 1000000007LL;
ll pow_mod(ll a, ll b) {
ll ans = 1LL, res = a;
while(b) {
if(1LL & b) {
ans = (ans * res) % ha;
}
res = (res * res) % ha;
b /= 2LL;
}
return ans;
}
ll inv(ll x) {
return pow_mod(x, ha - 2LL);
}
ll f[maxn], g[maxn];
void process_f() {
int N = 1000000;
f[0] = 0LL; f[1] = 1LL;
for(int i = 2; i <= N; i ++) {
f[i] = (f[i - 1] + f[i - 2]) % ha;
}
std::fill(g + 1, g + 1 + N, 1LL);
for(int i = 1; i <= N; i ++) {
g[i] = (f[i] * inv(g[i])) % ha;
for(int j = i * 2; j <= N; j += i) {
g[j] = (g[j] * g[i]) % ha;
}
}
#ifdef LOCAL
puts("g has been processed!");
#endif
}
bool vis[maxn];
void process_vis() {
int N = 1000000;
for(int i = 1; i <= N; i ++) {
for(int j = i * 2; j <= N; j += i) {
vis[i] = (vis[i] || vis[j]);
}
}
}
int main() {
process_f();
int n; scanf("%d", &n);
for(int i = 1; i <= n; i ++) {
int a; scanf("%d", &a);
vis[a] = true;
}
process_vis();
ll ans = 1LL;
for(int i = 1; i <= 1000000; i ++) {
if(vis[i]) {
ans = (ans * g[i]) % ha;
}
}
printf("%lld\n", ans);
return 0;
}
[BZOJ 4259]残缺的字符串
终于肝掉了这道提……
考虑用编辑距离函数处理这个问题。考虑每个可能的结尾\(i\)，定义\(f_i\)（这里把\(*\)看做\(0\)把）：
\[f_i = \sum_{j = 0}^{m - 1} A_j B_{i - m + 1 + j}\,(A_j - B_{i - m + 1 + j})^2\]
考虑翻转\(A\)，式子变成：
\[f_i = \sum_{j = 0}^{m - 1} A_{m - 1 - j}\,B_{i - m + 1 + j}\,(A_j - B_{i - m + 1 + j})^2\]
式子给人一种卷积的即视感……再化简一下，得到：
\[f_i = \sum_{j = 0}^{m - 1} A_{m - 1 - j}^3\,B_{i - m + 1 + j}\,- 2A_{m - 1 - j}^2\,B_{i - m + 1 + j}^2\,+ A_{m - 1 - j}\,B_{i - m + 1 + j}^3\]
这是三个卷积……跑三遍FFT即可得到\(f_i\)，当且仅当\(f_i = 0\)时，以\(i\)结尾的长度为\(m\)的子串才能匹配。
顺便提个几个坑……这个题用long double + std::complex会因为常数巨大而T掉。并且注意eps不要调那么细啊……这个东西精度误差极大，我直接调成1了。
代码：
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cctype>
#include <cmath>
#include <algorithm>
#include <utility>
#include <queue>
#include <complex>
constexpr int maxn = 300005;
using R = double;
using C = std::complex<R>;
constexpr R pi = acos(-1);
inline int flip(int bi, int x) {
int ans = 0;
for(int i = 0; i < bi; i ++) {
if(x & (1 << i)) {
ans += (1 << (bi - i - 1));
}
}
return ans;
}
inline void FFT(C *A, int bi, int n, R flag = 1.00) {
for(int i = 0; i < n; i ++) {
int fl = flip(bi, i);
if(i < fl) std::swap(A[fl], A[i]);
}
for(int L = 1; L < n; L <<= 1) {
R rd = pi / (R(L));
C xi(cos(rd), sin(flag * rd));
for(int i = 0; i < n; i += (L << 1)) {
C w(1.00, 0.00);
for(int j = i; j < i + L; j ++) {
C v1 = A[j], v2 = A[j + L];
A[j] = v1 + v2 * w;
A[j + L] = v1 - v2 * w;
w = w * xi;
}
}
}
}
constexpr R eps = 1.00;
inline int sign(R x) {
if(fabs(x) < eps) {
return 0;
} else {
if(x < 0.00) {
return -1;
} else {
return 1;
}
}
}
int main() {
static char A[maxn], B[maxn];
static int V1[maxn << 2], V2[maxn << 2];
static C P1[maxn << 2], P2[maxn << 2], P3[maxn << 2];
static R F[maxn << 2];
int m, n; scanf("%d%d", &m, &n);
if(m > n) {
puts("0"); return 0;
}
scanf("%s%s", A, B);
std::reverse(A, A + m);
for(int i = 0; i < m; i ++) {
V1[i] = (A[i] == '*') ? 0 : (A[i] - 'a' + 1);
}
for(int i = 0; i < n; i ++) {
V2[i] = (B[i] == '*') ? 0 : (B[i] - 'a' + 1);
}
int sz = 1, bi = 0;
while(sz <= (m + n)) {
sz *= 2; bi ++;
}
std::fill(P1, P1 + sz, C(0, 0));
std::fill(P2, P2 + sz, C(0, 0));
for(int i = 0; i < m; i ++) {
P1[i] = V1[i] * V1[i] * V1[i];
}
for(int i = 0; i < n; i ++) {
P2[i] = V2[i];
}
FFT(P1, bi, sz); FFT(P2, bi, sz);
for(int i = 0; i < sz; i ++) {
P3[i] = P1[i] * P2[i];
}
FFT(P3, bi, sz, -1.00);
for(int i = 0; i < sz; i ++) {
F[i] += P3[i].real() / (R(sz));
}
std::fill(P1, P1 + sz, C(0, 0));
std::fill(P2, P2 + sz, C(0, 0));
for(int i = 0; i < m; i ++) {
P1[i] = V1[i] * V1[i];
}
for(int i = 0; i < n; i ++) {
P2[i] = V2[i] * V2[i];
}
FFT(P1, bi, sz); FFT(P2, bi, sz);
for(int i = 0; i < sz; i ++) {
P3[i] = P1[i] * P2[i];
}
FFT(P3, bi, sz, -1.00);
for(int i = 0; i < sz; i ++) {
F[i] += (R(-2)) * P3[i].real() / (R(sz));
}
std::fill(P1, P1 + sz, C(0, 0));
std::fill(P2, P2 + sz, C(0, 0));
for(int i = 0; i < m; i ++) {
P1[i] = V1[i];
}
for(int i = 0; i < n; i ++) {
P2[i] = V2[i] * V2[i] * V2[i];
}
FFT(P1, bi, sz); FFT(P2, bi, sz);
for(int i = 0; i < sz; i ++) {
P3[i] = P1[i] * P2[i];
}
FFT(P3, bi, sz, -1.00);
for(int i = 0; i < sz; i ++) {
F[i] += P3[i].real() / (R(sz));
}
static int ans[maxn]; int cnt = 0;
for(int i = 0; i < n; i ++) {
#ifdef LOCAL
printf("F[%d] : %.10Lf\n", i, F[i]);
#endif
if(sign(F[i]) == 0 && i >= m - 1) {
ans[cnt ++] = i - m + 1 + 1;
}
}
printf("%d\n", cnt);
for(int i = 0; i < cnt; i ++) {
if(i) putchar(' ');
printf("%d", ans[i]);
}
putchar('\n');
return 0;
}
[BZOJ 3622]已经没有什么好害怕的了
已经……没有什么好害怕的了
OI给了我足够的欢乐了，OI让我看到了新的世界，让我有了人生中最好的几位朋友，我没有什么可以遗憾的了……就算可能要面临最残酷的结局吧……
有点跑题了呢……还是说这个题吧
显然可以看出来这个题要求你的从大到小的匹配（姑且称之为好的匹配）恰好要有\(\frac{n + k}{2}\)个，剩下的全都不是好的匹配。
首先把糖果（记为\(A\)）和药片（记为\(B\)）分别排序，对于每一个\(A\)中元素就很容易得到\(B\)中有多少比它小的。考虑设计状态\(f[i][j]\)表示对（排序后的）\(A\)的前\(i\)项，已经确定了的好的匹配有至少\(j\)个。这个DP转移显然。
然后发现\(f[n][j]\times (n - j)!\)就是对于所有元素的完美匹配中，有至少\(j\)个是好的匹配的方案数。这正是广义容斥原理的用武之地。考虑记\(t_j = f[n][j]\times (n - j)!\)，那么答案就是（令\(\frac{n + k}{2} = l\)）：
\[\sum_{i = l}^n (-1)^{n - i}\binom{i}{l} t_i\]
代码：
/**************************************************************
Problem: 3622
User: danihao123
Language: C++
Result: Accepted
Time:1844 ms
Memory:63680 kb
****************************************************************/
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cctype>
#include <algorithm>
#include <utility>
typedef long long ll;
const ll ha = 1000000009LL;
const int maxn = 2005;
ll C[maxn][maxn], F[maxn];
inline void process() {
int N = 2000;
C[0][0] = 1LL;
for(int i = 1; i <= N; i ++) {
C[i][0] = C[i][i] = 1LL;
for(int j = 1; j < i; j ++) {
C[i][j] = (C[i - 1][j] + C[i - 1][j - 1]) % ha;
}
}
F[0] = 1LL;
for(int i = 1; i <= N; i ++) {
F[i] = (F[i - 1] * (ll(i))) % ha;
}
}
int A[maxn], B[maxn];
ll cp[maxn]; int n;
void gen_cp() {
std::sort(A + 1, A + 1 + n); std::sort(B + 1, B + 1 + n);
for(int i = 1; i <= n; i ++) {
cp[i] = n - (B + n + 1 - std::upper_bound(B + 1, B + 1 + n, A[i]));
#ifdef LOCAL
printf("cp[%d] : %lld\n", i, cp[i]);
#endif
}
}
ll f[maxn][maxn];
ll dp(int dec) {
process(); gen_cp();
f[0][0] = 1LL;
for(int i = 1; i <= n; i ++) {
for(int j = 0; j <= n; j ++) {
f[i][j] = f[i - 1][j];
if(j > 0) {
f[i][j] += (f[i - 1][j - 1] * std::max(0LL, cp[i] - ll(j - 1))) % ha;
if(f[i][j] >= ha) f[i][j] -= ha;
}
#ifdef LOCAL
printf("f[%d][%d] : %lld\n", i, j, f[i][j]);
#endif
}
}
for(int j = 0; j <= n; j ++) {
f[n][j] = (f[n][j] * F[n - j]) % ha;
}
ll ans = 0;
for(int j = dec; j <= n; j ++) {
ll a = ((j - dec) % 2 == 0) ? 1LL : (ha - 1LL);
ll delta = (C[j][dec] * f[n][j]) % ha;
ans += (a * delta) % ha;
if(ans >= ha) ans -= ha;
}
return ans;
}
int main() {
int k; scanf("%d%d", &n, &k);
if((n - k) & 1) {
puts("0"); return 0;
}
for(int i = 1; i <= n; i ++) {
scanf("%d", &A[i]);
}
for(int i = 1; i <= n; i ++) {
scanf("%d", &B[i]);
}
printf("%lld\n", dp((n + k) / 2));
return 0;
}
[HDU 4625]JZPTREE
第一次用到第二类斯特林数……
首先众所周知第二类斯特林数有一个很好康的式子：
\[x^k = \sum_{i = 0}^k S(k, i) x^{\underline{i}}\]
然后带进去试试？一番化简之后得到：
\[E_u = \sum_{i = 0}^k S(k, i) \sum_{v = 1}^n dist(u, v)^{\underline{i}}\]
那个下降幂很恶心……把它强行搞成组合数后发现：
\[E_u = \sum_{i = 0}^k S(k, i) i! \sum_{v = 1}^n \binom{dist(u, v)}{i}\]
然后考虑每个点预处理后面那个和式……显然可以树形DP苟。用经典的树形DP套路，设两个状态分别表示一个点的子树内部的贡献以及外面的点给他的贡献。至于转移，由于那是个组合数，所以可以考虑用Pascal定理苟……
代码：
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cctype>
#include <cmath>
#include <algorithm>
#include <utility>
const int maxn = 50005;
const int maxk = 505;
const int ha = 10007;
int first[maxn];
int next[maxn << 1], to[maxn << 1];
int cnt;
void add_edge(int u, int v) {
cnt ++;
next[cnt] = first[u]; first[u] = cnt;
to[cnt] = v;
}
int n, k;
int down[maxn][maxk];
void dfs_1(int x, int fa = 0) {
down[x][0] = 1;
for(int i = first[x]; i; i = next[i]) {
int v = to[i];
if(v != fa) {
dfs_1(v, x);
down[x][0] = (down[x][0] + down[v][0]) % ha;
for(int j = 1; j <= k; j ++) {
down[x][j] = (down[x][j] + down[v][j] + down[v][j - 1]) % ha;
}
}
}
}
int up[maxn][maxk];
void dfs_2(int x, int fa = 0) {
if(fa) {
up[x][0] = n - down[x][0];
for(int j = 1; j <= k; j ++) {
up[x][j] = (up[fa][j] + up[fa][j - 1] + down[fa][j] + down[fa][j - 1]) % ha;
up[x][j] = (up[x][j] - (2 * down[x][j - 1]) % ha + ha) % ha;
up[x][j] = (up[x][j] - down[x][j] + ha) % ha;
if(j > 1) up[x][j] = (up[x][j] - down[x][j - 2] + ha) % ha;
}
}
for(int i = first[x]; i; i = next[i]) {
int v = to[i];
if(v != fa) {
dfs_2(v, x);
}
}
}
int S[maxk][maxk];
int fac[maxn];
void process() {
S[1][1] = 1;
for(int i = 2; i <= k; i ++) {
S[i][0] = 0;
for(int j = 1; j < i; j ++) {
S[i][j] = ((j * S[i - 1][j]) % ha + S[i - 1][j - 1]) % ha;
}
S[i][i] = 1;
}
fac[0] = 1;
for(int i = 1; i <= n; i ++) {
fac[i] = (fac[i - 1] * (i % ha)) % ha;
}
}
int main() {
int T; scanf("%d", &T);
n = 50000; k = 500; process();
while(T --) {
scanf("%d%d", &n, &k);
cnt = 0; memset(first, 0, sizeof(first));
for(int i = 0; i < n - 1; i ++) {
int u, v; scanf("%d%d", &u, &v);
add_edge(u, v); add_edge(v, u);
}
memset(down, 0, sizeof(down));
memset(up, 0, sizeof(up));
dfs_1(1); dfs_2(1);
for(int i = 1; i <= n; i ++) {
int ans = 0;
for(int j = 1; j <= k; j ++) {
int f_v = (down[i][j] + up[i][j]) % ha;
ans = (ans + (S[k][j] * ((f_v * fac[j]) % ha)) % ha) % ha;
}
printf("%d\n", ans);
}
}
return 0;
}
