BZOJ 2005 能量采集 （欧拉函数）

转载至http://blog.csdn.net/clove_unique

题目链接

https://vjudge.net/problem/HYSBZ-2005

题解

首先证明对于某个点（x,y）,k=gcd(x,y)-1： 
设gcd(x,y)=t，令x=at，y=bt，那么在这条直线上的整数点可以表示为(a,b)(2a,2b)(3a,3b)……(x,y)，由于不算x,y，则答案为gcd(x,y)-1 
那么总损耗2k+1=2×gcd(x,y)-1。 
我们最终要求的式子为： 
∑i=1n∑j=1m(gcd(i,j)∗2−1) 
=2∗∑i=1n∑j=1mgcd(i,j)−n∗m 
那么我们只需要算出∑i=1n∑j=1mgcd(i,j)这个式子就可以了 
推导如下： 
∑i=1n∑j=1mgcd(i,j) 
=∑i=1n∑j=1m∑d|gcd(i,j)ϕ(d) 
=∑i=1n∑j=1m∑d=1n[d|i][d|j]ϕ(d) 
=∑d=1n∑i=1n[d|i]∑j=1m[d|j]ϕ(d) 
=∑d=1n⌊nd⌋⌊md⌋ϕ(d)

实际上⌊nd⌋⌊md⌋只有(n√+m−−√)个取值。 
可以用分块来求。 
需要预处理phi的前缀和。

#include<iostream>#include<cstring>#include<cstdio>using namespace std;#define LL long longconst int max_n=1e5+5;LL n,m,ans;LL p[max_n],phi[max_n],prime[max_n];inline void get_phi(){    phi[1]=1;    for (int i=2;i<=n;++i){        if (!p[i]){            prime[++prime[0]]=i;            phi[i]=i-1;        }        for (int j=1;j<=prime[0]&&i*prime[j]<=n;++j){            p[i*prime[j]]=1;            if (i%prime[j]==0){                phi[i*prime[j]]=phi[i]*prime[j];                break;            }            else              phi[i*prime[j]]=phi[i]*(prime[j]-1);        }        phi[i]+=phi[i-1];    }}int main(){    scanf("%lld%lld",&n,&m);    if (n>m) swap(n,m);    get_phi();    for (LL i=1,j;i<=n;i=j+1){        j=min(n/(n/i),m/(m/i));        ans+=(LL)(phi[j]-phi[i-1])*(n/i)*(m/i);    }    printf("%lld\n",ans*2-n*m);}

另外，本题还有一种程序更简便的解法。令f(k)为gcd(i, j) == k的数对的个数，令F(k)为gcd(i, j)为k的倍数的数对的个数。则F(k) = ∑f(i*k)，其中 i*k <= min(n, m)，另外也有F(k) = (n/k) * (m/k)，所以f(k) = (n/k)*(m/k) - f(2*k) - f(3*k) - f(4*k) ... ... ，逆推出来即可。

#include<cstdio>#include<cstring>#include<string>#include<cctype>#include<iostream>#include<set>#include<map>#include<cmath>#include<sstream>#include<vector>#include<stack>#include<queue>#include<algorithm>#define fin(a) freopen("a.txt","r",stdin)#define fout(a) freopen("a.txt","w",stdout)typedef long long LL;using namespace std;const int INF = 1e8 + 10;const int maxn = 1e5 + 10;LL f[maxn]; //f[i] : gcd(x, y) == i的数的个数int main() {   int n, m;   scanf("%d%d", &n, &m);   int t = min(n, m);   LL ans = 0;   for(int i = t; i >= 1; i--) {      f[i] = (LL)(m/i)*(n/i);      for(int j = i+i; j <= t; j += i)         f[i] -= f[j];      ans += f[i] * (2*i-1);   }   printf("%lld\n", ans);   return 0;}