Lucas&exLucas
lucas
long long f[10010],g[10010];
long long ksm(long long x,long long y,long long mo)
{
long long cnt=1;
while(y)
{
if (y&1) cnt=(cnt*x)%mo;
x=(x*x)%mo;
y>>=1;
}
return cnt;
}
long long NY(long long a,long long p) {return ksm(a,p-2,p);}
void getfg(long long n)
{
f[0]=g[0]=1;
for (long long i=1;i<=n;i++) (f[i]=i*f[i-1])%=p,g[i]=NY(f[i],p);
}
long long getC(long long x,long long y) {if (x<y) return 0;return (f[x]*g[x-y]%p*g[y])%p;}
long long lucas(long long x,long long y)
{
if (y == 0) return 1;
return (getC(x%p,y%p)*lucas(x/p,y/p))%p;
}
exlucas
msjing甚至贺了三份
第一份
点击查看代码
#include<bits/stdc++.h>
using namespace std;
constexpr int maxn=50010;
long long read() {long long x=0,f=1;char ch=getchar();while (ch<'0' || ch>'9'){if (ch=='-') f=-1;ch=getchar();}while (ch>='0' && ch<='9'){x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}return x*f;}
long long n,m,p;
long long exgcd(long long a,long long b,long long &x,long long &y)
{
if (!b) {x=1,y=0;return a;}
long long res=exgcd(b,a%b,x,y);
long long t=x;
x=y,y=t-a/b*y;
return res;
}
long long inv(long long a,long long p)
{
long long x,y;
exgcd(a,p,x,y);
return (x%p+p)%p;
}
long long ksm(long long x,long long y,long long p)
{
long long cnt=1;
while (y)
{
if (y&1) cnt=(cnt*x)%p;
x=(x*x)%p;
y>>=1;
}
return cnt;
}
long long F(long long n,long long p,long long kp)
{
if (!n) return 1;
long long ro=1,re=1;
for (long long i=1;i<=kp;i++) if (i%p) (ro*=i)%=kp;
ro=ksm(ro,n/kp,kp);
// printf("%lld\n",ro);
for (long long i=kp*(n/kp);i<=n;i++) if (i%p) (re*=(i%kp))%=kp;
return F(n/p,p,kp)*ro%kp*re%kp;
}
long long G(long long n,long long p)
{
if (n<p) return 0;
return G(n/p,p)+n/p;
}
long long C(long long n,long long m,long long p,long long kp)
{
long long fz=F(n,p,kp),fm1=inv(F(m,p,kp),kp),fm2=inv(F(n-m,p,kp),kp);
// printf("%lld,%lld,%lld\n",fz,fm1,fm2);
long long mi=ksm(p,G(n,p)-G(m,p)-G(n-m,p),kp);
return fz*fm1%kp*fm2%kp*mi%kp;
}
long long a[1010],b[1010];
long long exlucas(long long n,long long m,long long p)
{
long long ljc=p,tot=0;
for (long long tmp=2;tmp*tmp<=p;tmp++)
{
if (!(ljc%tmp))
{
long long kp=1;
while (!(ljc%tmp)) kp*=tmp,ljc/=tmp;
// printf("%lld,%lld\n",kp,ljc);
a[++tot]=kp,b[tot]=C(n,m,tmp,kp);
// printf("%lld\n",b[tot]);
}
}
if (ljc!=1) a[++tot]=ljc,b[tot]=C(n,m,ljc,ljc);
// for (int i=1;i<=tot;i++) printf("%lld,%lld\n",a[i],b[i]);
long long ans=0;
for (long long i=1;i<=tot;i++)
{
long long M=p/a[i],T=inv(M,a[i]);
(ans+=b[i]*M%p*T%p)%=p;
}
return ans;
}
int main()
{
n=read(),m=read(),p=read();
printf("%lld\n",exlucas(n,m,p));
return 0;
}
点击查看代码
#include<bits/stdc++.h>
using namespace std;
constexpr int maxn=50010;
long long read() {long long x=0,f=1;char ch=getchar();while (ch<'0' || ch>'9'){if (ch=='-') f=-1;ch=getchar();}while (ch>='0' && ch<='9'){x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}return x*f;}
long long n,m,p;
long long exgcd(long long a,long long b,long long &x,long long &y)
{
if (!b) {x=1,y=0;return a;}
long long res=exgcd(b,a%b,x,y);
long long t=x;
x=y,y=t-a/b*y;
return res;
}
long long inv(long long a,long long p)
{
long long x,y;
exgcd(a,p,x,y);
return (x%p+p)%p;
}
long long ksm(long long x,long long y,long long p)
{
long long cnt=1;
while (y)
{
if (y&1) cnt=(cnt*x)%p;
x=(x*x)%p;
y>>=1;
}
return cnt;
}
long long fac(long long n,long long pi,long long pk)
{
if (!n) return 1;
long long res=1;
for (long long i=2;i<=pk;i++) if (i%pi) (res*=i)%=pk;
res=ksm(res,n/pk,pk);
for (long long i=2;i<=n%pk;i++) if (i%pi) (res*=i)%=pk;
return res*fac(n/pi,pi,pk)%pk;
}
long long crt(long long b,long long mo) {return b*inv(p/mo,mo)%p*(p/mo)%p;}
long long C(long long n,long long m,long long pi,long long pk)
{
long long up=fac(n,pi,pk),d1=fac(m,pi,pk),d2=fac(n-m,pi,pk);
long long k=0;
for (long long i=n;i;i/=pi) k+=i/pi;
for (long long i=m;i;i/=pi) k-=i/pi;
for (long long i=n-m;i;i/=pi) k-=i/pi;
return up*inv(d1,pk)%pk*inv(d2,pk)%pk*ksm(pi,k,pk)%pk;
}
long long exlucas(long long n,long long m)
{
long long res=0,tmp=p,pk;
for (int i=2;i*i<=p;i++)
{
if (tmp%i == 0)
{
pk=1;
while (tmp%i == 0) pk*=i,tmp/=i;
(res+=crt(C(n,m,i,pk),pk))%=p;
}
}
if (tmp>1) (res+=crt(C(n,m,tmp,tmp),tmp))%=p;
return res;
}
int main()
{
n=read(),m=read(),p=read();
printf("%lld\n",exlucas(n,m));
return 0;
}
点击查看代码
#include<bits/stdc++.h>
using namespace std;
constexpr int maxn=10010;
long long read() {long long x=0,f=1;char ch=getchar();while (ch<'0' || ch>'9'){if (ch=='-') f=-1;ch=getchar();}while (ch>='0' && ch<='9'){x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}return x*f;}
long long n,m,P,tot;
struct _ {long long pos,val;}c[maxn];
void exgcd(long long a,long long b,long long &x,long long &y)
{
if (!b) {x=1,y=0;return;}
exgcd(b,a%b,x,y);
long long t=x;
x=y,y=t-a/b*y;
}
long long inv(long long a,long long p)
{
long long x,y;
exgcd(a,p,x,y);
return (x+p)%p;
}
long long ksm(long long x,long long y,long long p)
{
long long cnt=1;
while (y)
{
if (y&1) cnt=(cnt*x)%p;
x=(x*x)%p;
y>>=1;
}
return cnt;
}
long long rem[maxn];/*remainder*/
void bk(long long x)
{
for (long long i=2;i*i<=x;i++)
{
if (x%i == 0)
{
c[++tot]={i,1};
while (x%i == 0) c[tot].val*=i,x/=i;
}
}
if (x>1) c[++tot]={x,x};
}
long long fac(long long n,long long p,long long k)
{
if (!n) return 1;
long long res=1;
for (long long i=2;i<k;i++) if (i%p) (res*=i)%=k;
res=ksm(res,n/k,k);
for (long long i=2;i<=n%k;i++) if (i%p) (res*=i)%=k;
return res*fac(n/p,p,k)%k;
}
long long C(long long n,long long m,long long p,long long k)
{
if (n<m) return 0;
long long aa=fac(n,p,k),bb=fac(m,p,k),cc=fac(n-m,p,k);
long long cnt=0;
for (long long i=p;i<=n;i*=p) cnt+=n/i;
for (long long i=p;i<=m;i*=p) cnt-=m/i;
for (long long i=p;i<=n-m;i*=p) cnt-=(n-m)/i;
return aa*inv(bb,k)%k*inv(cc,k)%k*ksm(p,cnt,k)%k;
}
long long exlucas()
{
bk(P);
long long ans=0;
for (long long i=1;i<=tot;i++) rem[i]=C(n,m,c[i].pos,c[i].val);
for (long long i=1;i<=tot;i++)
{
long long mi=P/c[i].val;
ans+=rem[i]%P*(mi*inv(mi,c[i].val)%P)%P;
}
(ans+=P)%=P;
return ans;
}
int main()
{
n=read(),m=read(),P=read();
printf("%lld\n",exlucas());
return 0;
}
浙公网安备 33010602011771号