• 莫比乌斯反演入门


          转载自-----   http://blog.csdn.net/qw4990/article/details/14055183

    这个文章主要讲一下ACM中1个常用的莫比乌斯反演公式,看到很多博客上面公式是有,但是都没证明,《组合数学》上的证明又没看懂,

    就自己想了种证明方法,觉得比《组合数学》的证明简单些,就写一下,希望对初学莫比乌斯反演的同学有帮助。

    PS:下面公式出现的sigma是累加,另外建议大家看的时候 把公式在纸上写出来!

    一:什么是莫比乌斯反演

    简单点的说,就是先给出一个函数 F(n) ,然后再由 F(n)定义一个新函数 G(n)

    其中   G(n) = sigma(F(d)) (其中d被“包含”于n)  

    然后 现在我们不知道 F(n)的值 , 却知道 G(n), 接着我们就可以通过 反演由G(n)反向得到F(n)

    什么叫 (其中d被“包含”于n) ?以及怎么理解反演? 通过下面的几个例子说明

    例1:

    我们直接定义 G(n)=sigma(F(i)) (1<=i<=n)    {这里的每个F(i),相对于G(n)实际上就是一种包含关系了!!}

    然后我们现在已经知道 G(n)=n*(n+1)/2;

    接下来 我们要通过 G(n)反向得到F(n) 的过程,就是反演

    当然,这个问题很简单,很容易都可以看出来 F(n)=n ~~

    --------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------

    例2:

    我们先令  S,X 都表示集合 比如 S={1,4,6} X={2} 等  并令|S|表示 S中元素的个数

    接着定义 集合上的函数 F(S)   /*具体怎么定义不用管,我们只需要知道有这么一个关于集合的函数F就好了 :) */ 

    然后再定义 G(S)=sigma(F(X)) (其中X是S的子集)   {这里也是一种包含关系,集合的包含!!}

    接着我们不知道F(S),想通过G(S) 来得到 F(S)

    这个问题相对于例1就复杂多了,但实际上我们已经有现成的关于集合包含的莫比乌斯反演公式了  :)

    F(S)=sigma((-1)^(|S|-|X|) * G(X))  (其中X是S的子集)

    是不是感觉有点神奇?

    大家可以自己写个程序来验证一下。

    下面就是我的验证程序:

    我定义 F(S)=|S| , 然后先 计算出 F(S) ,接着 计算出 G(S) , 然后 比较由G(S)反演得到的 F(S)和 |S| 的大小

    下面是 我的程序

    #include <iostream>
    #include <math.h>
    using namespace std;
    
    #define base 10
    #define REP(i,n) for(int i=0;i<(n);i++)
    int F[1<<base],G[1<<base];
    // 集合用二进制表示 base表示集合最多10个元素
    
    int Cal(int x){ // 计算 |x|
        int sum=0;
        while(x) sum+=(x&1),x/=2;
        return sum;
    }
    
    int main(){
        REP(S,1<<base) F[S]=Cal(S); // 计算出最开始的F(S)
        REP(S,1<<base){   // 计算G(S)
            G[S]=0;
            for(int X=S;X;X=(X-1)&S) G[S]+=F[X]; //用X遍历S集合
        }
        REP(S,1<<base){     // 计算反演的 F(S)
            F[S]=0;
            for(int X=S;X;X=(X-1)&S)
                F[S]+=(int)pow(-1,Cal(S)-Cal(X))*G[X];
        }
        bool flag=1;    // 验证一下
        REP(S,1<<base)
            if(F[S]!=Cal(S)) flag=0;
        if(flag) cout<<"YES"<<endl;
        else cout<<"NO"<<endl;
    }

    最后得到的结果 当然是 YES 咯!:)

    关于这个 反演公式 的证明,先不要着急,看完文章过后,你自己都能摸索着证明了!!

    现在先大概理解反演是个什么就行了!!

    -----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------

    例3:

    先令  d|n 表示 d能整除n  比如 2|4 (=.=)

    定义 关于 整数 的函数 F(n)

    然后 定义 G(n)=sigma(F(d)) (其中d|n)  

    上面的这种包含关系就更复杂了,只有当d是n的因子的时候,F(d)才会被包含在G(n)中。

    不过这种 包含关系 在 ACM中遇到的最多,所以我会详细讲一下这种类型的 反演。

    相信明白了这个过后,例2的反演也能够自己证明了。

    具体的讲解见下一章节 :)

    二:一类反演

    这个一类反演就是例3中的那一类咯= =

    我先直接给出结论吧

    原式 :  G(n)=sigma(F(d))  (其中d|n)  

    反演公式:   F(n)=sigma(U(n/d)*G(d))   这里U是一个函数,他是每一项 G(d) 的系数,他的定义见下面

      (强烈建议关于U的定义这一段可以先跳过,先认为他是G的系数就行了,可以跳到下面红字位置)

        (1).U是一个关于整数的函数

        (2).U[x] = 1 当且仅当 x能够分解成偶数个不同质数的乘积  (其中1不能被分解,所以1的分解出的质数个数是0,所以U[1]=1)

        (3).U[x] = -1 当且仅当 x能够分解成奇数个不同质数的乘积

        (4).U[x] = 0  除开(2),(3)的其他情况

        看上面关于U的定义可能有点看晕了,通俗一点的说

    对于一个 x , 分解因式过后 有  x=(p1^e1)*(p2^e2)...*(pr^er)

    如果 ei中(1<=i<=r)有一个数ei大于1  那么  U[x] = 0;

    不然的话   U[x] = (-1)^r

    依旧来两个例子(我最喜欢举例子了 = =)

    U[1]=1;定义中的说明

    U[2]=-1;    分解式 2=2;

    U[6]=1;   分解式 6=2*3

    U[9]=0;   9=3^2; 出现了e>1

    U[12]=0;   12=2^2*3;

    跳到这里 :)

    上面就是关于这类反演公式的定义,不要头晕= =,继续往下看吧

    在证明之前,我们先想一下,为什么反演公式会是   F(n)=sigma(U(n/d)*G(d))   这样的型式?

    依旧通过例题来找规律 (^ ^)

    我们令 n=6;

    那么 在计算 F(6)的时候,我们会用到 G(1) G(2) G(3) G(6)

    我们考察者4个G

    G(1) = F(1)

    G(2) = F(1)+F(2)

    G(3) = F(1)+F(3)

    G(6) = F(1)+F(2)+F(3)+F(6)

    观察上面可以发现 每个 G(n)都是由一些F(d)累加得到的 

    当我们需要逆向有G得到F(n)时,只需要将一些 与 F(n) 有关的 G进行容斥!!!!! 最终组合得到F(n)!!!

    比如 F(6) = G(6)-G(2)-G(3)+G(1) !!!!

    有些神奇!! 不过这类莫比乌斯反演的实质也就是容斥原理的应用!!

    那么我们现在知道为什么 这类反演公式会是 这个形式了,而且对其原理也有了更深的理解,现在该想一想公式的细节了。

    既然我们知道要得到 F(n) ,只需要将与其相关的 G进行容斥就可以,那么剩下的问题就是每个G的系数!!!

    我们以 求解 F(6)为例子来说明 ,并定义一个系数函数 H(d,n).

    其中 H(d,n)表示 求解F(n)时,G(d)的系数  (其中d|n)

    所以可以得到这个式子 F(6) = H(6,6)*G(6)+H(2,6)*G(2)+H(3,6)*G(3)+H(1,6)*G(1)

    我们用 a,b,c,d分别替代 四个H(6,6),H(2,6),H(3,6),H(1,6),并且把对应的G用F表示出来,得到

    F(6)=a*(F(6)+F(3)+F(2)+F(1))+b*(F(2)+F(1))+c*(F(3)+F(1))+d*F(1),再变形一下,又有

    F(6)*(a-1)+F(3)*(a+c)+F(2)*(a+b)+F(1)*(a+b+c+d)=0,把F(6),F(3),F(2),F(1)当作不同的元,则得到了下面的方程组!!!

    a-1==0

    a+c==0

    a+b==0

    a+b+c+d==0

    由此发现,四个未知数,四个方程,只需要解出方程,就能知道对于G的系数。

    再深入的想一下,对于每个 F(n),假设他的因子数为,m,则通过这种方式,总能设出m个未知数,m个方程,

    这样总能找到解,而这也为莫比乌斯反演的可能性作出了解释!!

    现在我们要证明一个结论,即使H(a,b)==H(1,b/a)!!这个结论很重要,具体分析见下 :)

    我们以求解 F(8)为例子,与F(8)相关的 H 有 ,H(8,8),H(4,8),H(2,8),H(1,8)

    F(8)=H(8,8)*G(8)+H(4,8)*G(4)+H(2,8)*G(2)+H(1,8)*G(1)

    首先看 H(8,8),其值可以直接确定,因为把F(8)当作元的话,左边一个F(8),而在右边F(8)只在G(8)中出现,所以H(8,8)==1

    同理 对于 F(n),其G(n)的系数H(n,n)==1,所以H(8,8)==H(1,1)

    再来看H(4,8),,首先想,F(4)在哪些地方出现,发现 在G(8)和G(4)出现,因为左边不含F(4),而前面G(8)的系数又已经确定,

    所以这里H(4,8)*G(4)的作用就是为了抵消前面G(8)的代换中,出现的F(4),所以 H(4,8)==-H(8,8)==-H(2,2)==H(1,2),{H(1,2)==-H(1,1)请大家自己验证一下}

    同理对于H(2,8),他是为了抵消前面在G(8)和G(4)中出现的F(2),所以H(2,8)相当于受到H(4,4)和H(2,4)的影响(假设这个结论对n==4也成立,H(2,4)==H(1,2)),

    所以H(2,8)==H(1,4)

    找到规律过后,总结一下,假设n的因子有 d1,d2,d3...dm 其中 d1>d2>d3...>dm

    我们依次确定H(di,n)的值,当我们在确定H(di,n)的值时,前面的值已经确定,即H(dj,n)(j<i)的值已经确定,

    H(di,n)会受到前面一些H(dj,n)的影响,当且仅当 dj>di且 di|dj 。

    假设 H(a,b)==H(1,b/a)对前面的 H(dj,n)和 所有的H(k,m)其中m<n 已经成立(首先对于H(n,n)已经成立),那么有  H(dj,n)==H(1,n/dj)==H(dj/di,n/di)

    这样就把前面对H(di,n)造成影响的H由 H(dj,n)转为了 H(dj/di,n/di) ,所以H(di,n) == H(1,n/di)

    既然 H(a,b) 都可以 写成 H(1,b/a) , 于是我们把H的第一个元素略去,简写为 H(x)

    说到这里,就可以把H和U联系起来了,其实 U(x) = H(x) = H(1,x)  

    再来,我们就可以给U(x)赋予一个更具体的意义, U(x)表示在计算 F(x)时,G(1)的系数!!(因为U(x)==H(1,x))

    接下来,我们来尝试一下,如何用上面那个U(x)的新意义,来计算U(x)的值!!

    首先需要明确2点! 

    一是G(x)中,一定包含一个F(1),因为 1|x

    二是,F(1)==G(1)

    (0).如果 x==1

    因为 F(1)==G(1) 所以 U[1]=1;

    (1).假设 x 是一个 质数   

    F(x) = U(1)*G(x)+U(x)*G(1)

    带入U(1) == 1,  因为G(x)中含有一个F(1),而左边不含F(1),所以我们需要利用G(1)来消去F(1)

    所以得到 U(x)=-1

    (2).假设 x 可以写成2个不同质数的乘积  x=p*q

    那么  F(x)=U(1)*G(pq)+U(q)*G(p)+U(p)*G(q)+U(x)*G(1)

    这里 U(1),U(p),U(q) 就是前面2种情况

    带入系数,因为左边没有 F(1),所以为了抵消右边的F(1),我们需要令 U(x)=1;

    (3).假设 x 可以写成3个不同质数的乘积  x=p*p1*p2  我们令  z = p1*p2

    F(x) = U(1)*G(pz)+U(z)*G(p)+U(p)*G(z)+U(x)*G(1);

    其中 U(1),U(p),U(z) 分别为前面几种情况,带入过后 ,为抵消F(1)   得到 U(x)=-1

    由此可以相同的方式向下递推,得到第一条结论

    如果 x = p1*p2...*pr , 其中pi是互异的质数,那么 U[x] = (-1)^r  -----------------------   1!!

    (4).假设 x 可以写成一个质数的平方  x=p^2

    F(x) = U(1)*G(x)+U(p)*G(p)+U(x)*G(1)

    带入系数 得到 U(x)=0;

    (5).假设 x 可以写成一个质数的三次方  x=p^3

    F(x) = U(1)*G(x)+U(p)*G(p^2)+U(p^2)*G(p)+U(x)*G(1)

    带入系数后   U(x)=0;

    由此可用相同方式向下递推,得到第二条结论

    如果 x = p^e (e>1) U[x] = 0; -------------------------- 2!!

    (6).假设 x 可写成 x = p^e*q  其中p,q为不同质数,e>1

    F(x) = U(1)*G(x)+U(q)*G(p^e)+U(p^e)*G(q)+U(x)*G(1)

    带入系数后   U(x) = 0;

    由此可继续向下递推,得到第二条结论的加强版!!

    如果 x = p^e*z 其中p为质数, z为任意数,e>1  那么  U[x] = 0 ----------------------2!!

    由此,我们得到了 U[x] 的计算方法!!即是U定义中给出的那样!!(没看定义的同学此时再跳回去看吧)

    三:应用

    得到了公式,也知道了他是怎么来的,现在就用一个应用来加深理解吧  :)

    首先我们要给出 第二部分 中那个公式的另外一种形式 = =  我们把它称为形式二吧~

    原式 :  G(n)=sigma(F(d))  (其中n|d,d<=N)  

    反演公式:   F(n)=sigma(U(d/n)*G(d))  (其中n|d,d<=N) 这里U[x]的计算方式和上面的相同!!

    注意上面的 n|d 和 d/n  和上面是相反的

    证明方法和上面差不多,大致说一下

    还是先设置一个系数函数 H(d,n) 表示求解 F(n)时 G(d)出现的次数,

    接着 用与上面类似的方法变化H(d,n) 为 H (d/n,1) ---> H(x,1)

    则联系 U(x) == H(x,1)  表示  在计算 F(1)时,G(x)的系数

    以 x 为质数为例子,由于 G(1)=F(1)+F(2)....+F(N)

    F(1) = G(1)+U(2)G(2)...+U(x)G(x)...+U(N)G(N)

    因为 x 为质数 所以 F(x)这一项 只在G(1)里面出现了一次,而其他地方只会在 G(x)出现

    所以我们需要让 U(x)=-1 来抵消 F(x)

    剩下的步骤就和上面差不多了,分类讨论一下,就可以求出这种情况下的U的计算方式,和上面相同!!

    ---------------------------------------------------------------------------------------------------------------------------------------------------------------------

    接下来就真正的开始演示怎么用 莫比乌斯反演 简化计算了 !!

    看下面这个问题!

    给出a,b 其中 (1<=a,b<=10^6) 

    求满足条件的 x,y 的对数,使得 1<=x<=a,1<=y<=b,且gcd(x,y) == 1。

    其中 (2,3) (3,2) 算两对!

    直接暴力显然复杂度太大,我们用莫比乌斯反演来解决。

    令N = max(a,b)

    然后定义 F(n) 表示满足条件的 gcd(x,y)==n的 (x,y) 对数

    在定义 G(n) 表示满足 n | gcd(x,y) 的(x,y)对数  即 gcd(x,y)%n==0 的x,y对数

    那么根据定义,有 G(n) = sigma(F(d)) (n|d,d<=N) 

    于是我们需要求的就是 F(1)

    怎么解决?

    首先根据G(n)的定义,可以很容易发现 G(n) = (a/n)*(b/n)这里是向下整除 (提示:把n当成最小的元)

    然后 我们只直接计算 F(1) 即可

    带入 G(n) 的公式 有  F(1) = sigma(U[i]*(a/i)*(b/i)) (1<=i<=N)

    至于U[]的值,可以提前用筛法在O(n)的时间内处理出来,这样总的时间复杂度就是 O(n),问题得到解决!!

    下面附上我自己求U[]的代码 (效率并不是严格上的O(n),不过一般情况下已经足够)

    #include <iostream>
    #include <cstring>
    using namespace std;
    
    #define REP(i,n) for(int i=0;i<(n);i++)
    #define FOR(i,j,k) for(int i=j;i<=(k);i++)
    #define maxn 1000009
    
    int u[maxn];
    bool vis[maxn];
    void Work(){        // 计算U[]
        memset(vis,0,sizeof(vis));
        memset(u,0,sizeof(u));
        for(int i=2;i<maxn;i++) if(!vis[i])
            for(int j=i;j<maxn;j+=i){
                if(u[j]==-1) continue;
                if((j/i)%i==0)  u[j]=-1;
                else u[j] ++;
                vis[j]=1;
            }
        FOR(i,1,maxn-1)
            if(u[i]==-1) u[i]=0;
            else if(u[i]%2) u[i]=-1;
            else u[i]=1;
    }
    
    int main(){
        Work();
    }

    四:进阶

    在ACM中,可以利用 莫比乌斯反演 来求解很多关于 Gcd 的问题 

    推荐几道基础题: SPOJ 7001 , ZOJ 3435, HDU 1695.

    想做更多的题的话,自己去HUST OJ搜索吧 :) 

    最后再说一下上面的证明方法都是个人YY的,感觉比《组合数学》上的证明简单些(数学太渣orz...那个证明我是没看太懂),写下来

    给初学莫比乌斯反演的童鞋当个资料(= =)。关于上面的证明我暂时没发现什么错误,如果发现错误,请在回复里面指出!另外

    形式二的证明应该可以由形式一直接得到,不过我没想出什么好办法,知道的神牛也请在评论中说一下!

    ( ̄. ̄) 完~~
  • 相关阅读:
    阿里数据中台
    智慧公安及移动执法平台(中云微迅)
    平安智慧城市的智慧政务
    数据分析师干货-业务逻辑篇
    数据分析师常用商业模型(一)
    智慧政务大数据云平台解决方案
    说几个flink好做spark却不好做的场景
    spark比flink好用的点
    干货好文!自底向上——知识图谱构建技术初探
    用户分析模型
  • 原文地址:https://www.cnblogs.com/chenyang920/p/4811995.html
Copyright © 2020-2023  润新知