首先这个题面就很胃疼
而且感觉讲的不是很清楚
实际上是要求满足如下条件的树的个数：
对于每个非叶节点，有
c
+
2
c+2
c+2个儿子，其中有
c
∈
A
c\in A
c∈A个叶子节点和2个非叶节点
且点的编号满足父亲小于儿子
按照套路设 f [ i ] f[i] f[i]为 i i i个点的答案可以列出 d p dp dp式
f
[
i
]
=
1
2
∑
j
∑
k
[
i
−
j
−
k
−
1
∈
A
]
(
i
−
1
j
)
(
i
−
j
−
1
k
)
f
j
f
k
f[i]=\frac 1 2\sum_{j}\sum_{k}[i-j-k-1\in A]{i-1\choose j}{i-j-1\choose k}f_jf_k
f[i]=21∑j∑k[i−j−k−1∈A](ji−1)(ki−j−1)fjfk
中间
1
2
\frac 1 2
21是由于两个非叶儿子会把每种情况算2次
边界
f
[
1
]
=
1
f[1]=1
f[1]=1
设 A ( x ) = ∑ i ∈ A x i i ! , f ( x ) = ∑ i f [ i ] x i i ! A(x)=\sum_{i\in A}\frac{x^i}{i!},f(x)=\sum_{i}f[i]\frac{x^i}{i!} A(x)=∑i∈Ai!xi,f(x)=∑if[i]i!xi
那么有 ∑ i f [ i + 1 ] x i i + 1 = f ′ = 1 2 f 2 A + 1 \sum_{i}f[i+1]\frac{x^i}{i+1}=f'=\frac 1 2f^2A+1 ∑if[i+1]i+1xi=f′=21f2A+1
这个可以分治
N
T
T
NTT
NTT做到
O
(
n
l
o
g
2
n
)
O(nlog^2n)
O(nlog2n)而且跑的比
n
l
o
g
n
nlogn
nlogn快
不过注意如果
l
≠
1
l\not=1
l=1的话实际上情况只统计到了一次
而前面有一个
1
2
\frac 1 2
21,所以算的时候需要乘2
#include<bits/stdc++.h>
using namespace std;
#define cs const
#define pb push_back
#define re register
#define ll long long
#define pii pair<int,int>
#define fi first
#define bg begin
#define se second
#define poly vector<int>
cs int RLEN=1<<20|1;
inline char gc(){
static char ibuf[RLEN],*ib,*ob;
(ib==ob)&&(ob=(ib=ibuf)+fread(ibuf,1,RLEN,stdin));
return (ib==ob)?EOF:*ib++;
}
#define gc getchar
inline int read(){
char ch=gc();
int res=0;bool f=1;
while(!isdigit(ch))f^=ch=='-',ch=gc();
while(isdigit(ch))res=(res+(res<<2)<<1)+(ch^48),ch=gc();
return f?res:-res;
}
template<class tp>inline void chemx(tp &a,tp b){a<b?a=b:0;}
template<class tp>inline void chemn(tp &a,tp b){a>b?a=b:0;}
cs int mod=998244353,G=3;
inline int add(int a,int b){return (a+=b)>=mod?a-mod:a;}
inline void Add(int &a,int b){(a+=b)>=mod?(a-=mod):0;}
inline int dec(int a,int b){a-=b;return a+(a>>31&mod);}
inline void Dec(int &a,int b){a-=b,a+=a>>31&mod;}
inline int mul(int a,int b){static ll r;r=1ll*a*b;return r>=mod?r%mod:r;}
inline void Mul(int &a,int b){static ll r;r=1ll*a*b;a=r>=mod?r%mod:r;}
inline int ksm(int a,int b,int res=1){for(;b;b>>=1,Mul(a,a))(b&1)&&(Mul(res,a),1);return res;}
inline int Inv(int x){return ksm(x,mod-2);}
cs int C=19;
int *w[C+1],rev[(1<<C)|5];
inline void init_w(){
for(int i=1;i<=C;i++)w[i]=new int[1<<(i-1)];
int wn=ksm(G,(mod-1)/(1<<C));
w[C][0]=1;
for(int i=1;i<(1<<(C-1));i++)w[C][i]=mul(w[C][i-1],wn);
for(int i=C-1;i;i--)
for(int j=0;j<(1<<(i-1));j++)w[i][j]=w[i+1][j<<1];
}
inline void init_rev(int lim){
for(int i=0;i<lim;i++)rev[i]=(rev[i>>1]>>1)|((i&1)*(lim>>1));
}
inline void ntt(poly &f,int lim,int kd){
for(int i=0;i<lim;i++)if(i>rev[i])swap(f[i],f[rev[i]]);
for(int mid=1,l=1,a0,a1;mid<lim;mid<<=1,l++)
for(int i=0;i<lim;i+=mid<<1)
for(int j=0;j<mid;j++)
a0=f[i+j],a1=mul(f[i+j+mid],w[l][j]),f[i+j]=add(a0,a1),f[i+j+mid]=dec(a0,a1);
if(kd==-1){
reverse(f.bg()+1,f.bg()+lim);
for(int i=0,iv=Inv(lim);i<lim;i++)Mul(f[i],iv);
}
}
inline poly operator *(poly a,poly b){
int deg=a.size()+b.size()-1,lim=1;
if(deg<=32){
poly c(deg,0);
for(int i=0;i<a.size();i++)
for(int j=0;j<b.size();j++)
Add(c[i+j],mul(a[i],b[j]));
return c;
}
while(lim<deg)lim<<=1;
init_rev(lim);
a.resize(lim),b.resize(lim);
ntt(a,lim,1),ntt(b,lim,1);
for(int i=0;i<lim;i++)Mul(a[i],b[i]);
ntt(a,lim,-1),a.resize(deg);
return a;
}
cs int N=200005;
int inv[N<<1],fac[N<<1],ifac[N<<1];
inline void init_inv(){
cs int len=(N-5)<<1;
fac[0]=ifac[0]=inv[0]=inv[1]=1;
for(int i=2;i<=len;i++)inv[i]=mul(mod-mod/i,inv[mod%i]);
for(int i=1;i<=len;i++)fac[i]=mul(fac[i-1],i),ifac[i]=mul(ifac[i-1],inv[i]);
}
int a[N],g[N],f[N],n;
char s[N];
void cdq(int l,int r){
if(l==r)return;
int mid=(l+r)>>1,l1=r-l+1;
cdq(l,mid);
poly h,p;
for(int i=l;i<=mid;i++)h.pb(f[i]);
for(int i=1;i<=l1;i++)p.pb(f[i]);
h=h*p;int coef=1+(l>1);
for(int i=mid+1;i<=r;i++)Add(g[i],mul(coef,h[i-l-1]));
h.clear(),p.clear();
for(int i=l;i<=mid;i++)h.pb(g[i]);
for(int i=1;i<=l1;i++)p.pb(a[i]);
h=h*p;
for(int i=mid+1;i<=r;i++)Add(f[i],mul(inv[i*2],h[i-l-1]));
cdq(mid+1,r);
}
int main(){
#ifdef Stargazer
freopen("lx.cpp","r",stdin);
#endif
init_w(),init_inv();
n=read();
scanf("%s",s+1);
for(int i=1;i<=n;i++)a[i]=mul((s[i]=='1'),ifac[i-1]);
f[1]=1,cdq(1,n);
for(int i=1;i<=n;i++)cout<<mul(f[i],fac[i])<<'\n';
}
接下来是
O
(
n
l
o
g
n
)
O(nlogn)
O(nlogn)的做法
而且是解常微分方程的一个通用套路
先设
h
(
x
)
=
1
2
A
x
2
+
1
h(x)=\frac 1 2Ax^2+1
h(x)=21Ax2+1
那么实际要解的就是
f
′
(
x
)
=
h
(
f
(
x
)
)
f'(x)=h(f(x))
f′(x)=h(f(x))
按照牛顿迭代的通用方法
假设已经求出
f
0
′
(
x
)
≡
h
(
f
0
(
x
)
)
%
x
⌈
x
2
⌉
f_0'(x)\equiv h(f_0(x))\%x^{\lceil \frac x 2\rceil}
f0′(x)≡h(f0(x))%x⌈2x⌉
求
f
′
(
x
)
≡
h
(
f
(
x
)
)
%
x
n
f'(x)\equiv h(f(x))\% x^n
f′(x)≡h(f(x))%xn
那么
f
′
(
x
)
=
h
(
f
0
)
+
h
′
(
f
0
)
(
f
−
f
0
)
f'(x)=h(f_0)+h'(f_0)(f-f_0)
f′(x)=h(f0)+h′(f0)(f−f0)
由于这个东西无法方便的牛顿迭代
考虑构造
v
(
x
)
=
e
x
p
(
−
∫
h
′
(
f
0
(
x
)
)
d
x
)
v(x)=exp(-\int h'(f_0(x))\mathrm{dx})
v(x)=exp(−∫h′(f0(x))dx)
那么有
v
′
=
−
h
′
(
f
0
)
v
v'=-h'(f_0)v
v′=−h′(f0)v
等式两边乘上
v
v
v
有
f
′
v
=
v
(
h
(
f
0
)
+
h
′
(
f
0
)
(
f
−
f
0
)
)
f'v=v(h(f_0)+h'(f_0)(f-f_0))
f′v=v(h(f0)+h′(f0)(f−f0))
f
′
v
−
f
h
′
(
f
0
)
v
=
v
(
h
(
f
0
)
−
h
′
(
f
0
)
f
0
)
f'v-fh'(f_0)v=v(h(f_0)-h'(f_0)f_0)
f′v−fh′(f0)v=v(h(f0)−h′(f0)f0)
f
′
v
+
f
v
′
=
v
(
h
(
f
0
)
−
h
′
(
f
0
)
f
0
)
f'v+fv'=v(h(f_0)-h'(f_0)f_0)
f′v+fv′=v(h(f0)−h′(f0)f0)
(
f
v
)
′
=
v
(
h
(
f
0
)
−
h
′
(
f
0
)
f
0
)
(fv)'=v(h(f_0)-h'(f_0)f_0)
(fv)′=v(h(f0)−h′(f0)f0)
f
=
1
v
∫
(
v
(
h
(
f
0
)
−
h
′
(
f
0
)
f
0
)
)
d
x
f=\frac{1}{v}\int(v(h(f_0)-h'(f_0)f_0))\mathrm{dx}
f=v1∫(v(h(f0)−h′(f0)f0))dx
因为这个和具体 h h h没有关系所以可以通用
对于原来的 h h h
那么就有 f = 1 v ∫ ( v ( 1 − 1 2 A f 0 2 ) ) d x f=\frac 1 v\int(v(1-\frac 1 2Af_0^2))\mathrm{dx} f=v1∫(v(1−21Af02))dx
复杂度 O ( n l o g n ) O(nlogn) O(nlogn)跑的远没有 n l o g 2 n nlog^2n nlog2n快
这东西知道一下就好了
除了理论复杂度以外其他都比不过分治
n
t
t
ntt
ntt
#include<bits/stdc++.h>
using namespace std;
#define cs const
#define pb push_back
#define re register
#define ll long long
#define pii pair<int,int>
#define fi first
#define bg begin
#define se second
#define poly vector<int>
cs int RLEN=1<<20|1;
inline char gc(){
static char ibuf[RLEN],*ib,*ob;
(ib==ob)&&(ob=(ib=ibuf)+fread(ibuf,1,RLEN,stdin));
return (ib==ob)?EOF:*ib++;
}
#define gc getchar
inline int read(){
char ch=gc();
int res=0;bool f=1;
while(!isdigit(ch))f^=ch=='-',ch=gc();
while(isdigit(ch))res=(res+(res<<2)<<1)+(ch^48),ch=gc();
return f?res:-res;
}
template<class tp>inline void chemx(tp &a,tp b){a<b?a=b:0;}
template<class tp>inline void chemn(tp &a,tp b){a>b?a=b:0;}
cs int mod=998244353,G=3;
inline int add(int a,int b){return (a+=b)>=mod?a-mod:a;}
inline void Add(int &a,int b){(a+=b)>=mod?(a-=mod):0;}
inline int dec(int a,int b){a-=b;return a+(a>>31&mod);}
inline void Dec(int &a,int b){a-=b,a+=a>>31&mod;}
inline int mul(int a,int b){static ll r;r=1ll*a*b;return r>=mod?r%mod:r;}
inline void Mul(int &a,int b){static ll r;r=1ll*a*b;a=r>=mod?r%mod:r;}
inline int ksm(int a,int b,int res=1){for(;b;b>>=1,Mul(a,a))(b&1)&&(Mul(res,a),1);return res;}
inline int Inv(int x){return ksm(x,mod-2);}
cs int C=20;
int *w[C+1],rev[(1<<C)|5],inv[(1<<C)+5];
inline void init_w(){
for(int i=1;i<=C;i++)w[i]=new int[1<<(i-1)];
int wn=ksm(G,(mod-1)/(1<<C));
inv[0]=inv[1]=w[C][0]=1;
for(int i=1;i<(1<<(C-1));i++)w[C][i]=mul(w[C][i-1],wn);
for(int i=C-1;i;i--)
for(int j=0;j<(1<<(i-1));j++)w[i][j]=w[i+1][j<<1];
for(int i=2,lim=(1<<(C-1));i<lim;i++)inv[i]=mul(mod-mod/i,inv[mod%i]);
}
inline void init_rev(int lim){
for(int i=0;i<lim;i++)rev[i]=(rev[i>>1]>>1)|((i&1)*(lim>>1));
}
inline void ntt(poly &f,int lim,int kd){
for(int i=0;i<lim;i++)if(i>rev[i])swap(f[i],f[rev[i]]);
for(int mid=1,l=1,a0,a1;mid<lim;mid<<=1,l++)
for(int i=0;i<lim;i+=mid<<1)
for(int j=0;j<mid;j++)
a0=f[i+j],a1=mul(f[i+j+mid],w[l][j]),f[i+j]=add(a0,a1),f[i+j+mid]=dec(a0,a1);
if(kd==-1){
reverse(f.bg()+1,f.bg()+lim);
for(int i=0,iv=Inv(lim);i<lim;i++)Mul(f[i],iv);
}
}
inline poly operator *(poly a,poly b){
int deg=a.size()+b.size()-1,lim=1;
if(deg<=32){
poly c(deg,0);
for(int i=0;i<a.size();i++)
for(int j=0;j<b.size();j++)
Add(c[i+j],mul(a[i],b[j]));
return c;
}
while(lim<deg)lim<<=1;
init_rev(lim);
a.resize(lim),b.resize(lim);
ntt(a,lim,1),ntt(b,lim,1);
for(int i=0;i<lim;i++)Mul(a[i],b[i]);
ntt(a,lim,-1),a.resize(deg);
return a;
}
inline poly Inv(poly a,int deg){
poly b(1,Inv(a[0])),c;
for(int lim=4;lim<(deg<<2);lim<<=1){
c.resize(lim>>1);
for(int i=0;i<(lim>>1);i++)c[i]=(i<a.size())?a[i]:0;
init_rev(lim);
b.resize(lim),c.resize(lim);
ntt(c,lim,1),ntt(b,lim,1);
for(int i=0;i<lim;i++)Mul(b[i],dec(2,mul(c[i],b[i])));
ntt(b,lim,-1),b.resize(lim>>1);
}b.resize(deg);return b;
}
inline poly deriv(poly a){
for(int i=0;i+1<a.size();i++)a[i]=mul(a[i+1],i+1);
a.pop_back();return a;
}
inline poly integ(poly a){
a.pb(0);
for(int i=(int)a.size()-1;~i;i--)a[i]=mul(a[i-1],inv[i]);
a[0]=0;return a;
}
inline poly Ln(poly a,int deg){
a=integ(Inv(a,deg)*deriv(a)),a.resize(deg);return a;
}
inline poly exp(poly a,int deg){
poly b(1,1),c;
for(int lim=2;lim<(deg<<1);lim<<=1){
c=Ln(b,lim);
for(int i=0;i<lim;i++)c[i]=dec((i<a.size())?a[i]:0,c[i]);
Add(c[0],1);
b=b*c,b.resize(lim);
}b.resize(deg);
return b;
}
cs int N=200010;
int n,a[N],fac[N],ifac[N];
inline poly calc(int lim){
if(lim==1)return poly(1,0);
poly f0=calc((lim+1)>>1);
poly p,v;
for(int i=0;i<lim;i++)p.pb(a[i]);
p=p*f0,v=integ(p),v.resize(lim);
for(int i=0;i<lim;i++)v[i]=dec(0,v[i]);
v=exp(v,lim);
poly h=p*f0;h.resize(lim);
for(int i=0;i<lim;i++)h[i]=dec(0,mul(inv[2],h[i]));
Add(h[0],1),h=integ(v*h);
h=h*Inv(v,lim),h.resize(lim);
return h;
}
inline void init_inv(){
cs int len=N-5;
fac[0]=ifac[0]=1;
for(int i=1;i<=len;i++)fac[i]=mul(fac[i-1],i);
ifac[len]=Inv(fac[len]);
for(int i=len-1;i;i--)ifac[i]=mul(ifac[i+1],i+1);
}
char s[N];
int main(){
#ifdef Stargazer
freopen("lx.cpp","r",stdin);
#endif
init_w(),init_inv();
n=read();
scanf("%s",s);
for(int i=0;i<n;i++)a[i]=mul(ifac[i],s[i]=='1');
poly ans=calc(n+1);
for(int i=1;i<=n;i++)cout<<mul(fac[i],ans[i])<<'\n';
return 0;
}
304
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
