组合数学的实质还是DP,但是从通式角度处理的话有利于FFT等的实现。
首先推荐$Candy?$的球划分问题集合:
http://www.cnblogs.com/candy99/p/6400735.html
以下部分节选自
http://blog.csdn.net/sr_19930829/article/details/40888349
第一类Stirling数
定理:第一类Stirling数$s(p,k)$计数的是把p个对象排成k个非空循环排列的方法数。
证明:把上述定理叙述中的循环排列叫做圆圈。递推公式为:
$s(p,p)=1 (p>=0)$ 有p个人和p个圆圈,每个圆圈就只有一个人
$s(p,0)=0 (p>=1)$ 如果至少有1个人,那么任何的安排都至少包含一个圆圈$$s(p,k)=(p-1)*s(p-1,k)+s(p-1,k-1)$$ 设人被标上$1,2,\ldots,p$。将这p个人排成k个圆圈有两种情况。第一种排法是在一个圆圈里只有标号为$p$的人自己,排法有$s(p-1,k-1)$个。第二种排法中,$p$至少和另一个人在一个圆圈里。这些排法可以通过把$1,2,\ldots,p-1$排成k个圆圈再把p放在$1,2,\ldots,p-1$任何一人的左边得到,因此第二种类型的排法共有$(p-1)*s(p-1,k)$种排法。
在证明中我们所做的就是把$\{1,2,\ldots,p\}$划分到k个非空且不可区分的盒子,然后将每个盒子中的元素排成一个循环排列。
组合数通式:$$s(k,i)=\frac{C_n^k\ k!}{\sum_{i=0}^{k}(-1)^{i+k}n^i}$$
第一类斯特林数快速求法:
首先我们有:$n^{\underline{k}}=\sum\limits_{i=0}^{k}(-1)^{k-i}s(k,i)\cdot n^i$,$n^{\overline{k}}=\sum\limits_{i=0}^{k}s(k,i)\cdot n^i$
也就是$s(n,k)=[x^k]\prod\limits_{i=0}^{n-1}(x+i)$
于是一种显然的做法就是分治NTT,复杂度$O(n\log^2 n)$。
1 #include<cstdio> 2 #include<algorithm> 3 #define rep(i,l,r) for (int i=(l); i<=(r); i++) 4 using namespace std; 5 6 const int N=300010,mod=998244353; 7 int n,A,B,top,stk[52],rev[N],fac[N],inv[N],S[N],tmp[52][N]; 8 9 int ksm(int a,int b){ 10 int res=1; 11 for (; b; a=1ll*a*a%mod,b>>=1) 12 if (b & 1) res=1ll*res*a%mod; 13 return res; 14 } 15 16 void NTT(int a[],int n,int f){ 17 for (int i=0; i<n; i++) if (i<rev[i]) swap(a[i],a[rev[i]]); 18 for (int i=1; i<n; i<<=1){ 19 int wn=ksm(3,f ? (mod-1)/(i<<1) : (mod-1)-(mod-1)/(i<<1)); 20 for (int p=i<<1,j=0; j<n; j+=p){ 21 int w=1; 22 for (int k=0; k<i; k++,w=1ll*w*wn%mod){ 23 int x=a[j+k],y=1ll*w*a[i+j+k]%mod; 24 a[j+k]=(x+y)%mod; a[i+j+k]=(x-y+mod)%mod; 25 } 26 } 27 } 28 if (f) return; 29 int inv=ksm(n,mod-2); 30 for (int i=0; i<n; i++) a[i]=1ll*a[i]*inv%mod; 31 } 32 33 int solve(int l,int r,int a[]){ 34 if (l==r){ a[0]=l; a[1]=1; return 1; } 35 int ls=stk[top--],rs=stk[top--],mid=(l+r)>>1; 36 int l1=solve(l,mid,tmp[ls]),l2=solve(mid+1,r,tmp[rs]),n=1,L=0; 37 for (; n<=l1+l2; n<<=1) L++; 38 for (int i=0; i<n; i++) rev[i]=(rev[i>>1]>>1)|((i&1)<<(L-1)); 39 NTT(tmp[ls],n,1); NTT(tmp[rs],n,1); 40 for (int i=0; i<n; i++) a[i]=1ll*tmp[ls][i]*tmp[rs][i]%mod; 41 NTT(a,n,0); stk[++top]=ls; stk[++top]=rs; 42 for (int i=0; i<n; i++) tmp[ls][i]=tmp[rs][i]=0; 43 return l1+l2; 44 } 45 46 int main(){ 47 freopen("960g.in","r",stdin); 48 freopen("960g.out","w",stdout); 49 scanf("%d%d%d",&n,&A,&B); 50 if (n==1){ if (A==1 && B==1) puts("1"); else puts("0"); return 0; } 51 rep(i,1,50) stk[i]=i; top=50; 52 solve(0,n-2,S); 53 fac[0]=1; rep(i,1,A+B) fac[i]=1ll*fac[i-1]*i%mod; 54 inv[A+B]=ksm(fac[A+B],mod-2); 55 for (int i=A+B-1; ~i; i--) inv[i]=1ll*inv[i+1]*(i+1)%mod; 56 printf("%lld\n",1ll*S[A+B-2]*fac[A+B-2]%mod*inv[A-1]%mod*inv[B-1]%mod); 57 return 0; 58 }