Sandomhet

Möbius

Introduces Mobius function and mobius inversion formula

Thu, Feb 1, 2024 5 min read

Möbius

莫比乌斯函数

定义

μ(n)={1n=10n含有平方因子(−1)kk为本质不同的质因子个数\mu(n)= \begin{cases} 1 &n=1 \\ 0 &n含有平方因子 \\ (-1)^k &k为本质不同的质因子个数 \end{cases}

性质

莫比乌斯函数是积性函数,且有

∑d∣nμ(d)={1n=10n≠1\sum\limits_{d|n}\mu(d)= \begin{cases} 1 &n=1\\ 0 &n\neq 1\\ \end{cases}

即∑d∣nμ(d)=ε(n)\sum\limits_{d|n}\mu(d)=\varepsilon(n),卷积形式为μ∗1=ε\mu*1=\varepsilon

证明

设n=∏i=1kpici,n′=∏i=1kpin=\prod\limits_{i=1}^{k}p_i^{c_i},n'=\prod\limits_{i=1}^{k}p_i ,有莫比乌斯函数的定义可知∑d∣nμ(d)=∑d∣n′μ(d)=∑i=0kCki(−1)i\sum\limits_{d|n}\mu(d)=\sum\limits_{d|n'}\mu(d)=\sum\limits_{i=0}^{k}C_k^i(-1)^i ,(把枚举因数换了一种形式),根据二项式定理,∑i=0kCki(−1)i1k−i=(1−1)k\sum\limits_{i=0}^{k}C_k^i(-1)^i1^{k-i}=(1-1)^k,所以当k=0(n=1)k=0(n=1)时为1,其余为0。得证。

莫比乌斯变换

公式

1.若f(n)=∑d∣ng(d)f(n)=\sum\limits_{d|n}g(d),则g(n)=∑d∣nμ(d)f(nd)g(n)=\sum\limits_{d|n}\mu(d)f(\frac{n}{d})

2.若f(n)=∑n∣dg(d)f(n)=\sum\limits_{n|d}g(d),则g(n)=∑n∣dμ(d)f(dn)g(n)=\sum\limits_{n|d}\mu(d)f(\frac{d}{n})

证明

1.∑d∣nμ(d)f(nd)=∑d∣nμ(d)∑k∣ndg(k)=∑k∣ng(k)∑d∣nkμ(d)=∑k∣ng(k)ε(nk)=g(n)\sum\limits_{d|n}\mu(d)f(\frac{n}{d})=\sum\limits_{d|n}\mu(d)\sum\limits_{k|\frac{n}{d}}g(k)=\sum\limits_{k|n}g(k)\sum\limits_{d|\frac{n}{k}}\mu(d)=\sum\limits_{k|n}g(k)\varepsilon(\frac{n}{k})=g(n)

2.∑n∣dμ(d)f(dn)=∑k=1+∞μ(kn)f(k)=∑k=1+∞μ(kn)∑k∣dg(d)=∑n∣dg(d)∑kn∣dμ(kn)=∑n∣dg(d)∑k∣dnμ(k)=∑n∣dg(d)ε(dn)=g(n)\sum\limits_{n|d}\mu(d)f(\frac{d}{n})=\sum\limits_{k=1}^{+\infty}\mu(kn)f(k)=\sum\limits_{k=1}^{+\infty}\mu(kn)\sum\limits_{k|d}g(d)=\sum\limits_{n|d}g(d)\sum\limits_{kn|d}\mu(kn)=\sum\limits_{n|d}g(d)\sum\limits_{k|\frac{d}{n}}\mu(k)=\sum\limits_{n|d}g(d)\varepsilon(\frac{d}{n})=g(n)

常用套路(例题+解释)

  1. [gcd⁡(i,j)==1]=ε(gcd⁡(i,j))=∑d∣gcd⁡(i,j)μ(d)[\gcd(i, j) == 1]=\varepsilon(\gcd(i,j))=\sum\limits_{d|\gcd(i, j)}\mu(d)

  2. ∑i=1n∑j=1m[gcd⁡(i,j)==k] \sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}[\gcd(i,j)==k]

题目链接:YY的GCD、[HAOI2011]Problem b、[POI2007]ZAP-Queries

∑i=1n∑j=1m[gcd⁡(i,j)==k]=∑i=1⌊nk⌋∑j=1⌊mk⌋[gcd⁡(i,j)==1]=∑i=1⌊nk⌋∑j=1⌊mk⌋∑d∣gcd⁡(i,j)μ(d)=∑d=1nμ(d)∑i=1⌊nkd⌋∑j=1⌊mkd⌋=∑d=1nμ(d)⌊nkd⌋⌊mkd⌋\begin{aligned} \sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}[\gcd(i,j)==k] &= \sum\limits_{i=1}^{\lfloor\frac{n}{k}\rfloor}\sum\limits_{j=1}^{\lfloor\frac{m}{k}\rfloor}[\gcd(i,j)==1] \\ &=\sum\limits_{i=1}^{\lfloor\frac{n}{k}\rfloor}\sum\limits_{j=1}^{\lfloor\frac{m}{k}\rfloor}\sum\limits_{d|\gcd(i, j)}\mu(d) \\ &=\sum\limits_{d=1}^{n}\mu(d)\sum\limits_{i=1}^{\lfloor\frac{n}{kd}\rfloor}\sum\limits_{j=1}^{\lfloor\frac{m}{kd}\rfloor} \\ &=\sum\limits_{d=1}^{n}\mu(d)\lfloor\frac{n}{kd}\rfloor\lfloor\frac{m}{kd}\rfloor \\ \end{aligned}
  1. ∑i=1n∑j=1mgcd⁡(i,j) \sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\gcd(i,j)
∑i=1n∑j=1mgcd⁡(i,j)=∑k=1nk∑i=1n∑j=1m[gcd⁡(i,j)==k]=∑k=1nk∑i=1⌊nk⌋∑j=1⌊mk⌋[gcd⁡(i,j)==1]=∑k=1nk∑i=1⌊nk⌋∑j=1⌊mk⌋∑d∣gcd⁡(i,j)μ(d)=∑k=1nk∑d=1nμ(d)∑i=1,d∣i⌊nk⌋∑j=1,d∣j⌊mk⌋=∑k=1nk∑d=1nμ(d)∑i=1⌊nkd⌋∑j=1⌊mkd⌋=∑k=1n∑d=1nkμ(d)⌊nkd⌋⌊mkd⌋令t=kd,原式=∑t=1n∑d∣tμ(d)(td)⌊nt⌋⌊mt⌋\begin{aligned} \sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\gcd(i,j) &=\sum\limits_{k=1}^{n}k\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}[\gcd(i, j)==k] \\ &=\sum\limits_{k=1}^{n}k\sum\limits_{i=1}^{\lfloor \frac{n}{k} \rfloor}\sum\limits_{j=1}^{\lfloor \frac{m}{k} \rfloor}[\gcd(i, j)==1] \\ &=\sum\limits_{k=1}^{n}k\sum\limits_{i=1}^{\lfloor \frac{n}{k} \rfloor}\sum\limits_{j=1}^{\lfloor \frac{m}{k} \rfloor}\sum\limits_{d|\gcd(i, j)}\mu(d) \\ &=\sum\limits_{k=1}^{n}k\sum\limits_{d=1}^{n}\mu(d)\sum\limits_{i=1, d|i}^{\lfloor \frac{n}{k} \rfloor}\sum\limits_{j=1, d|j}^{\lfloor \frac{m}{k} \rfloor} \\ &=\sum\limits_{k=1}^{n}k\sum\limits_{d=1}^{n}\mu(d)\sum\limits_{i=1}^{\lfloor \frac{n}{kd} \rfloor}\sum\limits_{j=1}^{\lfloor \frac{m}{kd} \rfloor} \\ &=\sum\limits_{k=1}^{n}\sum\limits_{d=1}^{n}k\mu(d)\lfloor \frac{n}{kd} \rfloor\lfloor \frac{m}{kd} \rfloor \\ 令t=kd,原式&=\sum\limits_{t=1}^{n}\sum\limits_{d|t}\mu(d)(\frac{t}{d})\lfloor \frac{n}{t} \rfloor\lfloor \frac{m}{t} \rfloor \\ \end{aligned}

定义f(t)=∑d∣tμ(d)(td)f(t)=\sum\limits_{d|t}\mu(d)(\frac{t}{d})

原式=∑t=1nf(t)⌊nt⌋⌊mt⌋原式=\sum\limits_{t=1}^{n}f(t)\lfloor \frac{n}{t} \rfloor\lfloor \frac{m}{t} \rfloor
  1. ∑i=1n∑j=1mlcm(i,j) \sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}lcm(i,j)

题目链接:[国家集训队]Crash的数字表格/JZPTAB

∑i=1n∑j=1mlcm(i,j)=∑i=1n∑j=1mijgcd⁡(i,j)=∑k=1n∑i=1n∑j=1mijk[gcd⁡(i,j)==k]=∑k=1n∑i=1⌊nk⌋∑j=1⌊mk⌋ijk2k[gcd⁡(i,j)==1]=∑k=1nk∑i=1⌊nk⌋∑j=1⌊mk⌋ij∑d∣gcd⁡(i,j)μ(d)=∑k=1n∑d=1nkμ(d)d2∑i=1⌊nkd⌋∑j=1⌊mkd⌋ij令t=kd,原式=∑t=1nt∑d∣tdμ(d)∑i=1⌊nt⌋i∑j=1⌊mt⌋j\begin{aligned} \sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}lcm(i,j) &=\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\frac{ij}{\gcd(i,j)} \\ &=\sum\limits_{k=1}^{n}\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\frac{ij}{k}[\gcd(i,j)==k] \\ &=\sum\limits_{k=1}^{n}\sum\limits_{i=1}^{\lfloor \frac{n}{k} \rfloor}\sum\limits_{j=1}^{\lfloor \frac{m}{k} \rfloor}\frac{ijk^2}{k}[\gcd(i,j)==1] \\ &=\sum\limits_{k=1}^{n}k\sum\limits_{i=1}^{\lfloor \frac{n}{k} \rfloor}\sum\limits_{j=1}^{\lfloor \frac{m}{k} \rfloor}ij\sum\limits_{d|\gcd(i, j)}\mu(d) \\ &=\sum\limits_{k=1}^{n}\sum\limits_{d=1}^{n}k\mu(d)d^2\sum\limits_{i=1}^{\lfloor \frac{n}{kd} \rfloor}\sum\limits_{j=1}^{\lfloor \frac{m}{kd} \rfloor}ij \\ 令t=kd,原式&=\sum\limits_{t=1}^{n}t\sum\limits_{d|t}d\mu(d)\sum\limits_{i=1}^{\lfloor \frac{n}{t} \rfloor}i\sum\limits_{j=1}^{\lfloor \frac{m}{t} \rfloor}j \\ \end{aligned}

定义f(t)=∑d∣tdμ(d),g(n)=∑i=1nif(t)=\sum\limits_{d|t}d\mu(d),g(n)=\sum\limits_{i=1}^{n}i,ff显然可以筛出来,gg显然是等差数列通项公式。

原式=∑t=1ntf(t)g(⌊nt⌋)g(⌊mt⌋)原式=\sum\limits_{t=1}^{n}tf(t)g(\lfloor \frac{n}{t} \rfloor)g(\lfloor \frac{m}{t} \rfloor)
  1. ∑i=1n∑j=1md(ij) \sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}d(ij) (dd为约数个数和)

题目链接:[SDOI2015]约数个数和

非常重要的性质:d(ij)=∑x∣i∑y∣j[gcd⁡(x,y)==1]d(ij)=\sum\limits_{x|i}\sum\limits_{y|j}[\gcd(x,y)==1]

证明:考虑dd函数的定义——约数个数和,d(n)=∏p=1k(cp+1)=∑p∣nd(n)=\prod\limits_{p=1}^{k}(c_p+1)=\sum\limits_{p|n}\quad 因为每个质数相互独立,所以先只考虑质因子pp,假设ii的因子中有aa个pp,jj中有bb个pp,那么pp的贡献就是(a+b+1)(a+b+1)。x、yx、y 互质即是要求pp最多只能出现在一方。从组合角度来理解:1.ii中选pp,jj不选pp,aa种方案;2.jj中选pp,ii不选pp,bb 种方案;3.i、ji、j都不选pp,11种方案。根据乘法原理与乘法结合律,pp的贡献为(a+b+1)(a+b+1) 。同样的,也可以利用映射关系:1.对于次数小于等于aa的pp因子,强制让它从ii中选择;2.当大于时,让不够的那部分从jj 中选择。由于两者是一一映射,所以正确性可以保证。 (目前貌似无法用严格的数学公式推导,所以只能通过逻辑感性理解了。)

有了性质,我们就可以大力推式子了!!

∑i=1n∑j=1md(ij)=∑i=1n∑j=1m∑x∣i∑y∣j[gcd⁡(x,y)==1]=∑x=1n∑y=1m[gcd⁡(x,y)==1]∑i=1⌊nx⌋∑j=1⌊my⌋=∑d=1nμ(d)∑x=1⌊nd⌋∑y=1⌊md⌋⌊nxd⌋⌊myd⌋=∑d=1nμ(d)∑x=1⌊nd⌋⌊nxd⌋∑y=1⌊md⌋⌊myd⌋\begin{aligned} \sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}d(ij) &=\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\sum\limits_{x|i}\sum\limits_{y|j}[\gcd(x,y)==1] \\ &=\sum\limits_{x=1}^{n}\sum\limits_{y=1}^{m}[\gcd(x,y)==1]\sum\limits_{i=1}^{\lfloor\frac{n}{x}\rfloor}\sum\limits_{j=1}^{{\lfloor\frac{m}{y}\rfloor}} \\ &=\sum\limits_{d=1}^{n}\mu(d)\sum\limits_{x=1}^{\lfloor \frac{n}{d} \rfloor}\sum\limits_{y=1}^{\lfloor \frac{m}{d} \rfloor}\lfloor\frac{n}{xd}\rfloor \lfloor\frac{m}{yd}\rfloor \\ &=\sum\limits_{d=1}^{n}\mu(d)\sum\limits_{x=1}^{\lfloor \frac{n}{d} \rfloor}\lfloor\frac{n}{xd}\rfloor\sum\limits_{y=1}^{\lfloor \frac{m}{d} \rfloor} \lfloor\frac{m}{yd}\rfloor \\ \end{aligned}

定义f(n)=∑i=1n⌊ni⌋f(n)=\sum\limits_{i=1}^{n}\lfloor\frac{n}{i}\rfloor,可以通过整除分块预处理出来O(1)O(1)查询。

原式=∑d=1nμ(d)f(⌊nd⌋)f(⌊md⌋)原式=\sum\limits_{d=1}^{n}\mu(d)f(\lfloor \frac{n}{d} \rfloor) f(\lfloor \frac{m}{d} \rfloor)
代码实现
#define sandom signed
#include <cstdio>
#include <iostream>
#define int long long 
using namespace std;
const int Z = 5e4 + 10;
inline int read() { int x = 0; int f = 0; char c = getchar(); while (!isdigit(c)) f |= c == '-', c = getchar(); while (isdigit(c)) x = (x << 1) + (x << 3) + (c ^ 48), c = getchar(); return f ? -x : x; }
inline void write(int x) { static int wrt[20]; int TP = 0; if (x < 0) putchar('-'), x = -x; while (x >= 10) wrt[++TP] = x % 10, x /= 10; wrt[++TP] = x; while (TP) putchar(wrt[TP--] | 48); putchar('\n'); }

int n, m, s, ans;
bool prime[Z];
int ip[Z];
int miu[Z], f[Z];
void Mobius(int n) {
miu[1] = 1;
for (int i = 2; i <= n; i++) {
if (!prime[i]) ip[++ip[0]] = i, miu[i] = -1;
for (int j = 1; j <= ip[0]; j++) {
int k = i * ip[j];
if (k > n) break;
prime[k] = 1;
if (i % ip[j] == 0) { miu[k] = 0; break; }
else miu[k] = -miu[i];
}
}
for (int i = 1; i <= n; i++) miu[i] += miu[i - 1];//前缀和
}
void init(int n) { //整除分块预处理的函数
for (int k = 1; k <= n; k++)
for (int l = 1, r; l <= k; l = r + 1) {
r = k / (k / l);
f[k] += (k / l) * (r - l + 1);
}
}
inline int NTB(int n, int m) {
int k = min(n, m);
int ans = 0;
for (int l = 1, r; l <= k; l = r + 1) {
r = min(n / (n / l), m / (m / l));
ans += f[n / l] * f[m / l] * (miu[r] - miu[l - 1]);
}
return ans;
}

sandom main() {
Mobius(5e4); init(5e4);
int T = read();
while (T--) {
n = read(), m = read();
write(NTB(n, m));
}
return 0;
}
  1. ∑i=1n∑j=1mφ(ij) \sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\varphi(ij)

题目链接:毒瘤之神的考验

重要的性质:φ(ij)=φ(i)φ(j)gcd⁡(i,j)φ(gcd⁡(i,j))\varphi(ij)=\frac{\varphi(i)\varphi(j)\gcd(i,j)}{\varphi(\gcd(i,j))}

证明:从φ\varphi 的定义式入手,易得φ(ij)=ij∏p∣i or p∣j(1−1p),φ(gcd⁡(i,j))=gcd(i,j)∏p∣i and p∣j(1−1p),φ(i)φ(j)=ij∏p∣i(1−1p)∏p∣j(1−1p)\varphi(ij)=ij\prod\limits_{p|i\ or\ p|j}(1-\frac{1}{p}), \varphi(\gcd(i, j))=gcd(i,j)\prod\limits_{p|i\ and\ p|j}(1-\frac{1}{p}),\varphi(i)\varphi(j)=ij\prod\limits_{p|i}(1-\frac{1}{p})\prod\limits_{p|j}(1-\frac{1}{p}) ,根据容斥原理得证。

∑i=1n∑j=1mφ(ij)=∑i=1n∑j=1mφ(i)φ(j)gcd⁡(i,j)φ(gcd⁡(i,j))=∑k=1nkφ(k)∑i=1nφ(i)∑j=1mφ(j)[gcd⁡(i,j)==k]=∑k=1nkφ(k)∑i=1⌊nk⌋φ(ik)∑j=1⌊mk⌋φ(jk)∑d∣gcd⁡(i,j)μ(d)=∑k=1nkφ(k)∑d=1nμ(d)∑i=1⌊nkd⌋φ(ikd)∑j=1⌊mkd⌋φ(jkd)令t=kd,原式=∑t=1n∑k∣tkμ(tk)φ(k)∑i=1⌊nt⌋φ(it)∑j=1⌊mt⌋φ(jt)\begin{aligned} \sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\varphi(ij) &=\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\frac{\varphi(i)\varphi(j)\gcd(i,j)}{\varphi(\gcd(i,j))} \\ &=\sum\limits_{k=1}^{n}\frac{k}{\varphi(k)}\sum\limits_{i=1}^{n}\varphi(i)\sum\limits_{j=1}^{m}\varphi(j)[\gcd(i,j)==k] \\ &=\sum\limits_{k=1}^{n}\frac{k}{\varphi(k)}\sum\limits_{i=1}^{\lfloor \frac{n}{k} \rfloor}\varphi(ik)\sum\limits_{j=1}^{\lfloor \frac{m}{k} \rfloor}\varphi(jk)\sum\limits_{d|\gcd(i,j)}\mu(d) \\ &=\sum\limits_{k=1}^{n}\frac{k}{\varphi(k)}\sum\limits_{d=1}^{n}\mu(d)\sum\limits_{i=1}^{\lfloor \frac{n}{kd} \rfloor}\varphi(ikd)\sum\limits_{j=1}^{\lfloor \frac{m}{kd} \rfloor}\varphi(jkd) \\ 令t=kd,原式&=\sum\limits_{t=1}^{n}\sum\limits_{k|t}\frac{k\mu(\frac{t}{k})}{\varphi(k)}\sum\limits_{i=1}^{\lfloor \frac{n}{t} \rfloor}\varphi(it)\sum\limits_{j=1}^{\lfloor \frac{m}{t} \rfloor}\varphi(jt) \\ \end{aligned}

定义f(t)=∑k∣tkμ(tk)φ(k),g(t,n)=∑i=1nφ(it)f(t)=\sum\limits_{k|t}\frac{k\mu(\frac{t}{k})}{\varphi(k)},g(t, n)=\sum\limits_{i=1}^{n}\varphi(it),

原式=∑t=1nf(t)g(t,⌊nt⌋)g(t,⌊mt⌋)原式=\sum\limits_{t=1}^{n}f(t)g(t, \lfloor \frac{n}{t} \rfloor)g(t, \lfloor \frac{m}{t} \rfloor)

但是发现还是不能在O(n)O(\sqrt{n})时间内得出答案,继续定义h(a,b,n)=∑t=1nf(t)g(t,a)g(t,b)h(a, b, n)=\sum\limits_{t=1}^{n}f(t)g(t,a)g(t,b) ,可以通过这个函数进行类似于根号分治的做法来平衡复杂度。

首先ff函数使用埃氏筛即可,g、hg、h显然有递推式g(t,n)=g(t,n−1)+φ(nt),h(a,b,n)=h(a,b,n−1)+f(n)g(n,a)g(n,b)g(t,n)=g(t,n-1)+\varphi(nt),h(a,b,n)=h(a,b,n-1)+f(n)g(n,a)g(n,b)。

hh因为内存问题,显然不能全部预处理出来。而根据整除函数的性质,当tt比较小时,⌊nt⌋\lfloor \frac{n}{t} \rfloor 较大,且连续段很短,根本没有必要整除分块,直接暴算即可。对于后面的tt,⌊nt⌋\lfloor \frac{n}{t} \rfloor 较小,因为内存开的下,所以可以选择预处理。记块长为BB,取[35,50][35, 50]都行,来作为hh前两维所管辖的范围。

当然需要用到一些卡常技巧,比如说动态内存池。

代码实现
#define sandom signed
#include <bits/stdc++.h>
using namespace std; 
const int Z = 1e5 + 2; const int M = 50; const int mod = 998244353;
inline char getc() { static char buf[1 << 18], *p1, *p2; if (p1 == p2) { p1 = buf, p2 = buf + fread(buf, 1, 1 << 18, stdin); if (p1 == p2) return EOF; } return *p1++; }
inline int read() { int x = 0, f = 0; char c = getc(); while (!isdigit(c)) f = c == '-', c = getc(); while (isdigit(c)) x = (x << 1) + (x << 3) + (c ^ 48), c = getc(); return f ? -x : x; }
inline void write(int x) { static int wrt[20]; int TP = 0; if (x < 0) putchar('-'), x = -x; while (x >= 10) wrt[++TP] = x % 10, x /= 10; wrt[++TP] = x; while (TP) putchar(wrt[TP--] | 48); putchar('\n'); }
inline int max(int a, int b) { return a > b ? a : b; } inline int min(int a, int b) { return a < b ? a : b; }
inline int qpow(int a, int b, int c) { int res = 1; while (b) { if (b & 1) res = 1ll * res * a % c; a = 1ll * a * a % c; b >>= 1; } return res; }

int n, m, ans;
bool prime[Z];
int ip[Z], miu[Z], phi[Z], inv[Z];
int f[Z], *g[Z], *h[M + 2][M + 2];
void Linear(int n) { //线性筛莫比乌斯函数和欧拉函数
    miu[1] = phi[1] = 1;
    for (int i = 2; i <= n; ++i) {
        if (!prime[i]) {
            ip[++ip[0]] = i;
            miu[i] = -1; phi[i] = i - 1;
        }
        for (int j = 1; j <= ip[0]; ++j) {
            int k = i * ip[j];
            if (k > n) break;
            prime[k] = 1;
            if (i % ip[j] == 0) { miu[k] = 0, phi[k] = phi[i] * ip[j]; break; }
            else miu[k] = -miu[i], phi[k] = phi[i] * (ip[j] - 1);
        }
    } 
    for (int i = 1; i <= n; ++i) inv[i] = qpow(phi[i], mod - 2, mod);//计算欧拉函数的逆元
}
void init(int n, int m) { //递推预处理相关函数
    for (int j = 1; j <= n; ++j) for (int i = j; i <= n; i += j) (f[i] += 1ll * j * miu[i / j] * inv[j] % mod) %= mod;
    for (int j = 1; j <= n; ++j) {
        int tmp = n / j;
        g[j] = new int [tmp + 2](); //():清零,g[j][0] = 0
        for (int i = 1; i <= tmp; ++i) g[j][i] = (g[j][i - 1] + phi[i * j]) % mod;
    }
    for (int a = 1; a <= m; ++a) for (int b = 1; b <= m; ++b) {
        int tmp = n / max(a, b);
        h[a][b] = new int [tmp + 2](); 
        for (int i = 1; i <= tmp; ++i) h[a][b][i] = (h[a][b][i - 1] + 1ll * f[i] * g[i][a] % mod * g[i][b] % mod) % mod;
    }
}
int calc(int a, int b, int l, int r) { //小块暴力
    int res = 0;
    for (int i = l; i <= r; ++i) (res += 1ll * f[i] * g[i][a] % mod * g[i][b] % mod) %= mod;
    return res;
}
int solve(int n, int m) { //整除分块
    ans = 0;
    for (int l = 1, r; l <= n; l = r + 1) {
        int tmp1 = n / l, tmp2 = m / l;
        r = min(n / tmp1, m / tmp2);
        if (max(tmp1, tmp2) > M) (ans += calc(tmp1, tmp2, l, r)) %= mod;//直接暴力计算
        else (ans += (h[tmp1][tmp2][r] - h[tmp1][tmp2][l - 1]) % mod) %= mod;//预处理出来了
    }
    return ans;
}

sandom main() {
    int T = read();
    Linear(1e5); init(1e5, M);
    while (T--) {
        n = read(), m = read();
        if (n > m) swap(n, m);
        write((solve(n, m) + mod) % mod);
    }
    return 0;
}
Search

Loading search…