51nod1584：加权约数和 （莫比乌斯反演）

题目传送门：http://www.51nod.com/onlineJudge/questionCode.html#!problemId=1584

---

题目分析：这题应该算是[SDOI2015]约数个数和的加强版吧，如果还没有做过这道题的同学可以先想一想这题。至于它的解法你们可以自行百度（或看看我写的题解）
好吧，回到这题。
我们要搞清楚一个东西：σ1(ij)（即i*j的约数和）究竟等于什么？（以下部分是我拿到这题时自己YY出来的）

---

在[SDOI2015]约数个数和这题中，我们有一个很巧妙的公式：σ0(ij)=∑x|i∑y|j[(x,y)=1]，它主要是由k+q+1=∑ki=0∑qj=0[(pi,pj)=1]推来的。那我们来看看它的本质是什么，我们做出i=1−>n,j=1−>n时(pi,pj)的表格：

(pi,pj) i=0 i=1 …… i=k j=q 1 p …… pmin(k,q) …… …… …… …… …… j=1 1 p …… p j=0 1 1 …… 1

很明显，我们发现表格中1的个数刚好就是k+q+1，而表格中1的要求是i或j至少有一个是0。然而这种限制是不可能用于多个质数pi合并的，所以我们才写成了∑ki=0∑qj=0[(pi,pj)=1]这种形式，这样所有的质数合并到一起时，就会变成∑x|i∑y|j[(x,y)=1]这种方便处理的形式。
现在我们要求的是σ1(ab)。如果单独考虑一个质数p，令a=pk∗……,b=pq∗……，那么我们要在答案中贡献(1+p+p2+……+pk+q)这个因子。由于列出gcd和lcm的表格，格子中都不可能出现pk+q，于是我列出了pi∗pj的表格：

pi∗pj i=0 i=1 …… i=k j=q pq pq+1 …… pk+q …… …… …… …… …… j=1 p p2 …… pk+1 j=0 1 p …… pk

这样我们就看得很清楚了，我们要求的是∑ki=0∑qj=0pi+j[i=0或j=q]，也就是表格倒”L”字形的部分。然而这样还是没办法用于多个质数合并，于是我改了改，写成∑ki=0∑qj=0pi+j[(pi,pqpj)=1]，然后多个质数合并在一起之后，就变成了最终式子：

σ1(ij)=∑x|i∑y|jxy[(x,jy)=1]

---

以上就是我在草稿纸上自己YY的过程，后来我才知道就是这样做的。
然后我们把这个式子带进原题中，基本上就求解了。
说实话，在我写这篇总结之前，网上已经有大神给出了十分优秀的题解，但鉴于这提示本人自己推导的，我还是将自己的推导过程写在下面，以作纪念吧。

---

题目要求：

∑i=1n∑j=1nmax(i,j)σ1(ij)

我们发现这个max很讨厌，于是令j只枚举到i，变形成：

2∑i=1n∑j=1iiσ1(ij)−∑i=1nσ1(i2)

其中右边是要减去重复计算的ans(i,i)，这一部分用线性筛后再做一次前缀和即可O(1)得到。线性筛时若i%prime[j]!=0，直接用积性函数的性质σ1(i∗prime[j])=σ1(i)∗σ1(prime[j])即可。若i是prime[j]的倍数，则σ1(i∗prime[j])=σ1(i)+σ1(iprime[j]h(i))∗(prime[j]h(i)+2+prime[j]h(i)+1)，其中h(i)是i的最小质因数的最高次幂，这个也可以用线性筛维护（唉，老套路了）。

---

然后考虑左边2后面那一坨，转化成：

∑i=1ni∑j=1i∑x|i∑y|jxjy[(x,y)=1]

将x,y的枚举拉到i,j前面，重新令a=ix,b=jy，整理得到：

∑x=1nx2∑a=1⌊nx⌋a∑y=1ax∑b=1⌊axy⌋b∑k|x,k|yμ(k)

然后我们将k的枚举拉到最前面，重新令i=xk,j=yk，得到：

∑k=1nk2μ(k)∑i=1⌊nk⌋i2∑a=1⌊nik⌋a∑y=1ai∑b=1⌊aiy⌋b

我们记G(n)=∑ni=1∑⌊ni⌋d=1d（很明显这就是约数和的前缀和，可以用线性筛求），F(n)=nG(n)∑d|nd（很明显，这也是可以线性求的），这样原式子就变成了：

∑k=1nk2μ(k)∑i=1⌊nk⌋F(i)

对F函数做前缀和，就可以做到预处理O(n)，询问O(n√)，然而50000组数据，还是会超。

---

式子推到上面我就推不下去了，感觉没办法让时间更优。
后来上网看了看别人的题解，原来再换个枚举顺序就可以了……（我是智障）
记d=i*k，外层枚举d，得到：

∑d=1n∑k|dkμ(k)F(dk)

然后就变成了预处理O(n*ln(n))，询问O(1)。

---

CODE：

#include<iostream>#include<string>#include<cstring>#include<cmath>#include<cstdio>#include<cstdlib>#include<stdio.h>#include<algorithm>using namespace std;const int maxn=1000010;const long long M=1000000007;typedef long long LL;int pk[maxn];LL h[maxn];LL p[maxn];LL mu[maxn];bool vis[maxn];int prime[maxn];int cur=0;LL g[maxn];LL f[maxn];int t,n;void Preparation(){    h[1]=p[1]=mu[1]=1;    for (int i=2; i<maxn; i++)    {        if (!vis[i]) pk[i]=i,h[i]=i+1,p[i]=i,p[i]*=i,p[i]=(p[i]+i+1)%M,mu[i]=-1,prime[++cur]=i;        for (int j=1; j<=cur && i*prime[j]<maxn; j++)        {            int k=i*prime[j];            vis[k]=true;            if (i%prime[j]) pk[k]=prime[j],h[k]=h[i]*h[ pk[k] ],p[k]=p[i]*p[ pk[k] ]%M,mu[k]=-mu[i];            else            {                pk[k]=pk[i]*prime[j];                h[k]=h[i]+h[ i/pk[i] ]*pk[k];                p[k]=pk[i]+pk[k],p[k]=p[k]*pk[k]%M*p[ i/pk[i] ]%M,p[k]=(p[k]+p[i])%M;                mu[k]=0;                break;            }        }    }    for (int i=1; i<maxn; i++) p[i]=p[i]*i%M;    for (int i=2; i<maxn; i++) p[i]=(p[i]+p[i-1])%M;    for (int i=1; i<maxn; i++) mu[i]=(mu[i]+M)*i%M*i%M;    g[1]=1;    for (int i=2; i<maxn; i++) g[i]=(g[i-1]+h[i])%M;    for (int i=1; i<maxn; i++) g[i]=g[i]*i%M*h[i]%M;    for (int i=1; i<maxn; i++)        for (int j=1; i*j<maxn; j++)            f[i*j]=(f[i*j]+ mu[i]*g[j]%M )%M;    for (int i=2; i<maxn; i++) f[i]=(f[i]+f[i-1])%M;}int main(){    freopen("c.in","r",stdin);    freopen("c.out","w",stdout);    Preparation();    scanf("%d",&t);    for (int i=1; i<=t; i++)    {        scanf("%d",&n);        LL ans=(f[n]*2LL-p[n]+M)%M;        printf("Case #%d: %lld\n",i,ans);    }    return 0;}