51nod 1584 加权约数和

原题链接.

题解：

必知知识：

σ(n)是n的约数和。
若$n=\prod p_i^{q_i}$
则$σ(n)=\prod \sum_{j=0}^{q_i} p_i^j$
$=\prod { p_i^{q_i + 1} - 1\over p_i - 1}$
知道了这个就可以线性筛法了，如有不懂得见code。

---

必要结论：

$σ(i*j)=\sum_{p|i} \sum_{q|j} p *{ j \over q} (gcd(p, q) =1)$

证明：

$\sum_{p|i} \sum_{q|j} p *{ j \over q} (gcd(p, q) =1)$
$=\sum_{p|i} \sum_{q|j} p *q (gcd(p, { j \over q}) =1)$

对每一个质因子分类讨论：

设$i=\sum p_i^{u_i}$,$j=p_i^{v_i}$
$p=\sum p_i^{l_i}$,$q=\sum p_i^{r_i}$

1.若$u_i>=l_i>0$,则$r_i = v_i$,此时可以表示出$p_i^{(v_i+1)->(u_i+v_i)}$。
2.$l_i=0$，则$v_i>=r_i>=0$,此时可以表示出$p_i^{0->v_i}$.

于是可以得出$p*q$一定可以表示出每个质因子$p_i^{0->(u_i+v_i)}$，合起来就包括了所有i*j的所有约数。

∴$σ(i*j)=\sum_{p|i} \sum_{q|j} p *{ j \over q} (gcd(p, q) =1)$

---

有了gcd就可以反演了。

$Ans = \sum_{i=1}^n \sum_{j=1}^n σ(i*j)*max(i, j)$
$=\sum_{i=1}^n \sum_{j=1}^n\sum_{p|i} \sum_{q|j} p *{ j \over q}*max(i, j) (gcd(p, q) =1)$
$=\sum_{p =1}^n \sum_{q =1}^n \sum_{i = 1}^{{\lfloor {n \over p} \rfloor}} \sum_{j = 1}^{{\lfloor {n \over q} \rfloor}} p*j*max(pi,qj)(gcd(p, q)=1)$
$=\sum_{p =1}^n \sum_{q =1}^n \sum_{i = 1}^{{\lfloor {n \over p} \rfloor}} \sum_{j = 1}^{{\lfloor {n \over q} \rfloor}} p*j*max(pi,qj)\sum_{d | gcd(i, j)}μ(d)$
$=\sum_{d=1}^n μ(d) * d \sum_{p =1}^{{\lfloor {n \over d} \rfloor}} \sum_{q =1}^{{\lfloor {n \over d} \rfloor}} \sum_{i = 1}^{{\lfloor {n \over d*p} \rfloor}} \sum_{j = 1}^{{\lfloor {n \over d*q} \rfloor}} p*j*max(d*pi,d*qj)$
$=\sum_{d=1}^n μ(d) * d^2*\sum_{p =1}^{{\lfloor {n \over d} \rfloor}} \sum_{q =1}^{{\lfloor {n \over d} \rfloor}} \sum_{i = 1}^{{\lfloor {n \over d*p} \rfloor}} \sum_{j = 1}^{{\lfloor {n \over d*q} \rfloor}} p*j*max(pi,qj)$

注意到d后面是和n/d有关的一个函数。
设$g(n/d) =\sum_{p =1}^{{\lfloor {n \over d} \rfloor}} \sum_{q =1}^{{\lfloor {n \over d} \rfloor}} \sum_{i = 1}^{{\lfloor {n \over d*p} \rfloor}} \sum_{j = 1}^{{\lfloor {n \over d*q} \rfloor}} p*j*max(pi,qj)$

$g(n) =\sum_{p =1}^{n} \sum_{q =1}^{n} \sum_{i = 1}^{{\lfloor {n \over p} \rfloor}} \sum_{j = 1}^{{\lfloor {n \over q} \rfloor}} p*j*max(pi,qj)$
仔细观察你发现这就是一个σ.
$g(n) =\sum_{i=1}^n\sum_{j=1}^n σ(i)*σ(j)*max(i,j)$
$=(\sum_{i=1}^n σ(i)*i*\sum_{j=1}^i *σ(j))*2 - (\sum_{i=1}^n i *σ(i)^2)$
用线性筛法预处理1-n的σ，就可以预处理g(1-n)了。

转回去，现在给出n，求$\sum_{d=1}^n g(\lfloor {n\over d}\rfloor)*μ(d) * d^2$
设$M(d) = μ(d) * d^2$
则求$\sum_{d=1}^n g(\lfloor {n\over d}\rfloor)*M(d)$

1≤数据组数≤50000，1≤N≤1000000。

不要以为你的暴力分块能过。

用差分表O(n log n)预处理即可O(1)查询了。

Code：

#include<cstdio>
#include<cstring>
#define ll long long
#define fo(i, x, y) for(int i = x; i <= y; i ++)
#define min(a, b) ((a) < (b) ? (a) : (b))
#define max(a, b) ((a) > (b) ? (a) : (b))
using namespace std;

const ll N = 1000000;
const ll mo = 1e9 + 7, ni_2 = 5e8 + 4, ni_6 = 166666668;

bool bz[N + 5];
ll p[N];
ll u[N + 5], v[N + 5], s[N + 5];
ll mu[N + 5], s2[N + 5], s3[N + 5], s4[N + 5];

int T, n; ll ans;

int main() {
    mu[1] = u[1] = v[1] = s[1] = 1;
    fo(i, 2, N) {
        if(!bz[i]) p[++ p[0]] = i, mu[i] = -1, u[i] = v[i] = i, s[i] = i + 1;
        fo(j, 1, p[0]) {
            ll k = i * p[j];
            if(k > N) break;
            bz[k] = 1;
            if(i % p[j] == 0) {
                mu[k] = 0;
                u[k] = u[i]; v[k] = v[i] * u[k];
                s[k] = s[k / v[k]] * ((v[k] * u[k] - 1) / (u[k] - 1));
                break;
            }
            u[k] = v[k] = p[j];
            s[k] = s[i] * s[p[j]];
            mu[k] = -mu[i];
        }
    }
    fo(i, 1, N) mu[i] = mu[i] * i % mo * i % mo;
    fo(i, 1, N) s2[i] = (s2[i - 1] + s[i]) % mo;
    fo(i, 1, N) s3[i] = (s3[i - 1] + s[i] * i % mo * s2[i] % mo) % mo;
    fo(i, 1, N) s4[i] = (s4[i - 1] + s[i] * s[i] % mo * i % mo) % mo;
    fo(i, 1, N) s[i] = (s3[i] * 2 - s4[i] + mo) % mo;
    memset(s2, 0, sizeof(s2));
    fo(i, 1, N) fo(j, 1, N / i) {
        ll x = i * j, y = min(N, i * (j + 1) - 1);
        s2[x] += mu[i] * s[j] % mo; s2[y + 1] -= mu[i] * s[j] % mo;
    }
    fo(i, 1, N) s3[i] = ((s3[i - 1] + s2[i]) % mo + mo) % mo;
    scanf("%lld", &T);
    fo(ii, 1, T) {
        ans = 0;
        scanf("%d", &n);
        printf("Case #%d: %d\n", ii, s3[n]);
    }
}