zoukankan      html  css  js  c++  java
  • BZOJ3129/洛谷P3301方程(SDOI2013)容斥原理+扩展Lucas定理

    题意:给定方程x1+x2+....xn=m,每个x是正整数。但是对前n1个数做了限制x1<=a1,x2<=a2...xn1<=an1,同时对第n1+1到n1+n2个数也做了限制xn1+1>=an1+1....xn1+n2>=an1+n2,输出方程解个数。

    解法:首先如果对数字没有任何要求(应该是只要求是非负数)的话,答案就是C(n+m+1,m+1)原理是隔板法。但是此题有各种限制,我们想办法解决限制使得答案往无限制上面靠。

    首先是解决要正整数,那么每个数字减一即可,就是m-=n。

    然后对于n1+1到n1+n2的限制,大于等于的限制也简单,也是xi>=ai的话就使xi-=ai即可。

    好了到最麻烦的xi<=ai限制,这个限制我们没办法通过数字处理来解决。但是我们观察到这些限制个数小于等于8个,自然会想到容斥原理。

    那么容斥原理常规套路,无限制-至少一个突破限制(也就是xi>ai为突破限制)方案数+至少两个数字突破限制.........

    ok,解题思路就是上面这些。但是这题比较麻烦的在于:要计算的组合数非常大且模数不为质数!!!那么Lucas定理也不管用了,我们只能用到exLucas定理。

    推荐看这篇博客:https://www.cnblogs.com/candy99/p/6637629.html

    直接套上大佬的exLucas会TLE,有一些地方需要优化,①Lucas函数内每次都要找MOD的因子太慢了,先预处理出来。②Fac函数内每次暴力计算阶乘太慢了,需要在调用Fac之前即C函数处先预处理阶乘。

    此题代码(注释的是修改大佬模板地方,加上优化):

     1 #include <bits/stdc++.h>
     2 using namespace std;
     3 typedef long long ll;
     4 const int N=1e5+10;
     5 ll n,m,n1,n2,MOD,a[N];
     6 ll tot=0,f[N],c[N];
     7 
     8 ll fac[N];
     9 void Initfac(ll p,ll pr) {  //预处理阶乘 
    10     fac[0]=1;
    11     for(ll i=1;i<=pr;i++) 
    12         if(i%p) fac[i]=fac[i-1]*i%pr;
    13         else fac[i]=fac[i-1];
    14 }
    15 
    16 ll Pow(ll a,ll b,ll P){
    17     ll ans=1;
    18     for(;b;b>>=1,a=a*a%P)
    19         if(b&1) ans=ans*a%P;
    20     return ans;
    21 }
    22 void exgcd(ll a,ll b,ll &d,ll &x,ll &y){
    23     if(b==0) d=a,x=1,y=0;
    24     else exgcd(b,a%b,d,y,x),y-=(a/b)*x;
    25 }
    26 ll Inv(ll a,ll n){
    27     ll d,x,y;
    28     exgcd(a,n,d,x,y);
    29     return d==1?(x+n)%n:-1;
    30 }
    31 ll Fac(ll n,ll p,ll pr){
    32     if(n==0) return 1;
    33     ll re=1;
    34     //for(ll i=2;i<=pr;i++) if(i%p) re=re*i%pr;
    35     re=fac[pr];
    36     re=Pow(re,n/pr,pr);
    37     ll r=n%pr;
    38     //for(int i=2;i<=r;i++) if(i%p) re=re*i%pr;
    39     re=re*fac[r]%pr;
    40     return re*Fac(n/p,p,pr)%pr;
    41 }
    42 ll C(ll n,ll m,ll p,ll pr){
    43     if(n<m) return 0;
    44     Initfac(p,pr);  //预处理阶乘 
    45     ll x=Fac(n,p,pr),y=Fac(m,p,pr),z=Fac(n-m,p,pr);
    46     ll c=0;
    47     for(ll i=n;i;i/=p) c+=i/p;
    48     for(ll i=m;i;i/=p) c-=i/p;
    49     for(ll i=n-m;i;i/=p) c-=i/p;
    50     ll a=x*Inv(y,pr)%pr*Inv(z,pr)%pr*Pow(p,c,pr)%pr;
    51     return a*(MOD/pr)%MOD*Inv(MOD/pr,pr)%MOD;
    52 }
    53 ll Lucas(ll n,ll m){   //exLucas 
    54     ll x=MOD,re=0;
    55 //    for(ll i=2;i<=MOD;i++) if(x%i==0){
    56 //        ll pr=1;
    57 //        while(x%i==0) x/=i,pr*=i;
    58 //        re=(re+C(n,m,i,pr))%MOD;
    59 //    }
    60     for (int i=1;i<=tot;i++) re=(re+C(n,m,f[i],c[i]))%MOD;
    61     return re;
    62 }
    63 
    64 void prework(ll n) {
    65     ll t=n;
    66     for (int i=2;(ll)i*i<=n;i++) {
    67         if (t%i==0) {
    68             tot++; f[tot]=i; c[tot]=1;
    69             while (t%i==0) t/=i,c[tot]*=i;
    70         }
    71     }
    72     if (t>1) f[++tot]=t,c[tot]=t;
    73 }
    74 
    75 int main()
    76 {
    77     int T; scanf("%d%lld",&T,&MOD);
    78     prework(MOD);
    79     while(T--) {
    80         scanf("%lld%d%d%lld",&n,&n1,&n2,&m);
    81         for(int i=1;i<=n1+n2;i++) scanf("%lld",&a[i]);
    82         for(int i=1;i<=n2;i++) m-=a[n1+i]-1;
    83         m-=n;
    84 
    85         ll ans=0;
    86         for(int i=0;i<(1<<n1);i++) {
    87             int tot=0;
    88             ll x=m;
    89             for(int j=0;j<n1;j++)
    90                 if (i&(1<<j)) { tot++; x-=a[j+1]; }
    91             if (tot%2) ans=(ans-Lucas(x+n-1,n-1)+MOD)%MOD;
    92             else ans=(ans+Lucas(x+n-1,n-1))%MOD;
    93         }
    94         printf("%lld
    ",ans);
    95     }
    96     return 0;
    97 }
    View Code

    记录下大佬的exLucas定理模板:

     1 ll Pow(ll a,ll b,ll P){
     2     ll ans=1;
     3     for(;b;b>>=1,a=a*a%P)
     4         if(b&1) ans=ans*a%P;
     5     return ans;
     6 }
     7 void exgcd(ll a,ll b,ll &d,ll &x,ll &y){
     8     if(b==0) d=a,x=1,y=0;
     9     else exgcd(b,a%b,d,y,x),y-=(a/b)*x;
    10 }
    11 ll Inv(ll a,ll n){
    12     ll d,x,y;
    13     exgcd(a,n,d,x,y);
    14     return d==1?(x+n)%n:-1;
    15 }
    16 ll Fac(ll n,ll p,ll pr){
    17     if(n==0) return 1;
    18     ll re=1;
    19     for(ll i=2;i<=pr;i++) if(i%p) re=re*i%pr;
    20     re=Pow(re,n/pr,pr);
    21     ll r=n%pr;
    22     for(int i=2;i<=r;i++) if(i%p) re=re*i%pr;
    23     return re*Fac(n/p,p,pr)%pr;
    24 }
    25 ll C(ll n,ll m,ll p,ll pr){
    26     if(n<m) return 0;
    27     ll x=Fac(n,p,pr),y=Fac(m,p,pr),z=Fac(n-m,p,pr);
    28     ll c=0;
    29     for(ll i=n;i;i/=p) c+=i/p;
    30     for(ll i=m;i;i/=p) c-=i/p;
    31     for(ll i=n-m;i;i/=p) c-=i/p;
    32     ll a=x*Inv(y,pr)%pr*Inv(z,pr)%pr*Pow(p,c,pr)%pr;
    33     return a*(MOD/pr)%MOD*Inv(MOD/pr,pr)%MOD;
    34 }
    35 ll Lucas(ll n,ll m){
    36     ll x=MOD,re=0;
    37     for(ll i=2;i<=MOD;i++) if(x%i==0){
    38         ll pr=1;
    39         while(x%i==0) x/=i,pr*=i;
    40         re=(re+C(n,m,i,pr))%MOD;
    41     }
    42     return re;
    43 }
    View Code
  • 相关阅读:
    常用 SQL 语句使用的总结
    LC 583. Delete Operation for Two Strings
    LC 873. Length of Longest Fibonacci Subsequence
    LC 869. Reordered Power of 2
    LC 900. RLE Iterator
    LC 974. Subarray Sums Divisible by K
    LC 973. K Closest Points to Origin
    LC 975. Odd Even Jump
    LC 901. Online Stock Span
    LC 722. Remove Comments
  • 原文地址:https://www.cnblogs.com/clno1/p/11662402.html
Copyright © 2011-2022 走看看