数学的恶魔

莫比乌斯反演

前置知识：数论分块

基本内容

数论分块可以快速计算一些含有除法向下取整的和式（即形如 $\sum _{i=1}^{n}f(i)g(\lfloor \frac{n}{i}\rfloor )$ 的式子），只要可以在 $O(1)$ 内计算 $f(r)-f(l)$ 或已经预处理出 $f$ 的前缀和时，可以用数论分块在 $O(\sqrt{n})$ 的时间内计算答案。

它主要利用富比尼定理（Fubini’s theorem）（听起来很厉害但其实没啥）：我们注意到，有很多的 $i$ 它们的 $\lfloor \frac{n}{i}\rfloor$ 是相等的，这启发我们分段打包计算。

举个例子，如果我们要求 $\sum _{i=1}^n\lfloor \frac{n}{i}\rfloor$ ，我们发现很多 $i$ 都对应同一个 $\lfloor \frac{n}{i}\rfloor$ ，它们呈块状分布，在打个表仔细观察（证明），发现对于每一个块，它的最后一个数是 $n/(n/i)$ （ $c++$ 自带下取整，一般我们称其这个式子的值为 $g(i)$ ），如对于 $\sum _{i=1}^{10}\lfloor \frac{10}{i}\rfloor$ ，分成： $(1)(2)(3)(4,5)(6,7,8,9,10)$ ，其中 $1=10/(10/1),2=10/(10/2),3=10/(10/3),5=10/(10/5),10=10/(10/10)$ ，于是，求解上述问题的代码如下（常用）：

UVA11526 H(n)

for(int l=1,r;l<=n;l=r+1){
    r=n/(n/l);
    ans+=(r-l+1)*(n/l);
}

类似的，我们可以用 $O(\sqrt{n})$ 计算 $g(\lfloor \frac{n}{i}\rfloor )$ 那么只要 $f$ 好求我们就能算了！

一般而言，我们可以统计一个前缀和（积），因为每当我们使用分块跳过一个区间的时候，其所对应的函数值也跳过了一个区间。所以此时，就需要乘上那一个区间的函数值（可以自己举个例子理解一下）。

但要是不一般（比如 $O(n)$ 的前缀和都会TLE），那……去死吧！（dalao可以用杜教筛试一试）

引理一

$$\mathrm{\forall }a,b,c\in \mathbb{Z},\lfloor \frac{a}{bc}\rfloor =\lfloor \frac{\lfloor \frac{a}{b}\rfloor }{c}\rfloor$$

证明：
$$\frac{a}{b}=\lfloor \frac{a}{b}\rfloor +r(0\le r<1)\Rightarrow \lfloor \frac{a}{bc}\rfloor =\lfloor \frac{1}{c}(\lfloor \frac{a}{b}\rfloor +r)\rfloor =\lfloor \frac{\lfloor \frac{a}{b}\rfloor }{c}+\frac{r}{c}\rfloor =\lfloor \frac{\lfloor \frac{a}{b}\rfloor }{c}\rfloor$$

引理二

$$\mathrm{\forall }n\in {\mathbb{N}}_{\mathbb{+}},|\lfloor \frac{n}{d}\rfloor |d\in {\mathbb{N}}_{\mathbb{+}}\wedge d\le n|\le \lfloor 2\sqrt{n}\rfloor$$

其中 $|V|$ 表示集合 $V$ 的元素个数

该引理是时间复杂度的保证

证明：

对于 $d\le \lfloor \sqrt{n}\rfloor$ ，由于 $d$ 最多有 $\lfloor \sqrt{n}\rfloor$ 种不同取值， $\lfloor \frac{n}{d}\rfloor$ 也只有 $\lfloor \sqrt{n}\rfloor$ 种取值

对于 $d>\lfloor \sqrt{n}\rfloor$ ，则 $\lfloor \frac{n}{d}\rfloor \le \lfloor \sqrt{n}\rfloor$ ，也只有 $\lfloor \sqrt{n}\rfloor$ 种取值

前置知识：积性函数

若函数 $f(x)$ 为积性函数，则有 $\mathrm{\forall }gcd(a,b)=1,f(ab)=f(a)f(b)$

可以证明以下函数全是积性函数：

• $\phi (n)$ －欧拉函数
• $\mu (n)$ －莫比乌斯函数
• $gcd(n,k)$ －最大公因子，当k固定的情况
• $d(n)$ －正因子数目
• $\sigma (n)$ －所有正因子之和
• $\sigma k(n)$ －因子函数， $n$ 的所有正因子的 $k$ 次幂之和，当中 $k$ 可为任何复数
• $1(n)$ －不变的函数，定义为 $1(n)=1$ （完全积性）
• $Id(n)$ －单位函数，定义为 $Id(n)=n$ （完全积性）
• $Idk(n)$ －幂函数，对于任何复数、实数 $k$ ，定义为 $Idk(n)=n^k$ （完全积性）
• $\epsilon (n)$ －定义为：若 $n=1$ ， $\epsilon (n)=1$ ；若 $n>1$ ， $\epsilon (n)=0$ 。别称为“对于狄利克雷卷积的乘法单位”（完全积性）
• $\lambda (n)$ －刘维尔函数，关于能整除 $n$ 的质因子的数目
• $\gamma (n)$ ，定义为 $\gamma (n)=(-1{)}^{\omega (n)}$，在此加性函数 $\omega (n)$ 是不同能整除n的质数的数目
• 另外，所有狄利克雷特征均是完全积性的

积性函数一般可以用线性筛筛得

前置知识：狄利克雷卷积

定义

狄利克雷（ $Dirichlet$ ）卷积定义为：

对于两个数论函数 $f,g$ ，有：
$$(f\ast g)(n)=\sum _{d|n}f(d)g(\frac{n}{d})$$

性质

狄利克雷卷积满足以下运算律：

• 交换律： $(f\ast g)=(g\ast f)$
• 结合律： $((f\ast g)\ast h)=(f\ast (g\ast h))$
• 分配律： $(f\ast (g+h))=((f\ast g)+(f\ast h))$
• $(f\ast \epsilon )=f$ ，其中 $\epsilon$ 为狄利克雷卷积的单位元

for example：

• $\epsilon =(\mu \ast 1)\iff \epsilon (n)=\sum _{d|n}\mu (d)1(\frac{n}{d})=\sum _{d|n}\mu (d)$
• $d=(1\ast 1)\iff d(n)=\sum _{d|n}1(d)1(\frac{n}{d})=\sum _{d|n}1$
• $Id=(\phi \ast 1)\iff Id(n)=n=\sum _{d|n}\phi (n)$
• $\phi =(\mu \ast Id)\iff \phi (n)=\sum _{d|n}d\mu (\frac{n}{d})$

求法

对于卷积的第 $n$ 项可以直接枚举约数，在 $O(\sqrt{n})$ 的时间内计算

但是对于狄利克雷卷积的前 $n$ 项
$$h[i]=\sum _{d|i}f[d]g[i/d]$$
如果一项一项算的话就需要 $O(n\sqrt{n})$ ，但实际上可以优化

设 $x=d,y=\frac{i}{d}$ 分别枚举 $x,y$ 对于 $h[xy]+=f[x]g[y]$ 即可。
时间复杂度 $O(n\log n)$

for(int i=1;i<=n;++i) h[i]=0;
for(int i=1;i<=n;++i)
    for(int j=1;j*i<=n;++j) h[i*j]=(h[i*j]+f[i]*g[j]%P)%P;

对于 $k$ 次卷积，还可以快速幂优化

struct HanShu{
    int f[N];
    void inint(){
    	for(int i=1;i<=n;++i) f[i]=0;
    }
    HanShu operator*(const HanShu _t){
	    HanShu h;
		h.inint();
	    for(int i=1;i<=n;++i)
	        for(int j=1;j*i<=n;++j) h.f[i*j]=(h.f[i*j]+f[i]*_t.f[j]%P)%P;
	    return h;
	}
	void operator=(const HanShu _t){
	    for(int i=1;i<=n;++i) f[i]=_t.f[i];
	}
};
HanShu H_qpow(HanShu x,int y){
    HanShu res;
    int _t=0;
    while(y){
    	if(y&1){
    		++_t;
    		if(_t==1) res=x;
    		else res=res*x;
    	}
    	x=x*x;
    	y>>=1;
    }
    return res;
}

前置知识：莫比乌斯函数

定义

设正整数 $N$ 分解质因数为 $N=p_1^{c_1}{p}_2^{c_2}\dots p_m^{c_m}$ ，定义函数
$$\mu (N)=\{\begin{array}{ll}0& \mathrm{\exists }i\in [1,m],c_i>1\\ 1& m\equiv 0(mod2),\mathrm{\forall }i\in [1,m],c_i=1\\ -1& m\equiv 1(mod2),\mathrm{\forall }i\in [1,m],c_i=1\end{array}$$
我们称 $\mu$ 为莫比乌斯（ $M\ddot{o}bius$ ）函数。

性质和结论

• 莫比乌斯函数是积性函数

• $\sum _{d|n}\mu (d)=[n==1]$ ，其中 $[p]$ 表示 $p$ 为真时取 $1$ 其他时候取 $0$ ，该式也可表示为 $\sum _{d|n}\mu (d)=\epsilon (n)$ ，也就是 $\mu \ast 1=\epsilon$

  证明：

  设 $n=\prod _{i=1}^k{p}_i^{c_i},n^{\prime }=\prod _{i=1}^k{p}_i$

  那么 $\sum _{d|n}\mu (d)=\sum _{d|n^{\prime }}\mu (d)=\sum _{i=0}^k{C}_k^i(-1{)}^i$

  又因为二项式定理： $(a+b{)}^k=\sum _{i=0}^k{C}_k^i{a}^i{b}^{k-i}$

  有： $\sum _{i=0}^k{C}_k^i(-1{)}^i=\sum _{i=0}^k{C}_k^i(-1{)}^i{1}^{k-i}=((-1)+1{)}^k$

  当且仅当 $k=0$ ，即 $n=1$ 时原式为 $1$ ，其余时候都为 $0$

• $[gcd(i,j)==1]=\sum _{d|gcd(i,j)}\mu (d)$ 该结论在反演中十分常用，正确性显然

筛法

若只求一项莫比乌斯函数，分解质因数即可；

若要求多项可以在线性筛素数时一并求出：

int primes[N],mu[N],cnt=0;
bool v[N];
mu[1]=1; //1特判
for(int i=2;i<=n;++i){
    if(!v[i]) primes[++cnt]=i,mu[i]=-1; //质数
    for(int j=1;j<=cnt&&primes[j]*i<=n;++j){
		v[primes[j]*i]=true;
        if(i%primes[j]==0){ //若有一个质因数指数大于1
            mu[i*primes[j]]=0;
            break;
        }
        mu[i*primes[j]]=-mu[i];
    }
}

正题：莫比乌斯反演

定理

若 $F(n)$ 和 $f(n)$ 是定义在 $\mathbb{N}$ 上的函数，且有
$$F(n)=\sum _{d|n}f(d)$$
则：
$$f(n)=\sum _{d|n}\mu (d)F(\frac{n}{d})$$
我们称该定理为莫比乌斯反演定理（原来是用 $f$ 表示 $F$ ，现在反过来了）。

证明如下：
$$\begin{array}{}\text{(1)}& \sum _{d|n}\mu (d)F(\frac{n}{d})& =\sum _{d|n}\mu (d)\sum _{i|\frac{n}{d}}f(i)\text{(2)}& & =\sum _{i|n}f(i)\sum _{d|\frac{n}{i}}\mu (d)\text{(3)}& & =\sum _{i|n}f(i)[\frac{n}{i}==1]\text{(4)}& & =f(n)\end{array}$$
当然，不难发现，狄利克雷卷积也可以证明：

$$\begin{array}{}\text{(5)}& & F(n)=\sum _{d|n}f(d)\text{(6)}& & \Rightarrow F=f\ast 1\text{(7)}& & \Rightarrow F\ast \mu =f\ast 1\ast \mu \text{(8)}& \text{又}\because & \mu \ast 1=\epsilon \text{(9)}& & \Rightarrow F\ast \mu =f\ast (1\ast \mu )=f\text{(10)}& & \Rightarrow f=\mu \ast F\end{array}$$

然而，在莫比乌斯反演中更常用的是这个式子：

当 $F(n)$ 和 $f(n)$ 满足：
$$F(n)=\sum _{n|d}f(d)$$
则有：
$$f(n)=\sum _{n|d}\mu (\frac{d}{n})F(d)$$
其实就是枚举约数变成了枚举倍数，若用卷积证明，过程都一样，下面用定义证明：
$$\begin{array}{}\text{(11)}& \sum _{n|d}\mu (\frac{d}{n})F(d)& =\sum _{n|d}\mu (\frac{d}{n})\sum _{d|i}f(i)\text{(12)}& \text{不妨设}{d}^{\prime }=\frac{d}{n}\text{(13)}& \text{原式}& =\sum _{n|i}f(i)\sum _{d^{\prime }|\frac{i}{n}}\mu (d^{\prime })\text{(14)}& & =\sum _{n|i}f(i)[\frac{i}{n}==1]\text{(15)}& & =f(n)\end{array}$$
以上就是莫比乌斯反演的常用形式（要记住）

例题一

Problem b

分析：

若我们将数对 $(x,y)$ 对应到坐标系内，其实是要求在左下角为 $(a,c)$ 右上角为 $(b,d)$ 的矩形中有多少个点 $(x,y)$ 满足 $gcd(x,y)==k$

由矩形考虑二维前缀和，设 $S(x,y)$ 表示以 $(x,y)$ 为右上角， $(0,0)$ 为左下角的矩形中满足要求的点的个数，则 $Ans=S(b,d)-S(a-1,d)-S(b,c-1)+S(a-1,c-1)$

考虑用莫比乌斯反演，我们需要定义 $F,f$ 且让 $F(n)=\sum _{n|d}f(d)$ （即满足反演条件），要注意的是， $F$ 应该是比较好计算的函数（至少要比 $f$ 好计算，不然你反演用 $F$ 表示 $f$ 干嘛），还有就是 $f$ 要与答案有关系（不然求出 $f$ 干嘛）

那么，不妨令 $f(n)=\sum _{x=1}^N\sum _{y=1}^M[gcd(x,y)==n]$ 则我们要求的 $S(N,M)$ 就是 $f(k)$

又因为 $F(n)=\sum _{n|d}f(d)$ ，故定义 $F(n)=\sum _{x=1}^N\sum _{y=1}^M[n|gcd(x,y)]$

反演得： $f(n)=\sum _{n|d}\mu (\frac{d}{n})F(d)$

又因为由 $F$ 的定义， $n|gcd(x,y)\iff n|x\wedge n|y$ ，故 $F(d)=\lfloor \frac{N}{d}\rfloor \lfloor \frac{M}{d}\rfloor$

所以 $f(n)=\sum _{n|d}\mu (\frac{d}{n})\lfloor \frac{N}{d}\rfloor \lfloor \frac{M}{d}\rfloor$

枚举倍数比枚举因数要简单，所以我们令 $d^{\prime }=\frac{d}{n},N^{\prime }=\frac{N}{n},M^{\prime }=\frac{M}{n}$ ，原式化为 $f(n)=\sum _{d^{\prime }}\mu (d^{\prime })\lfloor \frac{N^{\prime }}{d}\rfloor \lfloor \frac{M^{\prime }}{d}\rfloor$ 预处理 $\mu$ 的前缀和后用数论分块可以求

代码：

#include<bits/stdc++.h>
#define LL long long
using namespace std;
const int N=50000+5;
int pri[N],cnt,mu[N];
bool v[N];
LL sum[N];
void pred(){
	mu[1]=1;
	for(int i=2;i<N;++i){
		if(!v[i]) pri[++cnt]=i,mu[i]=-1;
		for(int j=1;j<=cnt&&pri[j]*i<N;++j){
			v[pri[j]*i]=true;
			if(i%pri[j]==0){
				mu[i*pri[j]]=0;
				break;
			}
			mu[i*pri[j]]=-mu[i];
		}
	}
	sum[0]=0;
	for(int i=1;i<N;++i) sum[i]=sum[i-1]+mu[i];
}
int g(int k,int x){ //得到(k/x)中x所在块的右端点
	return k/(k/x);
}
LL S(int x,int y,int k){
	LL res=0;
	x=x/k,y=y/k;
	int _t=min(x,y);
	for(int l=1,r;l<=_t;l=r+1){
		r=min(_t,min(g(x,l),g(y,l))); //每次在两个取整中找小的,注意不能超过_t
		res+=(sum[r]-sum[l-1])*(x/l)*(y/l);
	}
	return res;
}
int main(){
	int n;
	int a,b,c,d,k;
	pred();
	scanf("%d",&n);
	while(n--){
		scanf("%d%d%d%d%d",&a,&b,&c,&d,&k);
		printf("%lld\n",S(b,d,k)-S(a-1,d,k)-S(b,c-1,k)+S(a-1,c-1,k));
	}
	return 0;
}

例题二

约数个数和

分析：

先看原题给我们的函数，由于 $d$ 是积性函数（忘了的上去翻），所以有 $d(ij)=\sum _{x|i}\sum _{y|j}[gcd(x,y)==1]$ （呃呃呃，鬼想的到啊！）

简单证明：

设 $i=p_1^{{\alpha }_1}{p}_2^{{\alpha }_2}\dots p_k^{{\alpha }_k},j=p_1^{{\beta }_1}{p}_2^{{\beta }_2}\dots p_k^{{\beta }_k}$ ，则 $ij=i=p_1^{{\alpha }_1+{\beta }_1}{p}_2^{{\alpha }_2+{\beta }_2}\dots p_k^{{\alpha }_k+{\beta }_k}$ ，由约数个数定理 $d(ij)=({\alpha }_1+{\beta }_1+1)({\alpha }_2+{\beta }_2+1)\dots ({\alpha }_k+{\beta }_k+1)$

再看等式右边， $x,y$ 的因数也一定是在 $p_1,p_2,\dots ,p_k$ 中，对于质因子 $p_1$ ，若 $x$ 含 $p_1^t,(t\ne 0)$ ，则有意义的 $y$ 一定不含 $p_1$ ，此时 $t$ 有 ${\alpha }_1$ 种取值，同理，若 $y$ 含 $p_1^t,(t\ne 0)$ 有 ${\beta }_1$ 种取值，加上一种 $x,y$ 都不含 $p_1$ ，共 ${\alpha }_1+{\beta }_1+1$ 种有意义的取值，类似的，讨论其它质因数，乘法原理得左式等于右式

QED

我们发现 $Ans=\sum _{i=1}^N\sum _{j=1}^M\sum _{x|i}\sum _{y|j}[gcd(x,y)==1]$ ，类似与上题，令： $f(n)=\sum _{i=1}^N\sum _{j=1}^M\sum _{x|i}\sum _{y|j}[gcd(x,y)==n],F(n)=\sum _{i=1}^N\sum _{j=1}^M\sum _{x|i}\sum _{y|j}[n|gcd(x,y)]$ 那么 $Ans=f(1)$

则 $F(n)=\sum _{n|d}f(d)$ ，反演得 $f(n)=\sum _{n|d}\mu (\frac{d}{n})F(d)$

故 $f(1)=\sum _{d=1}^N\mu (d)F(d)$ ，为了求出 $f(1)$ ，我们来看 $F(n)$ ，由于 $[n|gcd(x,y)]$ 与 $i,j$ 无关，我们可以枚举 $x,y$ 计算有多少 $i,j$ 符合要求， $F(n)=\sum _{x=1}^N\sum _{y=1}^M\lfloor \frac{N}{x}\rfloor \lfloor \frac{M}{y}\rfloor [n|gcd(x,y)]$

同上题一样，可设 $x^{\prime }=\frac{x}{n},y^{\prime }=\frac{y}{n},N^{\prime }=\frac{N}{n},M^{\prime }=\frac{M}{n}$ ，得 $F(n)=\sum _{x^{\prime }=1}^{N^{\prime }}\sum _{y^{\prime }=1}^{M^{\prime }}\lfloor \frac{N^{\prime }}{x^{\prime }}\rfloor \lfloor \frac{M^{\prime }}{y^{\prime }}\rfloor$

调换循环顺序得 $F(n)=(\sum _{x^{\prime }=1}^{N^{\prime }}\lfloor \frac{N^{\prime }}{x^{\prime }}\rfloor )(\sum _{y^{\prime }=1}^{M^{\prime }}\lfloor \frac{M^{\prime }}{y^{\prime }}\rfloor )$

令 $h(k)=\sum _{i=1}^k\lfloor \frac{k}{i}\rfloor$ ，有 $F(n)=h(N^{\prime })h(M^{\prime })$ ，故 $f(1)=\sum _{d=1}^N\mu (d)h(\frac{N}{d})h(\frac{M}{d})$ ，预处理所有的 $h(x)$ ，数论分块解决

代码：

#include<bits/stdc++.h>
#define LL long long
using namespace std;
const int N=50000+5;
int pri[N],cnt,mu[N];
LL sum[N],h[N];
bool v[N];
int g(int k,int x){
	return k/(k/x);
}
void pred(){
	mu[1]=1;
	for(int i=2;i<N;++i){
		if(!v[i]) pri[++cnt]=i,mu[i]=-1;
		for(int j=1;j<=cnt&&pri[j]*i<N;++j){
			v[pri[j]*i]=true;
			if(i%pri[j]==0){
				mu[i*pri[j]]=0;
				break;
			}
			mu[i*pri[j]]=-mu[i];
		}
	}
	sum[0]=0;
	for(int i=1;i<N;++i) sum[i]=sum[i-1]+mu[i];
	for(int i=1;i<N;++i)
		for(int l=1,r;l<=i;l=r+1){
			r=min(i,g(i,l));
			h[i]+=(r-l+1)*(i/l);
		}
}
int main(){
	pred();
	int T;
	scanf("%d",&T);
	while(T--){
		int n,m;
		scanf("%d%d",&n,&m);
		LL ans=0;
		int _t=min(n,m);
		for(int l=1,r;l<=_t;l=r+1){
			r=min(_t,min(g(n,l),g(m,l)));
			ans+=(sum[r]-sum[l-1])*h[n/l]*h[m/l];
		}
		printf("%lld\n",ans);
	}
	return 0;
}

例题三

龙哥的问题

本题比较简单，甚至不用莫反

考虑答案：
$$\begin{array}{}\text{(16)}& Ans& =\sum _{i=1}^{n}gcd(i,n)\text{(17)}& & =\sum _{d|n}(d\times \sum _{i=1}^n[gcd(i,n)==d])\text{(18)}& & =\sum _{d|n}(d\times \sum _{i=1}^n[gcd(\frac{i}{d},\frac{n}{d})==1])\text{(19)}& & =\sum _{d|n}d\times \phi (\frac{n}{d})\end{array}$$
对于每一个 $n$ 枚举约数（ $int$ 范围内最多的一个数约数只有 $1600$ 多个） ，再 $\sqrt{n}$ 求 $\phi$ ，本题的范围可以通过

那可否继续优化呢？令 $k=\frac{n}{d}$ ，再继续化一下：
$$\begin{array}{}\text{(20)}& Ans& =\sum _{d|n}d\times \phi (\frac{n}{d})\text{(21)}& & =\sum _{k|n}\frac{n}{k}\times \phi (k)\text{(22)}& & =\sum _{k|n}\frac{n}{k}\times k\prod (1-\frac{1}{p_i})\text{(23)}& & =n\sum _{k|n}\prod (1-\frac{1}{p_i})\end{array}$$
其中 $k_t=\prod p_i^{c_{i,t}}$ ，又因为 $k|n$ ，不妨设 $n=\prod p_i^{{\alpha }_i}$ ，下面的变化比较常用，要记住：
$$\begin{array}{}\text{(24)}& & k|d\text{(25)}& \Rightarrow & \mathrm{\forall }{k}_t,0\le c_{i,t}\le {\alpha }_i\text{(26)}& \Rightarrow & \sum _{k|n}k=\sum _t\prod _i{p}_i^{c_{i,t}}=\prod _t\sum _i{p}_i^{c_{i,t}}\text{(27)}& & =(p_1^0+p_1^2+\dots +p_1^{{\alpha }_1})(p_2^0+p_2^2+\dots +p_2^{{\alpha }_2})\dots \end{array}$$

如果不理解的话可以这样想：每个括号里取一项作 $p_i^{c_{i,t}}$ 乘起来刚好对应一个 $k_t$

那么考虑上式对答案的贡献，我们发现对于所有的 $k$ ，若取第 $i$ 个质因子的指数为 $0$ （即取 $p_i^0$ ）对答案贡献为 $1$ ，其他情况（即 $p_i^c,c\ne 0$ ）对答案的贡献都是 $1-\frac{1}{p_i}$ ，故答案可化为：
$$\begin{array}{}\text{(28)}& Ans& =n\prod _i(1+(1-\frac{1}{p_i})+(1-\frac{1}{p_i})+\dots +(1-\frac{1}{p_i}))\text{(29)}& & =n\prod _i(1+a(1-\frac{1}{p_i}))\end{array}$$
分解 $n$ 的质因子，可以在 $O(\sqrt{n})$ 时间求解，但注意，本题对 $longlong$ 卡的非常死，甚至在计算 $ans$ 时先乘再除都会炸 $longlong$

代码：

#include<bits/stdc++.h>
#define LL long long
using namespace std;

int main(){
	LL n;
	scanf("%lld",&n);
	LL ans=n; //最初的n 
	for(LL i=2;i*i<=n;++i)
		if(n%i==0){
			LL a=0;
			while(n%i==0) ++a,n/=i;
			ans/=i,ans*=i+a*i-a; //一定要先除 
		}
	if(n>1) ans/=n,ans*=n+n-1;
	printf("%lld\n",ans);
	return 0;
}

例题四

来看一道难题：

数字表格

这道题要点很杂很多，我们慢慢看：

首先明确
$$Ans=\prod _{i=1}^n\prod _{j=1}^{m}fb(gcd(i,j))$$
其中 $fb(x)$ 表示Fibonacci数列第 $x$ 项

由于 $\prod$ 不好处理，考虑枚举 $gcd$ ：
$$Ans=\prod _{d=1}^N[fb(d){]}^{f(d)}$$
其中 $N=min(n,m)$ ， $f(d)=\sum _{i=1}^m\sum _{j=1}^n[gcd(i,j)==d]$

发现这个 $f$ 的定义我们非常熟悉，于是一通莫反：
$$\begin{array}{}\text{(30)}& \text{令}& F(d)=\sum _{i=1}^m\sum _{j=1}^n[d|gcd(i,j)]\text{(31)}& \text{则有}& F(d)=\sum _{r|d}f(r)\text{(32)}& \text{反演得}& f(d)=\sum _{d|r}\mu (\frac{r}{d})F(r)\text{(33)}& \text{又有}& d|gcd(i,j)\iff d|i,d|j\text{(34)}& \text{故}& F(d)=\lfloor \frac{n}{d}\rfloor \lfloor \frac{m}{d}\rfloor \text{(35)}& \text{带入得}& f(d)=\sum _{d|r}\mu (\frac{r}{d})\lfloor \frac{n}{r}\rfloor \lfloor \frac{m}{r}\rfloor \text{(36)}& \text{令}w=\frac{r}{d},\text{有}& f(d)=\sum _{w=1}^{\lfloor \frac{N}{d}\rfloor }\mu (w)\lfloor \frac{n}{wd}\rfloor \lfloor \frac{m}{wd}\rfloor \end{array}$$
此时我们发现在莫反后， $f(d)$ 可以在 $O(\sqrt{n})$ 内计算，但是，我们观察 $Ans$ ，对于每个 $d$ ， $f(d)$ 都必须重新计算（因为 $f$ 与 $n,m,d$ 都有关），这样的时间负责度为 $O(Tn(\sqrt{n}+\log n))$ 直接爆炸好吧

必须考虑优化，但 $fb$ 和 $\prod$ 都太不好搞了， $f$ 也不好在变形（我在这里卡了好久），既然两个都不好搞，我们不妨把它们一起来看看：
$$Ans=\prod _{d=1}^N[fb(d){]}^{\sum _{w=1}^{\lfloor \frac{N}{d}\rfloor }\mu (w)\lfloor \frac{n}{wd}\rfloor \lfloor \frac{m}{wd}\rfloor }$$
我们发现一个关键点：每次询问只修改 $n,m$ 而 $d,w$ 对于每次询问来说是一样的，这启发我们把几个与 $n,m$ 无关的项拿出来预处理：

$$\begin{array}{}\text{(37)}& Ans=\prod _{d=1}^N[fb(d){]}^{\sum _{w=1}^{\lfloor \frac{N}{d}\rfloor }\mu (w)\lfloor \frac{n}{wd}\rfloor \lfloor \frac{m}{wd}\rfloor }\text{(38)}& =\prod _{d=1}^N{[fb(d){]}^{\sum _{w=1}^{\lfloor \frac{N}{d}\rfloor }\mu (w)}}^{\lfloor \frac{n}{wd}\rfloor \lfloor \frac{m}{wd}\rfloor }\text{(39)}& \text{令}p=wd,\text{有}Ans=\prod _{p=1}^N{[\prod _{d|p}fb(d){]}^{\mu (\frac{p}{d})}}^{\lfloor \frac{n}{p}\rfloor \lfloor \frac{m}{p}\rfloor }\text{(40)}& \text{令}h(p)=\prod _{d|p}fb(d{)}^{\mu (\frac{p}{d})},\text{有}Ans=\prod _{p=1}^{N}h(p{)}^{\lfloor \frac{n}{p}\rfloor \lfloor \frac{m}{p}\rfloor }\end{array}$$

发现指数位置的 $\lfloor \frac{n}{p}\rfloor \lfloor \frac{m}{p}\rfloor$ 可以数论分块，而 $h$ 可以提前预处理出前缀积，如此，回答询问的总时间变为 $O(T\sqrt{n}\log n)$

下面考虑如何计算 $h$ ，由于因数不好枚举，我们改为枚举倍数：
$$\begin{array}{}\text{(41)}& h(p)& =\prod _{d|p}fb(d{)}^{\mu (\frac{p}{d})}\text{(42)}& & =\prod _{i=1}^{N}fb(\frac{p}{i}{)}^{\mu (i)}\text{(43)}& \text{变形得}h(ij)& =\prod _{i=1}^N\prod _{j=1}^{\lfloor \frac{N}{i}\rfloor }fb(j{)}^{\mu (i)}\end{array}$$
预处理时间为 $O(n\log n)$

代码如下：

#include<bits/stdc++.h>
#define LL long long
using namespace std;
const int N=1000000+5,P=1000000000+7;
LL pri[N],cnt,mu[N],fb[N],h[N],inv_fb[N],mul[N];
bool v[N];
int g(int k,int x){
	return k/(k/x);
}
LL qpow(LL x,LL y){
	LL res=1;
	while(y){
		if(y&1) res=(res*x)%P;
		x=(x*x)%P;
		y>>=1;
	}
	return res;
}
void pred(){
	mu[1]=1;
	for(int i=2;i<N;++i){
		if(!v[i]) pri[++cnt]=i,mu[i]=-1;
		for(int j=1;j<=cnt&&pri[j]*i<N;++j){
			v[pri[j]*i]=true;
			if(i%pri[j]==0){
				mu[i*pri[j]]=0;
				break;
			}
			mu[i*pri[j]]=-mu[i];
		}
	}
	fb[0]=0;
	fb[1]=1;
	inv_fb[1]=1; //特殊处理
	for(int i=2;i<N;++i) fb[i]=(fb[i-1]+fb[i-2])%P,inv_fb[i]=qpow(fb[i],P-2); //求逆元
	for(int i=1;i<N;++i) h[i]=1;
    for(int i=1;i<N;++i){
        if(!mu[i])continue;
        for(int j=1;j*i<N;++j)
        	h[j*i]=h[j*i]*(mu[i]==1?fb[j]:inv_fb[j])%P; //注意mu可能为负一
    }
	mul[0]=1;
	for(int i=1;i<N;++i) mul[i]=(mul[i-1]*h[i])%P;
}
LL get_ans(int x,int y){
	LL res=1,as;
	int _t=min(x,y);
	for(int l=1,r;l<=_t;l=r+1){
		r=min(_t,min(g(x,l),g(y,l)));
		as=mul[r]*qpow(mul[l-1],P-2)%P; //前缀积得底数
		res=res*qpow(as,1ll*(x/l)*(y/l)%(P-1))%P;
	}
	return res;
}
int main(){
	pred();
	int T;
	scanf("%d",&T);
	int n,m;
	while(T--){
		scanf("%d%d",&n,&m);
		printf("%lld\n",get_ans(n,m));
	}
	return 0;
}

补充：变换方法

下面记录一些常见的变换方法：

• $$gcd(i,j)==1\iff \sum _{d|i\wedge d|j}\mu (d)$$

• 像例题四一样的完全展开预处理无关项

• 像例题二一样 $d(ij)=\sum _{x|i}\sum _{y|j}[gcd(x,y)==1]$

• 改变未知数使得两个 $\sum$ 分离，预处理其中一个（也可以看作将两数积设出，继续化简）

  如YY的GCD，莫反后得
  $$\sum _{k\in prime}\sum _{d=1}^{\lfloor \frac{min(N,M)}{k}\rfloor }\mu (d)\lfloor \frac{N}{kd}\rfloor \lfloor \frac{M}{kd}\rfloor$$
  若直接计算，对于每个询问，最坏时间复杂度为 $O(n)$ 会TLE，可以变化如下：

  $$\begin{array}{}\text{(44)}& \text{设}T& =kd\text{(45)}& \text{则原式}& =\sum _{k\in primes}\sum _{d=1}^{\lfloor \frac{min(N,M)}{k}\rfloor }\mu (d)\lfloor \frac{N}{T}\rfloor \lfloor \frac{M}{T}\rfloor \text{(46)}& & =\sum _{T=1}^{min(N,M)}\sum _{k|T\wedge k\in prime}\mu (\frac{T}{k})\lfloor \frac{N}{T}\rfloor \lfloor \frac{M}{T}\rfloor \text{(47)}& & =\sum _{T=1}^{min(N,M)}\lfloor \frac{N}{T}\rfloor \lfloor \frac{M}{T}\rfloor \sum _{k|T\wedge k\in prime}\mu (\frac{T}{k})\end{array}$$

  后面的 $\sum$ 与询问无关，可以预处理

• 注意卷积中的一些常用等式，如 $Id\ast \mu =\phi$ ：

  $$\begin{array}{}\text{(48)}& Ans& =\sum _{i=1}^n\sum _{j=1}^{n}gcd(i,j)\text{(49)}& \text{莫反得}Ans& =\sum _{r=1}^n\lfloor \frac{n}{r}\rfloor \lfloor \frac{n}{r}\rfloor \sum _{d|r}d\ast \mu (\frac{r}{d})\text{(50)}& & =\sum _{r=1}^n\lfloor \frac{n}{r}\rfloor \lfloor \frac{n}{r}\rfloor \phi (r)\end{array}$$