【数论-莫比乌斯入门】hdu 1695 GCD

莫比乌斯反演：

                                                                                      

莫比乌斯反演变形：

                               假设：n<m

                               :gcd(i,j)=d的倍数的数对的个数；

                             ：gcd(i,j)=d的数对的个数；

参考blog：https://blog.sengxian.com/algorithms/mobius-inversion-formula

kuangbin的模板代码不错！

*****************************************************************************************************************************************************************

hdu 1695 GCD

题意：给出两个区间x:(1,n)，y:(1,m);求gcd(x,y)=d的数对个数；比如d=1、(4，5)(5，4)算一对；

思路：

根据数学知识易知gcd(x,y)=d的倍数的数对个数为：

利用莫比乌斯反演，使用变形公式:

则，

那么计算f(d)的复杂度为O(n)，优化：可以对(n/dx)*(m/dx)进行分块：

找到某一段，这一段的(n/dx)都相等，(m/dx)都相等，求取μ的前缀和，批量计算这一段的答案。

复杂度变为O(sqrt(n)+sqrt(m))；

分块举例：例如，n=32,m=40,d=2； x：1~min(n,m)/d

分块代码：

last=min(n/(n/i),m/(m/i))

代码：

#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;

typedef long long ll;

const int N=100100;
bool check[N+10];
int prime[N+10];
int mu[N+10];

void Moblus()
{
    memset(check,false,sizeof(check));
    mu[1]=1;
    int tot=0;
    for(int i=2;i<=N;i++)
    {
        if(!check[i])
        {
            prime[tot++]=i;
            mu[i]=-1;
        }
        for(int j=0;j<tot;j++)
        {
            if(i*prime[j]>N) break;
            check[i*prime[j]]=true;
            if(i%prime[j]==0)
            {
                mu[i*prime[j]]=0;
                break;
            }
            else
            {
                mu[i*prime[j]]=-mu[i];
            }
        }
    }
}

int sum[N+10];
int main()
{
    int t;
    Moblus();
    sum[0]=0;
    for(int i=1;i<=N;i++)
        sum[i]=sum[i-1]+mu[i];

//    for(int i=1;i<=5;i++)
//        printf("%d\n",sum[i]);

    scanf("%d",&t);
    int co=0;
    while(t--)
    {
        int a,n,c,m,k;
        scanf("%d%d%d%d%d",&a,&n,&c,&m,&k);
        if(n>m) swap(n,m);

        if(k==0)
        {
            printf("Case %d: 0\n",++co);
            continue;
        }
        else
            n/=k,m/=k;

        ll ans1=0,ans2=0,ans=0;

        for(int i=1,la=1;i<=n;i=la+1)
        {
            la=n/(n/i);                   //分块；
            ans1+=(ll)(sum[la]-sum[i-1])*(n/i)*(n/i);
        }
        ans1/=2;
//        printf("%d\n",ans1);

        for(int i=1,la=1;i<=n;i=la+1)
        {
            la=min(n/(n/i),m/(m/i));      //分块；
            ans2+=(ll)(sum[la]-sum[i-1])*(n/i)*(m/i);
        }
//        printf("%d\n",ans2);
        ans=ans2-ans1;

        printf("Case %d: %I64d\n",++co,ans);
    }
    return 0;
}