【51Nod1584】加权约数和（数论）
题面
题解
要求的是\[\sum_{i=1}^n\sum_{j=1}^n max(i,j)\sigma(ij)\]
这个\(max\)太讨厌了，直接枚举一半乘个二。
\[2\sum_{i=1}^n\sum_{j=1}^{i}i\sigma(ij)-\sum_{i=1}^ni\sigma(i^2)\]
后面这一半可以直接预处理，只需要把\(i\)分解，可以做到调和级数的复杂度。
只考虑前面这一半，显然只需要考虑的是\(\sigma(ij)\)这个东西。
那么我们考虑在\(i\)中枚举一个约数，在\(j\)中枚举一个约数，然后把这两个约数合并一下，看看能不能让每个约数只被计算一次。
\[\sigma(ij)=\sum_{u|i}\sum_{v|j}[gcd(u,\frac{j}{v})=1]uv\]
证明的话，大概就是我们的目标是让每个约数只被计算一次，首先在\(i\)中枚举一个约数肯定没有问题，在\(j\)中枚举一个质因数也没有问题。对于\(uv\)这个数而言，我们把只在\(i\)中有的因子和只在\(j\)中有的因子给丢掉，只考虑在\(i,j\)中都含有的因子\(u',v'\)，对于一个数\(uv\)而言，可能算重的情况是\(u'\)从\(v'\)那里抢走了一个质因子，而此时\(\frac{j}{v}\)就会对应的乘上那个质因子，使得\(gcd\neq 1\)，所以每个数只会被计算一次。
有了这个式子就很好搞了，首先把这个式子换一个形式：
\[\sigma(ij)=\sum_{u|i}\sum_{v|j}[gcd(u,v)=1]\frac{uj}{v}\]
带回去得到：
\[\sum_{i=1}^n i\sum_{j=1}^i\sum_{u|i}\sum_{v|j}[gcd(u,v)=1]\frac{uj}{v}\]
考虑对于每一个\(i\)分别计算答案，所以我们设
\[\begin{aligned} f[n]&=n\sum_{j=1}^n\sum_{u|n}\sum_{v|j}[gcd(u,v)=1]\frac{uj}{v}\\ &=n\sum_{j=1}^n \sum_{u|n}\sum_{v|j}\frac{uj}{v}\sum_{k|u,k|v}\mu(k)\\ &=n\sum_{j=1}^n\sum_{k|n,k|j}\mu(k)\sum_{k|u,u|n}\sum_{k|v,v|j}\frac{uj}{v}\\ &=n\sum_{j=1}^n\sum_{k|n,k|j}\mu(k)(k\sigma(\frac{n}{k}))\sigma(\frac{j}{k})\\ &=n\sum_{k|n}\mu(k)k\sigma(\frac{n}{k})\sum_{i=1}^{n/k}\sigma(i) \end{aligned}\]
然后就是前缀和计算就行了。
所有东西可以线性筛，中间要求逆就直接快速幂了。。。。
#include<iostream>
#include<cstdio>
using namespace std;
#define ll long long
#define MOD 1000000007
#define MAX 1000100
inline int read()
{
int x=0;bool t=false;char ch=getchar();
while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();
if(ch=='-')t=true,ch=getchar();
while(ch<='9'&&ch>='0')x=x*10+ch-48,ch=getchar();
return t?-x:x;
}
int fpow(int a,int b){int s=1;while(b){if(b&1)s=1ll*s*a%MOD;a=1ll*a*a%MOD;b>>=1;}return s;}
bool zs[MAX];
int pri[MAX],tot;
int mu[MAX],sig[MAX],ssig[MAX],pw[MAX],spw[MAX],dpw[MAX],dspw[MAX],dsig[MAX];
int f[MAX],s[MAX];
void Sieve(int n)
{
mu[1]=1;dsig[1]=sig[1]=1;
for(int i=2;i<=n;++i)
{
if(!zs[i])
{
pri[++tot]=i,mu[i]=MOD-1;
sig[i]=i+1,pw[i]=i,spw[i]=i+1;
dsig[i]=dspw[i]=(1+i+1ll*i*i)%MOD;dpw[i]=1ll*i*i%MOD;
}
for(int j=1;j<=tot&&i*pri[j]<=n;++j)
{
zs[i*pri[j]]=true;
if(i%pri[j])
{
mu[i*pri[j]]=MOD-mu[i];
sig[i*pri[j]]=1ll*sig[i]*sig[pri[j]]%MOD;
pw[i*pri[j]]=pri[j],spw[i*pri[j]]=1+pri[j];
dsig[i*pri[j]]=1ll*dsig[i]*dsig[pri[j]]%MOD;
dpw[i*pri[j]]=dpw[pri[j]];
dspw[i*pri[j]]=dspw[pri[j]];
}
else
{
mu[i*pri[j]]=0;
sig[i*pri[j]]=1ll*sig[i]*fpow(spw[i],MOD-2)%MOD*(spw[i]+pw[i]*pri[j])%MOD;
pw[i*pri[j]]=pw[i]*pri[j];
spw[i*pri[j]]=(spw[i]+pw[i]*pri[j])%MOD;
dspw[i*pri[j]]=(dspw[i]+1ll*dpw[i]*pri[j]%MOD+1ll*dpw[i]*pri[j]%MOD*pri[j]%MOD)%MOD;
dsig[i*pri[j]]=1ll*dsig[i]*fpow(dspw[i],MOD-2)%MOD*dspw[i*pri[j]]%MOD;
dpw[i*pri[j]]=1ll*dpw[i]*pri[j]%MOD*pri[j]%MOD;
break;
}
}
}
for(int i=1;i<=n;++i)dsig[i]=1ll*dsig[i]%MOD*i%MOD;
for(int i=1;i<=n;++i)ssig[i]=(ssig[i-1]+sig[i])%MOD;
for(int i=1;i<=n;++i)
{
if(mu[i])
for(int j=i;j<=n;j+=i)
f[j]=(f[j]+1ll*mu[i]*i%MOD*sig[j/i]%MOD*ssig[j/i])%MOD;
f[i]=1ll*f[i]*i%MOD;s[i]=(s[i-1]+2ll*f[i]+MOD-dsig[i])%MOD;
}
}
int main()
{
Sieve(MAX-1);
int T=read();
for(int i=1;i<=T;++i)
printf("Case #%d: %d\n",i,s[read()]);
return 0;
}