orangejuice's blog

挂了请务必叫我

「GDOI2020模拟赛day1」Permutation

为了便于叙述,设原题中的 \(n\)\(N\)

题目分析

要求一个 \(1-n\) 的环排列,看成是一个环遍历

发现每条边的权值限制了遍历过程中穿过这条边的次数

\(1\) 为根,强制从 \(1\) 开始遍历,考虑以一个个自由段(即未确定前后连接关系的子段)的形式维护 \(u\) 子树中的序列段

那么,只需要满足 \(u\) 子树中自由段的个数为 \(\frac{w(u,fa_u)}{2}\) (每一个自由段的两端均对应一次跨越)即可

分析即可发现 \(w(u,fa_u)\) 显然是 \(O(size_u)\) 级别的,因此考虑树形背包

那么我们就需要支持合并两组自由段

\[ \ \]

\(O(N^3-N^2\log N)\)

设当前 \(1\ldots n\) 个自由端,合并上 \(m\) 个自由段(从子树合并上来是恰好 \(m\) 个)

注意这些自由段之间是无序的

先考虑一个简单的模型:

\(n\) 个无序自由段拼接成 \(m\) 个无序自由段,设方案数为 \(W(n,m)\) ,则考虑

先把 \(n\) 个自由段排列,然后选出 \(n-m\) 个间隔连接,然后除掉得到的 \(m\) 个段之间的排列

得到 \(\begin{aligned} W(n,m)=\frac{n!}{m!}\binom{n-1}{n-m} \end{aligned}\)

类似的,可以把 \(n+m\) 个自由段排成一排,合并成若干段

但是显然存在的问题就是:可能在两组自由段之间形成了连接,这样的连接是非法的

因此考虑容斥

\(\begin{aligned} dp'_{d}\longleftarrow\sum_{i=1}^ndp_{u,i}\sum_{j=1}^idp_v (-1)^{i-j}W(i,j) \sum_{k=1}^m (-1)^{m-k}W(m-1,k)\sum_{d=1}^{j+k}W(j+k,d)\end{aligned}\)

将上式分解为四步转移

\(n,i\rightarrow j,O(n^2)\)

\(m\rightarrow k,O(m)\)

\(j,k\rightarrow j+k,O(nm)\)

\(j+k\rightarrow d,O((n+m)^2)\)

其中 \(O(n^2,(n+m)^2)\) 的部分如果用卷积优化,即可做到 \(O(N^2\log N)\)

但是 \(O(N^3)\)\(pts75\) 了...

Tips:注意在将 \(1\) 号节点加入序列时,用上面的方法无法保证它在序列首

需要特殊处理,始终强制它在第一个

\(O(N^2)\)

当我发现这个做法不需要任何优化就可以做到 \(O(N^2)\) 的时候。。。。

把这个容斥的过程爆开

先对于每个儿子的 \(dp\) 值按照容斥系数进行上文中 \(m\rightarrow k\) 的变换,复杂度为 \(O(size_u)\)

然后进行背包合并,由树形背包的复杂度分析,总复杂度为 \(O(N^2)\)

最后发现其实我们只需要知道 \(dp_{u,w(u,fa_u)}\) ,因此这里也只需要 \(O(size_u)\)

综上,复杂度为 \(O(N^2)\)

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
#include<bits/stdc++.h>
using namespace std;

typedef long long ll;
typedef unsigned long long ull;
typedef double db;
typedef pair <int,int> Pii;
#define reg register
#define mp make_pair
#define pb push_back
#define Mod1(x) ((x>=P)&&(x-=P))
#define Mod2(x) ((x<0)&&(x+=P))
#define rep(i,a,b) for(int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(int i=a,i##end=b;i>=i##end;--i)
template <class T> inline void cmin(T& a,T b){ ((a>b)&&(a=b)); }
template <class T> inline void cmax(T& a,T b){ ((a<b)&&(a=b)); }

const int DEBUG=1;
#define Log(...) (DEBUG&&(fprintf(stderr,__VA_ARGS__)))

char IO;
int rd(){
int s=0,f=0;
while(!isdigit(IO=getchar())) if(IO=='-') f=1;
do s=(s<<1)+(s<<3)+(IO^'0');
while(isdigit(IO=getchar()));
return f?-s:s;
}

bool Mbe;
const int N=5010;

int n,P;
ll qpow(ll x,ll k=P-2){
ll res=1;
for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P;
return res;
}
int I[N],J[N],C[N][N],W[N][N];

struct Edge{
int to,nxt,w;
} e[N<<1];
int head[N],ecnt;
void AddEdge(int u,int v,int w){
e[++ecnt]=(Edge){v,head[u],w};
head[u]=ecnt;
}

int sz[N],w[N],dp[N][N],T[N];
int A[N],B[N];

void dfs(int u,int f){
sz[u]=0,dp[u][0]=1;
for(int i=head[u];i;i=e[i].nxt){
int v=e[i].to;
if(v==f) continue;
w[v]=e[i].w,dfs(v,u);
rep(a,0,sz[u]) T[a]=dp[u][a],dp[u][a]=0;
rep(a,0,sz[u]) rep(b,1,sz[v]) dp[u][a+b]=(dp[u][a+b]+1ll*T[a]*dp[v][b])%P;
sz[u]+=sz[v];
}
int s=0;
rep(i,w[u],sz[u]+1) {
int e=W[i][w[u]];
if(u==1) e=1ll*C[i-1][i-w[u]]*J[i-1]%P*I[w[u]-1]%P;
s=(s+1ll*e*dp[u][i-1])%P;
}
// 求出我们需要的点值

rep(i,w[u]+1,sz[u]) dp[u][i]=0;
if(!s) puts("0"),exit(0);
// 容斥系数变换
rep(i,1,w[u]) {
dp[u][i]=1ll*s*W[w[u]][i]%P;
if((w[u]-i)&1) dp[u][i]=P-dp[u][i];
}
sz[u]=w[u];
}

bool Med;
int main(){
Log("Memory taken %.2lf\n",(&Med-&Mbe)/1024.0/1024.0);
n=rd(),P=rd();
rep(i,J[0]=1,n) J[i]=1ll*J[i-1]*i%P;
I[n]=qpow(J[n]);
drep(i,n,1) I[i-1]=1ll*I[i]*i%P;
rep(i,0,n) rep(j,C[i][0]=1,i) C[i][j]=C[i-1][j]+C[i-1][j-1],Mod1(C[i][j]),W[i][j]=1ll*C[i-1][i-j]*J[i]%P*I[j]%P;
rep(i,2,n){
int u=rd(),v=rd(),w=rd();
if(w&1) return puts("0");
AddEdge(u,v,w/=2),AddEdge(v,u,w);
}
dfs(w[1]=1,0);
dp[1][1]=1ll*dp[1][1]*n%P;
printf("%d\n",dp[1][1]);
}

「TJOI / HEOI2016」求和

题目大意:

\(\displaystyle \sum_{i=0}^n\sum_{j=0}^i \begin{Bmatrix}i\\ j\end{Bmatrix}2^j\cdot j!\)

由于第二类斯特林数的生成函数 \(S_m(x)=\cfrac{1}{m!}(e^x-1)^m\)

所以求的东西就是 \(\displaystyle F(x)=\sum_{i=0} (2e^x-2)^i=\frac{1}{3-2e^x}\)\(n\) 项系数。

Read more »

「lych_cys模拟题2018」橘子树

你可能需要一点点简单知识

你可能需要一点点预先知识

设值域上限为 \(M=2\cdot 10^5,T=3\cdot 10^5\)

\(F(n)\)

筛去 \(n^{\frac{1}{3}}\) 以内的质因数,此时 \(n\) 剩下的质因数 \(>n^{\frac{1}{3}}\)

剩下的 \(n\) 中,可能出现贡献的只有 \(n\) 是完全平方数的情况, \(O(1)\) 判断即可

复杂度为 \(O(n^{\frac{1}{3}})\)

\[ \ \]

求出 \(i\) 结点长了 \(t\) 时间的时候被收掉的答案

预处理 \(a_i,b_i\) 的答案 \(c_i=F(a_i),d_i=F(b_i)\) ,复杂度为 \(O(n\pi(M^{\frac{1}{3}}))\)

那么此时查询的权值就是 \(a_i\cdot t\) ,如果 \(>b_i\) 特判掉

否则,显然 \(c_i\) 的贡献可以先分离掉考虑进答案

对于剩下的部分,不妨设 \(x=\frac{a_i}{c_i}\) ,显然 \(x\) 不包含平方因子

考虑 \(x\)\(t\) 的合并所产生的贡献,实际上就是 \(\gcd(x,t)\) (让 \(x\) 的单个因子与 \(t\) 匹配一下)

然后对于 \(t\) 剩下的部分 \(\frac{t}{\gcd(x,t)}\) ,由于 \(t\leq 3\cdot 10^5\) ,可以预处理一下答案

综合一下上面的部分,那么一个点的答案就是 \(c_i\gcd(x,t)^2 F(\frac{x}{\gcd(t,x)})\)

因此查询一个点的复杂度就是 \(\gcd\)\(O(\log T)\)

\[ \ \]

综合利用上面的方法,勉强可以拿到 \(nm\log T\) 的20分

剩下的部分,可以考虑把树树剖一下,然后每个查询可以化为在 \(\text{dfs}\) 序上的若干更新区间 \(L_i,R_i,t_i\)

由于题目考虑的实际上是查询总和,因此可以对于每个点计算答案

不妨从 \(1-n\) 扫描 \(\text{dfs}\) 序,对于每个区间在 \(L_i\) 插入 \(t_i\) ,在 \(R_i+1\) 删除 \(t_i\)

用一个set来维护 \(t_i\) 的顺序,复杂度为 \(O(n\log ^2n)\)

在插入/删除set元素的同时,维护每一个 \(\Delta_i=t_i-t_{i-1}\) ,也就是我们要查询的每个数

查询单点时,由于这样的 \(\sum \Delta_i\leq T\) ,所以最多包含 \(O(\sqrt T)\) 种不同的 \(\Delta_i\)

用另一个set之类的东西维护这些不同的位置,然后暴力查询,复杂度为 \(O(n\sqrt T\log T)\)

总复杂度为 \(O(n\pi(M^{\frac{1}{3}})+n\sqrt T\log T)\) 左右

\[ \ \]

\[ \ \]

\[ \ \]


后记

由于特(智)殊(力)的(缺)原(陷)因

我考场上写了一个树剖+分块+莫比乌斯反演的做法,而且还没调出来,代码是这个,复杂度大概差不多 \(O(n\sqrt n \log_n^{\frac{3}{2}})\)

不提了。。。

「余姚中学 2019 联测 Day 6」解码

先不考虑求 \(p,q\)

根据人人都知道的欧拉定理 \(x^c\equiv x^{c\mod \varphi(n)} (\mod n)\)

那么 \(\varphi(n)=(p-1)(q-1)\) ,而 \((c,\varphi(n))=1\)

所以求出 \(\frac{1}{c} \pmod {\varphi(n)}\)

带入原来的式子 \((x^c)^{\frac{1}{c}\pmod {\varphi(n)}}\equiv x^1\pmod n\)

\(x=m^{\frac{1}{c}\pmod {\varphi(n)}}\pmod n\)

求逆最好用扩展欧几里得算法,复杂度为 \(O(\log n)\)

那么直接快速幂即可,但是如果快速幂套快速乘复杂度为 \(O(\log ^2)\) ,实际常数极大,很有可能超时(如果用long double O(1)快速乘另谈。。。)

由于知道 \(p,q\) 可以分别求出 \(m^{\frac{1}{c}\pmod {\varphi(n)}}\pmod p\)\(m^{\frac{1}{c}\pmod {\varphi(n)}}\pmod q\)

然后中国剩余定理合并,即 \(O(\log n)\)


现在问题是求 \(p,q\)

对于前3档分,由于素数密度是 \(O(\log n)\) 的,所以 \(\sqrt n -p\) 期望只有 \(\log n\)

而对于最后一档分,考虑更好表示,枚举 \(p+q\) 解出答案,发现

\(4n\leq (p+q)^2= (q-p)^2+4pq\leq \lambda^2+4n\)

\(2\sqrt{n}\leq (p+q)\le \sqrt{\lambda^2+4n}\)

由于 \(n\ge p^2\ge 10^{18},\lambda \leq 3\cdot 10^5\) ,这个范围实际非常非常非常非常小,大概只有 \(22\)

tips:实际上 \(4n\) 可能爆long long

枚举 \(p+q\) 之后, \(O(1)\) 解出 \(p,q\) 即可

竟然有人问怎么解 \(p,q\) ,我震惊了

\(q-p=\sqrt{(p+q)^2-4n}\) ,然后判一下是不是整数就好了

写的时候害怕sqrt炸精度,很多奇怪的语句请忽略

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
#include<bits/stdc++.h>
using namespace std;

#define reg register
typedef long long ll;
typedef unsigned long long ull;
#define rep(i,a,b) for(int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(int i=a,i##end=b;i>=i##end;--i)

#define pb push_back
#define Mod1(x) ((x>=P)&&(x-=P))
#define Mod2(x) ((x<0)&&(x+=P))

template <class T> inline void cmin(T &a,T b){ if(a>b) a=b; }
template <class T> inline void cmax(T &a,T b){ if(a<b) a=b; }

char IO;
template<class T=int> T rd(){
T s=0; int f=0;
while(!isdigit(IO=getchar())) if(IO=='-') f=1;
do s=(s<<1)+(s<<3)+(IO^'0');
while(isdigit(IO=getchar()));
return f?-s:s;
}

const int N=3e5+10;

ll n,m,c;

ll qmul(ll x,ll k,ll P){
k=(k%P+P)%P;
ll res=0;
for(;k;k>>=1,x=(x+x)%P) if(k&1) res=(res+x)%P;
return res;
}

ll qpow(ll x,ll k,ll P) {
ll res=1;
for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P;
return res;
}

void Exgcd(ll a,ll b,ll &x,ll &y){
if(b==0) x=1,y=0;
else Exgcd(b,a%b,y,x),y-=a/b*x;
}

ll Inv(ll a,ll P){
ll x,y;
Exgcd(a,P,x,y);
return (x%P+P)%P;
}

int main(){
freopen("rsa.in","r",stdin),freopen("rsa.out","w",stdout);
rep(kase,1,rd()) {
n=rd<ll>(),m=rd<ll>(),c=rd<ll>();
ll T=sqrt(n),p=-1,q;
for(ll i=T;i>=T-100;--i) if(n%i==0) {
p=i,q=n/i;
break;
}
if(p==-1) {
// 2*sqrt(n) <= p+q <= sqrt(4*n+lambda * lambda)
//ll R=ceil(sqrt((long double)4*n+9e10)+1);
ll L=ceil(sqrt(n+0.5));
if((L-1)*(L-1)>=n) L--;

for(ll x=2*L;;++x) {
ull y=x*x-4*(ull)n;
// x=p+q
// t=q-p
ull t=sqrt(y);
if((t+1)*(t+1)==y) t++;
if((t-1)*(t-1)==y) t--;
if(t*t==y) {
p=(x-t)/2;
q=(x+t)/2;
break;
}
}
}

ll t=Inv(c,(p-1)*(q-1));
// t= 1/c (mod phi(n))

ll k1=p,b1=qpow(m%p,t,p);
ll k2=q,b2=qpow(m%q,t,q);
// k1 x + b1 = k2 y + b2
// k1 x = b2-b1 (mod k2)
// x= (b2-b1)/k1;
// x' = (b2-b1)/k1 (mod k2) * k1 + b1
ll inv=qpow(k1,k2-2,k2);
b1=(b2-b1)%k2*inv%k2 * k1+b1;
k1*=k2;
b1=(b1%k1+k1)%k1;
printf("%lld\n",b1);
}
}




「余姚中学 2019 联测 Day 4」随机除法

好题,难就难在转移的高位前缀和

首先是一个浅显的 \(\text{dp}\) 状态,令 \(n=\Pi prime_i^{c_i}\)

则状态只跟 \(\{c_i\}\) 有关,这是一个可重集合,强制定义 \(c_i\ge c_{i-1}\) 最小表示出所有不同状态

搜索一下 \(\text{dp}\) 状态,发现只有 \(170000\) 左右的状态数

直接枚举因数转移复杂度显然是升天的,直接枚举子集状态转移复杂度也很高,并且不好确定系数

所以用一个高位前缀和处理优化枚举因数的转移

再说一遍,是高位前缀和枚举因数的转移,不同的因数可能对应同一个状态

常见的高位前缀和是 \(dp_{i,S}=dp_{i-1,S}+dp_{i,S-\{S_i\}}\)

转移具有单调性,对于状态排序之后,定义辅助数组 \(dp_{i,j}\) 表示对于 \(i\) 这个状态

它的子集(注意这个子集是未排序的)中和它不同的最低位置 \(\ge j\) 的总和

计算高位前缀和时,每次转移只会对于一个位置改变

枚举状态时,取得位置是 \(j\) ,调用时需要排序

而排完序之后 \(j\) 可能会后移,所以需要定义成 \(\ge j\) 的,否则会算多

比如转移 \((1,1,1)\leftarrow (0,1,1),(1,0,1),(1,1,0)\)

如果定义成 \(\le j\) 的状态,三种状态转移之后都变成 \((1,1,0)\)

原先在这个状态里的三个位置编号是 \((0,1,2)\)

如果都去 \((1,1,0)\) 这个状态里转移过来,原先 \((0,1,2)\) 对应的下标位置改变,变成

\((1,2,0)\)

\((0,2,1)\)

\((0,1,2)\)

我们访问的时候访问的应该是子状态中不同位置 \(\le j\) 的总和

而从这个下标改变的状态里转移过来时,原先 \(>j\) 的下标被移移动进 \(\le j\) 的范围

再转移就错了

所以正确定义状态之后就可以高位前缀和了

存储和访问状态可以用 \(\text{Hash,Trie,int128}\) 三种方法存储, \(\text{int128}\) 真香啊

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
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef __int128 Node;
#define rep(i,a,b) for(int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(int i=a,i##end=b;i>=i##end;--i)
#define Mod1(x) ((x>=P)&&(x-=P))

const int N=180000,P=1e9+7;

ll qpow(ll x,ll k=P-2) {
ll res=1;
for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P;
return res;
}
Node Max,st[N];
int m,a[100],cnt,dp[N][20],F[N],pri[]={0,2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79};
char str[30];

int Find(Node s){ return lower_bound(st+1,st+cnt+1,s)-st; }
void Load(Node s){
rep(i,1,20) a[i]=0;
for(int i=1;i<=20 && s>1; ++i) while(s%pri[i]==0) s=s/pri[i],a[i]++;
}
void dfs(int p,int lst,Node s) {
st[++cnt]=s;
rep(i,1,lst) {
if((s*=pri[p])>Max) return;
dfs(p+1,i,s);
}
}


int main(){
freopen("div.in","r",stdin),freopen("div.out","w",stdout);
Max=1e12,Max=Max*Max;
dfs(1,100,1);
fprintf(stderr,"States Number = %d\n",cnt);
sort(st+1,st+cnt+1);
rep(i,2,cnt) {
Node now; Load(now=st[i]);
int c=1;
rep(i,1,18) c=1ll*c*(a[i]+1)%P;
drep(j,18,0) {
if(!a[j]) dp[i][j]=dp[i][j+1];
else {
int k=j;
while(k>1 && a[k-1]==a[k]) --k;
int p=Find(now/pri[j]);
drep(d,j,k) {
dp[i][d]=dp[i][d+1]+dp[p][d];
Mod1(dp[i][d]);
}
j=k;
}
}
F[i]=(dp[i][1]+c)*qpow(c-1)%P;
rep(j,1,18) dp[i][j]+=F[i],Mod1(dp[i][j]);
}
while(~scanf("%s%d",str,&m)) {
Node n=0;
for(int i=0;str[i];++i) n=n*10+(str[i]-'0');
rep(i,1,m) {
int x;scanf("%d",&x); a[i]=0;
while(n%x==0) n=n/x,a[i]++;
}
sort(a+1,a+m+1,greater <int> ());
n=1;
rep(i,1,m) rep(j,1,a[i]) n=n*pri[i];
printf("%d\n",F[Find(n)]);
}
}

「清华集训 2017」小 Y 和二叉树

原题数据好像没有卡这个情况

1
2
3
4
5
6
5
1 2
2 1 3
3 2 4 5
1 3
1 3

输出是

1
1 2 3 4 5

首先考虑一个 \(O(n^2)\) 的暴力:

枚举一个点为根,向下展开树,此时只需要决策左儿子和右儿子的顺序

当两个子树都存在时,由于两个子树包含的元素不同,所以可以直接把 两个子树序列首较小 (显然不会出现相同的情况) 的一个放在前面即可

实际上我们可以发现,这样得到的序列第一个元素必然是 编号最小的不同时包含左右儿子 的结点

不妨称固定根之后,这样的结点为叶子

\[ \ \]

显然的性质:任何一个度数 \(\leq 2\) 的点可以作为答案序列的第一个点

设原树上最小的 \(\leq 2\) 的点为 \(root\) ,接下来对于 \(root\) 的不同情况讨论,要在强制 \(root\) 为序列首的情况下,求得最小的序列

不妨先预处理出结点 \(u\) 子树里最小的叶子 \(mi_u\)

1.没有相邻结点,结束

2.有两个相邻结点,此时要使自己为序列首,必然有一个结点是自己的父亲,有一个结点是自己的右儿子

右儿子会被先遍历到,所以可以直接考虑比较两个相邻结点 作为 右儿子时谁的序列首 较小

即比较两个子树中最小的叶子即可

3.只有一个相邻结点,设其为 \(v\)

此时要使得自己为序列第一个,只有两种可能,此时同样可以考虑比较序列首元素

3-1.让相邻结点作为自己的父亲,此时下一个元素一定是 \(v\)

3-2.让相邻结点作为自己的右儿子,此时下一个元素一定是 \(mi_v\)

如果 \(v\ne mi_v\) ,显然好决策

\(v=mi_v\) 时,必然满足 \(v\) 是一个叶子,此时如果将 \(v\) 放在父亲上, \(v\) 的另一个相邻结点(如果存在)

可以放在 \(v\) 的父亲或者是 \(v\) 的右子树,如果放在 \(v\) 的右子树,那么这种情况与 \(v\) 被放在 \(u\) 的子树等价

也就是说,把 \(v\) 放在父亲可以决策出的序列情况,包含了把 \(v\) 放在右儿子的情况

所以这种情况也应当把 \(v\) 放在父亲上

实现上,不妨用两个 \(dfs\) 处理最后的决策,一个强制当前结点 \(u\) 为序列首,一个求出 \(u\) 子树的最优方案

在代码里就是 \(Solve,dfs\_get\)

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
#include<bits/stdc++.h>
using namespace std;

#define rep(i,a,b) for(int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(int i=a,i##end=b;i>=i##end;--i)
template <class T> inline void cmin(T &a,T b){ ((a>b)&&(a=b)); }

char IO;
int rd(){
int s=0,f=0;
while(!isdigit(IO=getchar())) if(IO=='-') f=1;
do s=(s<<1)+(s<<3)+(IO^'0');
while(isdigit(IO=getchar()));
return f?-s:s;
}

const int N=1e6+10;

int n;
int c[N],s[N][3];
int mi[N];

int rt=1e9;
void dfs(int u,int f) {
mi[u]=1e9;
int cnt=0;
rep(i,0,c[u]-1) if(s[u][i]!=f) {
int v=s[u][i]; cnt++;
dfs(v,u),cmin(mi[u],mi[v]);
}
if(cnt<=1) cmin(mi[u],u);
}

int vis[N];
void dfs_get(int u) {
vis[u]=1;
int a=-1,b=-1;
rep(i,0,c[u]-1) if(!vis[s[u][i]]) {
int v=s[u][i];
if(~a) b=v;
else a=v;
}
if(a==-1) printf("%d ",u);
else if(b==-1) {
if(mi[a]<u) dfs_get(a),printf("%d ",u);
else printf("%d ",u),dfs_get(a);
} else {
if(mi[a]>mi[b]) swap(a,b);
dfs_get(a),printf("%d ",u),dfs_get(b);
}
}

void Solve(int u) {
int cnt=0;
rep(i,0,c[u]-1) if(!vis[s[u][i]]) cnt++;
vis[u]=1,printf("%d ",u);
if(cnt==1) {
rep(i,0,c[u]-1) if(!vis[s[u][i]]) {
int v=s[u][i];
if(v>mi[v]) dfs_get(s[u][i]);
else Solve(v);
}
} else if(cnt==2) {
int a=-1,b=-1;
rep(i,0,c[u]-1) if(!vis[s[u][i]]) {
int v=s[u][i];
if(~a) b=v;
else a=v;
}
if(mi[a]>mi[b]) swap(a,b);
dfs_get(a),Solve(b);
}
}

int main(){
//freopen("binary.in","r",stdin),freopen("binary.out","w",stdout);
n=rd();
if(n==1) return puts("1"),0;
rep(i,1,n) {
c[i]=rd();
if(c[i]<=2) cmin(rt,i);
rep(j,0,c[i]-1) s[i][j]=rd();
}
dfs(rt,0);
Solve(rt);
}

「清华集训 2017」小 Y 和恐怖的奴隶主

是不是这题太水了都没人写啊

本题官方题解提供的做法实际上复杂度非常高

Part1

很显然本题的 \(\text{dp}\) 是存储每种血量的随从数量

设状态数量的上限是 \(S\)

\(m=3,k=8\) 时,这样的状态一共有 \(S=165+1\)

如果直接 \(dp\) ,每次转移是 \(O(1)\) 的,可以做到 \(O(n\cdot S)\) ,显然无法处理 \(n\) 较大的情况

用矩阵优化 \(\text{dp}\) 转移,朴素的实现可以做到 \(O(T\cdot \log n\cdot S^3)\)

如果预处理出转移矩阵的幂次,每次查询时只有列向量与方阵的乘法,所以复杂度是 \(O(\log n\cdot S^3+T\log n\cdot S^2)\)

实际极限的复杂度预估在 \(60\cdot 166^3+500\cdot 60\cdot 166^2\approx 11\cdot 10^8\)

官方题解提供的做法就是在在这个算法上进行常数优化,这wtm

\[ \ \]

Part2

对于已知 \(m,k\) 来说,设每个 \(n\) 构成的答案数列为 \(a_n\)

预先处理一部分的答案 \(a_1\cdots a_n\) ,使用 \(\text{Berlekamp-Massey }\) 算法求出序列的最短线性递推

发现总是在 \(O(S)\) 的长度以后,递推序列不再改变,即总能得到一个长度为 \(O(S)\) 的全局线性递推

即总可以在 \(O(S^2)\) 的时间内求出答案序列 \(a_n\)线性递推式

那么对于求得的线性递推式,问题转化为对于每个查询的 \(n\) ,求常系数线性递推数列的第 \(n\) 项答案

特征多项式的做法,可以做到单组查询 \(O(S\log n\log S)\) 的之间求出

那么总复杂度就是 \(O(S^2+T\cdot S \log S\log n)\)

理论上来说,复杂度上限应该只有 \(166^2+500\cdot 166\cdot 10\cdot 60\approx 0.5\cdot 10^8\)

理论上来说,这个复杂度无论是不是渐进意义下都比矩阵快

但是实际实践中,由于多项式运算的大常数,不优秀的实现下甚至可能超时

考虑到本题多查询的性质,我改变了倍增的基数 \(D\) ,并且预处理出倍增用到的 \(x^i\mod \lambda\)

预处理部分的复杂度是 \(O(\log_D^n \cdot (D-1)\cdot S\log S)\)

查询部分的复杂度是 \(O(\log _D^n S\log S)\)

由于我的多项式模板不够成熟

在LOJ上,经过这样的魔改,已经能跑得比矩阵块

但是在UOJ上,我同样的代码竟然慢了四倍,而这份代码在UOJ上的运行时间还略有提高

表示很是绝望。。

事实上,这种做法在 \(m,k\) 较大时同样可行,在 \(S\) 较大, \(T\) 较小的情况下,实际运行总能有更好的表现

「雅礼集训 2018 Day4」Magic(分治NTT)

题目的条件简直无法计算恰好为 \(k\) 的方案数,所以考虑计算 \(\ge k\) 的方案数

所以可以强制有 \(k\) 个相邻位置相同,但是不确定相同的是那些颜色

对每个颜色 \(a_i\) 考虑,设把 \(a_i\) 这个颜色分成了 \(b_i\) 个联通块(即强制了 \(a_i-b_i\) 个相邻位置相同)

那方案数就是 \(C(a_i-1,b_i-1)\) (是一个简单的插板问题)

得到每种颜色的联通块个数 \(b_i\) ,那么这些联通块之间排列的方案数就是 \(\frac{(\sum b_i)!}{b_i!}\)

容易得到 \(a_i\rightarrow b_i\) 的方案数,直接合并 \(b_i\) 的方案数,是一个背包问题,所以考虑用分治 \(NTT\) 快速合并

那么得到了 \(\sum b_i=n-k\) 的所有方案,复杂度 \(O(n\log n\log m)\)

最后的容斥比较明显的是一个二项式反演的形式,可以 \(O(n)\) 计算

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
#include<bits/stdc++.h>
using namespace std;

typedef long long ll;
#define rep(i,a,b) for(int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(int i=a,i##end=b;i>=i##end;--i)
#define pb push_back
#define Mod1(x) ((x>=P)&&(x-=P))
#define Mod2(x) ((x<0)&&(x+=P))

char IO;
int rd(){
int s=0;
while(!isdigit(IO=getchar()));
do s=(s<<1)+(s<<3)+(IO^'0');
while(isdigit(IO=getchar()));
return s;
}

const int N=1<<17,P=998244353;

int n,m,k;
int Inv[N+1],Fac[N+1],FInv[N+1];
ll C(int n,int m){ return n<0||m<0||n<m? 0 : 1ll*Fac[n]*FInv[m]%P*FInv[n-m]%P; }
ll qpow(ll x,ll k=P-2){
ll res=1;
for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P;
return res;
}
typedef vector <int> Poly;
int w[N|10],rev[N];
void Init(){
Inv[0]=Inv[1]=Fac[0]=Fac[1]=FInv[0]=FInv[1]=1;
rep(i,2,N){
Fac[i]=1ll*Fac[i-1]*i%P;
Inv[i]=1ll*(P-P/i)*Inv[P%i]%P;
FInv[i]=1ll*FInv[i-1]*Inv[i]%P;
}
w[N>>1]=1;
ll t=qpow(3,(P-1)/N);
rep(i,(N>>1)+1,N-1) w[i]=w[i-1]*t%P;
drep(i,(N>>1)-1,1) w[i]=w[i<<1];
}

int Init(int n){
int R=1,cc=-1;
while(R<n) R<<=1,cc++;
rep(i,1,R-1) rev[i]=(rev[i>>1]>>1)|((i&1)<<cc);
return R;
}
void NTT(int n,Poly &a,int f){
if((int)a.size()<n) a.resize(n);
rep(i,1,n-1) if(i<rev[i]) swap(a[i],a[rev[i]]);
for(int i=1;i<n;i<<=1){
int *e=w+i;
for(int l=0;l<n;l+=i*2){
for(int j=l;j<l+i;++j) {
int t=1ll*a[j+i]*e[j-l]%P;
a[j+i]=a[j]-t; Mod2(a[j+i]);
a[j]+=t; Mod1(a[j]);
}
}
}
if(f==-1){
reverse(a.begin()+1,a.end());
rep(i,0,n-1) a[i]=1ll*a[i]*Inv[n]%P;
}
}

Poly operator * (Poly a,Poly b){
int n=a.size()+b.size()-1,R=Init(n);
NTT(R,a,1),NTT(R,b,1);
rep(i,0,R-1) a[i]=1ll*a[i]*b[i]%P;
NTT(R,a,-1),a.resize(n);
return a;
}

Poly Solve(int l,int r){
if(l==r){
int x=rd();
Poly F(x+1);
rep(y,1,x) F[y]=C(x-1,y-1)*FInv[y]%P;
return F;
}
int mid=(l+r)>>1;
return Solve(l,mid)*Solve(mid+1,r);
}

int main(){
freopen("magic.in","r",stdin),freopen("magic.out","w",stdout);
Init(),n=rd(),m=rd(),k=rd();
Poly dp=Solve(1,n);
rep(i,n,m) dp[i]=1ll*dp[i]*Fac[i]%P;
int i=m-k;
rep(j,n,i-1) dp[i]=(dp[i]+(((i-j)&1)?-1:1)*dp[j]*C(m-j,m-i)%P+P)%P;
ll ans=(dp[i]%P+P)%P;
printf("%lld\n",ans);
}

「雅礼集训 2018 Day8」B

Solution1

设到达一个点的时间为 \(T_u\) ,从这个点出去的时间为 \(T_u'\)

那么显然满足 \(T_u\leq T_u'\leq T_u+t_u\) ,答案就是 \(\sum (t_u-(T'_u-T_u))\cdot c_u\)

对于一条边满足 \(T_v\ge T'_u\) ,二分答案之后,容易发现这是一个线性规划问题

可以暴力单纯形解决掉(当然是水的,但是好像还挺快。。)

Loj Submission

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
#include<bits/stdc++.h>
using namespace std;
typedef double db;
typedef pair <int,int> Pii;
#define reg register
#define mp make_pair
#define pb push_back
#define Mod1(x) ((x>=P)&&(x-=P))
#define Mod2(x) ((x<0)&&(x+=P))
#define rep(i,a,b) for(int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(int i=a,i##end=b;i>=i##end;--i)
template <class T> inline void cmin(T &a,T b){ ((a>b)&&(a=b)); }
template <class T> inline void cmax(T &a,T b){ ((a<b)&&(a=b)); }

char IO;
template <class T=int> T rd(){
T s=0; int f=0;
while(!isdigit(IO=getchar())) f|=IO=='-';
do s=(s<<1)+(s<<3)+(IO^'0');
while(isdigit(IO=getchar()));
return f?-s:s;
}

const int N=630,P=1e9+7;
const db eps=1e-9;

int n,m,W,t[N],c[N];
int U[N],V[N];
db A[N][N];


db Simplex(int n){
srand(time(NULL));
auto Pivot=[&](int x,int y) {
static int P[N],C;
rep(i,0,n) A[y][i]=-A[x][i]/A[x][y];
A[y][x]=1/A[x][y],A[y][y]=0;
rep(i,0,n) A[x][i]=0;
C=0;
rep(i,0,n) if(abs(A[y][i])>eps) P[++C]=i;
rep(i,0,n) if(abs(A[i][y])>eps) {
db t=A[i][y]; A[i][y]=0;
rep(j,1,C) A[i][P[j]]+=t*A[y][P[j]];
}
};
vector <int> V;
rep(i,1,n) V.pb(i);
//random_shuffle(V.begin(),V.end());
while(1) {
int u=0,v=0;
for(int i:V) if(!u || A[i][0]<A[u][0]) u=i;
if(A[u][0]>-eps) break;
for(int i:V) if(A[u][i]>eps) v=i;
if(!v) return puts("Infeasible"),0;
Pivot(u,v);
}

while(1) {
int u=0,v=0;
for(int i:V) if(!v || A[0][i]>A[0][v]) v=i;
if(A[0][v]<eps) break;
for(int i:V) if(A[i][v]<-eps) if(!u || (A[i][0]/A[i][v] > A[u][0]/A[u][v])) u=i;
if(!u) return puts("Unbounded"),0;
Pivot(u,v);
}
return A[0][0];
}
int outd[N];
int Check(int lim) {
memset(A,0,sizeof A);
int cnt=n*2;
rep(i,1,n) A[0][i+n]=c[i],A[0][i]=-c[i],A[0][0]-=t[i]*c[i];
rep(i,1,n) {
A[++cnt][i]=-1; A[cnt][i+n]=1;
A[++cnt][i]=1; A[cnt][i+n]=-1; A[cnt][0]=t[i];
}
rep(i,1,m) A[++cnt][U[i]+n]=-1,A[cnt][V[i]]=1;
rep(i,1,n) if(!outd[i]) A[++cnt][0]=lim,A[cnt][i+n]=-1;
db res=-Simplex(cnt);
return res<=W+eps;
}

int main(){
freopen("soft.in","r",stdin),freopen("soft.out","w",stdout);
n=rd(),m=rd(),W=rd();
int l=0,r=0,res=-1;
rep(i,1,n) t[i]=rd(),r+=t[i];
rep(i,1,n) c[i]=rd();
rep(i,1,m) U[i]=rd(),V[i]=rd(),outd[U[i]]++;
while(l<=r) {
int mid=(l+r)>>1;
if(Check(mid)) r=mid-1,res=mid;
else l=mid+1;
}
printf("%d\n",res);
}

\[ \ \]

Solution2

二分答案 \(\text{lim}\) ,问题转化为求最小花费

设每个点减少了 \(x_i\)

考虑限制有两种,一种是路径长度的限制,一种是每个点大小的限制

\(\text{minimize:} \sum x_i\cdot c_i\)

\(\displaystyle \forall p\in paths , \sum x_{p_i}\ge \sum t_{p_i}-\text{lim}\)

\(-x_i\ge -t_i\)

对偶一下,设对于路径 \(p\)\(\sum x_{p_i}\) 的对偶变量为 \(y_p\)\(-x_i\) 的对偶变量为 \(z_i\)

\(\text{maximize}:\sum y_p\cdot (\sum t_{p_i}-\text{lim})-z_i\cdot t_i\)

\(\displaystyle \forall i\in[1,n], \sum_{p\in paths,i\in p} y_p-z_i\leq c\)

考虑对偶变量 \(y_p\)\(z_i\) 有什么意义

此时,选择一条路径 \(y_p\) ,会使得 关于路径上的点的限制+1 , 使得答案增加 \(\sum t_{p_i}-\text{lim}\)

\(z_i\) 是关于每个单点的变量,可以用 \(t_i\) 代价使得每个 \(i\) 的限制-1

那么可以考虑转化为一个路径覆盖问题,选择一条路径覆盖路径上的点,且得到 \(\sum t_{p_i}-\text{lim}\) 的价值

限制式子转化为:每个点被覆盖次数大于 \(c\) 时,再选择就要付出 \(t_i\) 的代价令 \(z_i\) 加一

带权的路径覆盖容易转化为费用流模型,可以把每个点拆成入点出点,每个点被覆盖前 \(c_i\) 次,价值为 \(t_i\) ,之后就为0

因此每个点的入点向出点连 \((c_i,t_i),(\infty,0)\) 两条边即可,路径的 \(\text{-lim}\) 可以在源点前加入

求一次最大费用可行流,最终得到的答案是原问题的最小代价

是我EK写得太丑的说: Loj Submission

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
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
#define rep(i,a,b) for(int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(int i=a,i##end=b;i>=i##end;--i)
template <class T> inline void cmin(T &a,T b){ ((a>b)&&(a=b)); }
template <class T> inline void cmax(T &a,T b){ ((a<b)&&(a=b)); }

char IO;
template <class T=int> T rd(){
T s=0; int f=0;
while(!isdigit(IO=getchar())) f|=IO=='-';
do s=(s<<1)+(s<<3)+(IO^'0');
while(isdigit(IO=getchar()));
return f?-s:s;
}

enum{N=10101,INF=1<<30};

int n,m,W;
int c[N],t[N],X[N],Y[N];
struct Edge{
int to,nxt,w,c;
} e[N];
int head[N],ecnt=1;
int V,S,T;
void AddEdge(int u,int v,int w,int c){ e[++ecnt]=(Edge){v,head[u],w,c},head[u]=ecnt; }
void Link(int u,int v,int w,int c){ AddEdge(u,v,w,c),AddEdge(v,u,0,-c); }

int pre[N],dis[N],inq[N];
int SPFA(int lim){
rep(i,1,V) dis[i]=-INF;
static queue <int> que;
dis[S]=-lim,que.push(S);
while(!que.empty()) {
int u=que.front(); que.pop(),inq[u]=0;
for(int i=head[u];i;i=e[i].nxt) {
int v=e[i].to;
if(!e[i].w || dis[v]>=dis[u]+e[i].c) continue;
dis[v]=dis[u]+e[i].c,pre[v]=i;
if(!inq[v]) que.push(v),inq[v]=1;
}
}
return dis[T]>0;
}

int Check(int lim){
S=++V,T=++V;
rep(i,1,n) {
Link(S,++V,INF,0);
Link(V,V+1,c[i],t[i]);
Link(V,V+1,INF,0);
Link(++V,T,INF,0);
}
rep(i,1,m) Link(X[i]*2+2,Y[i]*2+1,INF,0);
int ans=0;
while(SPFA(lim)){
int w=INF;
for(int i=T;i!=S;i=e[pre[i]^1].to) cmin(w,e[pre[i]].w);
for(int i=T;i!=S;i=e[pre[i]^1].to) e[pre[i]].w-=w,e[pre[i]^1].w+=w;
ans+=dis[T]*w;
}
rep(i,1,V) head[i]=0;
ecnt=1,V=0;
return ans<=W;
}

int main(){
freopen("soft.in","r",stdin),freopen("soft.out","w",stdout);
n=rd(),m=rd(),W=rd();
int l=0,r=0,res=-1;
rep(i,1,n) r+=t[i]=rd();
rep(i,1,n) c[i]=rd();
rep(i,1,m) X[i]=rd(),Y[i]=rd();
for(int mid;l<=r;) Check(mid=(l+r)>>1)?r=mid-1,res=mid:l=mid+1;
printf("%d\n",res);
}

「517模拟赛1」ABC难题

Tips: 模数是 \(10^8+7\)

不妨令序列总长为 \(n\) ,包含 \(m\) 种不同的数字

发现存在 \(ABC\) 的子序列只需要满足最左边的 \(A\) 和最右边的 \(C\) 之间有一个 \(B\)

由于可能出现多个 \(ABC\) ,考虑计算不存在 \(ABC\) 的情况

那么可以枚举最左边的 \(A\) 和最右边的 \(C\) 分别为 \(X,Y\) ,那么需要不存在 \(ABC\) ,并且枚举的 \(X,Y\) 合法 的充要条件是

\([1,X-1]\) 中不存在 \(A\)\([X+1,Y-1]\) 中不存在 \(B\)\([Y+1,n]\) 中不存在 \(C\)

这三个限制,每一个限制会让一种颜色能够选择的字母数量-1

那么枚举 \(X\) ,扫描每一个 \(Y\) (注意不一定满足 \(X< Y\) ),同时记录每个数字是否出现在 \([1,X-1],[X+1,Y-1],[Y+1,n]\)

在每次扫描时改变限制,同步统计出受到 \(0,1,2,3\) 个限制的数字的个数,然后答案就是 \((3-i)^k\) 之积

然而这样还不够,因为可能根本不存在 \(A,C\) ,不存在 \(A\) 或者 \(C\) 的答案为 \(2^m\) ,同时不存在 \(A,C\) 的答案为 \(1^m\) ,容斥一下即可

0%