[Noi 2010][BZOJ 2005][JZOJ 2225]能量采集

description

栋栋有一块长方形的地，他在地上种了一种能量植物，这种植物可以采集太阳光的能量。在这些植物采集能量后，栋栋再使用一个能量汇集机器把这些植物采集到的能量汇集到一起。 栋栋的植物种得非常整齐，一共有n列，每列有m棵，植物的横竖间距都一样，因此对于每一棵植物，栋栋可以用一个坐标(x, y)来表示，其中x的范围是1至n，表示是在第x列，y的范围是1至m，表示是在第x列的第y棵。 由于能量汇集机器较大，不便移动，栋栋将它放在了一个角上，坐标正好是(0, 0)。 能量汇集机器在汇集的过程中有一定的能量损失。如果一棵植物与能量汇集机器连接而成的线段上有k棵植物，则能量的损失为2k + 1。例如，当能量汇集机器收集坐标为(2, 4)的植物时，由于连接线段上存在一棵植物(1, 2)，会产生3的能量损失。注意，如果一棵植物与能量汇集机器连接的线段上没有植物，则能量损失为1。现在要计算总的能量损失。 下面给出了一个能量采集的例子，其中n = 5，m = 4，一共有20棵植物，在每棵植物上标明了能量汇集机器收集它的能量时产生的能量损失。 在这个例子中，总共产生了36的能量损失。

Solution

显然的莫比乌斯反演题，
要求x,y的公共因子个数，等价于求$\gcd(x,y)-1$，
证明：
x,y的公共因子除了1以外的最小的两个就是$\frac{x}{\gcd(x,y)},\frac{y}{\gcd(x,y)}$，
第二小的就是$\frac{x}{\gcd(x,y)}*2,\frac{y}{\gcd(x,y)}*2$，
一直下去，最大的就是x,y，再减去自身
所以最多有$\gcd(x,y)-1$个，
那么题目就变成了：

$$\sum_{x=1}^{n}\sum_{y=1}^m 2*\gcd(x,y)-1$$
先求出 $ans=\sum_{x=1}^{n}\sum_{y=1}^m \gcd(x,y)$
用莫比乌斯反演的标准形式来写：
设f(d)为 $\gcd(x,y)=d$的个数，g(d)为 $d|\gcd(x,y)(i<=n)$的x,y对数，
有：
$$g(d)=\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}f(id)$$
繁衍：
$$f(d)=\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}g(id)*\mu(i)$$

我们发现，g(T)很好求：

$$g(T)=\lfloor\frac{n}{T}\rfloor*\lfloor\frac{m}{T}\rfloor$$
$$ans=\sum_{d=1}^nd*\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\mu(i)*g(id)$$
$$ans=\sum_{d=1}^nd*\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\mu(i)*\lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor$$
发现后面一块可以分块求，
又发现前面一块也可以分块求，
再仔细一想：
$n,m<=10^5$
好像暴力也可以过耶！

Code

暴力

#include<iostream>
#include<cstdio>
#include<cstdlib>
#define fo(i,a,b) for(LL i=a;i<=b;i++)
using namespace std;
typedef long long LL;
const int N=100500;
LL n,m;
LL pr[N],mu[N];
bool prz[N];
LL ans;
void pre(int n)
{
    mu[1]=1;
    fo(i,2,n)
    {
        if(!prz[i])pr[++pr[0]]=i,mu[i]=-1;
        fo(j,1,pr[0])
        {
            int t=i*pr[j];
            if(t>n)break;
            prz[t]=1;
            if(i%pr[j]==0)break;
            mu[t]=-mu[i];
        }
    }
    fo(i,1,n)mu[i]+=mu[i-1];
}
int main()
{
    scanf("%d%d",&n,&m);
    if(n>m)swap(n,m);
    pre(m);
    ans=0;
    fo(d,1,n)
    {
        LL ww=0;
        fo(i,1,n/d)ww+=(mu[i]-mu[i-1])*(n/(i*d))*(m/(i*d));
        ans+=d*ww;
    }
    printf("%lld\n",ans*2-m*n);
    return 0;
}

一个分块

#include<iostream>
#include<cstdio>
#include<cstdlib>
#define fo(i,a,b) for(LL i=a;i<=b;i++)
using namespace std;
typedef long long LL;
const int N=100500;
LL n,m;
LL pr[N],mu[N];
bool prz[N];
LL ans;
LL min(LL a,LL b){return a<b?a:b;}
void pre(int n)
{
    mu[1]=1;
    fo(i,2,n)
    {
        if(!prz[i])pr[++pr[0]]=i,mu[i]=-1;
        fo(j,1,pr[0])
        {
            int t=i*pr[j];
            if(t>n)break;
            prz[t]=1;
            if(i%pr[j]==0)break;
            mu[t]=-mu[i];
        }
    }
    fo(i,1,n)mu[i]+=mu[i-1];
}
LL fk(int d)
{
    LL ans=0;
    LL i=1,r;
    while(i<=n/d)
    {
        LL t=i*d;
        r=min(min(n/(n/t),m/(m/t))/d,n/d);
        ans+=(mu[r]-mu[i-1])*(n/t)*(m/t);
        i=r+1;
    }
    return ans;
}
int main()
{
    scanf("%d%d",&n,&m);
    if(n>m)swap(n,m);
    pre(m);
    ans=0;
    fo(d,1,n)ans+=d*fk(d);
    printf("%lld\n",ans*2-m*n);
    return 0;
}