首页 > 其他分享 >【学习笔记】狄利克雷卷积与高级筛法

【学习笔记】狄利克雷卷积与高级筛法

时间:2023-07-03 18:12:37浏览次数:44  
标签:lfloor right 筛法 狄利克 卷积 dfrac sum mathrm left

狄利克雷卷积

概念

对于数论函数 \(f,g\),定义其狄利克雷卷积 \(h=f*g\),满足:

\[h(n)=(f*g)(n)=\sum_{d\mid n} f(d)g\left(\dfrac{n}{d}\right) \]

运算律:

  • 满足交换律,显然具有对称性。

  • 满足结合律,等价于三个 \(d_i\) 贡献到 \(n\)。

  • 满足加法的分配率。

常见数论函数:

  • \(\mathrm{I}\),常函数,值恒为 \(1\)。

  • \(\mathrm{Id}_k\),幂函数,值为 \(n^k\)。

  • \(\epsilon\),狄利克雷卷积的单位元,除 \(\epsilon(1)=1\) 外其余均为 \(0\)。

  • \(\mu\),莫比乌斯函数,狄利克雷卷积的逆元。

  • \(\varphi\),欧拉函数。

  • \(\sigma_k\),除数函数,定义 \(\sigma_k(n)=\sum_{d\mid n}d^k\)。

常见数论函数卷积变换

  • \(\mu *\mathrm{I}=\epsilon\)

  • \(\varphi *\mathrm{I}=\mathrm{Id}\)

  • \(\mathrm{Id}_k *\mathrm{I}=\sigma_k\)

  • \(\mu* \mathrm{Id}=\varphi\)

  • \((\mu \cdot \mathrm{Id}_k)*\mathrm{Id}_k=\epsilon\)

  • \((\varphi\cdot \mathrm{Id}_k)*\mathrm{Id}_k=\mathrm{Id}_{k+1}\)

杜教筛

原理

求数论函数前缀和。

设 \(h=f*g,S=\sum f\),容易推得:

\[\begin{aligned} \sum_{i=1}^n h(i)&=\sum_{i=1}^n \sum_{j\mid i} f(j)g\left(\dfrac{i}{j}\right)\\ &=\sum_{j=1}^n g(j)\sum_{i=1}^{\left\lfloor n/i\right\rfloor} f(i)\\ &=\sum_{i=1}^n g(i)S_f\left(\left\lfloor \dfrac{n}{i}\right\rfloor\right) \end{aligned} \]

整理一下可以写成:

\[S_f(n)=\sum_{i=1}^n h(i)-\sum_{i=2}^{n} g(i)S_f\left(\left\lfloor \dfrac{n}{i}\right\rfloor\right) \]

如果构造出可以快速计算前缀和的函数 \(g,h\),就可以递归求解了。

时间复杂度

首先需要证明的是:

\[\left\lfloor\dfrac{\left\lfloor\tfrac{n}{a}\right\rfloor}{b}\right\rfloor=\left\lfloor\dfrac{n}{ab}\right\rfloor \]

设 \(\mathrm{RHS}=k\),则有:

\[k\times ab\le n<(k+1)\times ab \]

于是:

\[k\times b\le \dfrac{n}{a}<(k+1)\times b \]

从而:

\[k\times b\le \left\lfloor\dfrac{n}{a}\right\rfloor<(k+1)\times b \]

\[k\le \dfrac{\left\lfloor\tfrac{n}{a}\right\rfloor}{b}<k+1 \]

即:

\[\left\lfloor\dfrac{\left\lfloor\tfrac{n}{a}\right\rfloor}{b}\right\rfloor=\left\lfloor\dfrac{n}{ab}\right\rfloor \]


这说明当我们整除分块向下递归时,递归到 \(n'=\left\lfloor \dfrac{n}{d}\right\rfloor\),枚举 \(d'\) 得到 \(\left\lfloor \dfrac{n'}{d'}\right\rfloor=\left\lfloor \dfrac{\left\lfloor \tfrac{n}{d}\right\rfloor}{d'}\right\rfloor\),一定可以表示为 \(\left\lfloor \dfrac{n}{k}\right\rfloor\) 的形式,而这一定在第一层递归时处理到(不考虑先后顺序),于是在记忆化的加持下,计算复杂度只需要考虑第一层递归即可。

设当前块长 \(i\),对复杂度的贡献就是 \(O\left(\sqrt{\dfrac{n}{i}}\right)\)。

当 \(i> \sqrt{n}\) 时,第二层最大 \(O(n^{1/4})\),所以总复杂度 \(O(n^{3/4})\)。

当 \(i\le \sqrt{n}\) 时,就写成:

\[\begin{aligned} T(n)&=O\left(\sum_{i=1}^{\sqrt{n}} \sqrt{\dfrac{n}{i}}\right)\\ &=O\left(\sqrt{n}\int_{1}^n \dfrac{1}{\sqrt{x}} \mathrm{d}x\right)\\ &=O(n^{3/4}) \end{aligned} \]

但较小的前缀和是可以预处理出的,假设预处理出了 \([1,n^k]\) 的前缀和,那么也就是只有 \(\dfrac{n}{i}>n^k\) 的部分才需要递归,也就变成:

\[\begin{aligned} T(n)&=O\left(n^k+\sum_{i=1}^{n^{1-k}} \sqrt{\dfrac{n}{i}}\right)\\ &=O\left(n^k+\sqrt{n}\int_{1}^{n^{1-k}} \dfrac{1}{\sqrt{x}}\mathrm{d}x\right)\\ &=O(n^k+n^{1-k/2}) \end{aligned} \]

取 \(k=\dfrac{2}{3}\) 得到理论最优复杂度 \(O(n^{2/3})\)。

整除分块时(无论嵌套几层)使用杜教筛求数论函数前缀和,由于所有询问的区间右端点为 \(\left\lfloor \dfrac{n}{\left\lfloor \tfrac{n}{d}\right\rfloor}\right\rfloor\),可以表示为上述记忆化的形式,所以并不会增加复杂度。

点击查看代码

namespace Phi{
    int pr[maxn],mn[maxn],mnpw[maxn],phi[maxn];
    bool vis[maxn];
    ull s1[maxn],s2[maxn];
    bool mark[maxn];
    inline void linear_sieve(){
        mn[1]=1,phi[1]=1;
        for(int i=2;i<=lim;++i){
            if(!vis[i]) pr[++pr[0]]=i,mn[i]=i,mnpw[i]=i,phi[i]=i-1;
            for(int j=1;j<=pr[0]&&i*pr[j]<=lim;++j){
                vis[i*pr[j]]=1,mn[i*pr[j]]=pr[j];
                if(mn[i]==pr[j]) mnpw[i*pr[j]]=mnpw[i]*pr[j];
                else mnpw[i*pr[j]]=pr[j];
                if(i*pr[j]==mnpw[i*pr[j]]) phi[i*pr[j]]=phi[i]*pr[j];
                else phi[i*pr[j]]=phi[i*pr[j]/mnpw[i*pr[j]]]*phi[mnpw[i*pr[j]]];
                if(i%pr[j]==0) break;
            }
        }
        for(int i=1;i<=lim;++i) s1[i]=s1[i-1]+phi[i];
    }
    ull get_sum(int N){
        if(N<=lim) return s1[N];
        if(mark[n/N]) return s2[n/N];
        mark[n/N]=1;
        ull sum=0;
        for(unsigned int l=2,r;l<=N;l=r+1){
            r=N/(N/l);
            sum+=1ull*(r-l+1)*get_sum(N/l);
        }
        if(N&1) sum=1ull*(N+1)/2*N-sum;
        else sum=1ull*N/2*(N+1)-sum;
        return s2[n/N]=sum;
    }
    inline ull solve(){
        memset(mark,0,sizeof(mark));
        return get_sum(n);
    }
}

namespace Mu{
    int pr[maxn],mu[maxn];
    bool vis[maxn];
    int s1[maxn],s2[maxn];
    bool mark[maxn];
    inline void linear_sieve(){
        mu[1]=1;
        for(int i=2;i<=lim;++i){
            if(!vis[i]) pr[++pr[0]]=i,mu[i]=-1;
            for(int j=1;j<=pr[0]&&i*pr[j]<=lim;++j){
                vis[i*pr[j]]=1,mu[i*pr[j]]=-mu[i];
                if(i%pr[j]==0){
                    mu[i*pr[j]]=0;
                    break;
                }
            }
        }
        for(int i=1;i<=lim;++i) s1[i]=s1[i-1]+mu[i];
    }
    int get_sum(int N){
        if(N<=lim) return s1[N];
        if(mark[n/N]) return s2[n/N];
        mark[n/N]=1;
        ull sum=0;
        for(unsigned int l=2,r;l<=N;l=r+1){
            r=N/(N/l);
            sum+=(r-l+1)*get_sum(N/l);
        }
        sum=1-sum;
        return s2[n/N]=sum;
    }
    inline int solve(){
        memset(mark,0,sizeof(mark));
        return get_sum(n);
    }
}

常用构造

最常用的是 \(\mu *\mathrm{I}=\epsilon\) 以及 \(\varphi *\mathrm{I}=\mathrm{Id}\)。

在此基础上有:\((\mu\cdot \mathrm{Id}_k)*\mathrm{Id}_k=\epsilon\),\((\varphi\cdot \mathrm{Id}_k) * \mathrm{Id}_{k}=\mathrm{Id}_{k+1}\)。

另外根据 \(\sigma_k\) 的定义有 \(\mathrm{Id}_k*\mathrm{I}=\sigma_k\)。

例题

Luogu-P5218 无聊的水题 II

容易发现是要求选出一个集合满足 \(\gcd=1\)。

设 \(F(n)\) 为值域 \([1,n]\) 的方案数,那么 \(\gcd=k\) 就是 \(F\left(\left\lfloor \dfrac{n}{k}\right\rfloor\right)\),可以单步容斥。

\[F(n)=2^n-1-\sum_{k=2}^n F\left(\left\lfloor \dfrac{n}{k}\right\rfloor\right) \]

式子长得像杜教筛,复杂度 \(O(n^{3/4})\) 不太好过。

考虑设 \(f(d)\) 为 \(\gcd=d\) 的方案数,\(g(d)\) 为 \(d \mid \gcd\) 的方案数,则 \(g(d)=2^{\left\lfloor n/d\right\rfloor}-1\),且有:

\[g(d)=\sum_{d\mid n} f(n) \]

根据莫比乌斯反演有:

\[f(d)=\sum_{d\mid n} \mu\left(\dfrac{n}{d}\right) g(n) \]

所以:

\[f(1)=\sum_{i=1}^n \mu(i)g(i)=\sum_{i=1}^n \mu(i)(2^{\left\lfloor n/i\right\rfloor}-1) \]

杜教筛预处理+整除分块即可,使用光速幂复杂度 \(O(n^{2/3})\)。

LOJ-6229 这是一道简单的数学题

把 \(j\) 的上下界改成 \([1,n]\),处理一下 \(i=j\) 的情况即可。

简单反演得到:

\[\sum_{T=1}^n\sum_{d\mid T} \mu(d)d^2 S_1\left(\left\lfloor \dfrac{n}{d}\right\rfloor\right)^2 \]

\(S_1(x)\) 是 \(1\) 次自然数幂和。

于是要求 \(\sum_{d\mid n} \mu(d)d^2\) 的前缀和,不难发现这个函数是 \((\mu \cdot \mathrm{Id}_2)* \mathrm{I}\),根据狄利克雷卷积的交换律以及结合律,可以得到 \((\mu \cdot \mathrm{Id}_2)*\mathrm{I}*\mathrm{Id}_2=(\mu \cdot \mathrm{Id}_2)*\mathrm{Id}_2*\mathrm{I}=\epsilon*\mathrm{I}=\mathrm{I}\)。

于是有:

\[S_f(n)=n-\sum_{i=2}^{n} i^2S_f\left(\left\lfloor \dfrac{n}{i}\right\rfloor\right) \]

杜教筛处理即可。

LOJ-6491 XXOI2018 简单的最大公约数

和上面的题差不多,先设 \(f(m)\) 表示值域 \([1,m]\),\(n\) 个数 \(\gcd=1\) 的方案数,答案是 \(\sum_{k=1}^m k\times f\left(\left\lfloor \dfrac{m}{k}\right\rfloor\right)\),最终求解用整除分块。

\[f(m)=m^n-\sum_{i=2}^m f\left(\left\lfloor \dfrac{m}{i}\right\rfloor\right) \]

直接算还是 \(O(n^{3/4})\)。

希望预处理一些 \(f\)。

依旧是莫比乌斯反演,把 \(f(m,k)\) 扩充到值域 \([1,m]\),\(\gcd=k\) 的方案数,\(g(m,k)\) 定义为值域 \([1,m]\),\(k\mid \gcd\) 的方案数,于是反演得到:

\[g(m,k)=\left\lfloor \dfrac{m}{k}\right\rfloor^{n} \]

\[f(m,k)=\sum_{k\mid m'} \mu\left(\dfrac{m'}{k}\right) \left\lfloor \dfrac{m}{m'}\right\rfloor^{n} \]

只关心 \(k=1\),所以有:

\[f(m)=\sum_{k=1}^m \mu(k) \left\lfloor \dfrac{m}{k}\right\rfloor^{n} \]

这个直接求一行比较困难,考虑差分:

\[\begin{aligned} \Delta f(m)&=f(m)-f(m-1)\\ &=\sum_{k=1}^m \mu(k) \left(\left\lfloor \dfrac{m}{k}\right\rfloor^{n}-\left\lfloor \dfrac{m-1}{k}\right\rfloor^{n}\right) \end{aligned} \]

分析一下什么时候 \(\left\lfloor \dfrac{m}{k}\right\rfloor\neq \left\lfloor \dfrac{m-1}{k}\right\rfloor^{n}\),显然一定是 \(k\mid m\),即 \(m\) 恰好到了下一个块里。

于是这个预处理就是调和级数 \(O(n\log n)\) 的,由于 \(\mu\) 的性质,实际产生贡献的会更少。

这样大致处理 \(2\times 10^6\) 左右,再套类似杜教筛的东西就可以通过了。


另一个做法是使用欧拉反演,即利用 \(\varphi *\mathrm{I}=\mathrm{Id}\) 的性质,可以得到:

\[\begin{aligned} &\sum_{i_1=1}^m\sum_{i_2=1}^m\cdots\sum_{i_n=1}^m \gcd(i_1,i_2,\cdots,i_n)\\ =&\sum_{i_1=1}^m\sum_{i_2=1}^m\cdots\sum_{i_n=1}^m \sum_{d\mid \gcd(i_1,i_2,\cdots,i_n)} \varphi(d)\\ =&\sum_{d=1}^m \varphi(d) \left\lfloor\dfrac{m}{d}\right\rfloor^n \end{aligned} \]

杜教筛朴素计算即可。

PN 筛

Powerful Number

定义一个数 \(n\) 是 Powerful Number(PN) 当且仅当 \(n=\prod_{p\in \mathbf{P}} p_i^{c_i},c_i>1\)。

每个 PN 一定可以被唯一表示为 \(a^2b^3\) 的形式,次数为偶数的放在 \(a\) 部分,次数为奇数的将 \(3\) 次放在 \(b\) 部分,剩余放在 \(a\) 部分。

于是可以通过,,枚举 \(a\) 的值来计算 PN 的个数,大致是:

\[O\left(\sum_{i=1}^{\sqrt{n}}\sqrt[3]{\dfrac{n}{i^2}}\right)=O\left(\int_1^{\sqrt{n}} \sqrt[3]{\dfrac{n}{i^2}}\right)=O(\sqrt{n}) \]

所以可以通过 DFS 来得出所有 PN。

求积性函数前缀和

模板题为例。

给定积性函数 \(f\),求其前缀和。

考虑构造积性函数 \(g\),满足 \(g(p)=f(p)\),这里 \(g(p)=f(p)=p(p-1)\),可以构造 \(g=\varphi\cdot\mathrm{Id}\)。

再找到积性函数 \(h\),满足 \(f=g*h\),由于 \(f(p)=g(1)h(p)+g(p)h(1)\),积性函数 \(1\) 处点值一定为 \(1\),所以 \(f(p)=h(p)+g(p)\),即 \(h(p)=0\),这样所有可以被质因子筛到的点值处都是 \(0\),也就是所有非 PN 点值都是 \(0\)。

这样化简 \(f\) 前缀和的式子:

\[\begin{aligned} &\sum_{i=1}^n f(i)\\ =&\sum_{i=1}^n \sum_{j\mid i}g(i)h\left(\dfrac{i}{j}\right)\\ =&\sum_{i=1}^n h(i)S_g\left(\left\lfloor\dfrac{n}{i}\right\rfloor\right)\\ =&\sum_{\substack{i=1\\ \text{i is PN}}}^n h(i)S_g\left(\left\lfloor\dfrac{n}{i}\right\rfloor\right) \end{aligned}\]

\(g\) 的前缀和可以用杜教筛求出(有更优秀的也可以),计算答案在 DFS 求出所有 PN 时进行,其中求 \(h\) 可以求出所有 \(h(p^c)\) 再累乘。

求 \(h(p^c)\) 最常规的办法是利用 \(f=g*h\),移项得到 \(h(p^c)=f(p^c)-\sum_{i=0}^{c-1} h(p^i)g(p^{c-i})\)。

点击查看代码
ll n;
int pr[maptrn],phi[maptrn];
bool vis[maptrn];
int s1[maptrn],s2[maptrn];
bool mark[maptrn];
int pw[maptrn/10][40],h[maptrn/10][40];
inline void linear_sieve(){
    phi[1]=1;
    for(int i=2;i<=lim;++i){
        if(!vis[i]) pr[++pr[0]]=i,phi[i]=i-1;
        for(int j=1;j<=pr[0]&&i*pr[j]<=lim;++j){
            vis[i*pr[j]]=1,phi[i*pr[j]]=phi[i]*phi[pr[j]];
            if(i%pr[j]==0){
                phi[i*pr[j]]=phi[i]*pr[j];
                break;
            }        
        }
    }
    for(int i=1;i<=lim;++i) s1[i]=(s1[i-1]+1ll*i*phi[i]%mod)%mod;
    for(int i=1;i<=pr[0];++i){
        pw[i][0]=1,h[i][0]=1;
        for(ll j=pr[i],e=1;;j*=pr[i],++e){
            pw[i][e]=j%mod;
            int now=1ll*pw[i][e]*(pw[i][e]-1)%mod;
            for(int k=0;k<e;++k){
                now=(now-1ll*h[i][k]*(pr[i]-1)%mod*pw[i][e-k-1]%mod*pw[i][e-k]%mod+mod)%mod;
            }
            h[i][e]=now;
            if(j>n/pr[i]) break;
        }
    }
}
int get_sumg(ll N){
    if(N<=lim) return s1[N];
    if(mark[n/N]) return s2[n/N]; 
    mark[n/N]=1;
    int tmpn=N%mod;
    int res=1ll*tmpn*(tmpn+1)%mod*(2ll*tmpn%mod+1)%mod*inv6%mod;
    for(ll l=2,r;l<=N;l=r+1){
        r=N/(N/l);
        int tmpl=(l-1)%mod,tmpr=r%mod;
        int now=(1ll*tmpr*(tmpr+1)/2-1ll*tmpl*(tmpl+1)/2)%mod;
        res=(res-1ll*now*get_sumg(N/l)%mod+mod)%mod;
    }
    return s2[n/N]=res;
}
int ans;
void dfs(ll N,int ptr,int res){
    ans=(ans+1ll*res*get_sumg(n/N)%mod)%mod;
    if(ptr>pr[0]) return;
    for(int i=ptr;i<=pr[0];++i){
        if(N>n/(1ll*pr[i]*pr[i])) break;
        for(ll j=N*pr[i]*pr[i],e=2;;j*=pr[i],++e){
            dfs(j,i+1,1ll*res*h[i][e]%mod);
            if(j>n/pr[i]) break;
        }
    }
}

int main(){
    scanf("%lld",&n);
    linear_sieve();
    dfs(1,1,1);
    printf("%d\n",ans);
    return 0;
}

复杂度分析

求 \(g\) 前缀和复杂度认为是杜教筛复杂度 \(O(n^{2/3})\)。

只预处理 \([1,\sqrt{n}]\) 的质数,每个质数枚举 \(p^c\) 并计算 \(h\) 的复杂度是 \(O(\log^2 n)\),素数密度 \(\pi(n)=O\left(\dfrac{n}{\log n}\right)\),所以总复杂度是 \(O(\sqrt{n}\log n)\),实际这个上界比较松。

于是得到大致是 \(O(n^{2/3})\) 的复杂度。

算法流程总结

  1. 用积性函数拟合出合适的 \(g\)。

  2. 构造快速求 \(g\) 前缀和的方法。

  3. 预处理出 \(h(p^c)\)。

  4. DFS 出所有 PN 并计算答案。

例题

LOJ-6053 简单的函数

除了 \(2\) 以外的 \(p\) 处点值是 \(p-1\),先不考虑特例的话拟合成 \(\varphi\) 是非常好的,在此基础上 \(\varphi(2)=1\),所以将 \(2\) 的倍数处拟合成 \(3\varphi\)。(肯定不是 \(\varphi+2\) 吧)

\[g(n)=\begin{cases}3\varphi(n)&2\mid n\\\varphi(n)&2\nmid n\end{cases} \]

把 \(\varphi\) 提出一个会好一点,顺便把系数也带出来:

\[l(n)=\begin{cases}\varphi(n)&2\mid n\\0&2\nmid n\end{cases} \]

于是 \(g(n)=\varphi(n)+2l(n)\),也即 \(S_g(n)=S_\varphi(n)+2S_l(n)\)。

注意到 \(S_l(n)\) 有意义的 \(n\) 都是偶数,可以用一个 \(S(n)=\sum_{i=1}^n \varphi(2i)\) 来代替,即 \(S_l(n)=S\left(\left\lfloor\dfrac{n}{2}\right\rfloor\right)\)。

把 \(\varphi(2n)\) 详写成关于 \(n\) 的:

\[\varphi(2n)=\begin{cases}2\varphi(n)&2\mid n\\\varphi(n)&2\nmid n\end{cases} \]

参照上面的化简方法,得到:

\[\varphi(2n)-\varphi(n)=\begin{cases}\varphi(n)&2\mid n\\0&2\nmid n\end{cases} \]

这个和 \(l\) 一模一样,这里当然再把 \(S_l\) 换成 \(S\)。

\[S(n)=S_\varphi(n)+S\left(\left\lfloor\dfrac{n}{2}\right\rfloor\right) \]

放回最初的式子 \(S_g(n)=S_\varphi(n)+2S_h(n)=S_\varphi(n)+2S\left(\left\lfloor\dfrac{n}{2}\right\rfloor\right)\),这样一直把 \(S\) 展开,就是一个倍增的形式,单次计算是 \(O(\log n)\),所求点值都在杜教筛射程范围内。

构造的另一函数 \(h\) 就是朴素的 \(O(\sqrt{n}\log n)\) 预处理。

根据 PN 的数量级和预处理的复杂度,总体仍是 \(O(n^{2/3})\)。

Min_25 筛

Min_25 筛是一种基于埃筛的筛法,可以对数论函数 \(f(n)\) 求前缀和,要求 \(f(p)\) 是关于 \(p\) 的低次多项式,且 \(f(p^k)\) 可快速计算。

构造质数点值前缀和函数

在求解的开始,先用线性筛筛出小于等于 \(\sqrt{n}\) 的质数,定义 \(\mathrm{minp}(i)\) 表示 \(i\) 的最小质因子。

设 \(g(n)=\sum_{p\in \mathbf{P},p\le n} f(p)\),这个式子不好直接求,先设 \(g_k(n)=\sum_{p\in \mathbf{P},p\le n} p^k\),即 \(g(n)=\sum_{k} a_kg_k(n)\),其中 \(a_k\) 是 \(f(p)\) 这个关于 \(p\) 的低次多项式 \(k\) 次项系数。

设 \(g_k(n,j)=\sum_{i=1}^n i^k[i\in\mathbf{P}\lor\mathrm{minp}(i)>p_j]\),这个定义相当于埃筛筛掉了含小于等于 \(p_j\) 质因子的合数,容易发现找到最小的质数 \(p_x\) 满足 \(p_x\ge \sqrt{n}\),那么 \(g_k(n)=g_k(n,x)\),因为此时不再有最小值因子大于 \(p_x\) 的合数。

考虑递推,有:

\[g_k(n,j)= \begin{cases} g_k(n,j-1)&p_j^2\ge n\\ g_k(n,j-1)-p_j^k\left(g\left(\left\lfloor\dfrac{n}{p_j}\right\rfloor,j-1\right)-g(p_{j-1,j-1})\right)&p_j^2<n \end{cases} \]

第二个递推式实际上就是把 \(\mathrm{minp}(i)=p_j\) 的数拆成 \(p_j\) 和另一部分的乘积,这一部分就表现在 \(g\left(\left\lfloor\dfrac{n}{p_j}\right\rfloor,j-1\right)\) 中,但这里面还包含了质数,所以要减去。

于是只关心 \(O(\sqrt{n})\) 个块筛位置的点值,而 \(p_{j-1}\le \sqrt{n}\),所以一定可以被块筛筛到,这样计算复杂度也是经典方法了,具体是:

\[\begin{aligned} &O\left(\sum_{i=1}^{\sqrt{n}} \pi(\sqrt{i})+\sum_{i=1}^{\sqrt{n}}\pi\left(\sqrt{\dfrac{n}{i}}\right)\right)\\ =&O\left(\sum_{i=1}^{\sqrt{n}}\pi\left(\sqrt{\dfrac{n}{i}}\right)\right)\\ =&O\left(\int_1^{\sqrt{n}} \dfrac{\sqrt{\tfrac{n}{x}}}{\log \sqrt{\tfrac{n}{x}}} \mathrm{d}x\right)\\ =&O\left(\dfrac{\sqrt{n}\displaystyle\int_1^{\sqrt{n}} \dfrac{1}{\sqrt{x}}\mathrm{d}x}{\log n}\right)\\ =&O\left(\dfrac{n^{3/4}}{\log n}\right) \end{aligned}\]

构造前缀和函数

方法一:递归

设 \(S(n)=\sum_{i\le n}f(i)\),仿照上面,扩充 \(S\) 的定义为 \(S(n,j)=\sum_{i\le n}f(i)[\mathrm{minp}(i)>p_j]\),这样答案就是 \(S(n,0)+1\)。

直接递归,式子是:

\[S(n,j)=g(n)-g(p_j)+\sum_{\substack{p_j<p_k\le \sqrt{n}\\p_k^2\le p_k^{e+1}\le n}} f(p_k^e)S\left(\left\lfloor\dfrac{n}{p_k^e}\right\rfloor,k\right)+f(p_k^{e+1}) \]

就是分质数合数来计算,后半部分相当于枚举了一个最小值因子及幂次,然后递归求解,\(f(p_k^{e+1})\) 的加入源于 \(S\) 并不计算 \(1\) 的点值,因此要补充计算,而 \(f(p_k)\) 已经在质数处算过了。

复杂度证明见论文,结论是复杂度为 \(O(n^{1-\epsilon})\),但 \(n\le 10^{13}\) 时为 \(O\left(\dfrac{n^{3/4}}{\log n}\right)\)。

实现常数很小,但只能求单点。

点击查看代码
inline int S1(ll N){
    N%=mod;
    return N*(N+1)%mod*inv2%mod;
}
inline int S2(ll N){
    N%=mod;
    return N*(N+1)%mod*(2*N+1)%mod*inv6%mod;
}
inline int f(ll N){
    N%=mod;
    return N*(N-1+mod)%mod;
}

ll n;
int lim;
int pr[maxn],cntpr;
bool vis[maxn];
ll id[maxn<<1];
int id1[maxn],id2[maxn];
int g1[maxn<<1],g2[maxn<<1];

inline int get_id(ll N){
    if(N<=lim) return id1[N];
    else return id2[n/N];
}

inline void linear_sieve(){
    for(int i=2;i<=lim;++i){
        if(!vis[i]) pr[++cntpr]=i;
        for(int j=1;j<=cntpr&&i*pr[j]<=lim;++j){
            vis[i*pr[j]]=1;
            if(i%pr[j]==0) break;
        }
    }
}

int S(ll N,int j){
    if(pr[j]>=N) return 0;
    int res=0;
    res=(g1[get_id(N)]-g1[get_id(pr[j])]+mod)%mod;
    for(int k=j+1;k<=cntpr&&1ll*pr[k]*pr[k]<=N;++k){
        ll now=pr[k];
        for(int e=1;now<=N;++e,now*=pr[k]){
            res=(res+1ll*f(now)*(S(N/now,k)+(e>1))%mod)%mod;
        }
    }
    return res;
}

int main(){
    scanf("%lld",&n);
    lim=sqrt(n)+1;
    linear_sieve();
    for(ll l=1,r;l<=n;l=r+1){
        r=n/(n/l);
        id[++id[0]]=n/l;
        if(n/l<=lim) id1[n/l]=id[0];
        else id2[n/(n/l)]=id[0];
        g1[id[0]]=(S1(n/l)-1+mod)%mod,g2[id[0]]=(S2(n/l)-1+mod)%mod;
    }
    for(int j=1;j<=cntpr;++j){
        for(int i=1;i<=id[0]&&1ll*pr[j]*pr[j]<=id[i];++i){
            g1[i]=(g1[i]-1ll*pr[j]*(g1[get_id(id[i]/pr[j])]-g1[get_id(pr[j-1])]+mod)%mod+mod)%mod;
            g2[i]=(g2[i]-1ll*pr[j]*pr[j]%mod*(g2[get_id(id[i]/pr[j])]-g2[get_id(pr[j-1])]+mod)%mod+mod)%mod;
        }
    }
    for(int i=1;i<=id[0];++i) g1[i]=(g2[i]-g1[i]+mod)%mod;
    printf("%d\n",(S(n,0)+1)%mod);
    return 0;
}

方法二:递推

考虑把 \(S\) 的状态设计和 \(g\) 类比为 \(S(n,j)=\sum_{i\le n} f(i)[i\in\mathbf{P}\lor\mathrm{minp}(i)>p_j]\)。

这样目标是得到 \(S(n,0)+1\),初始 \(S(n,+\infty)=g(n)\)(实际上就是只有质数部分,第二维取第一个大于 \(\sqrt{n}\) 的质数编号即可)。

递推过程:

\[S(n,j)= \begin{cases} S(n,j+1)&p_{j+1}^2> n\\ S(n,j+1)+\displaystyle\sum_{p_{j+1}^2\le p_{j+1}^{e+1}\le n}f(p_{j+1}^e)\left(S\left(\left\lfloor\dfrac{n}{p_{j+1}^e}\right\rfloor\right)-S(p_{j+1})+f(p_{j+1}^{e+1})\right)&p_{j+1}^2\le n \end{cases} \]

这里对 \(p_{j+1}^{e+1}\) 做限制的原因是保证了 \(\left\lfloor\dfrac{n}{p_{j+1}^e}\right\rfloor\ge p_{j+1}\),从而消去 \(S(p_{j+1})\) 是正确的。

复杂度同上,但常数较大,优越性在于递推算法求出了块筛。

点击查看代码
inline int S1(ll N){
    N%=mod;
    return N*(N+1)%mod*inv2%mod;
}
inline int S2(ll N){
    N%=mod;
    return N*(N+1)%mod*(2*N+1)%mod*inv6%mod;
}
inline int f(ll N){
    N%=mod;
    return N*(N-1+mod)%mod;
}

ll n;
int lim;
int pr[maxn],cntpr;
bool vis[maxn];
ll id[maxn<<1];
int id1[maxn],id2[maxn];
int g1[maxn<<1],g2[maxn<<1];
int S[maxn<<1];

inline int get_id(ll N){
    if(N<=lim) return id1[N];
    else return id2[n/N];
}

inline void linear_sieve(){
    for(int i=2;i<=lim;++i){
        if(!vis[i]) pr[++cntpr]=i;
        for(int j=1;j<=cntpr&&i*pr[j]<=lim;++j){
            vis[i*pr[j]]=1;
            if(i%pr[j]==0) break;
        }
    }
}

int main(){
    scanf("%lld",&n);
    lim=sqrt(n)+1;
    linear_sieve();
    for(ll l=1,r;l<=n;l=r+1){
        r=n/(n/l);
        id[++id[0]]=n/l;
        if(n/l<=lim) id1[n/l]=id[0];
        else id2[n/(n/l)]=id[0];
        g1[id[0]]=(S1(n/l)-1+mod)%mod,g2[id[0]]=(S2(n/l)-1+mod)%mod;
    }
    for(int j=1;j<=cntpr;++j){
        for(int i=1;i<=id[0]&&1ll*pr[j]*pr[j]<=id[i];++i){
            g1[i]=(g1[i]-1ll*pr[j]*(g1[get_id(id[i]/pr[j])]-g1[get_id(pr[j-1])]+mod)%mod+mod)%mod;
            g2[i]=(g2[i]-1ll*pr[j]*pr[j]%mod*(g2[get_id(id[i]/pr[j])]-g2[get_id(pr[j-1])]+mod)%mod+mod)%mod;
        }
    }
    for(int i=1;i<=id[0];++i) S[i]=g1[i]=(g2[i]-g1[i]+mod)%mod;
    for(int j=cntpr;j>=1;--j){
        for(int i=1;i<=id[0]&&1ll*pr[j]*pr[j]<=id[i];++i){
            ll now=pr[j];
            for(int e=1;now*pr[j]<=id[i];++e,now*=pr[j]){
                S[i]=(S[i]+1ll*f(now)*(S[get_id(id[i]/now)]-S[get_id(pr[j])]+mod)%mod+f(now*pr[j]))%mod;
            }
        }
    }
    printf("%d\n",(S[1]+1)%mod);
    return 0;
}

例题

LOJ-6235 区间素数个数

只需要做 Min_25 的第一步,求 \(g\) 即可。

LOJ-6027 from CommonAnt 质数计数 I

就是求模 \(4\) 余 \(1\) 的质数个数,注意到求 \(g\) 部分是把 \(\mathrm{minp}(i)=p_j\) 分成了 \(p_j\) 和其他两部分,因此 \(g\left(\left\lfloor\dfrac{n}{p_j}\right\rfloor,j-1\right)\) 中每个数的余数乘上 \(p_j\) 的余数才是得到 \(\mathrm{minp}(i)=p_j\) 合数的余数。因此不能只求余数为 \(1\) 的,计算时按余数分组即可。

LOJ-6028 from CommonAnt 质数计数 II

与上一题类似

LOJ-6053 简单的函数

先按 \(f(p)=p-1\) 求 \(g\),第二步算之前给前缀和再加上 \(2\) 即可。

注意第一步只需要保证完全积性,第二步只需要保证积性,所以这样处理是正确的。

LOJ-6785 简单的函数 10^13 版

要求块筛,用递推方法,需要大量卡常。

LOJ-6181 某个套路求和题

\(f(n)\) 不是积性函数,但是很简洁,容易发现质数位置是 \(-1\),含平方因子位置是 \(0\)。

不含平方因子的合数位置,设共 \(k\) 个质因子,那么只关心 \(\mu(d)=-1\) 的位置,也就是 \(k\) 中选出奇数个,方案数 \(2^{k-1}\),因此不含平方因子的合数位置是 \(1\)。

只看合数位置可以拟合成 \(\mu^2(n)\),这样质数位置从 \(-1\) 变成了 \(1\),答案就应该是:

\[\sum_{i\le n}\mu^2(i)-2\sum_{i\le n}[i\in \mathrm{P}] \]

前一个式子直接 Min_25 前缀和即可,容易发现后一个式子等价于第一步的 \(g\) 函数。

参考资料

狄利克雷卷积

杜教筛

PN 筛

Min_25 筛

标签:lfloor,right,筛法,狄利克,卷积,dfrac,sum,mathrm,left
From: https://www.cnblogs.com/SoyTony/p/Learning_Notes_about_Dirichlet_Convolution_and_Senior_Si

相关文章

  • Tensorflow整理[6].卷积神经网络
    概述对CIFAR-10数据集的分类是机器学习中一个公开的基准测试问题,其任务是对一组32x32RGB的图像进行分类,这些图像涵盖了10个类别:飞机,汽车,鸟,猫,鹿,狗,青蛙,马,船以及卡车。重点CIFAR-10教程演示了在TensorFlow上构建更大更复杂模型的几个种重要内容:相关核心数学对象,如卷积、修正......
  • TensorFlow10.4 卷积神经网络-ResNet与DenseNet及ResNet实战
    1ResNet我们是实验发现在我们堆叠更多的网络结构的时候,我们并不能又一个很好的结果,就是它网络层次变多了之后他会产生一个多层的loss的堆叠,使得梯度爆炸,或者梯度弥散。然后我们想了一个办法,就是我们比如说设置了一个30层的神经网络,我们在差也不能比22层的差。就是我们设置了一......
  • 深度卷积神经网络(AlexNet)
    1.AlexNet\(2012\)年,\(AlexNet\)横空出世。使用\(8\)层卷积神经网络,赢得\(ImageNet\2012\)图像识别挑战赛。\(AlexNet\) 网络结构:1.1第一个卷积层卷积运算:原始数据为\(227\times227\times3\)的图像。卷积核尺寸\(11\times11\times3\),步长\(4\),每次......
  • U-Net: 专注生物医学分割的卷积神经网络(翻译)
    原文链接:https://arxiv.org/pdf/1505.04597.pdf摘要:普遍认为,优秀的深度神经网络离不开数千个标注训练样本。在本文中,我们提出了一种网络和训练策略:该策略通过使用大量数据增强,从而充分利用带标注的训练样本;该网络结构包括了用于捕获上下文的收缩路径和用于实现精确定位的对称扩......
  • TensorFlow10.4 卷积神经网络-batchnorm
    我们发现这个sigmoid函数在小于-4或者大于4的时候他的导数趋近于0。然后我们送进去的input的值在[-100,100]之间,这样很容易引起梯度弥散的现象。所以我们一般情况下使用ReLU函数,但是我们有时候又不得不使用sigmoid函数。这个时候我们在送到下一层的时候我们应该先经过Normalizatio......
  • TensorFlow10.3 卷积神经网络-经典卷积网络(VGG,GoogLeNet)
    LeNet-5这个是5层的,3个c+s,然后有两个全连接层。AlexNet这里有8(5+3)层。就是之前的技术没有现在的好,所以它用了两块GTX580,然后让你它的模型分成两块,然后在两块显卡中跑。很好的把显存给分开来了。VGG之前都是用\(11*11\)的窗口,然后它用了\(3*3\)的窗口,这个\(3*3\)的窗......
  • TensorFlow10.2 卷积神经网络-CIFAR100 实战
    ▪Loaddatasets▪BuildNetwork▪Train▪Test这里先是进行卷积然后再进行全连接Loaddatasetsdefpreprocess(x,y):#[0~1]x=tf.cast(x,dtype=tf.float32)/255.y=tf.cast(y,dtype=tf.int32)returnx,y(x,y),(x_test,y_test)=dat......
  • TensorFlow10.2 卷积神经网络-卷积神经网络池化层与采样
    ▪Pooling▪upsample▪ReLU我们看一下这个Subsampling层就是这个:这一层起到ReduceDim的作用。1Max/Avgpooling(下采样)keras.layers.MaxPooling2D(pool_size=,strides=,padding='valid',data_format=None)pool_size:池化窗口大小strides:池化步长,默认值等于p......
  • 垃圾识别系统Python+TensorFlow+Django+卷积神经网络算法【完整代码系统】
    一、介绍垃圾识别系统,使用Python作为主要开发语言,基于深度学习TensorFlow框架,搭建卷积神经网络算法。并通过对5种垃圾数据集进行训练,最后得到一个识别精度较高的模型。并基于Django,开发网页端操作平台,实现用户上传一张垃圾图片识别其名称。二、效果展示三、演示视频+代码视......
  • 交通标志识别系统Python+TensorFlow+Django+卷积神经网络算法实现【完整代码】
    一、介绍使用Python作为主要开发语言,基于深度学习TensorFlow框架,搭建卷积神经网络算法。并通过对数据集进行训练,最后得到一个识别精度较高的模型。并基于Django,开发网页端操作平台,实现用户上传一张图片识别其名称。二、效果展示三、演示视频视频+完整代码:https://www.yuque.......