首页 > 编程语言 > 详细

[SDOI2014] 数表 - 莫比乌斯反演,树状数组,整除分块

时间:2020-06-16 23:52:15      阅读:116      评论:0      收藏:0      [点我收藏+]

Description

有一张 \(n\times m\) 的数表,其第 \(i\) 行第 \(j\) 列(\(1\le i\le n\)\(1\le j\le m\))的数值为能同时整除 \(i\)\(j\) 的所有自然数之和。给定 \(a\),计算数表中不大于 \(a\) 的数之和。

Solution

先按套路进行一些推导

技术分享图片

我们用 BIT 来维护 \(h(x)\),将所有询问按照 \(a\) 升序排序,当 \(a\) 的容许范围扩大时,就将对应的所有 \(x\) 满足 \(f(x)=a\)\(h\) 的贡献进行修改,具体地,对于 \(x\),它会对 \(ix\) 产生 \(a \mu(i)\) 的贡献,这部分单点修改即可。在进行整除分块计算的过程中,在 BIT 上区间求和即可。

#include <bits/stdc++.h>
using namespace std;

#define int long long
const int N = 1000005;
const int MAXN = 1000005;
const int mod = 1ll<<31;

namespace bit {
    const int N = 1000000;
    int ar[N]; // index: 1 ~ N
    int lowbit(int t) { return t & (-t); }
    void add(int i, int v) {
        for (; i < N; ar[i] += v, i += lowbit(i));
    }
    int sum(int i) {
        int s = 0;
        for (; i > 0; s += ar[i], i -= lowbit(i));
        return s;
    }
    int query(int i,int j) {
        return sum(j) - sum(i-1);
    }
}

bool isNotPrime[MAXN + 1];
int mu[MAXN + 1], phi[MAXN + 1], primes[MAXN + 1], cnt, f[N], idx[N];

map<int,int> mp;
vector<int> vec[N];

int n=1e5;

inline void euler() {
    isNotPrime[0] = isNotPrime[1] = true;
    mu[1] = 1;
    for (int i = 2; i <= MAXN; i++) {
        if (!isNotPrime[i]) {
            primes[++cnt] = i;
            mu[i] = -1;
        }
        for (int j = 1; j <= cnt; j++) {
            int t = i * primes[j];
            if (t > MAXN) break;
            isNotPrime[t] = true;
            if (i % primes[j] == 0) {
                mu[t] = 0;
                break;
            } else {
                mu[t] = -mu[i];
            }
        }
    }
    for(int i=1;i<=n;i++) {
        int j=1;
        for(;j*j<i;j++) if(i%j==0) {
            f[i]+=j;
            f[i]+=i/j;
        }
        if(j*j==i) {
            f[i]+=j;
        }
        mp[f[i]]++;
    }
    int ind=0;
    for(auto &i:mp) {
        i.second=++ind;
        idx[ind]=i.first;
    }
    for(int i=1;i<=n;i++) {
        vec[mp[f[i]]].push_back(i);
    }
}

struct qry {
    int l,r,a,id;
    bool operator < (const qry & b) {
        return a<b.a;
    }
} qr[N];

void refresh(int x,int a) {
    for(int i=1;i*x<=n;i++) {
        bit::add(i*x,a*mu[i]);
    }
}

void refresh(int a) {
    int t=mp[a];
    for(int i:vec[t]) {
        refresh(i,a);
    }
}

int solve(signed n,signed m) {
    if(n==0 || m==0) return 0;
    int ans=0;
    signed l=1,r=0;
    if(n>m) swap(n,m);
    while(l<=n) {
        r=min(n/(n/l),m/(m/l));
        ans+=bit::query(l,r)*(n/l)*(m/l);
        l=r+1;
    }
    return ans;
}

int res[N];

signed main() {
    euler();
    int t,n,m;
    ios::sync_with_stdio(false);
    cin>>t;
    int cc=0;
    while(t--) {
        ++cc;
        int x,y,z;
        cin>>x>>y>>z;
        qr[cc]={x,y,z,cc};
    }
    sort(qr+1,qr+cc+1);
    int pos=0;
    for(int i=1;i<=cc;i++) {
        auto it=mp.upper_bound(qr[i].a);
        --it;
        int tmp=0;
        if(it!=(--mp.begin())) tmp=it->second;
        while(pos<tmp) {
            ++pos;
            refresh(idx[pos]);
        }
        res[qr[i].id]=solve(qr[i].l,qr[i].r);
    }
    for(int i=1;i<=cc;i++) cout<<(res[i]%mod+mod)%mod<<endl;
}

[SDOI2014] 数表 - 莫比乌斯反演,树状数组,整除分块

原文:https://www.cnblogs.com/mollnn/p/13149906.html

(0)
(0)
   
举报
评论 一句话评论(0
关于我们 - 联系我们 - 留言反馈 - 联系我们:wmxa8@hotmail.com
© 2014 bubuko.com 版权所有
打开技术之扣,分享程序人生!