zoukankan      html  css  js  c++  java
  • bzoj4176

    莫比乌斯反演

    根据约数和个数公式

    $ans = sum_{i=1}^{n}sum_{j=1}^{n}sum_{x|i}sum_{y|j}{[gcd(i, j)==1]}$

    交换枚举顺序

    $ans = sum_{x=1}^{n}sum_{y=1}^{n}{[frac{n}{x}][frac{n}{y}]*[gcd(x, y)==1]}$

    $=sum_{x=1}^{n}sum_{y=1}^{n}{[frac{n}{x}][frac{n}{y}]sum_{d|x,y}{mu(d)}}$

    再交换枚举顺序

    $=sum_{d=1}^{n}{mu(d)sum_{i=1}^{frac{n}{d}}{frac{n}{di}}sum_{j=1}^{frac{n}{d}}{frac{n}{dj}}}$

    设$f(n)=sum_{i=1}^{n}{frac{n}{i}}$

    那么$=sum_{d=1}^{n}{mu(d)f(frac{n}{d})^{2}}$

    这就可以分块求了,$mu$用杜教筛求,$f$用二次分块

    反演就是不断交换求和顺序或者改变枚举变量

    #include <cstdio>
    #include <cstring>
    #include <algorithm>
    #include <map>
    using namespace std;
    typedef long long ll;
    const int N = 1e7 + 5, P = 1e9 + 7;
    int n;
    ll ans;
    int p[N], mark[N], mu[N];
    map<int, ll> sum_1;
    void ini() {
        mu[1] = 1;
        for(int i = 2; i < N; ++i) {
            if(!mark[i]) {
                p[++p[0]] = i;
                mu[i] = -1;
            }
            for(int j = 1; j <= p[0] && i * p[j] < N; ++j) {
                mark[i * p[j]] = 1;
                if(i % p[j] == 0) {
                    mu[i * p[j]] = 0;
                    break;
                }
                mu[i * p[j]] = -mu[i];
            }
        }
        for(int i = 1; i < N; ++i) {
            mu[i] += mu[i - 1];
        }
    }
    ll dj_m(int n) {
        if(n < N) {
            return mu[n];
        }
        if(sum_1.find(n) != sum_1.end()) {
            return sum_1[n];
        }
        ll ret = 1;
        for(int i = 2, j = 0; i <= n; i = j + 1) {
            j = n / (n / i);
            ret = (ret - (ll)(j - i + 1) * dj_m(n / i) % P + P) % P;
        }
        return sum_1[n] = ret;
    }
    ll calc(int n) {
        ll ret = 0;
        for(int i = 1, j = 0; i <= n; i = j + 1) {
            j = n / (n / i);
            ret = (ret + (ll)(n / i) * (j - i + 1) % P) % P;
        }
        return ret;
    }
    int main() {
        ini();
        scanf("%d", &n);
        for(int i = 1, j = 0; i <= n; i = j + 1) {
            j = n / (n / i);
            ll a = ((dj_m(j) - dj_m(i - 1)) % P + P) % P, b = calc(n / i);
            ans = (ans + a * b % P * b % P) % P;
        }
        printf("%lld
    ", ans);
        return 0;
    }
    View Code
  • 相关阅读:
    卧槽!缓存的问题太多了(雪崩、击穿、穿透…)一个个解决!
    Java 命名规范(非常全面,可以收藏)
    一次接口超时排查,花费了我两个星期。。
    LiveGBS和海康威视
    SQLite文件存储和读取
    Vue页面刷新原理:Cesium刷新机制
    MBtiles格式数据
    gitee:403错误
    uniapp是什么?
    HBuilderx怎么运行代码
  • 原文地址:https://www.cnblogs.com/19992147orz/p/8465159.html
Copyright © 2011-2022 走看看