【HDU6134】Battlestation Operational (莫比乌斯反演)

多校第八场没想到碰到了一直想学却没有学的莫比乌斯反演，比赛时不出意外的打出了gg，便下定决心要好好会一会高深的莫比乌斯反演(在菜鸡的眼里是这样),下面是我的一些学习资料：

http://www.cnblogs.com/chenyang920/p/4811995.html

http://blog.csdn.net/acdreamers/article/details/8542292

https://wenku.baidu.com/view/fbec9c63ba1aa8114431d9ac.html

 

看完这些之后不得不佩服这些数学家的伟大，这个算法对于解决组合问题，包含gcd的条件的问题实在太厉害了。

 

然后就要开始讲这道题目了，这场的多校赛题目都是巨长，我就大概说一下这道题的题意：

题目给出一个n*n的矩阵对于矩阵中每一个点(x, y) 如果x >= y && gcd(x, y) = 1, 那么这个点的值为 x / y 向上取整， 否则就为0，要求这个矩阵的各个点的值的总和。即如公式所见：

 

分析：

首先设

因为gcd(i, j)等于1才有效，因此i / j 一定无法整除（除了j = 1），且结果就为 i / j + 1,因此有如下转化，其中φ(i)为欧拉函数

然后就是重要发现：

 

然后运用莫比乌斯反演

可得：

 

继续转化：

其中τ(n)表示n的约数个数，是积性函数也可以线性求出

其实后来发现：

F(n) = F（n - 1） + d(n - 1) + 1， 其中d(n)表示n的约数的个数

这样可以递推出每一项的F(n),

然后利用莫比乌斯反演公式就可以求出每一项f(n),注意这里求每一项的f(n)可以枚举j, 然后j的倍数的答案可以很容易求出，因为调和级数这样复杂度为nlog(n), f(n)的前缀和就是答案

 

#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <algorithm>
#include <cstring>
#include <map>
#include <set>
#include <vector>
#include <queue>
#include <cmath>
#include <bitset>

using namespace std;

const int M = 1e9 + 7;
const int MAXN = 1000005;
int g[MAXN], pre[MAXN];

bool check[MAXN+10];
int prime[MAXN+10];
int mu[MAXN+10];
void Moblus()
{
    memset(check,false,sizeof(check));
    mu[1] = 1;
    int tot = 0;
    for(int i = 2; i <= 1000000; i++)
    {
        if( !check[i] )
        {
            prime[tot++] = i;
            mu[i] = -1;
        }
        for(int j = 0; j < tot; j++)
        {
            if(i * prime[j] > MAXN) break;
            check[i * prime[j]] = true;
            if( i % prime[j] == 0)
            {
                mu[i * prime[j]] = 0;
                break;
            }
            else
            {
                mu[i * prime[j]] = -mu[i];
            }
        }
    }
}

int factor[100][2];
int getFactors(int x)
{
    int fatCnt = 0;
    int tmp = x;
    for(int i = 0; prime[i] <= tmp / prime[i];i++)
    {
        factor[fatCnt][1] = 0;
        if(tmp%prime[i] == 0)
        {
            factor[fatCnt][0] = prime[i];
            while(tmp%prime[i] == 0)
            {
                factor[fatCnt][1]++;
                tmp /= prime[i];
            }
            fatCnt++;
        }
    }
    if(tmp != 1)
    {
        factor[fatCnt][0] = tmp;
        factor[fatCnt++][1] = 1;
    }
    return fatCnt;
}

int main()
{
    Moblus();
    g[1] = 1;
    int cnt = 1;
    for(int i = 2; i <= 1000000; ++i)
    {
        g[i] = g[i - 1] + cnt + 1;
        int m = getFactors(i);
        cnt = 1;
        for(int j = 0; j < m; ++j)
             cnt *= (factor[j][1] + 1);
    }

    for(int j = 1; j <= 1000000; ++j)
        for(int k = j, t = 1; k <= 1000000; k += j, t++)
            pre[k] += mu[t] * g[j];
    for(int i = 2; i <= 1000000; ++i)
        pre[i] = (pre[i] + pre[i - 1]) % M;
    int n;
    while(scanf("%d", &n) != EOF)
    {
        printf("%d\n", pre[n]);
    }
    return 0;
}

 

 

 

 

 

然后这里还附一份我觉得更简洁的代码：

 

#include<cstdio>
#include<cstring>
#include<vector>
using namespace std;
typedef long long int ll;

const ll mod = 1e9 + 7;
const int A = 1e6 + 10;
int mu[A],pri[A],tot;
ll d[A],cnt[A],sum[A];
bool vis[A];

void init(){
    tot = 0;mu[1] = d[1] = 1;
    for(int i=2 ;i<A ;i++){
        if(!vis[i]){mu[i] = -1;pri[++tot] = i;d[i] = 2;cnt[i] = 1;}
        for(int j=1 ;j<=tot && pri[j]*i<A ;j++){
            vis[i*pri[j]] = 1;
            if(i%pri[j] == 0){
                d[i*pri[j]] = d[i]/(cnt[i]+1)*(cnt[i]+2);
                cnt[i*pri[j]] = cnt[i] + 1;
                mu[i*pri[j]] = 0;
                break;
            }
            d[i*pri[j]] = d[i]<<1;
            cnt[i*pri[j]] = 1;
            mu[i*pri[j]] = -mu[i];
        }
    }
    sum[1] = 1;
    for(int i=2 ;i<A ;i++){
        sum[i] = (sum[i-1] + d[i-1] + 1)%mod;
    }
    for(int i=1 ;i<A ;i++){
        sum[i] = (sum[i] + sum[i-1])%mod;
        mu[i] = (mu[i] + mu[i-1])%mod;
    }
}


int main(){
    init();
    int n;
    while(~scanf("%d",&n)){
        ll ans = 0;
        int last;
        for(int i=1 ;i<=n ;i=last + 1){
            last = n/(n/i);
            ans = (ans + (mu[last]-mu[i-1])%mod*(sum[n/i])%mod)%mod;
        }
        ans = (ans%mod + mod)%mod;
        printf("%I64d\n",ans);
    }
    return 0;
}