莫反学习笔记

· · 算法·理论

数竞朋友们一致认为反演是几何的。

下文中的函数默认为积性函数。

本质是 \sum_{d | n} \mu(d) = \varepsilon (n)，即 \text{I} * \mu = \varepsilon。

一般会用到的：若 f(n) = \sum_{d | n} g(d)，则 g(n) = \sum_{d | n} \mu(\frac{n}{d}) f(d)。

P3455 & P2522

后者只是在前者的基础上加入容斥，这里讲前者。

原问即求

\sum_{i=1}^a \sum_{j=1}^b [\gcd(i,j) = d]

首先先把 =d 搞掉。

\sum_{i=1}^{\lfloor\frac{a}{d}\rfloor} \sum_{j=1}^{\lfloor\frac{b}{d}\rfloor} [\gcd(i,j) = 1]

接着把 [] 换成 \varepsilon。

\sum_{i=1}^{\lfloor\frac{a}{d}\rfloor} \sum_{j=1}^{\lfloor\frac{b}{d}\rfloor} \varepsilon(\gcd(i,j))

把 \varepsilon 展开成 \mu。

\sum_{i=1}^{\lfloor\frac{a}{d}\rfloor} \sum_{j=1}^{\lfloor\frac{b}{d}\rfloor} \sum_{k | \gcd(i,j)} \mu(k)

把 \mu 换出来。

\sum_{k=1}^{+\infty} \mu(k) \sum_{i=1}^{\lfloor\frac{a}{d}\rfloor} \sum_{j=1}^{\lfloor\frac{b}{d}\rfloor} [k|i\And k|j]

把后面化简一下。

\sum_{k=1}^{\min(\lfloor\frac{a}{d}\rfloor,\lfloor\frac{b}{d}\rfloor)} \mu(k) \lfloor\frac{a}{kd}\rfloor\lfloor\frac{b}{kd}\rfloor

线性筛后整除分块计算即可。

时间复杂度 O(n + T\sqrt{n})。

const int N = 5e4 + 5;
int cnt, p[N], pri[N], mo[N];
void sieve()
{
    mo[1] = 1;
    for (int i = 2; i < N; i++)
    {   
        if (!p[i])
        {
            p[i] = i;
            pri[++cnt] = i;
            mo[i] = -1ll;
        }
        for (int j = 1; j <= cnt; j++)
        {
            if (i * pri[j] > N) break;
            p[i * pri[j]] = j;
            if (i % pri[j] == 0)
            {
                mo[i * pri[j]] = 0;
                break;
            }
            mo[i * pri[j]] = -mo[i];
        }
    }
    for (int i = 1; i < N; i++) mo[i] = mo[i - 1] + mo[i];
}
int t, a, b, d;
int solve(int n, int m)
{
    int ret = 0;
    for (int i = 1, j; i <= min(n, m); i = j + 1)
    {
        j = min(n / (n / i), m / (m / i));
        ret += (mo[j] - mo[i-1]) * (m / i) * (n / i);
    }
    return ret;
}
int main()
{
    ios::sync_with_stdio(false); cin.tie(0); cout.tie(0);
    sieve();
    cin >> t;
    while (t--)
    {
        cin >> a >> b >> d;
        cout << solve(a / d, b / d) << "\n";
    }
    return 0;
}

SP5971

直接看吧。

\sum_{i=1}^n \text{lcm}(i,n)\\ = \sum_{i=1}^n \frac{ni}{\gcd(n,i)} \\ = \frac{1}{2} (\sum_{i=1}^{n-1} \frac{ni}{\gcd(n,i)} + \sum_{i=n-1}^{1} \frac{ni}{\gcd(n,i)} ) + n \\ = \frac{1}{2} (\sum_{i=1}^{n-1} \frac{ni}{\gcd(n,i)} + \sum_{i=n-1}^{1} \frac{ni}{\gcd(n,n-i)} ) + n \\ = \frac{1}{2} \sum_{i=1}^{n-1} \frac{n^2}{\gcd(n,i)} + n \\ = \frac{1}{2} \sum_{i=1}^{n} \frac{n^2}{\gcd(n,i)} + \frac{n}{2} \\ = \frac{1}{2} n\sum_{i=1}^{n} \frac{n}{\gcd(n,i)} + \frac{n}{2} \\

易知 \gcd(n,i)=d 有 \varphi(\frac{n}{d}) 个，故

\frac{1}{2} n\sum_{i=1}^{n} \frac{n}{\gcd(n,i)} + \frac{n}{2} \\ = \frac{1}{2} n \cdot (\sum_{d|n} \frac{n\varphi(\frac{n}{d})}{d} + 1) \\ = \frac{1}{2} n \cdot (\sum_{d|n} \frac{n}{d}\varphi(\frac{n}{d}) + 1) \\

设 d'=\dfrac{n}{d}，可得到

\frac{1}{2} n \cdot (\sum_{d'|n} d'\varphi(d') + 1)

打包一下，得到 \dfrac{1}{2} n \cdot (g(n) + 1)，考虑如何在线性筛时筛出 g，显然 g 是积性的。

推一下柿子，可知 g(tp) = g(t) + (g(t) - g(\frac{t}{p})) \cdot p^2，这是 p|t 时的柿子。

本题不是莫反，但是可以练习推柿子。

const int N = 1e6 + 5;
int cnt, p[N], pri[N];
ll g[N];
void sieve()
{
    g[1] = 1;
    for (int i = 2; i < N; i++)
    {   
        if (!p[i])
        {
            p[i] = i;
            pri[++cnt] = i;
            g[i] = 1ll * i * (i - 1) + 1;
        }
        for (int j = 1; j <= cnt; j++)
        {
            if (i * pri[j] > N) break;
            p[i * pri[j]] = j;
            if (i % pri[j] == 0)
            {
                g[i * pri[j]] = g[i] + (g[i] - g[i / pri[j]]) * pri[j] * pri[j];
                break;
            }
            g[i * pri[j]] = g[i] * g[pri[j]];
        }
    }
}
int t, n;
int main()
{
    ios::sync_with_stdio(false); cin.tie(0); cout.tie(0);
    sieve();
    cin >> t;
    while (t--)
    {
        cin >> n;
        cout << (g[n] + 1ll) * n / 2 << "\n";
    }
    return 0;
}

P1829

求

\sum_{i=1}^n \sum_{j=1}^m \text{lcm}(i,j)

首先转化，假定 n \le m。

\sum_{i=1}^n \sum_{j=1}^m \frac{i\cdot j}{\gcd(i,j)}

考虑枚举最大公因数 \gcd(i,j) = d。

\sum_{i=1}^n \sum_{j=1}^m \sum_{d | i, d | j, \gcd(\frac{i}{d},\frac{j}{d}) = 1} \frac{i\cdot j}{d}

把 \gcd 拉出来。

\sum_{d=1}^n d \sum_{i=1}^{\lfloor \frac{n}{d} \rfloor} \sum_{j=1}^{\lfloor \frac{m}{d} \rfloor} [\gcd(i,j)=1]i \cdot j

把后半部分单独拉出来，记为 f(n,m)。

f(n,m) = \sum_{i=1}^{n} \sum_{j=1}^{m} [\gcd(i,j)=1]i \cdot j

仿照 P3455 转化。

f(n,m) = \sum_{d=1}^n\sum_{d|i}^n \sum_{d|j}^m \mu(d)\cdot i \cdot j

显然转化成 \gcd(i,j) = 1 的形式，记 i' = \dfrac{i}{d}，j'=\dfrac{j}{d}。

f(n,m) = \sum_{d=1}^n \mu(d)\cdot d^2 \sum_{i=1}^{\lfloor \frac{n}{d} \rfloor} \sum_{j=1}^{\lfloor \frac{m}{d} \rfloor} i'j'

后面的一堆可以记为 g(n,m) = \dfrac{n\cdot(n+1)\cdot m\cdot(m+1)}{4}，故

f(n,m) = \sum_{d=1}^n \mu(d)\cdot d^2 \cdot g(\lfloor \frac{n}{d} \rfloor,\lfloor \frac{m}{d} \rfloor)

这个东西就可以数论分块了。

带回原式即得

\sum_{d=1}^n d \cdot f(\lfloor \frac{n}{d} \rfloor,\lfloor \frac{m}{d} \rfloor)

显然也可以数论分块。

最终复杂度是线性筛的 O(\max\{n,m\})。

注意取模。

const int N = 1e7 + 5, mod = 20101009, inv2 = 10050505, inv4 = 15075757;
int cnt, p[N], pri[N], mo[N], sum[N];
void sieve()
{
    mo[1] = 1;
    for (int i = 2; i < N; i++)
    {   
        if (!p[i])
        {
            p[i] = i;
            pri[++cnt] = i;
            mo[i] = -1ll;
        }
        for (int j = 1; j <= cnt; j++)
        {
            if (i * pri[j] > N) break;
            p[i * pri[j]] = j;
            if (i % pri[j] == 0)
            {
                mo[i * pri[j]] = 0;
                break;
            }
            mo[i * pri[j]] = -mo[i];
        }
    }
    for (int i = 1; i < N; i++) 
        sum[i] = (sum[i - 1] + 1ll * i * i % mod * (mo[i] + mod) % mod) % mod;
}
int n, m;
int s(int l, int r) {return 1ll * (l + r) * (r - l + 1) % mod * inv2 % mod;}
int g(int n, int m) {return 1ll * n * (n + 1) % mod * m % mod * (m + 1) % mod * inv4 % mod;}
int f(int n, int m)
{
    int ret = 0;
    for (int i = 1, j; i <= min(n, m); i = j + 1)
    {
        j = min(n / (n / i), m / (m / i));
        ret = (ret + 1ll * (sum[j] - sum[i-1] + mod) % mod * g(n / i, m / i) % mod) % mod;
    }
    return ret;
}
int solve(int n, int m)
{
    int ret = 0;
    for (int i = 1, j; i <= min(n, m); i = j + 1)
    {
        j = min(n / (n / i), m / (m / i));
        ret = (ret + 1ll * s(i, j) * f(n / i, m / i) % mod) % mod;
    }
    return ret;
}
int main()
{
    ios::sync_with_stdio(false); cin.tie(0); cout.tie(0);
    sieve();
    cin >> n >> m;
    cout << solve(n, m) << "\n";
    return 0;
}

P3312

题目即求

\sum_{i=1}^n \sum_{j=1}^m [\sigma(\gcd(i,j))\le a]\sigma(\gcd(i,j))

首先考虑没有 a 怎么操作，先设 n \le m。

\sum_{i=1}^n \sum_{j=1}^m \sigma(\gcd(i,j)) \\ = \sum_{d=1}^n \sum_{i=1}^n \sum_{j=1}^m \sigma(d)[\gcd(i,j)=d] \\ = \sum_{d=1}^n \sigma(d) \sum_{i=1}^n \sum_{j=1}^m [\gcd(i,j)=d]

后面的部分就比较典了，就是例一的结果。

= \sum_{d=1}^n \sigma(d) \sum_{k=1}^{\lfloor \frac{n}{d}\rfloor} \mu(k) \lfloor \frac{n}{kd}\rfloor\lfloor \frac{m}{kd}\rfloor

考虑令 L = kd。

= \sum_{L=1}^n \lfloor \frac{n}{L}\rfloor\lfloor \frac{m}{L}\rfloor \sum_{d|L} \mu(\frac{L}{d}) \sigma(d)

做到这就可以在没有 a 限制的时候用 O(\sqrt{n}) 解决了。

考虑如何处理 a 的限制。

记 f(L) = \sum_{d|L} \mu(\frac{L}{d}) \sigma(d) [\sigma(d) \le a]，考虑离线，将询问的 a 从小到大排序，将所有 d 也按 \sigma(d) 从小到大排序，每次移动 a 的时候，将对应 \sigma(d) 范围内 d 的倍数 L' 的 f(L') 加上 \mu(\frac{L'}{d}) \sigma(d) 即可。

每次询问用数论分块 + BIT 维护 f(L) 区间和即可。

考虑复杂度。

单次询问是 O(\sqrt{n} \log n) 的，将所有 d 的倍数做单点加是 O(\sum_{i=1}^n \lfloor \frac{n}{i}\rfloor \times \log n) = O(n \log^2 n) 的（调和级数），总复杂度 O(Q\sqrt{n} \log n + n \log^2 n)。

至于怎么筛 \sigma，可以参考 link。

const int N = 1e5 + 5;
int cnt, p[N], pri[N], mo[N], g[N];
pair <int, int> sig[N];
void sieve()
{
    mo[1] = g[1] = 1; sig[1] = {1, 1};
    for (int i = 2; i < N; i++)
    {   
        if (!p[i])
        {
            p[i] = i;
            pri[++cnt] = i;
            mo[i] = -1ll;
            sig[i] = {i + 1, i};
            g[i] = i + 1;
        }
        for (int j = 1; j <= cnt; j++)
        {
            if (i * pri[j] >= N) break;
            p[i * pri[j]] = pri[j];
            if (i % pri[j] == 0)
            {
                mo[i * pri[j]] = 0;
                g[i * pri[j]] = g[i] * pri[j] + 1;
                sig[i * pri[j]] = {sig[i].first / g[i] * g[i * pri[j]], i * pri[j]};
                break;
            }
            mo[i * pri[j]] = -mo[i];
            sig[i * pri[j]] = {sig[i].first * sig[pri[j]].first, i * pri[j]};
            g[i * pri[j]] = pri[j] + 1;
        }
    }
}
struct _ {int n, m, a, id, ans;} s[N];
int q;
int f[N];
void modify(int x, int v) {for (; x < N; x += x & -x) f[x] += v;}
int query(int x) {int ret = 0; for (; x; x -= x & -x) ret += f[x]; return ret;}
int solve(int n, int m)
{
    int ret = 0;
    for (int i = 1, j; i <= min(n, m); i = j + 1)
    {
        j = min(n / (n / i), m / (m / i));
        ret += (query(j) - query(i - 1)) * (n / i) * (m / i);
    }
    return ret;
}
int main()
{
    ios::sync_with_stdio(false); cin.tie(0); cout.tie(0);
    cin >> q; sieve();
    for (int i = 1; i <= q; i++) cin >> s[i].n >> s[i].m >> s[i].a, s[i].id = i;
    sort(s + 1, s + q + 1, [](_ x, _ y){return x.a < y.a;});
    sort(sig + 1, sig + N);
    for (int i = 1, j = 1; i <= q; i++)
    {
        while (j < N && sig[j].first <= s[i].a) 
        {
            int o = sig[j].second;
            for (int p = o; p < N; p += o) modify(p, mo[p/o] * sig[j].first);
            j++;
        }
        s[i].ans = solve(s[i].n, s[i].m);
    }
    sort(s + 1, s + q + 1, [](_ x, _ y){return x.id < y.id;});
    for (int i = 1; i <= q; i++) cout << (s[i].ans & (~(1 << 31))) << "\n";
    return 0;
}

P3704

题目即求

\prod_{i=1}^n \prod_{j=1}^m f_{\gcd(i,j)}

先转换一波，依然设 n \le m。

= \prod_{d=1}^n \prod_{i=1}^n \prod_{j=1}^m [\gcd(i,j) = d]f_d \\ = \prod_{d=1}^n f_d^{\sum_{i=1}^{\lfloor \frac{n}{d}\rfloor} \sum_{j=1}^{\lfloor \frac{m}{d}\rfloor} [\gcd(i,j)=1]} \\ = \prod_{d=1}^n f_d^{\sum_{k=1}^{\lfloor \frac{n}{d}\rfloor} \mu(k) \lfloor \frac{n}{kd}\rfloor\lfloor \frac{m}{kd}\rfloor}

这时我们仿照例四，设 L = kd。

= \prod_{L=1}^n\prod_{d | L} f_d^{\mu(\frac{L}{d}) \lfloor \frac{n}{L}\rfloor\lfloor \frac{m}{L}\rfloor} \\ = \prod_{L=1}^n{\Large (}\prod_{d | L} f_d^{\mu(\frac{L}{d})}{\Large )}^{\lfloor \frac{n}{L}\rfloor\lfloor \frac{m}{L}\rfloor}

内层可以用调和级数复杂度预处理，外层数论分块即可。

时间复杂度 O(n \log n + T \sqrt{n} \log n)，快速幂会多带一只 \log。

const int N = 1e6 + 5, mod = 1e9 + 7;
int cnt, p[N], pri[N], mo[N], f[N], _f[N], F[N];
void add(int &x, int y) {x += y; if (x >= mod) x -= mod;}
void mul(int &x, int y) {x = 1ll * x * y % mod;}
int ksm(int a, int b)
{
    int ret = 1;
    while (b)
    {
        if (b & 1) mul(ret, a);
        mul(a, a);
        b >>= 1;
    }
    return ret;
}
int inv(int x) {return ksm(x, mod - 2);}
void sieve()
{
    mo[1] = 1; F[0] = F[1] = f[1] = _f[1] = 1;
    for (int i = 2; i < N; i++)
    {   
        f[i] = (f[i-1] + f[i-2]) % mod; _f[i] = inv(f[i]);
        F[i] = 1;
        if (!p[i])
        {
            p[i] = i;
            pri[++cnt] = i;
            mo[i] = -1ll;
        }
        for (int j = 1; j <= cnt; j++)
        {
            if (i * pri[j] >= N) break;
            p[i * pri[j]] = pri[j];
            if (i % pri[j] == 0)
            {
                mo[i * pri[j]] = 0;
                break;
            }
            mo[i * pri[j]] = -mo[i];
        }
    }
    for (int i = 2; i < N; i++)
        for (int j = i; j < N; j += i)
            mul(F[j], (mo[j / i] == 1 ? f[i] : (mo[j / i] == -1 ? _f[i] : 1)));
    for (int i = 2; i < N; i++) mul(F[i], F[i-1]);
}
int t, n, m;
int solve(int n, int m)
{
    int ret = 1;
    for (int i = 1, j; i <= min(n, m); i = j + 1)
    {
        j = min(n / (n / i), m / (m / i));
        mul(ret, ksm(1ll * F[j] * inv(F[i-1]) % mod, 1ll * (n / i) * (m / i) % (mod - 1)));
    }
    return ret;
}
int main()
{
    ios::sync_with_stdio(false); cin.tie(0); cout.tie(0);
    cin >> t; sieve();
    while (t--)
    {
        cin >> n >> m;
        cout << solve(n, m) << "\n";
    }
    return 0;
}

P3172

把 N,K,L,H 写成 n,k,a,b。

记 f(x) 表示 \gcd = x 的答案，发现这样无法求，考虑 F(x) 表示 x | \gcd 的方案，则

F(x) = \sum_{x|L} f(L) \\ f(x) = \sum_{x|L} \mu(\frac{L}{x})F(L)

显然只有当 x=1 时 f(x) 可以整除分块，故让 a \gets \lceil \dfrac{a}{k} \rceil，b \gets \lfloor \dfrac{b}{k} \rfloor，现在就是求 f(1) = f(x) = \sum_{L=1}^{+\infty} \mu(L)F(L) 了。

易知

F(x) = (\lfloor \dfrac{b}{x} \rfloor - \lfloor \dfrac{a-1}{x} \rfloor)^n

整除分块即可，

吗？

可以发现，b 可以高达 10^9，所以你不能线筛。

那就杜教筛。复杂度 O(n^{\frac{2}{3}})。

const int N = 5e6 + 5, mod = 1e9 + 7;
typedef long long ll;
ll cnt, p[N], pri[N], mo[N], mosum[N];
void sieve()
{
    mo[1] = 1;
    for (int i = 2; i < N; i++)
    {   
        if (!p[i])
        {
            p[i] = i;
            pri[++cnt] = i;
            mo[i] = -1ll;
        }
        for (int j = 1ll; j <= cnt; j++)
        {
            if (i * pri[j] > N) break;
            p[i * pri[j]] = j;
            if (i % pri[j] == 0ll)
            {
                mo[i * pri[j]] = 0ll;
                break;
            }
            mo[i * pri[j]] = -mo[i];
        }
    }
    for (int i = 1; i < N; i++) mosum[i] = mosum[i - 1] + mo[i];
}
unordered_map <int, int> m;
ll cmo(ll x)
{
    if (x < N) return mosum[x];
    if (m.count(x)) return m[x];
    ll l, r;
    ll ans = 1;
    for (l = 2ll; l <= x; l = r + 1ll)
    {   
        r = x / (x / l);
        ans -= (r - l + 1ll) * cmo(x / l);
    }
    m[x] = ans;
    return ans;
}
int n, k, a, b;
int ksm(int a, int b) 
{
    int ret = 1;
    while (b)
    {
        if (b & 1) ret = 1ll * ret * a % mod;
        a = 1ll * a * a % mod;
        b >>= 1;
    }
    return ret;
}
int solve(int a, int b)
{
    int ret = 0;
    for (int i = 1, j; i <= b; i = j + 1)
    {
        if (a < i) j = (b / (b / i));
        else j = min(a / (a / i), b / (b / i));
        ret = (ret + 1ll * ksm(b / i - a / i, n) * (cmo(j) - cmo(i - 1) + mod) % mod) % mod;
    }
    return ret;
}
int main()
{
    ios::sync_with_stdio(false); cin.tie(0); cout.tie(0);
    cin >> n >> k >> a >> b; sieve();
    b /= k; a = (a + k - 1) / k; a--;
    cout << solve(a, b);
    return 0;
}

QwQ。