[BJ United Round #3] 押韵
先%%%%%%%%%%%%%%%%% EI
下文默认模数为P
简要题意:求:用k种元素,每种元素使用d的倍数次,排成一个长度为nd的序列 的方案数
这个题目的设定就让人想到两个离不开的元素 : (模数暗示了?)
指数型生成函数 + 单位根反演
显然可以得到每一种元素的指数型生成函数为
EGF(Element)=F(x)=∑d|ixii!
带入单位根反演 [d|n]=∑d−10ωindd
即F(x)=1d⋅d−1∑i=0eωidx
而总的生成函数就是 G(x)=Fk(x)
即G(x)=1dk⋅(d−1∑i=0eωidx)k
其中的和式幂次展开会得到一个kd项的多项式,我们要求[xn]G(x),就需要展开得到每一项的幂系数
所以显然我们需要先合并同类项一下。。。
而幂系数是一个单位根之和的形式,这就需要我们寻找单位根之间的关系
这里得到一个思路:用d次单位根中的φ(d)个作为基底,以简单的 有理数/整系数 表示出所有的ωid
对于d=4的情况比较简单,φ(d)=2,可以得到四个单位根分别为1,ω,−1,−ω
可以枚举得到的和为x+yω,然后求系数
优先考虑组合意义,可以发现就是在平面上每次可以走四个方向,k步之后最终到达(x,y)的方案数
两个维度分立的情况,还需要枚举每个维度走了几步,所以用一种巧妙的转化两个维度联系在一起
将平面旋转π8,并且扩大√2倍,得到新的坐标为(x−y,x+y),新的行走方向是(+1,+1),(−1,−1),(−1,+1),(+1,−1)
这样以来,每次每个维度都有行走,可以确定每个维度+1和−1的次数,直接组合数排列即可得到答案
对于d=6,甚至是更一般的情况的情况
只在代数层面来看单位根似乎十分抽象,不如从复平面单位根上找一找灵感
下面是d=6的情形

φ(6)=2,假设以−−→OA,−−→OB作为基底,可以直观地得到基底表达
|
−−→OA=1 |
−−→OB=ω |
−−→OA=ω06 |
1 |
0 |
−−→OB=ω16 |
0 |
1 |
−−→OC=ω26 |
-1 |
1 |
−−→OD=ω36 |
-1 |
0 |
−−→OE=ω46 |
0 |
-1 |
−−→OF=ω56 |
1 |
-1 |
由此我们得到了一个φ(d)维数的表达方法
把每一维看做不同元,也就是说,得到了一个φ(d)维,O(1)次的多项式,需要我们求高维多项式幂次
令N=kφ(d)
直接压位暴力多项式复杂度为O(NlogN−Nlog2N),而且面临着模数难以处理,常数大的问题
所以EI又用出了一个巧妙的暴力方法解决这个问题,以d=6为例,先做一下处理,得到要求的多项式
似乎每次k次幂总是求导+递推?
f(x,y)=x2y+xy2+y2+y+x+x2,g(x,y)=fk(x,y)
g(x,y)对于x求偏导,得到g′(x,y)=kfk−1(x,y)f′(x,y)
即g′(x,y)f(x,y)=kg(x,y)f′(x,y)
f′(x,y)=2xy+2x+y2+1
然后我们要解这个方程,考虑乘积为[xnym]一项两边的系数
左边=[xn−2ym−1]+[xn−1ym−2]+[xnym−2]+[xnym−1]+[xn−1ym]+[xn−2ym]
换成g(x,y)的系数应该是
左边=(n−1)[xn−1ym−1]+n[xnym−2]+(n+1)[xn+1ym−2]+(n+1)[xn+1ym−1]+n[xnym]+(n−1)[xn−1ym]
右边=2k[xn−1ym−1]+2k[xn−1ym]+k[xnym−2]+k[xnym]
其中[xn+1ym−1]只出现了一次,按照先n递增再m递增的顺序进行递推,即
[xnym]=2k[xn−2ym]+2k[xn−2ym+1]+k[xn−1ym−1]+k[xn−1ym+1]n
−(n−2)[xn−2ym]+(n−1)[xn−1ym−1]+n[xnym−1]+(n−1)[xn−1ym+1]+(n−2)[xn−2ym+1]n
边界条件是 [xi]=[yi](i≥k)=(ki−k) (由系数x,x2或y,y2得到)
由此带入递推即可
综上,得到的每项的系数的复杂度为O(d⋅kφ(d)) ,其中d为递推每项需要的时间
由系数得到答案仍然需要一次快速幂,因此依然带一个logP
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
| #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 pb push_back #define mp make_pair #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=min(a,b); } template <class T> inline void cmax(T &a,T b){ a=max(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=2e3+10,P=1049874433,G=7;
int n,k,d; ll qpow(ll x,ll k=P-2) { k%=P-1; ll res=1; for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P; return res; } int w[100],C[N][N],I[N*2];
int main() { n=rd(),k=rd(),d=rd(); w[0]=1,w[1]=qpow(G,(P-1)/d); rep(i,2,90) w[i]=1ll*w[i-1]*w[1]%P; rep(i,0,k) rep(j,C[i][0]=1,i) C[i][j]=C[i-1][j]+C[i-1][j-1],Mod1(C[i][j]); I[0]=I[1]=1; rep(i,2,k*2) I[i]=1ll*(P-P/i)*I[P%i]%P; if(d==1){ int ans=qpow(k,1ll*d*n); printf("%d\n",ans); } else if(d==2) { int ans=0; rep(i,0,k) ans=(ans+qpow((1ll*w[0]*i+1ll*w[1]*(k-i))%P,1ll*d*n)*C[k][i])%P; ans=ans*qpow(qpow(d,k))%P; printf("%d\n",ans); } else if(d==3) { int ans=0; rep(i,0,k) rep(j,0,k-i) ans=(ans+qpow((1ll*w[0]*i+1ll*w[1]*j+1ll*w[2]*(k-i-j))%P,1ll*d*n)*C[i+j][i]%P*C[k][i+j])%P; ans=ans*qpow(qpow(d,k))%P; printf("%d\n",ans); } else if(d==4) { int ans=0; rep(i,-k,k) rep(j,-k,k) if(abs(i)+abs(j)<=k && (k-i-j)%2==0) { ll x=qpow((1ll*w[0]*i+1ll*w[1]*j)%P,1ll*d*n); ll y=1ll*C[k][(abs(i-j)+k)/2]*C[k][(k+abs(i+j))/2]%P; ans=(ans+x*y)%P; } ans=(ans+P)*qpow(qpow(d,k))%P; printf("%d\n",ans); } else { static int F[N*2][N*2]; int ans=0; rep(i,0,k*2) rep(j,max(k-i,0),min(2*k,3*k-i)) { if(i==0) F[i][j]=C[k][j-k]; else if(j==0) F[i][j]=C[k][i-k]; else { int s=(2ll*k*(i>1?F[i-2][j]+F[i-2][j+1]:0)+1ll*k*(F[i-1][j-1]+F[i-1][j+1]))%P; int t=((i>1?1ll*(i-2)*(F[i-2][j]+F[i-2][j+1]):0)+ 1ll*(i-1)*F[i-1][j-1]+1ll*i*F[i][j-1]+1ll*(i-1)*F[i-1][j+1])%P; F[i][j]=1ll*(s-t+P)*I[i]%P; } ans=(ans+qpow((1ll*w[0]*(i-k)+1ll*w[1]*(j-k))%P,1ll*d*n)*F[i][j])%P; } ans=(ans+P)*qpow(qpow(d,k))%P; printf("%d\n",ans); } }
|