P4720 【模板】扩展卢卡斯定理/exLucas

题面传送门

感觉除了复杂度和卢卡斯定理没有任何关系。

首先我们考虑先把\(p\)唯一分解成\(\prod\limits_{i} p_i^{c_i}\)，然后对于每个\(p_i^{c_i}\)求出\(C_{n}^{m}\bmod p_i^{c_i}\)的值，然后CRT合并即可。

\(C_n^m\bmod p_{i}^{c_i}=\frac{n!}{m!(n-m)!}\bmod p_i^{c_i}\)。

一个基本的思路是，我们将每个阶乘拆成\(\frac{n!}{p^x}\)的形式，然后就可以exgcd求逆元了。

设\(f(n)=\frac{n!}{p^x}\)，则显然的，\(n!\)中至少有\(\frac{n}{p}\)个数有一个\(p\)，则可以写成\(f(n)=(\frac{n}{p})!(\prod\limits_{i=1,i\bmod p\ne1}^{p_i^{c_i}}{i})^{\frac{n}{p_i^{c_i}}}\prod \limits_{i=\lfloor \frac{n}{p_i^{c_i}}\rfloor}^{n}{i}\)

但是这样显然也是不对的，因为\((\frac{n}{p})!\)还有\(p\)的因数，因此写成递归形式，将\((\frac{n}{p})!\)换成\(f(\frac{n}{p})\)即可。

时间复杂度大概\(O(p\log p)\)。

code:

#include<bits/stdc++.h>
#define I inline
#define max(a,b) ((a)>(b)?(a):(b))
#define min(a,b) ((a)<(b)?(a):(b))
#define abs(x) ((x)>0?(x):-(x))
#define ll long long
#define db double
#define lb long db
#define N (100000+5)
#define M (6000000+5)
#define K (1500+5)
#define mod 998244353
#define Mod (mod-1)
#define eps (1e-5)
#define ull unsigned ll
#define it iterator
#define Gc() getchar() 
#define Me(x,y) memset(x,y,sizeof(x))
#define Mc(x,y) memcpy(x,y,sizeof(x))
#define d(x,y) ((k+1)*(x)+(y))
#define R(n) (1ll*rand()*rand()%(n)+1)
#define Pc(x) putchar(x)
#define LB lower_bound
#define UB upper_bound
#define PB push_back
using namespace std;
int p,k,H,B[N],C[N];ll n,m,ToT;
I ll mpow(ll x,ll y,int p){ll Ans=1;while(y) y&1&&(Ans=Ans*x%p),y>>=1,x=x*x%p;return Ans;}
I void EG(ll x,int p,ll &a,ll &b){if(!p){a=1;b=0;return;}EG(p,x%p,a,b);swap(a,b);b-=x/p*a;}
I ll GI(ll x,int p){ll y,z;EG(x,p,y,z);return (y%p+p)%p;}I ll GG(ll x,int p){if(!x) return 0;return x/p+GG(x/p,p);}
I ll GF(ll n,int p,int pk){if(!n) return 1;ll ToT=1;for(int i=0;i<=pk;i++) if(i%p) ToT=ToT*i%pk;ToT=mpow(ToT,n/pk,pk);for(ll i=n/pk*pk;i<=n;i++) if(i%p) ToT=ToT*(i%pk)%pk;return ToT*GF(n/p,p,pk)%pk;}
I ll calc(ll n,ll m,int p,int pk){return GF(n,p,pk)*GI(GF(m,p,pk),pk)%pk*GI(GF(n-m,p,pk),pk)%pk*mpow(p,GG(n,p)-GG(m,p)-GG(n-m,p),pk)%pk;}
int main(){
	freopen("1.in","r",stdin);
	int i,j;scanf("%lld%lld%d",&n,&m,&p);k=p;for(i=2;i*i<=k;i++) {if(k%i==0){B[++H]=i;C[H]=1;while(k%i==0) C[H]*=i,k/=i;}}k^1&&(B[++H]=k,C[H]=k);
	for(i=1;i<=H;i++) ToT+=p/C[i]*GI(p/C[i],C[i])%p*calc(n,m,B[i],C[i])%p;printf("%lld\n",ToT%p);
}