莫比乌斯反演笔记

Posted by

前置知识

正因子求和

\sum_{d|n}表示对n的所有正因子求和,例如\sum_{d|8}=1^2+2+2+4^2+8^2

积性函数

对于所有的正整数中选取x,y,使得x,y互素,并且满足f(x)f(y)=f(xy)的函数f,就是积性函数。

如果当x,y不互素时,仍然满足上述条件,那么f是完全积性函数。

一些积性函数可以参见例子

莫比乌斯函数

莫比乌斯函数

筛法实现模板

void seive()
{
    mu[1] = 1;
    for ( int i = 2; i <= 10010; i ++ )
    {
        if ( !flag[i] ) prime[++cnt] = i, mu[i] = -1;
        for ( int j = 1; j <= cnt && i * prime[j] <= 10010; j ++ )
        {
            flag[i*prime[j]] = 1;
            if ( i % prime[j] == 0 ) break;
            mu[i*prime[j]] = -mu[i];
        }
    }
    for ( int i = 1; i <= 10010; i ++ )
        sum[i] = sum[i-1] + mu[i];
}

莫比乌斯反演

常用的两种形式
如果F(n)=\sum_{d|n}f(d),那么f(n)=\sum_{d|n}\mu(d)F(\frac{n}{d})
如果F(n)=\sum_{n|d}f(d),那么f(n)=\sum_{n|d}\mu(\frac{d}{n})F(d)

具体数学中定义的形式
d \backslash n表示d能被n整除
如果 F(n)=\sum_{d\backslash n} f(d),那么f(n)=\sum_{d\backslash n} \mu(d)F(\frac{n}{d})=\sum_{d\backslash n} \mu(\frac{n}{d})F(d)

操作步骤

  • 构造f函数,通常是设满足题目条件的个数为f(d)
  • 构造F函数,通常是即F(n)=\sum_{d|n}f(d)
  • 得到反演公式f(n)=\sum_{d|n}\mu(d)F(\frac{n}{d})f(n)=\sum_{n|d}\mu(\frac{d}{n})F(d)
  • F(d)进行计算,即满足题目条件的那些d的整数倍,那么F(d)就可以变成另外一种形式.
  • (可选)分块处理,对f(n)求和计算出答案.

例题

POJ Sky Code为例
题目给出n和n个正整数,问能找出多少个不同的四元组(a,b,c,d)使得该四元组的最大公因数为1.

f(x)是满足gcd(a,b,c,d)=d的个数
F(n)为有多少个四元组满足gcd(a,b,c,d)=d的整数倍,那么有F(n)=\sum_{d|n}f(d),通过莫比乌斯反演得到f(n)=\sum_{n|d}\mu(\frac{d}{n})F(d),根据题意求f(1),那么也就是求f(1)=\sum_{1|d}\mu(\frac{d}{1})F(d)=\sum_{d=1}^{\infty}\mu(d)F(d)

考虑F(d)为有多少个四元组满足gcd(a,b,c,d)=d的整数倍,对原数组中计算到能被d整除的数的个数m[i],然后再\binom{m[i]}{4},就可以求出F(d)了。
最后全部累加一遍,就可以得到答案。

//#include <bits/stdc++.h>
#include <cstdio>
#include <math.h>
#include <cstring>

using namespace std;

typedef long long ll;

int n;
int a[10020];
ll mu[10020], flag[10020];
ll prime[10020], sum[10020];
ll my_count[10020], c[10020];
int cnt;

void seive()
{
    mu[1] = 1;
    for ( int i = 2; i <= 10010; i ++ )
    {
        if ( !flag[i] ) prime[++cnt] = i, mu[i] = -1;
        for ( int j = 1; j <= cnt && i * prime[j] <= 10010; j ++ )
        {
            flag[i*prime[j]] = 1;
            if ( i % prime[j] == 0 ) break;
            mu[i*prime[j]] = -mu[i];
        }
    }
    for ( int i = 1; i <= 10010; i ++ )
        sum[i] = sum[i-1] + mu[i];
}

void init_cnt()
{
    memset( my_count, 0, sizeof(my_count) );
    for ( int i  = 1; i <= n; i ++ )
    {
        int x = a[i];
        int m = sqrt(x);
        for ( int j = 1; j <= m; j ++ )
        {
            if ( x % j == 0 ) my_count[j] ++, my_count[x/j] ++;
        }
        if ( m * m == a[i] ) my_count[m] --;
    }
}

void init_c()
{
    c[1] = 1;
    for ( int i = 2; i <= 10010; i ++ )
    {
        c[i] = c[i-1] * i;
    }
}

ll F(int m)
{
    if (m == 0)return 0;
    return (ll)m*(m - 1)*(m - 2)*(m - 3) / 24;
}

int main()
{
    init_c();
    seive();
    while ( ~scanf("%d", &n)  )
    {
        for ( int i = 1; i <= n; i ++ ) scanf("%d", &a[i]);
        init_cnt();
        ll ans = 0;
        for ( int i = 1; i <= 10010; i ++ )
        {
            ans += (ll)(F(my_count[i]) * mu[i]);
        }
        printf("%lld\n", ans);
    }

    return 0;
}

Leave a Reply

电子邮件地址不会被公开。 必填项已用*标注