首页 > 代码库 > HDU 1695

HDU 1695

2024-07-22 13:56:20   218人阅读

看见别人的用的莫比乌斯来做，我看了好久也没明白，实在佩服，看到是组合数学的内容，只好先留着，待我学了组合数学后再用莫比乌斯来写。

 

求GCD(X,Y)=K.其实即是在[1,X/K]和[1,Y/K]的区间内求GCD(X,Y)=1的对数。这样，假设X/K<=Y/K，在[1,X/K]的区间内可以用欧拉函数求出对数，这是显然的。

而在[X/K+1,Y/K]的区间内，则可以用容斥原理求出不与区间内的数互质的对数，再后来减去即可。

怎么样去求区间内不互质的数呢，设一个P在[X/K+1,Y/K]，分解质因数为P1*P2*P3.....*PS，然后用对[1,X/K]应用容斥原理即可。

#include <iostream>#include <cstdio>#include <algorithm>#include <cstring>#include <cmath>#define LL __int64using namespace std;const int Max=1000005;bool isprime[Max];int prime[Max],np;int phi[Max];int num[Max],no;LL all,ans;void exchange(int &a,int &b){    if(a>b){        int tmp=a;        a=b; b=tmp;    }}void do_prime(){    np=0;    memset(isprime,true, sizeof(isprime));    isprime[1]=false;    for(LL i=2;i<Max;i++){        if(isprime[i]){            prime[np++]=i;            for(LL j=i*i;j<Max;j+=i){                isprime[j]=false;            }        }    }}void do_phi(){    for(int i=1;i<Max;i++) phi[i]=i;    for(int i=2;i<Max;i+=2) phi[i]/=2;    for(int i=3;i<Max;i+=2){        if(phi[i]==i){            for(int j=i;j<Max;j+=i)            phi[j]=phi[j]/i*(i-1);        }    }}void initial(){    do_prime();    do_phi();}void dfs(int i,int nu,int x,int mu,int b){    if(nu==x){        all+=(LL)(b/mu); return;    }    if(i==no) return ;    dfs(i+1,nu+1,x,mu*num[i],b);    dfs(i+1,nu,x,mu,b);}int Rong(int b,int d){    int s=0;    for(int i=1;i<=no;i++){        all=0;        dfs(0,0,i,1,b);        if(i&1) s+=(int)all;        else s-=(int)all;    }    return b-s;}int main(){    int T,a,b,c,d,k;    scanf("%d",&T);    initial();    int kase=0;    while(T--){        scanf("%d%d%d%d%d",&a,&b,&c,&d,&k);        if(k>b||k>d||k==0){            printf("Case %d: 0\n",++kase);             continue;        }        exchange(b,d);        b/=k; d/=k;        all=ans=0;        for(int i=1;i<=b;i++)        ans+=(LL)phi[i];        for(int i=b+1;i<=d;i++){            no=0;            int tmp=i;            for(int j=0;(LL)prime[j]*(LL)prime[j]<=(LL)tmp&&j<np;j++){                if(tmp%prime[j]==0){                    num[no++]=prime[j];                    while(tmp%prime[j]==0){                        tmp/=prime[j];                    }                }            }            if(tmp>1)            num[no++]=tmp;            ans+=(LL)Rong(b,d);        }        printf("Case %d: %I64d\n",++kase,ans);    }    return 0;}

　　

HDU 1695


联系
我们