2019年ICPC南京网络赛E-K Sum 莫比乌斯反演+杜教筛+欧拉降幂

题意

$定义f_n(k)=\sum_{l_1=1}^n\sum_{l_2=1}^n…\sum_{l_k=1}^n[gcd(l_1,l_2,…,l_k)]^2$

$$求解\sum_{i=2}^kf_{n(i)}$$

思路

$$f_n(k)=\sum_{l_1=1}^n\sum_{l_2=1}^n…\sum_{l_k=1}^n[gcd(l_1,l_2,…,l_k)]^2$$

$$=\sum_{t=1}^n\sum_{l_1=1}^n\sum_{l_2=1}^n…\sum_{l_k=1}^n[gcd(l_1,l_2,…,l_k)=t]*t^2$$

$$=\sum_{t=1}^{n}t^2\sum_{l_1=1}^{\left \lfloor \frac{n}{t} \right \rfloor }\sum_{l_2=1}^{\left \lfloor \frac{n}{t} \right \rfloor }…\sum_{l_k=1}^{\left \lfloor \frac{n}{t} \right \rfloor }[gcd(l_1,l_2,…,l_k)=1]$$

$$=\sum_{t=1}^{n}t^2\sum_{l_1=1}^{\left \lfloor \frac{n}{t} \right \rfloor }\sum_{l_2=1}^{\left \lfloor \frac{n}{t} \right \rfloor }…\sum_{l_k=1}^{\left \lfloor \frac{n}{t} \right \rfloor }\sum_{dl_1\;…dl_k}\mu(d)$$

$$=\sum_{t=1}^nt^2\sum_{d=1}^{\left \lfloor \frac{n}{t} \right \rfloor}\mu(d)*\left \lfloor \frac{n}{dt} \right \rfloor^k$$

$令T=dt，则：$

$$f_n(k)=\sum_{T=1}^n\sum_{dT}\mu(d)\left \lfloor \frac{T}{d} \right \rfloor^2\left \lfloor \frac{n}{T} \right \rfloor^k$$

$$\sum_{i=2}^kf_n(i)=\sum_{i=2}^k\sum_{T=1}^n\sum_{dT}\mu(d)\left \lfloor \frac{T}{d} \right \rfloor^2\left \lfloor \frac{n}{T} \right \rfloor^i$$

$$\sum_{T=1}^n\sum_{dT}\mu(d)\left \lfloor \frac{T}{d} \right \rfloor^2\sum_{i=2}^k\left \lfloor \frac{n}{T} \right \rfloor^i$$

$令h(x)=\sum_{dx}\mu(d)\left \lfloor \frac{x}{d} \right \rfloor^2，设S(n)=\sum_{i=1}^nh(i)，用杜教筛求出h的前缀和。$

$根据杜教筛求前缀和公式：$

$$S(n)=\frac{\sum_{i=1}^n(f*g)i-\sum_{i=2}^ng[i]S(\left \lfloor \frac{n}{i} \right \rfloor)}{g[1]}$$

$为了简便公式，设g=I,则(h*g)i=(I*\mu *id^2)i=i^2$

$则$ $$S(n)=\frac{n(n+1)(2n+1)}{6}-\sum_{i=2}^nS(\left \lfloor \frac{n}{i} \right \rfloor)$$

$对于\sum_{i=2}^k\left \lfloor \frac{n}{T} \right \rfloor^i，利用等比公式，p=\left \lfloor \frac{n}{T} \right \rfloor.$

$$\sum_{i=2}^k\left \lfloor \frac{n}{T} \right \rfloor^i=\frac{(\frac{n}{T})^2*[1-(\frac{n}{T})^{k-1}]}{1-\frac{n}{T}}$$

$注意\frac{n}{T}=1的情况，答案为k-1。$

$因为k很大，所以可以用欧拉降幂预先处理再作等比，再配合杜教筛求前缀和。$

$$\sum_{i=2}^nf_n(i)=\sum_{T=1}^nh(i)*\frac{(\frac{n}{T})^2*[1-(\frac{n}{T})^{k-1}]}{1-\frac{n}{T}}$$

$总时间复杂度为O(T(n^{\frac{2}{3}}+\sqrt n\log n))$

Code(MS)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
##include "bits/stdc++.h"
using namespace std;

typedef long long ll;

##define endl "\n"
const ll mod = 1e9 + 7;

const int N = 1e5 + 10;
bool is_prime[N];
int prime[N], cnt;
ll h[N], s[N];
ll inv6, k1, k2;

void sieve() {
    h[1] = 1;
    is_prime[0] = is_prime[1] = true;
    for(int i = 2;i < N; i++) {
        if(!is_prime[i]) {
            prime[++cnt] = i;
            h[i] = 1ll * i * i % mod - 1;
        }
        for(int j = 1;j <= cnt && i * prime[j] < N; j++) {
            is_prime[i * prime[j]] = true;
            if(i % prime[j] == 0) {
                h[i * prime[j]] = h[i] * prime[j] % mod * prime[j] % mod;
                break;
            }
            else {
                h[i * prime[j]] = h[i] * h[prime[j]] % mod;
            }
        }
    }
    for(int i = 1;i < N; i++) s[i] = (s[i - 1] + h[i]) % 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;
}

map<ll, ll> mp;
ll S(ll x) {
    if(x < N) return s[x];
    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 / l), x);
        ans = (ans - 1ll * (r - l + 1) * S(x / l) % mod) % mod;
    }
    return mp[x] = (ans % mod + mod) % mod;
}

ll Sum(ll p) {
    if(p == 1)
        return k2 - 1;
    else
        return p * p % mod * (quick_pow(p, k1 - 1) - 1) % mod * quick_pow(p - 1, mod - 2) % mod;
}

void init() {
    sieve();
    inv6 = quick_pow(6, mod - 2);
}

void solve() {
    init();
    int _; cin >> _;
    while(_--) {
        k1 = k2 = 0;
        int n; cin >> n;
        string K; cin >> K;
        for(int i = 0;i < K.length(); i++) {
            k1 = (k1 * 10 + K[i] - '0') % (mod - 1);
            k2 = (k2 * 10 + K[i] - '0') % mod;
        }
        k1 += mod - 1;
        ll ans = 0;
        ll pre = 0, now = 0;
        for(int l = 1, r;l <= n; l = r + 1) {
            r = min(n / (n / l), n);
            now = S(r);
            ans = (ans + (now - pre + mod) % mod * Sum(n / l) % mod) % mod;
            pre = now;
        }
        cout << (ans % mod + mod) % mod << endl;
    }
}

signed main() {
    solve();
}

本文作者：jujimeizuo
本文地址： https://blog.jujimeizuo.cn/2021/03/04/2019-icpc-online-nanjing-e/
本博客所有文章除特别声明外，均采用 CC BY-SA 3.0 协议。转载请注明出处！