数论柿子题的发现

~~鉴于博客园惨淡现状，搬点自己博客的文章~~

写在前面：下面除法默认下取整，为了方便不写分数，没有说明 $(i,j)$ 都表示取 $\gcd$。

希望下面的每一题都能让你们学到一个 $\texttt{trick}$，本文不适合入门，建议先学习一些推式子知识以及筛发（杜教筛，$\min25$ 筛）再来看毒瘤的推柿子。

两个前置知识（推柿子通法）

定义 $\epsilon(n)=[n=1]$。我们有：$n=\sum \limits_{d\mid n}\varphi(d)$ 和 $\epsilon(n)=\sum\limits_{d\mid n}\mu(d)$ 。

对于每个数 $x$ 它和 $n$ 的 $\gcd$ 是一定的。于是有 $n=\sum\limits_{i=1}^n\sum\limits_{j=1}^n[(j,n)=i]=\sum\limits_{i\mid n}\sum\limits_{i\mid j}[(j/i,n/i)=1]=\sum\limits_{i\mid n}\sum\limits_{j=1}^{n/i}[(j,n/i)=1]=\sum\limits_{i\mid n}\varphi(n/i)=\sum\limits_{i\mid n}\varphi(i)$。

当 $n=1$ 是，$\epsilon(1)=1=\mu(1)$。

设 $n>1$ 有 $k$ 个素因子，则 $\sum\limits_{d\mid n}\mu(d)=\sum\limits_{i=0}^k \dbinom{k}{i}(-1)^i=(-1+1)^k=0=\epsilon(n)$。

优越性：当 $n=\gcd(i,j)$ 时可以利用 $\sum\limits_{d\mid \gcd(i,j)}=\sum\limits_{d\mid i,d\mid j}$ 来消去 $\gcd$ 从而推柿子。

尝逝练习：P3768，loj 6629（都要杜教筛）。

一种类型 $\gcd$ 求和通法

设 $f$ 为一个定义域为 $1\sim n$ ，能 $O(n)$ 预处理 $O(1)$ 单点查询的函数。

求 $\sum\limits_{i=1}^n\sum\limits_{j=1}^n f((i,j))$。

$\sum\limits_{i=1}^n\sum\limits_{j=1}^n f((i,j))=\sum\limits_{d=1}^nf(d)\sum\limits_{i=1}^{[n/d]}\sum\limits_{j=1}^{[n/d]} [(i,j)=1]=\sum\limits_{d=1}^nf(d)(2\sum\limits_{i=1}^{[n/d]}\varphi(n)-1)$。

等式的推导可以看P2398的第二篇题解。

$O(n)$，如果 $f$ 能快速求前缀和，设求前缀和的复杂度为 $g(x)$，则总复杂度可以通过整除分块和杜教筛做到 $O(g(x)n^{2/3})$，一般 $g(x)=O(1)$。

应用：P2398。

$\gcd$ 卷积 $O(n\ln n)$

有两个长度为 $n$ 的数组 $f,g$，求数组 $h$ 满足 $h_k=\sum\limits_{\gcd(i,j)=k} f_ig_j$ 。对 $998244353$ 取模。

由于这个数组难求，于是考虑求数组 $H$ 满足 $H_k=\sum\limits_{k\mid\gcd(i,j)} f_ig_j=\sum\limits_{k\mid i}f_i\sum\limits_{k\mid j} g_j$ 。显然能 $O(n\ln n)$ 求出。

由于 $H_k=\sum\limits_{k\mid d} h_d$ 。容斥一下得：$h_k=\sum\limits_{k|d} \mu(d)H_d$，这可以算是莫反的另一种形式，也可以理解为用 $\mu$ 代替 $-1$ 的次方进行容斥。

于是这里也能 $O(n\ln n)$ 。总复杂度 $O(n\ln n)$ 。这可以称作 $\gcd$ 卷积。

$\text{lcm}$ 卷积：其他条件同上，求 $h_k=\sum\limits_{\text{lcm}(i,j)=k} f_ig_j$ 。

考虑求数组 $H$ 满足 $H_k=\sum\limits_{\text{lcm}(i,j)\mid k} f_ig_j=\sum\limits_{i\mid k}f_i\sum\limits_{j\mid k} g_j$ 。就类似上面了。

应用：板子，CF1043F（参考这篇题解），SP31893。

快速 $\gcd$ 卷积 $O(n\log\log n)$ 和狄利克雷卷积相关加速

前置知识（非必要）：$\text{FWT}$，最好了解思想。

狄利克雷前/后缀和，前/后缀差分

详细的可以看这里（差分就是这里说的逆）。我只讲狄利克雷前缀和（因为其他几乎类似）。

注：由莫比乌斯反演，一个函数和 $\mu$ 做狄利克雷卷积即是对该函数做狄利克雷前缀差分。

板子题：P5495。

$\text{FWT}$ 中是后面几位相同，把首位为 $0$ 的贡献给首位为 $1$ 的。

这里变成把 $i$ 贡献到 $ip$ 上，$ip\le n$，$p$ 是当前枚举的素数。

注意一下实现的时候要把 $a$ 复制一遍（和 $\text{FWT}$ 一样）。注意：一定要注意下文代码中的枚举顺序，建议写的同时自己思考一下。

//线性筛素数,cnt是素数个数,pr是素数数组
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\to b:b_k=\sum\limits_{k\mid i} a_k$ ，其中 $1\le i,k\le n$。

要知 $a$ 求 $b$ 和知 $b$ 求 $a$。就是狄利克雷后缀和与狄利克雷后缀差分，看看博客里的板子即可。

大致代码：

//线性筛素数,cnt为素数个数,pr是素数数组
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);

$\text{lcm}$ 卷积由 $2$ 知道是用狄利克雷前缀和与狄利克雷前缀差分，也能一样求。

注意到在两个数组长度不一样的时候也是能一样做的。

由埃式筛的复杂度得出这是 $O(n\log\log n)$ 的。应用同 $2$。

自创：$2\times 10^7\ge n\ge m$，给定 $a_{1,2,\dots,m},b_{1,2,\dots,n},c_{1,2,\dots,m}$，求 $\sum\limits_{i=1}^n\sum\limits_{j=1}^m a_{(i,j)}b_ic_j$ 对 $998244353$ 取模。

$\sum\limits_{i=1}^n\sum\limits_{j=1}^m a_{(i,j)}b_ic_j=\sum\limits_{d=1}^ma_d\sum\limits_{i=1}^n\sum\limits_{j=1}^m [(i,j)=d]b_ic_j$。

后面就是 $\gcd$ 卷积的形式了。可以应用到P4449。

$f$ 为积性快速狄利克雷卷积 $O(n\log\log n)$

设 $f$ 是积性函数，$g$ 是任意函数，求他们的狄利克雷卷积函数 $h$ 满足：$h=f*g$。其中 $f,g,h$ 定义域为 $1\sim n$。

设 $p_{1,2,\dots,cnt}$ 是 $n$ 以内的所有素数，$F_i:p_i^\alpha\to \mathbb{Z}$ ，且 $F_i(p_i^\alpha)=f(p_i^\alpha)$ 。

考虑讲素数幂依次卷上去：$f*g=g*f=g\prod\limits_{i=1}^{cnt} F_i$，其中连乘是狄利克雷卷积乘法。

直接对每个素数幂贡献即可，复杂度为 $\sum\limits_{i=1}^{cnt}\sum\limits_{j\ge 1,p_i^j\le n} \lfloor n/p_i^j \rfloor\le n\sum\limits_{i=1}^{cnt}\sum\limits_{j\ge 1} 1/p_i^j=n\sum\limits_{i=1}^{cnt}\dfrac{1}{p_i-1}< n\sum\limits_{i=1}^{cnt}\dfrac{1}{p_i}$。

于是复杂度同埃式筛为 $O(n\log\log n)$。

大致代码：

//线性筛素数,cnt为素数个数,pr是素数数组
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(n\log\log n)$

设 $f$ 是积性函数，求 $g$ 满足 $f*g=\epsilon$。

同样是对素数幂处求逆，最后所有素数卷起来。素数幂求逆的过程你可以看做解方程，自己推一推就可以发现复杂度和卷积部分是一样的。

狄利克雷生成函数相关

懒得写啦，见博客，有狄利克雷生成函数，狄利克雷下的 $\ln n$ 定义，狄利克雷求导积分，$\ln,\exp$，快速幂。

狄利克雷卷积相关的积性函数判定

定理 $1$，若 $f,g$ 均为积性函数，则它们的点积为积性函数，即令 $h(n)=f(n)g(n)$，则 $h$ 为积性函数。
定理 $2$，若 $f,g$ 均为积性函数，则它们的狄利克雷卷积 $f\times g$ 为积性函数。

$\forall n,m$ 满足 $(n,m)=1$，$(f*g)(n)\times (f*g)(m)=\sum\limits_{d_1\mid n} f(d_1)g(n/d_1)\sum\limits_{d_2\mid m} f(d_2)g(m/d_2)=\sum\limits_{d_1\mid n,d_2\mid m} (f(d_1)f(d_2))\times (g(n/d_1)g(m/d_2))$

由于 $(n,m)=1$，于是 $(d_1,d_2)=(n/d_1,m/d_2)=1$。于是：

$(f*g)(n)\times (f*g)(m)=\sum\limits_{d_1\mid n,d_2\mid m} f(d_1d_2)\times g(nm/(d_1d_2))=\sum\limits_{D\mid nm} f(D)g(nm/D)=(f*g)(nm)$

定理 $3$，若 $f$ 均为积性函数，则它的狄利克雷前缀和，后缀和，前缀差分，后缀差分为积性函数。读者不难自证。

根据这几个定理我们能看出绝大多数复杂的式子是不是积性函数。

应用：P10117，P4464。

数组 $\gcd$ 和

给定一个 $n$ 项正整数数列 $a$，求 $\sum\limits_{i=1}^n\sum\limits_{j=1}^n \gcd(a_i,a_j)$ ，对 $998244353$ 取模。$n\le 10^5,a_i\le 10^9$ 。

这不就是 $\gcd$ 卷积吗？但是 $a_i\le 10^9$ ，设值域为 $V$，则复杂度就是 $V\log\log V$，$\texttt{TLE}$。

做如下变换：$\sum\limits_{i=1}^n\sum\limits_{j=1}^n \gcd(a_i,a_j)=\sum\limits_{i=1}^n\sum\limits_{j=1}^n\sum\limits_{d\mid\gcd(a_i,a_j)} \varphi(d)=\sum\limits_{d=1}^V \varphi(d)(\sum\limits_{i=1}^n\left[d\mid a_i\right])^2$。

注意到每个 $a_i$ 只会对 $(\sum\limits_{i=1}^n\left[d\mid a_i\right])^2$ 造成 $d(a_i)$ 次贡献。定义 $cnt_d$ 表示 $a_1,..,a_n$ 中有因子 $d$ 的个数，则答案就是 $\sum\limits_{d=1}^V \varphi(d)(cnt_d)^2$ 。设 $\max\{d(1),d(2),...,d(V)\}=d(v)$ ，我们注意到不为 $0$ 的 $cnt$ 最多只有 $n\times d(v)$ 个。

找出来计算即可，用个 $\text{hash}$ 映射。用 $\text{pollard-rho}$，$\sqrt[4]{V}$ 筛出素因子，然后 $O(d(v))$ 用 $\texttt{dfs}$ 求出所有因子，在求所有因子的时候顺便计算所有因子的 $\varphi$ 值，用 $\varphi$ 的积性即可。然后对于所有的因子算贡献即可，注意不要算重。

复杂度 $O(n\sqrt[4]{V}+nd(v))$。不过这里是考虑比较极端的情况，大部分时候 $V\le 10^5$ 或 $a_i$ 可以快速分解质因数（参考下面应用），于是一般复杂度能更优秀。欢迎来爆踩。

变式：loj 6539，CF1575G

给定一个 $n$ 项正整数排列 $a$，求 $\sum\limits_{i=1}^n\sum\limits_{j=1}^n (i,j)\times \gcd(a_i,a_j)$ ，对 $10^9+7$ 取模。$n\le 10^5$。

类似的，有 $\sum\limits_{i=1}^n\sum\limits_{j=1}^n (i,j)\times \gcd(a_i,a_j)=\sum\limits_{i=1}^n\sum\limits_{j=1}^n \sum\limits_{k\mid (i,j) }\varphi(k)\sum\limits_{l\mid (a_i,a_j) }\varphi(l)=\sum\limits_{k=1}^n\varphi(k)\sum\limits_{l=1}^n\varphi(l)(\sum \limits_{i=1}^{[n/k]}[l\mid a_{ik}])^2$

枚举每个 $k$ ，逐个增加 $k$ 的倍数，维护所有 $l$ 的 $\varphi(l)(\sum \limits_{i=1}^{[n/k]}[l\mid a_{ik}])^2$ ，就是直接在总和上增减，还是不懂看代码。

复杂度是 $O(\sum\limits_{i=1}^n d(i)^2)$ 的。考虑：

$\sum\limits_{i=1}^n d(i)^2=\sum\limits_{i=1}^n d(i)d(i)=\sum\limits_{i=1}^n \sum\limits_{\text{lcm}(j,k)\mid i}1=\sum\limits_{j=1}^n \sum\limits_{k=1}^n [n/\text{lcm}(j,k)]\le\sum\limits_{j=1}^n \sum\limits_{k=1}^n n/\text{lcm}(j,k)=\sum\limits_{d=1}^nd\sum\limits_{j=1}^{[n/d]} \sum\limits_{k=1}^{[n/d]} n/ij[(i,j)=1]\le \sum\limits_{d=1}^nn/d\sum\limits_{j=1}^{[n/d]} \sum\limits_{k=1}^{[n/d]} 1/ij=\sum\limits_{d=1}^nn/d\ln (n/d)^2\le n\ln ^3n$

于是复杂度是 $O(n\ln ^3n)$。

特殊方程求解

引子：uoj 62。建议先阅读题解，由于题解写的太好了我就不在博客里重复一遍了，写点自己的东西。

先放下参考代码。

//uoj 62
//https://uoj.ac/problem/62
#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;
}

考虑题解完成了这样一件事情。对于矩阵 $T_{i,j}=f((i,j))g(i)h(j)$，给定向量 $\vec{b}$，求解了 $T\vec{a}=\vec{b}$ 的解 $\vec{a}$。

下文欲探究矩阵 $T$ 逆的形式与可逆的条件。

一般地，不妨设 $f,g,h$ 均为正。暂时保证 $T$ 可逆。

考虑 $T^{-1}\vec{b}=\vec{a}$，于是 $T^{-1}$ 相当于上文反解出 $a$ 各个运算的矩阵成起来，具体地：

首先我们对 $f$ 莫反，找出了 $F$ 满足 $f(n)=\sum\limits_{d\mid n}F(d)$。
然后，我们对 $\vec{b}$ 的每个位置点乘了 $g^{-1}$，此时设 $A_{i,j}=g(i)^{-1}$，即进行了运算 $A\vec{b}$。
接着，我们对此时的 $\vec{b}$ 做了莫反，由博客知乘上了莫反矩阵 $B_{i,j}=[j\mid i]\mu(i/j)$。
而后，我们对 $\vec{b}$ 的每个位置点乘了 $F^{-1}$，即乘上 $C:C_{i,j}=F(i)^{-1}$。
此时，$f_d = \sum\limits_{j=1}^{n} [d \mid j] g(j)$，要知 $f$ 求 $g$。此时由莫反的转置形式，即乘上矩阵 $D_{i,j}=[i\mid j]\mu(j/i)$。
最后，对 $\vec{b}$ 点成 $h^{-1}$，即乘上 $E:E_{i,j}=h(i)^{-1}$。

考虑反解的过程，有：$T^{-1}\vec{b}=E\cdot D\cdot C\cdot B\cdot A\vec{b}$。可能能把这段乘积整理成更优美的性质，笔者推不动了。

由题设，$A,B,D,E$ 矩阵均存在，只有 $F$ 可能有某些项为 $0$。

有线性代数结论：$\prod A_{i}$ 若除了 $A_{k}$ 均可逆，则其矩阵乘积可逆当且仅当 $A_k$ 可逆。

于是 $T$ 可逆等价于 $F$ 的每项均非 $0$，就找到了可逆条件。

但注意：解方程不一定要满足 $T$ 可逆，具体参考题解那个博客中的无解情况判断。

von Mangoldt function

定义 von Mangoldt 函数 $\Lambda(n)=\begin{cases}\ln p(n=p^k,p\in\text{Prime})\\0(\texttt{otherwise})\end{cases}$。

性质 $1$：$\ln n=\sum\limits_{d\mid n} \Lambda(d)$，读者不难自证。

性质 $2$：对性质 $1$ 进行莫比乌斯反演，有了新的定义式：$\Lambda(n)=\sum\limits_{d\mid n} \mu(n/d)\ln d$。

结论 $1$：对性质 $2$ 两边取 $\exp$ 得：$\prod\limits_{d\mid n}d^{\mu(n/d)}=\begin{cases} p(n=p^k,p\in\text{Prime})\\1(\texttt{otherwise})\end{cases}$。

应用：P7360。

一种特殊情形下的 PN 筛优化积性函数求和

设 $f(n)$ 为积性函数，且 $\forall p\in \text{Prime},f(p)=1$，且 $f(p^k)$ 能 $O(1)$ 求得。试快速求 $S(n)=\sum\limits_{i=1}^n f(i)$。

类似 PN 筛 的，考虑构造 $g=f*\mu \Rightarrow g*1=f$，此时有：$S(n)=\sum\limits_{i=1}^n \sum\limits_{d\mid i}g(d)=\sum\limits_{d=1}^n g(d)\lfloor n/d\rfloor$。

并且我们注意到：$g(p^k)=f(p^k)-f(p^{k-1})$，于是 $g(n)$ 有值当且仅当 $n$ 没有次数为 $1$ 的素因子，即 $n$ 为 PN 数。

注意到 $1\sim n$ 的 PN 数只有 $O(\sqrt n)$ 个，于是我们枚举所有 PN 数，计算其 $g(d)\lfloor n/d\rfloor$ 加起来即可。

复杂度 $O(\sqrt n)$，远优于直接上任何一个筛法！

注：经研究，对于平凡的 $f(p)=c(c\neq 1)$ 的情况，取 $U(n)=\mu(n)c^{\omega(n)},g=f*U\Rightarrow g*U^{-1}=f$。

容易解得 $U^{-1}(n)=c^{\Omega(n)}$，其中 $\Omega(n)$ 表示 $n$ 的可重素因子个数。

此时 $S(n)=\sum\limits_{d=1}^n g(d)\sum\limits_{i=1}^{\lfloor n/d\rfloor} c^{\Omega(i)}$。

观察到后半部分那个求和无法快速求得，故平凡（$c\neq 1$）情况下，无法用此方法进行复杂度优化。

masterhuang的博客

博客