BZOJ 3129 [Sdoi2013]方程 不定方程解的个数+组合数取模
2015-08-27 09:22
441 查看
题意:链接
方法:不定方程解的个数+组合数取模
解析:
先看n1与n2的部分的限制。
对于后半部分的限制来说,我们直接减去An1+i−1A_{n1+i}-1就可以转化一下求正整数解。
但是前半部分呢?
跟上一道猴子那个很像。
所以我们容斥搞就行了。
但是这道题好像不好写的地方不在这?
这题TMD不就是礼物吗!
大组合数取模如何取?
请参见我《BZOJ 礼物》的题解。
另外吐槽题干
明明是X1+X2+…+Xn=m
并不是小于等于
代码:
#include <cstdio> #include <cstring> #include <iostream> #include <algorithm> #define N 21 #define pa pair<int,int> #define mp make_pair #define fi first #define se second using namespace std; typedef long long ll; ll a ; ll t,p,tot; ll prime ,mod ,num ,ans ; void exgcd(ll a,ll b,ll &x,ll &y,ll &gcd) { if(!b) { x=1,y=0,gcd=a; return; } exgcd(b,a%b,y,x,gcd); y=y-a/b*x; } ll get_inv(ll n,ll m) { ll x,y,gcd; exgcd(n,m,x,y,gcd); //(x/gcd%(m/gcd)+m/gcd)%(m/gcd) //gcd==1 return (x%m+m)%m; } ll quick_my(ll x,ll y,ll MOD) { ll ret=1; while(y) { if(y&1)ret=(ret*x)%MOD; x=(x*x)%MOD; y>>=1; } return ret; } void getp() { int x=p; for(int i=2;i*i<=x;i++) { if(x%i==0) { prime[++tot]=i,mod[tot]=1,num[tot]=0; while(x%i==0) { x/=i,mod[tot]*=i,num[tot]++; } } } if(x!=1) { prime[++tot]=x,mod[tot]=x,num[tot]=1; } } ll crt() { ll ret=0; for(int i=1;i<=tot;i++) { ret=(ret+p/mod[i]*get_inv(p/mod[i],mod[i])%p*ans[i])%p; } return ret; } pa get_factor(ll x,ll pri,ll Mod) { ll ans=1; if(x==0||x==1)return mp(0,1); ll cnt_circle=x/Mod,cnt_pri=x/pri; if(cnt_circle!=0) { for(int i=2;i<Mod;i++) { if(i%pri!=0) ans=(ans*(ll)i)%Mod; } ans=quick_my(ans,cnt_circle,Mod); } for(int i=cnt_circle*Mod+1;i<=x;i++) { if(i%pri!=0) ans=(ans*(ll)i)%Mod; } pa tmp=get_factor(cnt_pri,pri,Mod); return mp(cnt_pri+tmp.fi,ans*tmp.se%Mod); } ll get_C(ll n,ll m,int k) { if(m==0)return 1; pa tmp=get_factor(n,prime[k],mod[k]),tmp1=get_factor(m,prime[k],mod[k]),tmp2=get_factor(n-m,prime[k],mod[k]); return quick_my(prime[k],tmp.fi-tmp1.fi-tmp2.fi,mod[k])*tmp.se%mod[k]*get_inv(tmp1.se,mod[k])%mod[k]*get_inv(tmp2.se,mod[k])%mod[k]; } ll calc(ll n,ll m) { if(n==m||m==0)return 1; if(n<m)return 0; for(int i=1;i<=tot;i++)ans[i]=get_C(n,m,i); return crt(); } ll ANS; ll n,n1,n2,m; void dfs(int now,ll cnt,ll sum) { if(now==n1+1) { if(cnt%2==0)ANS=(ANS+calc(m-sum-1,n-1))%p; else ANS=((ANS-calc(m-sum-1,n-1))%p+p)%p; return; } dfs(now+1,cnt,sum); dfs(now+1,cnt+1,sum+a[now]); } int main() { scanf("%lld%lld",&t,&p); getp(); while(t--) { scanf("%lld%lld%lld%lld",&n,&n1,&n2,&m); for(int i=1;i<=n1;i++)scanf("%lld",&a[i]); for(int i=1;i<=n2;i++) { ll x; scanf("%lld",&x); m=m-x+1; } if(m<=0){puts("0");continue;} ANS=0; dfs(1,0,0); printf("%lld\n",ANS); } }
相关文章推荐
- JavaScript中如何将html字符串转化为Jquery对象或者Dom对象
- 微软Nokia 222:可拍照可上网 售价37美元 32GB的microSD卡扩展
- hdu5418
- MySQL数据库连接超时(wait_timeout)问题的处理
- (转)Do not use "using" for WCF Clients - 不要将WCF Client 放在 ‘Using’ 代码块中
- c++中的static关键字
- 1高级js总结-------函数匿名自调用
- centos 查找不到网卡解决方法。
- 单例设计模式
- HDU 2767 强连通分量
- JVM学习笔记
- Eclipse中最有用的快捷键组合
- Android学习—RecyclerView的使用(1)
- CentOS安装memcached及配置php的memcache扩展
- 设计解析-封面/banner
- Java六功能开发者应避免
- 为什么要圆通?
- MySQL字符集知识点
- HTML5 CSS3:诱人的3D旋转木马效果相册实例
- button以及Imagebutton的使用