2 条题解
-
0
by hansang:
/* by;hansang 求的是 C nm %P (P不一定是质数),把 C nm 设为 x 把 P 化为 p1^k1*p2^k2*... (pi是质数) 列出式子: x % (p1^k1) = a1 x % (p2^k2) = a2 ...... 用中国剩余定理合并可得答案 */ #include<bits/stdc++.h> using namespace std; typedef long long LL; const int N=1e6+10; LL a[N], c[N]; int cnt; //a 数组表示 x (C nm ) % (pi^ki) 的值,c 数组表示 pi^ki,cnt 表示 a c 数组长度 LL q_pow(LL a, LL b, LL P) //求 a ^b %P { LL ans=1%P; a%=P; //防止 P 为 1,a 先%P for(; b; b>>=1) //枚举 b 的二进制位 { if(b&1) ans=ans*a%P; //b 这个位为 1,ans 就乘这个位的 a a=a*a%P; // a=a^2 } return ans; } void exgcd(LL a, LL b, LL &d, LL &x, LL &y) //扩展欧几里得 { if(b==0) {d=a; x=1; y=0;} //找到最大公约数了,d=a else { exgcd(b, a%b, d, y, x); y-=(a/b)*x; /* 设新的 x=dx,y=dy 则有; ax +by = b*tx +(a%b)*ty 已知; a%b = a-[a/b]*b ([]为下取整) 则; ax +by = b*tx +(a%b)*ty = b*tx +(a-[a/b]*b)*ty = b*tx +a*ty -[a/b] *b*ty = a*ty +b*(tx-[a/b]*ty) 所以; ty=x,tx=(y+[a/b]*x) 倒推回来 y 就要减 (a/b)*x */ } } // exgcd 和 q_pow 都为很重要的板子,请熟背 LL inv(LL a, LL P) //逆元 (运用 exgcd) { LL A=a, B=P, K=1, d, x, y; // ax=1 (%P) = ax+P*?=1 = Ax+By=K exgcd(A, B, d, x, y); //算的是 Ax+By=d (d为 gcd(a, b)) x=x*(K/d); LL dx=abs(B/d); x=(x%dx+dx)%dx; //把 x变成最小正整数解 // x*(K/d)变成 Ax+By=K 的解,dx等于当前 x 与下一组解 x 的差 (下一组解x+abs(B/d),y-abs(A/d)) // x=(x%dx+dx)%dx 是常用的把 x变成最小正整数解的方式,如果还不懂可以手模一下 return x; } LL fac(LL n, LL P, LL Pk) // 算 n! 不含 P 的部分 (例如 P=3,则 fac(9)=1*2*(1)*4*5*(2)*7*8*(1)) { LL ans=1; if(n==0) return 1; //ans初始化为 1,n=0 值为 1 for(LL i=1; i<Pk; i++) if(i%P!=0) ans=(ans*i)%Pk; //ans等于 1 ~ Pk 里不含 P 的数的积 ans=q_pow(ans, n/Pk, Pk); // ans^(n/Pk) % Pk, ans等于 1 ~ n/Pk*Pk 里不含 P 的数的积 for(LL i=1; i<=n%Pk; i++) if(i%P!=0) ans=(ans*i)%Pk; //ans等于 -
0
by hansang:
/* by;hansang 求的是 C nm %P (P不一定是质数),把 C nm 设为 x 把 P 化为 p1^k1*p2^k2*... (pi是质数) 列出式子: x % (p1^k1) = a1 x % (p2^k2) = a2 ...... 用中国剩余定理合并可得答案 */ #include<bits/stdc++.h> using namespace std; typedef long long LL; const int N=1e6+10; LL a[N], c[N]; int cnt; //a 数组表示 x (C nm ) % (pi^ki) 的值,c 数组表示 pi^ki,cnt 表示 a c 数组长度 LL q_pow(LL a, LL b, LL P) //求 a ^b %P { LL ans=1%P; a%=P; //防止 P 为 1,a 先%P for(; b; b>>=1) //枚举 b 的二进制位 { if(b&1) ans=ans*a%P; //b 这个位为 1,ans 就乘这个位的 a a=a*a%P; // a=a^2 } return ans; } void exgcd(LL a, LL b, LL &d, LL &x, LL &y) //扩展欧几里得 { if(b==0) {d=a; x=1; y=0;} //找到最大公约数了,d=a else { exgcd(b, a%b, d, y, x); y-=(a/b)*x; /* 设新的 x=dx,y=dy 则有; ax +by = b*tx +(a%b)*ty 已知; a%b = a-[a/b]*b ([]为下取整) 则; ax +by = b*tx +(a%b)*ty = b*tx +(a-[a/b]*b)*ty = b*tx +a*ty -[a/b] *b*ty = a*ty +b*(tx-[a/b]*ty) 所以; ty=x,tx=(y+[a/b]*x) 倒推回来 y 就要减 (a/b)*x */ } } // exgcd 和 q_pow 都为很重要的板子,请熟背 LL inv(LL a, LL P) //逆元 (运用 exgcd) { LL A=a, B=P, K=1, d, x, y; // ax=1 (%P) = ax+P*?=1 = Ax+By=K exgcd(A, B, d, x, y); //算的是 Ax+By=d (d为 gcd(a, b)) x=x*(K/d); LL dx=abs(B/d); x=(x%dx+dx)%dx; //把 x变成最小正整数解 // x*(K/d)变成 Ax+By=K 的解,dx等于当前 x 与下一组解 x 的差 (下一组解x+abs(B/d),y-abs(A/d)) // x=(x%dx+dx)%dx 是常用的把 x变成最小正整数解的方式,如果还不懂可以手模一下 return x; } LL fac(LL n, LL P, LL Pk) // 算 n! 不含 P 的部分 (例如 P=3,则 fac(9)=1*2*(1)*4*5*(2)*7*8*(1)) { LL ans=1; if(n==0) return 1; //ans初始化为 1,n=0 值为 1 for(LL i=1; i<Pk; i++) if(i%P!=0) ans=(ans*i)%Pk; //ans等于 1 ~ Pk 里不含 P 的数的积 ans=q_pow(ans, n/Pk, Pk); // ans^(n/Pk) % Pk, ans等于 1 ~ n/Pk*Pk 里不含 P 的数的积 for(LL i=1; i<=n%Pk; i++) if(i%P!=0) ans=(ans*i)%Pk; //ans等于 1 ~ n 里不含 P 的数的积 return ans*fac(n/P, P, Pk)%Pk; //继续递归,因为有些数 = x* P^i,x 也要算进来 } LL C(LL n, LL m, LL P, LL Pk) // C n m { if(n<m) return 0; // n<m 值为 0 LL f1=fac(n, P, Pk), f2=fac(m, P, Pk), f3=fac(n-m, P, Pk), sum=0; // f1表示 n!中不含 P的部分, f2表示 m!中不含 P的部分,f3表示 (n-m)!中不含 P的部分, for(LL i=n; i; i/=P) sum+=i/P; //sum加上 n 里面有多少个 P for(LL i=m; i; i/=P) sum-=i/P; //sum减去 m 里面有多少个 P for(LL i=n-m; i; i/=P) sum-=i/P; //sum减去 (n-m) 里面有多少个 P return f1*inv(f2, Pk)%Pk*inv(f3, Pk)%Pk*q_pow(P, sum, Pk)%Pk; // f1乘上 f2的逆元 (%Pk) 乘上 f3的逆元 (%Pk) 再乘上 P^sum (%Pk) 就是 C n m } LL CRT() //中国剩余定理 { LL m=1, ans=0; // m 表示所有 (pi^ki) 的最小公倍数 (就是 P),ans初始化为 0 for(int i=1; i<=cnt; i++) m*=c[i]; for(int i=1; i<=cnt; i++) { ans=(ans+a[i]*(m/c[i])%m*inv(m/c[i], c[i])%m)%m; // (m/c[i]) 是除 c[i] 其它所有 (pi^ki) 的公倍数,加上对其它 c的取模结果不影响 // inv(m/c[i], c[i]) 是 (m/c[i]) 的逆元 (%c[i]),乘以 (m/c[i]) %c[i] 等于 1 // 再乘以 a[i] 就达成了不改变其它 c的取模结果,让 %c[i] 的结果等于 a[i] // %m 防止溢出,加减 m 对所有 c的取模结果没影响 (因为是公倍数) } return ans; } LL exLucas(LL n, LL m, LL P) //拓展卢卡斯 { cnt=0; //cnt 初始化 for(int i=2; i*i<=P; i++) if(P%i==0) //找出组成 P 的质数 { LL tmp=1; while(P%i==0) {tmp*=i; P/=i;} //把 P中的 i提取到 tmp里 cnt++; a[cnt]=C(n, m, i, tmp); c[cnt]=tmp; //存到 a c 数组里 } if(P>1) {cnt++; a[cnt]=C(n, m, P, P); c[cnt]=P;} //当 P 还剩有质数时 return CRT(); } int main() { LL n, m, P; scanf("%lld%lld%lld", &n, &m, &P); printf("%lld\n", exLucas(n, m, P)); return 0; }
by scy:#include<bits/stdc++.h> using namespace std; typedef long long LL; const int N=110; LL cnt,a[N],p[N],pk[N]; void exgcd(LL a,LL b,LL &d,LL &x,LL &y) { if(b==0){d=a;x=1;y=0;} else { exgcd(b,a%b,d,y,x); y-=(a/b)*x; } } LL qpow(LL a,LL b,LL P) { LL ans=1%P;a%=P; for(;b;b>>=1) { if(b&1)ans=ans*a%P; a=a*a%P; } return ans; } LL inv(LL a,LL P) { LL A=a,B=P,x,y,d;exgcd(A,B,d,x,y); LL dx=abs(B/d); x=(x%dx+dx)%dx; return x; } LL fac(LL n,LL p,LL pk) { if(n==0) return 1; LL ans=1;for(LL i=1;i<pk;i++)if(i%p!=0) ans=ans*i%pk; ans=qpow(ans,n/pk,pk); for(int i=1;i<=n%pk;i++)if(i%p!=0) ans=ans*i%pk; return ans*fac(n/p,p,pk) %pk; } LL C(LL n,LL m,LL p,LL pk) { if(n<m) return 0; LL f1=fac(n,p,pk),f2=fac(m,p,pk),f3=fac(n-m,p,pk); LL xyz=0; for(LL i=n;i;i/=p)xyz+=i/p; for(LL i=m;i;i/=p)xyz-=i/p; for(LL i=n-m;i;i/=p)xyz-=i/p; return f1*inv(f2,pk)%pk*inv(f3,pk)%pk * qpow(p,xyz,pk) %pk; } LL exlucas(LL n,LL m,LL P) { cnt=0;LL tp=P; for(LL i=2;i*i<=tp;i++) if(tp%i==0) { cnt++;p[cnt]=i;pk[cnt]=1;while(tp%i==0)tp/=i,pk[cnt]*=i; } if(tp>1){cnt++;p[cnt]=pk[cnt]=tp;};</p>for(int i=1;i<=cnt;i++) a[i]=C(n,m,p[i],pk[i]); LL x=0; for(int i=1;i<=cnt;i++)x=( x+ a[i] * (P/pk[i]) % P * inv(P/pk[i], pk[i]) % P )%P; return x;} int main() { LL n,m,p;scanf("%lld%lld%lld",&n,&m,&p); printf("%lld\n",exlucas(n,m,p)); }
- 1
信息
- ID
- 379
- 时间
- 5000ms
- 内存
- 512MiB
- 难度
- 7
- 标签
- 递交数
- 17
- 已通过
- 8
- 上传者