gcdlcm[组合数学]

来源:互联网 发布:长胶大联盟淘宝店铺 编辑:程序博客网 时间:2024/04/29 18:00

gcdlcm[组合数学]:

描述:
问有多少个 k 元组(a1,a2,a3,...,ak)(ai>=1)满足 gcd(a1,a2,...,ak)=D,
并且 lcm(a1,a2,a3,...,ak)=L。
输入:
第一行有三个整数 k,D,L 分别为题面所描述的。
输出:
输出答案 mod(1e9+7)
输入输出样例:
gcdlcm.in gcdlcm.out
2 1 6 4
数据范围:
对于 30%数据 2<=k<=100,D<=100,L<=100
对于 100%数据 2<=k<=1000000000 D<=1000000000 L<=1000000000


数学弱成渣了,毫无办法- -

---------------------------------------------------------------------------------------------------------------------------------------------------------

首先,D|L.  否则不合法,输出0.

然后,将D,L分解为D=p1^b1*p2^b2*p3^b3*```*pq^bq  ,L=p1^c1*p2^c2*p3^c3*````*pq^cq

然后选k个数出来,使得k个数分解出来的cnt[i][pj]满足:

cj>=(cnt[i][pj])>=bj        ,求方案数

考虑对于每个因子pi   我们要在[bi,ci]里面选k个,并且满足min(pij)=bi,max(pij)=ci.得到一个方案gi

然后把所有gi乘起来!


问题就转化为求gi了!

这里要用到数学的一个容斥:(数学太弱的我表示不会想到)

合法条件在pi中一定至少有一个选了bi,ci个.  不合法的就是 没有选bi个+没有选ci个-(没有选bi&&ci个)

我们考虑用总的减去不合法的

总的就是在[bi,ci]中选k个,为(ci-bi+1)^k

不选bi或者不选ci都是(ci-bi)^k

gi=(ci-bi+1)^k-2*(ci-bi)^k+(ci-bi-1)^k

这个用快速幂求一下就好了


复杂度其实不算很大,不过需要考虑分解出来的pq可能会非常大(如果>sqrt(n)就单独记录,因为数组开不下)

代码:

#include<cstdio>#include<iostream>#include<queue>#include<cmath>#include<algorithm>#include<cstring>#include<cstdlib>#define LL long longusing namespace std;int k;int D,L;const LL M=1e9+7;/*gcd(a1,a2,a3---an)=D=p1^b1 * p2^b2 * p3^b3//lcm(a2,a2,a3---an)=L=p1^c1 * p2^c2 * p3^c3                  pk//选法就是考虑每个因子pi 从[b1,c1]中选k个 得到方案数 gi   ans=π gi                                                            //i=1 对于每个gi,我们考虑从[b1,c1]中选k个,方案数 (c1-b1+1)^k合法状态一定有至少一个b1,c1那么不合法状态,  不选b1  ---(b1,c1]选k个,方案数  (c1-b1)^k                 不选c1  ---[b1,c1)选k个,方案数  (c1-b1)^k多减去了(b1,c1)的方案数  (c1-b1-1)^kgi=(c1-b1+1)^k-2*(c1-b1)^k+(c1-b1-1)^k; 复杂度 log max(D,L) * */ int b[1000000],c[1000000];int xx=-1;int ml=0;void fenjie(int x,int *a)//把x拆分到a数组 {memset(a,0,sizeof(a));int t=(int)sqrt(1.0*x);for(int i=2;i<=t;i++){if(x%i==0){ml=max(ml,i);while(x%i==0){x/=i;a[i]++;}}}if(x>1){if(x>1000000)xx=x;else {ml=max(ml,x);a[x]++;}}//}LL qk(LL x,LL y)//快速计算x^y mod M {LL res=1;while(y){if(y%2==1)res=((res%M)*(x%M))%M;y/=2;x=((x%M)*(x%M))%M;}return res;}LL ans;int main(){freopen("gcdlcm.in","r",stdin);freopen("gcdlcm.out","w",stdout);scanf("%d%d%d",&k,&D,&L);if(L%D)printf("0\n");else{fenjie(D,b);    fenjie(L,c);    LL g=0;    ans=1;    for(int i=2;i<=ml;i++) {if(!c[i])continue;LL d=(LL)c[i]-(LL)b[i];if(!d)g=1;else g=qk(d+1,k)%M-2*qk(d,k)%M+qk(d-1,k);g%=M;ans=((ans%M)*g)%M; }    if(xx!=-1)    {    if(D%xx!=0)    {    LL d=1;        g=qk(d+1,k)%M-2*qk(d,k)%M+qk(d-1,k);        ans=((ans%M)*(g%M))%M;    }        }    ans=(ans+M)%M;    printf("%I64d\n",ans);}return 0;}



1 0
原创粉丝点击