HDU 4746 Mophues(有趣的前缀和/莫比乌斯反演)

题目链接:
HDU 4746 Mophues
题意:
gcd(a,b)=ddkkCkCluckyn,m,pa[1,n],b[1,m]使gcd(a,b)plucky(a,b)
分析:
d=gcd(i,j) ,我们枚举 d ,对于每个d计算 i=ni=1j=mj=1(gcd(i,j)==d) ,也就是 gcd(i,j)=d 的数字对数,这个就是莫比乌斯反演搞搞了。
但是因为 d 需要满足d质因子分解后质因子数字个数 p ,如果这样子直接枚举并且采用前缀和加速的话,因为 d 不一定时连续的,所以就不行了。我们可以考虑将前缀和记成二维的形式,sum[i][j]表示 gcd 为i时质因子分解数字个数 j μ 贡献之和。因为最大公约数是 i ,我们在累加的时候需要加的是 μ(ji),其中 j i的倍数。还有一点,考虑 5105 范围质因子分解的数字最大约为 20 ,因为 220 = 10242 ,也就是 106 数据大概有 20 个。那么当 20p 时, [1,n],[1,m] 中的每对数字都满足 gcd(i,j) 质因子分解数字个数 <p <script type="math/tex" id="MathJax-Element-871"> nm.


#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
#include <bitset>
#include <string>
using namespace std;
typedef long long ll;
const int MAX_N = 500010;
const int SZ = 25;

int prime_cnt, prime[MAX_N], factor_num[MAX_N];
ll mu[MAX_N], sum[MAX_N][SZ];
bitset<MAX_N> bs;

void GetMu()
{
    prime_cnt = 0;
    bs.set();
    mu[1] = 1;
    for(int i = 2; i < MAX_N; ++i) {
        if(bs[i]) {
            prime[prime_cnt++] = i;
            mu[i] = -1;
        }
        for(int j = 0; j < prime_cnt && i * prime[j] < MAX_N; ++j) {
            bs[i * prime[j]] = 0;
            if(i % prime[j]) {
                mu[i * prime[j]] = -mu[i];
            }else {
                mu[i * prime[j]] = 0;
                break;
            }
        }
    }
    int MAX_SZ = 0;
    for(int i = 0; i < prime_cnt; ++i ) {
        for(int j = prime[i]; j < MAX_N; j += prime[i]) {
            int x = j, tmp = 0;
            while(x % prime[i] == 0) {
                x /= prime[i];
                tmp++;
            }
            factor_num[j] += tmp;
            MAX_SZ = max(MAX_SZ, factor_num[j]);
        }
    }
    //MAX_SZ = 18
    for(int i  = 1; i < MAX_N; ++i) {
        for(int j = i; j < MAX_N; j += i){
            sum[j][factor_num[i]] += mu[j / i];
        }
    }
    for(int i = 1; i < MAX_N; ++i) {
        sum[i][0] +=sum[i - 1][0];
        for(int j = 1; j < SZ; ++j) {
            sum[i][j] += (sum[i - 1][j] + sum[i][j- 1] - sum[i - 1][j - 1]);    
        }   
    }
}

inline ll solve(int n, int m, int p)
{
    ll ans = 0;
    int top = min(n, m), last;
    for(int i = 1; i <= top; i = last + 1){
        last = min(n / (n / i), m / (m / i));
        ans += (sum[last][p] - sum[i - 1][p]) * (n / i) * (m / i);
    }
    return ans;
}

int main()
{
    GetMu();
    int T, n, m, p;
    scanf("%d", &T);
    while(T--) {
        scanf("%d%d%d", &n, &m, &p);
        if(p >= 20) printf("%lld\n", (ll)n * m);
        else printf("%lld\n", solve(n, m, p));
    }
    return 0;
}
  • 0
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值