51nod-1227-平均最小公倍数

题意

定义 \(n\) 的平均最小公倍数：

\[A(n)=\frac{1}{n}\sum _{i=1}^n\text{lcm}(n,i) \]

求

\[\sum _{i=L}^RA(i) \]

\(n\le 10^9\) 。

分析

有趣的题，学到了一些东西。

我最开始不知道怎么都枚举gcd的时候是整除枚举，然后怎么都做不了。改求和指标为gcd的时候，直接从 1 到 \(n\) 枚举不是很正常的做法吗？于是就开始推………经过很久，答案变成了这个样子：

\[\sum _{i=1}^ng(i)f(\lfloor\frac{n}{i}\rfloor) \\ g(n)=\sum _{e|n}e\mu (e) \\ f(n)=\frac{n(n+1)(n+2)}{6} \]

由于取底全部都是从 \(n\) 出发的，而且根据底合并的性质，可以用杜教筛来快速计算 \(g\) 的前缀和，分段计算，复杂度为 \(O(n^\frac{3}{4})\) 。

然后TLE了。

网上有一篇题解说可以直接从 \(L\) 到 \(R\) 计算，不需要算两次，我不是很懂怎么做；用杜教筛只能从前缀和出发吧 。

肯定是做法有点问题。于是去找题解。这题非常简单嘛！不用像我那样。

枚举一对互质的数，再取它们的倍数，比值是定值！ 非常妙的方法啊！！！

\[\begin{aligned} \sum _{i=1}^n\sum _{j=1}^i\frac{j}{\gcd(i,j)}&=\sum _{i=1}^n\sum _{j=1}^i[\gcd(i,j)=1]\sum _{k=1}^{\lfloor\frac{n}{i}\rfloor}j \\ &=\sum _{i=1}^n\lfloor\frac{n}{i}\rfloor\sum _{j=1}^i[\gcd(i,j)=1]j \end{aligned} \]

右边的那个东西……不就是与 \(i\) 互质的数的和！由对称性，可以得到：

\[\sum _{i=1}^n[\gcd(i,n)=1]i=\frac{\varphi(n)*n+[n=1]}{2} \]

于是就可以简单地用杜教筛求 \(\varphi(n)*n+[n==1]\) 的前缀和，分段计算啦！通过预处理，复杂度为 \(O(n^\frac{2}{3})\) 。

#include<bits/stdc++.h>
#include<tr1/unordered_map>
using namespace std;
using namespace std::tr1;
typedef long long giant;
const int q=1e9+7;
inline int Plus(int x,int y) {return ((giant)x+(giant)y)%q;}
inline int Sub(int x,int y) {return Plus(x,q-y);}
inline int Multi(int x,int y) {return (giant)x*y%q;}
inline int mi(int x,int y) {
	int ret=1;
	for (;y;y>>=1,x=Multi(x,x)) if (y&1) ret=Multi(ret,x);
	return ret;
}
inline int inv(int x) {return mi(x,q-2);}
const int maxn=2e6+1;
const int i2=inv(2);
const int i6=inv(6);
bool np[maxn];
int p[maxn],ps=0,phi[maxn];
typedef unordered_map<int,int> Map;
Map PHI;
inline int sum(int n) {return Multi(Multi(n,n+1),i2);}
inline int f(int n) {return Multi(Multi(n,n+1),Multi((2ll*n+1)%q,i6));}
int g(int n) {
	if (n<maxn) return phi[n];
	Map::iterator it=PHI.find(n);
	if (it!=PHI.end()) return it->second;
	int &ret=PHI[n]=0;
	for (int i=2,j;i<=n;i=j+1) {
		j=n/(n/i);
		ret=Plus(ret,Multi(Sub(sum(j),sum(i-1)),g(n/i)));
	}
	return ret=Sub(f(n),ret);
}
int calc(int n) {
	int ret=0,last=0;
	for (int i=1,j;i<=n;i=j+1) {
		j=n/(n/i);
		int x=Plus(g(j),1);
		ret=Plus(ret,Multi(n/i,Sub(x,last)));
		last=x;
	}
	return Multi(ret,i2);
}
int main() {
#ifndef ONLINE_JUDGE
	freopen("test.in","r",stdin);
#endif
	phi[1]=1;
	for (int i=2;i<maxn;++i) {
		if (!np[i]) p[++ps]=i,phi[i]=i-1;
		for (int j=1,tmp;j<=ps && (tmp=i*p[j])<maxn;++j) {
			np[tmp]=true;
			if (i%p[j]==0) {
				phi[tmp]=Multi(phi[i],p[j]);
				break;
			}
			phi[tmp]=Multi(phi[i],phi[p[j]]);
		}
	}
	for (int i=2;i<maxn;++i) phi[i]=Plus(phi[i-1],Multi(phi[i],i));
	int L,R;
	cin>>L>>R;
	cout<<Sub(calc(R),calc(L-1))<<endl;
	return 0;
}