51nod 1238 杜教筛

传送门：51nod 1238

---

题意

求

$$G(N) = \sum_{i = 1}^N\sum_{j=1}^Nlcm(i, j)$$

---

题解

首先

$$G(N) = \sum_{i = 1}^N\sum_{j=1}^Nlcm(i, j)=2 \sum_{i = 1}^N\sum_{j=1}^ilcm(i, j) - \sum_{i=1}^Nlcm(i, i) \\= 2\sum_{i = 1}^Ni\sum_{d|i}\sum_{u=1}^{\frac{i}{d}}u[gcd(u,\frac{i}{d}) = 1] -\frac{N(N+1)}{2} \\=2\sum_{i = 1}^Ni\sum_{d|i}\frac{\frac{i}{d}\phi(\frac{i}{d}) + [\frac{i}{d} = 1]}{2} - \frac{N(N+1)}{2}\\=\sum_{i = 1}^Ni\sum_{d|i}\frac{i}{d}\phi(\frac{i}{d}) + [\frac{i}{d} = 1] -\frac{N(N+1)}{2}\\=\sum_{i = 1}^Ni(1+\sum_{d|i}\frac{i}{d}\phi(\frac{i}{d})) - \frac{N(N+1)}{2}\\=\sum_{i = 1}^Ni\sum_{d|i}\frac{i}{d}\phi(\frac{i}{d})\\=\sum_{i = 1}^Ni\sum_{d|i}d\phi(d)$$

---

其次枚举因子对i的贡献:

$$G(N) = \sum_{i = 1}^Ni\sum_{d|i}d\phi(d)\\= \sum_{i=1}^N\sum_{d=1}^{\lfloor\frac{N}{i}\rfloor}d\phi(d)$$

---

可用分块:
令

$$F(x)=\sum_{i=1}^xi\phi(i)$$
预处理 $N^{\frac{2}{3}}$以内的 $F(i)$
对于 $x > N^{\frac{2}{3}}$
杜教筛求解:

$$F(x)=\sum_{i=1}^xi\phi(i)=\sum_{i=1}^xi(i-\sum_{d|i,d<i}\phi(d))\\= \sum_{i=1}^xi^2-\sum_{i=1}^xi\sum_{d|i,d<i}\phi(d))\\= \frac{x(2x+1)(x+1)}{6}-\sum_{i=2}^xi\sum_{d=1}^{\lfloor\frac{x}{i}\rfloor}d\phi(d)\\=\frac{x(2x+1)(x+1)}{6}-\sum_{i=2}^xiF(\lfloor\frac{x}{i}\rfloor)$$

分块求解：

---

code:

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;

const ll mod = 1e9 + 7;
const ll inv = (mod + 1) / 2;
const ll inv2 = inv * (mod + 1) / 3 % mod;
const int mo = 2333333;
const int N = 1000001;

bool isPrime[N];
int prime[N];
ll phi[N];
ll sum[N];

void init()
{
    int cnt = 0;
    memset(isPrime, true, sizeof isPrime);
    isPrime[1] = false;
    phi[1] = 1;
    sum[1] = 1;
    for(int i = 2; i < N; ++i)
    {
        if(isPrime[i])
        {
            prime[++cnt] = i;
            phi[i] = i - 1;
        }
        for(int j = 1; j <= cnt && i * prime[j] < N; ++j)
        {
            isPrime[i * prime[j]] = false;
            if(i % prime[j] == 0)
            {
                phi[i * prime[j]] = phi[i] * prime[j];
                break;
            }
            phi[i * prime[j]] = phi[i] * (prime[j] - 1);
        }
    }
    for(ll i = 2; i < N; ++i)
    {
        sum[i] = (sum[i - 1] + i * i * phi[i] % mod) % mod;
    }
}


ll Sum(ll x)
{
    return x % mod * ((x * 2 + 1) % mod) % mod * ((x + 1) % mod) % mod * inv2 % mod;
}
ll getSum(ll R, ll L)
{
    return ((Sum(R) - Sum(L)) % mod + mod) % mod;
}
ll n;

int next[mo], last[mo];
ll t[mo], v[mo];
int l;
void add(int x, ll y, ll z)
{
    t[++l] = y;
    next[l] = last[x];
    last[x] = l;
    v[l] = z;
}
ll cal(ll x)
{
    if(x < N) return sum[x];
    int k = x % mo;
    for(int i = last[k]; i; i = next[i])
        if(t[i] == x) return v[i];
    ll res = x % mod * ((x + 1) % mod) % mod * inv % mod;
    res = res * res % mod;
   // if(x == 10000000000) cout << res << endl;
    ll r;
    for(ll i = 2; i <= x; i = r + 1)
    {
        ll tmp = x / i;
        r = x / tmp;
        ll mul = getSum(r, i - 1);
        //if(x == 10000000000)cout << mul << " ";
        ll ans = cal(tmp) * mul % mod;
       // if(x == 10000000000)cout << ans << " ";
        res = ((res - ans) % mod + mod) % mod;
       // if(x == 10000000000) cout << res << endl;
    }
    add(k, x, res);
    return res;
}


ll Calc(ll x)
{
    ll res = 0;
    ll r;
    for(ll i = 1; i <= x; i = r + 1)
    {
        ll tmp = x / i;
        r = x / tmp;
        ll mul = ((r - i + 1) % mod) % mod * ((r + i) % mod) % mod * inv % mod;
       // cout << mul << " ";
        res = (res + mul * cal(tmp)  % mod) % mod;
        //cout << res << endl;
    }
    r = 0;
  /**  r = x * inv % mod;
    ll mul = 1;
    for(ll i = 1; i * i <= x; ++i)
    {
        if(x % i == 0)
        {
            mul = (mul + i * getPhi(i)) % mod;
            if(i * i != x)  mul = (mul + x / i * getPhi(x / i)) % mod;
        }
    }
    r = r * mul % mod;*/
    return ((res - r) % mod + mod) % mod;
}
int main()
{

    init();
    cin >> n;
 // cout << inv << " " <<inv2 << endl;
    cout << Calc(n) << endl;

    return 0;
}