?
?題目大意
定義一個$S-$四面體表示六條邊由$S$根不同的木棍組成,定義一種染色方法合法當且僅當至少有$S$根木棍被染色且與每個頂點相鄰的三根木棍中至多有一根被染色,求有$N$個$S=1,2...N$四面體,求至少染$K$個的方案數。
?
題解
單獨考慮$S=1$四面體,染它有$9$中方案,否則將與每個頂點相鄰的$12$條邊單獨拿出來計算。
考慮其余四面體被染色的方案數:他們中有$k$個頂點$(k=0,1,2,3,4)$相鄰的木棍三條中有一條被染色,方案是$C_4^k\times 3^k$,剩余的邊則有$\sum\limits_{i=S-k}^{6S-12}C_{6S-12}^{i}$,用乘法原理乘一下即可。
由于$k$很小,可以只求$\sum\limits_{i=S}^{6S-12}C_{6S-12}^{i}$,其余加上組合數即可。而這個式子是楊輝三角上的某一行的一個后綴,所以是可以遞推的,考慮上一行答案對下一行答案的貢獻,只需要乘二再加上上一行后綴左側的那個組合數即可。
我們現在已經解決了每一個$S$的方案數,接下來就是求染出至少$K$個的方案數,也只需要遞推,設$G_x$表示染$x$個的方案數$V$表示當染前四邊形的方案數,$G_x=G_x+V\times G_{x-1}$。這是$N$個簡單最高次項次數是一次的多項式的卷積,用分治$FFT$解決即可。
有兩個細節,由于模數只有$10^5+3$,所以算組合數要用$lucas$定理。并且,$FFT$中每一位可能達到$10^{15}$級別,精度可能會爆炸,所以要優化或者使用$long\space double$。
#include<algorithm>
#include<iostream>
#include<cstring>
#include<cstdio>
#include<cmath>
#define LL long long
#define LD long double
#define mod 100003
#define M 100020
using namespace std;
int read(){int nm=0,fh=1; int cw=getchar();for(;!isdigit(cw);cw=getchar()) if(cw=='-') fh=-fh;for(;isdigit(cw);cw=getchar()) nm=nm*10+(cw-'0');return nm*fh;
}
const int P[]={1,12,54,108,81};
const LD PI=3.141592653589793238462643383;
int mul(int x,int y){return (LL)x*(LL)y%mod;}
int add(int x,int y){return (x+y)>=mod?x+y-mod:x+y;}
int mus(int x,int y){return (x-y)<0?x-y+mod:x-y;}
void upd(int &x,int y){x=add(x,y);}
int qpow(int x,int sq){int res=1;for(;sq;sq>>=1,x=mul(x,x)) if(sq&1) res=mul(res,x);return res;
}
int n,m,fac[mod],ifac[mod],G[M<<2],lg[M<<2];
int C(int tot,int tk){if(tot<0||tk<0||tot<tk)return 0;return mul(fac[tot],mul(ifac[tot-tk],ifac[tk]));}
int lucas(int tot,int tk){if(!tk) return 1;return mul(C(tot%mod,tk%mod),lucas(tot/mod,tk/mod));}
int rev[M<<2],F[M<<2],S[M];
struct comp{LD r,d; comp(){r=d=0.0;}comp(LD _r,LD _d){r=_r,d=_d;}comp operator +(const comp&k)const{return comp(r+k.r,d+k.d);}comp operator -(const comp&k)const{return comp(r-k.r,d-k.d);}comp operator *(const comp&k)const{return comp(r*k.r-d*k.d,r*k.d+d*k.r);}
}A[M<<2],B[M<<2];
void FFT(comp *x,int len,LD kd){for(int i=1;i<len;i++) if(i<rev[i]) swap(x[i],x[rev[i]]);for(int tt=1;tt<len;tt<<=1){comp unit(cos(PI*kd/(tt*1.0)),sin(PI*kd/(tt*1.0)));for(int st=0;st<len;st+=(tt<<1)){comp now(1.0,0.0);for(int pos=st;pos<st+tt;pos++,now=now*unit){comp t1=x[pos],t2=x[pos+tt]*now;x[pos]=t1+t2,x[pos+tt]=t1-t2;}}}if(kd<0.0){for(int i=0;i<len;i++) x[i].r/=(len*1.0);}
}
void tms(int *x,int *x1,int *x2,int n1,int n2){if(n1+n2<=120){memset(S,0,sizeof(int)*(n1+n2+1));for(int i=0;i<=n1;i++){for(int j=0;j<=n2;j++) S[i+j]+=mul(x1[i],x2[j]);}for(int i=0;i<=n1+n2;i++) x[i]=S[i]%mod;}else{int len=1,nw=-1;for(;len<=n1+n2+1;len<<=1,nw++);for(int i=1;i<len;i++) rev[i]=(rev[i>>1]>>1)|((i&1)<<nw);for(int i=0;i<=n1;i++) A[i]=comp(x1[i]*1.0,0.0);for(int i=0;i<=n2;i++) B[i]=comp(x2[i]*1.0,0.0);for(int i=n1+1;i<=len;i++) A[i]=comp(0.0,0.0);for(int i=n2+1;i<=len;i++) B[i]=comp(0.0,0.0);FFT(A,len,1.0),FFT(B,len,1.0);for(int i=0;i<len;i++) A[i]=A[i]*B[i]; FFT(A,len,-1.0);for(int i=0;i<=n1+n2;i++) x[i]=llround(A[i].r)%mod;}
}
void solve(int *x,int L,int R){if(L==R){x[0]=1,x[1]=G[L];return;}int mid=((L+R)>>1),ls,rs; ls=mid-L+1,rs=R-mid;solve(x,L,mid),solve(x+ls+1,mid+1,R);tms(x,x,x+ls+1,ls,rs);
}
int main(){fac[0]=ifac[0]=1,G[1]=9,G[2]=243,G[3]=16224,G[4]=46489;for(int i=1;i<mod;i++) fac[i]=mul(fac[i-1],i),ifac[i]=qpow(fac[i],mod-2);for(int i=5,K,pre,last=3797,rem;i<M;i++){for(K=(i-2)*6,pre=K-6;pre<K;pre++) last=add(add(last,last),lucas(pre,i-2));last=mus(last,lucas(K,i-1)),rem=last;for(int k=0;k<=4;k++) upd(G[i],mul(P[k],rem)),upd(rem,lucas(K,i-k-1));}for(int T=read(),ans=0;T;--T,ans=0){n=read(),m=read(),solve(F,1,n);for(int i=n;i>=m;i--) upd(ans,F[i]);printf("%d\n",ans);}return 0;
}