• 杜教筛入门


    以下主要的话都用无序列表表示。

    诶,是不是应该先讲背景

    有什么好讲的?

    问一个积性函数的前缀和,项数到1e10。

    前置知识

    线性筛积性函数

    正文

    钦定你已经可以再(O(sqrt{n}))的复杂度内求出:

    [sum_{i=1}^n{lfloor frac{n}{i} floor} ]

    [sum_{i=1}^n{i imes lfloor frac{n}{i} floor} ]

    对于第一个,先枚举小于(sqrt n)的i,得出这段的值;又因为(i)在一段区间内(lfloor frac{n}{i} floor)都是(i_0leq sqrt n),因此可以求前式和后式。

    • 形式化地讲,若(sum_{i=1}^n{f(i) imes g(lfloor frac{n}{i} floor)})中,若能方便的求得(f(i))的前缀和以及(g(x)),就能方便的求得原式。

    然后对于积性函数(f(i)),我们想求(S_n=sum_{i=1}^n{f(i)})

    那么,我们找一个积性函数(g),令(T_n=sum_{dmid n}g(d) imes f(frac{n}{d}))(就是狄利克雷卷积)。则有:

    [egin{aligned} sum_{i=1}^n{T_i} & =sum_{i=1}^n sum_{dmid i}{g(d) imes f(frac{i}{d})} \ & =sum_{d=1}^n{g(d) imes sum_{dmid i, ileq n}{f(frac{i}{d})}} \ & =sum_{d=1}^n{g(d) imes S_{lfloor frac{n}{d} floor}} end{aligned} ]

    然后,钦定(g(1)=1),那么就有

    [egin{aligned} S_n & =sum_{d=1}^n{g(d) imes S_{lfloor frac{n}{d} floor}}-sum_{d=2}^n{g(d) imes S_{lfloor frac{n}{d} floor}} \ & =sum_{i=1}^n{T_i}-sum_{d=2}^n{g(d) imes S_{lfloor frac{n}{d} floor}} end{aligned} ]

    如果不管怎么求(S_{lfloor frac{n}{d} floor})的话,发现满足以下两个条件就可以求(S_n)了。

    • 可以求g的前缀和
    • 可以求T的前缀和

    考虑怎么求上式后面的(S_{lfloor frac{n}{d} floor}),由于f是积性函数,必定可以用线性筛筛出前n项,于是可以用线性筛筛了f求出(S_{1..sqrt n}),至于大于(sqrt{n})的下标,可以记搜:因为(lfloor frac{lfloor frac{n}{d} floor}{e} floor=lfloor frac{n}{d imes e} floor),于是可以记录n/d的d,当(d>sqrt{n})时直接返回结果就行。

    • 于是,求前缀和成功转化成用人类智慧求一个g。

    模板题

    对于(phi),g为1(常值函数),T为i(自然数序列),用到的结论是(x=sum_{dmid x}{phi_d})

    证明有这样几个方向:

    1. 证明(forall imid n,1leq j<i,(i,j)=1),都有(i_1 imes j_1 eq i_2 imes j_2)。显然是要被叉翻的(其实是我在用3时不知道干嘛了)。
    2. 直接通过(phi)的计算公式和积性通过一波推理得到一些式子,化简得到n。(巨佬做法)
    3. 证明(forall imid n,1leq j<i,(i,j)=1),都有(n/i_1 imes j_1 eq n/i_2 imes j_2)。显然,因为互质,所以两个都是最简分数,于是当i,j不同时,分数不可能相等,于是证明了任意一个(xin [1,n])都只能由一个(i,j)对转移而来,即一一对应,证毕。

    以上2的证明:

    [egin{aligned} sum_{dmid n}{phi(d)} & =sum{phi(prod_{i=1}^m{p_i^j})} \ & =sum{prod_{i=1}^m{phi(p_i^j)}} \ & =prod_{i=1}^m{sum_{j=0}^{a_i}{phi(p_i^j)}} \ & =prod_{i=1}^m{(sum_{j=1}^{a_i}{(p_i^j-p_i^{j-1})}+1)} \ & =prod_{i=1}^m{p_i^{a_i}}=n end{aligned} ]

    然后写出来就可以了。

    以下是一份跑的非常慢的代码模板(洛谷模板题的关键代码)。

    const int N = 3000005, nn = 3000000;
    struct getSum{
        ll presum[N], aftersum[N];
        bool calced[N];
        inline ll sum(int n,int d,ll sT(int),ll sg(int),ll g(int)){
            if (n / d <= nn) return presum[n / d];
            if (calced[d]) return aftersum[d];
            int nn = n / d;
            ll ans = sT(nn);
            ans -= (sg(nn) - sg(nn / 2)) * presum[1];
            REP(i, 2, floor(sqrt(nn))){
                ans -= g(i) * sum(n, d * i, T, b, g);
                ans -= (sg(nn / i) - sg(Max(nn / (i + 1), i))) * presum[i];
            }
            calced[d] = 1;
            aftersum[d] = ans;
            return ans;
        }
    };
    int T, n;
    int b[N];
    ll phi[N], miu[N];
    int temp[N / 10], top;
    int main(){
        read(T);
        phi[1] = miu[1] = 1;
        REP(i, 2, nn){
            if (!b[i]){
                temp[++top] = i;
                phi[i] = i - 1;
                miu[i] = -1;
            }
            for (int j = 1; j <= top && i * temp[j] <= nn; ++j){
                b[i * temp[j]] = 1;
                if (i % temp[j]){
                    phi[i * temp[j]] = phi[i] * (temp[j] - 1);
                    miu[i * temp[j]] = -miu[i];
                }
                else{
                    phi[i * temp[j]] = phi[i] * temp[j];
                    miu[i * temp[j]] = 0;
                    break;
                }
            }
        }
        REP(i, 2, nn){
            phi[i] = phi[i] + phi[i - 1];
            miu[i] = miu[i] + miu[i - 1];
        }
        memcpy(van.presum, phi, 8 * (nn + 1));
        memcpy(deep.presum, miu, 8 * (nn + 1));
        while (T--){
            read(n);
            mem(van.calced);
            printf("%lld ", van.sum(n, 1, [](int n) { return 1LL * n * (n + 1) / 2; }, [](int n) { return 1LL * n; }, [](int n) { return 1LL; }));
            mem(deep.calced);
            printf("%lld
    ", deep.sum(n, 1, [](int n) { return 1LL; }, [](int n) { return 1LL * n; }, [](int n) { return 1LL; }));
        }
        return 0;
    }
    
  • 相关阅读:
    RocketMQ系列(一)基本概念
    怎样实现登录?| Cookie or JWT
    Hotspot GC研发工程师也许漏掉了一块逻辑
    初级Java工程师也能轻松进行JVM调优了
    自动化不知如何参数化(二)?xlrd来帮你解决
    自动化不知如何参数化(一)?xlrd来帮你解决
    SpringCloud系列之API网关(Gateway)服务Zuul
    SpringCloud系列之客户端负载均衡Netflix Ribbon
    SpringCloud系列之使用Feign进行服务调用
    Spring Security系列之极速入门与实践教程
  • 原文地址:https://www.cnblogs.com/pupuvovovovo/p/11679401.html
Copyright © 2020-2023  润新知