莫比乌斯反演例题集---（自用）

参考：大佬1 大佬2 大佬3 大佬4

P3455 [POI2007]ZAP-Queries

$$求解\sum _{i=1}^n\sum _{j=1}^m[gcd(i,j)=k]$$
$反演过程：$
$$\sum _{i=1}^n\sum _{j=1}^m[gcd(i,j)=k]$$

$$\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{k}\rfloor }\epsilon [gcd(i,j)=1]$$

$枚举d：$

$$\sum _{d=1}^{\lfloor min(\lfloor \frac{n}{k}\rfloor ,\lfloor \frac{m}{k}\rfloor )\rfloor }\mu (d)\lfloor \frac{n}{kd}\rfloor \lfloor \frac{m}{kd}\rfloor$$

$先预处理\mu ，然后整数分块做，复杂度为O(n+T\sqrt{n})$

Code

#include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

#define INF 0x7f7f7f
#define mem(a, b) memset(a , b , sizeof(a))
#define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const int N = 5e5 + 10;

int mu[N]; // 莫比乌斯函数
bool is_prime[N];
int prime[N];
int cnt;

ll sum[N]; // 记录莫比乌斯函数前缀和


void Mobi() // 莫比乌斯函数初始化
{
    mu[1] = 1;
    is_prime[0] = is_prime[1] = true;
    for(int i = 2;i < N; i++) {
        if (!is_prime[i]) {
            mu[i] = -1;
            prime[++cnt] = i;
        }
        for (int j = 1; j <= cnt && i * prime[j] < N; j++) {
            is_prime[i * prime[j]] = true;
            if (i % prime[j] == 0) {
                mu[i * prime[j]] = 0;
                break;
            }
            mu[i * prime[j]] = -mu[i];
        }
    }

    for(int i = 1;i < N; i++) {
        sum[i] = sum[i - 1] + mu[i];
    }
}

ll k;

ll Ans(ll n, ll m) {
    ll ans = 0;
    n /= k, m /= k;
    int mx = min(n, m);
    for(int l = 1, r;l <= mx; l = r + 1) {
        r = min(n / (n / l), m / (m / l));
        ans += (n / l) * (m / l) * (sum[r] - sum[l - 1]);
    }
    return ans;
}

void solve() {
    Mobi();
    int T;
    cin >> T;
    while(T--) {
        ll n, m;
        cin >> n >> m >> k;
        cout << Ans(n, m) << endl;
    }
}


signed main() {
    ios_base::sync_with_stdio(false);
    //cin.tie(nullptr);
    //cout.tie(nullptr);
#ifdef FZT_ACM_LOCAL
    freopen("in.txt", "r", stdin);
    freopen("out.txt", "w", stdout);
    signed test_index_for_debug = 1;
    char acm_local_for_debug = 0;
    do {
        if (acm_local_for_debug == '$') exit(0);
        if (test_index_for_debug > 20)
            throw runtime_error("Check the stdin!!!");
        auto start_clock_for_debug = clock();
        solve();
        auto end_clock_for_debug = clock();
        cout << "Test " << test_index_for_debug << " successful" << endl;
        cerr << "Test " << test_index_for_debug++ << " Run Time: "
             << double(end_clock_for_debug - start_clock_for_debug) / CLOCKS_PER_SEC << "s" << endl;
        cout << "--------------------------------------------------" << endl;
    } while (cin >> acm_local_for_debug && cin.putback(acm_local_for_debug));
#else
    solve();
#endif
    return 0;
}

P2257 YY的GCD

$$求解\sum _{p\in prime}\sum _{i=1}^n\sum _{j=1}^m[gcd(i,j)=p]$$
$反演过程和上面几乎一样，就多了一维p为质数的情况，直接化简为：$

$$\sum _{p\in prime}\sum _{d=1}^{\lfloor min(\lfloor \frac{n}{p}\rfloor ,\lfloor \frac{m}{p}\rfloor )\rfloor }\mu (d)\lfloor \frac{n}{pd}\rfloor \lfloor \frac{m}{pd}\rfloor$$

$这里的复杂度为O(n+T\sqrt{n}\sqrt{n})，但明显不够，所以需要继续化简，这里的优化可以用T替换pd，则d=\frac{T}{p}，替换后的式子为：$

$$\sum _{T=1}^{min(n,m)})\lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor \sum _{p|Tp\in prime}\mu (\frac{T}{p})$$

$$令F(x)=\sum _{p|xp\in prime}\mu (\frac{x}{p})$$

$得：$

$$\sum _{p\in prime}\sum _{i=1}^n\sum _{j=1}^m[gcd(i,j)=p]=\sum _{T=1}^{min(n,m)}\lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor F(T)$$

$O(n)处理\mu ，O(nlogn)处理F，O(\sqrt{n})询问分块，即复杂度为O(nlogn+T\sqrt{n})$

Code

#include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

#define INF 0x7f7f7f
#define mem(a, b) memset(a , b , sizeof(a))
#define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const int N = 1e7 + 10;

int mu[N]; // 莫比乌斯函数
bool is_prime[N];
int prime[N];
int cnt;

ll sum[N]; // 记录莫比乌斯函数前缀和
ll F[N];


void Mobi() // 莫比乌斯函数初始化
{
    mu[1] = 1;
    is_prime[0] = is_prime[1] = true;
    for(int i = 2;i < N; i++) {
        if (!is_prime[i]) {
            mu[i] = -1;
            prime[++cnt] = i;
        }
        for (int j = 1; j <= cnt && i * prime[j] < N; j++) {
            is_prime[i * prime[j]] = true;
            if (i % prime[j] == 0) {
                mu[i * prime[j]] = 0;
                break;
            }
            mu[i * prime[j]] = -mu[i];
        }
    }

    for(int i = 1;i <= cnt; i++) {
        for(int j = 1;j * prime[i] < N; j++) {
            F[j * prime[i]] += mu[j];
        }
    }

    for(int i = 1;i < N; i++) {
        F[i] += F[i - 1];
    }
}

ll k;

ll Ans(ll N, ll M) {
    ll ans = 0;
    ll n = N, m = M;
    int mx = min(n, m);
    for (int l = 1, r; l <= mx; l = r + 1) {
        r = min(n / (n / l), m / (m / l));
        ans += (n / l) * (m / l) * (F[r] - F[l - 1]);
    }
    return ans;
}

void solve() {
    Mobi();
    int T;
    cin >> T;
    while(T--) {
        ll n, m;
        cin >> n >> m;
        cout << Ans(n, m) << endl;
    }
}


signed main() {
    ios_base::sync_with_stdio(false);
    //cin.tie(nullptr);
    //cout.tie(nullptr);
#ifdef FZT_ACM_LOCAL
    freopen("in.txt", "r", stdin);
    freopen("out.txt", "w", stdout);
    signed test_index_for_debug = 1;
    char acm_local_for_debug = 0;
    do {
        if (acm_local_for_debug == '$') exit(0);
        if (test_index_for_debug > 20)
            throw runtime_error("Check the stdin!!!");
        auto start_clock_for_debug = clock();
        solve();
        auto end_clock_for_debug = clock();
        cout << "Test " << test_index_for_debug << " successful" << endl;
        cerr << "Test " << test_index_for_debug++ << " Run Time: "
             << double(end_clock_for_debug - start_clock_for_debug) / CLOCKS_PER_SEC << "s" << endl;
        cout << "--------------------------------------------------" << endl;
    } while (cin >> acm_local_for_debug && cin.putback(acm_local_for_debug));
#else
    solve();
#endif
    return 0;
}

P1390 公约数的和

$$求解\sum _{i=1}^n\sum _{j=i+1}^{n}gcd(i,j)$$

$先来看看基本形式：$
$$\sum _{i=1}^n\sum _{j=1}^{n}gcd(i,j)$$
$则ans=\frac{Calc(n,n)-\frac{n(n+1)}{2}}{2}，首先减去gcd(i,i)的个数\frac{n(n+1)}{2}，然后减去gcd(i,j)=gcd(j,i)的个数，直接除2即可，则：$
$$\sum _{i=1}^n\sum _{j=1}^n\sum _{k=1}^{n}k[gcd(i,j)=k]$$

$枚举k：$

$$\sum _{k=1}^{n}k\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }\sum _{j=\lfloor \frac{i+1}{k}\rfloor }^{\lfloor \frac{n}{k}\rfloor }[gcd(i,j)=1]$$

$$\sum _{k=1}^{n}k\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{k}\rfloor }\sum _{d|id|j}\mu (d)$$

$$\sum _{k=1}^{n}k\sum _{d=1}^{\lfloor \frac{n}{k}\rfloor }\mu (d)\sum _{i=1}^{\lfloor \frac{n}{kd}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{kd}\rfloor }1$$

$$\sum _{k=1}^{n}k\sum _{d=1}^{\lfloor \frac{n}{k}\rfloor }\mu (d)(\lfloor \frac{n}{kd}\rfloor {)}^2$$

$令T=kd并枚举T，得：$

$$\sum _{T=1}^n\sum _{d|T}\mu (d)\frac{T}{d}(\lfloor \frac{n}{T}\rfloor {)}^2$$

$中间部分有{\sum }_{d|T}\mu (d)\frac{T}{d}=\phi (T)，因为\phi =\mu \ast id，即得：$

$$\sum _{T=1}^n\phi (T)(\lfloor \frac{n}{T}\rfloor {)}^2$$

$则最终答案为：$

$$ans=\frac{{\sum }_{T=1}^n\phi (T)(\lfloor \frac{n}{T}\rfloor {)}^2-\frac{n(n+1)}{2}}{2}$$

$O(n)预处理\phi ，\sqrt{n}询问分块，总时间复杂度为O(n+\sqrt{n})$

Code

#include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

#define INF 0x7f7f7f
#define mem(a, b) memset(a , b , sizeof(a))
#define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const ll P = 1e9 + 7;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const int N = 2e6 + 10;

int phi[N];
bool is_prime[N];
int prime[N];
int tot;

ll f[N];

void Euler()
{
    phi[1] = 1;
    is_prime[1] = true;
    for(int i = 2;i < N; i++){
        if(!is_prime[i]) {
            phi[i] = i - 1;
            prime[++tot] = i;
        }
        for(int j = 1;j <= tot && i * prime[j] < N; j++){
            is_prime[i * prime[j]] = true;
            if(i % prime[j]) {
                phi[i * prime[j]] = phi[i] * (prime[j] - 1);
            }
            else{
                phi[i * prime[j]] = phi[i] * prime[j];
                break;
            }
        }
    }

    for(int i = 1;i < N; i++) {
        f[i] = f[i - 1] + phi[i];
    }
}

ll Calc(ll n) {
    ll ans = 0;

    for(ll l = 1, r;l <= n;l = r + 1) {
        r = min(n, n / (n / l));
        ans += (f[r] - f[l - 1]) * (n / l) * (n / l);
    }
    return ans;
}

void solve() {
    Euler();
    ll n;
    cin >> n;
    cout << (Calc(n) - n * (n + 1) / 2) / 2 << endl;
}

signed main() {
    ios_base::sync_with_stdio(false);
    //cin.tie(nullptr);
    //cout.tie(nullptr);
#ifdef FZT_ACM_LOCAL
    freopen("in.txt", "r", stdin);
    freopen("out.txt", "w", stdout);
    signed test_index_for_debug = 1;
    char acm_local_for_debug = 0;
    do {
        if (acm_local_for_debug == '$') exit(0);
        if (test_index_for_debug > 20)
            throw runtime_error("Check the stdin!!!");
        auto start_clock_for_debug = clock();
        solve();
        auto end_clock_for_debug = clock();
        cout << "Test " << test_index_for_debug << " successful" << endl;
        cerr << "Test " << test_index_for_debug++ << " Run Time: "
             << double(end_clock_for_debug - start_clock_for_debug) / CLOCKS_PER_SEC << "s" << endl;
        cout << "--------------------------------------------------" << endl;
    } while (cin >> acm_local_for_debug && cin.putback(acm_local_for_debug));
#else
    solve();
#endif
    return 0;
}

P1829 [国家集训队]Crash的数字表格 / JZPTAB

$$求解\sum _{i=1}^n\sum _{j=1}^{m}lcm(i,j)$$

$由gcd和lcm之间的关系得：$

$$\sum _{i=1}^n\sum _{j=1}^m\frac{ij}{gcd(i,j)}$$

$$\sum _{k=1}^{min(n,m)}\sum _{i=1}^n\sum _{j=1}^m\frac{ij}{k}[gcd(i,j)=k]$$

$$\sum _{k=1}^{min(n,m)}\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{k}\rfloor }\frac{ij{k}^2}{k}[gcd(i,j)=1]$$

$$\sum _{k=1}^{min(n,m)}k\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{k}\rfloor }ij\sum _{d|id|j}\mu (d)$$

$枚举d得：$

$$\sum _{k=1}^{min(n,m)}k\sum _{d=1}^{min(\lfloor \frac{n}{k}\rfloor ,\lfloor \frac{m}{k}\rfloor )}\mu (d)\sum _{i=1}^{\lfloor \frac{n}{kd}\rfloor }id\sum _{j=1}^{\lfloor \frac{m}{kd}\rfloor }jd$$

$$\sum _{k=1}^{min(n,m)}k\sum _{d=1}^{min(\lfloor \frac{n}{k}\rfloor ,\lfloor \frac{m}{k}\rfloor )}\mu (d)d^2\sum _{i=1}^{\lfloor \frac{n}{kd}\rfloor }i\sum _{j=1}^{\lfloor \frac{m}{kd}\rfloor }j$$

$令f(x)={\sum }_{i=1}^{x}i，F(n,m)={\sum }_{d=1}^{min(n,m)}\mu (d)d^{2}f(\frac{n}{d})f(\frac{m}{d})，得：$

$$\sum _{k=1}^{min(n,m)}kF(\lfloor \frac{n}{k}\rfloor ,\lfloor \frac{m}{k}\rfloor )$$

$O(n)处理\mu ,O(1)处理f（等差数列，在计算F的过程中直接直接计算），O(\sqrt{\sqrt{n}})分块询问F，答案分块询问O(\sqrt{n})，总复杂度为O(n+n^{\frac{3}{4}})$

Code

#include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

#define INF 0x7f7f7f
#define mem(a, b) memset(a , b , sizeof(a))
#define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const int N = 1e7 + 10;

const ll mod = 20101009;

int mu[N]; // 莫比乌斯函数
bool is_prime[N];
int prime[N];
int cnt;

ll sum[N];

void Mobi() // 莫比乌斯函数初始化
{
    mu[1] = 1;
    is_prime[0] = is_prime[1] = true;
    for(int i = 2;i < N; i++) {
        if (!is_prime[i]) {
            mu[i] = -1;
            prime[++cnt] = i;
        }
        for (int j = 1; j <= cnt && i * prime[j] < N; j++) {
            is_prime[i * prime[j]] = true;
            if (i % prime[j] == 0) {
                mu[i * prime[j]] = 0;
                break;
            }
            mu[i * prime[j]] = -mu[i];
        }
    }

    for(int i = 1;i < N; i++) {
        sum[i] = (sum[i - 1] + mu[i] * i % mod * i % mod) % mod;
    }
}

ll Ans(ll n, ll m)
{
    ll inv2 = (mod + 1) / 2;
    ll ans = 0;
    ll p = min(n, m);
    for(int k = 1;k <= p; k++) {
        int x = n / k, y = m / k, P = min(x, y);
        ll Sum = 0;
        for(int l = 1, r;l <= P; l = r + 1) {
            r = min(x / (x / l), (y / (y / l)));
            Sum = (Sum + (sum[r] - sum[l - 1] + mod) % mod * (1 + x / l) % mod * (x / l) % mod * inv2 % mod * (1 + y / l) % mod * (y / l) % mod * inv2 % mod) % mod;
        }
        ans = (ans + Sum * k % mod) % mod;
    }
    return ans;
}

void solve() {
    Mobi();
    ll n, m;
    cin >> n >> m;
    cout << Ans(n, m) << endl;

}

signed main() {
    ios_base::sync_with_stdio(false);
    //cin.tie(nullptr);
    //cout.tie(nullptr);
#ifdef FZT_ACM_LOCAL
    freopen("in.txt", "r", stdin);
    freopen("out.txt", "w", stdout);
    signed test_index_for_debug = 1;
    char acm_local_for_debug = 0;
    do {
        if (acm_local_for_debug == '$') exit(0);
        if (test_index_for_debug > 20)
            throw runtime_error("Check the stdin!!!");
        auto start_clock_for_debug = clock();
        solve();
        auto end_clock_for_debug = clock();
        cout << "Test " << test_index_for_debug << " successful" << endl;
        cerr << "Test " << test_index_for_debug++ << " Run Time: "
             << double(end_clock_for_debug - start_clock_for_debug) / CLOCKS_PER_SEC << "s" << endl;
        cout << "--------------------------------------------------" << endl;
    } while (cin >> acm_local_for_debug && cin.putback(acm_local_for_debug));
#else
    solve();
#endif
    return 0;
}

P3768 简单的数学题

$$求解\sum _{i=1}^n\sum _{j=1}^{n}ij(gcd(i,j))modp$$
$$\sum _{k=1}^n\sum _{i=1}^n\sum _{j=1}^{n}ijk[gcd(i,j)=k]$$

$$\sum _{k=1}^{n}k\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{k}\rfloor }ij{k}^2[gcd(i,j)=1]$$

$$\sum _{k=1}^n{k}^3\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{k}\rfloor }ij\sum _{d|id|j}\mu (d)$$

$$\sum _{k=1}^n{k}^3\sum _{d=1}^{\lfloor \frac{n}{k}\rfloor }\mu (d)\sum _{i=1}^{\lfloor \frac{n}{kd}\rfloor }id\sum _{j=1}^{\lfloor \frac{n}{kd}\rfloor }jd$$

$$\sum _{k=1}^n{k}^3\sum _{d=1}^{\lfloor \frac{n}{k}\rfloor }\mu (d)d^2\sum _{i=1}^{\lfloor \frac{n}{kd}\rfloor }i\sum _{j=1}^{\lfloor \frac{n}{kd}\rfloor }j$$

$令f(x)={\sum }_{i=1}^{x}i，则式子变为：$

$$\sum _{k=1}^n{k}^3\sum _{d=1}^{\lfloor \frac{n}{k}\rfloor }\mu (d)d^2{f}^2(\frac{n}{kd})$$

$到这里，时间复杂度为O(n+n\sqrt{n})，我没尝试过，但是还可以继续化简，化简方法有替换，用T替换kd，则：$

$$\sum _{T=1}^n{f}^2(\frac{n}{T})T^2\sum _{k|T}k\ast \mu (\frac{T}{k})$$

$$\sum _{T=1}^n{f}^2(\frac{n}{T})T^2(id\ast \mu )T$$

$$\sum _{T=1}^n{f}^2(\frac{n}{T})T^2\phi (T)$$

$可以杜教筛筛T^2\phi (T)的前缀和，因为n最高有10^{10}，明显数组存不下，所以对于5e6可以用数组直接过度求前缀和，对于后面就可以用杜教筛求前缀和，然后用map存下（STL的好处），过程如下：$

$参考：$杜教筛

$先看杜教筛的基本形式(与上面无关)：设S(n)={\sum }_{i=1}^{n}f(i)，定义一个新函数为g，则有$

$$\sum _{i=1}^n(f\ast g)i$$

$$\sum _{i=1}^n\sum _{xy=i}f(x)g(y)$$

$$\sum _{y=1}^{n}g(y)\sum _{xy\le n}f(x)$$

$$\sum _{y=1}^{n}g(y)\sum _{x=1}^{\lfloor \frac{n}{y}\rfloor }f(x)$$

$$\sum _{y=1}^{n}g(y)S(\lfloor \frac{n}{y}\rfloor )$$

$则有：$

$$\sum _{i=1}^n(f\ast g)i=g(1)S(n)+\sum _{y=2}^{n}g(y)S(\lfloor \frac{n}{y}\rfloor )$$

$$S(n)=\frac{{\sum }_{i=1}^n(f\ast g)i-{\sum }_{y=2}^{n}g(y)S(\lfloor \frac{n}{y}\rfloor )}{g(1)}$$

$对于本题，要找到合适的g函数，因为我们求的前缀和为S(n)={\sum }_{i=1}^n{i}^2\phi (i),f(i)=i^2\phi (i)，所以设g(n)=n^2。$

${\sum }_{i=1}^{n}g(i)=\frac{x(x+1)(2x+1)}{6}$

$(f\ast g)x={\sum }_{d|x}f(d)g(\frac{x}{d})={\sum }_{d|x}{d}^2\phi (d)\frac{x^2}{d^2}=x^2{\sum }_{d|x}\phi (d)=x^{2}id(x)=x^3$

${\sum }_{i=1}^n(f\ast g)i={\sum }_{i=1}^n{i}^3=\frac{x^2(x+1{)}^2}{4}$

$$S(n)=\frac{n^2(n+1{)}^2}{4}-\sum _{y=2}^{n}g(y)S(\lfloor \frac{n}{y}\rfloor )$$

$O(n^{\frac{2}{3}})处理T^2\phi (T)的前缀和，f可以O(1)算出，O(\sqrt{n})询问分块，总时间复杂度为O(n^{\frac{2}{3}})$

Code

#include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

#define INF 0x7f7f7f
#define mem(a, b) memset(a , b , sizeof(a))
#define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const int N = 5e6 + 10;

ll n, mod, inv4, inv6;
ll phi[N];
bool is_prime[N];
int prime[N];
int tot;

ll sum[N];

void Get_phi() {
    phi[1] = 1;
    is_prime[1] = true;
    for(int i = 2;i <= N; i++){
        if(!is_prime[i]) {
            phi[i] = i - 1;
            prime[++tot] = i;
        }
        for(int j = 1;j <= tot && i * prime[j] < N; j++){
            is_prime[i * prime[j]] = true;
            if(i % prime[j]) {
                phi[i * prime[j]] = phi[i] * (prime[j] - 1);
            }
            else{
                phi[i * prime[j]] = phi[i] * prime[j];
                break;
            }
        }
    }
    for(int i = 1;i < N; i++) {
        sum[i] = (sum[i - 1] + (ll)phi[i] * i % mod * i % mod) % mod;
    }
}

ll quick_pow(ll a, ll b) {
    ll ans = 1;
    while(b){
        if(b & 1) ans = ans * a % mod;
        a = a * a % mod;
        b >>= 1;
    }
    return ans;
}

ll sum1(ll x) {x %= mod; return (1 + x) % mod * x % mod * (1 + x) % mod * x % mod * inv4 % mod;}
ll sum2(ll x) {x %= mod; return x % mod * (1 + x) % mod * (2 * x + 1) % mod * inv6 % mod;}

map<ll, ll> mp;

ll Calc(ll x) {
    if(x < N) return sum[x];
    if(mp[x]) return mp[x];
    ll num = sum1(x);
    for(ll l = 2, r;l <= x;l = r + 1) {
        r = x / (x / l);
        num = (num - Calc(x / l) % mod * (sum2(r) - sum2(l - 1) + mod) % mod + mod) % mod;
    }
    return mp[x] = num;
}

ll Ans() {
    ll ans = 0;
    for(ll l = 1, r;l <= n; l = r + 1) {
        r = n / (n / l);
        ans = (ans + sum1(n / l) % mod * (Calc(r) - Calc(l - 1) % mod + mod) % mod) % mod;
    }
    return ans;
}

void solve() {
    cin >> mod >> n;
    Get_phi();
    inv4 = quick_pow(4, mod - 2);
    inv6 = quick_pow(6, mod - 2);
    cout << Ans() << endl;
}

signed main() {
    ios_base::sync_with_stdio(false);
    //cin.tie(nullptr);
    //cout.tie(nullptr);
#ifdef FZT_ACM_LOCAL
    freopen("in.txt", "r", stdin);
    freopen("out.txt", "w", stdout);
    signed test_index_for_debug = 1;
    char acm_local_for_debug = 0;
    do {
        if (acm_local_for_debug == '$') exit(0);
        if (test_index_for_debug > 20)
            throw runtime_error("Check the stdin!!!");
        auto start_clock_for_debug = clock();
        solve();
        auto end_clock_for_debug = clock();
        cout << "Test " << test_index_for_debug << " successful" << endl;
        cerr << "Test " << test_index_for_debug++ << " Run Time: "
             << double(end_clock_for_debug - start_clock_for_debug) / CLOCKS_PER_SEC << "s" << endl;
        cout << "--------------------------------------------------" << endl;
    } while (cin >> acm_local_for_debug && cin.putback(acm_local_for_debug));
#else
    solve();
#endif
    return 0;
}

P3327 [SDOI2015]约数个数和

$$求解\sum _{i=1}^n\sum _{j=1}^{m}d(ij)$$

$$约数的重要性质:(ij)=\sum _{x|i}\sum _{y|j}[gcd(x,y)=1]$$

https://www.cnblogs.com/sun123zxy/p/12295533.html

$得：$
$$\sum _{i=1}^n\sum _{j=1}^m\sum _{x|i}\sum _{y|j}[gcd(x,y)=1]$$

$改变枚举顺序，枚举x和y:$

$$\sum _{x=1}^n\sum _{y=1}^m\lfloor \frac{n}{x}\rfloor \lfloor \frac{n}{y}\rfloor [gcd(x,y)=1]$$

$$\sum _{i=1}^n\sum _{j=1}^m\lfloor \frac{n}{i}\rfloor \lfloor \frac{n}{j}\rfloor [gcd(i,j)=1]$$

$后面套用前面得：$

$$\sum _{i=1}^n\sum _{j=1}^m\lfloor \frac{n}{i}\rfloor \lfloor \frac{n}{j}\rfloor \sum _{d|id|j}\mu (d)$$

$枚举d得：$

$$\sum _{d=1}^{min(n,m)}\mu (d)\sum _{i=1}^n\sum _{j=1}^m\lfloor \frac{n}{i}\rfloor \lfloor \frac{n}{j}\rfloor [d|gcd(i,j)]$$

$由于[d|gcd(i,j)]成立的条件为d是gcd(i,j)的约数，所以可以把枚举i,j变换为枚举di,dj，从而[1|gcd(i,j)]可以省略。得：$

$$\sum _{d=1}^{min(n,m)}\mu (d)\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }\lfloor \frac{n}{id}\rfloor \lfloor \frac{n}{jd}\rfloor$$

$$\sum _{d=1}^{min(n,m)}\mu (d)\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\lfloor \frac{n}{id}\rfloor \sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }\lfloor \frac{n}{jd}\rfloor$$

$令f(x)={\sum }_{i=1}^x\lfloor \frac{x}{i}\rfloor ，则：$

$$\sum _{d=1}^{min(n,m)}\mu (d)f(\frac{n}{d})f(\frac{m}{d})$$

$O(n)处理\mu 和其前缀和，O(n\sqrt{n})处理f，O(\sqrt{n})询问分块，总复杂度为O(n\sqrt{n}+T\sqrt{n})$

Code

#include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

#define INF 0x7f7f7f
#define mem(a, b) memset(a , b , sizeof(a))
#define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const int N = 5e4 + 10;

ll mu[N]; // 莫比乌斯函数
bool is_prime[N];
int prime[N];
int cnt;

ll f[N];

void Mobi() // 莫比乌斯函数初始化
{
    mu[1] = 1;
    is_prime[0] = is_prime[1] = true;
    for(int i = 2;i < N; i++) {
        if (!is_prime[i]) {
            mu[i] = -1;
            prime[++cnt] = i;
        }
        for (int j = 1; j <= cnt && i * prime[j] < N; j++) {
            is_prime[i * prime[j]] = true;
            if (i % prime[j] == 0) {
                mu[i * prime[j]] = 0;
                break;
            }
            mu[i * prime[j]] = -mu[i];
        }
    }

    for(int i = 1;i < N; i++) {
        mu[i] += mu[i - 1];
    }

    for(int i = 1;i < N; i++) {
        ll ans = 0;
        for(int l = 1, r;l <= i; l = r + 1) {
            r = i / (i / l);
            ans += (r - l + 1) * (i / l);
        }
        f[i] = ans;
    }
}

ll Ans(ll n, ll m)
{
    ll ans = 0;
    ll k = min(n, m);

    for(ll l = 1, r;l <= k; l = r + 1) {
        r = min(n / (n / l), m / (m / l));
        ans += (mu[r] - mu[l - 1]) * f[m / l] * f[n / l];
    }
    return ans;
}


void solve() {
    int T;
    Mobi();
    cin >> T;
    while(T--) {
        ll n, m;
        cin >> n >> m;
        cout << Ans(n, m) << endl;
    }
}

signed main() {
    ios_base::sync_with_stdio(false);
    //cin.tie(nullptr);
    //cout.tie(nullptr);
#ifdef FZT_ACM_LOCAL
    freopen("in.txt", "r", stdin);
    freopen("out.txt", "w", stdout);
    signed test_index_for_debug = 1;
    char acm_local_for_debug = 0;
    do {
        if (acm_local_for_debug == '$') exit(0);
        if (test_index_for_debug > 20)
            throw runtime_error("Check the stdin!!!");
        auto start_clock_for_debug = clock();
        solve();
        auto end_clock_for_debug = clock();
        cout << "Test " << test_index_for_debug << " successful" << endl;
        cerr << "Test " << test_index_for_debug++ << " Run Time: "
             << double(end_clock_for_debug - start_clock_for_debug) / CLOCKS_PER_SEC << "s" << endl;
        cout << "--------------------------------------------------" << endl;
    } while (cin >> acm_local_for_debug && cin.putback(acm_local_for_debug));
#else
    solve();
#endif
    return 0;
}

P3312 [SDOI2014]数表

$$求解\sum _{i=1}^n\sum _{j=1}^m\sum _{d|id|j}d[\sum _{d|id|j}d\le a]$$

$令F(x)={\sum }_{i|x}i，则得：$

$$\sum _{i=1}^n\sum _{j=1}^{m}F(gcd(i,j))[F(gcd(i,j))\le a]$$

$$\sum _{d=1}^{min(n,m)}\sum _{i=1}^n\sum _{j=1}^{m}F(d)[gcd(i,j)=d][F(d)\le a]$$

$$\sum _{d=1}^{min(n,m)}F(d)[F(d)\le a]\sum _{i=1}^n\sum _{j=1}^m[gcd(i,j)=d]$$

$后半部分就套用前面得：$

$$\sum _{d=1}^{min(n,m)}F(d)[F(d)\le a]\sum _{d^{{}^{\prime }}=1}^{min(\lfloor \frac{n}{d}\rfloor ,\lfloor \frac{m}{d}\rfloor )}\mu (d^{{}^{\prime }})\lfloor \frac{n}{d{d}^{{}^{\prime }}}\rfloor \lfloor \frac{m}{d{d}^{{}^{\prime }}}\rfloor$$

$令T替换d{d}^{{}^{\prime }}得：$

$$\sum _{T=1}^{min(n,m)}\sum _{d|T}F(d)\mu (\frac{T}{d})\lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor [F(d)\le a]$$

$$\sum _{T=1}^{min(n,m)}\lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor \sum _{d|T}F(d)\mu (\frac{T}{d})[F(d)\le a]$$

$令f(x)={\sum }_{d|x}F(d)\mu (\frac{x}{d})[F(d)\le a]得：$

$$\sum _{T=1}^{min(n,m)}\lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor f(T)$$

$O(n)处理\mu ，O(nlogn)处理F，对于每次询问的a值可以离线排序，树状数组维护f，处理时间为O(nlo{g}^{2}n)，O(\sqrt{n}logn)询问分块，总时间复杂度为O(nlo{g}^{2}n)+T\sqrt{n}logn$

Code

#include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

#define INF 0x3f3f3f3f
#define lc u << 1
#define rc u << 1 | 1
#define m (l + r) / 2
#define mid (t[u].l + t[u].r) / 2
#define lowbit(x) x & (-x)
#define mem(a, b) memset(a , b , sizeof(a))
#define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const ll mod = 1e9 + 7;
// const double eps = 1e-6;
// const double PI = acos(-1);
// const double R = 0.57721566490153286060651209;

const int N = 1e5 + 1;
ll mod = (1ll * 1) << 31;

struct node {
    int id;
    ll data;
}F[N];
struct Node {
    int nn, mm, id;
    ll a;
}q[N];

int cnt, prime[N];
bool is_prime[N];
int mu[N];

bool cmp_F(node a, node b) {
    return a.data < b.data;
}

bool cmp_Q(Node a, Node b) {
    return a.a < b.a;
}

inline void init()
{
    mu[1] = 1;
    is_prime[0] = is_prime[1] = true;
    for(int i = 2;i < N; i++) {
        if (!is_prime[i]) {
            mu[i] = -1;
            prime[++cnt] = i;
        }
        for (int j = 1; j <= cnt && i * prime[j] < N; j++) {
            is_prime[i * prime[j]] = true;
            if (i % prime[j] == 0) {
                mu[i * prime[j]] = 0;
                break;
            }
            mu[i * prime[j]] = -mu[i];
        }
    }

    for(int i = 1;i < N; i++) {
        F[i].id = i;
        for(int j = 1;j * i < N; j++) {
            F[i * j].data += i;
        }
    }
    sort(F + 1, F + N, cmp_F);

}

ll t[N];
ll ans[N];

inline void add(int x, int v) {
    while(x < N) {
        t[x] = (t[x] + v) % mod;
        if(t[x] > mod)
            t[x] %= mod;
        x += lowbit(x);
    }
}

inline ll query(int x) {
    ll ans = 0;
    while(x) {
        ans = (ans + t[x]) % mod;
        if(ans > mod)
            ans %= mod;
        x -= lowbit(x);
    }
    return ans;
}

inline ll Calc(int nn, int mm) {
    if(nn > mm)
        swap(nn, mm);

    ll ans = 0;
    for(int l = 1, r;l <= nn; l = r + 1) {
        r = min(nn / (nn / l), mm / (mm / l));
        ans = (ans + 1ll * (nn / l) * (mm / l) % mod * (query(r) - query(l - 1)) % mod) % mod;
    }
    return ans;
}

//struct SegMent {
//    int l, r;
//    int sum;
//}t[N << 2];
//
//void push_up(int u) {
//    t[u].sum = t[lc].sum + t[rc].sum;
//}
//
//void build(int u, int l, int r) {
//    t[u].l = l;
//    t[u].r = r;
//    if(l == r) {
//        t[u].sum = a[l];
//        return ;
//    }
//    build(lc, l, m);
//    build(rc, m + 1, r);
//    push_up(u);
//}
//
//void update(int u, int p, int v) {
//    if(t[u].l == t[u].r) {
//        t[u].sum += v;
//        return ;
//    }
//    if(p <= mid) update(lc, p, v);
//    else update(rc, p, v);
//    push_up(u);
//}
//
//int query(int u, int ql, int qr) {
//    if(ql <= t[u].l && t[u].r <= qr) {
//        return t[u].sum;
//    }
//    int ans = 0;
//    if(ql <= mid)
//        ans += query(lc, ql, qr);
//    if(qr > mid)
//        ans += query(rc, ql, qr);
//    return ans;
//}

void solve() {
    init();
    int T;
    scanf("%d",&T);
    for(int i = 1;i <= T; i++) {
        scanf("%d%d%lld",&q[i].nn,&q[i].mm, &q[i].a);
        q[i].id = i;
    }
    sort(q + 1, q + T + 1, cmp_Q);

    int now = 1;
    for(int i = 1;i <= T; i++) {
        while(now < N && F[now].data <= q[i].a) {
            for(int j = 1;j * F[now].id < N; j++) {
                add(j * F[now].id, mu[j] * F[now].data);
            }
            now++;
        }
        ans[q[i].id] = (Calc(q[i].nn, q[i].mm) % mod + mod) % mod;
    }
    for(int i = 1;i <= T; i++) {
        printf("%lld\n",ans[i]);
    }
}

signed main() {
    ios_base::sync_with_stdio(false);
    //cin.tie(nullptr);
    //cout.tie(nullptr);
#ifdef FZT_ACM_LOCAL
    freopen("in.txt", "r", stdin);
    freopen("out.txt", "w", stdout);
    signed test_index_for_debug = 1;
    char acm_local_for_debug = 0;
    do {
        if (acm_local_for_debug == '$') exit(0);
        if (test_index_for_debug > 20)
            throw runtime_error("Check the stdin!!!");
        auto start_clock_for_debug = clock();
        solve();
        auto end_clock_for_debug = clock();
        cout << "Test " << test_index_for_debug << " successful" << endl;
        cerr << "Test " << test_index_for_debug++ << " Run Time: "
             << double(end_clock_for_debug - start_clock_for_debug) / CLOCKS_PER_SEC << "s" << endl;
        cout << "--------------------------------------------------" << endl;
    } while (cin >> acm_local_for_debug && cin.putback(acm_local_for_debug));
#else
    solve();
#endif
    return 0;
}

P3172 [CQOI2015]选数

$$求解\sum _{i=L}^H\sum _{j=L}^H...\sum _{x=L}^H[gcd(i,j...x)=k]$$

$本题求n个数的gcd为k的个数，根据上面套路，将H/k,L/k，这样就只用求n个数的gcd为1的个数了，注意：如果L整除k，则L=L/k，否则L=L/k+1，就是向上取整，那么有一种办法非常好，就让L/k+(\frac{k-1}{k})，即：$

$$\sum _{i=\lfloor \frac{L-1}{k}\rfloor +1}^{\lfloor \frac{H}{k}\rfloor }\sum _{j=\lfloor \frac{L-1}{k}\rfloor +1}^{\lfloor \frac{H}{k}\rfloor }...\sum _{x=\lfloor \frac{L-1}{k}\rfloor +1}^{\lfloor \frac{H}{k}\rfloor }[gcd(i,j...x)=1]$$

$$\sum _{i=\lfloor \frac{L-1}{k}\rfloor +1}^{\lfloor \frac{H}{k}\rfloor }\sum _{j=\lfloor \frac{L-1}{k}\rfloor +1}^{\lfloor \frac{H}{k}\rfloor }...\sum _{x=\lfloor \frac{L-1}{k}\rfloor +1}^{\lfloor \frac{H}{k}\rfloor }\sum _{d|id|j..d|x}\mu (d)$$

$枚举d，令l=\lfloor \frac{L-1}{k}\rfloor +1,r=\lfloor \frac{H}{k}\rfloor ：$

$$\sum _{d=1}^r\mu (d)\sum _{i=l}^r[d|i]\sum _{j=l}^r[d|j]...\sum _{x=l}^r[d|x]$$

$只有当i，j...x都是d的倍数的时候，最终答案+1，根据容斥原理，那么在[l,r]中d倍数有\frac{r}{d}-\frac{l-1}{d}，对于所有的i，j...x，他们当中所有[l,r]的d的倍数的个数有(\frac{r}{d}-\frac{l-1}{d}{)}^n，即为：$

$$\sum _{d=1}^r\mu (d)\ast (\frac{r}{d}-\frac{l-1}{d}{)}^n$$

$这里的H非常大，如果用线性筛去处理\mu ，会直接Wa掉，因为没法处理到10^9次方那么大的前缀和，数组也开不了那么大。后果：$

$这里的AC可能就是一小部分的运气罢了。如果用线性筛，RE的RE，WA的WA，所以推荐用杜教筛取筛。$

$设f(d)=\mu (d),S(r)={\sum }_{d=1}^{n}f(d)，这里的S就是要求的前缀和。$
$我们选择g=I,则{\sum }_{i=1}^r(f\ast g)i={\sum }_{i=1}^r(\mu \ast I)i={\sum }_{i=1}^r\epsilon (i)=1。$

$枚举y：$

$$\sum _{y=1}^{r}g(y)\sum _{x=1}^{\lfloor \frac{r}{y}\rfloor )}f(x)$$

$$\sum _{y=1}^{r}g(y)S(\lfloor \frac{r}{y}\rfloor )$$

$$\sum _{i=1}^r(f\ast g)i=g(1)S(r)+\sum _{y=2}^{r}g(y)S(\lfloor \frac{r}{y}\rfloor )$$

$$S(r)=\frac{{\sum }_{i=1}^r(f\ast g)i-{\sum }_{y=2}^{r}g(y)S(\lfloor \frac{r}{y}\rfloor )}{g(1)}$$

$当g=I时，即：$

$$S(r)=1-\sum _{y=2}^{r}S(\lfloor \frac{r}{y}\rfloor )$$

$杜教筛的复杂度为O(n\frac{2}{3})，比线性筛还要快一点，总时间复杂度为O(n^{\frac{2}{3}})。$

Code

#include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

#define INF 0x7f7f7f
#define mem(a, b) memset(a , b , sizeof(a))
#define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const ll mod = 1e9 + 7;

const int N = 1e6 + 10;

int mu[N]; // 莫比乌斯函数
bool is_prime[N];
int prime[N];
int cnt;
ll sum[N];

ll n, L, H, k;

ll quick_pow(ll a, ll b) {
    ll ans = 1;
    while(b) {
        if(b & 1) ans = ans * a % mod;
        a = a * a % mod;
        b >>= 1;
    }
    return ans;
}

void Mobi() // 莫比乌斯函数初始化
{
    mu[1] = 1;
    is_prime[0] = is_prime[1] = true;
    for(int i = 2;i < N; i++) {
        if (!is_prime[i]) {
            mu[i] = -1;
            prime[++cnt] = i;
        }
        for (int j = 1; j <= cnt && i * prime[j] < N; j++) {
            is_prime[i * prime[j]] = true;
            if (i % prime[j] == 0) {
                mu[i * prime[j]] = 0;
                break;
            }
            mu[i * prime[j]] = -mu[i];
        }
    }
    for(int i = 1;i < N; i++) {
        sum[i] = (sum[i - 1] + mu[i]) % mod;
    }
}

map<ll , ll> mp;

ll Calc(ll x) {
    if(x < N)
        return sum[x];
    if(mp[x])
        return mp[x];
    ll num = 1;
    for(int l = 2, r;l <= x; l = r + 1) {
        r = x / (x / l);
        num = (num - Calc(x / l) * (r - l + 1) % mod + mod) % mod;
    }
    return mp[x] = num;
}

ll Ans() {
    ll ans = 0;
    for(ll l = 1, r;l <= H;l = r + 1) {
        r = min(H / (H / l), (L / l) ? L / (L / l) : H + 2);
        ans = (ans + (Calc(r) - Calc(l - 1) + mod) % mod * quick_pow((H / l) - (L / l), n) % mod + mod) % mod;
    }
    return (ans + mod) % mod;
}

void solve() {
    cin >> n >> k >> L >> H;
    L = (L - 1) / k; // 因为后面需要L-1，所以这里干脆就不+1了，但是含义不同
    H = H / k;
    Mobi();
    cout << Ans() << endl;
}

signed main() {
    ios_base::sync_with_stdio(false);
    //cin.tie(nullptr);
    //cout.tie(nullptr);
#ifdef FZT_ACM_LOCAL
    freopen("in.txt", "r", stdin);
    freopen("out.txt", "w", stdout);
    signed test_index_for_debug = 1;
    char acm_local_for_debug = 0;
    do {
        if (acm_local_for_debug == '$') exit(0);
        if (test_index_for_debug > 20)
            throw runtime_error("Check the stdin!!!");
        auto start_clock_for_debug = clock();
        solve();
        auto end_clock_for_debug = clock();
        cout << "Test " << test_index_for_debug << " successful" << endl;
        cerr << "Test " << test_index_for_debug++ << " Run Time: "
             << double(end_clock_for_debug - start_clock_for_debug) / CLOCKS_PER_SEC << "s" << endl;
        cout << "--------------------------------------------------" << endl;
    } while (cin >> acm_local_for_debug && cin.putback(acm_local_for_debug));
#else
    solve();
#endif
    return 0;
}

AT5200 [AGC038C] LCMs

$$求解\sum _{i=1}^N\sum _{j=1}^{N}lcm(A_i,A_j)$$

$莫比乌斯反演处理的一般都是枚举变量之间的关系，所以我们要让A_i和A_j变为i和j的形式。$

$想要用i和j表示A_i和A_j并且不遗漏任何的A，所以让A等于某个i，就是让cn{t}_i={\sum }_{d=1}^N[A_d=i]，统计所有A等于某个i的个数，即cn{t}_i，j也同样如此，如下所示：$

$$\sum _{i=1}^N\sum _{j=1}^{N}cn{t}_{i}cn{t}_{j}lcm(i,j)$$

$这样还不够，因为i只能是[1,N]的，当Ai中有些值大于N的话，再枚举i从[1,N]就会漏掉大于N的Ai的情况，所以想要统计所有情况，即让N=Ai中最大的数即可，即M={\max }_{i=1}^N{A}_i。$
$所以我们需要反演的式子为：$

$$\sum _{i=1}^M\sum _{j=1}^{M}cn{t}_i\cdot cn{t}_j\cdot lcm(i,j)$$

$$\sum _{i=1}^M\sum _{j=1}^{M}cn{t}_i\cdot cn{t}_j\cdot \frac{ij}{gcd(i,j)}$$

$$\sum _{k=1}^M\sum _{i=1}^M\sum _{j=1}^{M}cn{t}_i\cdot cn{t}_j\cdot \frac{ij}{k[gcd(i,j)=k]}$$

$$\sum _{k=1}^M\sum _{i=1}^{\lfloor \frac{M}{k}\rfloor }\sum _{j=1}^{\lfloor \frac{M}{k}\rfloor }cn{t}_{ik}\cdot cn{t}_{jk}\cdot \frac{ij{k}^2}{k}[gcd(i,j)=1]$$

$$\sum _{k=1}^{M}k\sum _{i=1}^{\lfloor \frac{M}{k}\rfloor }\sum _{j=1}^{\lfloor \frac{M}{k}\rfloor }cn{t}_{ik}\cdot cn{t}_{jk}\cdot i\cdot j\sum _{d|id|j}\mu (d)$$

$$\sum _{k=1}^{M}k\sum _{d=1}^{\lfloor \frac{M}{k}\rfloor }\mu (d)\sum _{i=1}^{\lfloor \frac{M}{kd}\rfloor }cn{t}_{ikd}\cdot id\sum _{j=1}^{\lfloor \frac{M}{kd}\rfloor }cn{t}_{jkd}\cdot jd$$

$$\sum _{k=1}^{M}k\sum _{d=1}^{\lfloor \frac{M}{k}\rfloor }\mu (d)d^2(\sum _{i=1}^{\lfloor \frac{M}{kd}\rfloor }cn{t}_{ikd}\cdot i{)}^2$$

$令T=kd则：$

$$\sum _{T=1}^M(\sum _{i=1}^{\lfloor \frac{M}{T}\rfloor }cn{t}_{iT}\cdot i{)}^2\sum _{k|T}\mu (k)k^2\frac{T}{k}$$

$$\sum _{T=1}^{M}T(\sum _{i=1}^{\lfloor \frac{M}{T}\rfloor }cn{t}_{iT}\cdot i{)}^2\sum _{k|T}\mu (k)k$$

$令f(x)={\sum }_{i|x}\mu (i)i,g(x)={\sum }_{i=1}^{\lfloor \frac{M}{x}\rfloor }cn{t}_{ix}\cdot i=\frac{1}{x}{\sum }_{x|t}tcn{t}_t,则最后的式子为：$

$$\sum _{T=1}^{M}T{g}^2(T)f(T)$$

$O(n)处理A的次数，O(nlogn)处理f和g，总时间复杂度为O(nlogn)$

Code

#include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

#define INF 0x7f7f7f
#define mem(a, b) memset(a , b , sizeof(a))
#define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const int N = 5e5 + 10;

int mu[N]; // 莫比乌斯函数
bool is_prime[N];
int prime[N];
int tot;
ll f[N], g[N];

int M, n;
int cnt[N];

void Mobi() // 莫比乌斯函数初始化
{
    mu[1] = 1;
    is_prime[0] = is_prime[1] = true;
    for(int i = 2;i < N; i++) {
        if (!is_prime[i]) {
            mu[i] = -1;
            prime[++tot] = i;
        }
        for (int j = 1; j <= tot && i * prime[j] < N; j++) {
            is_prime[i * prime[j]] = true;
            if (i % prime[j] == 0) {
                mu[i * prime[j]] = 0;
                break;
            }
            mu[i * prime[j]] = -mu[i];
        }
    }
    for(int i = 1;i <= M; i++) {
        for(int x = i;x <= M; x += i) {
            f[x] += mu[i] * i;
        }
    }

    for(int x = 1;x <= M; x++) {
        for(int t = x;t <= M; t += x) {
            g[x] += t * cnt[t];
        }
        g[x] /= x;
    }
}

void solve() {
    cin >> n;
    M = n;
    for(int i = 1;i <= n; i++) {
        int x;
        cin >> x;
        cnt[x]++;
        M = max(M, x);
    }
    Mobi();

    ll ans = 0;

    for(int i = 1;i <= M; i++) {
        ans += i * g[i] * g[i] * f[i];
    }

    cout << ans << endl;
}

signed main() {
    ios_base::sync_with_stdio(false);
    //cin.tie(nullptr);
    //cout.tie(nullptr);
#ifdef FZT_ACM_LOCAL
    freopen("in.txt", "r", stdin);
    freopen("out.txt", "w", stdout);
    signed test_index_for_debug = 1;
    char acm_local_for_debug = 0;
    do {
        if (acm_local_for_debug == '$') exit(0);
        if (test_index_for_debug > 20)
            throw runtime_error("Check the stdin!!!");
        auto start_clock_for_debug = clock();
        solve();
        auto end_clock_for_debug = clock();
        cout << "Test " << test_index_for_debug << " successful" << endl;
        cerr << "Test " << test_index_for_debug++ << " Run Time: "
             << double(end_clock_for_debug - start_clock_for_debug) / CLOCKS_PER_SEC << "s" << endl;
        cout << "--------------------------------------------------" << endl;
    } while (cin >> acm_local_for_debug && cin.putback(acm_local_for_debug));
#else
    solve();
#endif
    return 0;
}

SP5971 LCMSUM - LCM Sum

$$求解\sum _{i=1}^{n}lcm(i,n)$$

$\lfloor \frac{n}{k}\rfloor$

$$\sum _{i=1}^n\frac{in}{gcd(i,n)}$$

$$\sum _{k=1}^n\sum _{i=1}^n\frac{in}{k[gcd(i,n)=k]}$$

$$\sum _{k=1}^n\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }in[gcd(i,\frac{n}{k})=1]$$

$$n\sum _{k=1}^n\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }i[gcd(i,\frac{n}{k})=1]$$

$$n\sum _{k|n}\sum _{i=1}^{k}i[gcd(i,k)=1]$$

$$n\sum _{k|n}\sum _{i=1}^{k}i\sum _{d|id|k}\mu (d)$$

$枚举d：$

$$n\sum _{k|n}\sum _{d|k}\mu (d)\sum _{i=1}^{\lfloor \frac{k}{d}\rfloor }id$$

$$n\sum _{k|n}\sum _{d|k}\mu (d)d\sum _{i=1}^{\lfloor \frac{k}{d}\rfloor }i$$

$$n\sum _{k|n}\sum _{d|k}\mu (d)d\frac{\lfloor \frac{k}{d}\rfloor \ast (\lfloor \frac{k}{d}\rfloor +1)}{2}$$

n 2 ∑ k ∣ n { ∑ d ∣ k μ ( d ) ∗ k 2 d + ∑ d ∣ k μ ( d ) ∗ k } \frac{n}{2}\sum_{k|n}\{\sum_{d|k}\mu(d)*\frac{k^2}{d}+\sum_{d|k}\mu(d)*k\} 2n​k∣n∑​{d∣k∑​μ(d)∗dk2​+d∣k∑​μ(d)∗k}

n 2 ∑ k ∣ n k { ∑ d ∣ k μ ( d ) ∗ k d + ∑ d ∣ k μ ( d ) } \frac{n}{2}\sum_{k|n}k\{\sum_{d|k}\mu(d)*\frac{k}{d}+\sum_{d|k}\mu(d)\} 2n​k∣n∑​k{d∣k∑​μ(d)∗dk​+d∣k∑​μ(d)}

$这时候发现两个狄利克雷卷积，\phi (n)=\mu \ast id={\sum }_{d|n}\mu (d)\ast \frac{n}{d},\epsilon (n)=\mu \ast I={\sum }_{d|n}\mu (d)。所以我们将其替换得：$

$$\frac{n}{2}\sum _{k|n}k[\phi (k)+\epsilon (k)]$$

$O(n)处理\phi ，\sqrt{n}询问答案，总复杂度为O(n+T\sqrt{n})。$

Code

#include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

#define INF 0x7f7f7f
#define mem(a, b) memset(a , b , sizeof(a))
#define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const int N = 1e6 + 10;

int phi[N];
bool is_prime[N];
int prime[N];
int tot;

ll f[N];

void Euler()
{
    phi[1] = 1;
    is_prime[1] = true;
    for(int i = 2;i < N; i++){
        if(!is_prime[i]) {
            phi[i] = i - 1;
            prime[++tot] = i;
        }
        for(int j = 1;j <= tot && i * prime[j] < N; j++){
            is_prime[i * prime[j]] = true;
            if(i % prime[j]) {
                phi[i * prime[j]] = phi[i] * (prime[j] - 1);
            }
            else{
                phi[i * prime[j]] = phi[i] * prime[j];
                break;
            }
        }
    }
    for(int i = 1;i < N; i++) {
        for(int j = i;j < N; j += i) {
            f[j]  += 1ll * phi[i] * i + (i == 1);
        }
    }
}

void solve() {
    Euler();
    int T;
    cin >> T;
    while(T--) {
        int n;
        cin >> n;
        cout << f[n] * n / 2 << endl;
    }
}

signed main() {
    ios_base::sync_with_stdio(false);
    //cin.tie(nullptr);
    //cout.tie(nullptr);
#ifdef FZT_ACM_LOCAL
    freopen("in.txt", "r", stdin);
    freopen("out.txt", "w", stdout);
    signed test_index_for_debug = 1;
    char acm_local_for_debug = 0;
    do {
        if (acm_local_for_debug == '$') exit(0);
        if (test_index_for_debug > 20)
            throw runtime_error("Check the stdin!!!");
        auto start_clock_for_debug = clock();
        solve();
        auto end_clock_for_debug = clock();
        cout << "Test " << test_index_for_debug << " successful" << endl;
        cerr << "Test " << test_index_for_debug++ << " Run Time: "
             << double(end_clock_for_debug - start_clock_for_debug) / CLOCKS_PER_SEC << "s" << endl;
        cout << "--------------------------------------------------" << endl;
    } while (cin >> acm_local_for_debug && cin.putback(acm_local_for_debug));
#else
    solve();
#endif
    return 0;
}

P5221 Product

$$求解\prod _{i=1}^n\prod _{j=1}^n\frac{lcm(i,j)}{gcd(i,j)}$$

$首先还是把lcm(i,j)化简得：$

$$\prod _{i=1}^n\prod _{j=1}^n\frac{ij}{[gcd(i,j){]}^2}$$

$$\prod _{i=1}^n\prod _{j=1}^{n}ij\ast \prod _{i=1}^n\prod _{j=1}^n\frac{1}{[gcd(i,j){]}^2}$$

$我们很容易推出前半部分的数值为(n!{)}^{2n}，后半部分的值就是{\prod }_{i=1}^n{\prod }_{j=1}^n[gcd(i,j){]}^2的逆元的平方。$

$所以现在求解：$

$$\prod _{i=1}^n\prod _{j=1}^{n}gcd(i,j)$$

$$\prod _{k=1}^n\prod _{i=1}^n\prod _{j=1}^{n}k[gcd(i,j)=k]$$

$$\prod _{k=1}^n{k}^{{\sum }_{i=1}^n{\sum }_{j=1}^n[gcd(i,j)=k]}$$

$这个式子的变换，本来是[1,n]中所有的gcd累乘，然后我们单独把gcd(i,j)=k的拿出来，k会有[1,n]的范围，然后我们再将这些k累乘就会达到之前的效果。$
$不过，[1,n]中的gcd不单单只有一个，可能会有多个相同的gcd，所以变换之后的式子里还需要统计gcd=k出现的次数，那只需要将他们累加即可。$

$所以没毛病。继续：$

$$\prod _{k=1}^n{k}^{{\sum }_{i=1}^{\lfloor \frac{n}{k}\rfloor }{\sum }_{j=1}^{\lfloor \frac{n}{k}\rfloor }[gcd(i,j)=1]}$$

$$\prod _{k=1}^n{k}^{{\sum }_{i=1}^{\lfloor \frac{n}{k}\rfloor }{\sum }_{j=1}^{\lfloor \frac{n}{k}\rfloor }{\sum }_{d|id|j}\mu (d)}$$

$枚举d：$

$$\prod _{k=1}^n{k}^{{\sum }_{d=1}^{\lfloor \frac{n}{k}\rfloor }\mu (d){\sum }_{i=1}^{\lfloor \frac{n}{kd}\rfloor }{\sum }_{j=1}^{\lfloor \frac{n}{kd}\rfloor }1}$$

$$\prod _{k=1}^n{k}^{{\sum }_{d=1}^{\lfloor \frac{n}{k}\rfloor }\mu (d)(\lfloor \frac{n}{kd}\rfloor {)}^2}$$

$最终答案为：$

$$(n!{)}^{2n}\prod _{k=1}^n{k}^{{\sum }_{d=1}^{\lfloor \frac{n}{k}\rfloor }\mu (d)(\lfloor \frac{n}{kd}\rfloor {)}^2}$$

$O(n)处理\mu 的前缀和,\sqrt{n}询问指数的分块，可以用欧拉降幂优化，总时间复杂度为O(n+nlo{g}_{2}n$

Code

#include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

#define INF 0x7f7f7f
#define mem(a, b) memset(a , b , sizeof(a))
#define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const ll P = 1e9 + 7;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const int N = 1000010;
const ll mod = 104857601;

bool is_prime[N];
int prime[N];
int mu[N]; // 莫比乌斯函数
int cnt;

int n;
ll fac = 1;

ll quick_pow(ll a, ll b)
{
    ll ans = 1;
    while(b) {
        if(b & 1) ans = ans * a % mod;
        a = a * a % mod;
        b >>= 1;
    }
    return ans;
}
void Mobi() // 莫比乌斯函数初始化
{
    mu[1] = 1;
    is_prime[0] = is_prime[1] = true;
    for(int i = 2;i <= n; i++) {
        if (!is_prime[i]) {
            mu[i] = -1;
            prime[++cnt] = i;
        }
        for (int j = 1; j <= cnt && i * prime[j] <= n; j++) {
            is_prime[i * prime[j]] = true;
            if (i % prime[j] == 0) {
                mu[i * prime[j]] = 0;
                break;
            }
            mu[i * prime[j]] = -mu[i];
        }
        mu[i] = mu[i - 1] + mu[i];
        fac = fac * i % mod;
    }
}

ll Ans() {
    ll ans = 1;
    for(int k = 1;k <= n; k++) {
        ll p = 0;
        int t = n / k;
        for(int l = 1, r;l <= t; l = r + 1) {
            r = min(t, t / (t / l));
            p = (p + (mu[r] - mu[l - 1] + mod - 1) * (t / l) % (mod - 1) * (t / l) % (mod - 1)) % (mod - 1);
        }
        ans = ans * quick_pow(k, p) % mod;
    }
    ll k = quick_pow(ans, mod - 2);
    fac = quick_pow(fac, n);
    return quick_pow(k, 2) * fac % mod * fac % mod;
}

void solve() {
    cin >> n;
    Mobi();
    cout << Ans() << endl;
}

signed main() {
    ios_base::sync_with_stdio(false);
    //cin.tie(nullptr);
    //cout.tie(nullptr);
#ifdef FZT_ACM_LOCAL
    freopen("in.txt", "r", stdin);
    freopen("out.txt", "w", stdout);
    signed test_index_for_debug = 1;
    char acm_local_for_debug = 0;
    do {
        if (acm_local_for_debug == '$') exit(0);
        if (test_index_for_debug > 20)
            throw runtime_error("Check the stdin!!!");
        auto start_clock_for_debug = clock();
        solve();
        auto end_clock_for_debug = clock();
        cout << "Test " << test_index_for_debug << " successful" << endl;
        cerr << "Test " << test_index_for_debug++ << " Run Time: "
             << double(end_clock_for_debug - start_clock_for_debug) / CLOCKS_PER_SEC << "s" << endl;
        cout << "--------------------------------------------------" << endl;
    } while (cin >> acm_local_for_debug && cin.putback(acm_local_for_debug));
#else
    solve();
#endif
    return 0;
}

HDU 6588 Function

$$求解\sum _{i=1}^{n}gcd(\lfloor \sqrt[3]{i}\rfloor ,i)$$

$我们先将这个一重和式转化为二重和式，得：$

$$\sum _{a=1}^{\lfloor \sqrt[3]{n}\rfloor }\sum _{i=1}^{n}gcd(a,i)[\lfloor \sqrt[3]{i}\rfloor =a]$$

$我们来看这个式子，只有当\lfloor \sqrt[3]{i}\rfloor =a时，才会有对答案的贡献为gcd(a,i)，所以a\le \lfloor \sqrt[3]{i}\rfloor <a+1，a^3\le i\le (a+1{)}^3-1，所以上式转化为：$

$$\sum _{a=1}^{\lfloor \sqrt[3]{n}\rfloor }\sum _{i=a^3}^{min(n,(a+1{)}^3-1)}gcd(a,i)$$

$然后我们将这个二重和式拆开得：$

$$\sum _{a=1}^{\lfloor \sqrt[3]{n}\rfloor -1}\sum _{i=a^3}^{min(n,(a+1{)}^3-1)}gcd(a,i)+\sum _{i=(\lfloor \sqrt[3]{n}\rfloor {)}^3}^{n}gcd(\lfloor \sqrt[3]{n}\rfloor ,i)$$

$令\sqrt[3]{n}=N，并使a=i,i=j（好看一点）得：$

$$\sum _{i=1}^{N-1}\sum _{j=i^3}^{(i+1{)}^3-1}gcd(i,j)+\sum _{j=N^3}^{n}gcd(N,j)$$

$先计算前半部分：$

$|--------------------------------|$

HDU 6833 A Very Easy Math Problem

$$求解\sum _{a_1=1}^n\sum _{a_2=1}^n...\sum _{a_x=1}^n\left(\prod _{j=1}^x{a}_j^k\right)f(gcd(a_1,a_2...a_x))gcd(a_1,a_2...a_x)$$
$令d=gcd(a_1,a_2...a_x)得：$

$$\sum _{a_1=1}^n\sum _{a_2=1}^n...\sum _{a_x=1}^n\left(\prod _{j=1}^x{a}_j^k\right)f(d)d$$

$枚举d得：$

$$\sum _{d=1}^{n}df(d)\sum _{a_1=1}^n\sum _{a_2=1}^n...\sum _{a_x=1}^n\left(\prod _{j=1}^x{a}_j^k\right)[gcd(a_1,a_2...a_x)=d]$$

$$\sum _{d=1}^n{d}^{kx+1}f(d)\sum _{a_1=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{a_2=1}^{\lfloor \frac{n}{d}\rfloor }...\sum _{a_x=1}^{\lfloor \frac{n}{d}\rfloor }\left(\prod _{j=1}^x{a}_j^k\right)[gcd(a_1,a_2...a_x)=1]$$

$$\sum _{d=1}^n{d}^{kx+1}f(d)\sum _{a_1=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{a_2=1}^{\lfloor \frac{n}{d}\rfloor }...\sum _{a_x=1}^{\lfloor \frac{n}{d}\rfloor }\left(\prod _{j=1}^x{a}_j^k\right)\sum _{t|a_{1}t|a2...t|a_x}\mu (t)$$

$$\sum _{d=1}^n{d}^{kx+1}f(d){\left(\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }{i}^k\right)}^x\sum _{t|a_{1}t|a2...t|a_x}\mu (t)$$

$枚举t得：$

$$\sum _{d=1}^n{d}^{kx+1}f(d)\sum _{t=1}^{\lfloor \frac{n}{d}\rfloor }\mu (t){\left(\sum _{i=1}^{\lfloor \frac{n}{dt}\rfloor }(it{)}^k\right)}^x$$

$$\sum _{d=1}^n{d}^{kx+1}f(d)\sum _{t=1}^{\lfloor \frac{n}{d}\rfloor }\mu (t)t^{kx}{\left(\sum _{i=1}^{\lfloor \frac{n}{dt}\rfloor }{i}^k\right)}^x$$

$令T=dt并枚举T得：$

$$\sum _{T=1}^n{\left(\sum _{i=1}^{\lfloor \frac{n}{dt}\rfloor }{i}^k\right)}^x{T}^{kx}\sum _{d|T}df(d)\mu (\frac{T}{d})$$

$令g(T)={\sum }_{d|T}df(d)\mu (\frac{T}{d})，则$

$$\sum _{T=1}^n{\left(\sum _{i=1}^{\lfloor \frac{n}{T}\rfloor }{i}^k\right)}^x{T}^{kx}g(T)$$

$O(nlogn)筛f和g，并对T^{kx}g(T)和{\sum }_{i=1}^{\lfloor \frac{n}{T}\rfloor }{i}^k做前缀和，最后\sqrt{n}分块计算，复杂度为O(nlogn+T\sqrt{n})$

Code

#include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

#define INF 0x3f3f3f3f
#define lowbit(x) x & (-x)
#define mem(a, b) memset(a , b , sizeof(a))
#define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
const ll mod = 1e9 + 7;
// const double eps = 1e-6;
// const double pi = acos(-1);

const int N = 2e5 + 10;

ll T, k, x;

int mu[N]; // 莫比乌斯函数
bool is_prime[N];
int prime[N];
int cnt;
ll g[N], f[N];
ll sumG[N], sumi[N];

ll quick_pow(ll a, ll b) {
    ll ans = 1;
    while(b) {
        if(b & 1) ans = ans * a % mod;
        a = a * a % mod;
        b >>= 1;
    }
    return ans % mod;
}

void Mobi() // 莫比乌斯函数初始化
{
    f[1] = mu[1] = 1;
    is_prime[0] = is_prime[1] = true;
    for(int i = 2;i < N; i++) {
        f[i] = 1;
        if (!is_prime[i]) {
            mu[i] = -1;
            prime[++cnt] = i;
        }
        for (int j = 1; j <= cnt && i * prime[j] < N; j++) {
            is_prime[i * prime[j]] = true;
            if (i % prime[j] == 0) {
                mu[i * prime[j]] = 0;
                break;
            }
            mu[i * prime[j]] = -mu[i];
        }
    }
    // 筛f函数
    for(int d = 2;d * d < N; d++) {
        for(int i = d * d;i < N; i += d * d) {
            f[i] = 0;
        }
    }
    // 筛g函数
    for(int d = 1;d < N; d++) {
        for(int i = d;i < N; i += d) {
            g[i] = (g[i] + 1ll * d * f[d] % mod * mu[i / d] % mod + mod) % mod;
        }
    }
    // 处理sumG函数和sumi函数前缀和
    for(int i = 1;i < N; i++) {
        ll t = quick_pow(i, k);
        sumi[i] = (sumi[i - 1] + t) % mod;
        sumG[i] = (sumG[i - 1] + quick_pow(t, x) * g[i] % mod) % mod;
    }
}

void solve() {
    cin >> T >> k >> x;
    Mobi();
    while(T--) {
        int n;
        cin >> n;
        ll ans = 0;
        for(int l = 1, r;l <= n; l = r + 1) {
            r = min(n, n / (n / l));
            ans = (ans + (sumG[r] - sumG[l - 1] + mod) % mod * quick_pow(sumi[n / l], x) % mod) % mod;
        }

        cout << ans << endl;
    }
}

signed main() {
    ios_base::sync_with_stdio(false);
    //cin.tie(nullptr);
    //cout.tie(nullptr);
#ifdef FZT_ACM_LOCAL
    freopen("in.txt", "r", stdin);
    freopen("out.txt", "w", stdout);
    signed test_index_for_debug = 1;
    char acm_local_for_debug = 0;
    do {
        if (acm_local_for_debug == '$') exit(0);
        if (test_index_for_debug > 20)
            throw runtime_error("Check the stdin!!!");
        auto start_clock_for_debug = clock();
        solve();
        auto end_clock_for_debug = clock();
        cout << "Test " << test_index_for_debug << " successful" << endl;
        cerr << "Test " << test_index_for_debug++ << " Run Time: "
             << double(end_clock_for_debug - start_clock_for_debug) / CLOCKS_PER_SEC << "s" << endl;
        cout << "--------------------------------------------------" << endl;
    } while (cin >> acm_local_for_debug && cin.putback(acm_local_for_debug));
#else
    solve();
#endif
    return 0;
}

huntian oy

$$求解\sum _{i=1}^n\sum _{j=1}^{i}gcd(i^a-j^a,i^b-j^b)[gcd(i,j)=1]$$

$$gcd(i^a-j^a,i^b-j^b)=i^{gcd(a,b)}-j^{gcd(a,b)}$$

$题目说gcd(a,b)=1，所以那个复杂的式子直接变成i-j.$

$$\sum _{i=1}^n\sum _{j=1}^i(i-j)[gcd(i,j)=1]$$

$$\sum _{i=1}^n\sum _{j=1}^{i}i[gcd(i,j)=1]-\sum _{i=1}^n\sum _{j=1}^{i}j[gcd(i,j)=1]$$

$左边为与i互质的个数，右边为比i小并且互质的数字和，即为：$

$$\sum _{i=1}^{n}i\phi (i)-\frac{{\sum }_{i=1}^{n}i\phi (i)+1}{2}$$

$$\frac{{\sum }_{i=1}^{n}i\phi (i)-1}{2}$$

$因为n有1e9，所以需要杜教筛优化。$

$$f\ast g(n)=\sum _{d|n}f(d)g(\frac{n}{d})=\sum _{d|n}id(d)\ast \phi (d)\ast g(\frac{n}{d})$$

$很显然，设g=id，则：$

$$f\ast g(n)=\sum _{d|n}\phi (d)id(n)=n\sum _{d|n}\phi (d)=id(n)(\phi \ast I)=i{d}^2(n)=n^2$$

$$所以得：S(n)=\sum _{i=1}^n(f\ast g)i-\sum _{i=2}^{n}g(i)S(\frac{n}{i})=\frac{n(n+1)(2n+1)}{6}-\sum _{i=2}^{n}iS(\frac{n}{i})$$

Code

#include "bits/stdc++.h"
using namespace std;

typedef long long ll;
typedef pair<int, int> pii;
typedef pair<double, double> pdd;
typedef pair<ll, ll> pll;

#define endl "\n"
#define eb emplace_back
#define mem(a, b) memset(a , b , sizeof(a))

const ll INF = 0x3f3f3f3f;
// const ll mod = 998244353;
const ll mod = 1e9 + 7;
const double eps = 1e-6;
const double PI = acos(-1);
const double R = 0.57721566490153286060651209;

template<typename T>
inline void read(T &a) {char c = getchar();T x = 0, f = 1;
    while (!isdigit(c)) {if (c == '-')f = -1;c = getchar();}
    while (isdigit(c)) {x = (x << 1) + (x << 3) + c - '0';c = getchar();}a = f * x;
}

const int N = 1e6 + 10;

bool is_prime[N];
int prime[N], cnt;
int phi[N];
ll sum_i_phi[N];

void init() {
    phi[1] = 1;
    for(int i = 2;i < N; i++) {
        if(!is_prime[i]) prime[++cnt] = i, phi[i] = i - 1;
        for(int j = 1;j <= cnt && i * prime[j] < N; j++) {
            is_prime[i * prime[j]] = 1;
            if(i % prime[j] == 0) {
                phi[i * prime[j]] = phi[i] * prime[j];
                break;
            }
            else phi[i * prime[j]] = phi[i] * phi[prime[j]];
        }
    }
    for(int i = 1;i < N; i++) {
        sum_i_phi[i] = (sum_i_phi[i - 1] + 1ll * i * phi[i] % mod) % mod;
    }
}

ll quick_pow(ll a, ll b) {
    ll ans = 1;
    while(b) {
        if(b & 1) ans = ans * a % mod;
         a = a * a % mod;
         b >>= 1;
    }
    return ans % mod;
}

ll inv2;
ll inv6;

map<ll, ll> mp;

ll S(ll x) {
    if(x < N) return sum_i_phi[x];
    else if(mp[x]) return mp[x];
    ll ans = x * (x + 1) % mod * (2 * x + 1) % mod * inv6 % mod;
    for(int l = 2, r;l <= x; l = r + 1) {
        r = min(x, x / (x / l));
        ans = (ans - 1ll * (r + l) % mod * (r - l + 1) % mod * inv2 % mod * S(x / l)) % mod;
    }
    return mp[x] = ans;
}

void solve() {
    init();
    inv2 = quick_pow(2, mod - 2);
    inv6 = quick_pow(6, mod - 2);
    int _; scanf("%d",&_);
    while(_--) {
        ll n, a, b; scanf("%lld%lld%lld",&n,&a,&b);
        ll ans = (S(n) - 1) % mod;
        ans = ans * inv2 % mod;
        printf("%lld\n",(ans % mod + mod) % mod);
    }
}

signed main() {
    ios_base::sync_with_stdio(false);
    // cin.tie(nullptr);
    // cout.tie(nullptr);
#ifdef FZT_ACM_LOCAL
    freopen("in.txt", "r", stdin);
    freopen("out.txt", "w", stdout);
    signed test_index_for_debug = 1;
    char acm_local_for_debug = 0;
    do {
        if (acm_local_for_debug == '$') exit(0);
        if (test_index_for_debug > 20)
            throw runtime_error("Check the stdin!!!");
        auto start_clock_for_debug = clock();
        solve();
        auto end_clock_for_debug = clock();
        cout << "Test " << test_index_for_debug << " successful" << endl;
        cerr << "Test " << test_index_for_debug++ << " Run Time: "
             << double(end_clock_for_debug - start_clock_for_debug) / CLOCKS_PER_SEC << "s" << endl;
        cout << "--------------------------------------------------" << endl;
    } while (cin >> acm_local_for_debug && cin.putback(acm_local_for_debug));
#else
    solve();
#endif
    return 0;
}

HDU 6428 Problem C. Calculate

$$求解\sum _{i=1}^A\sum _{j=1}^B\sum _{k=1}^C\varphi (gcd(i,j^2,k^3))mod{2}^{30}$$

$首先要把\varphi 和gcd分开，即\varphi (n)={\sum }_{d|n}(\varphi \ast \mu )(d)，证明如下：$
$$\varphi (n)=(\mu \ast id)(n)$$

$$\sum _{d|n}\mu (d)\frac{n}{d}$$

$$\sum _{d|n}\mu (d)\sum _{d^{\prime }|\frac{n}{d}}\varphi (d^{\prime })$$

$$\sum _{d{d}^{\prime }|n}\mu (d)\varphi (d^{\prime })$$

$$\sum _{d|n}\sum _{d^{\prime }|d}\mu (d^{\prime })\varphi (\frac{d}{d^{\prime }})$$

$$\sum _{d|n}(\varphi \ast \mu )(d)$$

$所以原式变为：$

$$\sum _{i=1}^A\sum _{j=1}^B\sum _{k=1}^C\sum _{d|id|j^{2}d|k^3}(\varphi \ast \mu )(d)$$

$$\sum _{d=1}^A(\varphi \ast \mu )(d)\sum _{i=1d|i}^A\sum _{j=1d|j^2}^B\sum _{k=1d|k^3}^{C}1$$
\left \lceil \right \rceil
$观察后面式子，对于一个x^k，若d|x^k，先把d分解为\prod p_i^{a_i}.$
$则\prod p_i^{a_i}|x^k，得\prod p_i^{\lceil \frac{a_i}{k}\rceil }|x，所以设$

$$f_k(n)=\prod p_i^{\lceil \frac{a_i}{k}\rceil }$$

$则：$

$$ans=\sum _{d=1}^A(\varphi \ast \mu )(d)\frac{A}{f_1(d)}\frac{B}{f_2(d)}\frac{C}{f_3(d)}$$

$对于f_k(n)，类似分解n的形式。$
$分解质因子的过程中，记录质因子的指数，每次质因子+1时，若\%k=1时，说明该向上取整了，于是f_k(n)\ast =该质因子.$

$由于\varphi 和\mu 都是积性函数，卷积之后还是积性函数，设g(d)=(\varphi \ast \mu )(d)，则$

$$\varphi (n)=\sum _{d|n}g(d)$$

$那么可以得：$

$$\varphi (p^k)=\varphi (p^{k-1})+g(p^k)$$

$$g(p^k)=\varphi (p^k)-\varphi (p^{k-1})$$

$$g(p^k)=(p-1)p^{k-1}-(p-1)p^{k-2}$$

$$g(p^k)=(p-1{)}^2{p}^{k-2}$$

$当我们欧拉筛的过程中：$

$$g(1)=1$$

$$g(p)=p-2$$

$$g(p^k)=(p-1{)}^2{p}^{k-2}$$

$$g(p_1^{k_1}{p}_2^{k_2})=g(p_1^{k_1})g(p_2^{k_2})$$

$$......$$

Code

#include "bits/stdc++.h"
using namespace std;

typedef long long ll;

const ll mod = (1ll << 30);

const int N = 1e7 + 10;
bool is_prime[N];
int prime[N], cnt;
int g[N];
int f1[N], f2[N], f3[N];
int deg[N];

void init() {
    f1[1] = f2[1] = f3[1] = g[1] = 1;
    for(int i = 2;i < N; i++) {
        f1[i] = i;
        if(!is_prime[i]) prime[++cnt] = f2[i] = f3[i] = i, g[i] = i - 2, deg[i] = 1;
        for(int j = 1;j <= cnt && i * prime[j] < N; j++) {
            int now = i * prime[j];
            is_prime[now] = 1;
            if(i % prime[j] == 0) {
                deg[now] = deg[i] + 1;
                int num = 1, tmp = i;
                while(num <= 3 && tmp % prime[j] == 0) num++, tmp /= prime[j];
                if(num == 1) g[now] = g[i] * g[prime[j]];
                else if(num == 2) g[now] = g[i / prime[j]] * (prime[j] - 1) * (prime[j] - 1);
                else g[now] = g[i] * prime[j];
                f2[now] = f2[i] * (deg[now] % 2 == 1 ? prime[j] : 1);
                f3[now] = f3[i] * (deg[now] % 3 == 1 ? prime[j] : 1);
                break;
            }
            else {
                deg[now] = 1;
                f2[now] = f2[i] * f2[prime[j]];
                f3[now] = f3[i] * f3[prime[j]];
                g[now] = g[i] * g[prime[j]];
            }
        }
    }
}

void solve() {
    init();
    int _; scanf("%d",&_);
    while(_--) {
        int A, B, C; scanf("%d%d%d",&A,&B,&C);
        ll ans = 0;
        for(int d = 1;d <= A; d++) {
            ans = (ans + 1ll * g[d] * (A / f1[d]) % mod * (B / f2[d]) % mod * (C / f3[d]) % mod) % mod;
        }
        printf("%lld\n",(ans % mod + mod) % mod);
    }
}

signed main() {
    ios_base::sync_with_stdio(false);
    // cin.tie(nullptr);
    // cout.tie(nullptr);
#ifdef FZT_ACM_LOCAL
    freopen("in.txt", "r", stdin);
    freopen("out.txt", "w", stdout);
    signed test_index_for_debug = 1;
    char acm_local_for_debug = 0;
    do {
        if (acm_local_for_debug == '$') exit(0);
        if (test_index_for_debug > 20)
            throw runtime_error("Check the stdin!!!");
        auto start_clock_for_debug = clock();
        solve();
        auto end_clock_for_debug = clock();
        cout << "Test " << test_index_for_debug << " successful" << endl;
        cerr << "Test " << test_index_for_debug++ << " Run Time: "
             << double(end_clock_for_debug - start_clock_for_debug) / CLOCKS_PER_SEC << "s" << endl;
        cout << "--------------------------------------------------" << endl;
    } while (cin >> acm_local_for_debug && cin.putback(acm_local_for_debug));
#else
    solve();
#endif
    return 0;
}