51nod 1584 加权约数和

题意

求：(sumlimits_{i=1}^nsumlimits_{j=1}^nmax(i,j)*f(ij))，(f(ij))表示(ij)的约数和。

先把(max)去掉：
(2*sumlimits_{i=1}^{n}sumlimits_{j=1}^ii*f(ij)-sumlimits_{i=1}^ni*f(i^2))
后面那个是积性函数，显然可以线性筛，接下来考虑前面那个。
(sumlimits_{i=1}^{n}sumlimits_{j=1}^{i}i*f(ij))
首先考虑(f(ij))
(f(ij)=sumlimits_{x|i}sumlimits_{y|j}xy*[gcd(x,frac{j}{y})=1])
证明：
考虑从(i,j)中分别枚举个因子乘起来，考虑怎么不算重。对于(i,j)都含有的质因子(p)，假设约数中为(p^k)，我们让(y)先有，(y)不够再考虑(x)，也就是说要么(x)没有(p)，要么(y)中的(p)已经取完了，这样就可以保证不重。
代回去：
(sumlimits_{i=1}^{n}isumlimits_{j=1}^{i}sumlimits_{x|i}sumlimits_{y|j}frac{xj}{y}[gcd(x,y)=1])
设(f(n)=nsumlimits_{i=1}^nsumlimits_{x|n}sumlimits_{y|i}frac{xi}{y}[gcd(x,y)=1])
对它进行化简：
(nsumlimits_{i=1}^nsumlimits_{x|n}sumlimits_{y|i}frac{xi}{y}[gcd(x,y)=1])
(nsumlimits_{i=1}^nsumlimits_{x|n}sumlimits_{y|i}frac{xi}{y}sumlimits_{k|n,k|i}mu(k))
(nsumlimits_{i=1}^nsumlimits_{k|n,k|i}mu(k)sumlimits_{k|x,x|n}sumlimits_{k|y,y|i}frac{xi}{y})
(nsumlimits_{i=1}^nsumlimits_{k|n,k|i}mu(k)sumlimits_{x|frac{n}{k}}sumlimits_{y|frac{i}{k}}frac{xi}{y})
(nsumlimits_{i=1}^nsumlimits_{k|n,k|i}mu(k)f(frac{n}{k})sumlimits_{y|frac{i}{k}}frac{i}{y})
(nsumlimits_{i=1}^nsumlimits_{k|n,k|i}mu(k)f(frac{n}{k})ksumlimits_{y|frac{i}{k}}frac{i}{yk})
(nsumlimits_{i=1}^nsumlimits_{k|n,k|i}mu(k)f(frac{n}{k})kf(frac{i}{k}))
(nsumlimits_{k|n}mu(k)kf(frac{n}{k})sumlimits_{i=1}^{frac{n}{k}}f(i))
现在就都能算了。
code:

#include<bits/stdc++.h>
using namespace std;
#define int long long
const int maxn=1000010;
const int mod=1000000007;
int T,n;
int ans[maxn],mu[maxn],f[maxn],g[maxn],vf[maxn],vg[maxn],sumf[maxn],sumvf[maxn],sumvg[maxn],h[maxn];
bool vis[maxn];
vector<int>prime;
inline int read()
{
	char c=getchar();int res=0,f=1;
	while(c<'0'||c>'9'){if(c=='-')f=-1;c=getchar();}
	while(c>='0'&&c<='9')res=res*10+c-'0',c=getchar();
	return res*f;
}
inline int power(int x,int k,int mod)
{
	int res=1;
	while(k)
	{
		if(k&1)res=res*x%mod;
		x=x*x%mod;k>>=1;
	}
	return res;
}
inline void pre_work(int n)
{
	f[1]=g[1]=mu[1]=1;vis[1]=1;
	for(int i=2;i<=n;i++)
	{
		if(!vis[i])
		{
			prime.push_back(i);
			mu[i]=mod-1;
			f[i]=i+1,vf[i]=i,sumvf[i]=i+1;
			g[i]=(1+i+i*i%mod)%mod,vg[i]=i*i%mod,sumvg[i]=(1+i+i*i%mod)%mod;			
		}
		for(unsigned int j=0;j<prime.size()&&i*prime[j]<=n;j++)
		{
			vis[i*prime[j]]=1;
			if(i%prime[j]==0)
			{
				mu[i*prime[j]]=0;
				vf[i*prime[j]]=vf[i]*prime[j]%mod,sumvf[i*prime[j]]=(sumvf[i]+vf[i]*prime[j]%mod)%mod;
				f[i*prime[j]]=f[i]*power(sumvf[i],mod-2,mod)%mod*sumvf[i*prime[j]]%mod;
				vg[i*prime[j]]=vg[i]*prime[j]%mod*prime[j]%mod;
				sumvg[i*prime[j]]=(sumvg[i]+vg[i]*prime[j]%mod+vg[i]*prime[j]%mod*prime[j]%mod)%mod;
				g[i*prime[j]]=g[i]*power(sumvg[i],mod-2,mod)%mod*sumvg[i*prime[j]]%mod;
				break;
			}
			mu[i*prime[j]]=mod-mu[i];
			vf[i*prime[j]]=prime[j],sumvf[i*prime[j]]=prime[j]+1;
			f[i*prime[j]]=f[i]*f[prime[j]];
			vg[i*prime[j]]=vg[prime[j]]%mod;
			sumvg[i*prime[j]]=sumvg[prime[j]]%mod;
			g[i*prime[j]]=g[i]*g[prime[j]]%mod;
		}
	}
	for(int i=1;i<=n;i++)sumf[i]=(sumf[i-1]+f[i])%mod;
	for(int i=1;i<=n;i++)g[i]=g[i]*i%mod;
	for(int i=1;i<=n;i++)
	{
		if(mu[i])
			for(int j=1;j*i<=n;j++)
				h[i*j]=(h[i*j]+mu[i]*i%mod*f[j]%mod*sumf[j]%mod)%mod;
		h[i]=h[i]*i%mod;ans[i]=(ans[i-1]+2*h[i]%mod+mod-g[i])%mod;
	}
}
signed main()
{
	pre_work(1000000);
	scanf("%lld",&T);
    for(int i=1;i<=T;i++)
        printf("Case #%lld: %lld
",i,ans[read()]);
	return 0;
}