題目傳送門

matthew99神犇的題解講得非常清楚明白，跪爛Orzzzzzzzzzzzzz

總結一下，本題有很多重要的突破口

1.Lucas定理

看到n,m特別大但模數特別小時，容易想到$lucas$定理

$C_{n}^{m}=C_{n/p}^{m/p}\cdot C_{n\;mod\;p}^{m\;mod\;p}\;(mod\;p)$

但普通的$lucas$顯然不適用於多次計算，我們可以把$lucas$定理展開

我們把$n$和$m$都看成兩個$p$進制數$a$和$b$

$C_{n}^{m}=\pi C_{a_{i}}^{b_{i}}\;(mod\;p)$

這個展開顯然成立

2.數位$DP$

想到了上一條進制轉化，很容易就聯想到數位$DP$

定義$f[i][j]$表示枚舉到第$i$位，余數是$j$的方案數

轉移十分經典，設現在要填上的數是$x$，$0$表示沒達到上界，$1$達到上界，設$z=C_{x}^{b_{i+1}}$

$f0[i][j\cdot z\;mod\;p]+=f0[i][j]$

$x<a_{i+1}$時，$f0[i][j\cdot z\;mod\;p]+=f1[i][j]$

$x=a_{i+1}$時，$f1[i][j\cdot z\;mod\;p]+=f1[i][j]$

總復雜度$O(p^2log_{p}n)$

3.原根優化與多項式乘法

上面的$p^{2}$轉移咋這麽無腦呢，是不是有啥優化啊？是的

由於$p$是一個質數，它一定存在一個原根$g$

這就要涉及到離散對數了。其實這是一個映射

$g^{ind(x)}\equiv x\;(mod\;p)$

$ind(x)\;(ind(x)\in[0,p-1))$與$x\;(x\in [1,p))$是一組一一映射

那麽$j\cdot z\;mod\;p$

$=g^{ind(j)+ind(z)\;(mod\;p-1)}\;mod\;p$

我們利用映射關系，把底數化成指數

這樣轉移變成了卷積的形式，用多項式乘法每次$O(plogp)$計算

計算出結果後，再逆映射回來得到實際的答案

而離散對數的映射並不能處理余數等於$0$的情況，我們每次$O(p)$單獨討論即可

總時間$O(plogplog_{p}n)$

  1 #include <bits/stdc++.h>
  2 #define N1 (1<<16)+10
  3 #define dd double
  4 #define ld long double
  5 #define ll long long
  6 #define uint unsigned int 
  7 #define i128 __int128
  8 using namespace std;
  9 
 10 const int inf=0x3f3f3f3f;
 11 i128 gi128()
 12 {
 13     i128 ret=0;int fh=1;char c=getchar();
 14     while(c<‘0‘||c>‘9‘){if(c==‘-‘)fh=-1;c=getchar();}
 15     while(c>=‘0‘&&c<=‘9‘){ret=ret*10+c-‘0‘;c=getchar();}
 16     return ret*fh;
 17 }
 18 ll qpow(ll x,ll y,const int &p)
 19 {
 20     ll ans=1;
 21     for(;y;x=x*x%p,y>>=1) if(y&1) ans=ans*x%p;
 22     return ans;
 23 }
 24 const ll p=998244353;
 25 namespace NTT{
 26 ll a[N1],b[N1],c[N1],Wn[N1],_Wn[N1];
 27 int r[N1];
 28 ll qpow(ll x,ll y)
 29 {
 30     ll ans=1;
 31     for(;y;x=x*x%p,y>>=1) if(y&1) ans=ans*x%p;
 32     return ans;
 33 }
 34 void NTT(ll *s,int len,int type)
 35 {
 36     int i,j,k; ll wn,w,t;
 37     for(i=0;i<len;i++) if(i<r[i]) swap(s[i],s[r[i]]);
 38     for(k=2;k<=len;k<<=1)
 39     {
 40         wn=(type>0)?Wn[k]:_Wn[k];
 41         for(i=0;i<len;i+=k)
 42         {
 43             for(j=0,w=1;j<(k>>1);j++,w=w*wn%p)
 44             {
 45                 t=w*s[i+j+(k>>1)]%p;
 46                 s[i+j+(k>>1)]=(s[i+j]+p-t)%p;
 47                 s[i+j]=(s[i+j]+t)%p;
 48             }
 49         }
 50     }
 51 }
 52 void Pre(int len,int L)
 53 {
 54     int i;
 55     for(i=0;i<len;i++) r[i]=(r[i>>1]>>1)|((i&1)<<(L-1));
 56     for(i=1;i<=len;i<<=1) Wn[i]=qpow(3,(p-1)/i), _Wn[i]=qpow(Wn[i],p-2);
 57 }
 58 void Main(int len)
 59 {
 60     int i,invl=qpow(len,p-2);
 61     NTT(a,len,1); NTT(b,len,1);
 62     for(i=0;i<len;i++) c[i]=a[i]*b[i]%p;
 63     NTT(c,len,-1);
 64     for(i=0;i<len;i++) c[i]=c[i]*invl%p;
 65     memset(a,0,sizeof(a)); memset(b,0,sizeof(b)); 
 66 }
 67 };
 68 
 69 int T,m,G;
 70 i128 n,l,r;
 71 int a[N1],b[N1],tmp[N1];
 72 ll f0[2][N1],f1[2][N1],mul[N1],_mul[N1],ans[N1];
 73 
 74 inline ll C(int N,int M)
 75 { 
 76     if(M>N) return 0;
 77     return mul[N]*_mul[M]%m*_mul[N-M]%m;
 78 }
 79 
 80 int pr[N1],use[N1],ind[N1],_ind[N1],son[N1];
 81 void get_ind()
 82 {
 83     int i,j,ns=0,flag,x,cnt=0,sz=0;
 84     for(i=2;i<=m-1;i++)
 85     {
 86         if(!use[i]) pr[++cnt]=i;
 87         for(j=1;j<=cnt&&i*pr[j]<=m-1;j++)
 88         {
 89             use[i*pr[j]]=1;
 90             if(i%pr[j]==0) break;
 91         }
 92     }
 93     for(j=1;j<=cnt;j++) if((m-1)%pr[j]==0) son[++sz]=(m-1)/pr[j];
 94     for(i=2;i<=m-2;i++) 
 95     {
 96         flag=1;
 97         for(j=1;j<=sz;j++) if(qpow(i,son[j],m)==1){ flag=0; break;}
 98         if(flag){ G=i; break; }
 99     }
100     ind[1]=0; _ind[0]=1; ind[0]=m-1;
101     for(i=1,x=1;i<=m-2;i++) x=x*G%m, ind[x]=i, _ind[i]=x; 
102 }
103 //using NTT::a; using NTT::b; using NTT::c; 
104 void solve(i128 w,int type)
105 {
106     int i,j,k,nw=0,nn=0,len,L,now=0,pst=1;
107     if(w<n){ ans[0]=((w+1)*type%p+ans[0]+p)%p; return; }
108     while(w>0) nw++,tmp[nw]=w%m,w/=m;
109     for(i=1;i<=nw;i++) a[nw-i+1]=tmp[i];
110     i128 N=n;
111     while(N>0) nn++,tmp[nn]=N%m,N/=m;
112     for(i=1;i<=nn;i++) b[nw-i+1]=tmp[i];
113     
114     for(len=1,L=0;len<m+m-1;len<<=1,L++);
115     NTT::Pre(len,L);
116     
117     memset(f0,0,sizeof(f0)); memset(f1,0,sizeof(f1));
118     for(j=0;j<a[1];j++) f0[1][C(j,b[1])]++;
119     f1[1][C(a[1],b[1])]++;
120     for(i=1;i<nw;i++)
121     {
122         memset(f0[now],0,sizeof(f0[now])); memset(f1[now],0,sizeof(f1[now])); 
123         
124         for(j=1;j<m;j++) NTT::a[ind[j]]=f0[pst][j];
125         for(j=0;j<m;j++) NTT::b[ind[C(j,b[i+1])]]++;
126         for(j=1;j<m;j++) (f0[now][0]+=f0[pst][j]*NTT::b[m-1])%=p;
127         (f0[now][0]+=f0[pst][0]*m)%=p; NTT::b[m-1]=0;
128         NTT::Main(len);
129         for(j=0;j<len;j++) (f0[now][_ind[j%(m-1)]]+=NTT::c[j])%=p;
130         
131         for(j=1;j<m;j++) NTT::a[ind[j]]=f1[pst][j];
132         for(j=0;j<a[i+1];j++) NTT::b[ind[C(j,b[i+1])]]++;
133         for(j=1;j<m;j++) (f0[now][0]+=f1[pst][j]*NTT::b[m-1])%=p;
134         (f0[now][0]+=f1[pst][0]*a[i+1])%=p; NTT::b[m-1]=0;
135         NTT::Main(len);
136         for(j=0;j<len;j++) (f0[now][_ind[j%(m-1)]]+=NTT::c[j])%=p;
137         
138         for(k=0;k<m;k++)
139         (f1[now][k*C(a[i+1],b[i+1])%m]+=f1[pst][k])%=p;
140         
141         swap(now,pst);
142     }
143     for(i=0;i<m;i++) (ans[i]+=(f0[pst][i]+f1[pst][i])*type%p+p)%=p;
144 }
145 
146 int main()
147 {
148     int i,j,x,cnt=0;
149     scanf("%d",&m); n=gi128(); l=gi128(); r=gi128(); l--;
150     mul[0]=mul[1]=_mul[0]=_mul[1]=1;
151     for(i=2;i<m;i++) mul[i]=mul[i-1]*i%m, _mul[i]=qpow(mul[i],m-2,m);
152     get_ind();
153     solve(r,1); 
154     solve(l,-1);
155     for(i=0;i<m;i++) printf("%lld\n",ans[i]);
156     return 0;
157 }  

UOJ #86 mx的組合數 (數位DP+NTT+原根優化)