狄利克雷生成函数学习笔记 gxy001

26 Jan 2021 11378字 38分 186次 数学 and 数论
如果这篇博客帮助到你，可以请我喝一杯咖啡~
CC BY-NC-SA 4.0 （除特别声明或转载文章外）

狄利克雷生成函数（DGF\text{DGF}）浅谈

前置知识：无。

数论函数

定义域为正整数，陪域为复数的函数。

狄利克雷卷积

对于数论函数 f,gf,g，定义它们的狄利克雷卷积为 (f∗g)(x)=∑d∣xf(d)g(xd)(f*g)(x)=\sum\limits_{d\mid x}f(d)g(\frac{x}{d})。

换种容易理解的表示方式即 (f∗g)(x)=∑a×b=xf(a)g(b)(f*g)(x)=\sum\limits_{a\times b=x}f(a)g(b)。

狄利克雷生成函数（DGF\text{DGF}）

数论函数 ff 在 ii 处的点值表示为 fif_i，则 ff 的狄利克雷生成函数 F(x)=∑i=1∞fiixF(x)=\sum\limits_{i=1}^\infin \frac{f_i}{i^x}。

黎曼 ζ\zeta 函数

ζ(x)=∑i=1∞1ix\zeta(x)=\sum\limits_{i=1}^\infin \frac{1}{i^x}，容易发现，这是常数函数 I(x)=1I(x)=1 的 DGF\text{DGF}。

这个函数与 OGF\text{OGF} 的 11−x\frac{1}{1-x}，EGF\text{EGF} 的 exe^x 一样重要。

积性函数

若 ∀gcd⁡(i,j)=1,f(i×j)=f(i)×f(j)\forall \gcd(i,j)=1,f(i\times j)=f(i)\times f(j)，则称数论函数 ff 为积性函数。

积性函数与 DGF\text{DGF}

考虑将每个质数的贡献分开计算：

下文中 Prime\mathrm{Prime} 指全体质数集合。

单位函数 ϵ\epsilon

11

常数函数 II

ζ(x)=∏p∈Prime∑i=0∞1pix=∏p∈Prime11−p−x\begin{aligned} \zeta(x)=&\prod\limits_{p\in \mathrm{Prime}}\sum\limits_{i=0}^\infin \frac{1}{p^{ix}}\\ =&\prod\limits_{p\in \mathrm{Prime}}\frac{1}{1-p^{-x}} \end{aligned}

莫比乌斯函数 μ\mu

∏p∈Prime(1+−1px)=∏p∈Prime(1−p−x)=1ζ(x)\begin{aligned} &\prod\limits_{p\in\mathrm{Prime}}(1+\frac{-1}{p^x})\\ =&\prod\limits_{p\in\mathrm{Prime}}(1-p^{-x})\\ =&\frac{1}{\zeta(x)} \end{aligned}

刘维尔函数 λ\lambda

∏p∈Prime∑i=0∞(−1)ipix=∏p∈Prime11+p−x=ζ(2x)ζ(x)\begin{aligned} &\prod\limits_{p\in\mathrm{Prime}}\sum\limits_{i=0}^\infin\frac{(-1)^i}{p^{ix}}\\ =&\prod\limits_{p\in\mathrm{Prime}}\frac{1}{1+p^{-x}}\\ =&\frac{\zeta(2x)}{\zeta(x)} \end{aligned}

μ2\mu^2（点积）

∏p∈Prime(1+1px)=∏p∈Prime(1+p−x)=ζ(x)ζ(2x)\begin{aligned} &\prod\limits_{p\in\mathrm{Prime}}(1+\frac{1}{p^x})\\ =&\prod\limits_{p\in\mathrm{Prime}}(1+p^{-x})\\ =&\frac{\zeta(x)}{\zeta(2x)} \end{aligned}

恒等函数 idkid_k

∑i=1∞ikix=∑i=1∞1ix−k=ζ(x−k)\begin{aligned} &\sum\limits_{i=1}^\infin \frac{i^k}{i^x}\\ =&\sum\limits_{i=1}^\infin \frac{1}{i^{x-k}}\\ =&\zeta(x-k) \end{aligned}

欧拉函数 φ\varphi

∏p∈Prime(1+∑i=1∞pi−pi−1pix)=∏p∈Prime(1+∑i=1∞pi−ix−∑i=1∞pi−ix−1)=∏p∈Prime(1+p1−x1−p1−x−p−x1−p1−x)=∏p∈Prime1−p−x1−p1−x=ζ(x−1)ζ(x)\begin{aligned} &\prod\limits_{p\in\mathrm{Prime}}\left(1+\sum\limits_{i=1}^\infin\frac{p^i-p^{i-1}}{p^{ix}}\right)\\ =&\prod\limits_{p\in\mathrm{Prime}}\left(1+\sum\limits_{i=1}^\infin p^{i-ix}-\sum\limits_{i=1}^\infin p^{i-ix-1}\right)\\ =&\prod\limits_{p\in\mathrm{Prime}}\left(1+\frac{p^{1-x}}{1-p^{1-x}}-\frac{p^{-x}}{1-p^{1-x}}\right)\\ =&\prod\limits_{p\in\mathrm{Prime}}\frac{1-p^{-x}}{1-p^{1-x}}\\ =&\frac{\zeta(x-1)}{\zeta(x)}\\ \end{aligned}

除数函数 σk\sigma_k

∏p∈Prime∑i=0∞∑j=0ipjkpix=ζ(x)ζ(x−k)\begin{aligned} &\prod\limits_{p\in\mathrm{Prime}}\sum\limits_{i=0}^\infin\frac{\sum\limits_{j=0}^ip^{jk}}{p^{ix}}\\ =&\zeta(x)\zeta(x-k) \end{aligned}

关于这个的证明，等一会再说。

一个有趣的事实

已知 DGF\text{DGF} F(x)F(x)，该数论函数点积 idkid_k 的 DGF\text{DGF} 为 F(x−k)F(x-k)，证明留给读者做练习。

DGF\text{DGF} 的乘法

显然有，(F∗G)(x)=∑i=1∞∑d∣ifdgidix(F*G)(x)=\sum\limits_{i=1}^\infin\frac{\sum\limits_{d\mid i}f_{d}g_{\frac{i}{d}}}{i^x}，两个 DGF\text{DGF} 的乘积就是这两个数论函数狄利克雷卷积的 DGF\text{DGF}。

我们可以通过这个做很多事，例如除数函数 DGF\text{DGF} 是 ζ(x)ζ(x−k)\zeta(x)\zeta(x-k) 可以通过这个性质给出一个很简单的证明（idk∗I=σkid_k*I=\sigma_k）。

一些不那么显然的狄利克雷卷积可以被 DGF\text{DGF} 乘法轻松地证明：

注：d=σ0,σ=σ1d=\sigma_0,\sigma=\sigma_1。

ζ(x−1)ζ(x)×ζ(x)=ζ(x−1)⇒φ∗I=id\frac{\zeta(x-1)}{\zeta(x)}\times\zeta(x)=\zeta(x-1)\Rightarrow\varphi*I=id

ζ(x−1)ζ(x)×ζ(x)2=ζ(x)ζ(x−1)⇒φ∗d=σ\frac{\zeta(x-1)}{\zeta(x)}\times\zeta(x)^2=\zeta(x)\zeta(x-1)\Rightarrow\varphi*d=\sigma

直接卷积的复杂度为 O(nlog⁡n)\mathrm O(n\log n)。

void mul(int const *f,int const *g,int *ans,int n){
    for(int i=1;i<=n;i++)
        for(int j=i;j<=n;j+=i)
            ans[j]+=f[i]*g[j/i];
}

DGF\text{DGF} 的除法

已知数论函数 F,GF,G，求 HH 使得 F=H∗GF=H*G。 Fn=∑d∣nHdGndHnG1=Fn−∑d∣n∧d≠nHdGndF_n=\sum\limits_{d\mid n}H_dG_{\frac{n}{d}}\\ H_nG_1=F_n-\sum\limits_{d\mid n\land d\ne n}H_dG_{\frac{n}{d}} 求出 HnH_n 后就枚举倍数更新即可。

void div(int const *f,int const *g,int *ans,int n){
    for(int i=1;i<=n;i++) ans[i]=f[i];
    for(int i=1;i<=n;i++){
        ans[i]/=g[1];
        for(int j=i+i;j<=n;j+=i)
            ans[j]-=ans[i]*g[j/i];
    }
}

DGF\text{DGF} 的求导与积分

dfnnxdx=−ln⁡nfnnx∫fnnxdx=−1ln⁡nfnnx+C\frac{d\frac{f_n}{n^x}}{dx}=-\ln n\frac{f_n}{n^x}\\ \int\frac{f_n}{n^x}dx=-\frac{1}{\ln n}\frac{f_n}{n^x}+C

n=1n=1 时特殊处理一下，这里使用质因子次数和的相反数作为 ln⁡n\ln n，它与 ln⁡n\ln n 有着类似的性质，我们并不在意 ln⁡n\ln n 的真实值，ln⁡n\ln n 一定会被消掉，所以直接用质因子次数和的相反数代替即可。

void der(int const *f,int *ans,int n){
    for(int i=1;i<=n;i++) ans[i]=f[i]*c[i];
}
void inte(int const *f,int *ans,int n){
    for(int i=1;i<=n;i++) ans[i]=f[i]/c[i];
}

DGF\text{DGF} 的对数函数

ln⁡(F(x))=∫F′(x)F(x)dx\ln(F(x))=\int\frac{F'(x)}{F(x)}dx
void ln(int const *f,int *ans,int n){
    static int tmp[maxn];
    der(f,tmp,n);
    div(tmp,f,ans,n);
    inte(ans,ans,n);
}

DGF\text{DGF} 的指数函数

exp⁡(F(x))=G(x)G′(x)=F′(x)exp(F(x))=F′(x)G(x)gnln⁡n=∑d∣ngndfdln⁡d\exp(F(x))=G(x)\\ G'(x)=F'(x)exp(F(x))=F'(x)G(x)\\ g_n\ln n=\sum\limits_{d\mid n}g_{\frac{n}{d}}f_d\ln d\\
void exp(int const *f,int *ans,int n){
    static int tmp[maxn];
	der(f,tmp,n);
	for(int i=1;i<=n;i++)ans[i]=0;
	ans[1]=1;
    for(int i=1;i<=n;i++){
        if(i!=1)ans[i]=1ll*ans[i]*pow(c[i],mod-2)%mod;
        for(int j=i+i;j<=n;j+=i)
            ans[j]=(ans[j]+1ll*ans[i]*tmp[j/i])%mod;
    }
}

狄利克雷 k 次根 加强版

已知数论函数 gg，求 ff 使得 fk=gf^k=g。 f=gkF=exp⁡(ln⁡(G)k)f=\sqrt[k]{g}\\ F=\exp(\frac{\ln(G)}{k}) 代码：

#include<cstdio>
int const mod=998244353;
int n,k,f[1000010],g[1000010],c[1000010];
int pow(int x,int y){
	int res=1;
	while(y){
		if(y&1) res=1ll*x*res%mod;
		x=1ll*x*x%mod;
		y>>=1;
	}
	return res;
}
void div(int const *f,int const *g,int *ans,int n){
    for(int i=1;i<=n;i++) ans[i]=f[i];
    int inv=pow(g[1],mod-2);
    for(int i=1;i<=n;i++){
        ans[i]=1ll*ans[i]*inv%mod;
        for(int j=i+i;j<=n;j+=i)
            ans[j]=(ans[j]-1ll*ans[i]*g[j/i]%mod+mod)%mod;
    }
}
void der(int const *f,int *ans,int n){
    for(int i=1;i<=n;i++) ans[i]=1ll*f[i]*c[i]%mod;
}
void inte(int const *f,int *ans,int n){
    for(int i=1;i<=n;i++) ans[i]=1ll*f[i]*pow(c[i],mod-2)%mod;
}
void ln(int const *f,int *ans,int n){
    static int tmp[1000010];
    der(f,tmp,n);
    div(tmp,f,ans,n);
    inte(ans,ans,n);
}
void exp(int const *f,int *ans,int n){
    static int tmp[1000010];
	der(f,tmp,n);
	for(int i=1;i<=n;i++)ans[i]=0;
	ans[1]=1;
    for(int i=1;i<=n;i++){
        if(i!=1)ans[i]=1ll*ans[i]*pow(c[i],mod-2)%mod;
        for(int j=i+i;j<=n;j+=i)
            ans[j]=(ans[j]+1ll*ans[i]*tmp[j/i])%mod;
    }
}
int np[1000010],p[1000010],cnt;
int main(){
	scanf("%d%d",&n,&k);
	for(int i=2;i<=n;i++){
		if(!np[i])p[++cnt]=i,c[i]=1;
		for(int j=1;j<=cnt&&i*p[j]<=n;j++){
			np[i*p[j]]=1,c[i*p[j]]=c[i]+1;
			if(i%p[j]==0)break;
		}
	}
	k=pow(k,mod-2);
	for(int i=1;i<=n;i++)scanf("%d",f+i);
	ln(f,g,n);
	for(int i=1;i<=n;i++)g[i]=1ll*g[i]*k%mod;
	exp(g,f,n);
	for(int i=1;i<=n;i++)printf("%d ",f[i]);
	return 0;
}

下面将介绍 DGF\text{DGF} 的用处。

构造杜教筛

先从几个简单的例子开始说起：

【模板】杜教筛

φ=ζ(x−1)ζ(x)ζ(x−1)ζ(x)×ζ(x)=ζ(x−1)⇒φ∗I=id\varphi=\frac{\zeta(x-1)}{\zeta(x)}\\ \frac{\zeta(x-1)}{\zeta(x)}\times\zeta(x)=\zeta(x-1)\\ \Rightarrow \varphi*I=id μ=1ζ(x)1ζ(x)×ζ(x)=1⇒μ∗I=ϵ\mu=\frac{1}{\zeta(x)}\\ \frac{1}{\zeta(x)}\times\zeta(x)=1\\ \Rightarrow \mu*I=\epsilon

看完这两个就应该能理解如何构造一个好的杜教筛了。

来点有难度的。

【模板】Min_25筛

我们先尝试把这个东西的 DGF\text{DGF} 用 ζ\zeta 表示。

∏p∈Prime(1+∑i=1∞pi(pi−1)pix)=∏p∈Prime(1+∑i=1∞p2ipix−∑i=1∞pipix)=∏p∈Prime(1+∑i=1∞pi(2−x)−∑i=1∞pi(1−x))=∏p∈Prime(1+p2−x1−p2−x−p1−x1−p1−x)=∏p∈Prime(1−p2−x)(1−p1−x)+p2−x(1−p1−x)−p1−x(1−p2−x)(1−p2−x)(1−p1−x)=∏p∈Prime1−p1−x−p2−x+p3−2x+p2−x−p3−2x−p1−x+p3−2x1−p1−x−p2−x+p3−2x=∏p∈Prime1−2p1−x+p3−2x1−p1−x−p2−x+p3−2x=ζ(x−1)ζ(x−2)∏p∈Prime(1−2p1−x+p3−2x)\begin{aligned} &\prod\limits_{p\in\mathrm{Prime}}\left(1+\sum\limits_{i=1}^\infin\frac{p^i(p^i-1)}{p^{ix}}\right)\\ =&\prod\limits_{p\in\mathrm{Prime}}\left(1+\sum\limits_{i=1}^\infin\frac{p^{2i}}{p^{ix}}-\sum\limits_{i=1}^\infin\frac{p^i}{p^{ix}}\right)\\ =&\prod\limits_{p\in\mathrm{Prime}}\left(1+\sum\limits_{i=1}^\infin p^{i(2-x)}-\sum\limits_{i=1}^\infin p^{i(1-x)}\right)\\ =&\prod\limits_{p\in\mathrm{Prime}}\left(1+\frac{p^{2-x}}{1-p^{2-x}}-\frac{p^{1-x}}{1-p^{1-x}}\right)\\ =&\prod\limits_{p\in\mathrm{Prime}}\frac{(1-p^{2-x})(1-p^{1-x})+p^{2-x}(1-p^{1-x})-p^{1-x}(1-p^{2-x})}{(1-p^{2-x})(1-p^{1-x})}\\ =&\prod\limits_{p\in\mathrm{Prime}}\frac{1-p^{1-x}-p^{2-x}+p^{3-2x}+p^{2-x}-p^{3-2x}-p^{1-x}+p^{3-2x}}{1-p^{1-x}-p^{2-x}+p^{3-2x}}\\ =&\prod\limits_{p\in\mathrm{Prime}}\frac{1-2p^{1-x}+p^{3-2x}}{1-p^{1-x}-p^{2-x}+p^{3-2x}}\\ =&\zeta(x-1)\zeta(x-2)\prod\limits_{p\in\mathrm{Prime}}(1-2p^{1-x}+p^{3-2x})\\ \end{aligned}

我们发现后面那个东西似乎不能很好的用 ζ\zeta 表示，考虑换个方向推。

如果有一个积性函数的 DGF\text{DGF} 为 ∏p∈Prime(1+∑i=2∞f(pi)pix)\prod\limits_{p\in\mathrm{Prime}}(1+\sum\limits_{i=2}^\infin\frac{f(p^i)}{p^{ix}})，那么我们是可以在 O(n)\mathrm O(\sqrt n) 内求出它的前缀和的，因为 [1,n][1,n] 之间的 Powerful Number\text{Powerful Number}（所有质因子次数都大于 11 的数）只有 O(n)\mathrm O(\sqrt n) 个。

证明：显然所有 Powerful Number\text{Powerful Number} 都可以表示成 a2b3a^2b^3 的形式。

O(∑a=1n(na2)13)=O(∫1n(nx2)13dx)=O(n)\begin{aligned} &\mathrm O\left(\sum_{a=1}^{\sqrt n}\left(\frac{n}{a^2}\right)^{\frac{1}{3}}\right)\\ =&\mathrm O\left(\int_{1}^{\sqrt n}\left(\frac{n}{x^2}\right)^{\frac{1}{3}}dx\right)\\ =&\mathrm O(\sqrt n)\\ \end{aligned}

我们现在想要把后面那个东西推成 ∏p∈Prime(1+∑i=2∞f(pi)pix)\prod\limits_{p\in\mathrm{Prime}}(1+\sum\limits_{i=2}^\infin\frac{f(p^i)}{p^{ix}}) 的形式。

ζ(x−1)ζ(x−2)∏p∈Prime(1−2p1−x+p3−2x)=ζ(x−1)ζ(x−2)∏p∈Prime(1−2p1−x+p2(1−x)+1)=ζ(x−1)ζ(x−2)∏p∈Prime((1−p1−x)2+p2(1−x)+1−p2(1−x))=ζ(x−2)ζ(x−1)∏p∈Prime(1−p1−x)2+p2(1−x)+1−p2(1−x)(1−p1−x)2=ζ(x−2)ζ(x−1)∏p∈Prime(1+p2(1−x)+1(1−p1−x)2−p2(1−x)(1−p1−x)2)=ζ(x−2)ζ(x−1)∏p∈Prime(1+∑i=2∞(i−1)pi(1−x)+1−∑i=2∞(i−1)pi(1−x))=ζ(x−2)ζ(x−1)∏p∈Prime(1+∑i=2∞ipi+1−pi+1−ipi+pipix)\begin{aligned} &\zeta(x-1)\zeta(x-2)\prod\limits_{p\in\mathrm{Prime}}(1-2p^{1-x}+p^{3-2x})\\ =&\zeta(x-1)\zeta(x-2)\prod_{p\in\mathrm{Prime}}(1-2p^{1-x}+p^{2(1-x)+1})\\ =&\zeta(x-1)\zeta(x-2)\prod_{p\in\mathrm{Prime}}((1-p^{1-x})^2+p^{2(1-x)+1}-p^{2(1-x)})\\ =&\frac{\zeta(x-2)}{\zeta(x-1)}\prod_{p\in\mathrm{Prime}}\frac{(1-p^{1-x})^2+p^{2(1-x)+1}-p^{2(1-x)}}{(1-p^{1-x})^2}\\ =&\frac{\zeta(x-2)}{\zeta(x-1)}\prod_{p\in\mathrm{Prime}}\left(1+\frac{p^{2(1-x)+1}}{(1-p^{1-x})^2}-\frac{p^{2(1-x)}}{(1-p^{1-x})^2}\right)\\ =&\frac{\zeta(x-2)}{\zeta(x-1)}\prod_{p\in\mathrm{Prime}}\left(1+\sum_{i=2}^\infin (i-1)p^{i(1-x)+1}-\sum_{i=2}^\infin (i-1)p^{i(1-x)}\right)\\ =&\frac{\zeta(x-2)}{\zeta(x-1)}\prod_{p\in\mathrm{Prime}}\left(1+\sum_{i=2}^\infin\frac{ip^{i+1}-p^{i+1}-ip^{i}+p^i}{p^{ix}}\right)\\ \end{aligned}

设要求的函数为 ff，ζ(x−2)ζ(x−1)\frac{\zeta(x-2)}{\zeta(x-1)} 对应的函数为 gg，∏p∈Prime(1+∑i=2∞ipi+1−pi+1−ipi+pipix)\prod\limits_{p\in\mathrm{Prime}}(1+\sum_{i=2}^\infin\frac{ip^{i+1}-p^{i+1}-ip^{i}+p^i}{p^{ix}}) 对应的函数为 hh。

则有 f=g∗hf=g*h，gg 可以直接杜教筛，至于杜教筛的构造，留给读者作为练习，hh 显然只在 Powerful Number\text{Powerful Number} 处有值，直接爆搜所有小于 n\sqrt n 的质数即可，算到每一个 Powerful Number\text{Powerful Number} 处统计 h(x)(∑i=1⌊nx⌋g(i))h(x)(\sum\limits_{i=1}^{\lfloor\frac{n}{x}\rfloor} g(i)) 即可求出 ∑i=1nf(i)\sum\limits_{i=1}^n f(i)。

代码：

#include<cstdio>
#include<cmath>
int const mod=1000000007,inv6=1000000008/6,inv2=1000000008/2;
int np[4641600],lim,cnt,p[464160],g[4641600],g2[4641600];
long long n;
int getsum(long long x){
	if(x<=lim) return g[x];
	if(g2[n/x]) return g2[n/x];
	int ans=x%mod*((x+1)%mod)%mod*((2*x+1)%mod)%mod*inv6%mod;
	for(long long l=2,r,d;l<=x;l=r+1){
		r=x/(x/l);
		ans=(ans-(l+r)%mod*(r-l+1)%mod*inv2%mod*getsum(x/l)%mod+mod)%mod;
	}
	return g2[n/x]=ans;
}
int ans;
void dfs(int k,long long m,int h){
	if(k>cnt||m*p[k]>n||m*p[k]*p[k]>n){
		long long const &p=n/m;
		if(p<=lim)ans=(ans+1ll*h*g[p])%mod;
		else ans=(ans+1ll*h*g2[n/p])%mod;
		return;
	}
	long long p=1ll*::p[k]*::p[k];
	dfs(k+1,m,h);
	for(int e=2;m*p<=n;p*=::p[k],++e)dfs(k+1,m*p,p%mod*(::p[k]-1)%mod*(e-1)%mod*h%mod);
}
int main(){
	scanf("%lld",&n);
	lim=pow(n,2.0/3.0);if(!lim)++lim;
	g[1]=1;
	for(int i=2;i<=lim;i++){
		if(!np[i])p[++cnt]=i,g[i]=i-1;
		for(int j=1,tmp;j<=cnt&&(tmp=i*p[j])<=lim;j++){
			np[tmp]=1;
			if(i%p[j]==0){
				g[tmp]=g[i]*p[j];
				break;
			}else g[tmp]=g[i]*g[p[j]];
		}
	}
	for(int i=1;i<=lim;i++)g[i]=(g[i-1]+1ll*i*g[i])%mod;
	getsum(n);
	dfs(1,1,1);
	printf("%d\n",ans);
	return 0;
}

推通项

DETER2

容易发现，这道题在高斯消元完成后 m’i,j={0i∤jf(i)i∣jm’_{i,j}=\begin{cases}0 & i\nmid j \\ f(i) & i\mid j\end{cases}。

原因是只有 ii 的因数行会影响第 ii 行的值，而对于 i∣ji\mid j 的 jj，这些行的对应位置都相等。

而对于 i∤ji\nmid j 的位置，因为 gcd⁡(i,j)∣i，gcd⁡(gcd⁡(i,j),j)=gcd⁡(i,j)\gcd(i,j)\mid i，\gcd(\gcd(i,j),j)=\gcd(i,j)，所以在消元到第 gcd⁡(i,j)\gcd(i,j) 行时，mi,j=mgcd⁡(i,j),jm_{i,j}=m_{\gcd(i,j),j}，这样一减就变成 00 了。

我们有 mi,i=gcd⁡(i,i)k−∑d∣i,d≠imd,im_{i,i}=\gcd(i,i)^k-\sum\limits_{d\mid i,d\ne i}m_{d,i} 即 f(i)=ik−∑d∣if(d)f(i)=i^k-\sum\limits_{d\mid i}f(d)，也就是 ∑d∣if(i)=idk(i)\sum\limits_{d\mid i}f(i)=id_k(i)，到这一步写个狄利克雷除法就可以 O(nlog⁡n)O(n\log n) 的做了。

但根据数据范围看，我们应该需要 O(n)O(n) 的解法，所以继续推，设 ff 的 DGF 为 F(x)F(x)，我们已知 11 的 DGF 为 ζ(x)\zeta(x)，idkid_k 的 DGF 为 ζ(x−k)\zeta(x-k)，所以 F(x)ζ(x)=ζ(x−k)F(x)\zeta(x)=\zeta(x-k)，F(x)=ζ(x−k)ζ(x)F(x)=\frac{\zeta(x-k)}{\zeta(x)}，事实上，k=1k=1 时 ff 就是 φ\varphi。

F(x)=ζ(x−k)ζ(x)=∏p∈Prime1−p−x1−p−x+k=∏p∈Prime(1+∑i=1∞pik−p(i−1)kpix)\begin{aligned} F(x)=&\frac{\zeta(x-k)}{\zeta(x)}\\ =&\prod\limits_{p\in \mathrm{Prime}} \frac{1-p^{-x}}{1-p^{-x+k}}\\ =&\prod\limits_{p\in \mathrm{Prime}}\left(1+\sum\limits_{i=1}^\infin\frac{p^{ik}-p^{(i-1)k}}{p^{ix}}\right) \end{aligned}

所以 ff 在素数幂 pip^i 处的取值为 pik−p(i−1)kp^{ik}-p^{(i-1)k}，线性筛即可。

答案即为 ∏i=1nf(i)\prod\limits_{i=1}^nf(i)。

#include<cstdio>
int pw[1000010],f[1000010],t,n,k,np[1000010],p[1000010];
int const mod=1e6+3;
int fpow(int x,int y){
	int res=1;
	while(y){
		if(y&1)res=1ll*res*x%mod;
		x=1ll*x*x%mod;
		y>>=1;
	}
	return res;
}
int main(){
	scanf("%d",&t);
	while(t--){
		scanf("%d%d",&n,&k);
		f[1]=1;
		int ans=1,cnt=0;
		for(int i=2;i<=n;i++){
			if(!np[i])p[++cnt]=i,f[i]=(mod+(pw[i]=fpow(i,k))-1)%mod;
			for(int j=1;j<=cnt&&p[j]*i<=n;j++){
				np[i*p[j]]=1;
				if(i%p[j]) f[i*p[j]]=1ll*f[i]*f[p[j]]%mod;
				else{f[i*p[j]]=1ll*f[i]*pw[p[j]]%mod;break;}
			}
			ans=1ll*ans*f[i]%mod;
		}
		printf("%d\n",ans);
	}
	return 0;
}
3 评论
NFLSCode Chrome 91.0.4472.114 Windows 7
2021-07-07回复

tql

JaCky080624 Chrome 91.0.4472.114 Windows 7
2021-07-02回复

qp
感觉挺有用的，以后请f们看一下blackness

Acc Robin Chrome 88.0.4324.190 Windows 7
2021-03-04回复

%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Powered By Valine
v1.4.14
欢迎来到本站！