写在前面:下面除法默认下取整,为了方便不写分数,没有说明 (i,j) 都表示取 gcd。
希望下面的每一题都能让你们学到一个 trick,本文不适合入门,建议先学习一些推式子知识以及筛发(杜教筛,min25 筛)再来看毒瘤的推柿子。
两个前置知识(推柿子通法)#
定义 ϵ(n)=[n=1]。我们有:n=∑d∣nφ(d) 和 ϵ(n)=∑d∣nμ(d) 。
对于每个数 x 它和 n 的 gcd 是一定的。于是有 n=n∑i=1n∑j=1[(j,n)=i]=∑i∣n∑i∣j[(j/i,n/i)=1]=∑i∣nn/i∑j=1[(j,n/i)=1]=∑i∣nφ(n/i)=∑i∣nφ(i)。
当 n=1 是,ϵ(1)=1=μ(1)。
设 n>1 有 k 个素因子,则 ∑d∣nμ(d)=k∑i=0(ki)(−1)i=(−1+1)k=0=ϵ(n)。
优越性:当 n=gcd(i,j) 时可以利用 ∑d∣gcd(i,j)=∑d∣i,d∣j 来消去 gcd 从而推柿子。
尝逝练习:P3768,loj 6629(都要杜教筛)。
一种类型 gcd 求和通法#
设 f 为一个定义域为 1∼n ,能 O(n) 预处理 O(1) 单点查询的函数。
求 n∑i=1n∑j=1f((i,j))。
n∑i=1n∑j=1f((i,j))=n∑d=1f(d)[n/d]∑i=1[n/d]∑j=1[(i,j)=1]=n∑d=1f(d)(2[n/d]∑i=1φ(n)−1)。
等式的推导可以看P2398的第二篇题解。
O(n),如果 f 能快速求前缀和,设求前缀和的复杂度为 g(x),则总复杂度可以通过整除分块和杜教筛做到 O(g(x)n2/3),一般 g(x)=O(1)。
应用:P2398。
gcd 卷积 O(nlnn)#
有两个长度为 n 的数组 f,g,求数组 h 满足 hk=∑gcd(i,j)=kfigj 。对 998244353 取模。
由于这个数组难求,于是考虑求数组 H 满足 Hk=∑k∣gcd(i,j)figj=∑k∣ifi∑k∣jgj 。显然能 O(nlnn) 求出。
由于 Hk=∑k∣dhd 。容斥一下得:hk=∑k|dμ(d)Hd,这可以算是莫反的另一种形式,也可以理解为用 μ 代替 −1 的次方进行容斥。
于是这里也能 O(nlnn) 。总复杂度 O(nlnn) 。这可以称作 gcd 卷积。
lcm 卷积:其他条件同上,求 hk=∑lcm(i,j)=kfigj 。
考虑求数组 H 满足 Hk=∑lcm(i,j)∣kfigj=∑i∣kfi∑j∣kgj 。就类似上面了。
应用:板子,CF1043F(参考这篇题解),SP31893。
快速 gcd 卷积 O(nloglogn) 和狄利克雷卷积相关加速#
前置知识(非必要):FWT,最好了解思想。
狄利克雷前/后缀和,前/后缀差分#
详细的可以看这里(差分就是这里说的逆)。我只讲狄利克雷前缀和(因为其他几乎类似)。
注:由莫比乌斯反演,一个函数和 μ 做狄利克雷卷积即是对该函数做狄利克雷前缀差分。
板子题:P5495。
FWT 中是后面几位相同,把首位为 0 的贡献给首位为 1 的。
这里变成把 i 贡献到 ip 上,ip≤n,p 是当前枚举的素数。
注意一下实现的时候要把 a 复制一遍(和 FWT 一样)。注意:一定要注意下文代码中的枚举顺序,建议写的同时自己思考一下。
for(int i=1;i<=cnt;i++) for(int j=1;j*pr[i]<=n;j++) a[j*pr[i]]+=a[j];
快速 gcd 卷积#
和 1 类似的,但这里考虑快速变换 a→b:bk=∑k∣iak ,其中 1≤i,k≤n。
要知 a 求 b 和知 b 求 a。就是狄利克雷后缀和与狄利克雷后缀差分,看看博客里的板子即可。
大致代码:
inline void FGT(int *a,int n){for(int i=1;i<=cnt;i++) for(int j=n/pr[i];j>=1;j--) a[j]+=a[pr[i]*j];}
inline void IFGT(int *a,int n){for(int i=cnt;i;i--) for(int j=1;pr[i]*j<=n;j++) a[j]-=a[pr[i]*j];}
FGT(f,n);FGT(g,n);for(int i=1;i<=n;i++) f[i]*=g[i];IFGT(f,n);
lcm 卷积由 2 知道是用狄利克雷前缀和与狄利克雷前缀差分,也能一样求。
注意到在两个数组长度不一样的时候也是能一样做的。
由埃式筛的复杂度得出这是 O(nloglogn) 的。应用同 2。
自创:2×107≥n≥m,给定 a1,2,…,m,b1,2,…,n,c1,2,…,m,求 n∑i=1m∑j=1a(i,j)bicj 对 998244353 取模。
n∑i=1m∑j=1a(i,j)bicj=m∑d=1adn∑i=1m∑j=1[(i,j)=d]bicj。
后面就是 gcd 卷积的形式了。可以应用到P4449。
f 为积性快速狄利克雷卷积 O(nloglogn)#
设 f 是积性函数,g 是任意函数,求他们的狄利克雷卷积函数 h 满足:h=f∗g。其中 f,g,h 定义域为 1∼n。
设 p1,2,…,cnt 是 n 以内的所有素数,Fi:pαi→Z ,且 Fi(pαi)=f(pαi) 。
考虑讲素数幂依次卷上去:f∗g=g∗f=gcnt∏i=1Fi,其中连乘是狄利克雷卷积乘法。
直接对每个素数幂贡献即可,复杂度为 cnt∑i=1∑j≥1,pji≤n⌊n/pji⌋≤ncnt∑i=1∑j≥11/pji=ncnt∑i=11pi−1<ncnt∑i=11pi。
于是复杂度同埃式筛为 O(nloglogn)。
大致代码:
for(int i=1;i<=cnt;i++) for(int j=n/pr[i];j;j--)
for(LL k=pr[i];j*k<=n;k*=pr[i]) g[j*k]=(g[j*k]+1ll*g[j]*f[k])%mod;
应用:P7580。
f 为积性快速狄利克雷逆 O(nloglogn)#
设 f 是积性函数,求 g 满足 f∗g=ϵ。
同样是对素数幂处求逆,最后所有素数卷起来。素数幂求逆的过程你可以看做解方程,自己推一推就可以发现复杂度和卷积部分是一样的。
狄利克雷生成函数相关#
懒得写啦,见博客,有狄利克雷生成函数,狄利克雷下的 lnn 定义,狄利克雷求导积分,ln,exp,快速幂。
狄利克雷卷积相关的积性函数判定#
-
定理 1,若 f,g 均为积性函数,则它们的点积为积性函数,即令 h(n)=f(n)g(n),则 h 为积性函数。
-
定理 2,若 f,g 均为积性函数,则它们的狄利克雷卷积 f×g 为积性函数。
∀n,m 满足 (n,m)=1,(f∗g)(n)×(f∗g)(m)=∑d1∣nf(d1)g(n/d1)∑d2∣mf(d2)g(m/d2)=∑d1∣n,d2∣m(f(d1)f(d2))×(g(n/d1)g(m/d2))
由于 (n,m)=1,于是 (d1,d2)=(n/d1,m/d2)=1。于是:
(f∗g)(n)×(f∗g)(m)=∑d1∣n,d2∣mf(d1d2)×g(nm/(d1d2))=∑D∣nmf(D)g(nm/D)=(f∗g)(nm)
- 定理 3,若 f 均为积性函数,则它的狄利克雷前缀和,后缀和,前缀差分,后缀差分为积性函数。读者不难自证。
根据这几个定理我们能看出绝大多数复杂的式子是不是积性函数。
应用:P10117,P4464。
数组 gcd 和#
给定一个 n 项正整数数列 a,求 n∑i=1n∑j=1gcd(ai,aj) ,对 998244353 取模。n≤105,ai≤109 。
这不就是 gcd 卷积吗?但是 ai≤109 ,设值域为 V,则复杂度就是 VloglogV,TLE。
做如下变换:n∑i=1n∑j=1gcd(ai,aj)=n∑i=1n∑j=1∑d∣gcd(ai,aj)φ(d)=V∑d=1φ(d)(n∑i=1[d∣ai])2。
注意到每个 ai 只会对 (n∑i=1[d∣ai])2 造成 d(ai) 次贡献。定义 cntd 表示 a1,..,an 中有因子 d 的个数,则答案就是 V∑d=1φ(d)(cntd)2 。设 max{d(1),d(2),...,d(V)}=d(v) ,我们注意到不为 0 的 cnt 最多只有 n×d(v) 个。
找出来计算即可,用个 hash 映射。用 pollard-rho,4√V 筛出素因子,然后 O(d(v)) 用 dfs 求出所有因子,在求所有因子的时候顺便计算所有因子的 φ 值,用 φ 的积性即可。然后对于所有的因子算贡献即可,注意不要算重。
复杂度 O(n4√V+nd(v))。不过这里是考虑比较极端的情况,大部分时候 V≤105 或 ai 可以快速分解质因数(参考下面应用),于是一般复杂度能更优秀 。欢迎来爆踩。
给定一个 n 项正整数排列 a,求 n∑i=1n∑j=1(i,j)×gcd(ai,aj) ,对 109+7 取模。n≤105。
类似的,有 n∑i=1n∑j=1(i,j)×gcd(ai,aj)=n∑i=1n∑j=1∑k∣(i,j)φ(k)∑l∣(ai,aj)φ(l)=n∑k=1φ(k)n∑l=1φ(l)([n/k]∑i=1[l∣aik])2
枚举每个 k ,逐个增加 k 的倍数,维护所有 l 的 φ(l)([n/k]∑i=1[l∣aik])2 ,就是直接在总和上增减,还是不懂看代码。
复杂度是 O(n∑i=1d(i)2) 的。考虑:
n∑i=1d(i)2=n∑i=1d(i)d(i)=n∑i=1∑lcm(j,k)∣i1=n∑j=1n∑k=1[n/lcm(j,k)]≤n∑j=1n∑k=1n/lcm(j,k)=n∑d=1d[n/d]∑j=1[n/d]∑k=1n/ij[(i,j)=1]≤n∑d=1n/d[n/d]∑j=1[n/d]∑k=11/ij=n∑d=1n/dln(n/d)2≤nln3n
于是复杂度是 O(nln3n)。
特殊方程求解#
引子:uoj 62。建议先阅读题解,由于题解写的太好了我就不在博客里重复一遍了,写点自己的东西。
先放下参考代码。
code
#include<bits/stdc++.h>
#define LL long long
#define fr(x) freopen(#x".in","r",stdin);freopen(#x".out","w",stdout);
using namespace std;
const int N=3e5+5,mod=998244353;
int n,c,d,q,a[N],f[N],I[N];
inline int md(int x){return x>=mod?x-mod:x;}
inline void ad(int &x,int y){x=md(x+y);}
inline int ksm(int x,int p){int s=1;for(;p;(p&1)&&(s=1ll*s*x%mod),x=1ll*x*x%mod,p>>=1);return s;}
inline void IFGT(int *a,int n){for(int i=1;i<=n;i++) for(int j=i+i;j<=n;j+=i) ad(a[j],mod-a[i]);}
inline void sol()
{
for(int i=1;i<=n;i++)
{
if(!f[i]){if(a[i]) return cout<<"-1\n",void();}
else a[i]=1ll*a[i]*ksm(f[i],mod-2)%mod;
}
for(int i=n;i;i--) for(int j=i+i;j<=n;j+=i) ad(a[i],mod-a[j]);
for(int i=1;i<=n;i++) cout<<(1ll*a[i]*I[i]%mod)<<" ";cout<<"\n";
}
int main()
{
ios::sync_with_stdio(0);cin.tie(0);cout.tie(0);cin>>n>>c>>d>>q;
d%=(mod-1);c=(c-d+mod-1)%(mod-1);
for(int i=1;i<=n;i++) f[i]=ksm(i,c),I[i]=ksm(i,mod-1-d);IFGT(f,n);
while(q--)
{
for(int i=1;i<=n;i++) cin>>a[i],a[i]=1ll*a[i]*I[i]%mod;
IFGT(a,n);sol();
}
return 0;
}
考虑题解完成了这样一件事情。对于矩阵 Ti,j=f((i,j))g(i)h(j),给定向量 →b,求解了 T→a=→b 的解 →a。
下文欲探究矩阵 T 逆的形式与可逆的条件。
一般地,不妨设 f,g,h 均为正。暂时保证 T 可逆。
考虑 T−1→b=→a,于是 T−1 相当于上文反解出 a 各个运算的矩阵成起来,具体地:
-
首先我们对 f 莫反,找出了 F 满足 f(n)=∑d∣nF(d)。
-
然后,我们对 →b 的每个位置点乘了 g−1,此时设 Ai,j=g(i)−1,即进行了运算 A→b。
-
接着,我们对此时的 →b 做了莫反,由博客知乘上了莫反矩阵 Bi,j=[j∣i]μ(i/j)。
-
而后,我们对 →b 的每个位置点乘了 F−1,即乘上 C:Ci,j=F(i)−1。
-
此时,fd=n∑j=1[d∣j]g(j),要知 f 求 g。此时由莫反的转置形式,即乘上矩阵 Di,j=[i∣j]μ(j/i)。
-
最后,对 →b 点成 h−1,即乘上 E:Ei,j=h(i)−1。
考虑反解的过程,有:T−1→b=E⋅D⋅C⋅B⋅A→b。可能能把这段乘积整理成更优美的性质,笔者推不动了。
由题设,A,B,D,E 矩阵均存在,只有 F 可能有某些项为 0。
有线性代数结论:∏Ai 若除了 Ak 均可逆,则其矩阵乘积可逆当且仅当 Ak 可逆。
于是 T 可逆等价于 F 的每项均非 0,就找到了可逆条件。
但注意:解方程不一定要满足 T 可逆,具体参考题解那个博客中的无解情况判断。
von Mangoldt function#
定义 von Mangoldt 函数 Λ(n)={lnp(n=pk,p∈Prime)0(otherwise)。
性质 1:lnn=∑d∣nΛ(d),读者不难自证。
性质 2:对性质 1 进行莫比乌斯反演,有了新的定义式:Λ(n)=∑d∣nμ(n/d)lnd。
结论 1:对性质 2 两边取 exp 得:∏d∣ndμ(n/d)={p(n=pk,p∈Prime)1(otherwise)。
应用:P7360。
一种特殊情形下的 PN 筛优化积性函数求和#
设 f(n) 为积性函数,且 ∀p∈Prime,f(p)=1,且 f(pk) 能 O(1) 求得。试快速求 S(n)=n∑i=1f(i)。
类似 PN 筛 的,考虑构造 g=f∗μ⇒g∗1=f,此时有:S(n)=n∑i=1∑d∣ig(d)=n∑d=1g(d)⌊n/d⌋。
并且我们注意到:g(pk)=f(pk)−f(pk−1),于是 g(n) 有值当且仅当 n 没有次数为 1 的素因子,即 n 为 PN 数。
注意到 1∼n 的 PN 数只有 O(√n) 个,于是我们枚举所有 PN 数,计算其 g(d)⌊n/d⌋ 加起来即可。
复杂度 O(√n),远优于直接上任何一个筛法!
应用:填数。
注:经研究,对于平凡的 f(p)=c(c≠1) 的情况,取 U(n)=μ(n)cω(n),g=f∗U⇒g∗U−1=f。
容易解得 U−1(n)=cΩ(n),其中 Ω(n) 表示 n 的可重素因子个数。
此时 S(n)=n∑d=1g(d)⌊n/d⌋∑i=1cΩ(i)。
观察到后半部分那个求和无法快速求得,故平凡(c≠1)情况下,无法用此方法进行复杂度优化。
【推荐】国内首个AI IDE,深度理解中文开发场景,立即下载体验Trae
【推荐】编程新体验,更懂你的AI,立即体验豆包MarsCode编程助手
【推荐】抖音旗下AI助手豆包,你的智能百科全书,全免费不限次数
【推荐】轻量又高性能的 SSH 工具 IShell:AI 加持,快人一步
· 物流快递公司核心技术能力-地址解析分单基础技术分享
· 单线程的Redis速度为什么快?
· 展开说说关于C#中ORM框架的用法!
· Pantheons:用 TypeScript 打造主流大模型对话的一站式集成库
· SQL Server 2025 AI相关能力初探