2021-07-09
查拉图斯特拉如是说
Description
Solution
有多项式快速多点求值的算法，我写的是另外一种更加简单的 (Theta(nlog n)) 的算法。
我们可以发现，我们就是要求：
[sum_{k=0}^{n} inom{n}{k}sum_{i=0}^{m} a_i imes k^i
]
[sum_{i=0}^{m} a_i imes sum_{k=0}^{n} inom{n}{k} imes k^i
]
也就是说，我们只需要快速知道：
[F(x)=sum_{i=0}^{infty} frac{x^i}{i!}sum_{k=0}^{n} inom{n}{k} imes k^i
]
[=sum_{k=0}^{n} inom{n}{k} sum_{i=0}^{infty} frac{(xk)^i}{i!}
]
[=sum_{k=0}^{n} inom{n}{k} (e^x)^k
]
[=(e^x+1)^n
]
这个直接多项式 (ln, ext{exp}) 即可。
Code
#include <bits/stdc++.h>
using namespace std;
#define Int register int
#define MAXN 800005
template <typename T> inline void read (T &t){t = 0;char c = getchar();int f = 1;while (c < '0' || c > '9'){if (c == '-') f = -f;c = getchar();}while (c >= '0' && c <= '9'){t = (t << 3) + (t << 1) + c - '0';c = getchar();} t *= f;}
template <typename T,typename ... Args> inline void read (T &t,Args&... args){read (t);read (args...);}
template <typename T> inline void write (T x){if (x < 0){x = -x;putchar ('-');}if (x > 9) write (x / 10);putchar (x % 10 + '0');}
#define mod 998244353
#define Gi 332748118
#define G 3
int mul (int a,int b){return 1ll * a * b % mod;}
int dec (int a,int b){return a >= b ? a - b : a + mod - b;}
int add (int a,int b){return a + b >= mod ? a + b - mod : a + b;}
int qkpow (int a,int b){
b %= mod - 1;
if (b < 0) b += mod - 1;
int res = 1;for (;b;b >>= 1,a = mul (a,a)) if (b & 1) res = mul (res,a);
return res;
}
void Add (int &a,int b){a = add (a,b);}
typedef vector<int> poly;
int w[MAXN],rev[MAXN];
#define SZ(A) ((int)A.size())
void init_ntt (){
int lim = 1 << 18;
for (Int i = 0;i < lim;++ i) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << 17);
int Wn = qkpow (3,(mod - 1) / lim);w[lim >> 1] = 1;
for (Int i = lim / 2 + 1;i < lim;++ i) w[i] = mul (w[i - 1],Wn);
for (Int i = lim / 2 - 1;i;-- i) w[i] = w[i << 1];
}
void ntt (poly &a,int type){
#define G 3
#define Gi 332748118
static int d[MAXN];int lim = a.size();
for (Int i = 0,z = 18 - __builtin_ctz(lim);i < lim;++ i) d[rev[i] >> z] = a[i];
for (Int i = 1;i < lim;i <<= 1)
for (Int j = 0;j < lim;j += i << 1)
for (Int k = 0;k < i;++ k){
int x = mul (w[i + k],d[i + j + k]);
d[i + j + k] = dec (d[j + k],x),d[j + k] = add (d[j + k],x);
}
for (Int i = 0;i < lim;++ i) a[i] = d[i] % mod;
if (type == -1){
reverse (a.begin() + 1,a.begin() + lim);
for (Int i = 0,Inv = qkpow (lim,mod - 2);i < lim;++ i) a[i] = mul (a[i],Inv);
}
#undef G
#undef Gi
}
poly operator * (poly A,poly B){
int lim = 1,l = 0,len = SZ(A) + SZ(B) - 1;
while (lim < SZ(A) + SZ(B)) lim <<= 1,++ l;
A.resize (lim),B.resize (lim),ntt (A,1),ntt (B,1);
for (Int i = 0;i < lim;++ i) A[i] = mul (A[i],B[i]);
ntt (A,-1),A.resize (len);
return A;
}
poly operator + (poly A,poly B){
int len = max (SZ(A),SZ(B));A.resize (len),B.resize (len);
for (Int i = 0;i < len;++ i) Add (A[i],B[i]);
return A;
}
poly operator - (poly A,poly B){
int len = max (SZ(A),SZ(B));A.resize (len),B.resize (len);
for (Int i = 0;i < len;++ i) A[i] = dec (A[i],B[i]);
return A;
}
int sqri,Ansf;
struct node{
int real,imag;
node (int _real,int _imag) : real (_real) , imag (_imag) {}
node operator * (node b){return node (add (1ll * real * b.real % mod,1ll * imag * b.imag % mod * sqri % mod),add (1ll * real * b.imag % mod,1ll * b.real * imag % mod));}
node operator + (node b){return node (add (real,b.real),add (imag,b.imag));}
node operator ^ (int b){node a = *this,res = node (1,0);for (;b;b >>= 1,a = a * a) if (b & 1) res = res * a;return res;}
};
bool check_if_qs (int n){//判断是否是一个二次剩余
return qkpow (n,(mod - 1) >> 1) == 1;
}
bool tryit (int n){
if (n == 0){
puts ("0");
return 1;
}
int a = rand ();sqri = dec (mul (a,a),n);
if (!a || !check_if_qs (sqri)){
int x0 = (node (a,1) ^ (mod + 1 >> 1)).real;
Ansf = min (x0,mod - x0);
return 1;
}
else return 0;
}
int Work (int n){
while (1) if (tryit (n)) break;
return Ansf;
}
poly operator * (poly A,int B){
for (Int i = 0;i < SZ(A);++ i) A[i] = mul (A[i],B);
return A;
}
poly operator * (int B,poly A){
for (Int i = 0;i < SZ(A);++ i) A[i] = mul (A[i],B);
return A;
}
poly inv (poly A,int n){
if (n == 1) return poly(1,qkpow (A[0],mod - 2));
poly F1,F0 = inv (A,(n + 1) >> 1);poly tmp;tmp.resize (n);for (Int i = 0;i < n;++ i) tmp[i] = A[i];
F1 = F0 * 2 - F0 * F0 * A,F1.resize (n);
return F1;
}
poly inv(poly A){return inv(A,SZ(A));}
#define inv2 (mod+1)/2
poly Sqrt(poly A,int n){
if (n == 1) return poly(1,Work(A[0]));
poly F1,F0 = Sqrt(A,(n + 1) >> 1);poly tmp;tmp.resize (n);for (Int i = 0;i < n;++ i) tmp[i] = A[i];
F1 = (A * inv(F0,n) + F0) * (inv2),F1.resize (n);
return F1;
}
poly Sqrt (poly A){return Sqrt(A,SZ(A));}
poly der (poly a){
for (Int i = 0;i < SZ (a) - 1;++ i) a[i] = mul (a[i + 1],i + 1);
a.pop_back ();return a;
}
poly ine (poly a){
a.push_back (0);
for (Int i = SZ (a) - 1;i;-- i) a[i] = mul (a[i - 1],qkpow (i,mod - 2));
a[0] = 0;return a;
}
poly ln (poly A,int n){
A = ine (der(A) * inv(A,n)),A.resize (n);
return A;
}
poly ln (poly A){return ln(A,SZ(A));}
void Putout (poly A){
cout << SZ(A) << ": ";
for (Int i = 0;i < SZ(A);++ i) cout << A[i] << " ";cout << endl;
}
poly exp (poly A,int n){
if (n == 1) return poly (1,1);
poly F1,F0 = exp (A,n + 1 >> 1);poly tmp;tmp.resize (n);for (Int i = 0;i < n;++ i) tmp[i] = A[i];
F1 = (poly(1,1) - ln(F0,n) + tmp) * F0,F1.resize (n);
return F1;
}
poly exp (poly A){return exp(A,SZ(A));}
poly A;
int n,m,a[MAXN],fac[MAXN],ifac[MAXN];
signed main(){
freopen ("number.in","r",stdin);
freopen ("number.out","w",stdout);
init_ntt(),read (n,m),A.resize (m + 1);
fac[0] = 1;for (Int i = 1;i <= m;++ i) fac[i] = mul (fac[i - 1],i);
ifac[m] = qkpow (fac[m],mod - 2);for (Int i = m;i;-- i) ifac[i - 1] = mul (ifac[i],i);
for (Int i = 0;i <= m;++ i) read (a[i]),A[i] = mul (inv2,ifac[i] + (i == 0));
A = ln (A);for (Int i = 0;i <= m;++ i) A[i] = mul (A[i],n);A = exp (A);
for (Int i = 0,v = qkpow (2,n);i <= m;++ i) A[i] = mul (A[i],v);
int ans = 0;for (Int i = 0;i <= m;++ i) Add (ans,mul (a[i],mul (fac[i],A[i])));
write (ans),putchar ('
');
return 0;
}
橡树上的逃亡
Description
Solution
假设设 (s_x) 表示树上以 (x) 为根的子树内 ([L,R]) 的叶子数量，那么一个点 (u) 的父亲边的贡献就是：
[1-(frac{s_u}{s_1})^k-(frac{s_1-s_u}{s_1})^k
]
下面这个 (s_1) 可以提出来，也就是说我们只需要求 (sum s_u^k+sum (s_1-s_u)^k) 。
然后你发现这个 (kle 10)，你就可以把后面这个二项式展开，也就是只需要求 (sum s_u^i) 。
这个可以树上做个前缀和和 ( ext{dfs}) 序上做前缀和解决。
复杂度 (Theta(nk^2))。
Code
#include <bits/stdc++.h>
using namespace std;
#define Int register int
#define mod 998244353
#define MAXN 200005
#define MAXM 15
template <typename T> inline void read (T &t){t = 0;char c = getchar();int f = 1;while (c < '0' || c > '9'){if (c == '-') f = -f;c = getchar();}while (c >= '0' && c <= '9'){t = (t << 3) + (t << 1) + c - '0';c = getchar();} t *= f;}
template <typename T,typename ... Args> inline void read (T &t,Args&... args){read (t);read (args...);}
template <typename T> inline void write (T x){if (x < 0){x = -x;putchar ('-');}if (x > 9) write (x / 10);putchar (x % 10 + '0');}
template <typename T> void chkmax (T &a,T b){a = max (a,b);}
template <typename T> void chkmin (T &a,T b){a = min (a,b);}
int mul (int a,int b){return 1ll * a * b % mod;}
int dec (int a,int b){return a >= b ? a - b : a + mod - b;}
int add (int a,int b){return a + b >= mod ? a + b - mod : a + b;}
int qkpow (int a,int b){
int res = 1;for (;b;b >>= 1,a = mul (a,a)) if (b & 1) res = mul (res,a);
return res;
}
void Add (int &a,int b){a = add (a,b);}
void Sub (int &a,int b){a = dec (a,b);}
int type,n,q;
vector <int> g[MAXN],lef;
int C[MAXM][MAXM],fa[MAXN],tot[MAXN],siz[MAXN],tur[MAXN],rig[MAXN],lft[MAXN],dep[MAXN],st[MAXN][21],sum[MAXN][MAXM],lst1[MAXN][MAXM],lst2[MAXN][MAXM],lst3[MAXN][MAXM];
void dfs (int u){
tot[u] = 1,lft[u] = n + 1;
for (Int i = 1;i <= 20;++ i) st[u][i] = st[st[u][i - 1]][i - 1];
if (!g[u].size()) siz[u] = 1,rig[u] = u,lft[u] = u,lef.push_back (u),tur[u] = lef.size();
for (Int v : g[u]) dep[v] = dep[u] + 1,st[v][0] = u,dfs (v),tot[u] += tot[v],siz[u] += siz[v],chkmax (rig[u],rig[v]),chkmin (lft[u],lft[v]);
if (u > 1) for (Int i = 0,v = 1;i <= 10;++ i) sum[u][i] = v,v = mul (v,siz[u]);
}
int getlca (int u,int v){
if (dep[u] < dep[v]) swap (u,v);
for (Int i = 20,dis = dep[u] - dep[v];~i;-- i) if (dis >> i & 1) u = st[u][i];
if (u == v) return u;
else{
for (Int i = 20;~i;-- i) if (st[u][i] ^ st[v][i]) u = st[u][i],v = st[v][i];
return st[u][0];
}
}
int getit (int lp,int rp,int lca,int K){
int con1 = dec (dec (sum[rp][K],sum[lp][K]),dec (lst1[rp][K],lst1[lca][K])),con2 = 0;
for (Int i = 0,v = 1;i <= K;++ i,v = mul (v,tur[lp] - 1)){
int fuc = dec (lst2[lp][K - i],lst2[lca][K - i]),tmp = mul (fuc,mul (C[K][i],v));
if (i & 1) Sub (con2,tmp);else Add (con2,tmp);
}
for (Int i = 0,v = 1;i <= K;++ i,v = mul (v,tur[rp] + 1)){
int fuc = dec (lst3[rp][K - i],lst3[lca][K - i]),tmp = mul (fuc,mul (C[K][i],v));
if (K - i & 1) Sub (con2,tmp);else Add (con2,tmp);
}
return add (con1,con2);
}
signed main(){
freopen ("tree.in","r",stdin);
freopen ("tree.out","w",stdout);
read (type,n,q);
for (Int i = 2;i <= n;++ i) read (fa[i]),g[fa[i]].push_back (i);
dfs (1);C[0][0] = 1;
for (Int i = 1;i <= 10;++ i){
C[i][0] = 1;
for (Int j = 1;j <= i;++ j) C[i][j] = add (C[i - 1][j],C[i - 1][j - 1]);
}
for (Int u = 2;u <= n;++ u){
for (Int i = 0;i <= 10;++ i) Add (sum[u][i],sum[u - 1][i]);
for (Int i = 0,v1 = 1,v2 = 1,v3 = 1;i <= 10;++ i,v1 = mul (v1,siz[u]),v2 = mul (v2,tur[rig[u]]),v3 = mul (v3,tur[lft[u]]))
lst1[u][i] = add (lst1[fa[u]][i],v1),lst2[u][i] = add (lst2[fa[u]][i],v2),lst3[u][i] = add (lst3[fa[u]][i],v3);
}
int lstans = 0;
while (q --> 0){
int L,R,K;read (L,R,K);
L ^= (type * lstans),R ^= (type * lstans);
if (lef[0] > R || lef[lef.size() - 1] < L) write (lstans = 0),putchar ('
');
else{
int lp = lower_bound (lef.begin(),lef.end(),L) - lef.begin(),rp = upper_bound (lef.begin(),lef.end(),R) - lef.begin() - 1,lca,all = rp - lp + 1;lp = lef[lp],rp = lef[rp],lca = getlca (lp,rp);
int iv = qkpow (all,mod - 2);iv = qkpow (iv,K);
if (all <= 1) write (lstans = 0),putchar ('
');
else{
int tmp = getit (lp,rp,lca,K);
for (Int i = 0,v = 1;i <= K;++ i,v = mul (v,all)){
int fuc = mul (C[K][i],mul (v,getit (lp,rp,lca,K - i)));
if (K - i & 1) Sub (tmp,fuc);else Add (tmp,fuc);
}
int nsiz = rp - lp + dep[lp] - dep[lca];
int ans = dec (n - 1,mul (iv,add (tmp,mul (n - 1 - nsiz,qkpow (all,K)))));
write (lstans = ans),putchar ('
');
}
}
}
return 0;
}
2021-07-11
悄悄话
Description
Solution
水题
数学考试
Description
Solution
考虑没有限制的时候，因为 (L_i,R_i) 很小，所以我们可以网络流，即原点往汇点拉一条链，每条边的边权为 (f(i,j),jin [L_i,R_i])，然后做最小割。
考虑限制，因为 (x_uge x_v+DRightarrow x_vle x_u-D)，也就是说，我们可以 ((u,i)) 往 ((v,i-D)) 连一条边权为 (infty) 的边。
判断是否有解是否判断是否割掉了 (infty) 。
Code
#include <bits/stdc++.h>
using namespace std;
#define Int register int
#define int long long
#define MAXN 205
template <typename T> inline void read (T &t){t = 0;char c = getchar();int f = 1;while (c < '0' || c > '9'){if (c != ' ' && c != '
') f = -f;c = getchar();}while (c >= '0' && c <= '9'){t = (t << 3) + (t << 1) + c - '0';c = getchar();} t *= f;}
template <typename T,typename ... Args> inline void read (T &t,Args&... args){read (t);read (args...);}
template <typename T> inline void write (T x){if (x < 0){x = -x;putchar ('-');}if (x > 9) write (x / 10);putchar (x % 10 + '0');}
template <typename T> void chkmax (T &a,T b){a = max (a,b);}
template <typename T> void chkmin (T &a,T b){a = min (a,b);}
int T,n,m,ind[MAXN][MAXN];
struct node{
int a,b,c,d,L,R;
}s[MAXN];
int f (int i,int x){return s[i].a * x * x * x + s[i].b * x * x + s[i].c * x + s[i].d;}
#define MAXM 200005
#define MAXX 40005
struct edge{
int v,w,nxt;
}e[MAXM << 1];
int top = 1,tot;
int head[MAXX];
void link (int u,int v,int w){
e[++ top] = edge {v,w,head[u]},head[u] = top;
e[++ top] = edge {u,0,head[v]},head[v] = top;
}
int dep[MAXX];
bool BFS (int S,int T){
queue <int> q;q.push (S);
for (Int i = 0;i <= tot;++ i) dep[i] = -1;dep[S] = 0;
while (!q.empty()){
int u = q.front();q.pop ();
for (Int i = head[u];i;i = e[i].nxt){
int v = e[i].v,w = e[i].w;
if (w && dep[v] == -1) dep[v] = dep[u] + 1,q.push (v);
}
}
return dep[T] != -1;
}
int dfs (int s,int t,int flow){
if (s == t) return flow;int ans = 0;
for (Int i = head[s];i;i = e[i].nxt){
if (flow == 0) break;
int v = e[i].v,w = e[i].w;
if (dep[v] == dep[s] + 1 && w){
int del = dfs (v,t,min (flow,w));
e[i].w -= del,e[i ^ 1].w += del;
ans += del,flow -= del;
}
}
if (ans == 0) dep[s] = -1;
return ans;
}
#define inf 1e15
#define mxv 1e8
int MaxFlow (int s,int t){
int ans = 0;
while (BFS (s,t)) ans += dfs (s,t,inf);
return ans;
}
int led[MAXN];
void Work (){
read (n,m),top = 1,memset (head,0,sizeof (head));
for (Int i = 1;i <= n;++ i) read (s[i].a,s[i].b,s[i].c,s[i].d,s[i].L,s[i].R);
int S = 0,T = 1;tot = 1;
for (Int i = 1;i <= n;++ i){
link (S,led[i] = ++ tot,inf);
for (Int j = s[i].L;j <= s[i].R;++ j,++ tot) link (tot,tot + 1,mxv - f (i,j));
link (tot,T,inf);
}
for (Int i = 1;i <= m;++ i){
int u,v,D;read (u,v,D);
for (Int j = s[u].L;j <= s[u].R;++ j)
if (s[v].L <= j - D && j - D <= s[v].R)
link (led[u] + j - s[u].L,led[v] + j - D - s[v].L,inf);
else if (s[v].R < j - D) link (led[u] + j - s[u].L,T,inf);
}
int tmp = MaxFlow (S,T);
if (tmp > inf) puts ("mei ji ge");
else write ((int)mxv * n - tmp),putchar ('
');
}
signed main(){
freopen ("sleep.in","r",stdin);
freopen ("sleep.out","w",stdout);
int T;read (T);
while (T --> 0) Work ();
return 0;
}