莫比乌斯反演

反演

​ 举个例子：

$$已知g(n)的前缀和f(n)=\sum_{i=1}^ng(i)，则可以通过f反求g =>g(n)=f(n)-f(n-1)$$

​ 像这样的过程就是反演，换言之，就是一个函数$f$由$g$推导而来，已知$f$,求$g$

1、提出一个问题

​

$$已知g(n)的因数和f(n)=\sum_{d|n}g(d),如何通过f反求g$$

1.1、先模拟出f

$$\begin{cases} f(1)=g(1)\\ f(2)=g(1)+g(2)\\ f(3)=g(1)+g(3)\\ f(4)=g(1)+g(2)+g(4)\\ f(5)=g(1)+g(5)\\ f(6)=g(1)+g(2)+g(3)+g(6)\\ f(7)=g(1)+g(7)\\ ..... \end{cases}$$

​ 通过以上的等式，可以将相应的g求出：

$$\begin{cases} g(1)=f(1)\\ g(2)=f(2)-f(1)\\ g(3)=f(3)-f(1)\\ g(4)=f(4)-f(2)\\ g(5)=f(5)-f(1)\\ g(6)=f(6)-f(3)-f(2)+f(1)\\ g(7)=f(7)-f(1)\\ ..... \end{cases}$$

​ 通过这个模拟可以看出一丝小小的规律：

​ 将$n$先改写成$p_1^{\alpha_1}p_2^{\alpha_2}p_3^{\alpha_3}….p_k^{\alpha_k}$，然后代入$g(n)$,根据上面的模拟将$g(n)$展开

$$g(n)=g(p_1^{\alpha_1}p_2^{\alpha_2}p_3^{\alpha_3}....p_k^{\alpha_k})\\ =f(p_1^{\alpha_1}p_2^{\alpha_2}p_3^{\alpha_3}....p_k^{\alpha_k})-f(p_1^{\alpha_1-1}p_2^{\alpha_2}p_3^{\alpha_3}....p_k^{\alpha_k})-f(p_1^{\alpha_1}p_2^{\alpha_2-1}p_3^{\alpha_3}....p_k^{\alpha_k}) -f(p_1^{\alpha_1}p_2^{\alpha_2}p_3^{\alpha_3-1}....p_k^{\alpha_k}) \\-f(p_1^{\alpha_1}p_2^{\alpha_2}p_3^{\alpha_3}....p_k^{\alpha_k-1}) \\+f(p_1^{\alpha_1-1}p_2^{\alpha_2-1}p_3^{\alpha_3}....p_k^{\alpha_k})+f(p_1^{\alpha_1}p_2^{\alpha_2}p_3^{\alpha_3-1}....p_k^{\alpha_k})+..... +f(.....p_{i}^{\alpha_i-1}....p_j^{\alpha_j-1}.....)$$

1.2、浓缩展开式

​ 我们通过上面的模拟可以确定最终反演的$g$，一定是长成上面那个复杂的等式的样子，但是这样一个等式不怎么好记，因此我们将上述展开式进行抽象：

$$g(p_1^{\alpha_1}p_2^{\alpha_2}p_3^{\alpha_3}....p_k^{\alpha_k})=\sum_{s\subseteq\{1,2,3,..k\}}(-1)^{|k|}f(p_1^{\alpha_1-i_1}p_2^{\alpha_2-i_2}p_3^{\alpha_3-i_3}....p_k^{\alpha_k-i_k}) \\ 其中 \begin{cases} i_j=1,j\in S \\ i_j=0,其他 \end{cases}$$

碰巧！！！在数论有这样一个函数——莫比乌斯函数$\mu$

1.3、莫比乌斯函数

​

$$\mu(n)= \begin{cases} 1,n=1\\ (-1)^k,n=p_1^{\alpha_1}p_2^{\alpha_2}p_3^{\alpha_3}....p_k^{\alpha_k}且\alpha_i=1 \\ 0，其他 \end{cases}$$

​ 这个时候我们便可以将$(-1)^{|k|}$用$\mu（n）$进行替换

​

$$g(p_1^{\alpha_1}p_2^{\alpha_2}p_3^{\alpha_3}....p_k^{\alpha_k})=\sum_{s\subseteq\{1,2,3,..k\}}\mu(p_1^{i_1}p_2^{i_2}p_3^{i_3}....p_k^{i_k}) f(p_1^{\alpha_1-i_1}p_2^{\alpha_2-i_2}p_3^{\alpha_3-i_3}....p_k^{\alpha_k-i_k}) \\ 再化简：\\ g(n)=\sum_{d|n}\mu(\frac{n}{d})f(d)$$

​ 举个例子便可知道是否正确了：

​

$$例子：g(6)=\mu(1)f(6)+\mu(2)f(3)+\mu(3)f(2)+\mu(6)f(1)\\ =f(6)-f(3)-f(2)+f(1)$$

2、莫比乌斯反演

​ 通过上述过程我们便证明了定理1：

2.1、定理一：

$$已知g(n)的因数和f(n)=\sum_{d|n}g(d),则g(n)=\sum_{d|n}\mu(\frac{n}{d})f(d)$$

​ 除此之外，还有一个常见的

2.2、定理二：

$$有两个函数f,g,且存在正整数N，对于所有n>N,f(n)=g(n)=0,则: \\ f(n)=\sum_{n|m,m\leq N}g(m) \iff g(n)=\sum_{n|m,m\leq N}\mu(\frac{m}{n})f(m)$$

​ 证明的方式和定理一异曲同工。

3、莫反好题：

3.1、牛客14648

题目连接：https://ac.nowcoder.com/acm/problem/14648

​

$$题目大意：有两个序列a,b，问：有多少个序列<x,y>满足gcd(x,y)=1 且 a_{b_x}=b_{a_y}$$ $$这道题的意思等价于求\sum_{1\leq x\leq n,1\leq y\leq n}[a_{b_x}=b_{a_y}][gcd(x,y)=1] \\ 那么设f(D)=\sum_{1\leq x\leq n,1\leq y\leq n}[a_{b_x}=b_{a_y}][gcd(x,y)=D] \\ 那么根据莫比乌斯反演定理二有： g(d)=\sum_{d|D}f(D)\\ 由于 \begin{cases} D|x \\ D|y \end{cases} \\所以有 \begin{cases} d|x \\ d|y \end{cases} \\所以g(d)=\sum_{ d|x, d|y }[a_{b_x}=b_{a_y}] \\ 这样一来g函数就很好算（时间复杂度O(nlogn)），然后利用反演即可求出f函数 \\ 由定理二可得 f(D)=\sum_{D|d,d\leq n}\mu(\frac{d}{D})g(d) \\由于我们所求为f(1)，因此将D=1代入即可 \\ 所以f(1)=\sum_{1|d,d\leq n}\mu(d)g(d) \\然可以利用\mu函数是积性函数的特性将莫比乌斯函数打表打出来(时间复杂度O(n)) \\所以最后的时间复杂度为O(nlogn+n)$$

#include <iostream>
#include <algorithm>
#include <cmath>
#include <vector>
#include <string.h>
using namespace std;
#define int long long
const int maxn = 1e5 + 10;
vector<int> primes;
bool isnp[maxn];
int mu[maxn];
int a[maxn];
int b[maxn];
void EularMu(int n)//求莫比乌斯函数
{
    mu[1] = 1;
    for (int i = 2; i <= n; i++)
    {
        if (!isnp[i])
        {
            primes.push_back(i);
            mu[i] = -1;
        }
        for (auto p : primes)
        {
            if (i * p > n)
                break;
            isnp[p * i] = 1;
            if (i % p == 0)
            {
                mu[i*p]=0;
                break;
            }
            mu[i * p] = -mu[i] ;
        }
    }
}
int num[maxn];
int getg(int d, int n)//求g函数
{
    int ans = 0;
    for (int x = d; x <= n; x += d)
    {
        num[a[b[x]]]++;
    }
    for (int y = d; y <= n; y += d)
    {
        ans += num[b[a[y]]];
    }
    for (int x = d; x <= n; x += d)
    {
        num[a[b[x]]]--;
    }
    return ans;
}
signed main()
{
    ios::sync_with_stdio(false);
    cin.tie(0);
    cout.tie(0);
    EularMu(maxn);
    int n;
    cin >> n;
    for (int i = 1; i <= n; i++)
    {
        cin >> a[i];
    }
    for (int i = 1; i <= n; i++)
    {
        cin >> b[i];
    }
    int ans = 0;
    for (int i = 1; i <= n; i++)
    {
        ans += mu[i] * getg(i, n);
    }
    cout << ans << endl;
    return 0;
}

3.2、AtCoder Grand Contest 038C

题目连接：https://atcoder.jp/contests/agc038/tasks/agc038_c

$$给一个长度为 n的序列 a_0,a_1,a_2,...,a_{n−1}\\ 求∑_{i=0}^{n−2}∑_{j=i+1}^{n−1}lcm(a_i,a_j)\ mod\ 998244353$$ $$解：首先利用和式的性质：\sum_{i=0}^{n-1}\sum_{j=0}^{n-1}a_{ij}=2\sum_{i=0}^{n-2}\sum_{j=i+1}^{n-1}a_{ij}+\sum_{i=0}^{n-1}a_{ii} \\这个性质可以想象成一个大矩阵a_{ij}求和，这个和式等价于<整个矩阵的和= 右上角+左下角+对角线> \\ ∴∑_{i=0}^{n−2}∑_{j=i+1}^{n−1}lcm(a_i,a_j)=\frac{∑_{i=0}^{n−1}∑_{j=0}^{n−1}lcm(a_i,a_j)-\sum_{i=0}^{n-1}a_i}{2} \\ 由于\sum_{i=0}^{n-1}a_i很好算，所以这个问题等价于如何求∑_{i=0}^{n−1}∑_{j=0}^{n−1}lcm(a_i,a_j) \\又∵lcm(a_i,a_j)=\frac{a_i*a_j}{gcd(a_i,a_j)} \\∴∑_{i=0}^{n−1}∑_{j=0}^{n−1}lcm(a_i,a_j)=∑_{i=0}^{n−1}∑_{j=0}^{n−1}\frac{a_i*a_j}{gcd(a_i,a_j)} \\∴∑_{i=0}^{n−1}∑_{j=0}^{n−1}lcm(a_i,a_j)=\sum_{d=1}^{max(a_0,a_2..a_{n-1})}\frac{1}{d}∑_{i=0}^{n−1}∑_{j=0}^{n−1}a_i*a_j[gcd(a_i,a_j)=d] \\这样一来，这个问题感觉就变得套路化了！ \\上套路： \\设f(d)=∑_{i=0}^{n−1}∑_{j=0}^{n−1}a_i*a_j[gcd(a_i,a_j)=d] \\做一个倍数和函数g\\ 则g(d)=\sum_{d|D}f(D)=\sum_{d|D}∑_{i=0}^{n−1}∑_{j=0}^{n−1}a_i*a_j[gcd(a_i,a_j)=D] \\ ∵D|a_i,D|a_j \\ ∴d|a_i,d|a_j \\ ∴g(d)=\sum_{d|D}∑_{i=0}^{n−1}a_i[d|a_i]∑_{j=0}^{n−1}a_j[d|a_j] \\ ∴g(d)=\sum_{d|D}(∑_{i=0}^{n−1}a_i[d|a_i])^2 \\这样一来通过枚举d就可以在O((max(a_i))log(max(a_i)))的时间复杂度下通过打表打出g函数 \\利用莫比乌斯反演定理二可得： \\f(d)=\sum_{d|D}\mu(\frac{D}{d})g(D) \\∴∑_{i=0}^{n−1}∑_{j=0}^{n−1}lcm(a_i,a_j)=\sum_{d=1}^{max(a_0,a_2..a_{n-1})}\frac{1}{d}\sum_{d|D}\mu(\frac{D}{d})g(D) \\ ∴∑_{i=0}^{n−2}∑_{j=i+1}^{n−1}lcm(a_i,a_j)=\frac{\sum_{d=1}^{max(a_0,a_2..a_{n-1})}\frac{1}{d}\sum_{d|D}\mu(\frac{D}{d})g(D)-\sum_{i=0}^{n-1}a_i}{2} \\ 积性函数\mu可以用欧拉筛在线性的时间o(n)内打表打出来，g函数也是在O(nlogn)的时间复杂度打表预处理 \\同理最后这个\sum_{d=1}^{max(a_0,a_2..a_{n-1})}\frac{1}{d}\sum_{d|D}\mu(\frac{D}{d})g(D)也可以在o(nlogn)的时间复杂度求出来 \\因此这个复杂的玩意儿可以在O(nlogn)的时间复杂度下解决！！！ \\最后还需要注意一下，此题还有有一个取模的操作 \\对于此题分数取模，最好预先打表处理逆元，用一个inv数组保存，不然会超时$$

下面贴上代码（加了注释）：

#include <iostream>
#include <algorithm>
#include <cmath>
#include <vector>
using namespace std;
#define int long long
const int maxn = 1e6 + 10;
const int mod = 998244353;
vector<int> primes;
bool isnp[maxn];
int mu[maxn];
int a[maxn];
int gf[maxn];
int cnt[maxn];
int inv[maxn];
void Eular_mu(int n)//欧拉筛求莫比乌斯函数
{
    mu[1] = 1;
    for (int i = 2; i <= n; i++)
    {
        if (!isnp[i])
        {
            primes.push_back(i);
            mu[i] = -1;
        }
        for (auto p : primes)
        {
            if (i * p > n)
                break;
            isnp[i * p] = 1;
            if (i % p == 0)
            {
                mu[i * p] = 0;
                break;
            }
            mu[i * p] = -mu[i];
        }
    }
}
int qp(int a, int x)//用快速幂
{
    int ans = 1;
    while (x)
    {
        if (x & 1)
        {
            ans = ans * a % mod;
        }
        a = a * a % mod;
        x >>= 1;
    }
    return ans;
}
signed main()
{
    ios::sync_with_stdio(false);
    cin.tie(0);
    cout.tie(0);
    int n;
    cin >> n;
    int maxx = 0;
    int sum = 0;
    for (int i = 0; i < n; i++)
    {
        cin >> a[i];
        maxx = max(a[i], maxx);
        sum += a[i];
        cnt[a[i]]++;//记录a[i]的数量
    }
    Eular_mu(maxn);
    for(int i=1;i<=maxn;i++){//快速幂套费马小定理初始化逆元
        inv[i]=qp(i,mod-2);
    }
    for (int d = 1; d <= maxx; d++)//求g函数
    {
        for (int i = d; i <= maxx; i += d)
        {
            gf[d] = (gf[d] + i * cnt[i] % mod) % mod;
//由于是枚举d的倍数,如果符合条件一定存在a[i]和此时d的倍数相等，这个时候cnt的作用体现出来
//利用cnt便可以直接求和
        }
        gf[d] = gf[d] * gf[d] % mod;
    }
    int ans = 0;
    for (int d = 1; d <= maxx; d++)
    {
        for (int i = d; i <= maxx; i += d)
        {
            ans = (ans + mu[i / d] * gf[i] %mod* inv[d] % mod) % mod;
//同理枚举d的倍数，原理和g函数相似，最后直接套推出来的公式
        }
    }
    cout << (ans - sum%mod+mod ) % mod * inv[2] % mod << endl;
//sum容易吧mod爆了，一定要取模
    return 0;
}